A protocol for ultra-high field laminar fMRI in the human brain

Summary Ultra-high field (UHF) neuroimaging affords the sub-millimeter resolution that allows researchers to interrogate brain computations at a finer scale than that afforded by standard fMRI techniques. Here, we present a step-by-step protocol for using UHF imaging (Siemens Terra 7T scanner) to measure activity in the human brain. We outline how to preprocess the data using a pipeline that combines tools from SPM, FreeSurfer, ITK-SNAP, and BrainVoyager and correct for vasculature-related confounders to improve the spatial accuracy of the fMRI signal. For complete details on the use and execution of this protocol, please refer to Jia et al. (2020) and Zamboni et al. (2020).

SUMMARY Ultra-high field (UHF) neuroimaging affords the sub-millimeter resolution that allows researchers to interrogate brain computations at a finer scale than that afforded by standard fMRI techniques. Here, we present a step-by-step protocol for using UHF imaging (Siemens Terra 7T scanner) to measure activity in the human brain. We outline how to preprocess the data using a pipeline that combines tools from SPM, FreeSurfer, ITK-SNAP, and BrainVoyager and correct for vasculature-related confounders to improve the spatial accuracy of the fMRI signal. For complete details on the use and execution of this protocol, please refer to Jia et al. (2020) and Zamboni et al. (2020).

