React-CR fMRI Analysis Pipeline

🗃️ Phase 1: Pull Data and BIDS Standardize

See here for information on pulling data from CBSCentral

See here for information on BIDS standardization

Step 1.1: Accessing CBSCentral and Remote Desktop

Set up VPN access:
- Use Cisco AnyConnect to access NCF resources:
  - Open Cisco AnyConnect and enter vpn.rc.fas.harvard.edu
  - Submit credentials:
    - RC Username: e.g. elaurent@fasse
    - Password: Your unique string
    - Second Password: 6-digit DuoMobile passcode
Log into CBSCentral:
- Use the following portals:
  - CBSCentral (Old Database)
  - CBSCentral02 (New Database)
- Open a Remote Desktop session to access files and data.
Set up your workspace:
- Open a new terminal window and navigate to project scripts directory:
```
cd /ncf/phelps_lab/users/elaurent/react-CR/scripts
```
- Load and open Sublime Text (for editing scripts):
```
module load ncf sublime_text/3.1.1-ncf
sublime_text
```

Step 1.2: Confirming scan numbers and configurations

Verify scan numbers:
- Confirm scan numbers using the Scan Number Information spreadsheet
  - If no information yet, add participant to spreadsheet referencing the two below:
    - fMRI data log (“scan log” tab)
    - Comprehensive log (“participant visit log” tab)
  - Note: For CBSCentral02, pull the _ND_RMS T1w file, which includes Root-Mean-Square adjustments but does* not*** apply gradient moment nulling (Non-Distorted)
- Update fMRI data log (pull data tab) with relevant information
If needed, create new reactCR_exceptions.yaml (i.e., for sessions that did not run as expected and have repeat or skipped scans)
- Note: I have been creating individualized .yaml exception files for documentation’s sake (e.g. reactCR_{SUBID}_{SESID}.yaml)

Step 1.3: Pulling data with `pull_data.sh`

Open pull_data.sh in Sublime Text and update variables:
- database=cbscentral02 (or cbscentral for data on old CBSCentral)
- ses={SESID}(e.g., 230913_ReactCR_CN021_ses01, as found on CBSCentral)
- yaml_file=reactCR_CBS2.yaml or reactCR_exceptions.yaml (or reactCR_CBS1.yaml for data on old CBSCentral)
Run script in terminal:
```
sh pull_data.sh
```

Step 1.4: Updating filenames with `fixarcget`

Overview: In CBS1, the session names will be saved in long format. For CBS2, the subject name will also be in long format. This will be reflected in the JSON files’ IntendedFor fields as well. Finally, the fieldmap filenames will need to be updated only if they contain “magnitude” instead of “epi”. This step uses the fixarcget step in iteachmri.py to fix the outlined discrepancies.

⚠️ Note:fixarcget can’t update fieldmaps in the same run as JSON session name fields, so if changing both, run in a separate call

Run fixarcget step:

python iteachmri.py --json data_iteachmri.json --steps fixarcget --sub {NEW_SUB} --ses {NEW_SES}

Update flags:
- {NEW_SUB} = short target (e.g., CN021)
- {NEW_SES} = 01, 02, or 03

Enter a response to each prompt:
1. Q: “Do you want to fix the subject name? [y/n]”
  
  A:y
  
  Q: “What is the old subject name?”
  
  A: {OLD_SUB}(e.g., ReactCRCN021)
  
  Q: “What is the new subject name?”
  
  A: {NEW_SUB} (e.g., CN021)
2. Q: “Do you want to fix the session name? [y/n]”
  
  A: y
  
  Q: “What is the old session name?”
  
  A: {OLD_SES}(e.g., 230913ReactCRCN021ses01, which is the CBSCentral session label without any underscores)
  
  Q: “What is the new session name?”
  
  A: {NEW_SES} (e.g., 01, 02 or 03)
3. Q: “Do you want to fix the fieldmap filenames? [y/n]”
  
  A: y or n
  
  ⚠️ Note: Replaces “magnitude” with “epi”. Only to be done if fieldmaps have “magnitude” in filename instead of “epi”
