Classical and Bayesian inference in neuroimaging: theory

Stochastic models of neuronal dynamics

Cortical activity is the product of interactions among neuronal populations. Macroscopic electrophysiological phenomena are generated by these interactions. In principle, the mechanisms of these interactions afford constraints on biologically plausible models of electrophysiological responses. In other words, the macroscopic features of cortical activity can be modelled in terms of the microscopic behaviour of neurons. An evoked response potential (ERP) is the mean electrical potential measured from an electrode on the scalp, in response to some event. The purpose of this paper is to outline a population density approach to modelling ERPs. We propose a biologically plausible model of neuronal activity that enables the estimation of physiologically meaningful parameters from electrophysiological data. The model encompasses four basic characteristics of neuronal activity and organization: (i) neurons are dynamic units, (ii) driven by stochastic forces, (iii) organized into populations with similar biophysical properties and response characteristics and (iv) multiple populations interact to form functional networks. This leads to a formulation of population dynamics in terms of the Fokker-Planck equation. The solution of this equation is the temporal evolution of a probability density over state-space, representing the distribution of an ensemble of trajectories. Each trajectory corresponds to the changing state of a neuron. Measurements can be modelled by taking expectations over this density, e.g. mean membrane potential, firing rate or energy consumption per neuron. The key motivation behind our approach is that ERPs represent an average response over many neurons. This means it is sufficient to model the probability density over neurons, because this implicitly models their average state. Although the dynamics of each neuron can be highly stochastic, the dynamics of the density is not. This means we can use Bayesian inference and estimation tools that have already been established for deterministic systems. The potential importance of modelling density dynamics (as opposed to more conventional neural mass models) is that they include interactions among the moments of neuronal states (e.g. the mean depolarization may depend on the variance of synaptic currents through nonlinear mechanisms).Here, we formulate a population model, based on biologically informed model-neurons with spike-rate adaptation and synaptic dynamics. Neuronal sub-populations are coupled to form an observation model, with the aim of estimating and making inferences about coupling among sub-populations using real data. We approximate the time-dependent solution of the system using a bi-orthogonal set and first-order perturbation expansion. For didactic purposes, the model is developed first in the context of deterministic input, and then extended to include stochastic effects. The approach is demonstrated using synthetic data, where model parameters are identified using a Bayesian estimation scheme we have described previously.

Spatio-temporal modeling of localized brain activity

Functional neuroimaging, including positron emission tomography (PET) and functional magnetic resonance imaging (fMRI), plays an important role in identifying specific brain regions associated with experimental stimuli or psychiatric disorders such as schizophrenia. PET and fMRI produce massive data sets that contain both temporal correlations from repeated scans and complex spatial correlations. Several methods exist for handling temporal correlations, some of which rely on transforming the response data to induce either a known or an independence covariance structure. Despite the presence of spatial correlations between the volume elements (voxels) comprising a brain scan, conventional methods perform voxel-by-voxel analyses of measured brain activity. We propose a two-stage spatio-temporal model for the estimation and testing of localized activity. Our second-stage model specifies a spatial auto-regression, capturing correlations within neural processing clusters defined by a data-driven cluster analysis. We use maximum likelihood methods to estimate parameters from our spatial autoregressive model. Our model protects against type-I errors, enables the detection of both localized and regional activations (including volume of interest effects), provides information on functional connectivity in the brain, and establishes a framework to produce spatially smoothed maps of distributed brain activity for each individual. We illustrate the application of our model using PET data from a study of working memory in individuals with schizophrenia.

Modelling event-related responses in the brain

The aim of this work was to investigate the mechanisms that shape evoked electroencephalographic (EEG) and magneto-encephalographic (MEG) responses. We used a neuronally plausible model to characterise the dependency of response components on the models parameters. This generative model was a neural mass model of hierarchically arranged areas using three kinds of inter-area connections (forward, backward and lateral). We investigated how responses, at each level of a cortical hierarchy, depended on the strength of connections or coupling. Our strategy was to systematically add connections and examine the responses of each successive architecture. We did this in the context of deterministic responses and then with stochastic spontaneous activity. Our aim was to show, in a simple way, how event-related dynamics depend on extrinsic connectivity. To emphasise the importance of nonlinear interactions, we tried to disambiguate the components of event-related potentials (ERPs) or event-related fields (ERFs) that can be explained by a linear superposition of trial-specific responses and those engendered nonlinearly (e.g., by phase-resetting). Our key conclusions were; (i) when forward connections, mediating bottom-up or extrinsic inputs, are sufficiently strong, nonlinear mechanisms cause a saturation of excitatory interneuron responses. This endows the system with an inherent stability that precludes nondissipative population dynamics. (ii) The duration of evoked transients increases with the hierarchical depth or level of processing. (iii) When backward connections are added, evoked transients become more protracted, exhibiting damped oscillations. These are formally identical to late or endogenous components seen empirically. This suggests that late components are mediated by reentrant dynamics within cortical hierarchies. (iv) Bilateral connections produce similar effects to backward connections but can also mediate zero-lag phase-locking among areas. (v) Finally, with spontaneous activity, ERPs/ERFs can arise from two distinct mechanisms: For low levels of (stimulus related and ongoing) activity, the systems response conforms to a quasi-linear superposition of separable responses to the fixed and stochastic inputs. This is consistent with classical assumptions that motivate trial averaging to suppress spontaneous activity and disclose the ERP/ERF. However, when activity is sufficiently high, there are nonlinear interactions between the fixed and stochastic inputs. This interaction is expressed as a phase-resetting and represents a qualitatively different explanation for the ERP/ERF.

