The beta component of gamma-band auditory steady-state responses in patients with schizophrenia

The mechanisms underlying circuit dysfunctions in schizophrenia (SCZ) remain poorly understood. Auditory steady-state responses (ASSRs), especially in the gamma and beta band, have been suggested as a potential biomarker for SCZ. While the reduction of 40 Hz power for 40 Hz drive has been well established and replicated in SCZ patients, studies are inconclusive when it comes to an increase in 20 Hz power during 40 Hz drive. There might be several factors explaining the inconsistencies, including differences in the sensitivity of the recording modality (EEG vs MEG), differences in stimuli (click-trains vs amplitude-modulated tones) and large differences in the amplitude of the stimuli. Here, we used a computational model of ASSR deficits in SCZ and explored the effect of three SCZ-associated microcircuit alterations: reduced GABA activity, increased GABA decay times and NMDA receptor hypofunction. We investigated the effect of input strength on gamma (40 Hz) and beta (20 Hz) band power during gamma ASSR stimulation and saw that the pronounced increase in beta power during gamma stimulation seen experimentally could only be reproduced in the model when GABA decay times were increased and only for a specific range of input strengths. More specifically, when the input was in this specific range, the rhythmic drive at 40 Hz produced a strong 40 Hz rhythm in the control network; however, in the ‘SCZ-like’ network, the prolonged inhibition led to a so-called ‘beat-skipping’, where the network would only strongly respond to every other input. This mechanism was responsible for the emergence of the pronounced 20 Hz beta peak in the power spectrum. The other two microcircuit alterations were not able to produce a substantial 20 Hz component but they further narrowed the input strength range for which the network produced a beta component when combined with increased GABAergic decay times. Our finding that the beta component only existed for a specific range of input strengths might explain the seemingly inconsistent reporting in experimental studies and suggests that future ASSR studies should systematically explore different amplitudes of their stimuli. Furthermore, we provide a mechanistic link between a microcircuit alteration and an electrophysiological marker in schizophrenia and argue that more complex ASSR stimuli are needed to disentangle the nonlinear interactions of microcircuit alterations. The computational modelling approach put forward here is ideally suited to facilitate the development of such stimuli in a theory-based fashion.

