Cours et Atelier sur Méthodes Mathématiques en Cartographie Cérébrale
Course and Workshop on Mathematical Methods in Brain Mapping

Cours: du 5-8 décembre 2000 Atelier: du 10-11 décembre 2000 Course: December 5-8, 2000 Workshop: December 10-11, 2000 Organizateur / Organizer: Keith Worsley (McGill University) Mis a jour / Updated: 1/11/00.

Tuesday December 5, 2000: Course

Wednesday December 6, 2000: Course

Thursday December 7, 2000: Course

Friday December 8, 2000: Course

Saturday December 9, 2000

20:00 WELCOMING PARTY: `Le Crocodile', the Université de Montréal student bar, 5414 Gatineau, corner Lacombe (tel 733-2125), just 200 metres from the Hotel Terrasse-Royale (ask at the hotel, or see map )

Sunday December 10, 2000: Workshop

Session I: fMRI; Chair:

12:00-14:00 Lunch

Session II: Fusion of fMRI and EEG; Chair:

17:15 Reception hosted by the Centre de Recherches Mathématiques

Monday, December 11, 2000: Workshop

Session III: Structure; Chair:

12:00-14:00 Lunch

Session IV: Cortical Surface Mapping; Chair:

18:00 Beer in Thompson House, McGill University, or "Le Crocodile"

Random fields and their geometry

1. THE BASICS OF GAUSSIAN RANDOM FIELDS: One of the advantages that Gaussian processes have over their Markov counterparts is that the ``general theory'' is almost independent of the structure of the parameter space. I shall show why this is the case, by first introducing this family of stochastic process, and then looking at the issue of sample path regularity from the viewpoint of ``metric entropy''. In doing so we shall be able to treat processes defined on the line, random fields defined in Euclidean space, and even set and function indexed processes all the same way.

2. MAXIMA OF GAUSSIAN FIELDS: From the point of view of applications, one of the most important issues in the study of any class of real valued stochastic processes is the precise determination of the exceedence probabilities P{sup_{t in T} X_t > u}, where T is a parameter set and u is some level. This lecture will concentrate on a number of ways of attacking this problem for Gaussian fields, as well as looking at the structure of Gaussian fields at high levels.

3. RANDOM GEOMETRY AND EULER CHARACTERISTICS: This lecture will centre around the ``excursion sets'' of Gaussian fields; viz. the sets {t in T: X_t > u} for some u. I shall discuss their geometric structure, and some of the ways this can be quantified. In the end, we shall see that probabilistic considerations make the ``Euler characteristic'' (which I shall describe in detail) the ``right way'' to quantify excursion sets, and I shall describe what can be done with this.

4. EXTENSIONS: This lecture will be less structured, and may change depending on what we manage to cover in the first three and where audience interest lies. Topics that could be included are: how to handle non-smooth (fractal) Gaussian fields, Gaussian fields on manifolds, non-Gaussian processes, the Gaussian-Markov ``isomorphism theorem''.

Some Thoughts on Probability Distributions for Shape

In order to estimate a deformation field, a smoothness constraint needs to be incorporated into the registration model. Some warping methods constrain the deformations by parameterizing them in terms of spatial basis functions, whereby the prior variability of the deformations is assumed to be infinite for the spatial frequencies modeled, and zero for those higher frequencies not captured by the basis functions. Other authors use a form of linear regularization, such as minimizing bending energy or elastic energy, which is generally fast to compute, but does not implicitly preserve topology. These deformations can be explicitly constrained to be one-to-one, but this results in probability distributions with artificially sharp edges. What is needed is a formulation where a warp becomes increasingly less likely as it approaches singularity.

Normally distributed variables are allowed to be negative, but if a deformation field is one to one, the relative lengths, areas and volumes encoded by the field all have to be positive. Similarly, if the length of an object is normally distributed, then it is unlikely that its volume will also be normally distributed. This motivates the use of log-normal distributions. If logs of lengths are normally distributed, then logs of areas and volumes are also normally distributed. This choice of log- normal distributions can also be justified in terms of symmetry. Given a pair of randomly chosen brain images, the size of a structure in one brain is just as likely to be n times the size of the equivalent structure in the other brain, as it is to be 1/n times the size. In other words, a spatial transformation should be considered just as likely as its inverse.

Another assumption is that the probability distributions should be independent of pose, so that all rotations and translations are equally likely. The Jacobian matrix (partial derivatives) at each point in a deformation field encodes relative sizes, but also the relative orientation of the images. By using Singular Value Decomposition, a Jacobian matrix can be decomposed into the product of a rotation matrix, a set of orthogonal zooms and a second rotation matrix. It is the orthogonal zooms that are assumed to be log-normally distributed.

We present an implementation of this Bayesian approach to deforming images using the problem of spatial normalization, an essential component of both functional neuroimaging analyses and computational neuroanatomy.

Nonlinear Regression for Dynamic PET Ligand Studies

The method estimates the variance of the parameters using nonlinear least squares theory. A number of underlying assumptions must be assessed before this can be applied, most notably no autocorrelation between the residuals and a good agreement between the Basis Function fitting technique (Gunn et al, 1997) and nonlinear least squares fitting.

Fixed effects analysis can be carried out using the parameters and their associated variances. A further extension, using the methods established for fMRI to estimate the variance ratio between fixed and random effects (Worsley et al, 2000), will be considered, to allow the method more general application.

Gunn RN, Lammertsma AA, Hume SP and Cunningham VJ. Parametric Imaging of Ligand-Receptor Binding in PET using a Simplified Reference Tissue Model. Neuroimage 6,279-287 (1997).

Worsley KJ, Liao C, Grabove M, Petre V, Ha B, Evans AC. A general statistical analysis for fMRI data. NeuroImage, 11:S648 (2000).

Ballistocardiogram Removal and Motion Correction for EEG in the Magnet

(This is joint work with Giorigio Bonmassar, Patrick Purdon, Iiro P. Jaaskelainen, Victor Solo and Jack Belliveau)

Diffusion Smoothing on the Cortical Surface via the Laplace-Beltrami operator.

We present two explicit methods of the diffusion smoothing on a triangulated brain surface using the Laplace-Beltrami operator. The first method called the 'parametric method' uses quadratic polynomials for local surface parameterization. Then using a conformal coordinate transform, the Laplace-Beltrami operator is reduced to the planar Laplacian. The second technique is based on the 'finite element method' and so it has the advantage of avoiding the local surface parameterization.

As an illustration, we demonstrate the diffusion smoothing of the mean curvatures on the triangulated cortical surface consisting of 81920 triangles and show how the smoothing incorporates the geodesic curvature information.

References:

A. Andrade, F. Kherif, J. Mangin, K.J. Worsley, O. Simon, S. Dehaene, D. L. Bihan and J. Poline, "Detection of fMRI Activation on the Cortical Surface", NeuroImage 2000, in press.

Moo K. Chung, J. Taylor, K. J. Worsley, J. O. Ramsay, S. Robbins, A. Evans, "Diffusion Smoothing on the Cortical Surface via the Laplace-Beltrami Operator", in preparation (http://www.math.mcgill.ca/chung/diffusion/diffusion.pdf)

Spatiotemporal Imaging of Brain Activity: From the Microscopic to the Systems Level

Dynamic causal modelling of fMRI time-series

Wiener Deconvolution of fMRI Impulse Reponse

Despite the nonlinearity, Wiener deconvolution can be used to deblur the response, with effectiveness that depends on the noise characteristics. It is found that the impulse response function is subject-dependant. Therefore, the impulse response function is measured for each subject using a short-stimulus paradigm in order to generate the inverse deconvolution filter.

It is suggested that such deconvolution methods may be effective in diminishing the hemodynamically-imposed temporal blurring, and may have potential applications in quantitating responses in event-related fMRI.

Spatial mixture modeling of fMRI data

I will present a pragmatic alternative, namely locally specified models, which are analytically tractable. While global features of the activation pattern cannot be studied in this framework, local properties may be included in the analysis in a parametric form. A mixture distribution is assumed for a region of voxels, where the components correspond to the different activation states, and the posterior probability that a voxel is active is derived in closed form. The parameters of the model are estimated directly from the data. The approach is compared to existing methods using data from statistical image analysis, synthetic fMRI data and visual stimulation data.

Quasi-Conformal Maps of Brain Surfaces

It is impossible to flatten a surface with intrinsic curvature (such as the brain) without introducing linear and areal distortion, but the Riemann Mapping Theorem proves that it is possible to preserve angular (conformal) information under flattening. I will describe a novel computational method which uses the mathematical theory of circle packings to create quasi-conformal flat maps of the cortical surface. Cortical maps obtained from this approach will be presented. These maps exhibit conformal behavior in that angular distortion is controlled and can be created in the Euclidean and hyperbolic planes and on a sphere. Möbius transformations can be used to interactively change the map focus and no extraneous cuts in the original cortical surface are required. A canonical coordinate system can be imposed on these maps. In addition, these maps are mathematically unique.

Spatial modelling of functional imaging data using anatomically informed basis functions

The application of this spatial modelling is illustrated using as anatomical information the reconstructed grey matter surface derived from high-resolution T1-weighted magnetic resonance images (MRI). The spatial components specified in the model are of low spatial frequency, following the grey matter surface. When extracting these smooth components, by fitting them to the model, one efficiently captures spatially smooth components localized close or within the grey matter sheet. Effectively, the method implements a spatially variable anatomically informed smoothing with anatomically informed basis functions (AIBF). AIBF can be used for the analysis of any functional imaging modality. We have applied it to simulated and real functional MRI (fMRI) and positron emission tomography (PET) data. Amongst the various applications are high-resolution modelling of single-subject data (e.g. fMRI), spatial deconvolution (PET) and the analysis of multiple subject data by using canonical anatomical information.

Introduction to Brain Imaging

A structural alternative to the deformable brain atlas paradigm

D. Rivière, J.-F. Mangin, D. Papadopoulos, J.-M. Martinez, V. Frouin and J. Régis, Automatic recognition of cortical sulci using a congregation of neural networks , MICCAI'2000, Pittsburgh, LNCS-1935, Springer Verlag, pp 40-49.

An innovation approach to dynamic brain imaging

Cortical Surface Mapping

Functional and dynamic Magnet Resonance Imaging using vector adaptive weights smoothing

A generalization of the Adaptive Weights Smoothing (AWS) method allows to handle time series of images which typically occur in functional and dynamic MRI. It is shown how signal detection in functional MRI and classification in dynamic MRI can benefit from spatially adaptive smoothing. The performance of the procedure is illustrated by applications to real and simulated data.

Literature:

J. Polzehl and V. Spokoiny (2000). Adaptive weights smoothing with applications to image restoration. J. of the Royal Statist. Soc, Ser B, 62:335-354.

J. Polzehl and V. Spokoiny (1999) Vector adaptive weights smoothing with applications to MRI. WIAS Preprint 519.

A Natural RKHS formulation of the EEG/MEG forward and inverse problems

The EEG and MEG data was represented by multipole expansions in the spatial frequency domain. This was achieved using: a) expansions of the electric and magnetic lead fields in terms of sets of orthogonal eigenfunctions of the Laplace operator; and b) a basis for the reproducing kernel Hilbert space of brain electrical sources. The use of this representation has allowed: 1) to characterize the null spaces of the electric and magnetic lead fields; and 2) to obtain closed expressions for the combined EEG and MEG general smoothing spline inverse solutions for the real case of discrete measurements due to a finite set of sensors.

In this work, the feasibility of the combination of fMRI data with the EEG and MEG is also discussed on the basis of this mathematical contruct.

Complementary Methods in FMRI Analysis

Firstly, our GLM method is outlined. This uses local robust autocorrelation estimation (with non-linear spatial averaging to increase estimation robustness) to drive time series prewhitening on a per-voxel level, before the GLM is applied. This results in an increase in activation estimation efficiency of up to a factor of 2, depending on the experimental paradigm. This work is being extended to allow dynamic modelling of the haemodynamic response to the input function and replacement of spatial smoothing and (arbitrarily-thresholded) clustering with a full spatio-temporal noise model and MRF-based Bayesian activation segmentation.

Next, our model-free method is described, using ICA to identify different spatial and temporal components in the data. This approach can be used to identify non-Gaussian (e.g. scanner-related, physiological and motion-related) artefacts and activation components without knowledge of the experimental design or the resulting physiological response.

Finally, the complementarity between the approaches is explored. ICA can be used to remove artefacts (effectively an active pre-filtering step) before GLM analysis, reducing unmodelled variance (noise) in the data. ICA can also aid the design of an optimal model for GLM analysis, or even be integrated into a hybrid approach where ICA components are directly used in the design matrix. Some aspects of GLM-style significance testing can be used in ICA to allow significance to be given to ICA results; for example one can use a GLM-style model to obtain significances of IC time series (after projection of the model in the subspace in which ICA is run).

Random fields without stationarity on surfaces and volumes

A fundamental assumption of this approach was that the random field was isotropic, i.e. invariant under the group of Euclidean rigid motions. However, when studying random fields on manifolds such as the cortical surface, there is no notion of isotropy because of the irregularity of the manifold. Further, the distribution of the maximum of the random field is invariant under diffeomorphisms of the parameter space, so that our inference should not be based on how we visualize the data, i.e. flattening the cortical surface to view data on the cortical surface should not affect the choice of threshold.

These ideas point to the fact that the relevant geometry in the study of such random fields is one that is intrinsic to the random field itself. In other words, the relevant geometry is derived from the (Riemannian) structure that the random field induces on the manifold and not (necessarily) from the space in which we visualize the data. In the case of an isotropic field on a Euclidean space, the assumption of isotropy restricts the geometry to be the geometry of the Euclidean space (modulo a constant). We describe how the EC approach can be extended to random fields on manifolds and how to estimate relevant geometric quantities that appear in the EC approximation of the tail of the distribution of the maximum.

Mathematical Methods for EEG/MEG and their Fusion with other Neuorimages

I. INTRODUCTION TO EEG/MEG TOMOGRAPHY: The nature of electrophysiological data. Promise of EEG/MEG Tomography. Tasks to be accomplished: i) Specification of the volume conductor model ii) Specification of the Current Distribution model iii) Methods for Inverse solution. Origin of the EEG and MEG. The EEG/MEG forward problem. Concept of Lead Field. Null spaces for the EEG,MEG and EEG/MEG lead field operators in spherical geometry. General considerations on solving the EEG/MEG inverse problems.

II. EEG/MEG TOMOGRAPHIC METHODS: Deterministic criteria for measuring the performance of inverse soluitons: localization error, blurring and visibility. A general bayesian framwork for EEG/MEG tomography. Overview of major methods for EEG/MEG tomography within this framework. First order bayesian methods: diopole fits, linera distributed inverse sollutions. Nonlinear first order methods. Second order methods. Variable Resolution EEG/MEG inverse soluion. Spatio-temproal modeling.

III. STATISTICAL ISSUES IN EEG/MEG TOMOGRAPHY: Detectability as a further criteria for measuring the preformance of an inverse solution. EEG inverse solutions in the time domain, frequency domain and time-frequency domain. Pecularities of applying random field theory to EEG/MEG tomography. Use of randomization methods. Measures of connectivity and causality. Some examples of experimental and clinical applications.

IV. FUSION OF EEG/MEG TOMGRAPHY WITH OTHER IMAGING MODALITIES: a general Bayesian framework for neuroimge fusion. Example of the specific case of fusion of MEG data with fMRI data.