Motion or activity: their role in intra- and inter-subject variation in fMRI

Functional MRI (fMRI) carries the potential for non-invasive measurements of brain activity. Typically, what are referred to as activation images are actually thresholded statistical parametric maps. These maps possess large inter-session variability. This is especially problematic when applying fMRI to pre-surgical planning because of a higher requirement for intra-subject precision. The purpose of this study was to investigate the impact of residual movement artefacts on intra-subject and inter-subject variability in the observed fMRI activation. Ten subjects were examined using three different word-generation tasks. Two of the subjects were examined 10 times on 10 different days using the same paradigms. We systematically investigated one approach of correcting for residual movement effects: the inclusion of regressors describing movement-related effects in the design matrix of a General Linear Model (GLM). The data were analysed with and without modeling the residual movement artefacts and the impact on inter-session variance was assessed using F-contrasts. Inclusion of motion parameters in the analysis significantly reduced both the intra-subject as well as the inter-subject-variance.

An empirical Bayesian solution to the source reconstruction problem in EEG

Distributed linear solutions of the EEG source localisation problem are used routinely. In contrast to discrete dipole equivalent models, distributed linear solutions do not assume a fixed number of active sources and rest on a discretised fully 3D representation of the electrical activity of the brain. The ensuing inverse problem is underdetermined and constraints or priors are required to ensure the uniqueness of the solution. In a Bayesian framework, the conditional expectation of the source distribution, given the data, is attained by carefully balancing the minimisation of the residuals induced by noise and the improbability of the estimates as determined by their priors. This balance is specified by hyperparameters that control the relative importance of fitting and conforming to various constraints. Here we formulate the conventional "Weighted Minimum Norm" (WMN) solution in terms of hierarchical linear models. An "Expectation-Maximisation" (EM) algorithm is used to obtain a "Restricted Maximum Likelihood" (ReML) estimate of the hyperparameters, before estimating the "Maximum a Posteriori" solution itself. This procedure can be considered a generalisation of previous work that encompasses multiple constraints. Our approach was compared with the "classic" WMN and Maximum Smoothness solutions, using a simplified 2D source model with synthetic noisy data. The ReML solution was assessed with four types of source location priors: no priors, accurate priors, inaccurate priors, and both accurate and inaccurate priors. The ReML approach proved useful as: (1) The regularisation (or influence of the a priori source covariance) increased as the noise level increased. (2) The localisation error (LE) was negligible when accurate location priors were used. (3) When accurate and inaccurate location priors were used simultaneously, the solution was not influenced by the inaccurate priors. The ReML solution was then applied to real somatosensory-evoked responses to illustrate the application in an empirical setting.

Models of Brain Function in Neuroimaging

Inferences about brain function, using neuroimaging data, rest on models of how the data were caused. These models can be quite diverse, ranging from conceptual models of functional anatomy to nonlinear mathematical models of hemodynamics. However, they all have to be internally consistent because they model the same thing. This consistency encompasses many levels of description and places constraints on the statistical models, adopted for data analysis, and the experimental designs they embody. The aim of this review is to introduce the key models used in imaging neuroscience and how they relate to each other. We start with anatomical models of functional brain architectures, which motivate some of the fundaments of neuroimaging. We then turn to basic statistical models (e.g., the general linear model) used for making classical and Bayesian inferences about where neuronal responses are expressed. By incorporating biophysical constraints, these basic models can be finessed and, in a dynamic setting, rendered causal. This allows us to infer how interactions among brain regions are mediated.

Unsupervised learning and mapping of active brain functional MRI signals based on hidden semi-Markov event sequence models

In this paper, a novel functional magnetic resonance imaging (fMRI) brain mapping method is presented within the statistical modeling framework of hidden semi-Markov event sequence models (HSMESMs). Neural activation detection is formulated at the voxel level in terms of time coupling between the sequence of hemodynamic response onsets (HROs) observed in the fMRI signal, and an HSMESM of the hidden sequence of task-induced neural activations. The sequence of HRO events is derived from a continuous wavelet transform (CWT) of the fMRI signal. The brain activation HSMESM is built from the timing information of the input stimulation protocol. The rich mathematical framework of HSMESMs makes these models an effective and versatile approach for fMRI data analysis. Solving for the HSMESM Evaluation and Learning problems enables the model to automatically detect neural activation embedded in a given set of fMRI signals, without requiring any template basis function or prior shape assumption for the fMRI response. Solving for the HSMESM Decoding problem allows to enrich brain mapping with activation lag mapping, activation mode visualizing, and hemodynamic response function analysis. Activation detection results obtained on synthetic and real epoch-related fMRI data demonstrate the superiority of the HSMESM mapping method with respect to a real application case of the statistical parametric mapping (SPM) approach. In addition, the HSMESM mapping method appears clearly insensitive to timing variations of the hemodynamic response, and exhibits low sensitivity to fluctuations of its shape.

Wavelet variance components in image space for spatiotemporal neuroimaging data