4. Q: “Do you want to fix the field map .json IntendedFor fields? [y/n]?”
  
  A: y
  
  Q: “What is the old session name?”
  
  A: {OLD_SES}(e.g., 230913ReactCRCN021ses01, which is the CBSCentral session label without any underscores)
  
  ⚠️ Note that this won’t work in the same run as changing the fieldmap name, so if changing both fieldmap and .json session name, run in separate calls

Step 1.5: Adding JSON parameters for CBS2 participants

Ensure that PhaseEncodingDirection, TotalReadoutTime, and SliceTiming parameters are added to the .json files for CBS2 participants (Emma: ask Jenn about this step)

Step 1.6: Copying data for quality control (QC)

Copy raw BIDS-formatted data to group directory:

cp -r /ncf/phelps_lab/users/elaurent/react-CR/bids/sub-{SUBID} /ncf/phelps_lab/users/MRI_data/reactCR/raw_data

Update the React-CR fMRI data log (pull data tab)

IGNORE FOR NOW: Additional notes and reminders (ask Jenn)

Ensure NUMECHOS is correctly updated in tasktype_consolidated.csv for multi-echo (ME) acquisitions.
Use the following command to install tedana (requires Python 3.9+) for ME acquisitions:
```
module load ncf miniconda3/py310_22.11.1-1 pip install tedana --no-cache-dir
```
Always double-check for missing fields or metadata errors after data pulls to avoid downstream issues.

🧠 Phase 2: Preprocess Data with fMRIPrep

Step 2.1: Running preprocessing with `fmriprep`

Navigate to project scripts directory (if not already there):
```
cd /ncf/phelps_lab/users/elaurent/react-CR/scripts
```
Run fmriprep step of iteachmri.pyscript:
```
python iteachmri.py --json data_iteachmri.json --steps fmriprep --sub {SUBID}
```
- Update parameters:
  - {SUBID} = e.g., TS007
- ⚠️ Backup method (in event that the above command fails; skip otherwise):
  - Update SUBJ parameter in fmriprep_slurm.sh bash script (if it doesn’t exist, pull from react-ext/scripts/analysisScripts; e.g. SUBJ=CN021) – remember to save the modified file!
  - Submit batch job:
```
sbatch fmriprep_slurm.sh
```
Update fMRI data log (preproc tab) with JobID number
Check progress:
```
sh sbatch_commands.txt
```
Note: Jobs can be identified by their JobID number

⚠️ If you can’t find script in this folder, ask Emma or Jenn – it’s in both of their home directories!

Step 2.2: Reviewing fMRIPrep output reports

Check HTML reports generated by fMRIPrep for detailed summaries of preprocessing [Emma: TO ADD DETAILS]
- Update fMRI data log (preproc tab)

✂️ Phase 3: Prepare Data for SPM

Overview: Due to issues with Nipype’s handling of run concatenation (necessary for certain contrasts), the final stages of preprocessing and first-level analysis are performed by MATLAB. This stage involves converting the fMRIPrep outputs into a format compatible with SPM.

Step 3.1: Verifying the configuration file

Set (cd) working directory to /ncf/phelps_lab/users/elaurent/react-CR/scripts

Double-check the parameters in the data_iteachmri.json file to ensure they match your study design

See here for an example of the expected structure

{
  "baseDir": "/ncf/phelps_lab/users/elaurent/react-CR/",
  "nipypeContainer": "/ncf/phelps_lab/users/elaurent/nipype.sif",
  "fmriprep_spaces": "T1w MNI152NLin2009cAsym",
  "space": "T1w",
  "taskArray": ["MEMORY", "REST"],
  "runArray": [4, 2],
  "volArray": [78, 200],
  "TR": 1.5,
  "conds": [
    ["memory1", "memory2", "memory3", "memory4", "arrows"],
    []
  ],
  "contrasts": [
    [
      ["1: Emotional > Neutral", "T", ["memory1", "memory2", "memory3", "memory4"], [-0.5, -0.5, 0.5, 0.5]],
      ["2: Memory > Arrows", "T", ["memory1", "memory2", "memory3", "memory4", "arrows"], [0.25, 0.25, 0.25, 0.25, -1]]
    ],
    []
  ]
}

