W.D. Penny, J. Ashburner, S. Kiebel,
R. Henson, D.E. Glaser, C. Phillips and K. Friston.
Wellcome Department of Cognitive Neurology,
University College London.
Statistical Paramemetric Mapping (SPM) is a method for detecting regional, task-related changes in brain activity from neuroimaging data. In this paper we provide a guide to the methodological literature on SPM, highlighting key papers and those of historical interest. Complementary approaches, which look at distributed task-related changes, such as multivariate analysis and effective connectivity, are also covered.
The fundamental text on Statistical Parametric Mapping (SPM) is the Human Brain Function (HBF) book [1]. This comprises a basic overview [2], and then chapters on the spatial transformation of images [3], characterizing brain images with the general linear model [7], making statistical inferences [5], characterizing distributed functional systems [3], characterizing functional integration [2] and a chapter looking at the basic types of study design [2].
In addition there are a number of general introductory articles and commentaries. In [9], the distinction is made between (i) making maps of functional responses in the brain and (ii) discerning the principles underlying their organization, where both approaches are considered fruitful. A more detailed and up-to-date review of functional brain imaging appears in [10]. A comprehensive tutorial on Statistical Parametric Mapping (SPM), experimental design and functional and effective connectivity appears in [11]. Researchers new to the field are encouraged to read [12].
For each section, we provide a list of methods papers and describe the contribution of each. The lists are not exhaustive and some papers appear in more than one section. We emphasise that this bibliography is focussed on the methodological literature on SPM from the Wellcome Department of Cognitive Neurology and is by no means meant to be an exhaustive bibliography of the more general area of neuroimaging analysis.
Before subsequent processing, images must be realigned and spatially normalised to remove movement artifacts or to put data from different subjects into a common anatomical frame. These methods are reviewed in [1].
The idea of using structural information in functional data to register images was first presented in [2]. This approach provided a linear affine normalisation of PET images (an affine transformation preserves parallel lines).
In [3], the authors consider the registration of PET images using linear translations to correct for positional shifts and nonlinear resampling to account for morphological differences among subjects.
Friston et al. [4] propose an algorithm for removing movement-related artifacts in fMRI. The method is composed of two stages. In the first stage, each image is transformed via a 6-parameter rigid-body transformation so as to be as similar as possible to the first image in the time series. This corrects for changes in instantaneous position of the body in the scanner. A second stage, based on an autoregressive-moving average (ARMA) model, corrects for changes in the recorded signal due to the sequence of body positions. This second component can be attributed to the spin-excitation history of the object being scanned.
In [5], the authors describe a general framework for registering images that involves both spatial and intensity transformations using a Gauss-Newton optimisation method. This paper introduced spatial basis functions for defining warps and is the foundation for the nonlinear normalisation used in SPM.
In [6], the authors present a method for co-registration of brain images from different modalities. Instead of matching the images directly, one performs intermediate within modality registration to two template images that are already in register. This paper also describes a segmentation method that is refined in [7]. The main contribution of this work is that it uses spatial priors from previously segmented images. In [8], the authors consider the co-registration of functional PET images with structural fMRI images (using the above method). The algorithm is compared to the Automatic Image Realignment (AIR) method. Ashburner [9] et al. extend the standard affine registration method by incorporating information about the variability in the shape and size of heads in the form of Bayesian priors. The paper is is important because it brought the spatial basis function approach into a Bayesian framework. In [10] and [11], the authors present a nonlinear spatial registration algorithm which uses a set of low spatial frequency basis functions (discrete cosines), and a smoothness constraint which is implemented using a Bayesian Maximum a Posteriori (MAP) estimator. In this short review article [12], the authors explain the motivation behind the development of the above registration methods - highlighting the need to remove motion artifacts and to bring data from many subjects into the same space for anatomical localization and intersubject averaging.
In [13], the authors present a finite element (hence high dimensional) approach to nonlinear image registration. Symmetric priors are imposed to ensure smooth transforms and a MAP solution is found. The work in [14] extends the algorithm from 2 to 3-dimensional images.
Registration methods are also central to the field of computational neuroanatomy which looks at differences in the structure (rather than function) of different human brains.
In [15], the authors illustrate a method for identifying macroscopic anatomical differences among the brains of different populations of subjects. The method involves spatial normalization of structural MR images (affine transformation and a nonlinear deformation, based on a discrete cosine basis set) and a canonical variates analysis (see section 8) to assess significant differences in deformations between groups. This is known as Deformation-Based Morphometry (DBM).
A related procedure is known as Voxel-Based Morphometry (VBM) [7] and was first described in [16]. VBM throws away global differences and focuses on local differences. The normalised images are first classified as being either grey matter, white matter or Cerebro-Spinal Fluid (CSF), using a refined version of the algorithm described in [6]. Following this, the data are smoothed before performing voxel-wise statistical tests.
The General Linear Model (GLM) and the theory of Gaussian Random Fields (GRFs) (or, simply, `Random Field Theory') lie at the heart of statistical parametric mapping. SPM was first formally introduced in [1] which described a GLM that used global activity as a confounding covariate. The defining features of SPM were, at that time, the construction of maps of a statistical parameter obtained by treating each voxel separately in a mass univariate approach. In [2], the problem of multiple comparisons was first addressed using the theory of stochastic processes and estimates of the smoothness of the images (this work was later generalised - as described in section 4). The two elements of SPM, namely GLM and GRF, are described together in [3].
The validity of many of the assumptions underlying SPM (eg. that error variance changes over voxels, that voxel error terms are Gaussian-distributed), is addressed in [4].
In [5], Holmes considers a number of issues underlying the statistical analysis of PET data. This covers GLMs and GRFs and again addresses the issue of homoskedacity (equal error variance) across voxels. A Markov Random Field (MRF) model is also proposed, for segmenting active areas from non-active areas. Finally, a non-parametric approach is proposed for the analysis of simple activation studies. This removes the need to assume that voxels are Gaussian-distributed, and develops a `distribution-free' procedure which is always valid, albeit at the cost of increased computation. The non-parametric analysis is also published in [6].
An introduction to the GLM and its application to neuroimaging (including to fMRI) is given in [7].
In SPM, the GLM is used to relate the activity at each voxel to the experimental task. Significant task-related changes in voxel activity are assessed using F, t or z statistics. To guard against false positives, a correction for multiple comparisons (ie. comparisons over the voxels) must be made. This is achieved using Gaussian Random Field (GRF) theory.
The core paper in this area is [1]. This provides an estimate of the p-value for local maxima of Gaussian, t, and F fields over search regions of any shape or size in any number of dimensions. This unifies earlier results on 2D [2] and 3D [3] images.
The above analysis requires an estimate of the smoothness of the images (images have some inherent smoothness but are also smoothed using Gaussian kernels (i) to allow for variability in subjects neuroanatomy and (ii) to ensure the validity of GRF theory). In [4], Poline et al. estimate the dependence of the resulting SPMs on the estimate of this parameter. Whilst the applied smoothness is usually fixed, [5] propose a scale-space procedure for assessing significance of activations over a range of proposed smoothings. In [6], the authors implement an improved estimation procedure for estimating smoothness in Gaussianised t-fields.
Another approach to assessing significance is based, not on the height of a cluster of activity, but on its spatial extent [7].
In [8], the authors consider a hierarchy of tests that are regarded as voxel-level, cluster-level and set-level inferences. These inferences have increasing power but decreasing spatial localization. The notion of set-level inferences allows one to assess the significance of distributed (rather than local) activations.
If the approximate location of an activation can be specified in advance then the significance of the activation can be assessed using the spatial extent or volume of the nearest activated region [9]. This test is particularly elegant as it does not require a correction for multiple comparisons.
In more recent work [10], Kiebel et al. propose a model for the analysis of PET and fMRI data that allows incorporation of prior neuroanatomical knowledge. Specifically, a set of basis functions placed on the grey matter surface obtained from a T1-weighted MRI image, allows one to specify a spatial smoothing on the data that both varies over the image and is anatomically informed. This improves the spatial resolution and the sensitivity of the resulting SPM analysis.
The mathematical basis of GRF theory is described in a series of peer-reviewed articles in statisical journals [11, 12, 13, 14, 15].
There are a number of useful reviews of different approaches to experimental design in neuroimaging, the most introductory being [1], with [2] being the most tutorial, and [3] being the most up-to-date (including a discussion of event-related fMRI designs).
The first analyses of PET images used a `cognitive subtraction' paradigm which relies upon the assumption of pure insertion - this assumes that a new cognitive component can be purely inserted ie. without affecting the expression of previous ones. This assumption is often invalid, however, as shown in [4] where it was demonstrated that interaction terms (between the new and previous cognitive components) are often significant.
More powerful experimental designs, such as parametric designs (involving continuous dependent variables, such as time) and factorial designs (which explicitly look at interaction terms) can be facilitated with the GLM.
The first published study using a parametric analysis was by Price et al. [5] where activity in primary auditory cortex was found to be linearly related to the word presentation rate (this was in contrast to, for example, Wernicke's area which responded solely to the occurence of a word rather than its presentation rate).
In [6] Buchel et al. extend the above analysis by allowing stimulus or task parameters to be nonlinear functions of the hemodynamic response. This was implemented using second-order polynomial expansions and applied to PET data. In [7] the method was implemented in the context of the General Linear Model and SPM's, based on omnibus F-statistics (eg. to test for any significant linear and/or nonlinear effect), were used to test for local linear or nonlinear dependencies in fMRI data.
The first factorial analysis of PET data, for example, found a significant interaction between motor activation and time during a paced finger opposition task [8]. The first psycho-pharmaceutical study using a factorial design was [9].
Despite these advances, the experimental paradigms so far described are restricted in that the inferences made relate to the particular subjects scanned. This is overcome in the Random Effects (RFX) Analysis procedure where inferences can be made about the populations from which the subjects are sampled (eg. males/females). This can be implemented in SPM using a two-level analysis procedure described in [10].
Conjunction analysis [11] looks for areas of activation that are common to a number of tasks. For example, in a phonological retrieval task whether subjects named words, objects, letters or colors there is always activation of the left thalamus (and a number of other areas). In a further paper [12], the authors show that conjunction analysis can be applied to data from multiple subjects to make inferences about populations. This work relied on a technical development by Worsley [13].
The relative merits of RFX versus conjunction analysis and the larger issues are discussed in [14].
To apply SPM to fMRI time-series it was necessary to account for two fundamental characteristics of fMRI data: (i) that the Hemodynamic Reponse Function (HRF) is transient, delayed and dispersed in time and (ii) that the sampling interval is shorter than the time constant of the HRF, making the observed time-series highly correlated. These issues were addressed in [1] where the HRF was modelled using a Poisson function and the intrinsic correlation was measured after first removing the signal components that were phase-locked to the stimulus (ie. the extrinsic components).
In a later paper [2], the authors consider a temporal smoothing of the fMRI time-series (using a Gaussian filter) so as to induce a known autocorrelation structure. This is then used to adjust the degrees of freedom used when making inferences from the GLM. In the further development of this model [3], Worsley and Friston introduced a general procedure for unbiased estimation of the error variance term. This series of papers culminates in [4], where it was shown that whitening (estimating the autocorrelations and then removing them - a procedure favoured by a number of other researchers) can result in a large bias in the resulting statistical inferences (ie. many false positives).
At this point in time, the favoured approach for handling correlations in fMRI time seris was therefore `smoothing' rather than `whitening'. More recent developments, however, have changed this view somewhat. With the advent of Hierarchical Bayesian Models [5],[6] a new approach, which might be termed `Bayesian whitening' is now the preferred method. This allows for both the error variance term and the intrinsic correlation term to be estimated in an unbiased manner.
A second key feature of fMRI data, as opposed to PET data, is that data can be collected more than once from any given subject. This allows for a quantification of the `within-subject' variability, which is to be contrasted with the `between-subject' variability. Both sources of variability can be accounted for the random effects analysis procedure discussed in the previous section. In earlier work Buchel et al. [7] showed that fMRI data from multiple subjects could be pooled in a `fixed-effects' analysis.
Event-related fMRI (efMRI) is the use of fMRI to characterize and detect transient hemodynamic responses to brief stimuli or tasks, and is analagous to the the study of event-related potentials in electrophysiology. A review of the statistical issues underlying efMRI is presented in [1] with special emphasis on determining the optimal experimental design for testing a given hypothesis.
In [2], temporal basis functions were used to model the `early' and `late' components of the evoked response (in the context of a block-design). These took the form of exponentially modulated sine functions. The model can detect `conventional' activations, where both components are expressed to the same degree, and differential activations such as are involved in tasks that do not require sustained attention.
Josephs et al. [3] proposed staggering stimulus times relative to scan acquisition times so as to realize a higher effective sampling rate. This effectively allowed, for the first time, whole-brain EPI scans to be used in an efMRI context. In [4], the GLM is employed to detect responses that are different in both magnitude and latency.
Friston et al. [5] characterise the evoked hemodynamic response using Volterra kernels. This allows for nonlinear components of the response to be modelled, such as the saturation of responses at high presentation rates. The linear models described above may be viewed as a first order approximation to this nonlinear model.
In [6], the authors consider stochastic experimental designs where an event (eg. a stimulus) has a certain probability of occuring at any given time point. This is to be contrasted with deterministic designs in which the timing of events is fixed. They then make the distinction between stationary stochastic designs where the probability is fixed and nonstationary designs where this probability may vary over time. They conclude that block designs are generally the most efficient for detecting differential responses, whereas designs including null events are the most efficient for detecting transient responses.
In [7], Hopfinger et al. looked at the sensitivity of event-related fMRI analyses to voxel size, spatial and temporal smoothing parameters and the basis set used to characterise the HRF.
In neuroimaging, multivariate analysis is concerned with characterising the relations between different areas during brain function; ie. characterising distributed rather than local functional systems. The general approach and a number of instances of it are described in a tutorial in [2].
The simplest and most frequently applied multivariate procedure is Principal Component Analyis (PCA) which can be implemented using Singular Value Decomposition (SVD). As applied to neuroimaging data PCA finds a set of spatial modes that are mutually uncorrelated both spatially and temporally. The modes are also ordered according to the amount of variance they explain. In [2], the authors applied PCA to the analysis of PET data. By comparing the temporal expression of the first few spatial eigenmodes with the variation in experimental condition a distributed functional system associated with the activation could be identified.
A more sophisticated use of PCA occurs in the context of Generalised Eigenimage Analysis (GEA) [3] where the principal component is found which is maximally expressed in one experimental condition/group and minimally expressed in another (eg. control and patient groups).
Friston et al. [4] apply Multi-Dimensional Scaling (MDS) to data from schizophrenic and control subjects. The MDS procedure maps functionally connected areas into similar positions in a 2D or 3D map. This results in a `functional map' rather than an `anatomical map'. For the data studied, abnormal connectivity patterns were observed in a schizophrenic group. The MDS algorithm works by preserving pairwise distances between data in the original high dimensional space and data in the low dimensional (2D or 3D) visualisation space. If the distance metric is chosen to be the Euclidian distance then MDS turns out to be equivalent to projecting the data onto the first 2 or 3 principal components.
A second useful multivariate procedure is the Multivariate Linear Model (MLM) (see eg. [3]) which maps one set of variables to another, ie. there are multiple independent variables (inputs) and multiple dependent variables (outputs). This is to be contrasted with the Linear Model (LM) or General Linear Model (GLM) which has multiple independent variables but only a single dependent variable (eg. the model for a voxel in a standard SPM analysis). The multivariate analysis of covariance (ManCova), for example, may be implemented using an MLM - the MLM is to ManCova what the LM is to the analysis of variance (Anova).
In [5], the authors use a MLM in an analysis of PET data. Wilk's Lambda is used to assess the significance of the amount of signal variance explained (in proportion to the amount of variance unexplained - the `noise'). Canonical Variates Analysis (CVA) is then used to find the associated `canonical' images. These are modes that are expressed most in the signal and least in the noise - CVA is therefore a generalisation of GEA which allows for multiple groups and parametric designs. In [6], the approach is applied to event-related fMRI where task-dependent HRFs were identified. In [5], the approach is applied to the analysis of evoked responses in EEG and MEG. Worsley et al. [7] apply a MLM and CVA to PET and fMRI data. The paper extends previous work by allowing for temporal correlations - hence the application to fMRI.
In [8], the authors derive a nonlinear PCA algorithm, implemented in a neural network, which finds a number of sources which are orthogonal in time. These sources generate the data via first and second order spatial modes - these may be said to constitute the generative network. The sources themselves are derived from the data using a separate `recognition' network. Standard PCA is recovered in the absence of second order spatial modes - or when the generative network implements the inverse function of the recognition network.
In [9], Ashburner et al. apply a cluster analysis procedure to PET radiotracer time series. This identified a number of archetypal time series (tracer kinetics) and associated areas.
In [10], Friston offers a critique of Independent Component Analysis (ICA). As applied to neuroimaging data, ICA attempts to find a set of spatial modes that are spatially independent. In contrast to PCA, which finds spatially uncorrelated modes, the ICA modes are sparser and more spatially localised. As with many other multivariate procedures (eg. PCA, MDS, cluster analysis), ICA is data-driven rather than hypothesis-driven and is therefore best suited to the generation of new hypotheses rather than the testing of existing ones.
Effective connectivity is defined in the neuroimaging community [1] as the influence one neural system exerts over another, either at a synaptic or cortical level. It is defined both by a neuroanatomical model, defining which areas are connected, and a mathematical model detailing how they are connected (linearly, nonlinearly, with or without temporal evolution etc.). A tutorial on effective connectivity appears in [2]. Friston et al. [1] initially focussed on linear time-varying connections.
Effective connectivity is to be contrasted with functional connectivity. In [3], the authors present a synthesis of functional and effective connectivity, defining functional connectivity as the temporal correlation between spatially remote neurophysiological processes. Operationally, this is implemented using PCA (see previous section) where the different spatial modes together account for the observed correlations in the data. In contrast, effective connectivity posits directed (and therefore causal - rather than associative) relations between variables. The differences are illustrated on PET and fMRI data and nonlinear time-varying effective connectivity models are illustrated.
In [4], the interactions between V1 and V2 were characterized using a nonlinear model of effective connectivity. This extended the linear model by allowing for modulatary connections, ie. the activity in V1 is linearly dependent on activity in other areas and dependent on the product of activity from other areas and the intrinsic activity in V1 (this is the modulatory term).
To characterize the interactions between multiple areas a `Structural Equation Model (SEM)' or a `Path Model' is required. This is specified by a directed graph. In its simplest incarnation a number of path coefficients define the linear relation between nodes on the graph (ie. neuroanatomical areas). These coefficients are then set so as to explain the observed covariance among the variables, as described in [5]. SEM was introduced to neuroimaging by McIntosh and Gonzalez-Lima [6].
Buchel and Friston [7] propose an SEM with time-dependent connections which are estimated using a Kalman filter algorithm. The approach is reviewed in [8] which also has a general discussion on effective connectiviy.
In [9], Friston et al. introduced the concept of Psychophysiological Interactions (PPIs) to neuroimaging. This model explains the activation at any given voxel as an interaction between a psychological variable and a physiological variable (plus the main effects of each). This is a powerful cross-fertilisation of two related concepts (i) factorial designs - which look at the interactions between two psychological variables and (ii) effective connecticity - which looks at the interactions between two physiological variables (ie. the activations in different anatomical areas).
In a study of associative learning of objects and their positions [10], increases in effective connectivity between distinct cortical systems specialised for spatial and object processing were observed (in addition to the expected repetition supression effect).
In a study of visual attentional mechanisms [11] measures of effective connectivity based on a nonlinear system identification showed that the effective connectivity between V2 and V5/MT is dependent on activity in parietal cortex. This provides evidence for the role of top-down `backwards' connections in visual processing.
The table below shows a list of those papers that have been cited at least 50 times in scientific journals - naturally, older papers appear more often than newer ones. The information was collated from the Institute for Scientific Information (ISI) citation database. Note that this database does not include volumes 1 and 2 of the Human Brain Mapping journal.
A list of SPM methods papers with at least 50 citations. The number of citations are in the left column.
722 Comparing functional (PET) images - the assessment .. (1991)
456 The relationship between global and local changes .. (1990)
436 A 3-dimensional statistical analysis .. (1992)
278 Localization in PET images: direct fitting .. (1989)
205 Plastic transformation of PET images (1991)
169 Analysis of fMRI time-series revisited (1995)
165 Functional connectivity - the PCA of large .. (1993)
158 Analysis of fMRI time-series revisited - again (1995)
152 Spatial registration and normalisation of images (1995)
143 Movement-related effects in fMRI time-series (1996)
102 A unified statistical approach for determining .. (1996)
98 Cognitive conjunction: a new approach to brain .. (1997)
97 Functional neuroanatomy of the human brain .. (1994)
76 Event-related fMRI (1997)
68 Local maxima and the expected Euler .. (1994)
63 The trouble with cognitive subtraction (1996)
57 Nonlinear event-related responses in fMRI (1998)
56 Event-related fMRI: characterising differential responses (1998)
52 Statistical parametric mapping: ontology and current issues (1995)
51 Characterizing evoked hemodynamics with fMRI (1995)
51 Modulation of connectivity in visual pathways .. (1997)
50 A voxel-based method for the statistical ... (1995)
Statistical Parametric Mapping: An Annotated Bibliography
This document was generated using the LaTeX2HTML translator Version 96.1 (Feb 5, 1996) Copyright © 1993, 1994, 1995, 1996, Nikos Drakos, Computer Based Learning Unit, University of Leeds.
The command line arguments were:
latex2html -split 0 -show_section_numbers methods_refs.tex.
The translation was initiated by Will Penny on Tuesday May 22 16:29:28 BST 2001