Interestingly, 50% of PV + neurons in the dorsolateral prefrontal cortex of SCZ patients have very low levels of the GAD67 isoform of glutamate decarboxylase 12 . This reduced expression of GAD67 mRNA has been demonstrated to decrease the GABA synthesis in cortical GABAergic neurons, which would then result in smaller amplitudes of inhibitory postsynaptic currents (IPSCs) 8 . Additionally, PV + neurons show reduced levels of the plasma membrane GABA transporter GAT1 in SCZ patients 13 . A reduced concentration of GAT1 has been shown to increase the time GABA molecules reside at the receptor and thus increase IPSC durations 14 .
Administration of NMDAR antagonists leads to the emergence of schizophrenia-like symptoms, such as hallucinations, delusions and thought disorder, in healthy subjects 15 . Based on these findings it has been hypothesized that the reduced inhibition found in SCZ might not be a consequence of the changes to PV + neurons described above, but could be attributable to an NMDAR hypofunction. Dysfunction of NMDARs in SCZ is supported by several lines of evidence 10 . Specifically, Carlen et al. 16 found that targeted deletion of NMDARs from PV + interneurons led to increased spontaneous gamma oscillations and a deficit in gamma induction. Interestingly, they could reproduce these results in an established circuit model 17 when they implemented NMDAR hypofunction as an overall decrease in interneuron excitability.
While the reduction of 40 Hz power for 40 Hz drive has been well established and replicated in SCZ patients 5,6,17 , studies are inconclusive when it comes to an increase in 20 Hz power during 40 Hz drive 5,17,18 . There might be several factors explaining the inconsistencies, including differences in the sensitivity of the recording modality (EEG vs MEG), differences in stimuli (click-trains vs amplitude-modulated tones) and large differences in the amplitude of the stimuli. Kwon et al. 6 used a click-train paradigm to study ASSRs at 20, 30, and 40 Hz in SCZ patients using EEG and found a prominent reduction of power at the driving frequency for 40 Hz drive but no changes of power at the driving frequency for 30 Hz and 20 Hz. Although Figure 3 in Kwon et al. 6 seems to show an increase of the subharmonic 20 Hz component for 40 Hz drive, no statistical comparison is presented in the article. Vierling-Claassen et al. 17 reproduced this reduction of power at the driving frequency for 40 Hz drive using the same paradigm with MEG. Additionally, they found changes of power at subharmonic frequencies, namely, an increase of power at 20 Hz for 40 Hz drive. Furthermore, Spencer et al. 19 also found an increase of the 20 Hz component for 40 Hz drive. Importantly, both Vierling et al. and Spencer et al. used relatively low stimulus amplitudes of 65, and 55 dB, respectively (Kwon et al. do not report the amplitude of their stimuli). However, three studies did not find differences in the 20 Hz component for 40 Hz drive between SCZ patients and healthy controls 18,20,21 , all using relatively high stimulus amplitudes of 80, 76 and 80 dB, respectively.
In this study, we used an established network model of ASSR deficits in SCZ 17,22 , where SCZ-like behaviour is produced by an increase in IPSC decay times, to examine the dependence of 40 Hz ASSRs on the strength or amplitude of the inputs. In our model we could only reproduce the emergence of 20 Hz component during 40 Hz stimulation seen experimentally if the input strength was in a narrow range. More specifically, very weak input to the network did not result in a pronounced oscillatory rhythm. When the input was in a specific range, the 40 Hz stimulation entrained a pure 40 Hz oscillation in the control network, whereas in the 'SCZ-like' network, the changed IPSC time course caused a so-called 'beat-skipping' , where the network would only strongly respond to every other input. This resulted in significant increase in 20 Hz power. Ultimately, if the input became too strong the increased IPSC decay time was insufficient to suppress the very strong 40 Hz input. This was reflected in a single peak at 40 Hz in the power spectrum. We then extended the network model to include more cellular-level alterations such as reduced GABA levels and NMDAR hypofunction. We found that the addition of further alterations did not change the input strength dependence of the 20 Hz component but further limited the parameter range where the component would occur. Our finding that the beta component only existed for a specific range of input strengths might explain the seemingly inconsistent reporting in experimental studies and suggests that future ASSR studies should systematically explore different amplitudes of their stimuli.

Results
Replication of previous findings. First, we validated the IPSC-SCZ-like model against experimental observations 17 and replicated the findings from previous modelling studies with this model 17,22 . The control network strongly entrains to the driving stimulus, regardless of the specific driving frequency (Figs. 1 and 2 left columns), and shows stronger entrainment at 40 Hz than for 30 and 20 Hz, consistent with experiments 17 . Furthermore, the control model replicates another important feature seen in experiments and the previous models, a strong 40 Hz component for 20 Hz drive (Figs. 1 and 2 left columns, third rows). The 'SCZ-like' network, where 'SCZ-like' behaviour is achieved by an increase in the GABAergic IPSC decay time constant (from 8 to 28 ms as in earlier studies), also reproduces important characteristics from experiments and previous models: First, the 'IPSC-SCZ-like' network shows a marked reduction in 40 Hz power for 40 Hz drive (Fig. 2, right column, first row), as previously found in experiments (see 5 for a meta-analysis) and models 17,22 . Furthermore, this network shows an emergent 20 Hz component at 40 Hz drive (Fig. 2, right column, first row) as seen in 6,17,22 but not in other experimental studies 5 ; and we see an increase in 20 Hz power and a relative decrease in 40 Hz power for 20 Hz drive in this condition (Fig. 2, right column, third row).
Input strength dependence of the 20 Hz component. Next, we explored the input strength dependence of the IPSC-SCZ-like model response to 40 Hz drive by multiplying the input strength of the SCZ model by factors ranging from 0.1 to 1.5 in steps of 0.1. Figure 3a shows that the 40 Hz power increases with increasing input strength. Figure 3b shows that the 20 Hz component emerges at an input strength of 80% of the default strength, sharply increases for stronger inputs around the the default SCZ network strength and then sharply decreases again for inputs of 120% of the default strength and higher. Thus, in a narrow range between 80 and 120% of the default inhibitory input strength the network response exhibits a shift of power from the gamma (40 Hz) to the beta (20 Hz) band. Figure 3c, www.nature.com/scientificreports/ to force the network into a coherent rhythm and, therefore, the powers at 40 Hz and at 20 Hz are very low. For input strengths around the default SCZ network value, the input strongly drives the network and forces it into a rhythm. However, the increased IPSC decay times prevent the excitatory pyramidal neurons from responding to every 40 Hz cycle and only allow them to spike every other cycle (Fig. 3d,g). Thus, the network rhythm displays a so-called 'beat-skipping' behaviour. The power spectrum of the response therefore shows both a 40 Hz (which is substantially smaller then for the control network) and a prominent 20 Hz peak (which is not visible for the control network). If the input, however, exceeds 120% of the default SCZ network value, the rhythmic input becomes strong enough to overcome the prolonged inhibition and forces the network into a gamma oscillation at 40 Hz. Here, the excitatory cells fire during each cycle and the power spectrum only shows a large 40 Hz peak (Fig. 3e,h).