BEFORE YOU BEGIN Participant screening
Timing: 1 h per participant (1 day before scanning) 1. Participants complete a safety screening form including exclusion criteria related to safety in the MR scanner (~10 min). 2. Participants complete the consent form to ensure that they understand the motivation of the study and any risks/benefits associated with it (~5 min). Then, we introduce our experiment to the participants (~5 min), give them instructions on what they need to do during the scan (~5 min), and ask them to practice the experimental task (~20 min). 3. Participants complete a personal information sheet, including questions related to health (age, medication, visual acuity etc.,~10 min).
Scan setup for brain imaging at 7T (Siemens Terra scanner) Timing: 30 min per participant 4. Turn on the stimulus computer and check the resolution and frame rate are set appropriately based on the experimental design.
5. Check whether scanner triggers and button box responses are registered by the stimulus computer. 6. Check the projection screen is placed in the correct position based on the viewing distance. 7. Ask the participant to change into MR appropriate scrubs and check for items that are not MRcompatible and not safe to bring into the scanner room. Ask the participant to lie on the scanner table in a comfortable position. Participants will require earplugs, pillows for the legs, immobilization cushions, respiratory belt, and a pulse oximeter. Note, this is a list of items / equipment necessary for visual experiments and additional items may be needed for other experiments.
Remind participants not to cross their hands or legs. Explain to participants that they may squeeze the safety ball if they feel uncomfortable.
CRITICAL: Move the table into the scanner bore slowly to avoid the participants experiencing dizziness. If the participant feels dizzy, interrupt the table from moving into isocenter. This will allow the participant to recover from the rapid field change. Once the participant is happy to proceed, move the table towards isocenter.
CRITICAL: Ask the participant to stay as still as possible. Head movement may cause image distortions that cannot be corrected. Therefore, the researcher may communicate with participants using Yes or No questions shown on the screen. Ask the participants to answer these questions by pressing buttons. Remind the participant not to talk unless necessary.
CRITICAL: The Health and Safety Executive limit set for occupational exposure to sound is 85 dB without ear protection (see https://www.hse.gov.uk/noise/employers.htm for details). However, most sequences used in the high field scanner are very loud (around 110 dB). Hence all participants being scanned and those remaining inside the magnet room require ear protection. Earplugs offer up to 37 dB of attenuation with additional attenuation offered by the dense immobilization cushions. The researcher should carefully instruct the participants on the proper use of earplugs, and verify good fit and function of hearing protection in place prior to starting the scan.
8. Adjust the mirror system to check the visible area of the screen before the participant has moved into the scanner bore. Field of view may be limited when using the 32-channel phased array head coil (Nova Medical, Inc., Wilmington, MA, USA, Figure 1) for brain imaging on 7T scanners. If you notice that the participant cannot see enough of the screen, ask them to adjust their head so that they can see as much of the screen as possible, yet ensuring that also the participant's head is placed the furthest back possible inside the head coil.
CRITICAL: If stimuli are not fully visible on the screen, you may need to display the stimuli at a different location by adapting the stimulus presentation scripts. 2. Localizer. A scan is run at isocentre of the magnet in order to have a basic localization of the brain with respect to the scanner coordinates. 3. B1 mapping. The B1 map is acquired covering the whole brain and aligned to AC-PC. This will generate a filp angle map. To ensure the correct flip angle is achieved in the target area, load the flip angle map into the 3D tab of the Syngo software (Siemens Terra, VE system). Use the circle ROI tool (tools menu) to draw a ROI in the target area on sagittal image and estimate the following: (900 3 Transmit Reference Amplitude) / Mean average. This figure should not be over 300 V. If it is higher than 300 V, use 300 V). 4. Acquire anatomical images using MP2RAGE T1-weighted sequence (TR = 5000 ms, TE = 2.56 ms, FOV = 208 3 208 mm 2 , resolution 0.65 3 0.65 3 0.65 mm 3 , number of slices: 240, slice orientation: sagittal, see the MP2RAGE T1-weighted sequence in the key resources table). 5. Position the fMRI slab on the area of interest (2D Gradient Echo, Echo Planar Imaging (GE-EPI) sequence: TR = 2060 ms, TE = 26.4 ms, FOV = 148 3 148 mm 2 , flip angle: 70 , resolution 0.8 3 0.8 3 0.8 mm 3 , number of slices: 56, partial Fourier = 6/8, GRAPPA factor = 3, Multi-Band factor = 2, bandwidth = 1034 Hz/Pixel, echo spacing = 1.09 ms, see the EPI sequence in the key resources table). At this point the researcher will need to: a. Make sure that the area of interest is well contained within the acquisition slab; b. Check the phase encoding direction; c. Slide through the slices to make sure that the area of interest is fully covered by the slab; d. Make sure that there is no wrap-around, ghosting or dropouts in the images. 6. Run one volume of the above scan to check the position of the field of view (FOV). 7. If position is poor, reiterate steps 5 and 6 until good slab coverage has been achieved. 8. Adjust the shim volume (i.e., green box in Figure 2) of the fMRI slab to cover the whole imaging volume (i.e., yellow box in Figure 2). 10 mm were added in all directions of the shim volume to ensure stable numeric optimization and reduce sensitivity to head motion. For imaging at 7T (Siemens Terra scanner), perform B0 shimming as bellow: a. Select Options > Adjustments > Transmitter, in 'amplitude temp', type the value obtained from B1 mapping. Click on Apply. b. Manual Shimming (Siemens Terra, VE line): i. Click on Options>Adjustments. Navigate to the 'Frequency' tab. Click Go until convergence is reached. Click Apply. ii. Repeat this procedure 3 times: Navigate to the 3D Shim tab: measure, calculate and apply. iii. Navigate to the 'Inter. Shim' tab. Adjust the Z, Y, X, Z 2 shim currents as directed, trying to decrease the FWHM of the water spectrum. This value should be as close to 30 Hz as possible (see Figure 3 for good vs. poor shimming examples). Click Stop, then Apply.
CRITICAL: The researcher will need to ensure that the curves resulting from shimming resemble a Gaussian distribution. If the shimming value is not sufficiently low, the participant may need to be repositioned to ensure that the coil receives the signal properly.
9. Measure the B0 map. a. Copy adjustment volume and centre of slices from the fMRI slab (step 5). b. Move slices to obtain full brain coverage and run. 10. Acquire an inverted phase encoding direction fMRI scan (e.g., FH: phase encoding direction of functional data from foot to head) used to correct for distortions (5 volumes): a. Copy adjustment volume from the B0 map. b. Copy Centre of slices with phase encoding direction selected from the fMRI slab (step 5). c. Change the phase encoding direction to FH and the number of measurements to 5. 11. Acquire the experimental data.
a. Copy adjustment volume from the B0 map. b. Copy Centre of slices with phase encoding direction selected from the fMRI slab (step 5). c. Change the phase encoding direction to HF and the number of measurements based on the experimental design. 12. Run steps 10 and 11 as a pair to reduce the between-run head movement. 13. Repeat steps 10-12 for all the experimental runs.
CRITICAL: At this point the researcher will need to ask the participant to stay as still as possible. Head movement may cause image distortions that cannot be corrected.
14. Save all the data acquired.