Neuroimaging studies place great emphasis on not only the estimation but also the standard error estimates of underlying parameters derived from a temporal model. This allows inferences to be made about the signal estimates and resulting conclusions to be drawn about the underlying data. It can often be advantageous to interrogate temporal models after spatial transformation of the data into the wavelet domain. Wavelet bases provide a multiresolution decomposition of the spatial data dimension and an ensuing reduction in spatial correlation. However, widespread acceptance of these wavelet techniques has been hampered by the limited ability to reconstruct both parametric and error estimates into the image domain after analysis of temporal models in the wavelet domain. This paper introduces a derivation and a fast implementation of a method for the calculation of the variance of the parametric images obtained from wavelet filters. The technique is proposed for a class of estimators that have been shown to be useful in neuroimaging studies. The techniques are demonstrated for both functional magnetic resonance imaging (fMRI) and positron emission tomography (PET) data sets.

Bayesian population modeling of effective connectivity

A hierarchical model based on the Multivariate Autoregessive (MAR) process is proposed to jointly model neurological time-series collected from multiple subjects, and to characterize the distribution of MAR coefficients across the population from which those subjects were drawn. Thus, inference about effective connectivity between brain regions may be generalized beyond those subjects studied. The posterior on population- and subject-level connectivity parameters are estimated in a Variational Bayesian (VB) framework, and structural model parameters are chosen by the corresponding evidence criteria. The significance of resulting connectivity statistics are evaluated by permutation-based approximations to the null distribution. The method is demonstrated on simulated data and on actual multi-subject neurological time-series.

Mixed-effects and fMRI studies

This note concerns mixed-effect (MFX) analyses in multisession functional magnetic resonance imaging (fMRI) studies. It clarifies the relationship between mixed-effect analyses and the two-stage "summary statistics" procedure (Holmes, A.P., Friston, K.J., 1998. Generalisability, random effects and population inference. NeuroImage 7, S754) that has been adopted widely for analyses of fMRI data at the group level. We describe a simple procedure, based on restricted maximum likelihood (ReML) estimates of covariance components, that enables full mixed-effects analyses in the context of statistical parametric mapping. Using this procedure, we compare the results of a full mixed-effects analysis with those obtained from the simpler two-stage procedure and comment on the situations when the two approaches may give different results.

Prefrontal-hippocampal dynamics involved in learning regularities across episodes

Using functional magnetic resonance imaging, the neural correlates of context-specific memories and invariant memories about regularities across episodes were investigated. Volunteers had to learn conjunctions between objects and positions. In an invariant learning condition, positions were held constant, enabling subjects to learn regularities across trials. By contrast, in a context-specific condition object-position conjunctions were trial unique. Performance increase in the invariant learning condition was paralleled by a learning-related increase of inferior frontal gyrus activation and ventral striatal activation and a decrease of hippocampus activation. Conversely, in the context-specific condition hippocampal activation was constant across trials. We argue that the learning-related hippocampal activation pattern might be due to reduced relational binding requirements once regularities are extracted. Furthermore, we propose that the learning-related prefrontal modulation reflects the requirement to extract and maintain regularities across trials and the adjustment of object-position conjunctions on the basis of the extracted knowledge. Finally, our data suggest that the ventral striatum encodes the increased predictability of spatial features as a function of learning. Taken together, these results indicate a transition of the relative roles of distinct brain regions during learning regularities across multiple episodes: regularity learning is characterized by a shift from a hippocampal to a prefrontal-striatal brain system.

Multigrid priors for a Bayesian approach to fMRI

We introduce multigrid priors to construct a Bayesian-inspired method to asses brain activity in functional magnetic resonance imaging (fMRI). A sequence of different scale grids is constructed over the image. Starting from the finest scale, coarse grain data variables are sequentially defined for each scale. Then we move back to finer scales, determining for each coarse scale a set of posterior probabilities. The posterior on a coarse scale is used as the prior for activity at the next finer scale. To test the method, we use a linear model with a given hemodynamic response function to construct the likelihood. We apply the method both to real and simulated data of a boxcar experiment. To measure the number of errors, we impose a decision to determine activity by setting a threshold on the posterior. Receiver operating characteristic (ROC) curves are used to study the dependence on threshold and on a few hyperparameters in the relation between specificity and sensitivity. We also study the deterioration of the results for real data, under information loss. This is done by decreasing the number of images in each period and also by decreasing the signal to noise ratio and compare the robustness to other methods.

Comparing methods of analyzing fMRI statistical parametric maps

Approaches for the analysis of statistical parametric maps (SPMs) can be crudely grouped into three main categories in which different philosophies are applied to delineate activated regions. These being type I error control thresholding, false discovery rate (FDR) control thresholding and posterior probability thresholding. To better understand the properties of these main approaches, we carried out a simulation study to compare the approaches as they would be used on real data sets. Using default settings, we find that posterior probability thresholding is the most powerful approach, and type I error control thresholding provides the lowest levels of type I error. False discovery rate control thresholding performs in between the other approaches for both these criteria, although for some parameter settings this approach can approximate the performance of posterior probability thresholding. Based on these results, we discuss the relative merits of the three approaches in an attempt to decide upon an optimal approach. We conclude that viewing the problem of delineating areas of activation as a classification problem provides a highly interpretable framework for comparing the methods. Within this framework, we highlight the role of the loss function, which explicitly penalizes the types of errors that may occur in a given analysis.

