De Profundis

Since school started back up a little over a month ago, I've been teaching my own introductory psychology course - the first one I've ever taught. It's been great; when I do it, I feel happy, energized, useful. Not only do I get to teach a couple of hundred persons, some of whom really like the class, but I even get to write my own test questions. If ten years ago you told me I would be lecturing in front of college students and creating exams, I would've told you to go kick rocks. But oddly enough, I enjoy it. For instance, I get to come up with test questions like this:

Sigmund Freud was:

A) correct
B) completely insane
C) a cocaine abuser
D) all of the above


Understandably, the writing here took a hit, as I left the mistress of my youth, blogging, for the wife of mine age, teaching. I figured that whatever basics someone needed to know they could find somewhere in the archives. I thought that I probably was spending too much time on it anyway, feeling a distinct negative correlation between the amount of material that was posted on here and my actual productivity as a researcher. Maybe once in while, I thought, I would post about relationship advice, or books I was reading, or possibly a sequel to my post about erotic neuroimaging journal titles. In any case, though, I thought that I written most of what I had wanted to. Quod scripsi, scripsi.

However, I got a jolt a couple of days ago when the government shut down. Normally this wouldn't affect me all that much, except for the AFNI developers sending out a warning that the entire website would be shut down and the message boards closed until further notice. Apparently it's illegal for federal employees to use government emails for correspondence during a shutdown. No more AFNI support; no more new AFNI materials; no more infamous Bob Cox insults, put-downs, and slams. All gone.

Which made me realize, God forbid, that AFNI won't be around forever. I doubt that even fMRI will be all that popular twenty, thirty years from now, as technology progresses and we find some new way to image the brain yet still get significant results no matter how badly we screw up. But the need for understanding the processes will still be there; people will always want easy, accessible introductions for learning new tools and new methods. So I will continue to do that, no matter how incomplete my understanding, no matter how many miles I run, no matter how many jars of Nutella I devour. And that's a campaign promise.

At least one researcher out there has offered to write about topics that I am either entirely ignorant of, or only have a limited grasp. Starting soon, I will begin to interleave some of those posts with what I write to make the blog more centralized for information relating to fMRI. In the meantime, if anyone reading wants to contribute, just let me know. We'll see how it works out.

Real Ultimate Power

I once read an article by a man who guaranteed that by following a series of simple steps, your blog would become famous overnight and make you filthy rich. He pointed to himself as a case study, having founded a blog devoted to nutrition advice and dietary tips for Caucasian males. He called it White Man's Secret. Within weeks, he claimed, he had millions of unique hits every single day from all over the world. Tens of thousands of women left their husbands or lovers and became his willing slaves. He was elected president of his Rotary chapter and was voted Best Table Topics Speaker at his local Toastmasters five times in a row. And he had done it all by following a series of simple steps that anybody could do.

Most important, he said - or screamed, rather, since everything was in upper case capital letters - was to never apologize for failing to update regularly. People wouldn't even notice it, he said, unless you called attention to it. He recommended treating readers like any one of the members of his now-overflowing seraglio - by playing with their emotions through deliberate neglect. Then, when they were at their most desperate, you would give them some attention. Not much. Just enough for them to become dependent on you.

This, he said, was an expression of power. Real Ultimate Power.

________________________________________

My reason for neglect isn't nearly as insidious. It's because I've been running too much.

A few weeks ago when I visited home I sat down with my dad and a calendar and we made a training program leading up to the Indianapolis marathon in November. At the beginning of each week he wrote the number of weeks left until race day. So, "12," "11," "10," and so on, until "0." Then in each cell of the calendar we calculated how many workouts and how many miles I would need to run to fill a weekly quota. We determined this using Jack Daniels' Running Formula.

I'm not talking about the stuff created by yeast poop. Jack Daniels is a real person. He is a coach with a doctorate in exercise physiology and fifteen years ago he published a book called Daniels' Running Formula. He also popularized the use of a physiological measurement called VDOT, which measures the maximum amount of oxygen used by the body. I have had my VDOT measured a couple of times by exercise physiologists, and both times the process was extremely uncomfortable. They put you on a treadmill and stick a tube in your mouth and then pinch your nostrils closed with a clothespin. You then start running, and the treadmill speed and incline increases every minute or so until you can't run anymore.