EXPECTED OUTCOMES
At the end of the process, you will have collected anatomical ( Figure 4A, 50 MB per participant) and functional ( Figure 4B, the size of the EPI data is~2 MB per volume) data. For our experiments, we typically collect 10 runs (~250 volumes per run, that is~5 GB per participant) brain images at high resolution. The anatomical scans (acquired using MP2RAGE T1-weighted sequence) are used to define the cortical layers, while the functional scans (acquired using GE-EPI sequence) are used  to interrogate brain activity underlying cognitive processes (e.g., perceptual learning, adaptation) at a finer scale.

QUANTIFICATION AND STATISTICAL ANALYSIS
We used a pipeline that combines tools from SPM, FreeSurfer, ITK-SNAP and BrainVoyager (troubleshooting 1, also see Table 1 for the function of each file type) to preprocess the brain imaging data and to correct for vasculature-related confounds. Figure 5 provides an overview and flowchart of all analyses.
Note: Disk space: approx. 700 MB per participant White matter-gray matter (WM-GM) segmentation is achieved using FreeSurfer ( Figure 6, troubleshooting 2) and following the steps described below.
1. Convert DICOM files to Nifti file: a. From Matlab, launch SPM toolbox, click on DICOM Import, and select the DICOM files and output directory, click GO. b. From a terminal window, set FreeSurfer environment: i. Export FREESURFER_HOME=/Your Path to FreeSurfer ii. Source $FREESURFER_HOME/SetUpFreeSurfer.sh iii. Export SUBJECT_DIR=/Your Path to FreeSurfer/Your folder name c. Now run recon-all to obtain surfaces and segmentations: i. In the expert.opt file, please set the mris_inflate to 100.
ii. recon-all -all -s OUTPUTFOLDERNAME -hires -i INPUTFILE.nii -expert expert.opt d. Check results in freeview i. CD into the {Subject folder} and on the terminal run freeview -v OUTPUTFOLDERNAME/ mri/T1.mgz (or brainmask.mgz / wm.mgz) e. Convert the results into Nifti format for adjustments in the following steps. From the terminal (while in the OUTPUTFOLDERNAME/mri/ folder): i. mri_label2vol -seg ribbon.mgz -temp rawavg.mgz -o ribbon-in-rawavg.mgz -regheader ribbon.mgz ii. mri_convert ribbon-in-rawavg.mgz ribbon-in-rawavg.nii iii. mri_convert rawavg.mgz rawavg.nii CRITICAL: To save time, the researcher can run these steps in 4-8 participants in parallel. 2. Assign a name, an integer value and a color to each label. It is possible to save the labels and reload them in the next session. 3. Labeling Options a. Foreground label: the label that will be applied by the tool you use (paintbrush or polygon) b. Background label: the label that the foreground will be applied to (if you have foreground label 1 and background label 2, you will DRAW in label 1, but you will only draw over label 2). This prevents you from assigning the same voxel to multiple labels. 4. Select Paintbrush tool (red circle in Figure 7)    Figure 9B). l. Repeat steps a to k for the right hemisphere.

Preprocessing of anatomical data: Generate cortical depth layers using BrainVoyager
Timing: Approx. 1 h per participant.   Figure 10, troubleshooting 3) ii. Save the 3 layers as wholeBrain_layers.voi (in the scanID_VMR folder)