Statistical parametric mapping for event-related potentials: I. Generic considerations

In this paper, we frame the strategy and motivations behind developments in statistical parametric mapping (SPM) for the analysis of electroencephalogram (EEG) data. This work deals specifically with SPM procedures for the analysis of event-related potentials (ERP). We place these developments in the larger context of integrating electrophysiological and hemodynamic measurements of evoked brain responses through the fusion of EEG and fMRI data. In this paper, we consider some fundamental issues when selecting an appropriate statistical model that enables diverse questions to be asked of the data and at the same time retains maximum sensitivity. The three key issues addressed in this paper are as follows: (i) should multivariate or mass univariate analyses be adopted, (ii) should time be treated as an experimental factor or as a dimension of the measured response variable, and (iii) how to form appropriate explanatory variables in a hierarchical observation model. We review the relative merits of the different options and explain the rationale for our choices. In brief, we motivate a mass univariate approach in terms of sensitivity to region-specific responses. This involves modeling responses at each voxel or space bin separately. In contradistinction, we treat time as an experimental factor to enable inferences about temporally distributed responses that encompass multiple time bins. In a companion paper, we develop statistical models of ERPs in the time domain that follow from the heuristics established here and illustrate the approach using simulated and real data.

Evaluation of mixed effects in event-related fMRI studies: impact of first-level design and filtering

With the introduction of event-related designs in fMRI, it has become crucial to optimize design efficiency and temporal filtering to detect activations at the 1st level with high sensitivity. We investigate the relevance of these issues for fMRI population studies, that is, 2nd-level analysis, for a set of event-related fMRI (er-fMRI) designs with different 1st-level efficiencies, adopting three distinct 1st-level filtering strategies as implemented in SPM99, SPM2, and FSL3.0. By theory, experiments, and simulations using physiological fMRI noise, we show that both design and filtering impact the outcome of the statistical analysis, not only at the 1st but also at the 2nd level. There are several reasons behind this finding. First, sensitivity is affected by both design and filtering, since the scan-to-scan variance, that is the fixed effect, is not negligible with respect to the between-subject variance, that is the random effect, in er-fMRI population studies. The impact of the fixed effects error on the sensitivity of the mixed effects analysis can be mitigated by an optimal choice of er-fMRI design and filtering. Moreover, the accuracy of the 1st- and 2nd-level parameter estimates also depend on design and filtering; especially, we show that inaccuracies caused by the presence of residual noise autocorrelations can be constrained by designs that have hemodynamic responses with a Gaussian distribution. In conclusion, designs with both good efficiency and decorrelating properties, for example, such as the geometric or Latin square probability distributions, combined with the "whitening" filters of SPM2 and FSL3.0, give the best result, both for 1st- and 2nd-level analysis of er-fMRI studies.

Statistical parametric mapping for event-related potentials (II): a hierarchical temporal model

In this paper, we describe a temporal model for event-related potentials (ERP) in the context of statistical parametric mapping (SPM). In brief, we project channel data onto a two-dimensional scalp surface or into three-dimensional brain space using some appropriate inverse solution. We then treat the spatiotemporal data in a mass-univariate fashion. This implicitly factorises the model into spatial and temporal components. The key contribution of this paper is the use of observation models that afford an explicit distinction between observation error and variation in the expression of ERPs. This distinction is created by employing a two-level hierarchical model, in which the first level models the ERP effects within-subject and trial type, while the second models differences in ERP expression among trial types and subjects. By bringing the analysis of ERP data into a classical hierarchical (i.e., mixed effects) framework, many apparently disparate approaches (e.g., conventional P300 analyses and time-frequency analyses of stimulus-locked oscillations) can be reconciled within the same estimation and inference procedure. Inference proceeds in the normal way using t or F statistics to test for effects that are localised in peristimulus time or in some time-frequency window. The use of F statistics is an important generalisation of classical approaches, because it allows one to test for effects that lie in a multidimensional subspace (i.e., of unknown but constrained form). We describe the analysis procedures, the underlying theory and compare its performance to established techniques.

fMRI correlates of the episodic retrieval of emotional contexts

Functional neuroimaging studies reveal differences in neural correlates of the retrieval of emotional and nonemotional memories. In the present experiment, encoding of emotionally neutral pictures in association with positively, neutrally or negatively valenced background contexts led to differential modulation of neural activity elicited in a subsequent recognition memory test for these pictures. Recognition of stimuli previously studied in emotional compared to neutral contexts elicited enhanced activity in structures previously implicated in episodic memory, including the parahippocampal cortex, hippocampus and prefrontal cortex. In addition, there was engagement of structures linked more specifically to emotional processing, including the amygdala, orbitofrontal cortex and anterior cingulate cortex. These emotion-related effects displayed both valence-independent and valence-specific components. We discuss the findings in terms of current models of emotional memory retrieval.

Geometrical interpretation of fMRI-guided MEG/EEG inverse estimates

