Using voxel-specific hemodynamic response function in EEG-fMRI data analysis: An estimation and detection model

Research groups who study epileptic spikes with simultaneous EEG-fMRI have used mostly the general linear model (GLM). A shortcoming of the GLM is that the specification of a simple hemodynamic response function (HRF) may lead to biased results. Other methods, which predict the hemodynamic response from the measured data, have been termed "recognition models". The merit of recognition models lies in the power of estimating the region-specific or voxel-specific HRF. We propose an approach that merges these two models in a general framework: estimate the HRF on the training data sets, and applying the estimated HRF on the other part of the data sets. The merit of this framework is that it can utilize the advantages of both models. A comparison of performance is made between the GLM with three fixed HRFs and the new model with voxel-specific HRFs. The main results are as follows: (1) in 18 of the 21 patients, the new model has a higher adjusted coefficient of multiple determination than the GLM with fixed HRF; (2) in some subjects, with the new model, we found areas of activation that had not been detected with the three fixed HRFs at our threshold of significance. The results suggest that the new model can do better than the fixed HRF GLM for the analysis of epileptic activity with EEG-fMRI.     

Using patient-specific hemodynamic response functions in combined EEG-fMRI studies in epilepsy

Combining electroencephalogram (EEG) and functional MRI (fMRI) allows localization of brain regions activated as a result of epileptic spikes. The statistical analysis of fMRI data usually includes a standard model of the hemodynamic response function (HRF) but it is not known how well this fits the actual HRF of epileptic spikes. The objective of this exploratory study was to compare the activated areas and t-statistical scores obtained with a standard HRF to those obtained with a patient-specific HRF. Eight patients with focal epilepsy were studied. We obtained an estimate of the patient-specific HRFs for each patient at the local maximum of activation in the standard HRF analysis. The activated areas obtained with the patient-specific HRFs were larger or similar to the originally activated areas. Additional activated areas were seen in five patients, and most were compatible with the EEG and anatomical MRI localization of epileptogenic and lesional regions. Using patient-specific HRFs brings increased sensitivity to the analysis of epileptic spikes by EEG-fMRI.     

Using voxel-specific hemodynamic response function in EEG-fMRI data analysis

Most existing analytical techniques for EEG-fMRI data need specific assumptions about the hemodynamic response function (HRF). These assumptions may not be appropriate when the HRF varies from subject to subject or from region to region. In this article, we introduce a deconvolution method for EEG-fMRI activation detection, which can be implemented with voxel-specific HRFs. A comparison of performance is made between three fixed HRFs and the deconvolution method under the framework of the general linear model. The main results are as follows: (1) the volume of detected regions from the deconvolved HRFs is larger. (2) In some subjects, the deconvolution technique can find areas of activation that have not been detected with the three fixed HRFs at our threshold of significance. (3) Deconvolution obtained higher adjusted coefficients of multiple determination compared to those obtained with the three fixed HRFs. The results suggest that the fixed HRF methods may not be the most appropriate for the analysis of epileptic activity with EEG-fMRI, and the deconvolution method may be a better choice.     

Constrained linear basis sets for HRF modelling using Variational Bayes

FMRI modelling requires flexible haemodynamic response function (HRF) modelling, with the HRF being allowed to vary spatially and between subjects. To achieve this flexibility, voxelwise parameterised HRFs have been proposed; however, inference on such models is very slow. An alternative approach is to use basis functions allowing inference to proceed in the more manageable General Linear Model (GLM) framework. However, a large amount of the subspace spanned by the basis functions produces nonsensical HRF shapes. In this work we propose a technique for choosing a basis set, and then the means to constrain the subspace spanned by the basis set to only include sensible HRF shapes. Penny et al. showed how Variational Bayes can be used to infer on the GLM for FMRI. Here we extend the work of Penny et al. to give inference on the GLM with constrained HRF basis functions and with spatial Markov Random Fields on the autoregressive noise parameters. Constraining the subspace spanned by the basis set allows for far superior separation of activating voxels from nonactivating voxels in FMRI data. We use spatial mixture modelling to produce final probabilities of activation and demonstrate increased sensitivity on an FMRI dataset.     

