31 - Resting State & InstaCorr: Part 2 of 2
February 15, 2019
June 1, 2018
Rick Reynolds, NIMH
All Captioned Videos AFNI Training Bootcamp
Rick Reynolds, NIMH
For more information and course materials, please visit the workshop website:
We recommend viewing the videos at 1920 x 1080 (Full HD) resolution for the best experience. A lot of code and type is displayed that may be hard to read otherwise.
RICK REYNOLDS: Let's first of all look at how to do a resting state analysis in AFNI using afni_proc, and some choices you can make in how to apply your choices on the command line. And then we'll fire that off. And while that's pondering life, the universe, and everything, we can look at some data and play around a bit.
So let's CD back into our AFNI_data6/FT analysis directory. We've already done the task analysis. We've done surface analysis.
We did the PPI analysis. We even did an analysis with task that was in grid space. We've got a lot of garbage here already.
But now we'll look at some-- we've got a couple of resting state examples in here. Script 6 and script 7 are examples of resting state analysis. So I'll just look at those briefly. I'll try not to babble too much, because time is disappearing quickly. So let me cat S06 here, just dump it to the terminal window.
So this was example 10. I'll have to actually go back and verify if it still looks like example 10. We have basically most of the same options in the script that we had for a task. Except now, again, maybe we'll start with the d spike operation.
But the rest of these blocks are perhaps the same blocks that we used for the task-based analyses. And then we still give it the same anatomical data set. We still request it to remove the pre steady state time points from the EPI data.
Note that I'm using the same EPI data from task. I don't care that it's task. This is just going to show you the analysis. So we'll pretend this is resting state data.
But of course, we have big, bold signals in it. So this is just to look at how the analysis is run. We're still registering to the main outlier time point and go to standard space. But now we have a couple of other extra options here.
Now we're going to segment the anatomical data set and erode it with the purpose of doing this white matter eroded regression. So is this good? I don't know. But we're going to create a white matter mask and average signal over it, and then project out an average signal as a tissue-based regressor. This is how you would do it.
This WME, if you look in the afni_proc help for what tissue classes are created automatically that comes out of this segmentation step so you can look at the help. Again, you can give afni_proc any masks that you have already generated that are in alignment with your T1. So if you've made a white matter mask, if you've made a ventricle mask, whatever, you can give them to afni_proc
And it will bring the mask along with any anatomical transformations that you apply, including a nonlinear work to standard space. It will bring your mask along. And then you can still use this mask for tissue-based regression or principal component regression of CSF, whatever.
So afni_proc's segmentation, the segmentation in AFNI IS not fantastic. We can get white matter, CSF, and gray matter. So even the CSF is all CSF. It's not just ventricle.
But you can do various things. You can import. You can intersect your CSF mask with a ventricle template mask, and then the intersection will be a reasonable ventricle mask for that subject.
You can apply it that way. So there are a lot of options. Again, we try to overwhelm you with options. But if you have a masking technique that you'd like, you can just hand afni_proc the masks.
So we censor motion. But note that we're being a little more strict here. A motion makes a bigger difference in resting state analysis. So you'll typically censor at a more strict level than you would with task. And we're also going to censor outliers.
In this example, it's very common, we're going to do the bandpassing, even though I may whine about it a bit. So we'll bandpass between the 0.01 and the 0.1 hertz. So in this case, we will throw away 60% of our degrees of freedom.
Now I'll try to actually show you that. We're also regressing out the motion parameters and their first differences. Let's fire this off. And then we'll come back to it, and then we can be more clear. We can make it more clear exactly what is happening.
So I'm going to run this thing-- TCSHS06@ap.rest. Again, it gives us this long format command. This is the nice way to run the command. But again, most of you are using Bash. So you could just type this much in, or even drop the XEF stuff and just tcsh_proc.ft.
But make sure you get the rest one. We have a few proc files here. And if you get one that has already been run, it's going to whine at you.
So just to see TCSH and the proc.ft.rest. I'm going to run the whole thing, so I capture the text I'll put to this file. So I can look at it later if I choose to. So let's just fire that off, and we'll let it go. That's piping.
When programs dump stuff in your terminal window, it's often written to two file streams, even though they're both displayed-- standard out and standard error. Standard error is a separate file stream for special messages. So you can capture the errors differentially.
But it all appears here. So if you want to capture everything, you have to capture both streams. And that pipes the standard error and the standard out to the same stream, so that you can capture it all.
You can do that with Bash, too. And you can see the afni_proc help for the syntax. Or you can give a dash Bash option to afni_proc, and it will tell you how to do it. Let's let that do its thing for a while. And now let's go do some Insta stuff.
So for this, let's leap into AFNI_data5. Hopefully all of you have this. I don't know if you do. But let's go into the AFNI_data5 directory.
How many people don't have AFNI_data5? So in this directory, there are just a handful of data sets. You've got these resting state data sets. You've got, basically, an original resting state data set, one that's gone through alignment, and one that's been fully processed by, say, afni_proc. We'll focus on the first one, but let's just run AFNI here now.
So let's see. Let's set the underlay to be the T1 data set that was aligned to the EPI. And so we pretend that we know where we are in space. And now, to do this, we could display the EPI time series data. But let's not bother with that. We've looked at what that looks like.
Let's click on define overlay. But we're not going to do any colorization right now just yet. There's this Insta stuff menu up here-- another dropdown selector panel with InstaCorr, InstaCalc, 3dTstat, and Group InstaCorr.
I'll just give you a quick rundown of what these are. We actually have-- just so you see it. I will not get to this. But there is in AFNI handouts, AFNI-- what's the Insta stuff.
Handout 20 is Insta stuff. It goes through some of these things. It doesn't have the new Tstat Insta stuff. But it has-- you can see clustering. It has InstaCorr. The fun part, blurring [INAUDIBLE] and then Group InstaCorr.
So it has descriptions of these things. So if you want to do them on your own later, you can read through the slides. Let's just do it do this stuff live, though. But just so you know, we have the documentation if you want to go back and look at it.
So this menu has the InstaCorr, InstaCalc. That's like 3D calc, but it's in the GUI. So if you want to add two data sets or multiply some data set by 1.7. Whatever you want to do, you can do that here.
But let's just leave it on InstaCorr, the first one. And then let's click on the Setup InstaCorr button-- Set up ICorr. And this interface is basically allowing you to do more or less all the pre-processing, except for motion registration, volume registration. It doesn't do that, but it will do a lot of other pre-processing.
And the data set that we'll choose, if we click on Choose data set Now on the first row, we'll choose this S620rest run1 plus the ridge data set that has not been processed, except for presumably registration. And you don't have to do this. Let me just briefly show you what the time series looks like.
That looks pretty nice. That EPI data has nice contrast. This is nice-looking EPI data. If you look at the graph window, and we get out of the mid-sagittal plane, you notice-- I don't know if you can see. But the red dot at the beginning of the time series all the way at the left is way high up.
So we have pretty steady state data here. If you look in the orthogonal views, the coronal or sagittal views, you can see the striping. This striping, I guess, I'll mention briefly, is similar to a motion stripe.
Why do we get it? Well, remember, this magnetization for the gradient is being passed over the data. How does that happen? When it first starts off, you get a large MRI signal.
I'm not going to get into the physics much. But at the very beginning, you get a big signal due to this. But it's going to pass over slice 0, 2, 4, 6, 8, etc., and then 1, 3, 5, 7, 9.
But it bleeds. It's not perfectly within each slice. This bleeds across neighboring slices to some degree. That's why we used alternate interleaved slice acquisition because of that bleeding.
So the even slices get hit by this first. And then you go back and get the odd slices, the odd slices have already been magnetized a bit. They've already gotten some of this. This field map has already hit them to some degree. So you actually see a difference between the even and the odd slices here.
The pre-magnetization steady state effect is more clear in these views because you see between even and odd slices. For motion, it's different. For motion, the subject is moving differentially over space.
So you'll see a Venetian blind effect there, too. But that's slightly different, unless there's an inflow artifact that you want to talk about. But anyway--
RICK REYNOLDS: Yeah, many scanners have dummy scans that you don't even get your hands on. In general, we might suggest-- you don't need to let the scanner throw away anything. For some reason, if your data doesn't-- because you can throw it away whenever you choose to.
But for example, if for some reason, your EPI data doesn't have fantastic contrast, this EPI data is beautiful. But sometimes, you don't have such nice contrast. What does the contrast help you do? It helps you to register.
If you look at the EPI data, and you can hardly tell the difference between white and gray and CSF-- if you can hardly tell there, even with large voxels, if you can't tell the difference, you know the computer's not going to tell the difference that well. And then aligning the EPI to the [INAUDIBLE] well, it's just going to stick one blob over the other blob and say we're done. If you want the nice cortical contours being followed in the registration, you need some contrast in the data. So you care about that.
So backing up, why would you collect the pre-steady state volumes? Well, if you don't have such great contrast for [INAUDIBLE] to EPI registration and the EPI data, you can use the pre-steady state volume. The very first volume often has much better contrast than later ones. So you can use this for the anatomical alignment step.
If I go forward in time a bit here, you can see the contrast is very different. Time 0.0, it's fantastic. If I go to a few time points later, it's still very good, but it's not as good as it was before. This still will register very well.
But in some cases, you can look at your own data. In some cases, it's not as good. Let's just rip through some of the processing steps. Let's remove, I don't know, three time points at the beginning. I don't even know what the TR is, but we don't have to worry about that
Let's apply perhaps a 5 millimeter blur to the data. We can automask it. That's great. Fine, we'll do bandpassing since it's set on by default, even though I'll cry about it a little bit.
Despiking takes time. Let's not despike. Seed radius-- so if you click on the miscellaneous options, you can control the seed radius. So what seed radius means-- if you're going to choose a seed for correlations, you can grab just one voxel and correlate if you want to. But also, when you choose a voxel, you can dilate it a bit by averaging with the neighbors, and get a little ROI average time series and use that for the correlations.
So what we want to do is, at some seed voxel, we're going to get its time series and correlate the rest of the volume with that. That's the correlation analysis you do with resting state, or the seed-based correlation analysis. So that's what we'll do here.
And let's drop a seed radius of 4 millimeters, say. You can use 4.2 if you want to. These don't have to be integral.
Polar 2 for detrending, fine. That's good enough. You can do many things here, too. You can actually do temporal windowing if you want to in this interface. But I won't do that.
Any questions about the options? So we're essentially doing pre-processing steps that you might otherwise do with afni_proc. If you're really going to play with this, you might as well just pre-process it with afni_proc, and then come here and don't do anything. Or possibly use a seed radius.
That might apply in any case. But this data hasn't been processed. So let's hit Set up and Keep. And most of this is fairly fast. You can see in your terminal window what it's doing.
The blurring operation is perhaps the slowest one. But even that doesn't take too long. So let me just lower you. And in my main controller window now, you see it says ready underneath there.
So fantastic. I will raise-- I'm going to leave the EPI images up there in case we want to look at them for some reason. But let's do the correlations while we're looking at the anatomical data set. So you, maybe start down in the visual area, and maybe hunt around for the default mode network. We can play a bit.
So how do you do an instant correlation? Well, you can click somewhere in the brain. Let's just click down here. And then you can right-click and you can say InstaCorr set and with thresholding. Let's go back to the main interface and just drop the threshold.
It's at 0.5. We don't need to look at high correlations. Let's look at the whole volume. That's prettier. You can ponder threshold on your own.
So boom, we get a big correlation in the visual cortex, say. Everyone see that? Pretty fast.
So to make it slightly faster, well, that's instant. But jumping to a coordinate, right-clicking InstaCorr Set, I'm already tired. So we can be more instant than that. If you hold down Shift and Control-- so hold down on both Shift and Control at once, and just keep them down.
And now left-click somewhere. Left click, left click, left click. Now we're getting more Insta right here.
Now when we click, this is happening right away. So the data has been pre-processed. So producing a correlation at this point, it's probably just a dot product, assuming the data has been also scaled to be unit height.
So sum of squares equal to 1, you can just do a dot product all over the place. Very fast. So again, what we're doing is we click somewhere. We get the seed time series. And then we're correlating with every other point in the volume and producing this map.
RICK REYNOLDS: This is just one subject. This data has not been pre-processed, except for what we're doing. It is probably registered, but otherwise unprocessed, except for what we do up here. So we did bandpassing.
We did some pre-processing steps in the setup. And then it finished that. And now it had data ready for these correlations.
Still, that's not fast enough. [INAUDIBLE] side whine to Bob about that, too. I want it to be more Insta.
So now hold your Shift and Control keys down again, and you can click. But not just click. You can click and drag.
And oh, you could look at this for hours. Drag in here. We perhaps-- I don't know if we're getting into the default mode network here.
I didn't. I'm not going to worry too much about what networks you need to find here. But just browsing your data and seeing how networks change over space is very, very interesting.
And I'm going all the way up here, up in the right. And what is this? I drift up here. That does not look good.
These very high correlations and a big blob that's not following anatomical structure, it's going through all the slices. You can see in the orthogonal views. It even goes way down here in this coronal view. That's horrible.
What is going on? What does the EPI data look like there? So we're in this location. What does the EPI times series look like?
I don't know. It looks fine. I don't see any particular strangeness there.
Why would we get such a correlation? You might actually be able to see-- you don't have to do this again. But let's see if we find anything if we look at our voxels.
That's still not good enough. Let's just make them taller. And I guess in some cases, you can actually see a temporal pattern repeating over space by looking at a lot of voxels at once.
I'm not quite seeing it here. But what's actually going on here is, this is a scanner artifact. This is one of the coils going bad a little bit.
And if the signal drops down by a few percent for a few TRs. Or maybe for a lot of time, and goes up and down, who knows? When the scanner coils start going bad, at first, it will change a little bit. And then as they get worse and worse and worse, you get more and more problems over time.
And this is what you see. You notice the size of this. This is 20, 30 milliradius spherical-ish thing, depending on where that receiver coil is on the proximity to the subject. This is what you see.
So this is a hardware artifact. And this is one of the things that running in [INAUDIBLE]. So let's just jump briefly back to the InstaCorr plug-in and choose as data set the IRTS data set. This is the exact same data set, but run through afni_proc, including the inadequate step.
So that means that each voxel, you take a local white matter average and regress that out differentially for each location and space. You're doing a voxel-wise regression there-- one degree of freedom for this step. But it varies across space. But otherwise, it's the same data pre-process state in approximately the same way that you might do in the GUI.
So I'll hit Set here. Now we need to change a few things. Let's start at 0. Though if you start at 3, who cares?
Blur-- let's drop the blur to 0, because I think the pre-processing [AUDIO OUT] blurred. Automask, I guess yes. We don't have a mask.
Bandpassing, let's turn that off. It's already been processed, however it's processed. We can leave the seed radius in there. [INAUDIBLE], who cares? It's already been detrended, but that won't have any effect.
So now we've basically taken off all the pre-processing, except we're still going to use the 4.2 millimeter seed radius. So let's hit Set up and Keep. And the [AUDIO OUT] now it's ready. Now you've got the green ready up here.
So now if we hold the Control-Shift-- well, we could just right-click. But now look at it. Now it looks more normal.
And you see, it's [AUDIO OUT] symmetric now. This is the exact same area that had the horrible problem. But it looks more like what it should look.
RICK REYNOLDS: Reflect?
RICK REYNOLDS: Project? Yeah, you don't-- well you could this talking with SUMA and have this projected in SUMA. I think Daniel will actually do that with you today. But you can also do InstaCorr in SUMA.
We ran a surface-based analysis in here. You could take the IRTS time series, and do this InstaCorr in SUMA. It'll look fantastic.
I'll mention, you can also do this at the group level. So this is dragging all around and getting volumetric correlation images move over time. It's fantastic with single subject.
At the group level, it isn't quite as Insta. But you could drop-- that's what Group InstaCorr is. Group InstaCorr, we'll do the little bit of pre-processing ahead of time for your 50 subjects. But then, you have all the same controls to do across 50 subjects. You could even do two groups and get a two-sample t-test. What is it doing?
Suppose the time series has been set up so that you can just do a dot product for the correlation. So what happens now? At each voxel, you correlate with the seed time series and you get a map of correlations.
Those are converted to z-scores. Those are normalized z-scores. Now you've got the normalized z-scores. You can do your t-test on those normalized z-scores. That gives you a t stat across space for every one of your 50 subjects, if your computer has enough RAM to handle this.
And it doesn't take that much. It's not so incredible. I could do this [INAUDIBLE]. I don't know, 50.
But I could do a lot of subjects on a laptop. And then after it does the t-test, it converts that t stat to a normalized z-score as well. And that's what it projects in the image.
So it's not as fast as this. It probably takes a second or two. Obviously, it's going to depend on how many subjects you have. So you won't be able to necessarily drag around and do this.
I mean, look at how fast it does the correlations, how fast that image changes. You can't do that with 50 subjects. But still, you can go pretty quickly. You can just [AUDIO OUT]. That's probably how you would use it.
It might take a second or two. But still, that would be a group analysis right there. But that's just for looking at your data and getting a feel for things.
And when you actually want to run the analysis, of course, you would script it so that you have the record of exactly what you did when you need to talk to people about that. You can even script Group InstaCorr to do these things for you, and leap back to the terminal window. In the single subject analysis directory, that analysis finished long ago.
Let me just briefly look at the proc script and see what that looked like. I guess if I'm going to just look at the proc script, I didn't have to [AUDIO OUT]. But we wanted to know how do we handle things? Like, if we're not going to use a fast Fourier transform for bandpassing, how do we handle this in a regression model that's censoring and all that, doing that all that at once?
So look at that tier, [AUDIO OUT] going to happen after the vol rate block. So the vol rate block doesn't look special. The blur block doesn't look special. We know we want tissue-based regressors to be created from the vol rate.
But we haven't done that yet. Afni_proc knows what the vol rate output is called. So it hasn't worried about it.
In the mask block, we should have a 3D seg command in here. So here's where AFNI is computing the segmentation of brain matter in CSF. But I won't focus on that.
Here's the regress block. So here's where the action happened. So these are all basically the same things as before.
You demean the motion parameters, get first differences. You want to create a sensor time series for time points where you had more than 2 millimeter motion. You actually want to combine this sensor with the outlier censoring that was computed up above, since in this case, we did request censoring outliers. So we do all that stuff here.
Now it [AUDIO OUT] have 48 magic. So here, we actually create bandpass regressors. And we still give it the band range of 0.01 to 0.1 hertz.
But we use this 1 dB port that's, like, bandpass orthogonal terms. That ort is our shortcut for saying orthogonalization terms that we use with [INAUDIBLE] and stuff like that-- ortvex and whatnot. So this creates 1D files, which is to say text files of those sinusoids. And we write it in this bpass run whatever .1d.
Let me just briefly show that to you. You don't have to try to keep up with all the other stuff. But in the rest results-- I'm sorry, what was it called? bpass. In the rest results directory, now we should have these bpass files for the three runs.
So if we just look at those all that's just a mess. But [AUDIO OUT] a text file of numbers. That's our 1D format. So obviously, we won't-- obviously, we'd rather plot it. We can't even really see this, can we?
Let's see. So you can see it there a bit. But what are we looking at? You see the sinusoids in the left. You see the 0.01 is our high-pass.
Those are the low frequencies. We're going to regress out things slower than them. So we have a few high-pass terms with low frequencies on the bottom. But then we have our high-frequency [AUDIO OUT] in the top. And they're just the sinusoid pair-- sine and cosines at all of the frequencies that we want to check out.
So you can put these in a linear regression and get rid of them that way. And it's identical to doing bandpassing. We have three runs, so we do band [AUDIO OUT] per run.
The bandpass terms do not cross the run breaks. So we keep the band passing separate per run, because bold effects don't cross run break. So you want to do these separately across the runs. So these are the terms.
And in the full x matrix-- I don't even know if I want to show this here. The full x matrix is going to have all three runs of that, plus the terms. [INAUDIBLE] it's not even-- OK, fantastic.
It's so tight with regressors that we can't even look at it. But we've got a lot of regressors in there. So we've got 431 time points after censoring and 304 regressors. So we're removing a lot of degrees of freedom for the bandpassing operation.
It's not magic. It's not just [INAUDIBLE] Fourier. It's elegant.
It's just a quick step. It's a very costly step. You can do it, but it's costly.
And again, if you dump 60% your degrees of freedom for bandpassing, and what's a common censor limit? 10%, 15%, 20%? Sure. In the famous papers that tell you how to analyze your data, they went way above 20%. So let's say 20% for censoring.
And then you've got all these [AUDIO OUT] polar terms. You've got your motion terms, motion derivatives. Maybe people are telling you, oh, use the first in 24.
So now you've got twice as many terms. You've got all these regressors of no interest plus tissue-based regressors. How many degrees of freedom are you going to be left with? 50%, 20%, 10%, 5%.
Well, we're basically looking at-- it's as if we had acquired seven time points of data. And we'll do our correlation analysis on that. So that's what it can get down to.
Even if this doesn't feel, you may just have a few degrees of freedom left. And it's as if you had acquired seven time points. So you have to be careful with that, so that you're not wiping out all your degrees of freedom. if you look at the-- here's the output from running this analysis.
If we look at the end of that, it tells us we had 19 time points censored, a number of time points per run. So out of 450 total time points, censoring didn't take off too much-- 19. [AUDIO OUT] down to 431. We used a total of 304 degrees of freedom. So that extra 127 degrees of freedom was for the bandpassing, or and for motion and stuff.
Oh, I'm sorry, no. The 304 is the sum of all the bandpassing and the motion and all that. And it leaves us with 127 degrees of freedom. So when you run this in afni_proc, you should look at this.
For your data, depending on your run length and stuff, you should have some lower limit in mind for when should you just drop your subjects. If you're throwing away all the degrees of freedom, if you get down to 25, do you continue, 30? Again, a lot of the big important papers on telling you how to analyze your data really were in the negatives here.
I don't mean to be all doom and gloom, of course. You still get pretty pictures. I mean, you look at the stuff that we just ran with the InstaCorr-- fantastic. You get nice group results.
But you just want to be careful, so that you're not throwing away everything or putting your data in a bad place. And you get the QC information here. Again, remember, it's trivial to throw all this information right here into a spreadsheet for your subjects. And so you can sort them and see which subjects-- you can get an idea of where the subjects lie in terms of the degrees of freedom left.
Let me just briefly show one last thing, this script 7. Say, example 11, that's just a slightly different example for a resting state analysis. But I just want to show you some of the different options that are in this one.
[AUDIO OUT] is we specify some anatomical follower ROI data sets. In this case, they're showing, like, the FreeSurfer parcellation right here. So you can bring those along. Sorry, Christina-- ventricles derived from those data sets, the white matter averages.
And you can name them. You can say I want an anatomical resolution parcellation and EPI resolution parcellations that I'm going to name aesig, fsvent, fswe. Why am I naming these things? Well, because for example, I want to regress out the top three principle components from the FreeSurfer ventricle ROI for this data. So I'm doing PC regression to get rid of some of the physiological noise.
And this example is not doing bandpassing, say. We do the fast AnataCorr from the FreeSurfer white matter. You can give it your own white matter segmentation, whatever you want to So there are some options. So afni_proc is pulling these masks along.
Let me make a quick little comment about these. We're not using these parcellations in this analysis. We're using the ventricles and the white matter, but not the full parcellation.
So why would you bother doing this? Well, this puts the parcellation in standard space, even with your nonlinear work step. So we're doing nonlinear registration to the template, and it's pulling all these masks along.
Now you've got your FreeSurfer parcellation in standard space for all your subjects. It's trivial to make probability masks for all of the ROIs. You can see how well the ROIs overlap across your subjects.
So you can combine these, make your own atlas. You can make your own FreeSurfer atlas with a maximum probability map. Very easy to [AUDIO OUT].
That's essentially how the Desai atlases were made. So some options to ponder. Any last questions?