Advice on structural longitudinal pipeline usage

[annchenknodt]

I'm hoping to use the longitudinal structural pipeline to improve the reliability of 6-year structural change measurements in a midlife sample.

I've successfully run the newly-released pipeline using qunex, but so far the output isn't quite what I expected - primarily, when I compute the reliability of the change measurements in my N=20 test-retest subsample, I find that it is generally quite a bit worse compared to the same measurements extracted from a separate run of longitudinal freesurfer (v6.0). I don't see any red flags in the QC images for the most part, although there is a bit of a pattern on the myelin maps that looks potentially artifactual; it had not appeared when i ran cross-sectional HCP pipeline on the same data in the past.

This all has me suspecting that my inputs to the pipeline might be a little off, and the most likely candidate seems to be that we have supplemented our standard T1 MPRAGE scans at both timepoints with an additional 4 repeats of a compressed-sensing T1 MPRAGE scan at the second time point as an attempt at boosting reliability using a "cluster scanning" approach. My first pass run of the pipeline includes ALL SIX scans in the longitudinal step (so they all go into the template), and I'm concerned that small differences in the acquisition may create an issue. (FWIW, when i compare the extracted measures from my separate run of FS, i find that the grey matter area/volume/thickness measures from the CS MPRAGE scans tend to be slightly smaller than those from our standard MPRAGE scans, but when I rescale them to match and combine the five measurements, I do find increased reliability of the change measure.)

Is it ill-advised to combine slightly different acquisitions in the longitudinal pipeline? I am in the process of re-running without the 4 CS scans and testing some other things to see what I find, but wanted to seek feedback here as well. One of the questions I ran into while getting the pipeline set up with the 4 CS scans included was what to use for the t1samplespacing parameter, as I determined its value to be .00000797 for the standard T1 and 0000082 for the CS scans, but I wasn't sure whether that small difference was likely to affect anything. Happy to provide further info as needed - thanks!

[mharms]

Can you see if you can attach an couple examples of the "potentially artifactual pattern" in the myelin maps?

It might be "ill-advised to combine slightly different acquisitions in the longitudinal pipeline", and it might not. To me, it's an empirical question that will depend on whether the differences induced by the different acquisitions are sufficiently large so as to not be consistent with the assumptions inherent in the longitudinal pipeline.

The t1samplespacing parameter is only relevant if you are doing readout distortion correction (RDC). Are you?

Curious to hear what you learn from the other re-running/testing that you mentioned.

[glasserm]

This is a fairly long post, but I must admit I am not sure what the measurements are that are being compared longitudinally. If these are FreeSurfer measures, we didn't change anything about how FreeSurfer does its measurements. If these are HCP Pipelines measures, which ones are they? There is a small bug currently in the pipeline that the medial wall ROI is not the same on the native meshes across timepoints, but we are fixing this (it is actually a FreeSurfer issue though). One difference between the HCP Pipelines and FreeSurfer is that we use only a single individual to MNI space registration and a single surface registration across all timepoints (otherwise relying on FreeSurfer's rigid alignment of the timepoints). The reason for this is that non-rigid registrations can be affected by local minima related to noise in the features and we want to avoid that.

[annchenknodt]

Thanks to both of you for your comments!

@mharms hopefully this works for attaching an example of the myelin map pattern i'm seeing. i've drawn pink arrows pointing to some of the most prominent instances - they are the most obvious in the left hemisphere at time 1 for the longitudinal version (i see a similar trend in other individuals as well). on the left is the same scan from the run of the cross-sectional pipeline we did several years ago, which corresponds more to how i expect this map to look.

I am doing RDC, but i admit that i haven't done a thorough evaluation of how big of a difference it makes. I'll post an update about the differences I see running the longitudinal pipeline with the CS MPRAGEs omitted once I have that in hand.

@glasserm apologies for being unclear about the measures I'm looking at! Primarily, I am extracting global SA and CT values respectively from:

• /.long.base/T1w/fsaverage_LR32k/.long.base.L.midthickness.32k_fs_LR.surf.gii
• /.long.base/MNINonLinear/fsaverage_LR32k/.long.base.thickness.32k_fs_LR.dscalar.nii
  I'm also looking at subcortical volumes from <SUBJECT>/<SESSION>.long.base/T1w/<SESSION>.long.base/stats/aseg.stats

I certainly have an incomplete understanding of how the HCP pipeline interacts with Freesurfer, though I had thought (and observed) that the Freesurfer values output within the pipeline will differ from a separate run of just Freesurfer, due to the readout distortion correction etc in the preFreesurfer step. Your comments are helpful and I'm working on furthering my understanding so that I can make the most of our new longitudinal data. Thanks!

[glasserm]

Something looks like it is going wrong. We are not, in principle, changing any of FreeSurfers measures but we do changes its inputs so that they have better bias correction, have readout distortion correction (if turned on), are ACPC aligned, and have a better brain extraction.

The longitudinal HCP Pipelines are definitely not thoroughly run on a lot of individuals yet and we may find bugs. Perhaps you should post your calls to all the pipelines. It may be necessary to share us some results to review too.

[mmilch01]