Key parameters
- "taskArray": List of tasks included in the analysis.
- "runArray": Number of runs per task.
- "volArray": Number of volumes per run.
- "TR": Repetition time.
- "conds": Conditions for each task.
- "contrasts": Contrast definitions for first-level analysis.
⚠️ Note: Update "runArray" or "volArray" if there are subject-specific exceptions

Step 3.2: Unzipping and splitting files with `fmriprep2spmmatlab`

Run the fmriprep2spmmatlab step in iteachmri.py to convert the fMRIPrep .nii.gz 4D outputs to SPM-compatible .nii 3D formats:
```
python iteachmri.py --json data_iteachmri.json --steps fmriprep2spmmatlab --sub {SUBID} --ses {SESID}
```
- Update flags:
  - {SUBID} = Subject number in the BIDS-formatted volumes (e.g., ID010)
  - {SESID} = Day as it appears in the BIDS-formatted volumes (e.g., 01, 02, or 03)
- What does the fmriprep2spmmatlab step do?
  - Converts data format from compressed .nii.gz files (fMRIPrep output) to expanded .nii files (SPM-compatible)
  - Splits 4D .nii files into individual 3D volumes, one for each TR – this step is crucial for SPM, which requires each time point as a separate file for first-level analysis
Log progress in fMRI data log (first-level tab)

Step 3.3: Creating BIDS-formatted events (.tsv) files

🚨 Emma: [TO ADD] (see R code on Drive and review to automate this process)
Log progress in fMRI data log (first-level tab)

📊 Phase 4: Run First-Level Analysis with SPM

Overview: This stage involves setting up and running first-level fMRI analyses using SPM in MATLAB. The analysis is split into two parts: (1) creating design matrices, and (2) model and contrast estimation.

Step 4.1: Creating design matrices (`_part1`)

Set (cd) working directory to /ncf/phelps_lab/users/elaurent/react-CR/scripts
Load SPM and MATLAB modules
- Using aliases (if set in your .bashrc) run:
```
load_spm load_matlab matlab
```
  ⚠️ Note: Will only work if your .bashrc file includes:
  
  alias load_spm='module load ncf spm/12.7771-linux_x64-ncf'
  
  alias load_matlab='module load ncf matlab/R2022b-fasrc01'
- Alternatively:
```
module load ncf spm/12.7771-linux_x64-ncf matlab/R2022b-fasrc01 matlab
```
  Note*:**This should always work!*
Run Part 1 of the analysis in MATLAB (once per session per subject)
- Make sure MATLAB working directory is set to appropriate /scripts folder
- In the MATLAB Command Window, enter:
```
SPM_firstlevel_reactCR_concat_part1('{SUBID}', '{SESID}', 'sub-{SUBID}_ses-{SESID}_TASK_first-level')
```
  - Update flags:
    - {SUBID} = full subject ID (e.g., ID010)
    - {SESID} = 2-digit session ID (e.g., 03)
- What is SPM_firstlevel_reactCR_concat_part1.m doing?
  - This script sets up and runs the first-level fMRI analysis:
    - Reads event and regressor files
    - Prepares the data (ask Jenn)
    - Specifies the model
    - Runs the job to generate the design matrix
    - Adjusts onset times? (ask Jenn; why 0.5s?):
      - Adjust onset times to account for the run number, add 0.5 seconds, and then divide by TR to convert to scan units.
Update fMRI data log (first-level tab)

Step 4.2: Model and contrast estimation (`_part2`)

