Hi Linda,

I read through your question, and as far as I can tell with my somewhat 
limited
background in signal processing, most of your reasoning about sampling rates
and so on is correct, but I think that in some ways, steering the 
discussion in
a different direction somewhat might be more helpful.

In short, one of the things to remember is that much of the models that are
built rely on assumptions that were made from empirical data.  So, for 
example,
no one really knows what the high-frequency components of a "true" hrf 
actually
are.  (by "true" I mean to refer to what the hrf might look like if we could
record it in the brain from a single trial with zero noise and super-high
sampling rate) The only way we know the shape of the hrf is by doing
event-related fMRI studies, and since those studies tend not to go any faster
than having a 1 second TR, we lack the high-frequency information that might
alias into a 1Hz sampling of the hrf.  However, the hrf does look 
pretty smooth
if a large number of averages are taken, so it is extremely reasonable 
to assume
that it is bandwidth-limited.

(there are also good physiological reasons that the "true" hrf is probably
smooth, but for most people designing cognitive er-fmri studies, the "why" of
the smoothness of the hrf isn't important)

Another thing to point out is that the essential problem of analyzing 
ER-fMRI is
how to separate the signal of interest from noise, which is much larger 
than the
relatively weak functional BOLD signal.  If there were high frequency 
components
in the "true" hrf that aliased into the signal, they will be indistinguishable
from noise.

The hrf function provided by SPM is referred to as "canonical"-- that 
is, it may
not match the "true" hrf of your subjects perfectly, but it can 
probably capture
most of the signal, and if the time-derivative and dispersion derivatives are
modelled as well, you catch even more.  Just as important as 
considering how to
appropriately model the hrf is the consideration of how to appropriately model
the noise as well.  This is indeed the essential problem of imaging: 
the GLM is
used to attribute the variance in the measured signal to regressors of 
interest
(like hrf-convolved onset trains); regressors not of interest (like regressors
used for head motion); and random noise.

Below the signature, I list a few papers that I have found really 
helpful, and I
also comment on your original letter.

Ken

AM Dale, Human Brain Mapping 8:109-114, 1999.
This paper develops the GLM in a concise manner, and, the thing about it that
was helpful to me, explicitly connects it to the earlier framework of
time-locked selective averaging.

MA Burock, ... AM Dale. NeuroReport 9, 3735-3739, 1998.
Great article on why jitter is important and what it does, with very intuitive
figures.

Design Efficiency in fMRI by Rik Henson
http://imaging.mrc-cbu.cam.ac.uk/imaging/DesignEfficiency
Great article on jittering, uses somewhat of a signal-processing perspective.


R.N.A. Henson, M.D. Rugg, and K.J. Friston. The choice of basis functions in
event-related fMRI. NeuroImage, 13(6):127, June 2001. Note: Supplement 1.
Keyword(s): event-related.
Good paper on how well different choices of basis expansion model the hrf.


Quoting Linda Seltzer <[log in to unmask]>:

> Signal processing discussion:
> This is what I have been looking into so far.  Please let me know whether
> this is correct.  I have seen the functions which Dr. Fromm has mentioned.
>     Consider a function of the amplitude of one voxel vs. time.  That is a
> time domain function with a sampling rate of 1/TR.  Before considering any
> upsampling, the onset vector will be a stream of delta functions and zeros,
> at the sampling rate 1/TR.  The hemodynamic response function (hrf) must
> also be at the same sampling rate.  The hrf acts as a filter which is
> convolved with the delta stream.  Now consider passing just one delta
> function through the filter.  Then the hrf is an impulse response of a
> finite-length digital filter.  The output of the convolution has to be band
> limited to half the sampling rate, i.e., half of 1/TR.  Otherwise there
> will be aliasing. Since I don't see any low pass filters accomplishing
> this, then the filter hrf itself has to have a frequency response whose
> passband has energy only below half the sampling rate.  I tried the
> function spm_hrf, and the filter it makes does appear to have a frequency
> response with this property.

Yes, but spm_hrf is only a pretty-good guess of what your particular subject's
hrf might look like, so consider it a rough approximation.


I ran spm_hrf with an input RT and then freqz
> in Matlab with a sampling rate of 1/TR.  So it seems that the way the basis
> function combinations are designed for spm_hrf, they create a filter with
> the right properties.  If anyone has more details on how this filter design
> works, I would like to read about it, since it is impressive that the blood
> response and filtering properties can be taken care of in one filter
> design.  Next, I am considering the upsampling.  I tried the hrf obtained
> in spm_hrf before the downsampling occurs. This means that the hrf has been
> upsampled by L=16.  Now the repetition of the spectrum occurs at L*Fs,
> where Fs is the sampling rate.  I looked at the output in Matlab and the
> bandlimiting was such that there was only a tiny bit of energy out of the
> original band, so that there is still no aliasing except for some possible
> tiny amount of noise.

Convolution in the time domain is straight multiplication in the frequency
domain.  So, after the convolution, the resulting signal cannot have any
frequencies that are not present in _both_ input signals.


This is why the hrf can be downsampled at the end of
> spm_hrf without requiring a decimation filter.  The next question is the
> jitter.  spm_hrf allows for delaying the filter's response.

I don't follow here-- are you referring to the use of the 
time-derivative of the
canonical hrf when you mention "delaying", or modification of the hrf 
using the
"params" argument of spm_hrf(TR, params)?

My question,
> which I will also look into independently, is how this affects the spectrum
> of the output of the convolution, or does it affect the phase only and not
> the magnitude?  I would also like to know how this is related to the idea
> of parametric modulation, which I am going to begin studying.
> Linda Seltzer
> [log in to unmask]
>
>> SPM mostly does things on the time side, not the frequency side, so a lot
>> of what you're describing isn't directly there but rather implicit only.
>> The viewpoint is conceptually linear models (from statistics), rather than
>> dsp, though the usual mathematical equivalences still exist of course, and
>> in the literature people look at things both ways.
>>
>> Upsampling/downsampling occurs when the model is specified.  The canonical
>> hemodynamic response function and the onset times are specified at a
>> relatively fine time scale, and that scale is used in the convolution you
>> refer to, then downsampled to the true experimental sampling rate.  In
>> spm_fMRI_design.m (SPM2, version % @(#)spm_fMRI_design.m 2.34 Karl Friston
>> 03/01/30), the code is
>>
>> % create convolved stimulus functions or inputs
>> %=======================================================
>> % Get inputs, neuronal causes or stimulus functions U
>> %-------------------------------------------------------
>>         U = spm_get_ons(SPM,s);
>> % Convolve stimulus functions with basis functions
>> %-------------------------------------------------------
>> [X,Xn,Fc] = spm_Volterra(U,bf,V);
>> % Resample regressors at acquisition times (32 bin offset)
>> %-------------------------------------------------------
>> try
>>         X = X([0:(k - 1)]*fMRI_T + fMRI_T0 + 32,:);
>> end
>>
>> Best,
>>
>> Stephen J. Fromm, PhD
>> Contractor, NIMH/MAP
>> (301) 451--9265
>>
>>
>>
>>
>> From: Linda Seltzer
>> Sent: Thu 2008-01-17 6:56 PM
>> To: Stephen J. Fromm; [log in to unmask]
>> Subject: Re: [SPM] Question on sampling rates, aliasing, noise and jitter
>>
>>
>> I'm also talking about where in the code each step is handled in SPM.
>> > On Wed, 16 Jan 2008 12:31:57 -0800, Linda Seltzer
>> > <[log in to unmask]>
>> > wrote:
>> >
>> > >Does anyone written a paper or tutorial applying the well-known
>> concepts
>> >
>> > of
>> > >sampling, aliasing, Nyquist theorem, etc. to data analysis in SPM?
>> > >Specifically, TR and sampling rate; spectrum of voxel data vs time,
>> the
>> > >same with noise; how the concept of jitter relates to this; upsampling
>> > to
>> > >convolve the delta fns with hemodynamic response (for the regressors)
>> > and
>> > >then downsampling the signal back, when each step is done.  I would
>> like
>> >
>> > to
>> > >see a description from a signal processing engineering viewpoint of how
>>
>> > all
>> > >of the sampling rates and bandwidths are handled.  If there is nothing
>> > >formal, but someone has made some notes deriving all of this, I would
>> > like
>> > >to read them and learn more.  I would like to have an overview of all
>> of
>> > >the sampling rates and spectrums involved.
>> > >Linda Seltzer
>> > >[log in to unmask]
>> > >[log in to unmask]
>> >
>> >=========================================================================
>> >
>> > IMHO a useful way to look at this is to not restrict consideration to
>> SPM,
>> >
>> > but think about fMRI analysis in general (restricted to the case of
>> > univariate statistics, as opposed to multivariate methods).  Then look
>> for
>> >
>> > articles in the published literature via Google or PubMed.  There's
>> quite
>> >
>> > a bit out there.
>> >
>> > Just to take an interesting an relevant example, the problems of aliased
>>
>> > noise from cardiac and respiratory cycles is discussed quite
>> frequently.
>> >
>>
>



----------------------------------------------------
Ken Roberts
Woldorff Laboratory
Center for Cognitive Neuroscience, Duke University
(919) 668-1334
----------------------------------------------------