Dear Karl,
thank you for your answer.
> The model used for serial correlations
> in standard MRI timeseries is usually described as "an AR model plus
> white noise". Practically, this means there are two covariance
> components (C{1} and C{2} in your code below), where the first is
> the identity or white component and the second is an AR component.
> This is simply a covariance component based upon an AR(1) process
> with an AR coefficient of 0.2 (or user-specified).
I understand, but as far as I can tell, SPM12 doesn't seem to be actually doing this.
According to your description (consistent with the paper), the components should be the identity matrix (truncated to the first 6 x 6 entries)
>> eye(6)
ans =
1 0 0 0 0 0
0 1 0 0 0 0
0 0 1 0 0 0
0 0 0 1 0 0
0 0 0 0 1 0
0 0 0 0 0 1
and the autocorrelation matrix of an AR(1) process with parameter 0.2
>> toeplitz(0.2 .^ (0 : 5))
ans =
1.0000 0.2000 0.0400 0.0080 0.0016 0.0003
0.2000 1.0000 0.2000 0.0400 0.0080 0.0016
0.0400 0.2000 1.0000 0.2000 0.0400 0.0080
0.0080 0.0400 0.2000 1.0000 0.2000 0.0400
0.0016 0.0080 0.0400 0.2000 1.0000 0.2000
0.0003 0.0016 0.0080 0.0400 0.2000 1.0000
But when I set up a first-level model with standard parameters, including modeling of serial correlations, these are the resulting covariance components:
>> full(SPM.xVi.Vi{1}(1:6, 1:6))
ans =
1.0000 -0.0000 -0.0401 -0.0160 -0.0048 -0.0013
-0.0000 1.0000 -0.0000 -0.0401 -0.0160 -0.0048
-0.0401 -0.0000 1.0000 -0.0000 -0.0401 -0.0160
-0.0160 -0.0401 -0.0000 1.0000 -0.0000 -0.0401
-0.0048 -0.0160 -0.0401 -0.0000 1.0000 -0.0000
-0.0013 -0.0048 -0.0160 -0.0401 -0.0000 1.0000
>> full(SPM.xVi.Vi{2}(1:6, 1:6))
1.0000 0.4000 0.1201 0.0320 0.0080 0.0019
0.4000 1.0000 0.4000 0.1201 0.0320 0.0080
0.1201 0.4000 1.0000 0.4000 0.1201 0.0320
0.0320 0.1201 0.4000 1.0000 0.4000 0.1201
0.0080 0.0320 0.1201 0.4000 1.0000 0.4000
0.0019 0.0080 0.0320 0.1201 0.4000 1.0000
It turns out that the mean of these two matrices is actually the AR(1) matrix:
>> full(SPM.xVi.Vi{1}(1:6, 1:6) + SPM.xVi.Vi{2}(1:6, 1:6)) / 2
ans =
1.0000 0.2000 0.0400 0.0080 0.0016 0.0003
0.2000 1.0000 0.2000 0.0400 0.0080 0.0016
0.0400 0.2000 1.0000 0.2000 0.0400 0.0080
0.0080 0.0400 0.2000 1.0000 0.2000 0.0400
0.0016 0.0080 0.0400 0.2000 1.0000 0.2000
0.0003 0.0016 0.0080 0.0400 0.2000 1.0000
but the complement
>> full(SPM.xVi.Vi{2}(1:6, 1:6) - SPM.xVi.Vi{1}(1:6, 1:6)) / 2
ans =
0 0.2000 0.0801 0.0240 0.0064 0.0016
0.2000 0 0.2000 0.0801 0.0240 0.0064
0.0801 0.2000 0 0.2000 0.0801 0.0240
0.0240 0.0801 0.2000 0 0.2000 0.0801
0.0064 0.0240 0.0801 0.2000 0 0.2000
0.0016 0.0064 0.0240 0.0801 0.2000 0
is far from being the identity matrix. According to the comments in spm_Ce, this is the derivative of the AR(1) matrix with respect to a at a = 0.2:
% C{i} = AR(a) - a*dAR(a)/da;
% C{i + 1} = AR(a) + a*dAR(a)/da;
Could you shed light on this behavior?
Thank you,
Carsten
|