Run Part 2 of the analysis in MATLAB (once per session per subject)
- In the MATLAB Command Window, enter:
```
SPM_firstlevel_reactCR_concat_part2('{SUBID}', '{SESID}', 'sub-{SUBID}_ses-{SESID}_TASK_first-level', {CONCAT})
```
  - Update flags:
    - {SUBID} = full subject ID (e.g., ID010)
    - {SESID} = 2-digit session ID (e.g., 03)
    - {CONCAT} = Flag to indicate if runs should be concatenated (1 for yes, 0 for no)
      - ⚠️ WARNING: If re-running, always set the last argument to 0 – Otherwise, you will concatenate files which have already been concatenated!
  - Example:
```
SPM_firstlevel_reactCR_concat_part2('ZM001', '02', 'sub-ZM001_ses-02_TASK_first-level', 1)
```
- What is SPM_firstlevel_reactCR_concat_part2.m doing?
  - Completes the analysis by estimating the model and contrasts; saves the outputs to the outDir argument specified in the call (e.g., /sub-DT002_ses-01_TASK_first-level/)
Update fMRI data log (first-level tab)

Step 4.3: Reviewing the design matrix in SPM

Launch SPM in MATLAB
- In the MATLAB Command Window, enter:
```
spm fmri
```
Check design matrix
- Single click on “Review” in the SPM GUI
- Visually inspect design matrix:
  - If outlier regressors are not visible, check for any preprocessing issues.
  - 🚨 Emma: [TO ADD] (ask Jenn what to else to look for and create guidelines for visual inspection; maybe add example image?)
Update fMRI data log (first-level tab)

🎣 Phase 5: Create ROIs and Extract Mean Timecourses

Overview: This stage involves creating regions of interest (ROIs) based on anatomical masks and generating timecourses of mean activity for each ROI. Additionally, we’ll cover creating a specific vmPFC ROI using peak voxel coordinates from contrasts.

Step 5.1: Updating parameters in `ROI_fishing.sh`

Use the terminal to navigate into the createROIs scripts directory:
```
cd /ncf/phelps_lab/users/elaurent/react-ext/scripts/createROIs
```
Open ROI_fishing.sh in Sublime Text and update the following parameters:

⚠️ Reminder: Make sure to save the updated file!
1. LABEL_NAME and LABEL_NUM: Defining ROIs based on FreeSurfer recon-all outputs
  - Cortical ROIs
  - Left Caudal/Dorsal ACC (dACC_L) = 1002 (aparcaseg)
  - Right Caudal/Dorsal ACC (dACC_R) = 2002 (aparcaseg)
  - Left Rostal/Ventral ACC (vACC_L) = 1026 (aparcaseg)
  - Right Rostral/Ventral ACC (vACC_R) = 2026 (aparcaseg)
  - Left Medial OFC (mOFC_L) = 1014 (aparcaseg)
  - Right Medial OFC (mOFC_R) = 2014 (aparcaseg)
  - Subcortical ROIs
  - Left Amygdala (amyg_L) = 18 (aseg)
  - Right Amygdala (amyg_R) = 54 (aseg)
  - Left Hippocampus (HPC_L) = 17 (aseg)
  - Right Hippocampus (HPC_R) = 53 (aseg)
2. PARCELLATION: Mask type to be used
  - desc-aparcaseg_dseg for cortical ROIs
  - desc-aseg_dseg for subcortical ROIs
3. SUBJARRAY: Subjects you want to iterate through
4. RUNARRAY: Run labels you want to run through (Emma: [TO ADD] examples)

Step 5.2: Running `ROI_fishing_slurm.sh`

Execute the script (ultimately for all subjects, runs, and anatomical ROIs):
```
sbatch ROI_fishing_slurm.sh
```
- Details: This script calls on ROI_fishing.sh and submits it as a batch job to the cluster.
Check and log outputs in fMRI data log (ROI fish tab):
- ROI masks: react-CR/bids/derivatives/sub-{SUBID}/ses-{SESID}/func/
  - Review the ROI masks in fsleyes to confirm that they match expected anatomical regions [Emma: ADD DETAILS]
- Mean timecourses: react-CR/tmcrs_ROI/ (subject-, run-, and ROI-specific) (Emma: [TO ADD] - current script saves it to “ROIs”)
  - Review extracted timecourses to confirm they have the expected number of TRs for each run and are non-zero values