The BOLD response to interictal epileptiform discharges

We studied single-event and average BOLD responses to EEG interictal epileptic discharges (IEDs) in four patients with focal epilepsy, using continuous EEG-fMRI during 80-min sessions. The detection of activated areas was performed by comparing the BOLD signal at each voxel to a model of the expected signal. Since little is known about the BOLD response to IEDs, we modeled it with the response to brief auditory events (G. H., NeuroImage 9, 416-429). For each activated area, we then obtained the time course of the BOLD signal for the complete session and computed the actual average hemodynamic response function (HRF) to IEDs. In two of four patients, we observed clear BOLD responses to single IEDs. The average response was composed of a positive lobe peaking between 6 and 7 s in all patients and a negative undershoot in three patients. There were important variations in amplitude and shape between average HRFs across patients. The average HRF presented a wider positive lobe than the Glover model in three patients and a longer undershoot in two. There was a remarkable similarity in the shape of the HRF across areas for patients presenting multiple activation sites. There was no clear correlation between the amplitude of individual BOLD responses and the amplitude of the corresponding EEG spike. The possibility of a longer HRF could be used to improve statistical detection of activation in simultaneous EEG-fMRI. The variability in average HRFs across patients could reflect in part different pathophysiological mechanisms.     

Variability of the hemodynamic response as a function of age and frequency of epileptic discharge in children with epilepsy

EEG-fMRI is a non-invasive tool to investigate epileptogenic networks in patients with epilepsy. Different patterns of BOLD responses have been observed in children as compared to adults. A high intra- and intersubject variability of the hemodynamic response function (HRF) to epileptic discharges has been observed in adults. The actual HRF to epileptic discharges in children and its dependence on age are unknown. We analyzed 64 EEG-fMRI event types in 37 children (3 months to 18 years), 92% showing a significant BOLD response. HRFs were calculated for each BOLD cluster using a Fourier basis set. After excluding HRFs with a low signal-to-noise ratio, 126 positive and 98 negative HRFs were analyzed. We evaluated age-dependent changes as well as the effect of increasing numbers of spikes. Peak time, amplitude and signal-to-noise ratio of the HRF and the t-statistic score of the cluster were used as dependent variables. We observed significantly longer peak times of the HRF in the youngest children (0 to 2 years), suggesting that the use of multiple HRFs might be important in this group. A different coupling between neuronal activity and metabolism or blood flow in young children may cause this phenomenon. Even if the t-value increased with frequent spikes, the amplitude of the HRF decreased significantly with spike frequency. This reflects a violation of the assumptions of the General Linear Model and therefore the use of alternative analysis techniques may be more appropriate with high spiking rates, a common situation in children.     

Independent component analysis reveals dynamic ictal BOLD responses in EEG-fMRI data from focal epilepsy patients