The last time I did this, the researchers told me that my VDOT was seventy-point-three. Jack Daniels has a specific workout plan for nearly any VDOT number. And so seventy-point-three was the number my dad and I were using when we filled out the calendar.

________________________________________

Here is what Jack Daniels says I should do in a typical week:

Sunday: 13-15 miles at marathon pace (5:40/mile); 2 miles warmup, 2 miles cooldown.
Monday: 10 miles
Tuesday: 12 miles in the morning, 4 miles in the evening
Wednesday: 8 mile recovery run
Thursday: 2 mile warmup. 2x12 minutes of anaerobic threshold pace (5:20/mile). 4 mile recovery jog followed by another 12 minutes anaerobic threshold pace (5:20/mile).
Friday: 12 miles
Saturday: 12 miles in the morning, 6 miles in the evening


When I read this I thought that Jack Daniels was nuttier than squirrel shit. Obviously there are some people who disagree; there are some out there who think this is easy and no big deal. Fly them.

I've stuck with it nonetheless; and the odd thing is, it seems to be working. My runs are getting easier; workouts are getting better; I feel more prepared and more motivated. The only problem is that during the day I'm too tired to do anything else. I get up in the morning around six and am usually done around seven-thirty. For the next hour or so I stretch a little bit and drink chocolate soy milk straight out of the carton. I put a towel on the living room floor and lie there sweating, and sometimes put on music I checked out from the library. For the past month or so I've been listening primarily to jazz - Duke Ellington, Earl Klugh, Ella Fitzgerald, Bill Evans, Thelonious Monk, to name a few. I listen to some classical stuff as well. This morning as I laid on my sweat-soaked towel on the floor, I listened to a sonata by Ravel. It was heavenly.

Then I go to work and feel the slow burn in my legs for the rest of the day. In a way it feels almost pleasurable. And then when I get home sometimes I run again. Sometimes I will see other people out there running, and wonder about whether they are faster than I am.

Some would call this insecurity. They're probably right. In any case, I like to call it competitiveness.

________________________________________

Speaking of competitiveness: A couple of years ago a girl told me that one of her previous boyfriends had run a marathon in two hours and thirty-two minutes. At the time, that was faster than I had run a marathon, and so naturally I became insanely jealous. It didn't help that I was hopelessly in love with this girl. The fact that she had dated someone with a faster marathon time was an insult to me.

For the next nine months I trained like a banshee and that summer I ran a marathon in two hours and thirty minutes. Somewhere in the back of my mind I was convinced that when I told her my marathon time, she would rip off all her clothes with uncontrollable passion.

That never happened. The next time she saw me she didn't even take off her sunglasses. Things have been pretty much the same.

________________________________________

By the way: When someone wants to run with you, wait for a really cold morning and tell them that you want to run then. At the last minute they will come up with an excuse not to go.

You know why? Because running in the cold sucks.

Pro Tip: Changing SPM Output Files

A couple of weeks ago, alert reader Paula DiNoto asked about how to output SPM files in NIFTI format instead of ANALYZE (.hdr/.img). Apparently, for those of you bleeding-edge adrenaline junkies who just have to have the latest thing, SPM12 doesn't give you this option in the interface. SPM5, however - you know, the one that the rest of us barbarians use - still does. (I'll save my gloating for when I meet you face-to-face at the conferences.)

Getting back to the story, I had to admit to Paula that I was, for the first time in my life, unable to figure out how to troubleshoot the problem. However, this was merely a ruse to get Paula motivated; she ended up solving the problem by changing "spm_defaults.m" and changing the default.images.format field to 'nii'. Obviously I already knew this, but I was pleased to see her figure it out for herself.


Thanks again to Paula, whom I once called, in a moment of self-forgetfulness, the most wonderful person I had ever met.

Back Next Week

Fellow brainbloggers,