For comparison, I'm attaching QC snapshots for two subjects with three time points each, processed with the longitudinal version at WashU.
Subject1.zip
Subject2.zip

[mharms]

For comparison, I'm attaching QC snapshots for two subjects with three time points each, processed with the longitudinal version at WashU.

I don't see any obvious artifact in the native space ("scene4") QC snapshots of the longitudinal pipeline myelin maps relative to their cross-sectional version in those 2 subjects.

[annchenknodt]

Ah ok, it's good to know that you all agree that the myelin maps don't look quite right.

Here is (I think) all the relevant code I'm using to run the pipeline:

Qunex call to run_recipe:

qdir=/work/long/qunex
con=$qdir/qunex_suite-1.1.0.sif
qunex_container run_recipe \
        --container="${con}" \
	--bind="${qdir}" \
        --recipe_file="$qdir/$projectname/sessions/specs/recipe.yaml" \
        --recipe="hcp_longitudinal" \
	--scheduler="SLURM,jobname=hcp_long,cpus-per-task=2,time=96:00:00,mem-per-cpu=16000,partition=scavenger"

Contents of recipe.yaml:

global_parameters:
  sessionsfolder    : /work/long/qunex/PROJECTNAME/sessions
  sessions             : p45,p52mprage
  overwrite           : "yes"
  batchfile            : /work/long/qunex/PROJECTNAME/processing/batch.txt
  parsessions       : 2

recipes:

  hcp_longitudinal:
      parsessions: 2
      
      commands:
          - create_session_info:
              mapping: /work/long/qunex/PROJECTNAME/sessions/specs/hcp_mapping.txt
          - setup_hcp
          - create_batch:
              targetfile: /work/long/qunex/PROJECTNAME/processing/batch.txt
              paramfile : /work/long/qunex/PROJECTNAME/sessions/specs/parameters.txt              
          - hcp_pre_freesurfer
          - hcp_freesurfer
          - hcp_post_freesurfer
          - hcp_long_freesurfer
          - hcp_long_post_freesurfer

Active options in parameters.txt:

_hcp_suffix                       :
_hcp_brainsize                 : 150
_hcp_t2                            : t2                                            
_hcp_t1samplespacing    : 0.0000080               
_hcp_t2samplespacing    : 0.0000026
_hcp_unwarpdir               : z
_hcp_gdcoeffs                 : NONE
_hcp_avgrdcmethod       : SiemensFieldMap
_hcp_topupconfig           : b02b0.cnf                                     
_hcp_printcom                :
_hcp_sephaseneg           : NONE
_hcp_sephasepos           : NONE
_hcp_seechospacing      : NONE                                          
_hcp_seunwarpdir          : NONE
_hcp_freesurfer_home    : ${FREESURFER_HOME}
_hcp_regname                : MSMSulc
_hcp_grayordinatesres   : 2
_hcp_hiresmesh             : 164
_hcp_lowresmesh          : 32

_hcp_echodiff             : 2.46            # ... the delta in TE times for the hi-res fieldmap image [''].
_run                      : run             # ... Run type: run - do the task, test - perform checks.
_log                      : keep            # ... Whether to remove ('remove') the temporary logs once jobs are completed, keep them in the study level processing/logs/comlogs folder ('keep' or 'study') in the hcp folder ('hcp') or in a <session id>/logs/comlogs folder ('sessions'). Multiple options can be specified separated by '|'.

As my first pass, I've copied over the parameters that worked for us when we originally ran our first time point through the cross-sectional version, but perhaps I'm missing an adjustment I should have made in transitioning to running the longitudinal pipeline within the qunex framework.

Also, I finished running it without the additional CS mprage scans at the second timepoint, so using only the identical standard mprage scans we collected at each timepoint (the first of which we've processed with these parameters previously). This seems to have reduced the apparent "artifact" in a handful of the participants, but left more of them unchanged. Attaching example native space myelin maps for each, where the left column is from the cross-sectional step, middle column is long step with just 2 standard mprages in the template, and right column is long step with the additional CS mprages in the template (but looking at the maps for the standard mprages). The top row is the first time point (where "artifact" is worse fairly consistently) and the bottom row is the second - you can see the middle column has slight improvement in the first example but not the second. So the issue is clearly not limited to my inclusion of these additional scans... I'm also noticing that the pattern of apparent edges also shows up to a lesser extent in the cross-step for the first timepoint (top left images).

I appreciate any thoughts and am happy to share further info!

[glasserm]

Is sharing us a dataset with artifact feasible?

[annchenknodt]

yes it is - what would be the best way to do so privately? thanks!

[glasserm]

Send me an e-mail.

[annchenknodt]

I wanted to come back and post the resolution to this publicly, for anyone who comes across this thread in the future.

Matt took a look at my data and discovered that the issue is most likely due to my use of a FLAIR for surface refinement - since it hasn't been implemented in the longitudinal pipeline, it was getting used as a T2, with suboptimal results. The low resolution of my FLAIR was causing further problems.

Once we wrap up some other projects, I may come back and do a PR to add the FLAIR option myself (Matt said it shouldn't be too difficult), but in the mean time we will stick with just longitudinal freesurfer for now.

Thanks again!