IntroductionSeizures occur rarely during EEG-fMRI acquisitions of epilepsy patients, but can potentially offer a better estimation of the epileptogenic zone than interictal activity. Independent component analysis (ICA) is a data-driven method that imposes minimal constraints on the hemodynamic response function (HRF). In particular, the investigation of HRFs with clear peaks, but varying latency, may be used to differentiate the ictal focus from propagated activity.MethodsICA was applied on ictal EEG-fMRI data from 15 patients. Components related to seizures were identified by fitting an HRF to the component time courses at the time of the ictal EEG events. HRFs with a clear peak were used to derive maps of significant BOLD responses and their associated peak delay. The results were then compared with those obtained from a general linear model (GLM) method. Concordance with the presumed epileptogenic focus was also assessed.ResultsThe ICA maps were significantly correlated with the GLM maps for each patient (Spearman's test, p < 0.05). The ictal BOLD responses identified by ICA always included the presumed epileptogenic zone, but were also more widespread, accounting for 20.3% of the brain volume on average. The method provided a classification of the components as a function of peak delay. BOLD response clusters associated with early HRF peaks were concordant with the suspected epileptogenic focus, while subsequent HRF peaks may correspond to ictal propagation.ConclusionICA applied to EEG-fMRI can detect areas of significant BOLD response to ictal events without having to predefine an HRF. By estimating the HRF peak time in each identified region, the method could also potentially provide a dynamic analysis of ictal BOLD responses, distinguishing onset from propagated activity.     

The Variability of Human, BOLD Hemodynamic Responses

Cerebral hemodynamic responses to brief periods of neural activity are delayed and dispersed in time. The specific shape of these responses is of some importance to the design and analysis of blood oxygenation level-dependent (BOLD), functional magnetic resonance imaging (fMRI) experiments. Using fMRI scanning, we examine here the characteristics and variability of hemodynamic responses from the central sulcus in human subjects during an event-related, simple reaction time task. Specifically, we determine the contribution of subject, day, and scanning session (within a day) to variability in the shape of evoked hemodynamic response. We find that while there is significant and substantial variability in the shape of responses collected across subjects, responses collected during multiple scans within a single subject are less variable. The results are discussed in terms of the impact of response variability upon sensitivity and specificity of analyses of event-related fMRI designs.     

Human recombinant granulocyte-macrophage colony-stimulating factor and interleukin 3 cause basophil histamine release.

Histamine-releasing factors (HRFs) have been shown to be released from a variety of human cells, including T lymphocytes and alveolar macrophages. We considered the possibility that known cytokines might possess such activity on human basophils and/or mast cells and therefore tested preparations of human recombinant IL 3, IL 4, IL 5, granulocyte colony-stimulating factor (G-CSF) and granulocyte-macrophage colony-stimulating factor (GM-CSF) upon a panel of basophil donors. IL 3 and GM-CSF possessed significant histamine-releasing activity in 8 of 10 and 12 of 14 subjects, respectively. In each instance, a dose response could be demonstrated. IL 4 and G-CSF had no such activity, whereas IL 5 had activity in only 2 of 14 donors tested. We conclude that IL 3 and GM-CSF represent two effective HRFs, and suggest that HRF, as isolated based upon histamine-releasing activity, is likely to be heterogeneous in terms of molecular identity. Whether previously described HRFs relate specifically to IL 3 or GM-CSF must await primary sequence analysis of HRF and/or studies with monospecific antisera.     

Variability of the BOLD response over time: an examination of within-session differences

Model-based analysis methods for fMRI data assume a priori knowledge of the time course of the hemodynamic response (HR) in reaction to experimental stimuli or events. This knowledge is incorporated into the hemodynamic response function (HRF), which is a common model of the HR. Although it is already known that the HR varies across individuals and brain regions, few studies have investigated how variations within one session affect the results of statistical analysis using the general linear model (GLM). In this study, we formally tested for a possible variation of the BOLD response during prolonged functional measurement (120 min). To provoke performance of simple visual, motor, and cognitive tasks, we opted for a combination of a variant of the Stroop task and rotating L's. In selected regions of interest, time courses were extracted and compared with regard to mean and maximum amplitudes throughout the time of functional measurement. Additionally, parameter estimates derived from the GLM were tested for differences over time. Although differences between conditions were found to be significant, results did not show significant variance due to a within-factor time. Similarly, a temporal change in the relation between conditions, in terms of an interaction between the within-factor time and the within-factor condition, was not detectable by a repeated measures ANOVA. Similar results were obtained for analysis of mean and maximum amplitudes as well as for the analyses of parameter estimates.     

Model-free fMRI group analysis using FENICA

Exploratory analysis of functional MRI data allows activation to be detected even if the time course differs from that which is expected. Independent Component Analysis (ICA) has emerged as a powerful approach, but current extensions to the analysis of group studies suffer from a number of drawbacks: they can be computationally demanding, results are dominated by technical and motion artefacts, and some methods require that time courses be the same for all subjects or that templates be defined to identify common components. We have developed a group ICA (gICA) method which is based on single-subject ICA decompositions and the assumption that the spatial distribution of signal changes in components which reflect activation is similar between subjects. This approach, which we have called Fully Exploratory Network Independent Component Analysis (FENICA), identifies group activation in two stages. ICA is performed on the single-subject level, then consistent components are identified via spatial correlation. Group activation maps are generated in a second-level GLM analysis. FENICA is applied to data from three studies employing a wide range of stimulus and presentation designs. These are an event-related motor task, a block-design cognition task and an event-related chemosensory experiment. In all cases, the group maps identified by FENICA as being the most consistent over subjects correspond to task activation. There is good agreement between FENICA results and regions identified in prior GLM-based studies. In the chemosensory task, additional regions are identified by FENICA and temporal concatenation ICA that we show is related to the stimulus, but exhibit a delayed response. FENICA is a fully exploratory method that allows activation to be identified without assumptions about temporal evolution, and isolates activation from other sources of signal fluctuation in fMRI. It has the advantage over other gICA methods that it is computationally undemanding, spotlights components relating to activation rather than artefacts, allows the use of familiar statistical thresholding through deployment of a higher level GLM analysis and can be applied to studies where the paradigm is different for all subjects.     

An empirical investigation into the number of subjects required for an event-related fMRI study

Optimising the number of subjects required for an event-related functional imaging study is critical for ensuring sufficient statistical power. We report an empirical investigation of this issue by employing a resampling approach to the data of 58 subjects drawn from four previous GO/NOGO studies. Using voxelwise measures and setting the activation map from the complete sample to be a "gold standard", analyses revealed the statistical power to be surprisingly low at typical sample sizes (n = 20). However, voxels that were significantly active from smaller samples tended to be true positives, that is, they were typically active in the gold standard map and correlated well with the gold standard activation measure. The numerous false negatives that resulted from the lower SNR of the smaller samples drove the poor statistical power of those samples. Splitting the sample into two groups provided a test of the reproducibility of activation maps that was assessed using an alternative measure that quantified the distances between centres-of-mass of activated areas. These analyses revealed that although the voxelwise overlap may be poor, the locations of activated areas provide some optimism for studies with typical sample sizes. With n = 20 in each of two groups, it was found that the centres-of-mass for 80% of activated areas fell within 25 mm of each other. The reported analyses, by quantifying the spatial reproducibility for various sample sizes performing a typical event-related cognitive task, thus provide an empirical measure of the disparity to be expected in comparing activation maps.     

Are fMRI event-related response constant in time? A model selection answer

An accurate estimation of the hemodynamic response function (HRF) in functional magnetic resonance imaging (fMRI) is crucial for a precise spatial and temporal estimate of the underlying neuronal processes. Recent works have proposed non-parametric estimation of the HRF under the hypotheses of linearity and stationarity in time. Biological literature suggests, however, that response magnitude may vary with attention or ongoing activity. We therefore test a more flexible model that allows for the variation of the magnitude of the HRF with time in a maximum likelihood framework. Under this model, the magnitude of the HRF evoked by a single event may vary across occurrences of the same type of event. This model is tested against a simpler model with a fixed magnitude using information theory. We develop a standard EM algorithm to identify the event magnitudes and the HRF. We test this hypothesis on a series of 32 regions (4 ROIS on eight subjects) of interest and find that the more flexible model is better than the usual model in most cases. The important implications for the analysis of fMRI time series for event-related neuroimaging experiments are discussed.     

Late Onset of Anterior Prefrontal Activity during True and False Recognition: An Event-Related fMRI Study

Previous studies using PET and fMRI to examine memory retrieval have been limited by the requirement to test different types of items in separate blocks and to average data across items and response types within blocks. We used recently developed procedures for analyzing event-related mixed trial data from fMRI experiments to compare brain activity during true recognition of previously studied words and false recognition of semantic associates. A previous PET study using blocked testing procedures reported similarities and differences in rCBF patterns associated with true and false recognition (Schacteret al.,1996a). We examined brain activity during blocked testing of studied words and nonstudied semantic associates (similar to PET), and also during event-related mixed trials, where studied words and nonstudied semantic associates are intermixed. Six subjects initially heard lists of semantically related words and were later tested for old/new recognition with studied words and nonstudied semantic associates, either in separate blocks or intermixed randomly for the event-related analysis. Compared to a fixation control condition, a variety of regions previously reported in the PET study showed significant activation for both true and false recognition, including anterior prefrontal, frontal opercular, medial parietal, and visual cortex extending into hippocampal/parahippocampal regions. Differences across trial types were not clearly present. Event-related analyses of time course data show a relatively late onset and sustained duration for anterior prefrontal signal changes compared to signal changes in other activated regions. Further study is needed to resolve whether this late onset originates from variance in hemodynamic response properties or is attributable to delayed neural activity. The delayed onset is consistent with the idea that anterior prefrontal regions participate in postretrieval monitoring processes.     

BOLD changes occur prior to epileptic spikes seen on scalp EEG

This study examined BOLD changes prior to interictal discharges in the EEG of patients with epilepsy. From a database of 143 EEG-fMRI studies, we selected the 16 data sets that showed both strong fMRI activation in the original analysis and only a single spike type in the EEG. Scans were then analyzed using seven model HRFs, peaking 3 or 1 s before the event, or 1, 3, 5, 7, or 9 s after it. An HRF was calculated using a deconvolution method for all activations seen in each analysis. The results showed that seven data sets had HRFs that peaked 1 s after the event or earlier, indicating a BOLD change starting prior to the spike seen on the scalp EEG. This is surprising since the BOLD change is expected to result from the spike. For most of the data sets with early peaking HRFs, the maximum activation in all of the statistical maps was when the model HRF peaked 1 s after the event, suggesting that the early activation was at least as important as any later activation. We suggest that this early activity is the result of neuronal changes occurring several seconds prior to a surface EEG event, but that these changes are not visible on the scalp. This is the first report of a BOLD response occurring several seconds prior to an interictal event seen on the scalp and could have important implications for our understanding of the generation of epileptic discharges.     

Multiple imputation of missing fMRI data in whole brain analysis

Whole brain fMRI analyses rarely include the entire brain because of missing data that result from data acquisition limits and susceptibility artifact, in particular. This missing data problem is typically addressed by omitting voxels from analysis, which may exclude brain regions that are of theoretical interest and increase the potential for Type II error at cortical boundaries or Type I error when spatial thresholds are used to establish significance. Imputation could significantly expand statistical map coverage, increase power, and enhance interpretations of fMRI results. We examined multiple imputation for group level analyses of missing fMRI data using methods that leverage the spatial information in fMRI datasets for both real and simulated data. Available case analysis, neighbor replacement, and regression based imputation approaches were compared in a general linear model framework to determine the extent to which these methods quantitatively (effect size) and qualitatively (spatial coverage) increased the sensitivity of group analyses. In both real and simulated data analysis, multiple imputation provided 1) variance that was most similar to estimates for voxels with no missing data, 2) fewer false positive errors in comparison to mean replacement, and 3) fewer false negative errors in comparison to available case analysis. Compared to the standard analysis approach of omitting voxels with missing data, imputation methods increased brain coverage in this study by 35% (from 33,323 to 45,071 voxels). In addition, multiple imputation increased the size of significant clusters by 58% and number of significant clusters across statistical thresholds, compared to the standard voxel omission approach. While neighbor replacement produced similar results, we recommend multiple imputation because it uses an informed sampling distribution to deal with missing data across subjects that can include neighbor values and other predictors. Multiple imputation is anticipated to be particularly useful for 1) large fMRI data sets with inconsistent missing voxels across subjects and 2) addressing the problem of increased artifact at ultra-high field, which significantly limit the extent of whole brain coverage and interpretations of results.     

