Analyzing fMRI Data 3dDeconvolve <BODY background="00116emb.gif"> <center> <h1> Instructions for Analyzing a fMRI Data Set 3dDeconvolve (Updated August 5, 2005) </h1> </center> <ul> <li>To do this analysis you should have fMRI files with .BRIK, .HEAD and anatomy files <li>this example uses the WWH data set. There are 10 subjects (178 179 181 183 184 185 186 187 188 189),2 stimuli type (pictwhere wordwhere), and 2 stimuli tasks (upper and lower referring to visual fields) <li>You can run this analysis from the command line but to make it easier scripts have been created which use loops to run each command multiple times for each subject,stimuli,task combination required (for more information on how loops work see the Getting Looped tutorial) <li>for each script you will have to change the subject numbers (sub) and the stimulus (stim) to the ones for your data set and save the file <li>to change the file to an executable (commands in this tutorial that are to be run from the command line with begin with a command prompt ex.[fMRI@wilder]$)<br> <blockquote>[fMRI@wilder]$ chmod a+x script_name </blockquote> <li>you can use all the scripts by copying the text into gedit </ul> <h3> Deconvolve </h3> <ol> <li> Preparing the Stimuli Files <ul> <li> The purpose of this step to make text files that tell the deconvolve program when the desired stimuli task turns on and off ordered according to the time course of the experiment, for each subject/stimuli combination (<a href=183_pictwhere.1D>see example text file</a>) <li> The easiest way to do this is to open the E-prime script or E-prime data aid files for each subject and copy the information referring to when a stimulus is turned on/off in to an excel file (<a href=WHERE_PICT.xls> see example xls file</a>) and then copy the only two stimuli columns into a text file. Remember to exclude the first two relaxed state volumes (could be different than two depending on the experiment). <li> 1=ON 0=OFF <li> Stimuli are numbered down the column and each column represents a different stimulus (The first two time sets are normally discarded.) <li> In our example, there is an upper (column one )and a lower visual field condition (column two), so in our example text file <blockquote> 1 0 = The stimuli is in the upper visual field <br> 0 1 = The stimuli is in the lower visual field <br> 0 0 = The stimuli is in the middle visual field or basline readings are being taken </blockquote> <li> Rename the text file to include a .1D extension so that it can be used in the deconvolve script. <li> Please note to be very careful about the accuracy of your stim files. It can be very easy to make mistakes and mess up the whole analysis. Make sure you have the same number of volumes in your stim files as your fMRI data. This is done by going into your subject data and find your 3d time series data for your condition. If you choose overlay, they should be the ones that are named epan:3D+t:{volumenumber}. That volume number should match the number of row from your stim files. </ul> <li> Preparing the data for deconvolving <ul> <li>Run the deconvolve_rename_code.sh script to rename the .orig files and copy them into the deconvolve directory <blockquote> #!/bin/tcsh <br> foreach sub ( 183 184 185 186 187 188 ) <br> foreach stim ( pictwhere wordwhere ) <br> echo $sub $stim <br> cp ./${sub}/${stim}?_000??+orig.HEAD ./deconvolve <br> cp ./${sub}/${stim}?_000??+orig.BRIK ./deconvolve <br> mv ./deconvolve/${stim}?_000??+orig.BRIK ./deconvolve/${sub}_${stim}+orig.BRIK <br> mv ./deconvolve/${stim}?_000??+orig.HEAD ./deconvolve/${sub}_${stim}+orig.HEAD <br> end <br> end </blockquote> <li> Run the sag_rename_code.sh to rename and copy the anatomical files into the deconvolve directory <li> This is done to put all the files in one directory so a script can be run. The anatomical orig and talairach (both need to be there) files for each subject need to be copied in a deconvolve directory. <li> If the maps have not been talairached then you can still do them later and just copy the orig anatomicals. This script is placed in the root directory (where you see all the directories of you subjects). <blockquote> #!/bin/tcsh <br> set filetype=orig <br> set cut=sagittal <br> foreach sub ( 183 184 185 186 187 188 ) <br> echo $sub <br> cp ./${sub}/${cut}_seT1_00???+${filetype}.BRIK ./deconvolve <br> cp ./${sub}/${cut}_seT1_00???+${filetype}.HEAD ./deconvolve <br> mv ./deconvolve/${cut}_seT1_00???+${filetype}.BRIK ./deconvolve/${sub}_${cut}+${filetype}.BRIK <br> mv ./deconvolve/${cut}_seT1_00???+${filetype}.HEAD ./deconvolve/${sub}_${cut}+${filetype}.HEAD <br> end </blockquote> <li> Please note that to copy the talairach versions, put tlrc instead of orig in the filetype. To copy coronal, axials then set those in the cut variable. </ul> <li> Deconvolving the dataset <ul> <li> Run the 3ddeconvolve.sh script in the deconvolve directory <br> <li> [] refers to the column in your text (renamed .1D) file referring to a specific stimuli starting at [0] as column 1. <blockquote> #!/bin/tcsh <br> foreach sub ( 183 184 185 186 187 188 ) <br> foreach stim ( pictwhere wordwhere ) <br> 3dDeconvolve \ <br> -input ${sub}_${stim}+orig \ <br> -nfirst 6 \ <br> -polort 1 \ <br> -xout \ <br> -num_stimts 2 \ <br> -stim_file 1 "${stim}.1D[0]" -stim_file 2 "${stim}.1D[1]" \ <br> -stim_label 1 ${sub}_${stim}_upper_block -stim_label 2 ${sub}_${stim}_lower_block \ <br> -stim_maxlag 1 5 -stim_maxlag 2 5 \ <br> -iresp 1 ${sub}_${stim}_upper_block_irf -iresp 2 ${sub}_${stim}_lower_block_irf \ <br> -tshift \ <br> -sresp 1 ${sub}_${stim}_upper_block_irf_sd -sresp 2 ${sub}_${stim}_lower_block_irf_sd \ <br> -bucket ${sub}_${stim}_bucket -tout -fout -rout <br> end <br> end </blockquote> <dl> <dt> Where <dd> -nfirst ${volume} sets the first volume <dd> -input sets the input .orig file <dd> -polort 1 sets the function to a linear regression <dd> -xout checks the resulting output matrix, if it is 1x1 then the maxlag has to be changed <dd> -num_stimts ${number of stimuli tasks} sets how many stimuli will be analyzed <dd> -stim_file ${task} ${sub}_${stim}.1D[${column}] defines the text file created in Step 1 <br> use column 0 for stimulus 1 and column 1 for stimulus 2 <dd> -stim_label ${stimulitask} ${sub}_${stim}_upper_block defines the label for the output file <dd> -stim_maxlag ${stimulitask} ${length of function} defines the desired length of the impulse response function (details on the IRF in appendix in end) <dd> -iresp ${stimulitask} ${sub}_${stim}_upper_block_irf outputs the impulse response function <dd> -tshift sets time correction <dd> -sresp ${stimulitask} ${sub}_${stim}_upper_block_irf_sd outputs the standard deviation <dd> -bucket ${sub}_${stim}_bucket sets the filename for the output bucket that contains all of the output <dd> -tout outputs t stats <dd> -fout outputs the F stats <dd> -rout outputs least squares </dl> </ul> <br> <li> Find the highest point of the impulse response function <ul> <li> The highest intensity found in the impulse response function will be used in the difference calculation. You will need to do this for each subject/stimuli combination <ol type=a> <li> Run afni <li> Switch Underlay (also called anatomicals) to an impulse response function (${stimulitask} ${sub}_${stim}_${block}_irf) <br> Switch Overlay (also called function map) to the bucket for the same subject <li> Go to Define Datamode > Misc > Overlay info to view contents of bucket file <br> Record with subbricks contain the R^2 and F-stat for each visual field (upper and lower in this case) <li> Go to Define Overlay and set Overlay = the sub-brick of R^2 and Thr = sub-brick of F-stat <br> Open the axial view of the current map <br> Set See Overlay and set the exponential to 1 (drop down box below t-stat) <br> Raise the t-value to clean up the image until only a few active voxels are seen <li> Highlight an active voxel and open the graph page <br> This is the impulse response function, click on the point where the graph is at its highest. This is the highest impulse intensity<br> There should be x + 1 points on the graph where x is the length of the impulse response (specificed by -stim_maxlag) <br> Counting from left to right starting at 0 record the time lag (index = time lag ) <li> View the contents of the bucket file again and record the sub-brick that corresponds to the time lag of the highest intensity and the Fstat for each stimulus (in our case for upper_block and lower_block) <br> <li> Also you need to record F-values for the whatever p-values you will want to use as you won't be using ETA values but actually F-values in your analysis step. This is done by going to Define Overlay and moving the t value bar until the p-value (bottom number) is right. Record the value shown beside your bar (this is your F-value). <br> The time lag is listed after each sub-brick in []. <li> If you haven't already talairach the anatomicals then use AFNI to talairach the sag files <br> <li> Switch Underlay to a talairached sagittal file and the Function to a bucket for the same subject Go to Define Datamode > Write (many) and select the buckets for the subject whose anatomical you selected to write out tlrc bucket file <br> Do this for each subject/stimuli bucket. (note that in order for the writing to work both orig and tlrc of the sag must be in the deconvolve directory and remember to click on talairach view first before going to Define Datamode) </ol> </ul> <br> <li> Analyze the results (This ends the deconvolve step and starts the Analysis step) <ol type=a> <li> Note that the analysis proceed like a regular difference maps except that you have to specify certain sub-bricks for the noneg and Difference calculation scripts and that you use a F-value instead of a ETA value. These scripts also differ slightly to reflect the differences in file naming. <li> Specific details and explanations of each step can be found in the Analyzing fMRI Data using Logical AND, Logical OR and Exclusive OR or <a href ="logical_analysis.html">here</a>. <li> Remove the negative intensities with with volcode_noneg.sh script <li> The sub-brick number enclosed in [] is one recorded in the previous step the matches the time lag of the highest intensity (Step 4f) <li> This script also is set for 2 stimuli type and tasks. If you have more tasks then copy the set task and 3dmerge line and paste in their corresponding sections in the script while changing the number on the task variable. (And don't forget to change the sub-brick corresponding to that stimuli task). This is run in the deconvolve directory <blockquote> #!/bin/tcsh <br> set task1=UPPER <br> set task2=LOWER <br> foreach stim ( pictwhere wordwhere ) <br> foreach sub ( 183 184 185 186 187 188 ) <br> echo $stim $sub <br> 3dmerge -1noneg -prefix ${sub}_${stim}_${task1}_noneg "${sub}_${stim}_bucket+tlrc[8]" <br> 3dmerge -1noneg -prefix ${sub}_${stim}_${task2}_noneg "${sub}_${stim}_bucket+tlrc[22]" <br> end <br> end </blockquote> <li> Calculate the difference maps using the script below <li> Output will be stored in the directory specified by the directory variable and note that value a and b in the below calculation should correspond to the sub-brick of your F-stat brick of each stimuli. Change the F-value that you found for the corresponding p-value in Step 4g. When running the script, do not be alarmed that input c and d are not used as they aren't in the NOT maps. This script is again run in the deconvolve directory. <blockquote> #!/bin/tcsh <br> set Fvalue=3.046 <br> set directory=./calcsFvalue${Fvalue} <br> mkdir ${directory} <br> foreach stim ( pictwhere wordwhere ) <br> foreach sub ( 183 184 185 186 187 188 ) <br> <br> #Calculate the first set of subtraction and intersection maps <br> echo starting 3dcalc on subject $sub <br> set task1cap=UPPER <br> set task2cap=LOWER <br> echo Calculating LOGICAL OR MAP <br> <br> 3dcalc -a "${sub}_${stim}_bucket+tlrc[17]" -b "${sub}_${stim}_bucket+tlrc[31]" -c "${sub}_${stim}_${task1cap}_noneg+tlrc[0]" -d "${sub}_${stim}_${task2cap}_noneg+tlrc[0]" -prefix ${sub}_${stim}_${task1cap}subLogOr${task2cap}_Fvalue${Fvalue} -session ${directory} -expr "((step(a-${Fvalue})*c)-(step(b-${Fvalue})*d))" <br> <br> echo Calculating EXCLUSIVE OR MAP <br> <br> 3dcalc -a "${sub}_${stim}_bucket+tlrc[17]" -b "${sub}_${stim}_bucket+tlrc[31]" -c "${sub}_${stim}_${task1cap}_noneg+tlrc[0]" -d "${sub}_${stim}_${task2cap}_noneg+tlrc[0]" -prefix ${sub}_${stim}_${task1cap}subExcOr${task2cap}_Fvalue${Fvalue} -session ${directory} -expr "(((step(a-${Fvalue})*c)-(step(b-${Fvalue})*d))*(1-(step(a-${Fvalue})*step(b-${Fvalue}))))" <br> <br> echo Calculating INTERSECTION MAP <br> <br> 3dcalc -a "${sub}_${stim}_bucket+tlrc[17]" -b "${sub}_${stim}_bucket+tlrc[31]" -c "${sub}_${stim}_${task1cap}_noneg+tlrc[0]" -d "${sub}_${stim}_${task2cap}_noneg+tlrc[0]" -prefix ${sub}_${stim}_${task1cap}int${task2cap}_Fvalue${Fvalue} -session ${directory} -expr "((step(a-${Fvalue})*step(b-${Fvalue}))*((c+d)/2))" <br> <br> echo Calculating NOT MAP 1 <br> <br> 3dcalc -a "${sub}_${stim}_bucket+tlrc[17]" -b "${sub}_${stim}_bucket+tlrc[31]" -c "${sub}_${stim}_${task1cap}_noneg+tlrc[0]" -d "${sub}_${stim}_${task2cap}_noneg+tlrc[0]" -prefix ${sub}_${stim}_${task1cap}not${task2cap}_Fvalue${Fvalue} -session ${directory} -expr "((step(a-${Fvalue})-(step(a-${Fvalue})*step(b-${Fvalue})))*c)" <br> <br> echo Calculating NOT MAP 2 <br> <br> 3dcalc -a "${sub}_${stim}_bucket+tlrc[17]" -b "${sub}_${stim}_bucket+tlrc[31]" -c "${sub}_${stim}_${task1cap}_noneg+tlrc[0]" -d "${sub}_${stim}_${task2cap}_noneg+tlrc[0]" -prefix ${sub}_${stim}_${task2cap}not${task1cap}_Fvalue${Fvalue} -session ${directory} -expr "((step(b-${Fvalue})-(step(b-${Fvalue})*step(a-${Fvalue})))*d)" <br> end <br> end </blockquote> <dl> <dt> Where <dd> -a is the threshold for stimulus 1 (or upper in this example) (Step 4f Fstat sub-brick) <dd> -b is the threshold for stimulus 2 (or lower in this example) (Step 4f Fstat sub-brick) <dd> -c is the intensity for stimulus 1 that has had negative intensities removed <dd> -d is the intensity for stimulus 2 that has had negative intensities removed </dl> </ol> <br> <li> Proceed with analysis the same as in the Logical Analysis at step 3 of Creating Subtraction and Intersection Maps given <a href ="logical_analysis.html">here</a> <ol type=a> <li> Use the deconvolve_cluster.sh to cluster the data <li> Use the deconvolve_blur.sh to blur the data <li> Use the deconvolve_anovameans.sh to calculate the means <li> Use the deconvolve_dup_means.sh to duplicate the means <li> Use the deconvovle_float.sh to change the data into float format </ol> </ol> <h3> Appendix </h3> <ol> <li> Explanation of the IRF (quoted from Gord Sarty) <blockquote> The IRF is how a dynamical system (one that can be described by differential equations - which is everything physical) responds to an ideal instantaneous "impulse" input function. Knowing the IRF is mathematically equivalent to knowing the differential equations for the system if the system is linear. The brain is almost surely not linear but another mathematical fact is that for short intervals, every system is approximately linear (that's the basis of calculus). So we assume a linear model for the hemodynamic response and hope that it's close. The response, let's call it HRF, to a real input, let's call it F, is the convolution between the IRF and F. Using * to denote convolution this is <br><br> HRF = F*IRF<br><br> The fMRI time-series is the HRF and we have the input F (the stim file) so we can deconvolve to get the IRF. Deconvolution is mathematically tricky - AFNI uses a statistical method (within the class of general linear models - GLM) to do the deconvolution but we must nail done the length of the IRF to define the GLM. After that AFNI finds the best, in a least squares sense, IRF of the pregiven length that will reproduce the measured HRF (as per the above formula).<br><br> </blockquote> </ol> Back to <a href =" index.html">U of S fMRI web page</a>. </HTML>