I'll be going home for a few days - running a few races, seeing a few friends, eating a few pizzas, kicking back a few Shirley Temples. Getting that messed up on Sprite and grenadine means that updates may be scarce. When I return, though, I will be going to FMRI pound-town on those resting-state analyses and hit you with the full report as promised.


Eric my man, take us outta here.

Psychophysiological Interactions Explained

That last tutorial on psychophysiological interactions - all thirteen minutes of it - may have been a little too much to digest in one sitting. However, you also don't understand what it's like to be in the zone, with the AFNI commands darting from my fingertips and the FMRI information tumbling uncontrollably from my lips; sometimes I get so excited that I accidentally spit a little when I talk, mispronounce words like "particularly," and am unable to string two coherent thoughts together. It's what my dating coach calls the vibe.

With that in mind, I'd like to provide some more supporting commentary. Think of it as supplementary material - that everybody reads!

Psychophysiological interactions (PPIs) are correlations that change depending on a given condition, just as in statistics. (The concept is identical, actually.) For example: Will a drug affect your nervous system differently from a placebo depending on whether you have been drinking or not? Does the effectiveness of Old Spice Nocturnal Creatures Scent compared to Dr. Will Brown's Nutty Time Special Attar (tm) depend on whether your haberdasher is Saks Fifth Avenue, or flophouse? Does the transcendental experience of devouring jar after jar of Nutella depend on whether you are alone, or noshing on it with friends? (Answer: As long as there's enough to go around, let the good times roll, baby!)

So much for the interaction part, but what about this psychophysiological thing? The "psychophysiological" term is composed of two parts, namely - and I'm not trying to insult your intelligence here - a psychological component and a physiological component. The "psychological" part of a PPI refers to which condition your participant was exposed to; perhaps the condition where she was instructed to imagine a warm, fun, exciting environment, such as reading Andy's Brain Blog with all of her sorority sisters during a pajama party on a Saturday night and trying to figure out what, exactly, makes that guy tick. Every time she imagines this, we code it with a 1; every time she imagines a contrasting condition, such as reading some other blog, we code that with a -1. And everything else gets coded as a 0.

Here are pictures to make this more concrete:

Figure 1: Sorority girls reading Andy's Brain Blog (aka Condition 1)

Figure 2: AvsBcoding with 1's, -1's, and 0's

The physiological component, on the other hand, simply refers to the underlying neural activity. Remember, however, that what we see in a typical FMRI timeseries isn't the actual neural activity; it's neural activity that has been smeared (i.e., convolved) with the hemodynamic response function, or HRF. To convert this back to the neural timeseries, we de-smear (i.e., deconvolve) the timeseries by removing the HRF - similar to a Fourier analysis if you have done one of those. (You have? NERD!)

To summarize: we assume that the observed timeseries is a convolution of the underlying neural activity with the HRF. Using NA for neural activity, HRF for the gamma response function, and (x) for the Kroenecker product, symbolically this looks like:

TimeSeries = NA (x) HRF

To isolate the NA part of the equation, we extract the TimeSeries from a contrast or some anatomically defined region, either using the AFNI Clusterize GUI and selecting the concatenated functional runs as an Auxiliary TimeSeries, or through a command such as 3dmaskave or 3dmaskdump. (In either case, make sure that the text file is one column, not one row; if it is one row, use 1dtranspose to fix it.) Let's call this ideal time series Seed_ts.1D.

Next, we generate an HRF using the waver command and shuffle it into a text file:

waver -dt [TR] -GAM -inline 1@1 > GammaHR.1D

Check this step using 1dplot:

1dplot GammaHR.1D

And you should see something like this:

Now that we have both our TimeSeries and our HRF, we can calculate NA, the underlying neural activity, using 3dTfitter (the -Faltung option, by the way, is German for "deconvolution"):

3dTfitter -RHS Seed_ts.1D -FALTUNG GammaHR.1D Seed_Neur 012 0

This will tease apart GammaHR.1, from the TimeSeries Seed_ts.1D, and store it in a file called Seed_Neur.1D. This should look similar to the ideal time course we generated in the first step, except now it is in "Neural time":

1dplot Seed_Neur.1D

