Differences between revisions 11 and 12
Deletions are marked like this. Additions are marked like this.
Line 1: Line 1:
Please refer to this site/make edits here for the most updated information: https://partnershealthcare.sharepoint.com/sites/LCN-Exvivo/SitePages/Exvivo-Recon.aspx

Reconning Exvivo Data

Exvivo data cannot be pushed through the normal recon-all stream due to the changes in tissue contrast after fixation. Below are the steps the FreeSurfer Development team follows to recon their exvivo data. It assumes you are working with 1mm isotropic, multiecho flash scans (8 echoes) with flip angles of 30, 25, 20, 15, 10 and 5 degrees collected at 1.5T.

Please feel free to add to this wiki page if you have found other methods that work as our method. Additionally, it is intended to be run on single hemisphere data. Others have been able to make it work with whole brain data but we do not have those steps documented at this time.

Recon steps

Directory Structure

The exvivo recon procedure assumes the same directory structure as the invivo recon procedure. You can make all the relevant directories with the following command:

mksubjdirs  subjid

Processing Scans

  1. Convert dicoms for each echo of each flip to mgz format. For example:
    foreach e (0 1 2 3 4 5 6 7)
    mri_convert -nth ${e} <mef30 run1 dicom> ../mri/flash/mef30_run1_echo${e}.mgz
    mri_convert -nth ${e} <mef25 run1 dicom> ../mri/flash/mef25_run1_echo${e}.mgz
    mri_convert -nth ${e} <mef20 run1 dicom> ../mri/flash/mef20_run1_echo${e}.mgz
    mri_convert -nth ${e} <mef15 run1 dicom> ../mri/flash/mef15_run1_echo${e}.mgz
    mri_convert -nth ${e} <mef10 run1 dicom> ../mri/flash/mef10_run1_echo${e}.mgz
    mri_convert -nth ${e} <mef5 run1 dicom> ../mri/flash/mef5_run1_echo${e}.mgz
    mri_convert -nth ${e} <mef30 run2 dicom> ../mri/flash/mef30_run2_echo${e}.mgz
    mri_convert -nth ${e} <mef25 run2 dicom> ../mri/flash/mef25_run2_echo${e}.mgz
    mri_convert -nth ${e} <mef20 run2 dicom> ../mri/flash/mef20_run2_echo${e}.mgz
    mri_convert -nth ${e} <mef15 run2 dicom> ../mri/flash/mef15_run2_echo${e}.mgz
    mri_convert -nth ${e} <mef10 run2 dicom> ../mri/flash/mef10_run2_echo${e}.mgz
    mri_convert -nth ${e} <mef5 run2 dicom> ../mri/flash/mef5_run2_echo${e}.mgz
    The -nth flag will create individual files for each echo stored in the dicom.
  2. Average the echoes for each flip angle using mri_average. For example:
    foreach fa (5 10 15 20 25 30)
    mri_average -noconform mef${fa}_run*_echo*.mgz  mef${fa}_avg.mgz
    This will call all the echoes for all the runs that exist for a particular flip angle and average them.
  3. Synthesize T1, T2, and PD maps using mri_ms_fitparms:
     mri_ms_fitparms  input  outdir
    Input should be all the individual echoes of all the flip angles collected. The different flip angles allow fitparms to synthesize T1 contrast, etc. You may want to consider registering each echo to the first or an average if there is motion during the scan.

Choose Orig Volume

For the orig.mgz, you want to choose a volume with even wm intensities and good grey matter/white matter contrast (like an MPRAGE). Look at the parameter map volumes and flip angle averages you created to find a volume with uniform wm and good contrast. Make a copy of this volume called orig.mgz and put it in the mri dir.

Tip: You may want to try out a few different volumes in the kvlThresholdImage tool (described below) to see which combination of volumes will give the best wm segmentation. You can skip to those instructions but be sure to come back here to follow the steps on orienting and conforming the orig and any other volumes you are using in the threshold step. Then run the kvlThresholdImage tool on the conformed volumes and follow the rest of the Segment White Matter steps to get the conformed wm segmentation BEFORE making any edits to the wm.