Magneto- and electroencephalography (MEG/EEG) and functional magnetic resonance imaging (fMRI) provide complementary information about the functional organization of the human brain. An important advantage of MEG/EEG is the millisecond time resolution in detecting electrical activity in the cerebral cortex. The interpretation of MEG/EEG signals, however, is limited by the difficulty of determining the spatial distribution of the neural activity. Functional MRI can help in the MEG/EEG source analysis by suggesting likely locations of activity. We present a geometric interpretation of fMRI-guided inverse solutions in which the MEG/EEG source estimate minimizes a distance to a subspace defined by the fMRI data. In this subspace regularization (SSR) approach, the fMRI bias does not assume preferred amplitudes for MEG/EEG sources, only locations. Characteristic dependence of the source estimates on the regularization parameters is illustrated with simulations. When the fMRI locations match the true MEG/EEG source locations, they serve to bias the underdetermined MEG/EEG inverse solution toward the fMRI loci. Importantly, when the fMRI loci do not match the true MEG/EEG loci, the solution is insensitive to those fMRI loci.

Bayesian model averaging in EEG/MEG imaging

In this paper, the Bayesian Theory is used to formulate the Inverse Problem (IP) of the EEG/MEG. This formulation offers a comparison framework for the wide range of inverse methods available and allows us to address the problem of model uncertainty that arises when dealing with different solutions for a single data. In this case, each model is defined by the set of assumptions of the inverse method used, as well as by the functional dependence between the data and the Primary Current Density (PCD) inside the brain. The key point is that the Bayesian Theory not only provides for posterior estimates of the parameters of interest (the PCD) for a given model, but also gives the possibility of finding posterior expected utilities unconditional on the models assumed. In the present work, this is achieved by considering a third level of inference that has been systematically omitted by previous Bayesian formulations of the IP. This level is known as Bayesian model averaging (BMA). The new approach is illustrated in the case of considering different anatomical constraints for solving the IP of the EEG in the frequency domain. This methodology allows us to address two of the main problems that affect linear inverse solutions (LIS): (a) the existence of ghost sources and (b) the tendency to underestimate deep activity. Both simulated and real experimental data are used to demonstrate the capabilities of the BMA approach, and some of the results are compared with the solutions obtained using the popular low-resolution electromagnetic tomography (LORETA) and its anatomically constraint version (cLORETA).

Multilevel linear modelling for FMRI group analysis using Bayesian inference

Functional magnetic resonance imaging studies often involve the acquisition of data from multiple sessions and/or multiple subjects. A hierarchical approach can be taken to modelling such data with a general linear model (GLM) at each level of the hierarchy introducing different random effects variance components. Inferring on these models is nontrivial with frequentist solutions being unavailable. A solution is to use a Bayesian framework. One important ingredient in this is the choice of prior on the variance components and top-level regression parameters. Due to the typically small numbers of sessions or subjects in neuroimaging, the choice of prior is critical. To alleviate this problem, we introduce to neuroimage modelling the approach of reference priors, which drives the choice of prior such that it is noninformative in an information-theoretic sense. We propose two inference techniques at the top level for multilevel hierarchies (a fast approach and a slower more accurate approach). We also demonstrate that we can infer on the top level of multilevel hierarchies by inferring on the levels of the hierarchy separately and passing summary statistics of a noncentral multivariate t distribution between them.

Fully Bayesian spatio-temporal modeling of FMRI data

We present a fully Bayesian approach to modeling in functional magnetic resonance imaging (FMRI), incorporating spatio-temporal noise modeling and haemodynamic response function (HRF) modeling. A fully Bayesian approach allows for the uncertainties in the noise and signal modeling to be incorporated together to provide full posterior distributions of the HRF parameters. The noise modeling is achieved via a nonseparable space-time vector autoregressive process. Previous FMRI noise models have either been purely temporal, separable or modeling deterministic trends. The specific form of the noise process is determined using model selection techniques. Notably, this results in the need for a spatially nonstationary and temporally stationary spatial component. Within the same full model, we also investigate the variation of the HRF in different areas of the activation, and for different experimental stimuli. We propose a novel HRF model made up of half-cosines, which allows distinct combinations of parameters to represent characteristics of interest. In addition, to adaptively avoid over-fitting we propose the use of automatic relevance determination priors to force certain parameters in the model to zero with high precision if there is no evidence to support them in the data. We apply the model to three datasets and observe matter-type dependence of the spatial and temporal noise, and a negative correlation between activation height and HRF time to main peak (although we suggest that this apparent correlation may be due to a number of different effects).

A neural mass model for MEG/EEG: coupling and neuronal dynamics

Although MEG/EEG signals are highly variable, systematic changes in distinct frequency bands are commonly encountered. These frequency-specific changes represent robust neural correlates of cognitive or perceptual processes (for example, alpha rhythms emerge on closing the eyes). However, their functional significance remains a matter of debate. Some of the mechanisms that generate these signals are known at the cellular level and rest on a balance of excitatory and inhibitory interactions within and between populations of neurons. The kinetics of the ensuing population dynamics determine the frequency of oscillations. In this work we extended the classical nonlinear lumped-parameter model of alpha rhythms, initially developed by Lopes da Silva and colleagues [Kybernetik 15 (1974) 27], to generate more complex dynamics. We show that the whole spectrum of MEG/EEG signals can be reproduced within the oscillatory regime of this model by simply changing the population kinetics. We used the model to examine the influence of coupling strength and propagation delay on the rhythms generated by coupled cortical areas. The main findings were that (1) coupling induces phase-locked activity, with a phase shift of 0 or pi when the coupling is bidirectional, and (2) both coupling and propagation delay are critical determinants of the MEG/EEG spectrum. In forthcoming articles, we will use this model to (1) estimate how neuronal interactions are expressed in MEG/EEG oscillations and establish the construct validity of various indices of nonlinear coupling, and (2) generate event-related transients to derive physiologically informed basis functions for statistical modelling of average evoked responses.