Ausgezeichnet! All we have to do now is multiply this by the AvsBcoding we made before, and we will have our psychophysiological interaction!

1deval -a Seed_Neur.1D -b AvsBcoding.1D -expr 'a*b' -prefix Inter_neu.1D

This will generate a file, Inter_neu.1D, that is our interaction at the neural level:

1dplot Inter_neu.1D

And, finally, we take it back into "BOLD time" by convolving this timeseries with the canonical HRF:

waver -dt 2 -GAM -input Inter_neu.1D -numout [numTRs] > Inter_ts.1D

1dplot Inter_ts.1D

Let's pause and review what we've done so far:

1. First, we extracted a timeseries from a voxel or group of voxels somewhere in the brain, either defined by a contrast or an anatomical region (there are several other ways to do this, of course; these are just the two most popular);
2. We then created a list of 1's, -1's, and 0's for our condition, contrasting condition, and all the other stuff, with one number for each TR;
3. Next, we created a Gamma basis function using the waver command, simulating the canonical HRF;
4. We used 3dTfitter to decouple the HRF from the timeseries, moving from "BOLD time" to "Neural time";
5. We created our psychophysiological interaction by multiplying our AvsBcoding by the Neural TimeSeries; and
6. Convolved this interaction with the HRF in order to move back to "BOLD time".

Our last step is to enter this interaction and the region's timeseries into our model. Since we are entering this into the model directly as is, and since we are not convolving it with anything, we will use the -stim_file option. Also remember to not use the -stim_base option, which is usually paired with -stim_file; we are still estimating parameters for this regressor, not just including it in the baseline model.

Below, I have slightly modified the 3dDeconvolve script used in the AFNI_data6 directory. The only two lines that have been changed are the stim_file options for regressors 9 and 10:

 #!/bin/tcsh

set subj = "FT"

# run the regression analysis
3dDeconvolve -input pb04.$subj.r*.scale+orig.HEAD                       \
    -polort 3                                                           \
    -num_stimts 10                                                       \
    -stim_times 1 stimuli/AV1_vis.txt 'BLOCK(20,1)'                     \
    -stim_label 1 Vrel                                                  \
    -stim_times 2 stimuli/AV2_aud.txt 'BLOCK(20,1)'                     \
    -stim_label 2 Arel                                                  \
    -stim_file 3 motion_demean.1D'[0]' -stim_base 3 -stim_label 3 roll  \
    -stim_file 4 motion_demean.1D'[1]' -stim_base 4 -stim_label 4 pitch \
    -stim_file 5 motion_demean.1D'[2]' -stim_base 5 -stim_label 5 yaw   \
    -stim_file 6 motion_demean.1D'[3]' -stim_base 6 -stim_label 6 dS    \
    -stim_file 7 motion_demean.1D'[4]' -stim_base 7 -stim_label 7 dL    \
    -stim_file 8 motion_demean.1D'[5]' -stim_base 8 -stim_label 8 dP    \
    -stim_file 9 Seed_ts.1D -stim_label 9 Seed_ts    \
    -stim_file 10 Inter_ts.1D -stim_label 10 Inter_ts    \
    -gltsym 'SYM: Vrel -Arel' -x1D X.xmat.1D -tout -xjpeg X.jpg                             \
    -bucket stats_Inter.$subj -rout


This will generate a beta weights for the Interaction term, which can then be taken to a group-level analysis using your command of choice. You can also do this with the R-values, but be aware that 3dDeconvolve will generate R-squared values; you will need to take the square root of these and then apply a Fisher's R-to-Z transform before you do a group analysis. Because of the extra steps involved, and because interpreting beta coefficients is much more straightforward, I recommend sticking with the beta weights.

Unless you're a real man like my dating coach, who recommends alpha weights.

Context-Dependent Correlations (aka PsychoPhysiological Interactions)

As promised, here is the tutorial on context-dependent correlations. Just to confuse you a little, I also refer to them as psychophysiological interactions, or PPIs, because different people call it different things. Kind of like how some people call me Andy, but other people call me Andy Pandy. Specifically, my mom.