Orient the Orig

You may find it easiest to put the exvivo volume in the correct RAS coordinates if it is not already.

Copy RAS Coordinates to all Volumes

If orientation was necessary, you will need to copy the new orientation parameters to the other volumes you will use to segment the white matter (described below). You can use mri_copy_params, mri_convert, or mri_vol2vol to do this. Just be sure to check that all the volumes are aligned after doing so.

Conform all Volumes

In order for the exvivo data to run through the FreeSurfer pipeline, it must be conformed to a 256x256x256 matrix. You can do this with the following command:

mri_convert -c  input  output

You will want to do this to the orig.mgz and all other volumes being used for Segmenting the White Matter (described below) and the

Segment White Matter

You may have to spend a significant amount of time trying out the different methods described below with different volumes to see what will get you the best wm segmentation. Volumes can be improved by adding control points to the wm to get more uniformity. No segmentation will be perfect and will likely require a lot of manual intervention (you will have to fill in all the subcortical structures and will have to probably erase many voxels surrounding the wm that are included in error). You will want to find the right balance of not too many incorrect voxels being included in the segmentation and not too many correct voxels being excluded.

  • Threshold Out Noise
    The first step is to threshold out structures that are not of any interest (air, fluid, sponges, etc.) and thus create a mask consisting of only the tissue classes you want to segment. By loading multiple volumes, you can interactively define threshold values on one volume and apply the resulting mask on another volume. Choose the volume(s) where fluid and air are a different contrast from brain. It's okay if there is not a lot of contrast between wm and gm. In the past, flash30.mgz created by mri_ms_fitparms has worked well. The orig.mgz will be used as the volume to create the mask from (to be used for segmentation later).

    kvlThresholdImage  orig.mgz  /path/to/second_volume.mgz

    Where it says Image to mask, choose the orig volume. Where it says Image to threshold, choose the other volume you loaded. Use the scroll bar to change the Lower threshold and Upper threshold until you see a purple overlay covering all the stuff that you do not want (you do not want to cover gm or wm).

    • Use the View section to choose whether you see all planes or just one.

    • You can scroll through the slices using the scroll bars at the top right of the screen in order to see how much you are masking out on the other slices.
    • Check where it says Show inverse mask to see the actual mask that will be created. This will show everything that you will keep for the segmentation step.

    • You can change the Overlay opacity to see the volume underneath. The lower threshold is originally set to the minimum intensity value in the image (0) and the upper threshold defaults to the maximum intensity value. Whatever is covered by the overlay will be set to an intensity value of zero in the masked image.

    Pressing the Mask image button will write out a file called "*_masked.mgz" (in this case orig_masked.mgz). This file will be the one you want to segment with in the next step (kvlEMSegment). (You can also choose to press the Write mask button to get a mask.mgz volume of just the mask and not the volume with noise masked out.)

    The second step is to actually segment the wm from the gm. This will classify tissues using an automatically estimated Gaussian mixture and polynomial MR bias field model. You would use this command from the directory where the masked volume is:

     kvlEMSegment  orig_masked.mgz  numberOfGaussians  biasFieldOrder  downSamplingFactor
    • numberOfGaussians = number of tissue classes. (You can try 2 for wm and gm or 3 for wm, gm, and CSF/plp.)
    • biasFieldOrder = a number between 4 and 8. (The higher the #, the higher the spatial frequencies in MR inhomogeneities the method tries to compensate for.)
    • downSamplingFactor= a number between 1 and 4. (The higher the #, the less voxels the method uses to estimate the model parameters: this speeds up computation time but possibly at a loss of segmentation accuracy. For example, when set to 4 this tells the algorithm to only consider every 4th voxel in the x-, y-, and z-directions for calculating the model parameters, speeding the computations up by a factor of 4^3 = 64.)

    The suggested values to start with are:

     kvlEMSegment  orig_masked.mgz 2 4 4

    Note: Every case is different so it is worth playing around with the different values above if you are unhappy with the segmentation (definitely try using 3 different tissue classes and/or a down sampling factor of 1 to see if it improves your segmentation).

    View Segmentation Results
    Results are written out as crispPosterior_classX.mgz, one for each class (with 0 being wm and 1 being gm). To see the result of the automated segmentation overlaid on the original MR image, you can use this viewer:

     kvlViewImage  orig_masked.mgz  crispPosterior_class0.mgz
    Or if you have several wm segmentations you would like to view at the same time to see which is better, you can do this in freeview:
     freeview -v orig.mgz  wm_file1.mgz:colormap=heatscale  wm_file2.mgz:colormap=jet

    If you are not satisfied with the segmentation, you can tweak the thresholds in the kvlThresholdImage step or try out different parameters in the kvlEMSegment step. If you still do not feel you are getting a good segmentation, try out some of the optional tips for improving it described below in the Optional section.

When you are happy with the segmentation you have, rename crispPosterior_class0.mgz to wm_init.mgz. You can now move on to the Editing WM Volume step. (It may be worth running the pretess steps in the Optional section below even if you are happy with the segmentation - it will probably add in some wm voxels and cut down on some editing.)

Important: Before moving on to wm editing, be sure that the volumes used to create the wm segmentation were conformed to 256 x 256 x 256.


If you can't get a good segmentation, try out some of the suggestions below.

Add Control Points

If there are any areas in the white matter that are not uniform, open the orig volume in tkmedit or freeview and add control points throughout wm at least 1 voxel away from the gm border. Be careful not to put any control points in subcortical structures. The control points will be saved as tmp/control.dat if you are using tkmedit (make sure you have a tmp dir created). If you are using freeview, it will save them in the directory you are in - you should either move the control.dat into the tmp dir or change the commands below to reflect their actual location.

Tip: You can load the wm segmentation as a 2nd volume while you add control points so you can see what areas need control points the most.

Run the commands below to normalize the volume. This will create a brain.mgz, T1.mgz and brain_mask.mgz.

cd mri
mkdir nu

mri_convert orig.mgz nu/nu0.mnc
cd nu
nu_correct -clobber nu0.mnc nu1.mnc
nu_correct -clobber nu1.mnc nu2.mnc
nu_correct -clobber nu2.mnc nu3.mnc
nu_correct -clobber nu3.mnc nu4.mnc
mri_convert nu4.mnc nu.mgz

cd ../

foreach s (1 2 3 4 5 6 7 8 9 10)
mri_normalize -sigma 10 -fonly ../tmp/control.dat ../mri/nu/nu.mgz ../mri/brain.mgz
cp ../mri/brain.mgz ../mri/T1.mgz
cp ../mri/brain.mgz ../mri/brain_mask.mgz


After this finishes, try the Threshold and Segmentation steps again using the T1.mgz for the thresholding. If satisfied, continue to the Editing WM section.

Use mri_pretess

This step can be used to fix pathological configurations in the wm segmentation (such as wm voxels touching only at the corner). mri_pretess is run later on while making surfaces but it might be worth running it before then if you are having trouble getting a good segmentation. It will not likely fill in huge gaps but might fill in smaller ones. It's a good idea to make a back-up of your wm segmentation before running this.

cd mri

foreach x (1 2 3 4 5 6 7 8 9 10)
mri_pretess wm_init.mgz  wm  brain.mgz  wm_fixed.mgz

For every iteration, it adds voxels to areas where voxels are connected on a diagonal. It uses brain.mgz as reference and outputs wm_fixed.

If you do not have a brain.mgz from the control points step above, you can use the orig.mgz instead.

Edit WM Segmentation

tkmedit -f orig.mgz -aux wm_init.mgz


freeview -v orig.mgz  wm_init.mgz:colormap=heatscale

Note: Must change brush value to 255 when wm_init.mgz is highlighted in order to be able to add voxels to an overlay in freeview.

Fill in ventricles (look at orig volume to see the borders of these subcortical structures), accumbens, caudate, vessels, putamen, pallidum, thalamus, ventralDC, lesions, hypointensities, and basal ganglia. Connect the temporal lobe if necessary. Remove the cerebellum if there is any. Hippocampus, the ventricle above hippocampus, and amygdala should not be filled in (they are too much trouble and don't matter because they do not cause topological problems). For an example of what needs to be filled in, look at the in vivo example of a wm.mgz which comes with your FreeSurfer installation, bert. Check the different views (sagittal, coronal, transversal) to see where connections need to be made. You can more quickly fill in huge gaps by using Shift+Middle click in tkmedit to fill in anything that is within an enclosed circle. Erase extraneous voxels or gm voxels labeled as wm which are not wm.

When you are done making edits, copy wm_edit.mgz to wm.mgz which will be used to create surfaces.

Make Surfaces

Note: If you did not go through the control points step, you will need to create a brain.mgz before running the steps below. It can be created by making a symbolic link to the orig.mgz:

ln -s orig.mgz brain.mgz

Set your SUBJECTS_DIR to the dir above where this subject is located. Run the following steps to create surfaces (be sure to change the appropriate values if it is a left hemisphere; what is given below is for a right hemisphere):

*Fixes pathological configurations in the wm.mgz; output is wm_pretess.mgz:
mri_pretess mri/wm.mgz wm  mri/brain.mgz  mri/wm_pretess.mgz

*Output is filled.mgz. Change -findrhv to -findlhv if running this on a left hemisphere.
mri_fill -fillonly -findrhv  mri/wm_pretess.mgz mri/filled.mgz

*Change to 127 if rh; 255 if lh
mri_tessellate mri/filled.mgz 127 surf/rh.orig
cp surf/rh.orig surf/rh.orig.nofix
mris_smooth surf/rh.orig surf/rh.smoothwm

*Inflates surface.
mris_inflate surf/rh.smoothwm surf/rh.inflated
cp surf/rh.inflated surf/rh.inflated.nofix

*Output is spherical surface, ?h.qsphere
mris_sphere -w 0 -inflate -in 200 -q surf/rh.inflated surf/rh.qsphere
cp surf/rh.qsphere surf/rh.qsphere.nofix

*Specify lh or rh
mris_fix_topology -mgz  subjid  rh

*Counts number of errors after topology is fixed. Want it to be zero. Still goes on if not zero but Bruce may have problems later if it is not zero.
mris_euler_number surf/rh.orig

*Smoothes tessellation of surface. We run this again since mris_fix_topology has fixed the errors and corrected the surfaces
mris_smooth surf/rh.orig surf/rh.smoothwm

*we run this again after mris_fix_topology
mris_inflate surf/rh.smoothwm surf/rh.inflated

**Use either mris_make_surfaces or mris_exvivo_surfaces. Generally, we use mris_make surfaces. You can only use mris_exvivo_surfaces if you have a multiecho flash scan with flip angles 5 and 30 (you will need to conform these volumes first before running the command). If you don't have those, use mris_make_surfaces.
mris_make_surfaces -w 0 -mgz -noaparc -noaseg subjid rh
#mris_exvivo_surfaces -mgz -W 0 -flash30 mef30_avg -flash5 mef5_avg  subjid rh

*Creates rh.sphere; registers surfaces to atlas; makes cortical parcellation and stats
recon-all -cortribbon -sphere -surfreg -cortparc -parcstats -subjid  subjid -hemi rh -noaseg

Check the Surfaces

View your results to make sure the surfaces accurately follow the gm/wm border in the orig.mgz:

freeview -v orig.mgz wm.mgz -f ../surf/?h.white:edgecolor=blue ../surf/?h.pial:edgecolor=red

Replace ? with 'r' or 'l' depending on the hemisphere. Turn off the white surface (Main Surface)- the orig surface is usually more accurate than the white surface for labeling the wm.

You can also view it with tkmedit if you set your SUBJECTS_DIR.

setenv SUBJECTS_DIR /path/to/subject
tkmedit subjid orig.mgz ?h.white -aux wm.mgz

ExvivoRecon (last edited 2021-07-12 13:27:58 by DevaniCordero)