Autorecon2: Difference between revisions

From CCN Wiki
Jump to navigation Jump to search
 
(25 intermediate revisions by 5 users not shown)
Line 63: Line 63:
We will probably find ourselves getting more experience with the issues described on this Freesurfer tutorial as time goes on: [https://surfer.nmr.mgh.harvard.edu/fswiki/FreeviewGuide/FreeviewWorkingWithData/FreeviewEditingaRecon]
We will probably find ourselves getting more experience with the issues described on this Freesurfer tutorial as time goes on: [https://surfer.nmr.mgh.harvard.edu/fswiki/FreeviewGuide/FreeviewWorkingWithData/FreeviewEditingaRecon]


'''Regenerating the Surface'''
===Regenerating the Surface===
When you are finished editing the voxels, you will need to regenerate the surfaces. Since the white matter hasn't been changed, you don't need to resegment the volume. You can regenerate the pial surface with:
When you are finished editing the voxels, you will need to regenerate the surfaces. Since the white matter hasn't been changed, as we do not currently edit the white matter (can be done with control points, here: https://surfer.nmr.mgh.harvard.edu/fswiki/FsTutorial/ControlPoints_tktools), you don't need to resegment the entire volume (e.g. you don't need to completely re-run autorecon2.sh). You can simply regenerate the pial surface with:
 
recon-all -autorecon-pial -subjid <subject_id>
 
Or use the autorecon2_fix.sh script.
 
autorecon2_fix.sh FS_T1_501
 
== The Mysterious Case of the Missing Chunks ==
 
After running autorecon2, you may notice that your surfaces don't exclude some of the anatomical voxels. This really only can be seen when checking quality via tkmedit, so be sure to do that. We've been primarily seeing this around the anterior temporal lobes. The cause is unknown, but we suspect it is related to the larger head coil used with the Semantic Imagery experiment.
 
We have been experimenting with different methods to enhance the T1 MPRAGE image so that autorecon2 is better able to detect the white matter tracts and the gray matter boundaries and avoid missing large chunks of the anterior temporal lobes. The following technique seems to work, and it makes use of the segmentation files that SPM generates, and the intensity thresholding provided by MRICron. As of Feb 8, 2018, only ''Brocasarea'', ''Accumbens'' and ''Wernickesarea'' have SPM12 configured correctly, so you will need to copy the subject/mri folder to one of those computers and work there.
 
===Incorporating SPM Segmentation Files into Recon-All===
The first step will be to segment the original MPRAGE file as described in the section on SPM Segmentation under Brainmask Editing for [[Autorecon1#SPM_Segmentation|Autorecon1]]. If you haven't done so already, follow those instructions first for the problematic brain.
 
At the end of this process, you will have a series of .nii files, prefixed with C1 to C5. We will be recombining these files and using MRICron to modify the WM.mgz file output by recon-all during the autorecon2 step.
 
===Copy Files to Fix Directory===
The next bunch of work will take place in the fix directory that contains all your segmented files. We will need the wm.mgz generated by the autorecon2 step. You can find that in $SUBJECTS_DIR/$SUBJID/mri/wm.mgz. We can use mri_convert to copy it over as a .nii, which is the format we need the file in anyways:
SUBJID=FS_0238
FIXDIR=~/Fix/0238/mri
mri_convert ${SUBJECTS_DIR}/${SUBJID}/mri/wm.mgz ${FIXDIR}/wm.nii
mri_convert ${SUBJECTS_DIR}/${SUBJID}/mri/brainmask.mgz ${FIXDIR}/brainmask.nii
===Explore Masks in MRICron===
At this point, we can load the brainmask.nii into MRICron, and then load the various segmentation masks as overlays to see whether the FreeSurfer white matter segmentation missed some critical tissue in the missing regions. In MRICron, select '''Overlay>Add''' and select wm.nii, which should then appear as a semitransparent red mask. Do likewise for the c2001.nii (SPM white matter segmentation) and c1001.nii (SPM gray matter segmentation) files. As you load these overlays, you should get an idea of what tissues SPM identified as WM that do not appear in the wm.mgz file. Likewise, the SPM gray matter volume should help identify the extent of the pial surface.
===3dcalc===
There is no prescribed formula for this, but the information from the previous step can be used to guide what sort of operations you might want to do in 3dcalc to generate a fixed wm.mgz or enhanced T1.mgz file. Below are a few plausible scenarios:
====White Matter Voxels Missing From wm.mgz====
If the SPM WM segmentation looks reasonable, you might consider simply using that file in place of wm.mgz and/or filled.mgz (wm.mgz is used to generate filled.mgz). I will note that SPM includes the cerebellum in its WM segmentation, whereas wm.mgz does not. I'm not certain if it's important to remove it, but you might simply convert a rescaled version of c2001.nii into filled.mgz. The voxels in c2001.nii range from 0.0 to 1.00, and are thus SPM's probabalistic determination that voxel v contains white matter, whereas filled.mgz contains values ranging from 0 to 255, so we rescale accordingly. This file can be converted to a .mgz file and replace the one generated by FreeSurfer:
3dcalc -a c2001.nii -expr 'a*255' -prefix rescaled_spm_wm.nii
mri_convert rescaled_spm_wm.nii filled.mgz
 
====Gray Matter Too Noisy/Dark====
In a problem volume I'm looking at now, the FreeSurfer WM segmentation seemed to get *most* of what SPM found, so that may not entirely explain the missing anterior temporal tissue. What else might be the problem? I noted that the brainmask.mgz volume is pretty noisy in this area, so perhaps the apparently random distribution of higher/lower values within (what should be) the same tissue is causing a problem. Perhaps making the GM tissue a bit more homogeneous would help. I might start by selecting GM voxels below a certain value and increasing them to a threshold minimum (perhaps the average value of the rest of the GM volume). According to the [[https://surfer.nmr.mgh.harvard.edu/fswiki/ReconAllDevTable|Recon-All wiki]] the tesselation step (the part that turns voxels into triangles) takes filled.mgz and norm.mgz as input volumes. This tells us that we can try using the SPM gray matter mask to fix norm.mgz and then re-tesselate.
#Use mri_convert to copy norm.mgz to a .nii file in your fix directory, as described above
#Load norm.nii into MRICron
#From the Draw menu, select '''Draw>Intensity Filter...'''
#*A dialog box will appear, allowing you to select the minimum/maximum values of the underlying image to be included in your mask.
#*Adjust the min/max value until the masked area (blue) seems to adequately capture white matter
#*Save your mask VOI. We will use this to modify the SPM GM segmentation file and create a new GM mask. I called mine something like "intensity_filter_wm.voi". The name doesn't matter; just that it makes sense.
#Load the SPM GM segmentation file (c001.nii)
#Load the intensity filter mask you just created
#From the Draw menu, select '''Draw>Mask image with voi>Delete regions with VOI'''
#*This will remove from your GM mask those voxels that you identified as WM
#Save your modified GM mask as a new file (e.g., "GM_intensity.nii")
At this point you now have a new GM mask, with values ranging from 0.0 to 1.0. However this volume may not match the voxel dimensions of norm.mgz, as I just discovered (FreeSurfer resampled norm.mgz into a 256x256x256 volume, whereas SPM's volume had fewer slices in the Z-plane). Fortunately, this just adds an extra couple of steps:
#Load norm.nii into MRICron
#Load your GM mask as an ROI ('''Draw>Open Voi...''', select GM_intensity.nii, or whatever is appropriate)
#Delete all voxels in norm.nii that are '''NOT''' in the mask:
#*'''Draw>Mask image with voi>Preserve regions with VOI'''
#Save your masked image, which should have the same dimensions as norm.nii
#*'''File>Save as NIfTI...''' (give it an appriate name, e.g., "GM_intensity_mask.nii")
We can now use this mask to modify norm.nii with combination of functions in 3dcalc to selectively boost certain voxel values. Suppose, for example, we want to bring all the GM voxels that are below 50 up to 75 (remember, norm.nii voxels range from 0-255, with WM having an average value of 110, and I found that the maximum value for GM voxels in my mask was 80):
3dcalc -a GM_intensity_Mask.nii -expr 'astep(a,50)' -prefix binaryGM.nii #binary mask, setting voxels>50 to 1, otherwise 0
3dcalc -a norm.nii -b binaryGM.nii -expr 'max(a,b*75)' -prefix boosted.nii
#boosted contains either the original norm.nii value or [0,1]*75, whichever is greater
3dcalc -a norm.nii -b boosted.nii -expr '(a+b)/2' -prefix norm2.nii #average the boosted and norm volumes
#the averaged volume will now be less noisy and have slightly brighter GM voxels
mri_convert norm2.nii norm.mgz #create mgz file
The next step would be to copy any new enhanced mgz files back to the FreeSurfer subject directory and restart the tesselation process
===Retesselation===
After you have regenerated filled.mgz and/or norm.mgz, you can see if these modified files will tesselate any better. Recon-all can be started at any point in the process, so you will be picking up at the -tesselate step:
$SUBJID=FS_0238
recon-all -s $SUBJID -tesselate
(<2 minutes)
It then looks like you can do the next few steps in one fell swoop:
recon-all -s $SUBJID -smooth1 -inflate1 -qsphere
(~10 minutes)
And then a few more steps all at once:
recon-all -s $SUBJID -fix -white -smooth2 -inflate2
(~4.5 hours)
 
 
 
 
[[Category: FreeSurfer]]

Latest revision as of 20:50, 8 February 2018

The -autorecon2 directive tells Freesurfer to convert the cubic voxel-based T1 (high-resolution anatomical) image into a triangular mesh surface model. It generates a separate surface mesh for each hemisphere, and from this point on in the process, each hemisphere will have its own independent set of files.

If you checked your data carefully after autorecon1 has finished you should find that the second step goes fairly smoothly.

There is a script, autorecon2.sh that is a wrapper for the recon-all command with the -autorecon2 directive. It simply requires you to specify a valid subjectid within your $SUBJECTS_DIR. Note that a valid subjectid will be the name of the directory generated by Freesurfer when executing the -autorecon1 directive. For example, if my source data prior to running Freesurfer came from:

  • /home
    • /mydirectory
      • /fmri
        • /sub_001
          • /MPRAGE

My $SUBJECTS_DIR might be set to be /home/mydirectory/fmri. After running autorecon1.sh, depending on the value given for the -subjid flag, a new directory might now exist containing all the resulting outputs, e.g.:

  • /home
    • /mydirectory
      • /fmri
        • /FS_001

For all subsequent processing steps, the subject id is FS_001 for this data set, and not 001.

If you have an autorecon2.sh script in your ~/bin/ directory, inspect it first before running it, and ensure that you understand what it's doing. For reference, here's what a script of mine looks like:

#!/bin/bash
# Usage: autorecon2.sh <subjectid>

PROJECTROOT=~/ubfs/cpmcnorg/openfmri/booth/ #Change this: where are your data located?

OLD_SUBJECTS_DIR=${SUBJECTS_DIR} #Your subjects directory might be initially set to something else
SUBJECTS_DIR=${PROJECTROOT} #now tell freesurfer where the subject in question resides
export SUBJECTS_DIR 

echo SUBJECTS_DIR: ${SUBJECTS_DIR}

subjectID=$1
recon-all -autorecon2 \
  -subjid ${subjectID}

SUBJECTS_DIR=${OLD_SUBJECTS_DIR} #reset subjects_dir to whatever it was before
export SUBJECTS_DIR 

When I ran my autorecon1.sh script, it created a new folder called FS_T?_* for the individual (I was working on longitudinal data, where there was a T1 and T2 time point for each participant). When I ran it for the T1 data for subject 501, it created a folder in $SUBJECTS_DIR called FS_T1_501. After editing the brainmask.mgz file for this participant, I ran autorecon2.sh like this:

autorecon2.sh FS_T1_501

The -autorecon2 directive takes several hours to complete. Initial estimates were as long as 11 hours, though it recently took only about 5 or 6 hours to run on a different participant. Either way, you should plan to do some other work while you're waiting for this step to be completed.

N.B.: This script may report that it exited with errors, with the error log file indicating that a script was unable to find ribbon.mgz. As far as we can tell, this error is inconsequential and can be ignored.

Quality Control

Note: There is a systematic blip that occurs on the medial left hemisphere near the corpus callosum and ventricle areas. Don't worry about this, since it's not grey matter anyway.

After this step has completed, you can overlay the surface mesh vector on to the source T1 anatomical image in tkmedit:

tkmedit <subject_id> brainmask.mgz -aux T1.mgz -surfs

You will want to confirm that the pial boundaries, which define where the grey matter is located, correspond to the source image data. The following wiki page has an example of where it includes dura or skull: Pial Editing Another potential problem that we have seen is where the grey matter signal is too faint, causing this stage to omit these voxels in the mesh. In one example, the top of the head appeared "flattened" because of signal drop-off. I'll hopefully be able to find a solution to this problem with the Freesurfer people.

Next, load up the generated mesh in a program called tksurfer:

tksurfer <subject_id> <hemi> <surf>

Where subject_id is the name that you passed to autorecon2.sh (e.g., FS_T1_501), hemi is the hemisphere you want to load (one of either lh or rh), and the surf flag indicates which surface you want to view (typically I look at the 'inflated' surface. tksurfer FS_1001 lh inflated

You can rotate the surface around along all 3D axes and confirm that it looks like what you would expect a normal brain to look like. Errors in the surface construction would be evident by strangely-shaped features, like spikes or pits along the surface.

Another useful application of our surface mesh is that it lets us fix potentially catastrophic misalignments manually (to get the functional images in the right ballpark) as described in the FreeSurfer wiki

We will probably find ourselves getting more experience with the issues described on this Freesurfer tutorial as time goes on: [1]

Regenerating the Surface

When you are finished editing the voxels, you will need to regenerate the surfaces. Since the white matter hasn't been changed, as we do not currently edit the white matter (can be done with control points, here: https://surfer.nmr.mgh.harvard.edu/fswiki/FsTutorial/ControlPoints_tktools), you don't need to resegment the entire volume (e.g. you don't need to completely re-run autorecon2.sh). You can simply regenerate the pial surface with:

recon-all -autorecon-pial -subjid <subject_id>

Or use the autorecon2_fix.sh script.

autorecon2_fix.sh FS_T1_501

The Mysterious Case of the Missing Chunks

After running autorecon2, you may notice that your surfaces don't exclude some of the anatomical voxels. This really only can be seen when checking quality via tkmedit, so be sure to do that. We've been primarily seeing this around the anterior temporal lobes. The cause is unknown, but we suspect it is related to the larger head coil used with the Semantic Imagery experiment.

We have been experimenting with different methods to enhance the T1 MPRAGE image so that autorecon2 is better able to detect the white matter tracts and the gray matter boundaries and avoid missing large chunks of the anterior temporal lobes. The following technique seems to work, and it makes use of the segmentation files that SPM generates, and the intensity thresholding provided by MRICron. As of Feb 8, 2018, only Brocasarea, Accumbens and Wernickesarea have SPM12 configured correctly, so you will need to copy the subject/mri folder to one of those computers and work there.

Incorporating SPM Segmentation Files into Recon-All

The first step will be to segment the original MPRAGE file as described in the section on SPM Segmentation under Brainmask Editing for Autorecon1. If you haven't done so already, follow those instructions first for the problematic brain.

At the end of this process, you will have a series of .nii files, prefixed with C1 to C5. We will be recombining these files and using MRICron to modify the WM.mgz file output by recon-all during the autorecon2 step.

Copy Files to Fix Directory

The next bunch of work will take place in the fix directory that contains all your segmented files. We will need the wm.mgz generated by the autorecon2 step. You can find that in $SUBJECTS_DIR/$SUBJID/mri/wm.mgz. We can use mri_convert to copy it over as a .nii, which is the format we need the file in anyways:

SUBJID=FS_0238
FIXDIR=~/Fix/0238/mri
mri_convert ${SUBJECTS_DIR}/${SUBJID}/mri/wm.mgz ${FIXDIR}/wm.nii
mri_convert ${SUBJECTS_DIR}/${SUBJID}/mri/brainmask.mgz ${FIXDIR}/brainmask.nii

Explore Masks in MRICron

At this point, we can load the brainmask.nii into MRICron, and then load the various segmentation masks as overlays to see whether the FreeSurfer white matter segmentation missed some critical tissue in the missing regions. In MRICron, select Overlay>Add and select wm.nii, which should then appear as a semitransparent red mask. Do likewise for the c2001.nii (SPM white matter segmentation) and c1001.nii (SPM gray matter segmentation) files. As you load these overlays, you should get an idea of what tissues SPM identified as WM that do not appear in the wm.mgz file. Likewise, the SPM gray matter volume should help identify the extent of the pial surface.

3dcalc

There is no prescribed formula for this, but the information from the previous step can be used to guide what sort of operations you might want to do in 3dcalc to generate a fixed wm.mgz or enhanced T1.mgz file. Below are a few plausible scenarios:

White Matter Voxels Missing From wm.mgz

If the SPM WM segmentation looks reasonable, you might consider simply using that file in place of wm.mgz and/or filled.mgz (wm.mgz is used to generate filled.mgz). I will note that SPM includes the cerebellum in its WM segmentation, whereas wm.mgz does not. I'm not certain if it's important to remove it, but you might simply convert a rescaled version of c2001.nii into filled.mgz. The voxels in c2001.nii range from 0.0 to 1.00, and are thus SPM's probabalistic determination that voxel v contains white matter, whereas filled.mgz contains values ranging from 0 to 255, so we rescale accordingly. This file can be converted to a .mgz file and replace the one generated by FreeSurfer:

3dcalc -a c2001.nii -expr 'a*255' -prefix rescaled_spm_wm.nii
mri_convert rescaled_spm_wm.nii filled.mgz

Gray Matter Too Noisy/Dark

In a problem volume I'm looking at now, the FreeSurfer WM segmentation seemed to get *most* of what SPM found, so that may not entirely explain the missing anterior temporal tissue. What else might be the problem? I noted that the brainmask.mgz volume is pretty noisy in this area, so perhaps the apparently random distribution of higher/lower values within (what should be) the same tissue is causing a problem. Perhaps making the GM tissue a bit more homogeneous would help. I might start by selecting GM voxels below a certain value and increasing them to a threshold minimum (perhaps the average value of the rest of the GM volume). According to the [wiki] the tesselation step (the part that turns voxels into triangles) takes filled.mgz and norm.mgz as input volumes. This tells us that we can try using the SPM gray matter mask to fix norm.mgz and then re-tesselate.

  1. Use mri_convert to copy norm.mgz to a .nii file in your fix directory, as described above
  2. Load norm.nii into MRICron
  3. From the Draw menu, select Draw>Intensity Filter...
    • A dialog box will appear, allowing you to select the minimum/maximum values of the underlying image to be included in your mask.
    • Adjust the min/max value until the masked area (blue) seems to adequately capture white matter
    • Save your mask VOI. We will use this to modify the SPM GM segmentation file and create a new GM mask. I called mine something like "intensity_filter_wm.voi". The name doesn't matter; just that it makes sense.
  4. Load the SPM GM segmentation file (c001.nii)
  5. Load the intensity filter mask you just created
  6. From the Draw menu, select Draw>Mask image with voi>Delete regions with VOI
    • This will remove from your GM mask those voxels that you identified as WM
  7. Save your modified GM mask as a new file (e.g., "GM_intensity.nii")

At this point you now have a new GM mask, with values ranging from 0.0 to 1.0. However this volume may not match the voxel dimensions of norm.mgz, as I just discovered (FreeSurfer resampled norm.mgz into a 256x256x256 volume, whereas SPM's volume had fewer slices in the Z-plane). Fortunately, this just adds an extra couple of steps:

  1. Load norm.nii into MRICron
  2. Load your GM mask as an ROI (Draw>Open Voi..., select GM_intensity.nii, or whatever is appropriate)
  3. Delete all voxels in norm.nii that are NOT in the mask:
    • Draw>Mask image with voi>Preserve regions with VOI
  4. Save your masked image, which should have the same dimensions as norm.nii
    • File>Save as NIfTI... (give it an appriate name, e.g., "GM_intensity_mask.nii")

We can now use this mask to modify norm.nii with combination of functions in 3dcalc to selectively boost certain voxel values. Suppose, for example, we want to bring all the GM voxels that are below 50 up to 75 (remember, norm.nii voxels range from 0-255, with WM having an average value of 110, and I found that the maximum value for GM voxels in my mask was 80):

3dcalc -a GM_intensity_Mask.nii -expr 'astep(a,50)' -prefix binaryGM.nii #binary mask, setting voxels>50 to 1, otherwise 0
3dcalc -a norm.nii -b binaryGM.nii -expr 'max(a,b*75)' -prefix boosted.nii 
#boosted contains either the original norm.nii value or [0,1]*75, whichever is greater
3dcalc -a norm.nii -b boosted.nii -expr '(a+b)/2' -prefix norm2.nii #average the boosted and norm volumes
#the averaged volume will now be less noisy and have slightly brighter GM voxels
mri_convert norm2.nii norm.mgz #create mgz file

The next step would be to copy any new enhanced mgz files back to the FreeSurfer subject directory and restart the tesselation process

Retesselation

After you have regenerated filled.mgz and/or norm.mgz, you can see if these modified files will tesselate any better. Recon-all can be started at any point in the process, so you will be picking up at the -tesselate step:

$SUBJID=FS_0238
recon-all -s $SUBJID -tesselate

(<2 minutes) It then looks like you can do the next few steps in one fell swoop:

recon-all -s $SUBJID -smooth1 -inflate1 -qsphere 

(~10 minutes) And then a few more steps all at once:

recon-all -s $SUBJID -fix -white -smooth2 -inflate2

(~4.5 hours)