Combinations of alterations and their input-strength-dependence. As explained earlier it is
unlikely that the microscopic alterations associated with schizophrenia occur in isolation. Therefore, we added two more alterations to the IPSC-SCZ-like model: (1) A reduction in GABA levels, implemented as a reduction of the inhibitory weights; (2) A hypofunction of NMDARs at inhibitory interneurons, implemented as a reduction in interneuron excitability. We first added these modifications individually and combined them in a final set of simulations.  Lastly, the full model combining all three alterations, displayed similar behaviour as the previous models but an even more pronounced shift of the 40 Hz components to lower input strengths for higher levels of GABA (Fig. 6).

Beta component in response to 20 Hz drive.
After exploring the 20 Hz component in response to 40 Hz ASSR drive, we finally also investigated the differential effect of the microcircuit alterations on the 20 Hz beta component in response to 20 Hz drive. We simulated the same three models as in the section before, namely the IPSC+gGABA-SCZ-like, the IPSC+bInh-SCZ-like and the Full-SCZ-like models, however, with 20 Hz drive. The results are summarised in Fig. 7. Overall, we found that adding both reduced GABA levels and NMDAR hypofunction, alone or in combination, changed the response to 20 Hz drive in a similar way.

Discussion
Robust evidence for a reduction in gamma power during 40 Hz ASSR stimulation in patients with schizophrenia has emerged over the last decades 5 . Modelling studies show that there are many possible microcircuit alterations that can explain this reduction, such as as IPSC decay time increases (e.g. Vierling-Claassen et al. 17 , Metzner et al. 23 , and Metzner 22 ), reduced GABA levels (e.g. Vierling-Claassen et al. 17 , Metzner et al. 24 ), NMDAR hypofunction (e.g. Carlen et al. 16 ), dopamine deficits (e.g. Kömek et al. 25 ) and SCZ-associated genetic variants of ion channels (Metzner et al. 24 ). Therefore gamma power at 40 Hz ASSR stimulation cannot disentangle the contributions of these different microcircuit alterations. However, a few experimental 17 and modelling studies 17,23 have indicated that the failure to fully entrain to the 40 Hz stimulus train is also reflected in an increased 20 Hz component during 40 Hz stimulation. While the existence of this beta component has not been robustly demonstrated experimentally, modelling suggests that it might offer an opportunity to mechanistically link gamma entrainment deficits to microcircuit alterations 17,23 . Therefore, in the present computational modelling study, we explored the conditions leading to the emergence of such a 20 Hz component. The aim of this work was to confirm the hypothesis that a beta component during 40 Hz ASSR drive is mechanistically linked to an increase in IPSC decay times and explore the conditions under which this beta component emerges, with a focus on the strength of the stimulation input. We could demonstrate  In both plots the network model has increased IPSC decay times (from 8 to 28 ms) and the I-E and I-I synaptic strength ( g ie and g ii , respectively) is varied from 100% (black) to 10% (lightest grey) in steps of 5%.  In both plots the network model has increased IPSC decay times (from 8 to 28 ms) and now both the I-E and I-I synaptic strength ( g ie and g ii , respectively) is varied from 100% (black) to 10% (lightest grey) in steps of 10% and simultaneously the interneuron excitability b inh is varied from − 0.01 (black) to − 0.1 (darkest grey) and then in steps of − 0.05 to − 0.6 (lightest grey). However, there are several other potential factors that could contribute to these inconsistent experimental results. First, a difference in stimuli might play a role, since some studies use amplitude-modulated tones as opposed to the click-trains used in Vierling-Claassen et al. 17 . While click-trains have been shown to produce stronger and more reliable responses than other stimuli 20,26,27 , a systematic analysis of potential differences of ASSRs at stimulus frequency or at (sub-)harmonics has not been carried out so far. Nevertheless, Brenner et al. assumed that such differences are likely to exist 28 and Griskova-Bulanova et al. 29 could demonstrate that flutter In all plots the network model has increased IPSC decay times (from 8 to 28 ms). In (a) and (b) the I-E and I-I synaptic strength ( g ie and g ii , respectively) is varied from 100% (black) to 10% (lightest grey) in steps of 5%. In (c) and (d) the interneuron excitability b inh is varied from − 0.01 (black) to − 0.1 (darkest grey) and then in steps of − 0.05 to − 0.6 (lightest grey). In (e) and (f) now both the I-E and I-I synaptic strength ( g ie and g ii , respectively) is varied from 100% (black) to 10% (lightest grey) in steps of 10% and simultaneously the interneuron excitability b inh is varied from − 0.01 (black) to − 0.1 (darkest grey) and then in steps of − 0.05 to − 0.6 (lightest grey).  30 found that flutter amplitude modulated tones were more pleasant and less arousing than click-trains and, importantly, were not modulated by attention during a variety of tasks. Second, most ASSR studies are performed using EEG recordings and the lower sensitivity of EEG compared to MEG might explain why the more subtle effect of the beta component has not been detected with EEG. Furthermore, as already pointed out by Vierling-Claassen et al. 17 , averaging in time before the transformation into the frequency domain can potentially reduce the 20 Hz component considerably since the beat-skipping behaviour can vary from trial to trial. The skipped beats can be the 1st, 3rd, 5th,... in one trial while being the 2nd, 4th, 6th,.. in another trial, and thus would cancel out when averaged over time prior to the frequency transform. This cancelling effect would further exacerbate the potential problems of EEG studies detecting the effect compared to studies using MEG. Lastly, it seems likely, that the GAT1 deficits and therefore, the prolonged IPSC decay times, are heterogeneous in patients and therefore, that for some patients a beta component is detectable and for some not. Whether the existence of a beta component might be a marker for a specific stage or subgroup of the disorder is not clear and we are not aware of studies exploring the GAT1 expression over the course of the disorder or for different patient subpopulations. Exploring this in animal models of schizophrenia during the emergence of the disorder and over the course of illness would allow us to incorporate more information on the development into the model to answer these questions. Nevertheless, our modelling results offer a plausible and simple explanation of the inconsistent experimental results and suggest a careful choice of the input strengths in ASSR experiments. One direct implication of our work is that 40 Hz ASSR stimuli at lower amplitudes should be used in future studies to experimentally confirm the robustness of the beta component. However, the modelling approach could also be used in several other ways to inform new experimental designs: (a) More complex stimuli could be constructed, based on modelling results with the presented model, that are able to separate the different microcircuit alterations. In general, using the model one could potentially design a stimulus or several stimuli that would yield different results depending on the assumption alteration A is present in contrast to alteration B. To give a specific example, a stimulus that consists of a superposition of a 20 Hz and a 40 Hz oscillations is likely to induce different responses in a circuit with prolonged IPSC decay times (i.e. one that shows beta resonance) compared to one where this is not the case. (b) The model (or other models of ASSR deficits in SCZ) can be used to explore potential treatment targets (for an example of such an approach see Siekmeier and van Maanen 31,32 ). However, whether an intervention represents a potential drug (or more general treatment) target depends on which microcircuit alterations are mechanistically linked to ASSR measures. Therefore, to increase the likelihood that intervention targets found in such modelling studies really translate to clinical practice, a correspondence between alterations and measures as outlined in a) needs to be established first. c) Lastly, one could for example use the model to explore the effects the type of ASSR stimulus (click-train, amplitude modulated tones, chirps,...) has on the existence of the beta component discussed in our work and to try to identify the most suitable stimulus type, i.e. the type that produces a beta component with the highest probability, which could then inform future experimental studies to confirm the existence of a robust beta component.
We used a simple computational model consisting of an excitatory population, representing pyramidal cells, and an inhibitory population, representing PV + inhibitory interneurons. While most of the experimental evidence for a reduction in GAT1, which in turn would lead to an increase in IPSC decay times, points towards chandelier cells 33 , we have previously shown that at realistically low ratios of chandelier cells to basket cells in a microcircuit, gamma and beta range ASSR changes as seen in SCZ patients, are most likely due to an increase of IPSC decay times at basket cell synapses 34 . Our simplified model does not incorporate other types of inhibitory interneurons such as somatostatin-positive (SST + ) or vasoactive intestinal peptide-positive (VIP + ), although they have been shown to play important functional roles in cortical microcircuits 35 . Furthermore, there is recent experimental evidence of alterations to SST + interneurons in schizophrenia. Hashimoto et al. 36 found a reduced expression of GAD67 in SST + neurons, while no reduction of GAT1 was apparent. Furthermore, Morris et al. 37 observed that both the density of SST neurons and the expression of SST + per neuron was reduced in schizophrenia. These changes have been found in most cortical layers with varying strength 37 and can be observed throughout cortex 38 . This suggests, however, that IPSC decay times at SST + interneuron synapses, a necessary requirement for the emergence of a beta component in our model, should remain intact in patients with schizophrenia. Additionally, the generation and maintenance of fast cortical rhythms in the beta and gamma range has been mainly attributed to PV + neurons [39][40][41] , although SST + neurons have recently also been found to be involved 42 . These findings suggest that alterations of SST + interneurons should only play a minor role in the emergence of a beta component in gamma ASSR tasks, and they were therefore not considered in the present study. Nevertheless, an exploration of the effects of SST + alterations on cortical rhythms, especially for the lower frequency bands such as theta and alpha and for theta-gamma cross-frequency coupling, is warranted.
Beyond the question whether a beta component emerges in ASSR responses of patients with schizophrenia or not, our modelling work addresses a broader and more important issue. In general, it has proven extremely difficult to map schizophrenia-associated alterations of local microcircuits to specific neurocognitive or electrophysiological markers. Similar difficulties exist for other neuropsychiatric disorders such as autism spectrum disorder. Our work here shows that, while the robust deficit in the 40 Hz response to 40 Hz drive is not specific to a single microcircuit alteration (and neither are putative changes for the responses to 20 Hz ASSR drive), the emergence of the 20 Hz in response to 40 Hz drive is. The computational model presented here mechanistically links the microcircuit change to the electrophysiological marker, thus, demonstrating the usefulness of mechanistic computational models in advancing our understanding of the relationship between features of the microcircuitry and non-invasive biomarkers, as we have argued before 43 . Furthermore, our simulation results show that the standard 40 Hz ASSR measure is not specific enough to resolve the complex, nonlinear interactions Scientific Reports | (2021) 11:20387 | https://doi.org/10.1038/s41598-021-99793-w www.nature.com/scientificreports/ on the local circuit level and that more complex experimental designs are needed to disentangle them. This becomes especially important when considering that not only changes to the glutamatergic and GABAergic synapses considered in this work influence gamma ASSRs, but also neuromodulators such as dopamine 25 and cell-intrinsic changes to voltage-gated ionic channels 24 . This is further underpinned by the low specificity of the 40 Hz ASSR to schizophrenia, as for example similar ASSR deficits have been found in autism spectrum disorder 44 and bipolar disorder 45 .
In addition to the well established and robust deficits in the gamma band modelled in this study, lower frequency bands (delta and theta) have also been implicated in schizophrenia [46][47][48] . Including these lower frequency deficits into computational models such as the one used here and studying their interactions (e.g. theta-gamma coupling) might further constrain the models and help to disentangle the contributions of different microcircuit alterations. Here, again a special focus should lie on a comprehensive characterisation of the ASSRs, i.e. subharmonic and harmonic responses should not be neglected as they might contain important information. This might be especially important for the low frequency deficits because there is evidence of more complex ASSRs comprising several harmonics in healthy individuals 49,50 .
In summary, with this computational study we provide insights into the mechanistic generation of ASSR frequency components in schizophrenia beyond the traditional 40 Hz power at 40 Hz drive. Furthermore, we are able to explain seemingly conflicting experimental findings and suggest a more thorough and careful consideration of the effect of stimulus strength when designing ASSR experiments. Finally, we argue for a more complex and model-driven design of gamma and beta ASSR experiments in schizophrenia and for other neuropsychiatric disorders, which might be better suited to disentangle the nonlinear contributions of different microcircuit alterations found in these disorders.