Detection of time-varying signals in event-related fMRI designs

In neuroimaging research on attention, cognitive control, decision-making, and other areas where response time (RT) is a critical variable, the temporal variability associated with the decision is often assumed to be inconsequential to the hemodynamic response (HDR) in rapid event-related designs. On this basis, the majority of published studies model brain activity lasting less than 4 s with brief impulses representing the onset of neural or cognitive events, which are then convolved with the hemodynamic impulse response function (HRF). However, electrophysiological studies have shown that decision-related neuronal activity is not instantaneous, but in fact, often lasts until the motor response. It is therefore possible that small differences in neural processing durations, similar to human RTs, will produce noticeable changes in the HDR, and therefore in the results of regression analyses. In this study we compare the effectiveness of traditional models that assume no temporal variance with a model that explicitly accounts for the duration of very brief epochs of neural activity. Using both simulations and fMRI data, we show that brief differences in duration are detectable, making it possible to dissociate the effects of stimulus intensity from stimulus duration, and that optimizing the model for the type of activity being detected improves the statistical power, consistency, and interpretability of results.     

Enhancing reproducibility of fMRI statistical maps using generalized canonical correlation analysis in NPAIRS framework

Common fMRI data processing techniques usually minimize a temporal cost function or fit a temporal model to extract an activity map. Here, we focus on extracting a highly, spatially reproducible statistical parametric map (SPM) from fMRI data using a cost function that does not depend on a model of the subjects' temporal response. Based on a generalized version of canonical correlation analysis (gCCA), we propose a method to extract a highly reproducible map by maximizing the sum of pair-wise correlations between some maps. In a group analysis, each map is calculated from a linear combination of fMRI scans of a subset of subjects under study. The proposed method is applied to BOLD fMRI datasets without any spatial smoothing from 10 subjects performing a simple reaction time (RT) task. Using the NPAIRS split-half resampling framework with a reproducibility measure based on SPM correlations, we compare the proposed approach with canonical variate analysis (CVA) and a simple general linear model (GLM). gCCA provides statistical parametric maps with higher reproducibility than CVA and GLM with correlation reproducibilities across independent split-half SPMs of 0.78, 0.46, and 0.41, respectively. Our results show that gCCA is an efficient approach for extracting the default mode network, assessing brain connectivity, and processing event-related and resting-state datasets in which the temporal BOLD signal varies from subject to subject.     

