Retain ICA products in output directory

[tsalo]

Currently, fMRIPrep only supports one ICA-based denoising method (AROMA), in which components classified as noise are retained in the confounds tsv. However, there are also plans to support MEICA in fMRIPrep, which involves running ICA and classifying components as noise or signal (much like AROMA). To facilitate both AROMA and MEICA, along with any other ICA-based denoising algorithms one might want to apply, I propose that we output BEP012-compatible versions of the ICA outputs (i.e., the mixing matrix, component maps, and the "component table"/json). (EDIT) I realize the mixing matrix is already outputted, and I assume it's in BEP012 format, but the other elements are not.

This stems from nipreps/fmriprep#1010 (comment) and nipreps/fmriprep#1784 (comment). The general idea was accepted in those issues, but I'd like to discuss the exact settings that should be used, given that AROMA and tedana have slightly different standards for their ICAs.

Open issues:

1. Space
   • tedana assumes native space
   • AROMA assumes standard space
   • Either one could run ICA on another space, retain the component time series, and then estimate component maps in their own required space as needed (see Running ICA on smoothed data ME-ICA/tedana#565).
2. Smoothing
   • tedana assumes unsmoothed data
   • AROMA assumes smoothed data
   • Again, as long as the component time series are valid, either approach can estimate component maps using whatever imaging data they want.
3. Storage of component metrics and classifications
   • tedana stores metrics and classifications in the decomposition json file. We could probably include a prefix in those fields to delineate metric/classification source (e.g., tedana vs. AROMA)
   • I don't think AROMA retains metrics in its current form, but hopefully a pure Python version would be able to do that.

Pinging tedana, AROMA, and BEP012 folks who might be able to weigh in:

• tedana: @CesarCaballeroGaudes, @handwerkerd, @emdupre
• AROMA: @jdkent, @rciric
• BEP012: @effigies

[tsalo]

Relating to point 3 above, what do folks think of this for the decomposition json?

{
    "Method": "Independent components analysis with MELODIC.",
    "ica_00": {
        "Description": "ICA fit to dimensionally-reduced optimally combined data.",
        "AROMA__classification": "rejected",
        "AROMA__rationale": "High HFC",
        "AROMA__maxRPcorr": 0.05,
        "AROMA__HFC": 0.4,
        "AROMA__edgeFract": 0.05,
        "AROMA__csfFract": 0.05,
        "tedana__classification": "rejected",
        "tedana__rationale": "Noise F-value is higher than signal F-value and high variance explained",
        "tedana__countnoise": 1119,
        "tedana__countsigFR2": 6128,
        "tedana__countsigFS0": 3778,
        "tedana__d_table_score": 18.4,
        "tedana__d_table_score_scrub": "n/a",
        "tedana__dice_FR2": 0.5568785491848324,
        "tedana__dice_FS0": 0.243486809765689,
        "tedana__kappa": 35.24615192099909,
        "tedana__kappa_ratio": 2.408554759865709,
        "tedana__normalized_variance_explained": 0.017924586667012896,
        "tedana__rho": 18.026996948546344,
        "tedana__signal-noise_p": 0.28063083660822047,
        "tedana__signal-noise_t": -1.0960543047701905,
        "tedana__variance_explained": 9.670718133027448
    },
    "ica_01": {
        "Description": "ICA fit to dimensionally-reduced optimally combined data.",
        "AROMA__classification": "accepted",
        "AROMA__rationale": "",
        "AROMA__maxRPcorr": 0.05,
        "AROMA__HFC": 0.1,
        "AROMA__edgeFract": 0.1,
        "AROMA__csfFract": 0.1,
        "tedana__classification": "accepted",
        "tedana__rationale": "",
        "tedana__countnoise": 1277,
        "tedana__countsigFR2": 4097,
        "tedana__countsigFS0": 3496,
        "tedana__d_table_score": 14.8,
        "tedana__d_table_score_scrub": 11.8,
        "tedana__dice_FR2": 0.5252645699806231,
        "tedana__dice_FS0": 0.3942153186930905,
        "tedana__kappa": 33.957022187827526,
        "tedana__kappa_ratio": 0.13136512366806724,
        "tedana__normalized_variance_explained": 0.0009146448952444163,
        "tedana__rho": 33.23361928836013,
        "tedana__signal-noise_p": 0.7822498063437031,
        "tedana__signal-noise_t": 0.2769873372743333,
        "tedana__variance explained": 0.508159640579781
    },
}

The component classification algorithm (in this case AROMA or tedana) is used as a prefix before the classification or metric. Each metric has its own entry under each component.

[oesteban]

I realize the mixing matrix is already outputted, and I assume it's in BEP012 format, but the other elements are not.

It is written, but I doubt it's in BEP012 format...

• AROMA assumes standard space

We want to move the decomposition to native space (#1784) so I guess the problem is actually re-implementing the AROMA workflow.

• tedana assumes unsmoothed data

this should not be a problem, are you referring to the consistency between tedana and AROMA?

Relating to point 3 above, what do folks think of this for the decomposition json?

Isn't this more of a BIDS-Derivatives issue? - that said, I'll defer to experts @rciric for this particular proposal.

[CesarCaballeroGaudes]

I don't think re-implementing the AROMA to operate in the standard space should be a problem.

Among the four criteria used in the original ICA-AROMA, the only one that could be 'difficult' is computing the percentage of variance of the spatial map in the ventricles, which in the original code is done based on the template (MNI) space. IMHO, applying the inverse transformation to bring the ventricles template to the original space should do the trick, or even relying on the ventricles segmentation on the T1-weighted image (from freesurfer or FAST).

The other 3 criteria are temporal (based on the spectrum of the IC timecourse, and correlation with expansion of realignment parameters) and percentage of variance of IC map in edge-brain voxels, which only needs a brain mask computed in the original space.

Does this make sense? I haven't checked the one implemented in fMRIprep, is it different from the original ICA-AROMA code. could you point me to the lines of code where both ICA is done.

[tsalo]

AROMA assumes standard space

We want to move the decomposition to native space (#1784) so I guess the problem is actually re-implementing the AROMA workflow.

@oesteban, I've commented on that issue so hopefully we'll be able to get some momentum going there.

Relating to point 3 above, what do folks think of this for the decomposition json?

Isn't this more of a BIDS-Derivatives issue? - that said, I'll defer to experts @rciric for this particular proposal.

I can open an issue in bids-specification (or comment directly on nipreps/fmriprep#519 if @effigies would prefer) if that would be a better place to discuss the structure of the decomposition files.

Does this make sense? I haven't checked the one implemented in fMRIprep, is it different from the original ICA-AROMA code.

@CesarCaballeroGaudes It is the original ICA-AROMA code, for now. I agree that it should be straightforward to apply the methods to native structural-space data (especially with tissue-type information available). If everyone is comfortable pursuing that change to AROMA, then we can commit the fMRIPrep ICA to native structural space.

could you point me to the lines of code where both ICA is done.

I believe that the relevant lines are housed within the fMRIPrep AROMA workflow here:

https://github.com/poldracklab/fmriprep/blob/c717c6b216886be5de75627c067018700f0efd75/fmriprep/workflows/bold/confounds.py#L684-L701

[effigies]

I don't have personal opinions on the ICA. If it works for people who use ICA to use a common MELODIC decomposition, that's fine, but it's not clear that most people will be happy with that.

If we are going to reimplement ICA-AROMA, I would suggest that we do it as an independent workflow or collection of workflows (perhaps a niflow). This would allow the components to be isolated and independently reviewed, and simplify the process of quantifying the effects of the changes you're proposing.

Regarding the decomposition.json, it looks like a table serialized to JSON, so I wonder why we wouldn't target TSVs. The derivatives spec, while a useful tool for standardizing common things, is inherently incomplete, and tools are allowed to produce unspecified derivatives. It makes more sense to me to spend a year or two with things in an unstable but useful state than wedging data into a specified but awkward structure.

[emdupre]

Regarding the decomposition.json, it looks like a table serialized to JSON, so I wonder why we wouldn't target TSVs. The derivatives spec, while a useful tool for standardizing common things, is inherently incomplete, and tools are allowed to produce unspecified derivatives. It makes more sense to me to spend a year or two with things in an unstable but useful state than wedging data into a specified but awkward structure.

I think we've generally aimed to build off of the fMRI derivatives BEP, which currently recommends a decomposition.json file specifically. Though I may have misunderstood how that file was originally envisioned !

Personally, I like having the metadata file in a JSON rather than a TSV, but that's only because it feels more "BIDS-y" to me. Though I'm certainly not the expert in this discussion on that 😸 All that to say, and seconding @oesteban: Maybe we should have this conversation on that specific part of the derivatives BEP ?

[tsalo]

If we are going to reimplement ICA-AROMA, I would suggest that we do it as an independent workflow or collection of workflows (perhaps a niflow). This would allow the components to be isolated and independently reviewed, and simplify the process of quantifying the effects of the changes you're proposing.

Would the niflow wrap a package or would all of AROMA be implemented purely within the workflow?

The metric storage conversation can be continued at bids-standard/bids-specification#519 (comment).

[effigies]

Would the niflow wrap a package or would all of AROMA be implemented purely within the workflow?

My interpretation was that people were interested in rewriting AROMA from scratch as a nipype workflow, which would allow for more flexibility to make other processing decisions regarding space and smoothing. I would assume there would be a mix of FSL interfaces and some pure Python components reimplemented as interfaces, and bundled up into ~3 sub-workflows (pre, main, post). But this is far from a plan anybody needs to follow.

[tsalo]

The metric storage conversation can be continued at bids-standard/bids-specification#519 (comment).

Per maintainers call today, two possible routes forward for decomposition statistics are (1) house everything in jsons (as tedana currently does) or (2) add a new suffix for TSVs containing these stats (e.g., _decompStats.tsv). The former is completely valid under the current rules, and the latter can always be adopted later, so the consensus was to stick with a json file for the immediate future, for both fMRIPrep and the Functional Derivatives BEP.

We'll want to standardize how these jsons are structured though, and since that structure will be fMRIPrep-specific we can discuss it here.

My interpretation was that people were interested in rewriting AROMA from scratch as a nipype workflow, which would allow for more flexibility to make other processing decisions regarding space and smoothing. I would assume there would be a mix of FSL interfaces and some pure Python components reimplemented as interfaces, and bundled up into ~3 sub-workflows (pre, main, post). But this is far from a plan anybody needs to follow.

@effigies That makes sense. I think that would be a good plan going forward, although I would love to see a pure Python implementation (which I think is possible). Hopefully I'll be able to take a crack at this workflow at some point soon.

[tsalo]

This can be transferred to fmripost-aroma, I think.

[tsalo]

We have ICA outputs now. If anyone wants us to change them, please open a new issue.