Preprocessing of functional data: Data quality check using BrainVoyager and mcheck toolbox
Timing: Approx. 30 min per participant.
Assess data quality by loading FMRs for each run and visually inspect the images, checking for: 11. Head motion: select Options/Time Course Movie. This generates a movie of the timeseries that is displayed in a loop. Check for large jumps in the images throughout the run that indicate head movement and make a note of the volumes where this occurs. 12. Signal: if you are interested in a specific ROI, you can inspect signal strength and pattern by selecting an ROI around the area you are interested in on one of the slices. The timecourse of the selected area will be displayed. Check for artifacts (e.g., signal drifts, drop outs, spikes). 13. Distortions: to visually assess extent of distortions: load the FMR obtained during the inverse phase encoding scans and toggle between the EPI and IPE images.
Another tool for assessing data quality is the mcheck toolbox. This is a Matlab toolbox that loads the FMR files and evaluates similarity in time-series by computing correlations of image intensity values (Figure 11).   other FMR for intra-session alignment: choose the Single Band Reference of the first EPI run, click OK, and GO. c. Check the results of motion correction. BrainVoyager will output the result of the 3D motion correction for each run. This step aligns each acquired volume to the first volume of the first run, and this is in turn aligned to the Single Band reference of the first run. To assess the result of the motion correction, run the checkMotion.m function of the mcheck toolbox with the output files from the steps above. The result should present improved mean correlation (see Figure 12 for an example This process can be divided into two main steps: initial alignment and fine-tuning adjustment. The former is performed using the 3D Volume Tools / Coregistration process (FMR-VMR alignment). The Figure 11. Example output of the mcheck toolbox (A) and (B) show poor correlation across runs and high correlation within runs. We would expect these values to be within the 0.9-1 interval. (C) represents the results of a k-means clustering algorithm: this works on single volumes rather than whole runs and indicates volumes that are more alike based on image properties.

OPEN ACCESS
STAR Protocols 2, 100415, June 18, 2021 fine-tuning adjustment is performed using Boundary Based Registration (BBR). Note that in order to run BBR, you first need to generate the cortex mesh for the whole brain from the fixed segmentation.

Initial Alignment
Here, a volume (VMR) representation of the respective run of the functional data needs to be created in a space that is close to the native space VMR. This is achieved by running the header-based initial alignment (IA) routine. Besides creating the IA transformation matrix, the header-based step creates and saves a functional volume that is in space close enough to the native space VMR so that the BBR routine is able to find a fine-tuning adjustment (FA) transformation that optimally aligns the two data sets. Here we: a. Run standard header-based initial alignement (IA) from the FMR-VMR Coregistration dialog (3D Volume Tools / Coregistration) b. Select the SBRef of the first run. c. Ensure that Save VMR for BBR is checked in the Initial Alignment tab of the FMR-VMR Coregistration dialog d. As output, the IA will generate and save in the same folder of the SBRef a VMR: scanID_#_un-dist_SCCTBL_3DMCTS_THPGLMFc2-TO-scanID_final_IA_For-BBR.vmr 17. Generate the Cortex Mesh for BBR BBR uses a cortex mesh (SRF) as input which is used to align the respective anatomical data (VMR) with another data set. To obtain a suitable cortex mesh, the underlying anatomical data set (VMR) needs to be segmented. Since the FMR data are aligned usually to the intra-session anatomical data set, this VMR file can be used for the segmentation process. It is also possible to first transform the intra-session VMR into ACPC (or TAL) space, perform standard (advanced) segmentation and afterward transforming the segmented brain back into native space. Here, we compute the segmentation in the native space: a. After having manually adjusted the WM/GM segmentation using ITK-SNAP and saved the file as scanID_segmentation.vmr, load this in BV b. Generate a mask for the WM only and its corresponding Mesh c. From 3D Volume Tools/Segmentation, set Value Range: Min = 0, Max = 149, New = 0, Range; set Value Range: Min = 150, Max = 150, New = 240, Range; do Gaussian filtering Figure 12. Example output of the mcheck toolbox after motion correction After motion correction, the mean spatial correction across volumes and runs are higher than 0.9.
(to reduce noise), then set 18. Fine-Tuning Adjustment using Boundary Based Registration (troubleshooting 5) BBR will use the cortex mesh in native VMR space that will be aligned to the created functional volume VMR: a. Load the created functional data that should be in the same space as the anatomical (mesh) data set (found in the SBRef folder) b. Open the created cortex mesh in native space from the functional volume c. Open the BBR dialog (Meshes/Boundary-Based Registration), click GO to start BBR d. This will produce the spatial transformation FA TRF matrix file.

Correct for superficial bias
It has been shown that the BOLD signal measured using GE-EPI (i.e., T2* weighted) is confounded by macro-and micro-vasculature signals (Uǧurbil, et al., 2003;Uludaǧ et al., 2009;Yacoub et al., 2005). The macro-vasculature contribution is due to veins penetrating the gray matter and running through its thickness, as well as large pial veins situated along the surface of the gray matter (Duvernoy et al., 1981). This results in increased sensitivity (i.e., strong BOLD effect) but decreased spatial specificity of the measured signal. The latter can be understood by the mechanics of the draining veins carrying deoxygenated hemoglobin downstream from the true neuronal site of neural activation, leading to a response spatially biased towards the pial surface, an effect known as superficial bias.
Below we describe our approach in controlling for this superficial bias when acquiring data at high resolution using GE-EPI: 21. Exclude the vein voxels. a. For each voxel in the region of interest, compute the temporal signal to noise ratio (tSNR).
We use tSNR to identify voxels near large veins that are expected to have large variance and low intensity signal due to the local concentration of deoxygenated hemoglobin resulting in a short T2* decay time (i.e., dark intensity in a T2* weighted image). We identify voxels with low tSNR and check their correspondence with voxels of lower intensities on the T2* weighted images. b. For each voxel in the region of interest, perform a GLM (stimulus vs. fixation condition) and compute a t score. It has been shown that high t-values on a fMRI statistical map are likely to arise from large pial veins (Kashyap et al., 2018;Polimeni et al., 2010). c. Exclude voxels with low tSNR values or t-score values above the 90 th percentile of the t-score distribution obtained by the GLM described above from further analysis. 22. Spatial regression analysis (Kok et al., 2016;Koster et al., 2018) a. For each voxel in the superficial layers, find the nearest neighbor in the middle layer. b. Regress out the mean time course of these voxels assigned to middle layers from the time course of voxels assigned to superficial layers. c. Save the time course after the regression for further analyses.

LIMITATIONS
In this protocol, we detail a step-by-step method for using GE-EPI to measure the laminar brain activity in the human cortex. It is important to note that despite the advances afforded by laminar fMRI, GE-EPI is limited by vascular-related contributions to the BOLD signal at the cortical surface resulting in loss of spatial specificity (Kay et al., 2019). Here, we describe possible controls for these potential confounds. Further, CBV imaging using vascular space occupancy (VASO) (Huber et al., 2019) can be used instead of GE-EPI to enhance the spatial specificity of laminar imaging in the human brain.

TROUBLESHOOTING Problem 1: BrainVoyager
In this protocol, we used BrainVoyager for data preprocessing (Quantification and Statistical Analysis). BrainVoyager includes up-to-date algorithms specifically for 7T data analyses (e.g., boundary-based registration, equi-volume depth model). The preprocessed data can be read into Matlab for further analyses. Since BrainVoyager is a commercial software, it may not be available to the general public.

Potential solution
For each step we used BrainVoyager, we listed some free software that can be used instead (Table 2).
Problem 2: FreeSurfer segmentation For some participants, FreeSurfer cannot detect the brain automatically and you may encounter a skull strip error (Preprocessing of Anatomical data: Segmentation using FreeSurfer and SPM).

Potential solution
Use the proton images to generate the brain mask. Use this mask to exclude the noise background, and then run the automatic segmentation again.
Problem 3: Definition of cortical layers For some participants, the whole-mesh cortical depth sampling method may lead to suboptimal geometric sampling (Preprocessing of Anatomical data: Generate Cortical Depth Layers using BrainVoyager).

Potential solution
Try to use the 2D Grids tool of BrainVoyager. Please find more details at https://download. brainvoyager.com/bv/doc/UsersGuide/HighResDataAnalysis/RegularGridCorticalDepthSampling. html. You can also try to use LAYNII instead.

Problem 4: Distortion correction
For some participants, you may still encounter large residual distortions after running the distortion correction using COPE (Preprocessing of Functional data: Distortion correction using COPE).

Potential solution
Use AFNI or FSL TOPUP for distortion correction.
Problem 5: Alignment between functional and anatomical images For some participants, the alignment between functional and anatomical images is not satisfactory (Preprocessing of Functional data: Coregistration between the functional data and the anatomical data using BrainVoyager).

Potential solution
Truncate the functional data (i.e., set all the time series of voxels out of the ROI to 0) to improve the coregistration of the ROI. Use recursive BBR (Van Mourik et al., 2019) or Advanced Normalization Tools (ANTs) to run the alignment.

RESOURCE AVAILABILITY
Lead contact Further information and requests for resources should be directed to the lead contact, Zoe Kourtzi (zk240@cam.ac.uk).

Materials availability
This protocol did not generate any unique reagents.
Data and code availability Further information and requests for the raw datasets and code generated by this protocol should be directed to and will be fulfilled by the lead contact, Zoe Kourtzi (zk240@cam.ac.uk).