More details about PPIs - and my childhood - coming up soon, but first I have a level 2 scanning emergency to take care of tonight. And I mean that in all seriousness, not the way I usually toss that phrase around - to get out of awkward dates. And now you know, ladies!

Future Functional Connectivity Tutorials, and Other Updates

A few notes:

1) The previous functional connectivity posts and tutorials are cribbed from Gang Chen's homepage, which is available here. Kind of like the way I crib quotations and passages from authors that no one reads anymore, and then pass it off as my own style to boost my pathologically low self-esteem. Keep in mind that most of these demonstrations deal with a single subject and simplified situations that you probably will not encounter in your research. Given these contrived examples, most of the results generated in these demos are relatively meaningless; it's up to you to learn and understand the concepts, and then apply them to your own data and make your own inferences. My task which I am trying to achieve is, by the power of Youtube tutorials, to make you hear, to make you feel — it is, before all, to make you understand. That — and no more, and it is everything. (That was Conrad, by the way.)

2) A lot of you - I'm talking a LOT of you players - have been making requests for MELODIC tutorials and resting state analyses in FSL. All I can say is, we'll get there, in time. Before that, however, I believe AFNI is better suited for building up one's intuition, and so we will be working through a few more connectivity topics in AFNI - specifically, context-dependent correlations, beta series correlations, and resting state connectivity. After that we will again cover the same concepts, but applied in FSL - by which time, given my glacial pace, either FMRI will have become a passé technique or the Andromeda galaxy will have crashed into us.

3) Recently you may have noticed the "Donate" button on the right sidebar of the blog. This was done at the request of one reader who felt the powerful, irrational urge to loosen his purse-strings and give some alms out of the goodness of his heart, which is located somewhere way, way down there, somewhere nearabouts the cockles. Although I can't fully understand this behavior - even less than I can understand why there is someone who still has purse-strings, or what cockles are, exactly - nevertheless it helps satisfy my cupidity and strokes my ego. Furthermore, in addition to serving Mammon, these tokens of gratitude motivate me to regularly produce new material and, as a bonus, help me to continue procrastinating on my dissertation. Now that's what I call a win-win-win.

4) Also, at least one of you has mailed me a two-pack of Nutella. This has pleased me greatly. My brain needs hazelnut spread for fuel, and the more it has, the hotter and better it burns.

5) If everything goes according to plan, we should cover context-dependent correlations this weekend, beta series correlations next week, and resting-state connectivity the week after that.

Lunch in Paris, dinner in London, comrade.

Creating Functional Connectivity Maps in AFNI (Second Step on the Road Towards Glory)

Having used 3dSynthesize to create a dataset with the effects of no interest, we are now ready to subtract those effects from the dataset that went into 3dDeconvolve; for you players using the AFNI_data6 tutorial data, this is the all_runs.FT+orig dataset. We can use 3dcalc for this operation:

3dcalc -a all_runs.FT+orig -b effectsNoInterest+orig -expr 'a-b' -prefix cleanData+orig

This will create a time series with all the drift correction and motion parameter schmutz subtracted out and cast into the pit. What is left over will be signal untainted by motion or drift effects, which we can then use to generate a seed time series.

Before that, however, there is the additional step of warping our data to a template space (such as MNI or Talairach). Assuming you have already warped your anatomical dataset to a template space either manually or using @auto_tlrc, and assuming that there is already good alignment between your anatomical and functional data, you can use the command adwarp:

adwarp -apar anat_final.FT+tlrc -dpar cleanData+orig -dxyz 3 -prefix EPI_subj_01

(Note: There are better ways to do warping, such as the warping option built into align_epi_anat.py (which calls upon @auto_tlrc), but for pedagogical purposes we will stick with adwarp.)

Once this is done, all we need to do is select a seed region for our functional connectivity analysis. This could be a single voxel, a region of interest that averages over an entire area of voxels, or a region selected based on a contrast or an anatomical landmark; it depends on your research question, or whether you want to mess around and do some exploratory analyses. For now, we will use a single voxel using XYZ coordinates centered on the left motor cortex (note that, in this coordinate system, left is positive and posterior is positive):