Advances in functional and structural MR image analysis and implementation as FSL

The techniques available for the interrogation and analysis of neuroimaging data have a large influence in determining the flexibility, sensitivity, and scope of neuroimaging experiments. The development of such methodologies has allowed investigators to address scientific questions that could not previously be answered and, as such, has become an important research area in its own right. In this paper, we present a review of the research carried out by the Analysis Group at the Oxford Centre for Functional MRI of the Brain (FMRIB). This research has focussed on the development of new methodologies for the analysis of both structural and functional magnetic resonance imaging data. The majority of the research laid out in this paper has been implemented as freely available software tools within FMRIB's Software Library (FSL).

Motor Affordance and its Role for Visual Working Memory: Evidence from fMRI studies

Abstract. We examined the role of motor affordances of objects for working memory retention processes. Three experiments are reported in which participants passively viewed pictures of real world objects or had to retain the objects in working memory for a comparison with an S2 stimulus. Brain activation was recorded by means of functional magnetic resonance imaging (fMRI). Retaining information about objects for which hand actions could easily be retrieved (manipulable objects) in working memory activated the hand region of the ventral premotor cortex (PMC) contralateral to the dominant hand. Conversely, nonmanipulable objects activated the left inferior frontal gyrus. This suggests that working memory for objects with motor affordance is based on motor programs associated with their use. An additional study revealed that motor program activation can be modulated by task demands: Holding manipulable objects in working memory for an upcoming motor comparison task was associated with left ventral PMC activation. However, retaining the same objects for a subsequent size comparison task led to activation in posterior brain regions. This suggests that the activation of hand motor programs are under top down control. By this they can flexibly be adapted to various task demands. It is argued that hand motor programs may serve a similar working memory function as speech motor programs for verbalizable working memory contents, and that the premotor system mediates the temporal integration of motor representations with other task-relevant representations in support of goal oriented behavior.

A heuristic for the degrees of freedom of statistics based on multiple variance parameters

In neuroimaging, data are often modeled using general linear models. Here, we focus on GLMs with error covariances which are modeled as a linear combination of multiple variance/covariance components. Each of these components is weighted by one variance parameter. In many analyses variance parameters are estimated using restricted maximum likelihood (ReML). Most classical approaches assume the error covariance matrix can be factorized into a single variance parameter and a nonspherical correlation matrix. In this context, the F test based on a single variance parameter, with a suitable correction to the degrees of freedom, is the standard inference tool. This correction can also be adapted to models with multiple variance parameters. However, this extension overlooks the uncertainty about the variance parameter estimates and P values tend to be underestimated. Here, we show how one can overcome this problem to render the F test more exact. This issue is important, because serial correlations in fMRI time series are generally modeled using multiple variance parameters. Another application is to hierarchical linear models, which are used for modeling multisubject data. To illustrate our approach, we apply it to some typical modeling scenarios in fMRI data analysis.

Assessing brain activity through spatial bayesian variable selection

Statistical parametric mapping (SPM), relying on the general linear model and classical hypothesis testing, is a benchmark tool for assessing human brain activity using data from fMRI experiments. Friston et al. discuss some limitations of this frequentist approach and point out promising Bayesian perspectives. In particular, a Bayesian formulation allows explicit modeling and estimation of activation probabilities. In this study, we directly address this issue and develop a new regression based approach using spatial Bayesian variable selection. Our method has several advantages. First, spatial correlation is directly modeled for activation probabilities and indirectly for activation amplitudes. As a consequence, there is no need for spatial adjustment in a postprocessing step. Second, anatomical prior information, such as the distribution of grey matter or expert knowledge, can be included as part of the model. Third, the method has superior edge-preservation properties as well as being fast to compute. When applied to data from a simple visual experiment, the results demonstrate improved sensitivity for detecting activated cortical areas and for better preserving details of activated structures.

Bridging Functional MR Images and Scientific Inference: Reproducibility Maps

Historically, reproducibility has been the sine qua non of experimental findings that are considered to be scientifically useful. Typically, findings from functional magnetic resonance imaging (fMRI) studies are assessed with statistical parametric maps (SPMs) using a p value threshold. However, a smaller p value does not imply that the observed result will be reproducible. In this study, we suggest interpreting SPMs in conjunction with reproducibility evidence. Reproducibility is defined as the extent to which the active status of a voxel remains the same across replicates conducted under the same conditions. We propose a methodology for assessing reproducibility in functional MR images without conducting separate experiments. Our procedures include the empirical Bayes method for estimating effects due to experimental stimuli, the threshold optimization procedure for assigning voxels to the active status, and the construction of reproducibility maps. In an empirical example, we implemented the proposed methodology to construct reproducibility maps based on data from the study by Ishai et al. (2000). The original experiments involved 12 human subjects and investigated brain regions most responsive to visual presentation of 3 categories of objects: faces, houses, and chairs. The brain regions identified included occipital, temporal, and fusiform gyri. Using our reproducibility analysis, we found that subjects in one of the experiments exercised at least 2 mechanisms in responding to visual objects when performing alternately matching and passive tasks. One gave activation maps closer to those reported in Ishai et al., and the other had related regions in the precuneus and posterior cingulate. The patterns of activated regions are reproducible for at least 4 out of 6 subjects involved in the experiment. Empirical application of the proposed methodology suggests that human brains exhibit different strategies to accomplish experimental tasks when responding to stimuli. It is important to correlate activations to subjects' behavior such as reaction time and response accuracy. Also, the latency between the stimulus presentation and the peak of the hemodynamic response function varies considerably among individual subjects according to types of stimuli and experimental tasks. These variations per se also deserve scientific inquiries. We conclude by discussing research directions relevant to reproducibility evidence in fMRI.