Methods
The model proposed here is based on a recent reimplementation 22 of the simple model presented by Vierling-Claassen et al. 17 , which has been used in previous studies of ASSR deficits 34 , and which is integrated in the ASS-RUnit model database, a framework for automated testing of ASSR models against observations from empirical studies 51 .
Single cell model. Single cells are represented as theta neurons (see e.g. 52 for an in-depth analysis and description of the theta neuron model).
The kth neuron in a network is described by the variable θ k , which represents the neuron state, and which is governed by the following equation where b is an externally applied current, S is the total synaptic input to the cell and N(t) is a time-varying noise input. The total synaptic input to a cell in a network amounts to where n is the number of presynaptic neurons, α j controls excitation and inhibition, i.e. is +1 for excitatory synapses and −1 for inhibitory ones, g jk is the synaptic strength from cell j to cell k and s jk is the synaptic gating variable from cell j to cell k . Synaptic gating variables evolve according to where τ j is the synaptic decay time, τ R the synaptic rise time and η is a scaling parameter. A single pacemaker cell provides rhythmic ASSR drive to the network. Additionally, Poissonian noise input is also given to all cells in the network, where a noise spike at time t n elicits the following excitatory postsynaptic potential (EPSP) N(t) where Ag gmax is the noise strength, τ exc is the synaptic decay time, τ R the synaptic rise time, and H the Heaviside function.
Network. We combined 20 excitatory pyramidal cells together with 10 inhibitory cells into a network model, following the earlier work of 17,22 .
A schematic depiction of the network can be found in Fig. 8. Populations connect to each other and also to themselves. The connectivity between any two populations is all-to-all. All populations also have two sources of input, the oscillatory drive input and a background noise input. The drive input periodically sends spikes with a given drive frequency to all populations, mimicking the rhythmic ASSR input. An overview of the model parameters can be found in Table 1.
To evaluate the oscillatory entrainment we calculate simulated MEG signals by summing all excitatory synaptic variables over all pyramidal cells (as in 17,22,34 ). As the main measures for entrainment we perform a Fourier transform on these 'MEG' signals and extract the power at 40 Hz and at 20 Hz.
α j g jk s jk www.nature.com/scientificreports/ Implementation of schizophrenogenic microcircuit alterations. We implemented changes to the GABAergic and glutamatergic neurotransmitter systems that have been associated with schizophrenia.
GABAergic system. 50% of PV + neurons in the dorsolateral prefrontal cortex lack detectable levels of GAD67 12 .
It has been suggested that the reduced expression of GAD67 mRNA likely implies a reduction in GABA synthesis in cortical GABAergic neurons, which in turn would lead to smaller amplitude IPSCs at the postsynaptic site 8 . We implemented this change as a reduction of the weight of inhibitory connections (for both I-E and I-I connections). Furthermore, a reduction in the plasma membrane GABA transporter GAT1 14 has also been found in PV + interneurons in SCZ patients 13 . GAT1 is a major contributor to the specificity of synapses by preventing spillover   14 and a reduction in GAT1 results in prolonged IPSC durations 14 . Here, we realised this change as an increase of the IPSC decay time constant (as in previous studies 17,22 ).
Glutamatergic system. NMDAR antagonists, such as phencyclidine and ketamine, produce symptoms, which are very similar to key clinical features of SCZ 15 . Convergent lines of evidence underpin that NMDARs are dysfunctional in SCZ 10 . Examples of direct evidence in favour of this hypothesis are changes in NMDAR-associated protein levels at the postsynaptic site 53 , a reduction in NMDAR-mediated signalling following neuregulin 1 activation of ErbB4 receptors 54 , lower levels of glutathione (a modulator at the redox-sensitive site of NMDARs) 55 and a reduction of kynurenine 3-monooxygenase that might increase kynurenic acid (an NMDAR antagonist) levels 56 . Indirect evidence, such as findings that putative risk genes for SCZ can affect NMDAR function 57 and that substances enhancing NMDARs might reduce symptom severity in SCZ 58 , further underpin this idea. A potential hypofunction of NMDARs would lead to lower levels of excitation of PV + interneurons and we therefore implemented this alterations by decreasing the applied current b inh to inhibitory cells.
In this study, we considered four different 'SCZ-like' networks, which comprised the following combinations of changes to the GABAergic and glutamatergic system described above: • IPSC-SCZ-like network: For this 'SCZ-like' network, we only implemented the increase of the IPSC decay time constant (as in previous studies 17,22 ). • IPSC+gGABA-SCZ-like network: Here, additionally to the increase of IPSC decay times, we also reduced the weight of the inhibitory GABAergic connections (as in other previous studies 23 ). • IPSC+bInh-SCZ-like network: For this network, we decreased the applied current to the inhibitory cells together with the increase in IPSC decay times. • Full-SCZ-like network: Here, we combined all three alterations, i.e. we implemented an increase in IPSC decay times, a decrease in inhibitory weights and a decrease in applied currents to inhibitory cells.
Replication of previous findings. In a first step, we validated the control and the 'IPSC-SCZ-like' model against experimental data and previous modelling studies. Here, we stimulated the network, either in the control or the 'IPSC-SCZ-like' condition, with ASSR stimuli at drive frequencies of 20, 30, and 40 Hz. We then analysed the response by plotting the resulting simulated MEG signal and by transforming it into the frequency domain via a Fourier transform (see Figs. 1 and 2). We compared the responses to experimental studies 6,17 and previous modelling studies 17,22 . Input strength dependence of the 20 Hz component. Since ASSR deficits are most pronounced in the gamma band in schizophrenia 5 , we focused on 40 Hz stimulation. Specifically, we explored how 40 Hz and, importantly, 20 Hz power depend on the input strength for 40 Hz drive. Therefore, we stimulated the 'IPSC-SCZ-like' for input strengths varying from 0.1 to 1.5 times the default value from the section before. Again, we analysed the simulated MEG signal and the resulting power spectra with a focus on the 40 Hz and 20 Hz power (see Fig. 3).

