Dear Guillaume & Will,
I noticed a minor inaccuracy in my code transforming AR coefficients into autocorrelation functions. According to comparison with the analytic solution for specific coefficient values, I the corrected version has an absolute error on the order of 1e-16. The updated results have minor changes, but the basic observations are the same.
There are 13256 in-mask voxels. Of these, 16 voxels have estimated AR coefficients whose roots' magnitudes are not all below 1 (describing nonstationary AR processes).
The autocorrelation functions of the remaining 13240 voxels are shown superimposed in
https://drive.google.com/open?id=1y6veRsvCYcLnEe-nDpx8FUn7ziPGGW8P
As you notice, a non-neglibile number of them have autocorrelation functions that do not decay to 0 over the course of 128 scans (= 256 s).
The spatial distribution of autocorrelation over lags 0 to 127 in slice 13 (of 23) containing 1150 in-mask voxels is shown in
https://drive.google.com/open?id=1XB2q7tWqUjlHVUMJLjiiidobKkv6Y089
It is apparent that the voxels with long-range autocorrelation lie mainly at the edge of the brain mask.
Guillaume, you asked for time series. The image
https://drive.google.com/open?id=1swKtR64OI-QIY-Q8ZMcEZenGCYJPA7CP
shows the underlying BOLD measurements, in the top panel for the 16 voxels with nonstationary AR, and in the lower panel for the 71 voxels where the autocorrelation at lag 127 is larger than 0.1. I can't say that I notice anything special or common to these timeseries.
For the moment I will proceed with my data analysis after excluding the 16 + 71 = 87 "weird" voxels. But I think it would be useful if you could look into this further. In particular, in my opinion an AR estimation method should not produce coefficients describing a nonstationary process. Do you agree?
Thank you!
Best,
Carsten
----- Original Message -----
> From: "Guillaume Flandin" <[log in to unmask]>
> To: "Carsten Allefeld" <[log in to unmask]>
> Cc: [log in to unmask], "William Penny (PSY)" <[log in to unmask]>
> Sent: Wednesday, 28 November, 2018 11:35:53 AM
> Subject: Re: [SPM] Strange AR processes from Bayesian estimation
>
> Dear Carsten,
>
> All the voxels exhibiting pathological behaviour seem to be at the
> boundary of the brain mask so, as you say, simply discarding them is
> probably the easiest thing to do. To understand what is going on, it
> would be useful to extract and display their time series.
> I copy this email to Will so that he can add any further comments or
> advice on the spatial noise prior.
>
> Best regards,
> Guillaume.
>
>
> On 27/11/2018 19:20, Carsten Allefeld wrote:
> > Dear Guillaume,
> >
> > thanks for replying!
> >
> >> These results were obtained with which spatial noise prior option?
> >> The
> >> interface lists five options: UGL, GMRF, LORETA, Tissue-type and
> >> Robust.
> >
> > I used the default, "UGL". Which one would you recommend?
> >
> > Options in detail:
> > fmri_est.method.Bayesian.space.volume.block_type = 'Slices';
> > fmri_est.method.Bayesian.signal = 'UGL';
> > fmri_est.method.Bayesian.ARP = 6;
> > fmri_est.method.Bayesian.noise.UGL = 1;
> > fmri_est.method.Bayesian.LogEv = 'No';
> > fmri_est.method.Bayesian.anova.first = 'No';
> > fmri_est.method.Bayesian.anova.second = 'No';
> > fmri_est.method.Bayesian.gcon = struct('name', {}, 'convec', {});
> >
> >> Could you show a map of where the voxels you are concerned about
> >> are?
> >
> > Attached are plots of the absolute value of the autocorrelation at
> > lags 0 to 127 in a middle slice (#13).
> > Comparison with the coregistered T1 indicates that they are located
> > mainly at the outer edge of gray matter (maybe meninges), but also
> > frontally and posteriorly slightly into the longitudinal fissure.
> >
> > That suggests I should simply exclude these voxels from further
> > analysis.
> >
> > Do you have a suggestion which criterion to use?
> > The data-based threshold (2.02e-9) discards more than half of the
> > brain, and any other threshold seems arbitrary.
> >
> > Best,
> > Carsten
> >
> >
> >>
> >> Best regards,
> >> Guillaume.
> >>
> >>
> >> On 27/11/2018 17:11, Carsten Allefeld wrote:
> >>> Hello all,
> >>>
> >>> I'm interested in getting local estimates of temporal
> >>> autocorrelation in SPM, and for that purpose used Bayesian
> >>> 1st-level estimation.
> >>> The fMRI data I used to test that have 3 sessions of 128 scans at
> >>> a
> >>> TR of 2 s and 64x64x23 voxels of size 4x4x4 mm, unsmoothed, of
> >>> which approximately 13,000 are within brain.
> >>>
> >>> I then extracted the AR coefficients (order 6) for the first
> >>> session (Sess1_AR_0001.nii to Sess1_AR_0006.nii) and used the
> >>> Yule–Walker equations iteratively to obtain the corresponding
> >>> autocorrelation function across lags 0 to 128.
> >>>
> >>> The results are strange (see attached plot):
> >>> – In 16 voxels the AR coefficients describe a non-stationary
> >>> process. After excluding them:
> >>> – At lag 127, 7906 voxels have an autocorrelation > 1e-6, 1651
> >>> voxels > 1e-3, and 113 voxels > 0.1.
> >>> – The largest negative autocorrelation at lag 127 is -2.02e-9. If
> >>> I
> >>> take that as an indicator of numerical/estimation precision,
> >>> there
> >>> are 8185 voxels where the autocorrelation at lag 127 is different
> >>> from 0 (> +2.02e-9).
> >>>
> >>> This makes me suspect that the AR estimation in "Bayesian
> >>> 1st-level" is not very reliable. Is there something I might have
> >>> done wrong?
> >>>
> >>> Is there a recommended postprocessing for the AR coefficients or
> >>> autocorrelation functions?
> >>> I thought about tapering à la FSL, or clustering as a crude form
> >>> of
> >>> spatial regularization.
> >>> Or should I simply exclude voxels with unbelievably long-range
> >>> autocorrelation?
> >>>
> >>> Thank!
> >>>
> >>> Best,
> >>> Carsten
> >>>
> >>
> >> --
> >> Guillaume Flandin, PhD
> >> Wellcome Centre for Human Neuroimaging
> >> UCL Queen Square Institute of Neurology
> >> London WC1N 3BG
> >>
>
> --
> Guillaume Flandin, PhD
> Wellcome Centre for Human Neuroimaging
> UCL Queen Square Institute of Neurology
> London WC1N 3BG
>
|