Contents
- Introduction
- User Guide
- FAQ
Running Dual Regression
Run MELODIC on your group data in Concat-ICA mode ("Multi-session temporal concatenation"). Find the file containing the ICA spatial maps output by the group-ICA; this will be called something like melodic_IC.nii.gz and will be inside a something.ica MELODIC output directory.
Use Glm (or any other method) to create your multi-subject design matrix and contrast files (design.mat / design.con).
- Run dual_regression. Just type the script name to get the usage - should be mostly self-explanatory:
The 4D group spatial IC maps file will be something like somewhere.ica/melodic_IC
The des_norm option determines whether to variance-normalise the timecourses created by stage 1 of the dual regression; it is these that are used as the regressors in stage 2. If you don't normalise them, then you will only test for RSN "shape" in your cross-subject testing. If you do normalise them, you are testing for RSN "shape" and "amplitude".
One easy way to get the list of inputs (all subjects' standard-space 4D timeseries files) at the end of the command is to use the following (instead of listing the files explicitly, by hand), to get the list of files that was fed into your group-ICA: `cat somewhere.gica/.filelist`
Explanation of outputs
dr_stage1_subject[#SUB].txt - the timeseries outputs of stage 1 of the dual-regression. One text file per subject, each containing columns of timeseries - one timeseries per group-ICA component. These timeseries can be fed into further network modelling, e.g., taking the N timeseries and generating an NxN correlation matrix.
dr_stage2_subject[#SUB].nii.gz - the spatial maps outputs of stage 2 of the dual-regression. One 4D image file per subject, and within each, one timepoint (3D image) per original group-ICA component. These are the GLM "parameter estimate" (PE) images, i.e., are not normalised by the residual within-subject noise. By default we recommend that it is these that are fed into stage 3 (the final cross-subject modelling).
dr_stage2_subject[#SUB]_Z.nii.gz - the Z-stat version of the above, which could be fed into the cross-subject modelling, but in general does not seem to work as well as using the PEs.
dr_stage2_ic[#ICA].nii.gz - the same as the PE images described above, but reorganised into being one 4D image file per group-ICA component, and, within each, having one timepoint (3D image) per subject. This reorganisation is to allow stage 3, the cross-subject modelling for each group-ICA component - so it is these files that would normally be fed into randomise.
dr_stage3_ic[#ICA]_tstat[#CON].nii.gz - the output of "stage 3", i.e. files created by running randomise, doing cross-subject statistics separately for each group-ICA component. You'll get one set of statistical output files per group-ICA component, and, within that set of statistical output files, one t-stat (etc.) per contrast in the cross-subject contrast file (design.con). The corresponding corrected (1-p) p-value images are called *corrp*.
Multiple-comparison correction across all RSNs
The need for correction, and correction via Bonferroni
Warning: The corrected p-values output by the final randomise (*corrp*) are fully corrected for multiple comparisons across voxels, but only for each RSN in its own right, and only doing one-tailed testing (for t-contrasts specified in design.con). This means that if you test (with randomise) all components found by the initial group-ICA, and you do not have a prior reason for only considering one of them, you should correct your corrected p-values by a further factor. For example, let's say that your group-ICA found 30 components, and you decided to ignore 18 of them as being artefact. You therefore only considered 12 RSNs as being of potential interest, and looked at the outputs of randomise for these 12, with your model being a two-group test (controls and patients). However, you didn't know whether you were looking for increases or decreases in RSN connectivity, and so you ran the two-group contrast both ways for each RSN. In this case, instead of your corrected p-values needing to be <0.05 for full significance, they really need to be < 0.05 / (12 * 2) = 0.002 !
Alternative: correction via combining statistics
An alternative to this full Bonferroni correction, that might be more sensitive (while still giving full control of overall family-wise-error, is to do the following "combining" statistic approach. For each voxel, this just keeps the maximum test statistic value across all RSNs, and builds up a null distribution of that. This mini-script needs running in sh (or bash):