General multilevel linear modeling for group analysis in FMRI

This article discusses general modeling of multisubject and/or multisession FMRI data. In particular, we show that a two-level mixed-effects model (where parameters of interest at the group level are estimated from parameter and variance estimates from the single-session level) can be made equivalent to a single complete mixed-effects model (where parameters of interest at the group level are estimated directly from all of the original single sessions' time series data) if the (co-)variance at the second level is set equal to the sum of the (co-)variances in the single-level form, using the BLUE with known covariances. This result has significant implications for group studies in FMRI, since it shows that the group analysis requires only values of the parameter estimates and their (co-)variance from the first level, generalizing the well-established "summary statistics" approach in FMRI. The simple and generalized framework allows different prewhitening and different first-level regressors to be used for each subject. The framework incorporates multiple levels and cases such as repeated measures, paired or unpaired t tests and F tests at the group level; explicit examples of such models are given in the article. Using numerical simulations based on typical first-level covariance structures from real FMRI data we demonstrate that by taking into account lower-level covariances and heterogeneity a substantial increase in higher-level Z score is possible.

Bayesian second-level analysis of functional magnetic resonance images

We propose a new method for the second-level analysis of functional MRI data based on Bayesian statistics. Our method does not require a computationally costly Bayesian model on the first level of analysis. Rather, modeling for single subjects is realized by means of the commonly applied General Linear Model. On the basis of the resulting parameter estimates for single subjects we calculate posterior probability maps and maps of the effect size for effects of interest in groups of subjects. A comparison of this method with the conventional analysis based on t statistics shows that the new approach is more robust against outliers. Moreover, our method overcomes some of the severe problems of null hypothesis significance tests such as the need to correct for multiple comparisons and facilitates inferences which are hard to formulate in terms of classical inferences.

Posterior probability maps and SPMs

This technical note describes the construction of posterior probability maps that enable conditional or Bayesian inferences about regionally specific effects in neuroimaging. Posterior probability maps are images of the probability or confidence that an activation exceeds some specified threshold, given the data. Posterior probability maps (PPMs) represent a complementary alternative to statistical parametric maps (SPMs) that are used to make classical inferences. However, a key problem in Bayesian inference is the specification of appropriate priors. This problem can be finessed using empirical Bayes in which prior variances are estimated from the data, under some simple assumptions about their form. Empirical Bayes requires a hierarchical observation model, in which higher levels can be regarded as providing prior constraints on lower levels. In neuroimaging, observations of the same effect over voxels provide a natural, two-level hierarchy that enables an empirical Bayesian approach. In this note we present a brief motivation and the operational details of a simple empirical Bayesian method for computing posterior probability maps. We then compare Bayesian and classical inference through the equivalent PPMs and SPMs testing for the same effect in the same data.

Variational Bayesian inference for fMRI time series

We describe a Bayesian estimation and inference procedure for fMRI time series based on the use of General Linear Models with Autoregressive (AR) error processes. We make use of the Variational Bayesian (VB) framework which approximates the true posterior density with a factorised density. The fidelity of this approximation is verified via Gibbs sampling. The VB approach provides a natural extension to previous Bayesian analyses which have used Empirical Bayes. VB has the advantage of taking into account the variability of hyperparameter estimates with little additional computational effort. Further, VB allows for automatic selection of the order of the AR process. Results are shown on simulated data and on data from an event-related fMRI experiment.

Modeling regional and psychophysiologic interactions in fMRI: the importance of hemodynamic deconvolution

The analysis of functional magnetic resonance imaging (fMRI) time-series data can provide information not only about task-related activity, but also about the connectivity (functional or effective) among regions and the influences of behavioral or physiologic states on that connectivity. Similar analyses have been performed in other imaging modalities, such as positron emission tomography. However, fMRI is unique because the information about the underlying neuronal activity is filtered or convolved with a hemodynamic response function. Previous studies of regional connectivity in fMRI have overlooked this convolution and have assumed that the observed hemodynamic response approximates the neuronal response. In this article, this assumption is revisited using estimates of underlying neuronal activity. These estimates use a parametric empirical Bayes formulation for hemodynamic deconvolution.

Robust Bayesian estimation of the hemodynamic response function in event-related BOLD fMRI using basic physiological information

In BOLD fMRI data analysis, robust and accurate estimation of the Hemodynamic Response Function (HRF) is still under investigation. Parametric methods assume the shape of the HRF to be known and constant throughout the brain, whereas non-parametric methods mostly rely on artificially increasing the signal-to-noise ratio. We extend and develop a previously proposed method that makes use of basic yet relevant temporal information about the underlying physiological process of the brain BOLD response in order to infer the HRF in a Bayesian framework. A general hypothesis test is also proposed, allowing to take advantage of the knowledge gained regarding the HRF to perform activation detection. The performances of the method are then evaluated by simulation. Great improvement is shown compared to the Maximum-Likelihood estimate in terms of estimation error, variance, and bias. Robustness of the estimators with regard to the actual noise structure or level, as well as the stimulus sequence, is also proven. Lastly, fMRI data with an event-related paradigm are analyzed. As suspected, the regions selected from highly discriminating activation maps resulting from the method exhibit a certain inter-regional homogeneity in term of HRF shape, as well as noticeable inter-regional differences.

Detection of grey matter loss in mild Alzheimer's disease with voxel based morphometry

OBJECTIVES

To test the applicability of an automated method of magnetic resonance image analysis (voxel based morphometry) to detect presence and severity of regional grey matter density reduction-a proxy of atrophy-in Alzheimer's disease.

METHODS

Twenty nine probable Alzheimer's patients and 26 non-demented controls (mini-mental state examinations mean (SD) 21 (4) and 29 (1)) underwent high resolution 3D brain magnetic resonance imaging. Spatial normalisation to a stereotactic template, segmentation into grey matter, white matter, and cerebrospinal fluid, and smoothing of the grey matter were carried out based on statistical parametric mapping (SPM99) algorithms. Analyses were carried out: (a) contrasting all Alzheimer's patients with all controls (p<0.05 corrected for multiple comparisons); (b) contrasting the three Alzheimer's patients with mini-mental state of 26 and higher with all controls (p<0.0001 uncorrected); and (c) correlating grey matter density with mini-mental state score within the Alzheimer's group (p<0.0001 uncorrected).

RESULTS

When all Alzheimer's patients were compared with controls, the largest atrophic regions corresponded to the right and left hippocampal/amygdalar complex. All parts of the hippocampus (head, body, and tail) were affected. More localised atrophic regions were in the temporal and cingulate gyri, precuneus, insular cortex, caudate nucleus, and frontal cortex. When the mildest Alzheimer's patients were contrasted with controls, the hippocampal/amygdalar complex were again found significantly atrophic bilaterally. The mini-mental state score correlated with grey matter density reduction in the temporal and posterior cingulate gyri, and precuneus, mainly to the right.

CONCLUSIONS

Voxel based morphometry with statistical parametric mapping is sensitive to regional grey matter density reduction in mild Alzheimer's disease.

Classical and Bayesian inference in neuroimaging: applications

In Friston et al. ((2002) Neuroimage 16: 465-483) we introduced empirical Bayes as a potentially useful way to estimate and make inferences about effects in hierarchical models. In this paper we present a series of models that exemplify the diversity of problems that can be addressed within this framework. In hierarchical linear observation models, both classical and empirical Bayesian approaches can be framed in terms of covariance component estimation (e.g., variance partitioning). To illustrate the use of the expectation-maximization (EM) algorithm in covariance component estimation we focus first on two important problems in fMRI: nonsphericity induced by (i) serial or temporal correlations among errors and (ii) variance components caused by the hierarchical nature of multisubject studies. In hierarchical observation models, variance components at higher levels can be used as constraints on the parameter estimates of lower levels. This enables the use of parametric empirical Bayesian (PEB) estimators, as distinct from classical maximum likelihood (ML) estimates. We develop this distinction to address: (i) The difference between response estimates based on ML and the conditional means from a Bayesian approach and the implications for estimates of intersubject variability. (ii) The relationship between fixed- and random-effect analyses. (iii) The specificity and sensitivity of Bayesian inference and, finally, (iv) the relative importance of the number of scans and subjects. The forgoing is concerned with within- and between-subject variability in multisubject hierarchical fMRI studies. In the second half of this paper we turn to Bayesian inference at the first (within-voxel) level, using PET data to show how priors can be derived from the (between-voxel) distribution of activations over the brain. This application uses exactly the same ideas and formalism but, in this instance, the second level is provided by observations over voxels as opposed to subjects. The ensuing posterior probability maps (PPMs) have enhanced anatomical precision and greater face validity, in relation to underlying anatomy. Furthermore, in comparison to conventional SPMs they are not confounded by the multiple comparison problem that, in a classical context, dictates high thresholds and low sensitivity. We conclude with some general comments on Bayesian approaches to image analysis and on some unresolved issues.

Bayesian estimation of dynamical systems: an application to fMRI

This paper presents a method for estimating the conditional or posterior distribution of the parameters of deterministic dynamical systems. The procedure conforms to an EM implementation of a Gauss-Newton search for the maximum of the conditional or posterior density. The inclusion of priors in the estimation procedure ensures robust and rapid convergence and the resulting conditional densities enable Bayesian inference about the model parameters. The method is demonstrated using an input-state-output model of the hemodynamic coupling between experimentally designed causes or factors in fMRI studies and the ensuing BOLD response. This example represents a generalization of current fMRI analysis models that accommodates nonlinearities and in which the parameters have an explicit physical interpretation. Second, the approach extends classical inference, based on the likelihood of the data given a null hypothesis about the parameters, to more plausible inferences about the parameters of the model given the data. This inference provides for confidence intervals based on the conditional density.