[pymvpa] Accounting for the haemodynamic lag

[Brian Murphy]

Hello David,

you have a few options. You can aim to get a single aggregate volume for
each trial, by a) assuming a standard HRF corresponds well to the
responses seen in your data, and use that to do a weighted average of
your trial volumes (NiPy has a built in HRF); or b) do a simple block
average as you suggest (e.g. taking an offset of 2 TRs, and block
averaging the following couple).

You could also do either of those things with a cross-validated
parameter setting step (to decide on the optimal offset and block
length, or HRF parameters) - keeping in mind that the HRF can vary by
individual participant, task and brain area (the 'standard' HRFs are
based on low level visual cortex responding to flashes of light).

If your main interest is getting good classification results, you could
also throw in all the volumes together (each trial would be comprised by
V voxels x T TRs), and let the machine learning methods decide which
volumes (or weighting thereof) are the ones it wants to listen too. The
advantage there is that there you make no assumptions at all about the
timing and shape of the response, and you don't assume a uniform
response across people/locations/tasks. With that approach you'll need a
classifier that works well with large numbers of co-linear dimensions -
e.g. PLR or Random Forests.

We wrote a paper covering some of that ground which you can look at for
background:
http://www.frontiersin.org/Journal/10.3389/fninf.2012.00024/abstract
... and I can dig out the associated code if that is helpful,

best of luck,

Brian

On Mon, 2014-01-20 at 09:55 +0000, David Watson wrote:
> Dear All,
> 
> I was wondering if anyone could give me some advice on how best to
> account for the haemodynamic lag of the BOLD signal when performing a
> pattern analysis on 4D fMRI data?  This seems like a fairly basic
> issue, but I am struggling to find a clear answer on how best to deal
> with it. I have spent some time reading around (e.g. this page from
> the Princeton toolbox was quite informative:
> http://code.google.com/p/princeton-mvpa-toolbox/wiki/HowtosRegressors#How_can_I_take_the_haemodynamic_lag_into_account) and I get the impression that there are two main ways that people tend to do this:
> 
> 1. Offset the timeseries or the sample labels by a suitable number of
> TRs. For instance, my TR is 3 seconds, and the lag is estimated to be
> approximately 6 seconds for most subjects, so I could either remove
> the first 2 TRs of the timeseries, or increment my sample labels along
> 2 time points. I could easily enough do this myself within python once
> I've loaded in my sample attributes and dataset, although maybe PyMVPA
> already has some built in support for this function that I have
> missed. But I am a little concerned as to how accurate this is likely
> to be, e.g. the lag is unlikely to be precisely 6 seconds in all
> subjects.
> 
> 2. Convolve my model regressors with an HRF. This option seems like it
> might be preferable, and I can easily enough derive a gamma HRF (e.g.
> the nitime package seems to provide one), but I'm not sure how I would
> then apply this to a given model within PyMVPA. Or does PyMVPA already
> provide some functionality to let me do this?
> 
> As it happens I have a block design so perhaps I could get away with
> just offsetting the timeseries, although convolving an HRF might still
> be preferable. But if I ever wanted to do an event-related design
> where measuring timings precisely is more important then I'm not sure
> if simply offsetting the timeseries would still be considered
> acceptable. Also, are there any other commonly used methods of
> accounting for the lag that I have missed?
> 
> 
> Regards,
> 
> David
> 
> 

-- 
Dr. Brian Murphy
Lecturer (Assistant Professor)
Knowledge & Data Engineering (EEECS)
Queen's University Belfast
brian.murphy at qub.ac.uk