Bayesian inverse analysis of neuromagnetic data using cortically constrained multiple dipoles

Abstract

A recently introduced Bayesian model for magnetoencephalographic (MEG) data consistently localized multiple simulated dipoles with the help of marginalization of spatiotemporal background noise covariance structure in the analysis [Jun et al., (2005): Neuroimage 28:84-98]. Here, we elaborated this model to include subject's individual brain surface reconstructions with cortical location and orientation constraints. To enable efficient Markov chain Monte Carlo sampling of the dipole locations, we adopted a parametrization of the source space surfaces with two continuous variables (i.e., spherical angle coordinates). Prior to analysis, we simplified the likelihood by exploiting only a small set of independent measurement combinations obtained by singular value decomposition of the gain matrix, which also makes the sampler significantly faster. We analyzed both realistically simulated and empirical MEG data recorded during simple auditory and visual stimulation. The results show that our model produces reasonable solutions and adequate data fits without much manual interaction. However, the rigid cortical constraints seemed to make the utilized scheme challenging as the sampler did not switch modes of the dipoles efficiently. This is problematic in the presence of evidently highly multimodal posterior distribution, and especially in the relative quantitative comparison of the different modes. To overcome the difficulties with the present model, we propose the use of loose orientation constraints and combined model of prelocalization utilizing the hierarchical minimum-norm estimate and multiple dipole sampling scheme.

Figure 1

The full spatiotemporal noise covariance structure is identical for both simulated and empirical data and it was estimated from the real raw data for both cases.

Figure 2

Top: Simulated source locations include small patches on the primary auditory and visual areas. The auditory patches contained two and visual patches ten source space points on the densest 16,000 grid. Dark shades of gray denote the sulci and light shades of gray the gyri. Bottom: The simulated source timecourses for each stimulus category were acquired by averaging few sensors of MEG evoked fields and scaled to a suitable value regardless of the source extent.

Figure 3

On the left side the grid used for inverse analysis was the same as the one used for forward calculating the measurements (4,000) and on the right side the grid size used for inverse analysis was coarser (3,000).

Figure 4

Illustration of the energy of the posterior distribution at each step of our sampler. The chains converge rather rapidly to a local mode, but they get stuck easily and may jump to a mode having lower energy only after several thousands of samples.

Figure 5

The posterior samples of {θ,ϕ}‐coordinates of five different one dipole source configuration runs (A) and for all the runs of the auditory empirical dataset analyzed with grid size 4,000 (B). In each subfigure, the same color pair trendlines represent the ith {θ,ϕ}‐coordinate pair (i.e., one dipole).

Figure 6

Left side: The dipole locations at the end of the chains for different stimulus types of SIMULATED data using the 4,000 grid size. The number of linearly independent measurement combinations for the SVD speed‐up strategy is also shown for different cases. Yellow color depicts areas where many of the chains had a dipole at the end and reddish color areas with fewer occurrences. Note that the scales vary and that this does not represent the posterior distribution of the samples, but only how the locations of different chains are distributed. Right side: Same as left side but for 8,000 grid size.

Figure 7

Left side: The dipole locations at the end of the chains for different stimulus types of EMPIRICAL data using the 4,000 grid size. The number of linearly independent measurement combinations for the SVD speed‐up strategy is also shown for different cases. Yellow color depicts areas where many of the chains had a dipole at the end and reddish color areas with fewer occurrences. Note that the scales vary and that this does not represent the posterior distribution of the samples, but only how the locations of different chains are distributed. Right side: Same as left side but for 8,000 grid size.

Figure 8

The distribution of the number of dipoles at the end of the chains for simulated and empirical audiovisual data.

Figure 9

Representative best solution mode (lowest posterior energy) for simulated auditory source using 8,000 grid. The original source on the left hemisphere (magenta patch) is mislocalized (red dipole) to a neighboring deeper gyrus and, therefore, the orientation is opposite and the solution amplitude larger as more current is needed to produce the same fields from a deeper location. The standard deviation of the current samples (dashed lines) do not reveal the great difficulty in the multimodality of the posterior as it only gives information of this particular mode. The data fit to all 198 measured gradiometer channels and 31 timepoints produced by the solution is good given that only n _SVD = 23 measurement combinations were used in the inverse analysis.

Figure 10

Representative best solution mode (lowest posterior energy) for empirical auditory source using 8,000 grids. This solution mode suggests that there is an earlier more superficial source (green dipole) and a later deeper source (red dipole) on the left hemisphere and one earlier source (blue dipole) on the right hemisphere. The data fit is adequate although there are some differences between the measured and forward calculated fields.

Figure 11

Representative best solution mode (lowest posterior energy) for empirical visual source using 8,000 grids. In this solution mode, there are two solution dipoles, one of which is on the left V1 (red dipole), and the other is located on the right V1 (blue dipole). We obtained several other modes with the visual category yielding results having two or three dipoles located either as in this mode or bit more superficially, so that on one hemisphere there are two dipoles, which do not land directly on V1 or V2. In such case, the two dipoles probably just compensate for the large‐scale and widespread activation pattern of primary visual areas. This suggests that as several different modes produce similar and good data fits, the underlying source is far more complicated than two or three dipolar sources.