Nugent, Allison C. (NIH/NIMH) [E] wrote:
> % Set VResMS scalefactor as 1/trRV (raw voxel data is ResSS)
> %-----------------------------------------------------------
> VResMS.pinfo(1) = 1/xX.trRV;
> VResMS = spm_create_vol(VResMS,'noopen');
>
> If I use spm_vol to load in ResMS.img, the resulting image handle (I'll
> call it ResMSimage) has the field ResMSimage.pinfo = [1; 0; 0]. If I
> look at SPM.VResMS.pinfo, however, the result is [1/xX.trRV; 0; 0]. So,
> all data pulled from SPM.VResMS is scaled by that factor.
Ah, well spotted... I'm very sorry for suggesting you had the wrong
ResMS.img or the wrong voxel...
I've just had a look in the archives and found this (its ref'd posts
are informative too)
http://www.jiscmail.ac.uk/cgi-bin/webadmin?A2=ind0603&L=SPM&P=R50240
> In the spm_spm.m documentation, it states that ResSS/xX.trRV = ResMS.
> So, does that mean that ResMS.img should more appropriately be named
> ResSS.img?
You're right, unless you access ResMS.img only via SPM.VResMS in which
case it will be MS. Worth knowing...
Also, I've just checked under SPM5, and there the correct 1/trRV
scaling factor is saved in V.pinfo from V=spm_vol('ResMS.img') as well
as in SPM.VResMS.pinfo, so with SPM5, the saved ResMS.img really is
the mean square (at least for NIfTI-compatible software, including
e.g. FSL as well as spm_vol).
At risk of overstating the case, always use SPM.VResMS instead of
reading ResMS.img directly and correcting, since then any code you
write will work under SPM5 as well as SPM2.
> spm_spm.m also states that xX.trRV is the trace of R*V -
> excuse my lack of statistical knowledge, but what is that
> non-mathematically?
Well, mathematically (!) it's described on p37 of the appendix of ch7
of HBF2:
http://www.fil.ion.ucl.ac.uk/spm/doc/books/hbf2/pdfs/Ch7.pdf
Also by Worsley and Friston: http://dx.doi.org/10.1006/nimg.1995.1023
(and presumably at length in their reference Seber 1977...)
Non-mathematically, I don't think I can give you a better answer than
it's what you have to divide the sum-squared error by to get an
unbiased estimate of the variance.
As pointed out in HBF, Tr(RV) just gives the standard J-p (or n-k more
commonly...) degrees of freedom if V is an identity matrix. If V is a
more complicated matrix describing a non-iid error structure then
this is no longer the case.
Note also (personally I don't think this is too clear in HBF) that the
Satterthwaite approximation to the denominator degrees of freedom for
the actual test statistic (which isn't simply Tr(RV), but
Tr(RV)^2/Tr(RVRV) instead) also reduces to the standard error degrees
of freedom if V is an identity matrix, since Tr(RVRV)=Tr(RR)=Tr(R),
where the last bit follows because R, "the residual forming matrix",
projects onto a subspace, and projecting onto space for the second
time does nothing since the first result is already in that space. I'm
not making this too clear either, am I?
Roughly: Tr(RV) and Tr(RV)^2/Tr(RVRV) replace the familiar degrees of
freedom in models with more complicated error structures (for example,
a two-group t-test with unequal group variances).
Hope that makes some small amount of sense (which is all it does to me),
Ged.
|