|
Size: 12464
Comment:
|
← Revision 74 as of 2021-05-03 07:53:08 ⇥
Size: 20095
Comment:
|
| Deletions are marked like this. | Additions are marked like this. |
| Line 1: | Line 1: |
| Compared with cross-sectional studies, a longitudinal design can significantly reduce the confounding effect of inter-individual morphological variability by using each subject as his or her own control. As a result, longitudinal imaging studies are getting increased interest and popularity in various aspects of neuroscience. The standard FreeSurfer pipeline is designed for the processing of individual data set, and thus not optimal for the processing of longitudinal data series. How to take into account the subject-wise and temporal correlation in a longitudinal data series and get more robust and reliable cortical and subcortical morphological measurements is an active research area in NMR center. Although further improvements are constantly being sought, as the result of an early attempt, a longitudinal processing pipeline has been designed for the processing of longitudinal data and has shown proved success. The basic idea of this longitudinal scheme is to borrow the processing results of an early time point to help the processing of a later time point. As easily noticed, the FreeSurfer cortical and subcortical segmentation and parcellation procedure involves solving many complex nonlinear optimization problems, such as the deformable surface reconstruction, the nonlinear atlas-image registration, and the nonlinear spherical surface registration. These nonlinear optimization problems are usually solved using iterative methods, and the final results are known to be highly sensitive to the selection of a particular starting point (a.k.a. algorithm initialization). It's our belief that by initializing the processing of a new data set in a longitudinal series using the processed results from another time point, we can reduce the random variation in the processing procedure and improve the robustness and sensitivity of the overall longitudinal analysis. Such an initialization scheme makes sense also because a longitudinal design is often targeted at detecting small or subtle changes. The longitudinal processing scheme is coded in the script called "recon-all-long", as compared to the standard FreeSurfer pipeline as described in the "recon-all" script. In the following, we list the major differences between the two scripts or the two processing streams. The two scripts will ultimatedly be combined as a single one. In the following, we assume that the longitudinal series have two time-points, tp1 and tp2, and tp1 has already been processed using the standard FreeSurfer recon-all. To process tp2 in the longitudinal way, "recon-all-long" can be used and the command reads like "recon-all-long -long -tp1 tp1 -subjid tp2 ...". "..." represents standard options as used when processing tp1. Difference 1: Register tp2 with tp1. The longitudinal scheme requires aligning the image data of tp2 to tp1. Since both data sets come from the same subject, a linear registration with a dof equal 6 or 9 is enough to get a good alignment between them. The alignment/registration can be computed using the FSL FLIRT tool, and the registration result (a coordinate transformation matrix) is saved and will be used to transfer information from tp1 to tp2 in later steps. As can be seen in the recon-all-long script, the alignment is computed using the following command: |
#pragma section-numbers 2 = Longitudinal Processsing (FS 5.1, 5.2, 5.3 and 6.0) = The largest '''differences''' between this release and previous versions: * Individual time points are '''resampled''' to the base/template space for the -long runs. This reduces variability even further and simplifies many of the algorithms as all input data is in correspondence across time. * The workflow below has not changed! All changes are internal. * '''Manual edits''' are incorporated and in many cases transferred between the cross, base and long runs. * Subject template creation uses cubic spline interpolation in FS 5.2 and later. * Processing single time point data through the longitudinal stream is now possible (usefull for LME where subjects with only a single time point can be included) (FS 5.2 or later). Usefull links: * '''Cite''': [[FreeSurferMethodsCitation]] - please cite our papers if you use this software * '''Edits''': [[LongitudinalEdits]] - how to do edits for longitudinal processing. * '''Tutorial''': [[FsTutorial/LongitudinalTutorial]] - a tutorial on running, post-processing and editing. * '''Statistics''': [[LongitudinalStatistics]] - information on post-processing and analysis. * '''Changes''': [[LongitudinalChangeLog]] - a more detailed change log of the different versions. * '''Data''': [[LongitudinalData]] - links to longitudinal datasets <<TableOfContents(2)>> ---- == Background == Compared with cross-sectional studies, a longitudinal design can significantly reduce the confounding effect of inter-individual morphological variability by using each subject as his or her own control. As a result, longitudinal imaging studies are getting increased interest and popularity in various aspects of neuroscience. The default !FreeSurfer pipeline is designed for the processing of individual data sets (cross-sectionally), and thus not optimal for the processing of longitudinal data series. It is an active research area at the Martinos Center for Biomedical Imaging, how to obtain robust and more reliable cortical and subcortical morphological measurements by incorporating additional (temporal) information in a longitudinal data series. The longitudinal scheme is designed to be unbiased with respect to any time point (TP). Instead of initializing it with information from a specific time point, a template volume is created and run through !FreeSurfer. This template can be seen as an initial guess for the segmentation and surface reconstruction. The !FreeSurfer cortical and subcortical segmentation and parcellation procedure involves solving many complex nonlinear optimization problems, such as the deformable surface reconstruction, the nonlinear atlas-image registration, and the nonlinear spherical surface registration. These nonlinear optimization problems are usually solved using iterative methods, and the final results are known to be highly sensitive to the selection of a particular starting point (a.k.a. algorithm initialization). It is our belief that by initializing the processing of a new data set in a longitudinal series using the processed results from the unbiased template, we can reduce the random variation in the processing procedure and improve the robustness and sensitivity of the overall longitudinal analysis. Such an initialization scheme makes sense also because a longitudinal design is often targeted at detecting small or subtle changes. Additionally to the template new probabilistic methods (temporal fusion) were introduced to further reduce the variability of results across time points. For these algorithms it becomes necessary to process all time points cross-sectionally first. The longitudinal processing scheme is coded in the recon-all script via the flags "'''-base'''" (template creation) and "'''-long'''" (longitudinal runs). == References == '''Please cite at least Reuter et al., 2012 if the longitudinal stream was used!''' '''Example Text:''' <<BR>> ''"To extract reliable volume and thickness estimates, images were automatically processed with the longitudinal stream (Reuter et al., 2012) in !FreeSurfer. Specifically an unbiased within-subject template space and image (Reuter and Fischl, 2011) is created using robust, inverse consistent registration (Reuter et al., 2010). Several processing steps, such as skull stripping, Talairach transforms, atlas registration as well as spherical surface maps and parcellations are then initialized with common information from the within-subject template, significantly increasing reliability and statistical power (Reuter et al., 2012)."'' '''''Highly Accurate Inverse Consistent Registration: A Robust Approach,'''''<<BR>> M. Reuter, H.D. Rosas, B. Fischl.<<BR>> !NeuroImage 53(4), pp. 1181-1196, 2010. http://dx.doi.org/10.1016/j.neuroimage.2010.07.020 <<BR>> http://reuter.mit.edu/papers/reuter-robreg10.pdf '''''Avoiding Asymmetry-Induced Bias in Longitudinal Image Processing,'''''<<BR>> M. Reuter, B. Fischl.<<BR>> !NeuroImage 57(1), pp. 19-21, 2011. http://dx.doi.org/10.1016/j.neuroimage.2011.02.076 <<BR>> http://reuter.mit.edu/papers/reuter-bias11.pdf '''''Within-Subject Template Estimation for Unbiased Longitudinal Image Analysis'''''<<BR>> M. Reuter, N.J. Schmansky, H.D. Rosas, B. Fischl.<<BR>> !NeuroImage 61(4), pp. 1402-1418, 2012. http://dx.doi.org/10.1016/j.neuroimage.2012.02.084 <<BR>> http://reuter.mit.edu/papers/reuter-long12.pdf Also see [[FreeSurferMethodsCitation]] for more !FreeSurfer references. == Workflow Summary == Note: All subjects and all time points are set up in a single SUBJECTS_DIR, which will also contain the base and the longitudinal runs after processing. '''Step 1. cross-sectionally process all time points with the default workflow (tpN is one of the timepoints):''' {{{ recon-all -all -s <tpNid> -i path_to_tpN_dcm }}} '''Step 2. create an unbiased template from all time points for each subject and process it with recon-all:''' {{{ recon-all -base <templateid> -tp <tp1id> -tp <tp2id> ... -all }}} can be started once all norm.mgz files are available from the cross sectional processing of the individual timepoint (step 1). The <templateid> represents a new average template for this subject and the name needs to be chosen by the user. A directory with this name will be created in $SUBJECTS_DIR. This <templateid> needs to be passed to the longitudinal processing of each timepoint of this subject below (step 3). Note, the '...' means that more timepoints can be used as "-tp <tpNid>". It is possible to pass only a single time point (e.g. for subjects that dropped out, see below). '''Step 3. "-long" longitudinally process all timepoints:''' {{{ recon-all -long <tpNid> <templateid> -all }}} The longitudinal processing scheme is coded in the recon-all script via the "-long" flag. These runs will be much faster than the cross sectional or template runs above. This step produces output subject data containing '''<tpNid>.long.<templateid>''' in the name (to help distinguish from the default stream). '''Step 4. Compare results from step 3:''' <<BR>> E.g. calculate differences between <tp1id>.long.<templateid> and <tp2id>.long.<templateid>. Do not compare to the template, as that is merely a blurry guess of where things are located in the specific subject, used only for the initialization of the longitudinal runs. ---- In the following, we describe the two processing streams. We assume that the longitudinal series has time-points tpN (tp1, tp2, ...). All of these have already been processed cross sectionally using the standard !FreeSurfer recon-all (step 1). First we discuss the construction of the template (step 2) and then the longitudinal processing of tpN (step 3). == Creation of Template (recon-all -base) == In step 2 the unbiased template is created for each subject using information from all of its time points. It is possible to start the template construction once the norm.mgz of all time points are available from step 1. After the template is fully processed it is used as an initial guess to initialize many steps in the longitudinal runs (step 3). Some results will not be initialized, but directly copied from the subject-template (brain mask, linear Talairach map, estimated total intracranial volume [[eTIV]]). The assumption is that head size does not change across time. Special situations: * '''Pediatric data''' violates the fixed head size assumption. In the presence of substantial longitudinal growth, the current pipeline may fail. Yet we have found that it is very robust with respect to limited head size differences. In this type of data it is essential to check the template/base image for skull strip errors and meaningful surfaces. Also the eTIV is fixed across time in the final longitudinal directories. In order to get eTIV measures for each time point, read that data from the aseg.stats in the cross sectional directories. * '''Single time-point subjects''': Because of drop-outs it can happen that some individuals were scanned only at a single time point. Instead of removing this valuable information, it is possible (e.g. in the linear mixed effects models LME) to include this data. For that it, however, becomes necessary to run the single time points through the longitudinal stream, to ensure that also those images undergo the same processing steps. This can be done by simply passing a single time point in the base creation step. See "Statistical Analysis of Longitudinal Neuroimage Data with Linear Mixed Effects Models. J.L.Bernal-Rusiel, D.N.Greve, M.Reuter et al. NeuroImage 66:249-260, 2012" for a discussion. The following paragraphs will explain details on what happens inside the recon-all -base. Only the differences to the regular cross sectional recon_all are discussed. Unless you are interested in the detailed changes, you may stop reading here. === Template Initialization (-base-init) === The template is created in the -base-init block of recon-all. This is achieved by [[mri_robust_template]], which constructs an unbiased mean or median (default) norm_template.mgz volume together with the transforms that align each tp's norm.mgz volume with the template (see also [[mri_robust_register]] which is used by [[mri_robust_template]] to construct the maps). The longitudinal scheme later requires aligning the image data of tpN to the template, thus all time points are aligned to an unbiased common space with the command: {{{ mri_robust_template --mov <tpsvols> --lta <tpsltas> --template <templateid>/mri/norm_template.mgz }}} where the input --mov <tpsvols> is a list of the time point's norm.mgz files and the output --lta <tpsltas> a list of the LTA registrations files that take each tp to the template and --template ... norm_template.mgz the median image. The LTA maps are stored in <templateid>/mri/transforms/<tpNid>_to_<templateid>.lta. Also the inverse maps are needed and constructed by {{{ mri_concatenate_lta -invert1 <tpNid>_to_<templateid>.lta identity.nofile <templateid>_to_<tpNid>.lta }}} Since all data sets come from the same subject, these rigid registrations with 6 dof (translation,rotation) are sufficient to get a good alignment between the (intensity normalized) images (i.e. norm.mgz). The registrations and its inverse will be used to transfer information between time points mutually and between time points and the template in the longitudinal stream. The routine [[mri_robust_template]] uses robust statistics to automatically detect outliers and align the rest of the image in an optimal manner. The median template is unbiased with respect to any timepoint and therefore perfectly suited as a template to produce initializations for several steps in the longitudinal runs (see below). After the registrations and the norm_template.mgz volume are created, the orig.mgz images from all TPs are mapped to the template location and averaged (again the default is the median image) to produce the <templateid>/mri/orig/001.mgz image. This image is then processed cross-sectionally with the standard !FreeSurfer stream, mainly to obtain the talairach.xfm, nu.mgz and T1.mgz image, however very early we switch over to the norm_template.mgz. It is possible to check the registrations and template creation either using the norm.mgz or orig.mgz of all time points. These images are first mapped to and resampled in the space of the template with: {{{ mri_convert -rl <templateid>/mri/norm.mgz -at <templateid>/mri/transforms/<tpNid>_to_<templateid>.lta <tpNid>/mri/norm.mgz <tpNid>_to_<templateid>.norm.mgz }}} The resampled files <tpNid>_to_<templateid>.norm.mgz can than be compared with each other or with the template's norm_template.mgz using e.g. tkmedit -f file1 -aux file2. If you want to compare aseg.mgz files, make sure, you use ''nearest'' for the resample type: {{{ mri_convert -rl <templateid>/mri/aseg.mgz -rt nearest -at <templateid>/mri/transforms/<tpNid>_to_<templateid>.lta <tpNid>/mri/aseg.mgz <tpNid>_to_<templateid>.aseg.mgz }}} Except for the following steps, the template is processed with the default cross sectional stream: === Normalization (-normalize) === Output the control points ({{{ctrl_vol.mgz}}}) and bias field ({{{bias_vol.mgz}}}). === Skull Strip (-skullstrip) === The brainmask.mgz are mapped (with nearest neighbor) from the cross sectional runs and averaged. This is similar to a logical OR (union) as it will only exclude voxels that are non-brain in almost all TP. === EM (GCA) Registration (-gcareg) === Use the norm_template.mgz instead of the nu.mgz for the Talairach registration. The resulting registration {{{talairach.lta}}} should be checked as it will be used in the longitudinal stream and will influence all time points. === CA Normalize (-canorm) === Use the norm_template.mgz instead of the nu.mgz for the normalization. This ensures that the norm_template is correctly normalized. ''(internal note: this step might not be necessary, should be checked if norm_template.mgz is already sufficiently normalized)'' ---- |
| Line 19: | Line 144: |
| "fsl_rigid_register -cost corratio -i $tp1dir/mri/orig.mgz -r $subjdir/mri/orig.mgz -o $subjdir/mri/orig_${tp1id}_to_${subjid}.mgz -dof 6 -ltamat $tp1dir/mri/transforms/$regfile" "$regfile" stores the registration as a LTA file. The filename is automatically set inside the script, unless the user changes it using the "-reg" option. Difference 2: Talairach Registration. In the original scheme, the Talairach registration is computed using MINC tools. In the longitudinal scheme, the Talairach registration for tp2 is computed by simply composing the linear registration from tp1 to tp2 and the Talairach registration for tp1. In particular, the following command is used: "mri_concatenate_lta -invert1 -tal $tp1dir/mri/orig.mgz $talairach_template_file $tp1dir/mri/transforms/$regfile $tp1dir/mri/transforms/talairach.xfm $subjdir/mri/transforms/talairach.xfm". "mri_concatenate_lta" is a program for combining two consecutive linear transformation/registration into a single total transformation. "-invert1" indicates that the first transformation matrix need be inverted before composing with the second one. "-tal src_file trg_file" gives the source and target volume names for the first Talairach registration. "$talairach_template_file" is set to "${MINC_BIN_DIR}/../share/mni_autoreg/average_305.mnc" by default. Difference 3: Intensity Normalization. If "-usetp1ctrlvol" is set, the longitudinal scheme will first map the automatically computed control points for tp1 (computed inside mri_normalize) to the data space of tp2, and then use them to perform intensity normalization for tp2. By default, the control points for tp1 are not saved in the original script. Thus, recon-all-long will first regenerate the control points using "mri_normalize -mask $tp1dir/mri/brain.mgz -W $tp1dir/mri/ctrl_vol.mgz $tp1dir/mri/bias_vol.mgz $tp1dir/mri/nu.mgz $tp1dir/mri/T1_tmp.mgz". Then the intensity normalization using tp1's control points (saved as a volume file "ctrl_vol.mgz") is performed as follows: "mri_normalize_tp2 -T1 $tp1dir/mri/nu.mgz -ctrl $tp1dir/mri/ctrl_vol.mgz -mask1 $tp1dir/mri/brain.mgz -xform $tp1dir/mri/transforms/$regfile $subjdir/mri/nu.mgz $subjdir/mri/T1.mgz" We haven't investigated the mapping of manually picked control points from tp1 to tp2. If manual control points are necessary, it's better that the user pick them manually for the tp2 volume too. Difference 4: Skull-stripping. In the longitudinal scheme, the brain-mask for tp2 is borrowed from tp1 (with the tp1-to-tp2 registration transformation properly applied): "mri_mask -xform $tp1dir/mri/transforms/$regfile T1.mgz $tp1dir/mri/$BM $BMA" Difference 5: GCA Registration. Similar to the Talairach registration, in the longitudinal scheme, the linear subcortical atlas registration for tp2 is computed by composing the linear atlas registration for tp1 and the linear mapping from tp1 to tp2: "mri_concatenate_lta -invert1 $tp1dir/mri/transforms/$regfile $tp1dir/mri/transforms/talairach.lta transforms/talairach.lta" "-invert1" again indicates that the first LTA need be inverted before composing with the second. The rationale is easily understandable by drawing the source and target space for each LTA. Difference 6: Canonical Registration (nonlinear). In the standard recon-all, this nonlinear atlas registration is initialized by the previous linear registration result: "mri_ca_register $UseCAAlign -mask brainmask.mgz -T transforms/talairach.lta $xopts norm.mgz $FREESURFER_HOME/average/$GCA transforms/talairach.m3z". In the longitudinal scheme, the nonlinear registration for tp2 is initialized by the nonlinear registration result of tp1, and only local refinements are further performed to get the final registration: "mri_ca_register $UseCAAlign -mask brainmask.mgz -levels 2 -A 1 -l $tp1dir/mri/transforms/talairach.m3z $tp1dir/mri/transforms/$regfile -T transforms/talairach.lta norm.mgz $FREESURFER_HOME/average/$GCA transforms/talairach.m3z". The options "-levels 2 -A 1" tell the program to do only local refinement. "-l $tp1dir/mri/transforms/talairach.m3z $tp1dir/mri/transforms/$regfile" is the longitudinal option that specifies the initial nonlinear registration and the mapping from tp1 to tp2. "-T transforms/talairach.lta" gives the linear registration as in the original scheme. This linear registration will still be used when cross-sequence atlas renormalization is needed (as specified by the $UseCAAlign option). Difference 7: Subcortical Segmentation. In the longitudinal processing scheme, the final subcortical segmentation step again uses the segmentation result from tp1 to initialize the iterative labeling process: "mri_ca_label $UseCAAlign $xopts -l $tp1dir/mri/aseg.mgz $tp1dir/mri/transforms/$regfile norm.mgz transforms/talairach.m3z $FREESURFER_HOME/average/$GCA aseg.auto.mgz". Difference 8: WM surface tessellation and topology correction. The whole process for generating the orig surface that involves WM tessellation, topology correction, and etc is no longer needed in the longitudinal processing of tp2, since the final gray/white surface from tp1 will be used as the initial gray/white surface for tp2. Consequently, the original process is replaced by the following two commands just before running mris_make_surfaces: "mris_transform --dst $subjdir/mri/orig.mgz $tp1dir/surf/${hemi}.white $tp1dir/mri/transforms/$regfile $subjdir/surf/${hemi}.orig", and "mris_transform --dst $subjdir/mri/orig.mgz $tp1dir/surf/${hemi}.white $tp1dir/mri/transforms/$regfile $subjdir/surf/${hemi}.orig_white". That is, the final white surface from tp1 is mapped to the tp2 volume and is used as the orig surface for tp2 and also used as the "orig_white", i.e., the "initial white surface", for the later surface deformation. Similarly, the final pial surface from tp1 is mapped to tp2 and will be used as the initial pial surface for tp2: "mris_transform --dst $subjdir/mri/orig.mgz $tp1dir/surf/${hemi}.pial $tp1dir/mri/transforms/$regfile $subjdir/surf/${hemi}.orig_pial". Difference 9: Make Final Surfaces. In the standard FreeSurfer, the orig surface is deformed to get the white surface, and the white surface is further deformed to get the pial surface. In the longitudinal scheme, the white and pial surfaces are both initialized using tp1's surface deformation results as described above, and the command line becomes: "mris_make_surfaces -orig_white orig_white -orig_pial orig_pial -orig orig_white -long -max 3.5 -mgz -T1 brain.finalsurfs $xopts $subjid $hemi". The options "orig_white" and "orig_pial" specifies the initialization for the white and pial surface respectively. "-long" indicates the longitudinal scheme. Note that in the original scheme, the white and pial surfaces are guaranteed not to cross each other. This non-cross-intersection between the final white and pial surfaces won't be guaranteed in the new scheme since the initial pial surface is independent of the final white surface. If non-cross-intersection is important, the user can remove the "-orig_pial" option and its argument, such that the pial surface is still initialized by the final white surface. But using both "-orig_white" and "-orig_pial" gives better thickness reproducibility. Difference 10: Spherical Mapping. In the longitudinal scheme, the spherical map ($hemi.sphere) for tp2 is simply copied from tp1. It's not hard to see why this is reasonable. Difference 11: Spherical Registration. Like in subcortical segmentation, the nonlinear spherical registration for tp2 in the longitudinal scheme is initialized by the sphere.reg from tp1: "mris_register -curv -nosulc -norot $tp1dir/surf/$hemi.sphere.reg $AvgTif ../surf/$hemi.sphere.reg". "-nosulc" and "-norot" tell the program to ignore the initial linear alignment step as in standard spherical registration. (NOTE: the same principle may be applied for the contra-surface registration for tp2, but it's not currently implemented) Difference 12: Cortical Parcellation. Similar to the subcortical labeling step, the cortical parcellation for tp2 in the longitudinal scheme is again initialized with the parcellation result from tp1: "mris_ca_label $xopts -long -R $tp1dir/label/${hemi}.aparc.annot $subjid $hemi ../surf/$hemi.sphere.reg $CPAtlas $annot" This is done also for the second parcellation using another surface labeling scheme: "mris_ca_label $xopts -long -R $tp1dir/label/${hemi}.aparc.a2005s.annot $subjid $hemi ../surf/$hemi.sphere.reg $CPAtlas $annot" The above is a complete list of all the differences between the standard recon-all and the longitudinal "recon-all-long" processing stream. The new options added to "recon-all-long" include "-long", "-tp1", "-usetp1ctrlvol", "-no-orig-pial", and "-regfile". Their meanings should be clear from the above description. Note that the "recon-all-long" script is modified from a dev version of "recon-all", and may not be consistent with the most recent version of "recon-all". |
== Longitudinal Stream (recon-all -long) == In the -long stream we set {{{NoRandomness = 1}}} (always use same seed for RNG). Affects [[mris_smooth]], [[mris_fix_topology]], [[mris_topo_fixer]], [[mris_sphere]], and [[mris_ca_label]]) === Input (-i) === Copy the {{{orig/00?.mgz}}} from the cross sectionals. === Motion Corrections (-motioncor) === Map 00?.mgz to base space and average there in one step to create {{{orig.mgz}}}. The transforms are available from the cross sectional step (if not, because an older version or the fsl flag was used, they are recreated with mri_robust_register). === NU Intensity Correction (-nuintensitycor) === Same processing as in cross sectional stream. === Talairach (-talairach) === Copy the {{{talairach.xfm}}} from the template (keeps edits). === Normalization (-normalization) === Default: Map and use control points {{{control.dat}}} from the cross sectional, if the file exists (e.g. if manual edits were made to the cross). If the control points file exists in the long run it will not be overwritten (to preserve potential edits done to the long). === Skull Strip (-skullstrip) === Copy the {{{brainmask.mgz}}} from the template to the current TP. Use it to mask the T1.mgz to obtain the final brainmask. Manual edits should be done in the base/template. === EM (GCA) Registration (-gcareg) === Copy the {{{talairach.lta}}} from the template. Edits should be done in base/template. ''(internal note: in the future maybe use the nu.mgz's from all TPs to construct the registration simultaneously (in the -base))'' === CA Normalize (-canorm) === The normalization is initialized with the {{{aseg.mgz}}} of the template copied to the current TP. Thus all TPs use similar control points for the normalization. ''(internal note: in the future maybe find the control points by looking at all nu.mgz simultaneously (in the -base))'' === CA Nonlinear Registration (-careg) === Uses the {{{talairach.m3z}}} from the template as initialization. Also different flags are used: {{{ mri_ca_register -levels 2 -A 1 -l <templateid>/mri/transforms/talairach.m3z <tpNid>/mri/transforms/<templateid>totpN_regfile }}} === CA Nonlinear Registration Inverse (-careginv) === Same processing as in cross sectional stream. === Remove Neck (-rmneck) === Same processing as in cross sectional stream. === EM Registration (with skull but no neck) (-skull-lta) === Same processing as in cross sectional stream. === CA Label (-calabel) === Copy the linear transformation of this TP (cross) to base/template from the template into local transform directory. Then create {{{aseg.fused.mgz}}} by mapping and incorporating segmentation information of all TPs (probabilistic voting). Finally uses the fused aseg as initialization to mri_ca_label to construct the final labels. This will indirectly incorporate aseg edits from the cross sectionals. Furthermore, the intensity scaling factors are passed from the template. ''(internal note: In the future maybe labeling all TPs simultaneously?)'' === Normalization 2 (-normalization2) === Same processing as in cross sectional stream. === Mask Brain Final Surface (-maskbfs) === Changes only concern manual edits: If no brain.finalsurfs.manedit.mgz file exists, check if it exists in the cross sectional run of this time point and map/copy edits from there. === WM Segmentation (-segmentation) === Changes only concern manual edits: If not edited in -long, copy edits and deletions from cross (default) or template (if -uselongbasewmedits is specified). === Cut/Fill (-fill) === Same processing as in cross sectional stream. === Tessellate (-tesselate) === Skip, because we will take initial surface (orig) from template/base anyway. === Orig Surface Smoothing 1 (-smooth1) === Skipped. === Inflation 1 (-inflate1) === Skipped. === QSphere (-qsphere) === Skipped. === Automatic Topology Fixer (-fix) === Skipped. === Final Surfaces (-finalsurfs) === Copy and use {{{?h.white}}} and {{{?h.pial}}} from the template to initialize white, pial and orig surfaces in current TP. This also ensures that all surfaces are indirectly registered (the vertex numbers agree) across all time points. === Surface Volume (-surfvolume) === Same processing as in cross sectional stream. === Orig Surface Smoothing 2 (-smooth2) === Same processing as in cross sectional stream. === Inflate 2 (-inflate2) === Same processing as in cross sectional stream. === ASeg Stats (-segstats) === Same processing as in cross sectional stream. === Spherical Inflation (-sphere) === {{{?h.sphere}}} copied from the template. === Ipsilateral Surface Registation (Spherical Morph) (-surfreg) === Extra flags to [[mris_register]]: {{{ -nosulc -norot <templateid>/surf/$hemi.sphere.reg }}} === Jacobian (-jacobian_white) === Same processing as in cross sectional stream. === Average Curvature (-avgcurv) === Same processing as in cross sectional stream. === Cortical Parcellation (-cortparc) === Extra flags with [[mris_ca_label]]: {{{ -long -R <templateid>/surf/$hemi.aparc.annot }}} === Parcellation Statistics (-parcstats) === Same processing as in cross sectional stream. === Cortical Parcellation 2 (-cortparc) === Extra flags with [[mris_ca_label]]: {{{ -long -R <templateid>/surf/$hemi.aparc.a2005s.annot }}} === Parcellation Statistics (-parcstats) === Same processing as in cross sectional stream. === Cortical Ribbon Mask (-cortribbon) === Same processing as in cross sectional stream. === Add Parcellation to ASeg (-aparc2aseg) === Same processing as in cross sectional stream. === Update WMparc(-wmparc) === Same processing as in cross sectional stream. ---- MartinReuter |
Longitudinal Processsing (FS 5.1, 5.2, 5.3 and 6.0)
The largest differences between this release and previous versions:
Individual time points are resampled to the base/template space for the -long runs. This reduces variability even further and simplifies many of the algorithms as all input data is in correspondence across time.
- The workflow below has not changed! All changes are internal.
Manual edits are incorporated and in many cases transferred between the cross, base and long runs.
- Subject template creation uses cubic spline interpolation in FS 5.2 and later.
- Processing single time point data through the longitudinal stream is now possible (usefull for LME where subjects with only a single time point can be included) (FS 5.2 or later).
Usefull links:
Cite: FreeSurferMethodsCitation - please cite our papers if you use this software
Edits: LongitudinalEdits - how to do edits for longitudinal processing.
Tutorial: FsTutorial/LongitudinalTutorial - a tutorial on running, post-processing and editing.
Statistics: LongitudinalStatistics - information on post-processing and analysis.
Changes: LongitudinalChangeLog - a more detailed change log of the different versions.
Data: LongitudinalData - links to longitudinal datasets
Contents
1. Background
Compared with cross-sectional studies, a longitudinal design can significantly reduce the confounding effect of inter-individual morphological variability by using each subject as his or her own control. As a result, longitudinal imaging studies are getting increased interest and popularity in various aspects of neuroscience. The default FreeSurfer pipeline is designed for the processing of individual data sets (cross-sectionally), and thus not optimal for the processing of longitudinal data series. It is an active research area at the Martinos Center for Biomedical Imaging, how to obtain robust and more reliable cortical and subcortical morphological measurements by incorporating additional (temporal) information in a longitudinal data series.
The longitudinal scheme is designed to be unbiased with respect to any time point (TP). Instead of initializing it with information from a specific time point, a template volume is created and run through FreeSurfer. This template can be seen as an initial guess for the segmentation and surface reconstruction. The FreeSurfer cortical and subcortical segmentation and parcellation procedure involves solving many complex nonlinear optimization problems, such as the deformable surface reconstruction, the nonlinear atlas-image registration, and the nonlinear spherical surface registration. These nonlinear optimization problems are usually solved using iterative methods, and the final results are known to be highly sensitive to the selection of a particular starting point (a.k.a. algorithm initialization). It is our belief that by initializing the processing of a new data set in a longitudinal series using the processed results from the unbiased template, we can reduce the random variation in the processing procedure and improve the robustness and sensitivity of the overall longitudinal analysis. Such an initialization scheme makes sense also because a longitudinal design is often targeted at detecting small or subtle changes. Additionally to the template new probabilistic methods (temporal fusion) were introduced to further reduce the variability of results across time points. For these algorithms it becomes necessary to process all time points cross-sectionally first.
The longitudinal processing scheme is coded in the recon-all script via the flags "-base" (template creation) and "-long" (longitudinal runs).
2. References
Please cite at least Reuter et al., 2012 if the longitudinal stream was used!
Example Text:
"To extract reliable volume and thickness estimates, images were automatically processed with the longitudinal stream (Reuter et al., 2012) in FreeSurfer. Specifically an unbiased within-subject template space and image (Reuter and Fischl, 2011) is created using robust, inverse consistent registration (Reuter et al., 2010). Several processing steps, such as skull stripping, Talairach transforms, atlas registration as well as spherical surface maps and parcellations are then initialized with common information from the within-subject template, significantly increasing reliability and statistical power (Reuter et al., 2012)."
Highly Accurate Inverse Consistent Registration: A Robust Approach,
M. Reuter, H.D. Rosas, B. Fischl.
NeuroImage 53(4), pp. 1181-1196, 2010.
Avoiding Asymmetry-Induced Bias in Longitudinal Image Processing,
M. Reuter, B. Fischl.
NeuroImage 57(1), pp. 19-21, 2011.
Within-Subject Template Estimation for Unbiased Longitudinal Image Analysis
M. Reuter, N.J. Schmansky, H.D. Rosas, B. Fischl.
NeuroImage 61(4), pp. 1402-1418, 2012.
Also see FreeSurferMethodsCitation for more FreeSurfer references.
3. Workflow Summary
Note: All subjects and all time points are set up in a single SUBJECTS_DIR, which will also contain the base and the longitudinal runs after processing.
Step 1. cross-sectionally process all time points with the default workflow (tpN is one of the timepoints):
recon-all -all -s <tpNid> -i path_to_tpN_dcm
Step 2. create an unbiased template from all time points for each subject and process it with recon-all:
recon-all -base <templateid> -tp <tp1id> -tp <tp2id> ... -all
can be started once all norm.mgz files are available from the cross sectional processing of the individual timepoint (step 1). The <templateid> represents a new average template for this subject and the name needs to be chosen by the user. A directory with this name will be created in $SUBJECTS_DIR. This <templateid> needs to be passed to the longitudinal processing of each timepoint of this subject below (step 3). Note, the '...' means that more timepoints can be used as "-tp <tpNid>". It is possible to pass only a single time point (e.g. for subjects that dropped out, see below).
Step 3. "-long" longitudinally process all timepoints:
recon-all -long <tpNid> <templateid> -all
The longitudinal processing scheme is coded in the recon-all script via the "-long" flag. These runs will be much faster than the cross sectional or template runs above. This step produces output subject data containing <tpNid>.long.<templateid> in the name (to help distinguish from the default stream).
Step 4. Compare results from step 3:
E.g. calculate differences between <tp1id>.long.<templateid> and <tp2id>.long.<templateid>. Do not compare to the template, as that is merely a blurry guess of where things are located in the specific subject, used only for the initialization of the longitudinal runs.
In the following, we describe the two processing streams. We assume that the longitudinal series has time-points tpN (tp1, tp2, ...). All of these have already been processed cross sectionally using the standard FreeSurfer recon-all (step 1). First we discuss the construction of the template (step 2) and then the longitudinal processing of tpN (step 3).
4. Creation of Template (recon-all -base)
In step 2 the unbiased template is created for each subject using information from all of its time points. It is possible to start the template construction once the norm.mgz of all time points are available from step 1. After the template is fully processed it is used as an initial guess to initialize many steps in the longitudinal runs (step 3). Some results will not be initialized, but directly copied from the subject-template (brain mask, linear Talairach map, estimated total intracranial volume eTIV). The assumption is that head size does not change across time.
Special situations:
Pediatric data violates the fixed head size assumption. In the presence of substantial longitudinal growth, the current pipeline may fail. Yet we have found that it is very robust with respect to limited head size differences. In this type of data it is essential to check the template/base image for skull strip errors and meaningful surfaces. Also the eTIV is fixed across time in the final longitudinal directories. In order to get eTIV measures for each time point, read that data from the aseg.stats in the cross sectional directories.
Single time-point subjects: Because of drop-outs it can happen that some individuals were scanned only at a single time point. Instead of removing this valuable information, it is possible (e.g. in the linear mixed effects models LME) to include this data. For that it, however, becomes necessary to run the single time points through the longitudinal stream, to ensure that also those images undergo the same processing steps. This can be done by simply passing a single time point in the base creation step. See "Statistical Analysis of Longitudinal Neuroimage Data with Linear Mixed Effects Models. J.L.Bernal-Rusiel, D.N.Greve, M.Reuter et al. NeuroImage 66:249-260, 2012" for a discussion.
The following paragraphs will explain details on what happens inside the recon-all -base. Only the differences to the regular cross sectional recon_all are discussed. Unless you are interested in the detailed changes, you may stop reading here.
4.1. Template Initialization (-base-init)
The template is created in the -base-init block of recon-all. This is achieved by mri_robust_template, which constructs an unbiased mean or median (default) norm_template.mgz volume together with the transforms that align each tp's norm.mgz volume with the template (see also mri_robust_register which is used by mri_robust_template to construct the maps). The longitudinal scheme later requires aligning the image data of tpN to the template, thus all time points are aligned to an unbiased common space with the command:
mri_robust_template --mov <tpsvols> --lta <tpsltas> --template <templateid>/mri/norm_template.mgz
where the input --mov <tpsvols> is a list of the time point's norm.mgz files and the output --lta <tpsltas> a list of the LTA registrations files that take each tp to the template and --template ... norm_template.mgz the median image. The LTA maps are stored in <templateid>/mri/transforms/<tpNid>_to_<templateid>.lta. Also the inverse maps are needed and constructed by
mri_concatenate_lta -invert1 <tpNid>_to_<templateid>.lta identity.nofile <templateid>_to_<tpNid>.lta
Since all data sets come from the same subject, these rigid registrations with 6 dof (translation,rotation) are sufficient to get a good alignment between the (intensity normalized) images (i.e. norm.mgz). The registrations and its inverse will be used to transfer information between time points mutually and between time points and the template in the longitudinal stream. The routine mri_robust_template uses robust statistics to automatically detect outliers and align the rest of the image in an optimal manner. The median template is unbiased with respect to any timepoint and therefore perfectly suited as a template to produce initializations for several steps in the longitudinal runs (see below).
After the registrations and the norm_template.mgz volume are created, the orig.mgz images from all TPs are mapped to the template location and averaged (again the default is the median image) to produce the <templateid>/mri/orig/001.mgz image. This image is then processed cross-sectionally with the standard FreeSurfer stream, mainly to obtain the talairach.xfm, nu.mgz and T1.mgz image, however very early we switch over to the norm_template.mgz.
It is possible to check the registrations and template creation either using the norm.mgz or orig.mgz of all time points. These images are first mapped to and resampled in the space of the template with:
mri_convert -rl <templateid>/mri/norm.mgz -at <templateid>/mri/transforms/<tpNid>_to_<templateid>.lta <tpNid>/mri/norm.mgz <tpNid>_to_<templateid>.norm.mgz
The resampled files <tpNid>_to_<templateid>.norm.mgz can than be compared with each other or with the template's norm_template.mgz using e.g. tkmedit -f file1 -aux file2. If you want to compare aseg.mgz files, make sure, you use nearest for the resample type:
mri_convert -rl <templateid>/mri/aseg.mgz -rt nearest -at <templateid>/mri/transforms/<tpNid>_to_<templateid>.lta <tpNid>/mri/aseg.mgz <tpNid>_to_<templateid>.aseg.mgz
Except for the following steps, the template is processed with the default cross sectional stream:
4.2. Normalization (-normalize)
Output the control points (ctrl_vol.mgz) and bias field (bias_vol.mgz).
4.3. Skull Strip (-skullstrip)
The brainmask.mgz are mapped (with nearest neighbor) from the cross sectional runs and averaged. This is similar to a logical OR (union) as it will only exclude voxels that are non-brain in almost all TP.
4.4. EM (GCA) Registration (-gcareg)
Use the norm_template.mgz instead of the nu.mgz for the Talairach registration. The resulting registration talairach.lta should be checked as it will be used in the longitudinal stream and will influence all time points.
4.5. CA Normalize (-canorm)
Use the norm_template.mgz instead of the nu.mgz for the normalization. This ensures that the norm_template is correctly normalized.
(internal note: this step might not be necessary, should be checked if norm_template.mgz is already sufficiently normalized)
5. Longitudinal Stream (recon-all -long)
In the -long stream we set NoRandomness = 1 (always use same seed for RNG). Affects mris_smooth, mris_fix_topology, mris_topo_fixer, mris_sphere, and mris_ca_label)
5.1. Input (-i)
Copy the orig/00?.mgz from the cross sectionals.
5.2. Motion Corrections (-motioncor)
Map 00?.mgz to base space and average there in one step to create orig.mgz. The transforms are available from the cross sectional step (if not, because an older version or the fsl flag was used, they are recreated with mri_robust_register).
5.3. NU Intensity Correction (-nuintensitycor)
Same processing as in cross sectional stream.
5.4. Talairach (-talairach)
Copy the talairach.xfm from the template (keeps edits).
5.5. Normalization (-normalization)
Default: Map and use control points control.dat from the cross sectional, if the file exists (e.g. if manual edits were made to the cross). If the control points file exists in the long run it will not be overwritten (to preserve potential edits done to the long).
5.6. Skull Strip (-skullstrip)
Copy the brainmask.mgz from the template to the current TP. Use it to mask the T1.mgz to obtain the final brainmask. Manual edits should be done in the base/template.
5.7. EM (GCA) Registration (-gcareg)
Copy the talairach.lta from the template. Edits should be done in base/template.
(internal note: in the future maybe use the nu.mgz's from all TPs to construct the registration simultaneously (in the -base))
5.8. CA Normalize (-canorm)
The normalization is initialized with the aseg.mgz of the template copied to the current TP. Thus all TPs use similar control points for the normalization.
(internal note: in the future maybe find the control points by looking at all nu.mgz simultaneously (in the -base))
5.9. CA Nonlinear Registration (-careg)
Uses the talairach.m3z from the template as initialization. Also different flags are used:
mri_ca_register -levels 2 -A 1 -l <templateid>/mri/transforms/talairach.m3z <tpNid>/mri/transforms/<templateid>totpN_regfile
5.10. CA Nonlinear Registration Inverse (-careginv)
Same processing as in cross sectional stream.
5.11. Remove Neck (-rmneck)
Same processing as in cross sectional stream.
5.12. EM Registration (with skull but no neck) (-skull-lta)
Same processing as in cross sectional stream.
5.13. CA Label (-calabel)
Copy the linear transformation of this TP (cross) to base/template from the template into local transform directory. Then create aseg.fused.mgz by mapping and incorporating segmentation information of all TPs (probabilistic voting). Finally uses the fused aseg as initialization to mri_ca_label to construct the final labels. This will indirectly incorporate aseg edits from the cross sectionals. Furthermore, the intensity scaling factors are passed from the template.
(internal note: In the future maybe labeling all TPs simultaneously?)
5.14. Normalization 2 (-normalization2)
Same processing as in cross sectional stream.
5.15. Mask Brain Final Surface (-maskbfs)
Changes only concern manual edits: If no brain.finalsurfs.manedit.mgz file exists, check if it exists in the cross sectional run of this time point and map/copy edits from there.
5.16. WM Segmentation (-segmentation)
Changes only concern manual edits: If not edited in -long, copy edits and deletions from cross (default) or template (if -uselongbasewmedits is specified).
5.17. Cut/Fill (-fill)
Same processing as in cross sectional stream.
5.18. Tessellate (-tesselate)
Skip, because we will take initial surface (orig) from template/base anyway.
5.19. Orig Surface Smoothing 1 (-smooth1)
Skipped.
5.20. Inflation 1 (-inflate1)
Skipped.
5.21. QSphere (-qsphere)
Skipped.
5.22. Automatic Topology Fixer (-fix)
Skipped.
5.23. Final Surfaces (-finalsurfs)
Copy and use ?h.white and ?h.pial from the template to initialize white, pial and orig surfaces in current TP. This also ensures that all surfaces are indirectly registered (the vertex numbers agree) across all time points.
5.24. Surface Volume (-surfvolume)
Same processing as in cross sectional stream.
5.25. Orig Surface Smoothing 2 (-smooth2)
Same processing as in cross sectional stream.
5.26. Inflate 2 (-inflate2)
Same processing as in cross sectional stream.
5.27. ASeg Stats (-segstats)
Same processing as in cross sectional stream.
5.28. Spherical Inflation (-sphere)
?h.sphere copied from the template.
5.29. Ipsilateral Surface Registation (Spherical Morph) (-surfreg)
Extra flags to mris_register:
-nosulc -norot <templateid>/surf/$hemi.sphere.reg
5.30. Jacobian (-jacobian_white)
Same processing as in cross sectional stream.
5.31. Average Curvature (-avgcurv)
Same processing as in cross sectional stream.
5.32. Cortical Parcellation (-cortparc)
Extra flags with mris_ca_label:
-long -R <templateid>/surf/$hemi.aparc.annot
5.33. Parcellation Statistics (-parcstats)
Same processing as in cross sectional stream.
5.34. Cortical Parcellation 2 (-cortparc)
Extra flags with mris_ca_label:
-long -R <templateid>/surf/$hemi.aparc.a2005s.annot
5.35. Parcellation Statistics (-parcstats)
Same processing as in cross sectional stream.
5.36. Cortical Ribbon Mask (-cortribbon)
Same processing as in cross sectional stream.
5.37. Add Parcellation to ASeg (-aparc2aseg)
Same processing as in cross sectional stream.
5.38. Update WMparc(-wmparc)
Same processing as in cross sectional stream.
