Abstract
To bridge the gap between preclinical cellular models of disease and in vivo imaging of human cognitive network dynamics, there is a pressing need for informative biophysical models. Here we assess dynamic causal models (DCM) of cortical network responses, as generative models of magnetoencephalographic observations during an auditory oddball roving paradigm in healthy adults. This paradigm induces robust perturbations that permeate frontotemporal networks, including an evoked 'mismatch negativity' response and transiently induced oscillations. Here, we probe GABAergic influences in the networks using double-blind placebo-controlled randomized-crossover administration of the GABA reuptake inhibitor, tiagabine (oral, 10 mg) in healthy older adults. We demonstrate the facility of conductance-based neural mass mean-field models, incorporating local synaptic connectivity, to investigate laminar-specific and GABAergic mechanisms of the auditory response. The neuronal model accurately recapitulated the observed magnetoencephalographic data. Using parametric empirical Bayes for optimal model inversion across both drug sessions, we identify the effect of tiagabine on GABAergic modulation of deep pyramidal and interneuronal cell populations. We found a transition of the main GABAergic drug effects from auditory cortex in standard trials to prefrontal cortex in deviant trials. The successful integration of pharmaco- magnetoencephalography with dynamic causal models of frontotemporal networks provides a potential platform on which to evaluate the effects of disease and pharmacological interventions.
SIGNIFICANCE STATEMENT Understanding human brain function and developing new treatments require good models of brain function. We tested a detailed generative model of cortical microcircuits that accurately reproduced human magnetoencephalography, to quantify network dynamics and connectivity in frontotemporal cortex. This approach identified the effect of a test drug (GABA-reuptake inhibitor, tiagabine) on neuronal function (GABA-ergic dynamics), opening the way for psychopharmacological studies in health and disease with the mechanistic precision afforded by generative models of the brain.
Introduction
Biophysically informed models of cognition and cognitive disorders facilitate the effective translation of the mechanisms and treatments of disease. Recent progress toward detailed generative models that replicate neurophysiological correlates of cognition based on cellular and network dynamics, such as “dynamic causal models” (DCM), make predictions that approximate observations by functional magnetic resonance imaging or electro- and magneto-encephalography (MEG) (Moran et al., 2013). To be most useful, these models should incorporate laminar, cellular and synaptic functions (Bastos et al., 2012), and adhere to basic principles of cortical connectivity (Shipp, 2016), while also being sufficiently tractable and accurate to study cognition.
The DCM framework developed to meet these criteria, with applications in health and neurological disorders (Kiebel et al., 2008; Stephan et al., 2008; Boly et al., 2011; Marreiros et al., 2015). DCMs draw on empirical priors for synaptic time constants and conductances, together with a mean-field forward model. They are optimized to match the observed neurophysiological data. DCMs are supported by extensive data for face-validity (Stephan et al., 2008, 2015) and construct-validity (Razi et al., 2015), but they must also achieve predictive validity (Moran et al., 2014; Gilbert and Moran, 2016; Shaw et al., 2019).
We tested the ability of DCMs to identify the effect of a pharmacological intervention. The DCMs were designed to model human frontotemporal cortical networks during an auditory oddball paradigm, with characteristic MEG responses to standard and deviant tones (<300 ms). The differential response to these tones (the mismatch negativity, MMN) is abnormal in many neurological diseases (Boly et al., 2011; Näätänen et al., 2011; Hughes et al., 2013), reflecting a change in prediction errors in hierarchical frontotemporal networks (Garrido et al., 2009b; Phillips et al., 2015).
To examine laminar- and synaptic-dynamics in response to auditory stimuli we developed a new DCM with six cell populations, called “ext-DCM.” In six connected regions (locations from Phillips et al., 2015, 2016), we used a conductance-based mean-field cortical modeling scheme (cf. Moran et al., 2013; Marreiros et al., 2015). For auditory mismatch responses, both thalamocortical and corticocortical connections integrate feedforward sensory inputs and feedback expectations. The network architecture controls the flow and integration of information, via cell- and neurotransmitter-specific interactions. The ext-DCM introduces new cortico-thalamic burst-firing cells (“tp” in Fig. 1a) that enable the model to generate beta activity from deep-layers (Roopun et al., 2008a, 2010; Bordas et al., 2015; Michalareas et al., 2016). The ext-DCM also separates the inhibitory interneuronal populations for superficial and deep pyramidal cells (Jiang et al., 2015). These extensions improve the DCMs' functionality in terms of laminar dynamics. We tested the model's ability to accurately generate evoked magnetoencephalographic responses (i.e., event-related fields, ERFs), under placebo and drug conditions.
The neuronal model. a, Intrinsic connectivities found in all nodes between layer 4 stellates (ss), inhibitory interneurons (ii), superficial pyramidal modules (sp), and deep pyramidal modules (dp). b, All six nodes used are represented as a network on the left, showing the extrinsic connectivities (solid line, forward; dotted line, backward; dashed line, lateral). A left hemisphere representation of these bilateral nodes in primary auditory cortex, superior temporal gyrus, and inferior frontal gyrus is shown in light, medium, and dark gray, respectively. c, Detailed view of the extrinsic population connections for forward (solid lines) and backward (dotted lines) connections. d, Matrices of the extrinsic and intrinsic connectivity weights, all of which had a permitted variance of 1/16. e, A process flow describing the steps taken in the meta-analysis phase.
With the ext-DCM, we used the drug tiagabine to test how well the neurophysiological model could identify changes in the causes of observed neuronal dynamics. Tiagabine is a gamma-amino-butyric acid (GABA) reuptake inhibitor. GABA is critical for the generation of physiological responses and rhythms in local and global processing (Whittington et al., 2000). This pharmacological specificity provides a more controlled acute test of DCMs than autoimmune (Symmonds et al., 2018) and genetic channelopathies (Gilbert et al., 2016).
Using parametric empirical Bayes to optimize the model across participants and drug conditions we examined how modeled GABAergic dynamics are altered by tiagabine. Based on the hypothesis that prediction and prediction error depend on short-term GABAergic plasticity (Castro-Alamancos and Connors, 1996; Garrido et al., 2009a; Mongillo et al., 2018; Spriggs et al., 2018), we predicted that upper and lower hierarchical frontotemporal processing would be differentially affected by tiagabine during standard and deviant tones.
In summary, the study's principal aims were as follows: (1) to introduce and assess the ext-DCM for generating the event-related fields observed by MEG, (2) to identify receptor-specific changes that govern these dynamics, comparing tiagabine and placebo treatment conditions, and (3) to assess whether these pharmacological effects are expressed dynamically across trial types and regions with laminar specificity.
Materials and Methods
Experimental design
We undertook a randomized placebo-controlled double-blind crossover study of the effects of tiagabine in 20 healthy adults (aged 67.5 ± 4.2, 10 male). Participants had no neurological or psychiatric illness and were recruited from the MRC Cognition and Brain Sciences and Join Dementia Research volunteer panels. The study was approved by the Cambridge Research Ethics Committee and written informed consent was acquired, in keeping with the declaration of Helsinki.
Neurophysiological responses were measured in an auditory roving oddball paradigm (Garrido et al., 2008). Binaural sinusoidal tones were presented in phase via ear-pieces for 75 ms (with 7.5 ms ramp up and down at start and end of the tone), at 500 ms intervals. The frequency of the tone increased or decreased in steps of 50 Hz (range 400–800 Hz). The change of frequency occurred after between 3 and 10 repetitions, with a truncated exponential distribution that approximated a stable expectancy of change over time. Auditory thresholds were assessed in quiet at 500, 1000, and 1500 Hz. Tones were presented at 60 dB above the average threshold for a standard population through the earpieces in the MEG.
Each participant attended two MEG sessions with a minimum 2-week interval. They received either 10 mg oral tiagabine or a placebo in randomized order. Blood samples were taken 105 min later, immediately before MEG data acquisition, to coincide with peak plasma levels and CNS penetration (Nutt et al., 2015).
Data acquisition and preprocessing
MEG used a 306-channel Vectorview acquisition system (Elekta) in a light Elekta Neuromag magnetically shielded room. This consists of a pair of gradiometers and a magnetometer at each of 102 locations, sampled at 1000 Hz. Vertical and horizontal EOGs tracked eye movements and 5 head-position indicator coils tracked head position. A MEG-Compatible 70 channel EEG cap (Easycap) using Ag/AgCl electrodes positioned according to the 10–20 system was used concurrently. A 3D digitizer (Fastrak; Polhemus) was used to record >100 scalp data points, nasion and bilateral preauricular fiducials. Subjects also underwent T1-weighted structural magnetic resonance imaging (MPRAGE sequence, TE = 2.9 msTR = 2000 ms, 1.1 mm isotropic voxels) using a 3 T Siemens PRISMA scanner.
MEG data preprocessing included head position alignment and movement compensation using six head coils placed around the head on the EEG cap, and used the temporal extension of Signal Space Separation with MaxFilter v2.2 (Elekta). The auto-detection of bad channels was combined with manual input of any channels logged as bad during data acquisition. The Statistical Parametric Mapping toolbox (SPM12) (The Wellcome Trust Centre for Neuroimaging, UCL, UK) was used for further preprocessing and analysis, in conjunction with modified and custom MATLAB scripts (MATLAB 2017a, The MathWorks). Data were Butterworth filtered between 1 and 180 Hz, epoched from −100 to 400 ms relative to the auditory stimuli, and artifact rejected using EOG, EEG, and MEG channel thresholding. Spectral analyses were performed using a multitaper method. The deviant trial was taken as the first trial of a train, regardless of the frequency and the sixth trial of a train was modeled as “standard.”
Source reconstruction used a forward model estimated using the single shell cortical mesh from each individual's T1-weighted MR structural scan. After coregistration using the fiducials and head points, local fields (LFs) for 6 sources of interest were source-reconstructed using SPM “COH” method, a combination of LORETA and minimum norm (Pascual-Marqui et al., 1994; Heers et al., 2016). Sources of interest were (with MNI coordinates in standard space following inverse normalization): left auditory cortex (LAud; −42, −22, 7), left superior temporal gyrus (LSTG; −61 −32 8), left inferior frontal gyrus (LIFG; −46 20 8), right auditory cortex (RAud; 46, −14, 8), right superior temporal gyrus (RSTG; 59 −25 8) and right inferior frontal gyrus (RIFG; 46 20 8). To create images of induced power, SPM-LORETA was used for source localization of a 5 mm3 regular grid at the MMN (150–250 ms) time window (100 ms in width, regularization = 0.05).
Correlation coefficients for comparing the actual and predicted ERFs were calculated using the corrcoef function (Pearson correlation) in MATLAB 2017a for each individual, condition and node.
Time-frequency analysis was performed in SPM12 using a multitaper method with 100 ms windows overlapped by 5 ms and a bandwidth of 3. Frequency bands were split into alpha (8–13 Hz), beta (14–29 Hz), low gamma (30–48 Hz), and high gamma (52–80 Hz).
Neuronal modeling: an extended canonical microcircuit model
We used conductance-based canonical mean field (CMM) models for evoked responses (Kiebel et al., 2008) using canonical microcircuit models (SPM12, DCM10). This approach to neurophysiologically informed modeling using DCM goes beyond descriptive biomarkers by providing a mechanistic link to realistic microscopic processes. A common approach in DCM is to invert the neuronal and spatial forward model as a single generative model, to solve the source reconstruction and biophysical modeling problems jointly by fitting the DCM to sensor data. However, we modeled source specific responses to suppress conditional dependencies between the neuronal parameters and the parameters of a spatial forward model. This affords more efficient estimators of neuronal parameters, providing the source reconstruction is sufficiently precise given the spatial topography of the network of interest. This has the advantage of compatibility with multiple studies of this task (Muthukumaraswamy et al., 2015; Gilbert and Moran, 2016; Shaw et al., 2017, 2019), including MEG and electrocorticography studies; the chosen network was based on the published bilateral A1, STG, IFG networks associated with the generation of the MMN response. Since this spatial element of the inverse problem was constrained, it is computationally more appropriate to source localize using SPM with prior expected sources. The subsequent DCM was then run on these virtual electrodes.
Our DCM included a conductance-based neural-mass model at each of the six anatomical locations, as shown in Figure 1. We compared the default 4-cell conductance canonical-microcircuit model with the ext-DCM, comprising 6 cell modules: a superficial pyramidal module (sp), a deep corticocortical pyramidal module (dp), a thalamic-projection pyramidal module (tp), a granular stellate module (ss), and separate supragranular and infragranular interneuron populations (si and di). Excitatory autapses existed for all excitatory cell modules and all modules were also governed by an inhibitory self-gain function that provided tonic inhibition to each module. The ext-DCM was compared with the standard 4-cell model that is standard in SPM and is described in detail in Kiebel et al. (2008). In summary, the 4-cell DCM lacks thalamocortical connectivity and has a unitary inhibitory population interacting with all pyramidal and stellate cell populations.
The intrinsic connectivities are shown in Figure 1a; note the excitatory conductances based on AMPA and NMDA and inhibitory GABA-A and GABA-B conductances. The model is an extension of the SPM conductance-based CMM model (SPM12, 2013): inclusion of separate supragranular and infragranular interneuron populations creates a more biophysically realistic model that allows a greater flexibility of independence of deep and superficial activity than in previous work (Bhatt et al., 2016; Shaw et al., 2019; Spriggs et al., 2018). Additionally, the new 'tp' population expressed a hyperpolarization-activated cation current (H-current) and a non-inactivating potassium current (M-current) to provide surrogate intrinsic dynamics involved in the characteristic intrinsic bursting behavior of these cells. These two currents were fixed together with the reversal potential and the slope on the sigmoid convolution of in-activation for the H-current (details of which parameters had a permitted variance is given in Table 1). This, coupled with the cell capacitances, differentiates the intrinsic activation of the “tp” population from the “dp” population. The populations also differed in their extrinsic connectivities, with dp populations forming corticocortical connections and tp populations allowing for cortico-thalamocortical connections. The thalamus was modeled implicitly, by an 80 ms delay in connectivity with permitted variance.
Parameter values used by the neuronal model shown with their permitted variances
Extrinsic connectivity between the six nodes is shown in Figure 1b, with the detailed extrinsic population connections shown in Figure 1c. In keeping with the established principle of differential cortical laminar projections of feedforwards vs feedback connectivity (Bastos et al., 2012), backward connections are facilitated by the dp cells terminating on “sp” and “si” cells, while forward connections run from sp cells to “ss” cells. Cortico-thalamo-cortical connections originate from tp cells and terminate following a thalamic delay at layer 4 ss cells. The presence or absence of connections between nodes was based on the fully connected models from Phillips et al. (2015) and Shaw et al. (2019), which in turn were derived from Garrido et al. (2008). This was used for the basis of an iterative process to find the most likely reduced model (described below).
A Gaussian kernel (peak 60 ms, half-width 8 ms) represented auditory input to layer 4 stellates in bilateral auditory and inferior frontal cortex.
Bayesian modeling and statistical analysis
We used Bayesian model inversion (estimation) and Bayesian model comparison (selection) to identify the best explanation for subject-specific data, in terms of neuronal and biophysical parameters. Parametric Empirical Bayes (PEB) was used for group inferences and to examine drug effects, as described in Zeidman et al. (2019). By inverting a “full” DCM per subject at the first level, PEB avoids the problem of different first-level DCMs falling into different local optima, and allows subsequent comparison between conditions. At the second level, the parameters of interest were included in the PEB, namely the GABAA synaptic connections. This restricted set of second level parameters was oriented to our GABA-ergic hypothesis, and to improve stability of neural system identifiability.
The DCM was run for each subject. Data were filtered between 0 and 48 Hz and a Tukey window was applied that did not attenuate signals 50 ms before or 350 ms after stimuli. Inversion of the full model was run separately for the standard and deviant trials and the parameter distributions passed to second level parametric empirical Bayesian with contrasts for both trial types and drug conditions. All intrinsic and extrinsic AMPA, NMDA, and GABA-A conductance scalings could vary independently in a manner that assumed symmetry between the two hemispheres. The prior means and permitted variances are summarized in Table 1.
Variational Bayesian statistics using the Laplace approximation determined the probable parameter space given the neuronal model and the data (Friston et al., 2007). The full model parameter space was reduced by iteratively searching for dependencies in this parameter space and systematically removing parameters not contributing to the free energy of the system (Henson et al., 2011). The optimized reduced model comprises all those parameters and connections found to contribute significantly to the system temporal dynamics. The comparison of full and reduced models is conceptually analogous to F tests in classical statistics, but inferences are Bayesian. A second-level PEB was run, optimizing GABAA-ergic synaptic parameters (representing inhibitory gain). This second level PEB identifies parameter that are estimated to differ significantly between task conditions, or differ between drug-sessions, or for which there is a drug-by-condition interaction. The parameter distributions from this reduced model were used to create a Bayesian model average of parameters that differ significantly across the contrasts of trial types and drug conditions. The implementation of PEB for model optimization and contrast estimation is summarized in Figure 1e.
For other data types, Bayesian t tests reported in the main text used JASP (JASP Team 2019, version 0.10.2). Frequentist statistical methods reported in the main text used MATLAB 2017a.
Code accessibility
The custom neuronal model used to generate these results is available at https://gitlab.com/tallie/edcm and works in conjunction with SPM12.
Results
Event-related fields and induced spectral power
Event-related responses to standard and deviant trials were consistent with previous findings (Hughes and Rowe, 2013; Phillips et al., 2015, 2016) (Fig. 2a, first and second rows) and show the expected M100, the primary response after the onset of a tone (80–120 ms), a difference signal (MMN) between the standard and deviant trials (150–250 ms) and an M300 visible in frontal nodes (250–380 ms). The M100 was significantly reduced by tiagabine on standard and deviant trials, in left temporal nodes (A1, and STG p < 0.05, paired t test), whereas the later response leading into the M300 was significantly reduced only on deviant trials in L/R IFG (p < 0.05, Bonferroni corrected for 6 regions).
ERFs. a, Mean ERFs across all subjects for all six nodes for the standard and deviant trials from 0 to 380 ms. The difference wave (MMN) is also shown. ERFs from the placebo condition are shown in blue and from the tiagabine condition in red. Significant changes with time across the drug condition are shown as a thick black line within each axis (p < 0.05, Bonferroni corrected for 6 regions). Shaded areas represent SEM. b, Significant differences for induced spectra power were found in the alpha (α), beta (β), and lower and higher gamma bands (γ1 and γ2) (FWE cluster corrected at p < 0.001). Here they are shown as flat scalp maps (lower plots) with rostro-caudal activity versus time (upper plots). The time axis runs from 0 to 380 ms poststimulus. c, Source-reconstructed T-contrasts (p < 0.001) created for those frequency bands showing spatial changes across the drug condition in the 135–235 ms time window.
The difference waveform (i.e., the deviant-the standard) reveals a typical biphasic MMN between 150 and 250 ms, observed in primary auditory cortex and STG (Fig. 2a, third row). Tiagabine significantly reduced the second peak of the MMN (p < 0.05) with bilateral IFG nodes and RSTG showing reductions in the first peak of the mismatch response on tiagabine (p < 0.05). As with the deviant response, LIFG showed a significant reduction of the later MMN peak and the M300 on tiagabine (p < 0.05).
The temporal profile of spectral power differences (see Materials and Methods for time-frequency analysis) matched that of the ERFs, including spectral counterparts to M100, MMN, continuing through the M300 window (Fig. 2b,c). During the M100, alpha-power (8–12 Hz) decreases on tiagabine were localized to temporal cortex and beta (14–29 Hz) decreases more prominently to posterior temporal cortex. During the MMN, increases in low and high gamma (30–48 Hz and 52–80 Hz, respectively) were observed broadly across right frontal cortex, including IFG. Low gamma also showed increases in right temporal cortex.
Such changes in the observed spatiotemporal physiology on tiagabine will be dependent on changes in local and global network connectivity. The extended conductance-based dynamic causal model was therefore used to infer the causes of the observed physiological changes.
Dynamical causal model
The residuals (difference between the actual and generated ERFs) were greater (worse) for the 4-cell DCM than for the ext-DCM (Bayesian paired sample t test: BF = 8.5 × 1028) as shown in Figure 3a. Bayesian model comparison of the 4-cell versus ext-DCM confirmed that the ext-DCM performed better (i.e., was a more likely generator of the observed MEG) than the 4-cell DCM (BF = 40.6, Fig. 3b). Note that the model-evidences are corrected for differences in model complexity. Further analyses use the ext-DCM only.
Comparison between model and data. a, Residual differences between the observed and model-generated ERFs are shown for both the standard four-cell conductance-based DCM and the ext-DCM. ERFs from all nodes for every subject are concatenated along the y-axis. b, Bayesian model comparison of the four-cell conductance-based DCM and the ext-DCM favors the ext-DCM, plotted here in terms of the posterior model probability (RFX Bayes factor = 40.6). c, Predicted ERFs are shown for the standard and deviant conditions, along with the difference wave (Std-Dev). The placebo and tiagabine conditions are depicted in blue and red, respectively, with significant differences (p < 0.05, Bonferroni corrected for 6 regions) shown as a thick black line within each axis. d, Correlation coefficient between prediction and data for each node and each condition. Box plots represent the distribution over subjects with small dots representing outliers and larger black circles representing the correlation coefficient of the mean response of all subjects for each node and each condition.
Figure 3c demonstrates the evoked-response generated by the conductance-based dynamic causal model at each node, for both drug conditions, using the optimal ext-DCM model as determined by parametric empirical Bayes (see the Materials and Methods). Figure 3d shows the correlation between generated and observed data, for both standards' and deviants' responses, for both drug conditions at each node. Boxplots indicate the spread of single-subject correlations across the group (open circles are outliers), and black closed circles indicate the correlation of the mean response across all subjects for each condition and node. Note how the periods of difference between the placebo and drug conditions (black lines in Fig. 3c) are accurately generated (cf. 'predicted') by the model, with a high match to the observed data in Figure 2a.
The modeled responses are explained in terms of the parameters of the optimized model. Using parametric empirical Bayes, condition effects on model parameters (connection and synaptic parameters) were compared across the standard and deviant conditions, as well as across the placebo and tiagabine conditions. Figure 4 shows the effect of tiagabine on the intrinsic GABAergic connectivity, assuming symmetry (three bilateral averaged nodes are shown). We confirmed that tiagabine significantly increases tonic GABAergic inhibition (posterior probability given for each parameter in Fig. 4a). This was seen primarily in the deep layer pyramidal and interneuron populations in primary auditory cortex and STG (Fig. 4a). An interaction between drug and condition was found for the deep interneurons of auditory cortex (posterior p ≈ 1.0).
Prediction of hidden states. a, Significant differences in the modulation of GABA-A synaptic scaling for each of the three symmetric nodes. Green/red show significantly greater/lesser GABA-A synaptic scaling for tiagabine than the placebo. Posterior probability p-values are shown next to each connection. b, To explore the functional differentiation between regions during the task conditions with respect to tonic inhibition, tonic GABA-A scaling on deep interneurons in IFG, STG, and Aud, for each individual is plotted for the placebo and tiagabine conditions. The standard and deviant conditions are plotted separately in the left and right columns respectively. Pairwise Bayesian t test statistics are reported on each plot, showing the Bayes factor for each of the six comparisons. When there is evidence for a difference, or evidence for no difference, the Bayes factor is shown in green or blue respectively. c, Correlation demonstrating the dynamic balance that persists between phasic and tonic inhibition (see main text discussion). Linear fit with 95% confidence bounds for tonic GABA-A scaling on deep inhibitory neurons versus phasic GABA-A scaling from deep inhibitory neurons to thalamic projecting pyramidals is shown (Bayesian correlation pairs, Bayes factor = 398.43).
Figure 4b compares GABA-A conductance scaling on deep interneurons between placebo and tiagabine conditions, plotted for each individual. There was very strong evidence for differences between the two drug conditions in primary auditory areas for the standard condition (BF = 782356), and in IFG and STG for the deviant condition (BF = 3.58 × 107 and BF = 166 respectively). This difference between primary auditory cortex and association cortex in STG/IFG, is in keeping with the functional differentiation of upper versus lower levels in a hierarchical neural network with backwards prediction and forward prediction error. Conversely, there was evidence of no difference between the two drug conditions for the standard condition in IFG (BF = 0.274) and for the deviant condition in Aud (BF = 0.241).
The correlation between tonic and phasic inhibition was explored for each region and condition. In the frontal cortex, a strong negative relationship was found between the tonic inhibition of deep inhibitory cells and their phasic inhibition onto cortico-thalamic cells (Fig. 4c, Bayesian correlation pairs, BF = 398.43).
Discussion
The principal insights from this study are that an extended conductance-based canonical mean-field method of dynamic causal modeling succeeds in identifying the modulation of GABAergic dynamics by the GABA-reuptake inhibitor tiagabine and is tractable and an accurate generator of event-related fields that match those observed by magnetoencephalography, improving on an earlier four-cell model. Moreover, the ext-DCM suggests the effect of drug to be both laminar-specific and dynamically modulated in different regions according to task condition. This opens the way for psychopharmacological studies in health and disease with the mechanistic precision afforded by using ext-DCMs as generative models.
We demonstrate that the intrinsic connectivity within hierarchical brain networks changes between conditions in the mismatch task. The approach is of generalized relevance to hierarchical network models of cognition such as speech (Cope et al., 2017), semantic (Adams et al., 2019) and visual perception (Muthukumaraswamy et al., 2013a,b). Moreover, the laminar and pharmacological specificity provided by the ext-DCM has the potential to quantify neuropathology in dementia, developmental and psychiatric disorders (Duyckaerts et al., 1986; Kinoshita et al., 1996; Ferrer, 1999; Ji et al., 2018; Shaw et al., 2019).
Understanding the MMN in terms of short-term plasticity
Tiagabine modulated the GABA-egic dynamics across the trial types, implicating both local tonic- and phasic effects. Repetitive activation with the same stimulus attenuated the ERF (reduction in N1/N2 by sixth repetition, Fig. 2). The model indicated higher tonic inhibition in the deep layers. We interpret this as local short-term plastic changes in deep-layer inhibition (Knott et al., 2002; Hensch, 2005; Jääskeläinen et al., 2007), regulating salient information (Mongillo et al., 2018).
The model suggested that tiagabine-induced increases of extracellular GABA leads to greater tonic inhibition, consistent with overspill of GABA onto extra-synaptic receptors (Semyanov et al., 2004). The effect was modulated differently in primary and associative processing areas: for tonic inhibition of deep interneurons the drug's efficacy was highest in prefrontal cortex for deviant trials and in auditory cortex for standard trials. In other words, GABAergic effects are modulated differentially in upper and lower areas of the hierarchy dependent on the coding context. We speculate that this reflects differential emphasis on beliefs (and feedback predictions) versus feedforward sensory prediction errors in prefrontal versus primary auditory cortex; and that lower tonic inhibition at the presentation of a deviant tone relates to homeostatic competition between phasic and tonic inhibition (Wu et al., 2013). Increased phasic activation of deep-layer projections is necessary for feedback of top-down information on context, which in turn increases phasic (and decreases tonic) activation of deep interneurons. Decreasing tonic inhibition likely increases the interneuron population activation (Semyanov et al., 2004), leading to increased phasic inhibition onto deep pyramidal cells. This relationship was confirmed (Fig. 4c) between tonic inhibition of deep IFG interneurons and phasic inhibition of deep IFG thalamic-projection neurons. Figure 4b shows that whereas a drop in deep interneuron tonic inhibition was observed on deviant trials (vs standard), tiagabine abolished the effect. It is to be expected that increases in exogenous GABA would increase tonic GABAergic currents.
GABA-ergic modulation of evoked and induced responses
Tiagabine affects oscillatory dynamics, which may influence behavior (Coenen et al., 1995; Magazzini et al., 2016; Port et al., 2017; Wyss et al., 2017). It remains a challenge to relate systemic drug effects with local frequency-spectral phenomena. It has been proposed that beta-band activity is associated with infragranular cortical projection neurons with intrinsically bursting profiles (Groh et al., 2010; Roopun et al., 2010; Kim et al., 2015). We found that Tiagabine reduced the induced beta-band activity in temporal areas. The model suggests that tonic inhibition is increased on intrinsically bursting thalamic projection neurons in STG, which could increase rebound bursting via intrinsic M- and H-currents (Roopun et al., 2008a,b).
Conversely, it has been shown that gamma-band activity is dependent on the GABA-A receptor activation and the phasic interplay of interneuron-pyramidal cell networks, particularly in the superficial layers (Buffalo et al., 2011; Whittington et al., 2011). In the mismatch temporal window (Fig. 2b) peak gamma increased occurring at the start of the mismatch period. This is consistent with thalamic input (Di and Barth, 1992, 1993; Sukov and Barth, 2001) governing the envelope of gamma activity in the superficial layers (Metherate and Cruikshank, 1999).
Overall, the observed dynamics and the model posterior parameters are consistent with knowledge of network activation within the context of beta- and gamma-rhythm generation in cortex.
Generative models of drug effects on cognitive physiology
Tiagabine's effect was largely confined to deep layers. As we modeled evoked activity it is difficult to speculate on how this influences gamma activity across the network, however, a reduction in deep-layer influence may increase local cortical processing associated with gamma-band activity in the superficial layers. As GABA levels are typically lower in older versus younger adults, tiagabine may act “restoratively.” This is corroborated with lower-frequency responses that are dependent on GABA (Mathias et al., 2001). Finally, we speculate that the reduced M100 on tiagabine results from the widespread increased tonic inhibition represented in the model (Fig. 4), reducing local population activity.
Study limitations
Our study was motivated by the need for mechanistic studies of human cortical function, underlying cognition, disease and therapeutics. Despite support for our three principal hypotheses, and background validation studies (Moran et al., 2014), evidence from one study may not generalize to other tasks and populations. There are study-specific considerations that limit our inferences, in relation to our participants, our model, and drug of choice. For example, our participants were healthy, and therefore have normal age related variance in GABA (Gao et al., 2013; Eavri et al., 2018). They were older than those studied by Nutt et al. (2015), and age effects could interact with the effects of tiagabine (Nutt et al., 2015). Our study was not designed to examine the effect of age or aging, but rather to focus on the normal brain in midlife and later life. Further work would be required to examine the effects of aging on the ext-DCM.
Our model provides a simplified substrate for the neurophysiological processes. It is more detailed than previous canonical microcircuit convolution models (Moran et al., 2013) in an effort to improve the modeling of specific dynamics from distinct cell populations, their differing connectivities, synaptic time constants, and voltage-gated conductances. The extended model can produce a spectrum of fast and slow responses, with fast responses involved in local processing dominated by superficial layers and slower responses associated with feedback of information dominated by deep layers (Roopun et al., 2006; Kramer et al., 2008; Whittington et al., 2011). It can incorporate delayed activity associated with local, corticocortical and cortico-thalamo-cortical connections. Currently, this system is a simplified network acting as a neural mass, and can represent relevant cortical interactions involved in ERF generation in the context of this task and study. It does this by allowing forward and backward modulation of activity between deep and superficial layers, where synaptic time constants corroborate with standard GABA, NMDA and AMPA receptor decays. The six specified nodes are commonly cited in the literature in the context of this task (Garrido et al., 2009b; Phillips et al., 2015). Although they are not a complete representation of possible network configurations, they have been shown to capture critical aspects of cortical function: here the network has been supplemented with modeled exogenous and endogenous inputs via thalamus.
We emphasize Bayesian statistical analyses over classical frequentist methods. Where parameter estimates derived from earlier DCMs are used for frequentist statistical tests, they have excellent reliability across sessions, and similar power to fMRI and EEG studies (Rowe et al., 2010; Goulden et al., 2012; Bernal-Casas et al., 2013). Frequentist approaches are familiar to many readers, and have been the norm for comparison of ERFs, and we therefore include them selectively. Such a frequentist approach is surpassed by the direct inferences on posterior probability inherent in DCMs' Bayesian inference, including PEB.
Tiagabine is a relatively specific blocker of GAT-1 at the concentrations used, but does not distinguish between the mechanisms activated by GABA (Bowery et al., 1987; Mody and Pearce, 2004; Lee and Maguire, 2014). The timing of the MEG coincided with expected peak plasma levels, but levels may vary between individuals and future studies could include levels as a covariate of interest, or model time-varying responses in relation to drug levels (Muthukumaraswamy et al., 2013b).
In conclusion, we have used a conductance-based model of cortical neuronal dynamics to study GABA-ergic interactions and probe laminar-specific physiological responses to tiagabine. The model accurately generated physiological data that matched the MEG responses and confirmed the effect of tiagabine on tonic GABA-A inhibitory gain within frontal and temporal cortical circuits. Our data provide support for mechanistic studies of neurological disorders, including but not limited to GABAergic impairments (Murley and Rowe, 2018). They also point to new approaches for experimental medicine studies in humans that aim for the laminar, cellular, or synaptic precision made possible in new generations of dynamic causal models.
Footnotes
This work was supported by the Wellcome Trust (Grant 103838), the National Institute for Health Research Cambridge Biomedical Research Centre and the Medical Research Council (Grants MC_U105597119, MC_U_00005/12, and SUAG/004/91365), the Cambridge Centre for Parkinson-plus, and by a Holt Fellowship. A.D.S. is supported by a Wellcome Trust Strategic Award (104943/Z/14/Z).
The authors declare no competing financial interests.
- Correspondence should be addressed to James B. Rowe at james.rowe{at}mrc-cbu.cam.ac.uk
