assigning arbitrary values to surface parcels: method 2

The previous post describes a method for assigning arbitrary values to surface MMP parcels via GIfTI files. Tim Coalson kindly pointed me to another method, which I'll demonstrate here. Both methods work, but one or the other might be easier in particular situations.

In the previous post I used matlab and ran wb_command from the command prompt (on Windows, by opening a command window in the directory with wb_command.exe, then using full paths to the input and output files). Here, I use R, and call wb_command from within R using its system() function. You may need to adjust system for other operating systems, or simply replace it with print and copy-paste the full line to the command prompt.

 rm(list=ls());  # clear R's memory  
   
 setwd("d:/Workbench/workbench-v1.2.3/bin_windows64/"); # location of wb_command.exe, to call wb_command functions via system() within R  
   
 local.path <- "d:/temp/";  # local directory for reading and writing files  
 s1200.path <- "d:/Workbench/HCP_S1200_GroupAvg_v1/"; # HCP S1200 Group Average Data Release  
 template.fname <- "MMPtemplate.pscalar.nii";  # filename for the template we'll make  
   
 # make a template pscalar.nii from the MMP atlas  
 system(paste0("wb_command -cifti-parcellate ", s1200.path, "S1200.thickness_MSMAll.32k_fs_LR.dscalar.nii ",   
        s1200.path, "Q1-Q6_RelatedValidation210.CorticalAreas_dil_Final_Final_Areas_Group_Colors.32k_fs_LR.dlabel.nii COLUMN ",  
        local.path, template.fname));  
 if (!file.exists(paste0(local.path, template.fname))) { stop(paste0("missing: ", local.path, template.fname)); }  
   
 # you can make a text version of the template, which has 360 rows, but the values in each row aren't the parcel numbers.  
 # system(paste0("wb_command -cifti-convert -to-text ", local.path, template.fname, " ", local.path, "temp.txt"));  
   
 # the text file needs to be arranged with the 180 right hemisphere parcels in the first 180 rows (in order),   
 # then the 180 parcels for the left hemisphere.  
 # build a text file with the values to plot.   
 out.vec <- rep(0,360);  
 out.vec[c(1,10,15)] <- 1;  # right hemisphere parcels given value 1  
 out.vec[c(180+1,180+10,180+15)] <- 2;  # corresponding left hemisphere parcels given value 2  
 write.table(out.vec, paste0(local.path, "plotDemo.txt"), col.names=FALSE, row.names=FALSE);  
   
 # create a CIFTI from the text file for viewing in Workbench  
 system(paste0("wb_command -cifti-convert -from-text ", local.path, "plotDemo.txt ", local.path, template.fname, " ", local.path, "plotDemo.pscalar.nii"));  
 if (!file.exists(paste0(local.path, "plotDemo.pscalar.nii"))) { stop(paste("missing:", local.path, "plotDemo.pscalar.nii")); }  

And here is the resulting image, with parcels 1, 10, and 15 plotted in pink on the right hemisphere (1), and yellow on the left (2). Instructions for generating this view are below the jump.

The Workbench screenshot above has the borders of the MMP parcels, and a custom coloring to make parcels with the value 1 pink, and 2, yellow. If you've never used Workbench, start with this tutorial.You should already have the HCP 1200 release files locally, along with the plotDemo.pscalar.nii file made above.

First open Workbench without opening a spec file (click the Skip button). Choose File -> Open File, and navigate to HCP_S1200_GroupAvg_v1 (the HCP 1200 release files). Put the Files of type: selector box to Surface Files (*surf.gii) and pick S1200.R.midthickness_MSMAll.32k_fs_LR.surf.gii and S1200.L.midthickness_MSMAll.32k_fs_LR.surf.gii (or your favorite inflation). Workbench should now show blank brains.

Next, go to File -> Open File again, but set the Files of type: selector to "Connectivity - Dense Label Files (*.dlabel.nii)" and pick Q1-Q6_RelatedValidation210.CorticalAreas_dil_Final_Final_Areas_Group_Colors.32k_fs_LR.dlabel.nii. Workbench won't show the file yet, but it should have been loaded (and its name shown in the File boxes in the Overlay ToolBox).

Go to File -> Open File one more time, and navigate to the location of plotDemo.pscalar.nii and open that file (Files of type: "Connectivity - Parcel Scalar Files (*.pscalar.nii)").

Now that all the files are loaded, they can be displayed. To reproduce the screenshot above, first click the first two Layers in the Overlay Toolbox On (green arrow), and set the top File to the Q1-Q6... and the second to the plotDemo.pscalar.nii (blue arrows).

It's hard to see the plotDemo parcels on top of the MMP parcels (from the Q1-Q6... file), so switch the MMP parcels to just show the borders. Do this by clicking the little wrench icon in the top row (blue arrow), which brings up the Overlay and Map Settings window. On the Labels tab switch the Drawing Type to Outline Color and click Close.

Now we can see the plotDemo parcels, but not in the right colors. Click the wrench icon in the second row to bring up the Overlay and Map Settings window again, but with many more options (since plotDemo.pscalar.nii and the Q1-Q6 files are of different types). To reproduce the pink and yellow colors shown above, set the Palette to HSB8_clrmid (blue lines) and click Close. Finally, to show the color bar click the rainbow icon in the second row, marked with a red arrow.