The utility of near-infrared spectroscopy in the regression of low-frequency physiological noise from functional magnetic resonance imaging data

Near-infrared spectroscopy (NIRS) signals have been shown to correlate with resting-state BOLD-fMRI data across the whole brain volume, particularly at frequencies below 0.1 Hz. While the physiological origins of this correlation remain unclear, its existence may have a practical application in minimizing the background physiological noise present in BOLD-fMRI recordings. We performed simultaneous, resting-state fMRI and 28-channel NIRS in seven adult subjects in order to assess the utility of NIRS signals in the regression of physiological noise from fMRI data. We calculated the variance of the residual error in a general linear model of the baseline fMRI signal, and the reduction of this variance achieved by including NIRS signals in the model. In addition, we introduced a sequence of simulated hemodynamic response functions (HRFs) into the resting-state fMRI data of each subject in order to quantify the effectiveness of NIRS signals in optimizing the recovery of that HRF. For comparison, these calculations were also performed using a pulse and respiration RETROICOR model. Our results show that the use of 10 or more NIRS channels can reduce variance in the residual error by as much as 36% on average across the whole cortex. However the same number of low-pass filtered white noise regressors is shown to produce a reduction of 19%. The RETROICOR model obtained a variance reduction of 6.4%. Our HRF simulation showed that the mean-squared error (MSE) between the recovered and true HRFs is reduced by 21% on average when 10 NIRS channels are applied and by introducing an optimized time lag between the NIRS and fMRI time series, a single NIRS channel can provide an average MSE reduction of 14%. The RETROICOR model did not provide a significant change in MSE. By each of the metrics calculated, NIRS recording is shown to be of significant benefit to the regression of low-frequency physiological noise from fMRI data.     

Artifactual fMRI group and condition differences driven by performance confounds

Analysis techniques comparing groups or conditions that vary in performance are open to a possible confound driven by those performance differences, if these errors are ignored. Disproportionate numbers of errors may either introduce noise into the signal of interest or confound the signal of interest with additional signal associated with specific error-related processes. Two inhibitory task datasets were reanalysed, one comparing young and elderly groups, the other comparing high and low conflict conditions within the same group of subjects. The data were analysed twice using event-related techniques, one treating correct and error responses separately, the other treating error responses as if they were correct. It was found that the activation maps differed considerably, with the inclusion of errors leading to many false positive and false negative activation clusters. Using performance as a covariate, analyses of covariance (ANCOVA) were used to try to correct these differences without success. Data simulations that varied the number of errors included in the analyses found that surprisingly few errors could significantly alter activation maps. Consequently, brain-imaging investigations that do not accommodate error contributions to functional signals are at risk of misinterpreting activation patterns.