3dmaskdump -noijk -dbox 20 19 53 EPI_subj_01+tlrc > ideal_ts.1D

Because 3dmaskdump outputs a single row, we will need to transpose this into a column:

1dtranspose ideal_ts.1D LeftMotorIdeal.1D

After you have output your ideal timeseries, open up the AFNI viewer and make sure that what was output matches up with the actual timeseries. In this example, set EPI_subj_01+tlrc as an underlay and right click in the slices viewer to jump to XYZ coordinates of 20, 19, 53; open up the timeseries graph as well, and scroll to the first time point. Does the value there match up with the first entry of your ideal time file? If not, you may have mis-specified the voxel you wanted, or you are using a different XYZ grid than is displayed in the upper left corner of the AFNI viewer. If you run into this situation, consult your local AFNI guru; if you don't have one of those, well...

In any case, we now have all the ingredients for the functional connectivity analysis, and we're ready to pull the trigger. Assemble your materials using 3dfim+ (which, as far as I can understand, stands for "functional intensity map," or something like that):

3dfim+ -polort 3 -input EPI_subj_01+tlrc -ideal_file LeftMotorIdeal.1D -out Correlation -bucket Corr_subj01

You will now have a dataset, Corr_subj_01+tlrc, which is a functional connectivity map between your seed region and the rest of the brain. Note that the threshold slider will now be in R values, instead of beta coefficients or t-statistics; therefore, make sure the order of magnitude button (next to the '**' below the slider) is 0, to restrict the range between +1 and -1.

The last step is to convert these correlation values into z-scores, to reduce skewness and make the distribution more normal. This is done through Fisher's R-to-Z transformation, given by the formula z =  (1/2) * ln((1+r)/(1-r)). Using 3dcalc, we can apply this transformation:

3dcalc -a Corr_subj01+tlrc -expr 'log((1+a)/(1-a))/2' -prefix Corr_subj01_Z+tlrc

Those Z-maps can then be entered into a second-level design just like any other, using 3dttest+ or some other variation.

Sanity check: When viewing the correlation maps, jump to your seed region and observe the r-value. It should be 1, which makes sense given that this was the region we were correlating with every other part of the brain, and correlating it with itself should be perfect.

Like me.

Debussy: La Fille Aux Cheveux De Lin

FMRI analysis can be stressful; graduate school can be stressful; life and relationships can be stressful. Therefore, it can be salutary to take stock of everything once in a while, just to keep things in perspective; because all of those things listed above can, sometimes, if you pay attention - be very, very good.

Just kidding; you're probably going to have an extremely stressful week anyway. Regardless, here's some relaxing music to help you out.

AFNI Command: 3dSynthesize (First Step on the Road Towards Glory)

Before we get into full-blown functional connectivity mode, we need to get acquainted with a command called 3dSynthesize. This will create a fitted time series based on selected regressors from your model, and will be used for subtracting out effects of no interest - that is, things like drift correction and motion regressors - so that the correlation we observe is protected from systematic differences in things like head motion, which can be particularly problematic in patient populations and children.

3dSynthesize requires a bucket dataset only containing coefficients from your model; therefore, you may need to go back and rerun 3dDeconvolve with the -cbucket option to make that magic happen. Once done, feed 3dSynthesize the coefficient bucket, the X matrix from 3dDeconvolve, and a selection of columns that you want for a fitted time series, e.g.:

3dSynthesize -prefix effectsNoInterest -cbucket CStats.FT+orig -matrix X.xmat.1D -select 0 1 2 3 4 5 6 7 8 9 10 11 14 15 16 17 18 19

The columns I've selected here are based on the dataset in AFNI_data6; obviously, yours will probably differ. The columns can be selected based on number, as above, or the actual label names, which can be retrieved from the X matrix with the grep command:

grep ColumnLabels X.xmat.1D

In any case, check the X matrix with either aiv or 1dplot to make sure that the columns you are extracting are indeed the effects that you want.

Hit the video for a brief demonstration, as well as a 30-second monologue to get you pumped up to use 3dSynthesize.