Combinations of alterations and their input-strength-dependence.
Lastly, we investigated how the input-strength-dependence of the beta component of the 40 Hz ASSR changes when more than one microcircuit alteration were present. We thus explored the 'IPSC + gGABA-SCZ-like' , 'IPSC+bInh-SCZ-like' and 'Full-SCZ-like' model (see above), which included reduced GABA levels, NMDAR hypofunction, or both, in addition to the other SCZ-associated microcircuit alteration that was already present in the IPSC-like-SCZ model. We then stimulated these models with 40 Hz ASSR drive with varying input strengths as in the section before and analysed the simulated MEG signal and its 40 and 20 Hz power (see Figs. 4, 5 and 6).
Implementation details and code availability. The model was implemented using Python 2.7.9 and numpy 1.9.3. Analysis and visualisation of the model output was also done in Python using the numpy and matplotlib packages (matplotlib 1.4.3). Model equations were numerically solved using a simple forward Euler scheme. A single model run simulated a 500 ms trial and the time step was chosen such that this resulted in 2 13 = 8192 data points. The model output was unaffected by using a smaller time step.
Simulation results varied from trial to trial because of the stochastic nature of the background input. Therefore, we always performed 20 simulation trials, each with a different realisation of the noise process. We then averaged these trials in time to get an average simulated MEG signal and all analyses were based on this average signal.