Introduction

Stress is an important contributor to the onset, maintenance, and relapse of major depressive disorder (MDD), with estimates suggesting that up to 80% of first major depressive episodes (MDEs) are preceded by major life events [1]. Evaluating the neurobiological effects of stress is therefore critical for understanding the pathophysiology of MDD. Females are at twice the risk for developing MDD as males; thus identifying sex-dependent neural responses to negative stress in MDD is important [2, 3]. Acute stress recruits bottom-up and/or top-down regulatory brain pathways [4]. Specifically, the brainstem senses perturbations to homeostatic equilibrium and directly activates the hypothalamus–pituitary–adrenal (HPA) axis and autonomic nervous system through the hypothalamus, resulting in the secretion of corticotropin from the pituitary gland and subsequently cortisol from the adrenal glands. In healthy individuals, the limbic system (e.g., amygdala, hippocampus) is initially deactivated in response to stress [5, 6], which is correlated with cortisol secretion and thought to index a defense response. The prefrontal cortex is a key component of this response, providing top-down control to enable limbic deactivation [7]. This adaptive pattern is disrupted in MDD, with exaggerated limbic responses to stressors and reduced top-down control from frontal regions [8,9,10]. These neural effects are linked to HPA axis reactivity being altered in MDD [10], with studies reporting acute hyper- [11] and more chronically hypo-secretion of cortisol in response to stress, particularly in females (for a review, see ref. [12]). Such interplay between the HPA axis, cortical and limbic regions in stress circuitry call for a network approach to investigations of acute stress. Independent component analysis (ICA) is a data-driven technique that identifies intrinsically coupled functional networks and allows subsequent multivariate approaches [13] and network modeling [14]. However, the many networks uncovered by this data-driven approach raises potential multiple comparison issues.

To guide our thinking around the networks associated with stress and psychopathology we relied on Menon’s triple network model [15], which identifies networks important for understanding cognitive dysfunction across psychiatric disorders [16, 17]. These include the frontoparietal control network (FPN), important for cognitive control and adaptive behavioral regulation, the default mode network (DMN), implicated in self-referential mentation, and the salience network (SN), critical in the detection and mapping of salient external and internal inputs. The FPN and SN typically increase during stimulus-driven cognitive and affective information processing, while the DMN de-activates during task engagement. In terms of responses to acute stress, connectivity between the FPN↔SN has been highlighted as crucial to shifting between vigilance and recovery [18]. In MDD, deficits have manifested as increased DMN activation [16, 17, 19], associated with increased rumination [20], and decreased ability to deactivate the DMN during cognitive or affective tasks. Reduced FPN activation is also observed in MDD [21, 22], suggesting diminished top-down control over limbic regions and aberrant SN mapping, which ultimately impacts switching between the FPN and DMN and is associated with symptom severity in MDD [23]. In line with this, evidence points to decreased connectivity between the FPN↔DMN and increased connectivity between the DMN↔SN [23, 24] in MDD. However, a recent meta-analysis [17] suggests increased connectivity between the FPN↔DMN. Therefore, further investigations are needed.

Studies of acute negative stress responses in healthy controls (HC) implicate regions such as the anterior cingulate cortex (ACC), basal ganglia, hippocampus, insula, amygdala, and frontal regions [6]. We previously demonstrated that these stress responses were significantly associated with steroid hormone responses in females [25]. A network modeling approach found increased SN amplitude [26] and decreased SN↔DMN connectivity after stress [19]. Moreover, among cortisol responders, increases in amygdala connectivity with the DMN, mPFC, OFC, and hippocampus emerged [27, 28]. In addition, under acute stress, decreases in amygdala connectivity with the dorsolateral prefrontal cortex (dlPFC) and ACC have also been reported in HC [29]. These findings suggest that stress increases connectivity within limbic regions and decreases connectivity between limbic and prefrontal regions, resulting in decreased top-down control. Critically, in individuals with current or remitted MDD (rMDD), acute stress delivered through the Montreal Imaging Stress Test (MIST) reduced activation in the vmPFC [30]; moreover, among those with rMDD, stress-related activation of the dlPFC and striatum, as well as hyperconnectivity between the striatum and amygdala, emerged [30, 31]. Collectively, these findings highlight interactions between limbic and cortical regions subserving the experience of, and reaction to, acute stress that are impaired in a partially trait-like fashion in MDD. Given this circuitry is highly sexually dimorphic and females have higher rates of MDD, an understanding of the neural responses within females will contribute importantly to the overall understanding of MDD [32].

Growing evidence indicates that network interactions are affected by neuromodulators, including γ-aminobutyric acid (GABA)—the main inhibitory neurotransmitter [33, 34]. Reduced GABA levels have been observed in cerebrospinal fluid [35] and cortical brain tissues [36] of patients with MDD. Magnetic resonance spectroscopy (MRS) quantifies in vivo GABA in brain regions of interest non-invasively, with findings indicating reduced prefrontal and ACC GABAergic transmission in acute stress [37] and in MDD [38, 39], which can be reversed with treatment [40,41,42]. Associations between GABA and large-scale network connectivity have been found in healthy samples [43, 44] and in post-partum depression [45]. However, relationships between GABA and large-scale network connectivity in the response to stress in depression are poorly understood.

In sum, the present study takes a multi-pronged approach combining multimodal neuroimaging (GABA MRS and functional magnetic resonance imaging (fMRI)) with network modeling to assess the effects of acute stress in females with MDD and, to examine putative trait characteristics in rMDD. We hypothesized that in this female sample MDD would be associated with a blunted cortisol response to stress (although this can be dependent on chronicity), reduced rACC and/or left dlPFC GABA levels, and MDD and rMDD would be associated with aberrant changes in network amplitude and between-network connectivity in the FPN, SN, and DMN. We also predicted that those with rMDD would show greater increases in top-down control signals in the FPN compared with MDD alone. Further, we hypothesized that lower rACC and dlPFC GABA would be associated with deficits in engaging and disengaging these networks under stress. Finally, we anticipated that GABA in regions overlapping with the networks of interest (rACC in SN and dlPFC in FPN) would moderate network amplitude/connectivity through increased inhibition in these regions.

Materials, patients, and methods

Participants

Fifty-one unmedicated females (aged 18–25) were recruited from the community (see Supplemental Table S1). They provided written informed consent to a protocol approved by the Partners Human Research Committee. Participants were assessed by a clinician (at two study sites with high inter-rater reliability, see Supplemental Information) using the Structured Clinical Interview for the DSM-5 (SCID-5; [46]). Nineteen participants were assessed as having current MDE (MDD group), 15 were fully remitted from at least one prior MDE (rMDD), and 17 were HC. Exclusion criteria included other comorbid psychiatric disorders (except for simple phobia, social anxiety, and generalized anxiety disorder if secondary to MDD), use of psychoactive drugs, recent recreational drug use, past substance use disorder, and more than five alcohol-related blackouts.

Procedure

The imaging session took place in the early follicular phase and in the afternoon (to control for diurnal variability of cortisol response) [47]. After IV insertion, participants entered the 3 T scanner and underwent a single ~45-min GABA MRS scan before completing an acute laboratory stressor task during fMRI. Based on extensive piloting before study onset, to improve the potency of the stressor following piloting, this protocol combined the MIST [48] and the Maastricht Acute Stress Test (MAST) [49] into a single hybrid stressor and compared measurements to baseline. As depicted in Fig. 1, participants carried out four blocks of arithmetic problems, each lasting ~3.5 min. Participants used a three-button touchpad to navigate left or right to the correct digit (0 through 9) on a rotary dial and submitted their answer with a top button press. During the first block (pre-MAST), no time pressure was imposed, and the participants received trial-by-trial feedback of their performance (“correct” and “incorrect”), constituting a “no-stress” baseline condition. After the first block, the scanner table was brought out and the participant (whilst still lying on the scanner table) was asked to complete a 12-min MAST protocol: two experimenters acting as “doctors” (whom the participant had not met yet) entered the scanner suite and gave instructions for the MAST task, which involved interleaving blocks of mental arithmetic (counting backward from a four-digit number out loud in steps of 17) and immersing their hand in ice-cold (0–2 °C) water. The length of the blocks was determined by a computer algorithm (thus, introducing unpredictability and uncontrollability) and the “doctors” provided an evaluation of the arithmetic task for a socio-evaluative component. After the MAST, the “doctors” informed the participant that they would continue monitoring their performance from the scanner suite. The table was returned to the scanner and the participant completed the second block of the MIST (post-MAST 1). In this block the problems were also untimed, providing a direct comparison to the pre-MAST block. In the third block of the MIST (post-MAST 2), stress was increased. Specifically, time pressure was imposed on each problem (calibrated by participants’ no-stress responses), and performance monitoring compared to the purported average was shown on a mock performance bar. After the MAST 2 block, one “doctor” gave negative verbal feedback over the intercom, saying the performance was well below average and that the participant would need to improve in this final block to make the data useable. Finally, the participant completed another block (post-MAST 3), identical to post-MAST 2. Self-report affective ratings were collected pre- and post-stress.

Fig. 1: Combined stressor.
figure 1

Initial BOLD acquisition under baseline condition (pre-MAST, untimed problems), followed by an acute stressor (MAST [49]), followed by subsequent BOLD acquisition under stress conditions (post-MAST1—untimed problems, post-MAST2—timed problems and progress bar), followed by negative verbal feedback over scanner intercom, followed by final BOLD acquisition (post-MAST3, timed problems with progress bar).

Full size image

MRI data acquisition and preprocessing

A 3 T Siemens MAGNETOM Prisma scanner (Siemens Medical Systems, Iselin, NJ) equipped with a 64-channel head coil was used to acquire high-resolution functional and structural MRI data (see Supplement). Functional MRI data were pre-processed using fMRIPrep 1.5.8 ([50]; RRID:SCR_016216), which is based on Nipype 1.4.1 ([51]; RRID:SCR_002502).

GABA MRS acquisition and processing

The T1-weighted structural images were used to place two independent voxels in the bilateral rACC (17.5 ml; 35 × 20 × 25 mm3) and left dlPFC (18.75 ml; 25 × 30 × 25 mm3) for MRS data collection (see Supplement and Supplemental Figs. S1 and S2) [52]. Single measurement GABA+ concentrations are reported as GABA+/water.

Blood cortisol collection and analysis

We previously published [11] method for acquiring HPA-axis hormone level changes in response to acute negative stress “in real-time” using serial blood samples was applied (see Supplement). Serum cortisol changes from stress were quantified using the area under the curve with respect to ground (AUCg) and with respect to increasing (AUCi) calculations [53] to take into account all measurements and variability in the timing of cortisol peak following stress [54].

ICA of functional data

For a data-driven, circuit-level evaluation of stress-related neural changes, we ran a group ICA of the task fMRI data from all runs, and all participants temporally concatenated together using MELODIC in FSL. The model order (number of components estimated by ICA) was set between 10 and 40 for fine-tuning to estimate the model order that provided good separation between networks and between networks and noise. The resulting independent component maps were thresholded using Gaussian-gamma mixture modeling with a p value of (p = 0.05) to identify the signals in each component. Task activation networks were identified that correspond [55, 56] to resting-state networks highlighted in the triple network model. A dual regression approach [13] was used to extract network timecourses for each run for each participant that was then used to estimate the connectivity between each network pair and the amplitude (or task activation) of each network.

Model comparison approach

We tested for effects of stress and group using FSLNETS v0.6 [57] to implement general linear models with non-parametric permutation testing for estimation and inference (see Supplement). Our main outcome variables were partial correlations between network time-series (network pairs) and individual network amplitudes. GABA and cortisol AUC outliers (>1.5× outside the interquartile range) were first removed (cortisol: 2× HC, 1× MDD, dlPFC GABA: 1× MDD) and linear regression evaluated group differences using the lm function in R [58]. For repeated measures (main outcome network variables, self-report), mixed-effects regressions were conducted using the lmer4 package in R [58], a method that efficiently handles non-independence of repeated measures. We first specified an intercept-only null model to estimate the total systematic variance in the outcome variable where the outcome was regressed on a random effect variable for participants (random intercept). Next, we specified a three-step mixed-effects model with REML estimation. Firstly, we added the main predictor—stress (4 levels/timepoints: pre-MAST, post-MAST1, post-MAST2, and post-MAST3) and secondly—group (3 levels: HC, MDD, and rMDD). Thirdly, we added the two-way interaction of stress and group. Model selection was identified using likelihood ratio tests [59]. We compared the null model (only random effects) with three models:

Null model: value ~ (1|participant)

Model 1 (main effect of stress): value ~ stress + (1|participant)

Model 2 (additive model): value ~ stress + group + (1|participant)

Model 3 (interactive model): value ~ stress + group + stress * group + (1|participant)

95% confidence intervals are reported. Power to achieve an R2 of 0.3 using a random-effects model with two predictors (group and stress) was calculated using G*Power 3.1. For a one-tailed test, the required sample was 44 and for a two-tailed test, the required sample size was 50.

To leverage our multimodal data but also deal with possibly missing values more efficiently, mixed-effects regression in R was used to ascertain the effects of stress, group, GABA, cortisol AUC, and affective ratings on our main outcome variables of between-network covariance and network amplitude [58]. Specifically, we regressed connectivity and amplitude values on group, stress, rACC GABA, dlPFC GABA, cortisol, and affective ratings with a random effect for participants (intercept). Interactions between other predictors and stress were tested by comparing models with and without each interaction term. As the full sample did not have all of the measures (Supplemental Table S2), model selection was carried out using Akaike Information Criterion (AIC), a difference of ten between AICs being necessary for a model to be considered significantly improved [60].

Results

Network maps

Visual inspection of the group ICA spatial maps in consultation with a neuroimaging biostatistician (L.N.) resulted in the selection of the model order 40 group ICA results for analyses. Based on Menon’s triple network model [15] and additional regions of interest for stress tasks, a priori task networks of interest were as follows: (1) right FPN (Supplemental Fig. S3A); (2) ventromedial prefrontal, striatal, ACC network (vmPFC–Str–ACC; Fig. S3B); (3) SN (Fig. S3C); (4) a network including temporal regions, the insula, and amygdala (Temp–Ins–Amyg; Fig. S3D); and (5) the DMN (Fig. S3E). A maximum of five networks was selected to reduce multiple comparisons.

Stress results

A one-sample t test (against zero) on the AUCg and AUCi revealed that cortisol concentrations were significantly increased across all participants (AUCg: t(33) = 15.96, p < 0.001; AUCi: t(34) = 2.60, p < 0.007) indicating that the stressor elicited the intended effect. AUCg and AUCi were then regressed on a between-subjects factor of the group (HC, MDD, and rMDD). This model significantly predicted AUCg (F(2, 31) = 3.73, p = 0.04, R2 = 0.19). Follow-up tests showed that MDD vs. rMDDs had lower AUCg (t(31) = 2.71, p = 0.03) (Supplemental Table S3 and Supplemental Fig. S4), but HC did not differ from the others (all p > 0.3). Similarly, for all affective rating measures, significant main effects of stress and MDD group emerged (Supplemental Fig. S5), indicating that the stressor increased negative affect and decreased positive affect across participants as expected [49] and that the MDD group had increased negative affect compared to the HC group.

GABA MRS results

GABA MRS values from the rACC and dlPFC were regressed on the group (HC, MDD, and rMDD). The model significantly predicted rACC GABA (F(2, 40) = 3.23, p = 0.05, R2 = 0.14). Follow-up tests showed that, as hypothesized, the MDD group had significantly lower levels of rACC GABA than HC (t(40) = 2.49, p = 0.04) (Fig. 2), whereas rMDDs did not (t(40) = 1.74, p = 0.20). A post hoc analysis showed that gray matter volume in the rACC did not differ between the HC and MDD groups (t(40) = 0.22, p = 0.97), although the rMDD group showed trends for lower gray matter volume in the rACC than the HCs (t(40) = 2.31, p = 0.07) and the MDDs (t(40) = 2.23, p = 0.08).

Fig. 2: Decreased rostral Anterior Cingulate cortex GABA in current major depressive disorder vs. healthy controls.
figure 2

The dark line inside the box represents the median. The top of the box is 75th percentile and the bottom of the box is the 25th percentile. The endpoints of the lines (aka whiskers) are at a distance of 1.5× interquartile range (the distance between 25th and 75th percentiles). HC (N = 13) mean (M) = 1.26, standard deviation (SD) = 0.20; MDD (N = 17) M = 1.10, SD = 0.16; rMDD (N = 14) M = 1.14, SD = 0.19.

Full size image

FMRI results: effects of stress and diagnosis

Network amplitude in the FPN was best explained by a model including stress + group (model 2 above), indicating a decrease in activation in this network from stress at post-MAST2 (β = −0.13, 95% CI: [−0.24 to −0.02]) and post-MAST3 (β = −0.22, 95% CI: [−0.33 to −0.11]), and an overall decrease in network amplitude in MDD (β = −0.26, 95% CI: [−0.12 to −0.040]) and rMDD (β = −0.23, 95% CI: [−0.37 to −0.08]) compared to HC (Fig. 3). Similarly, network amplitude in the SN was best explained by the stress + group model (model 2), indicating a decrease in activation for the MDD (β = −0.28, 95% CI: [−0.046 to −0.09]) and rMDD (β = −0.19, 95% CI: [−0.38 to −0.00]) groups compared to HC but the individual predictors were not significant.

Network amplitude in the vmPFC–Str–ACC network was best explained by a model including stress only (model 1), owing to increases in activation at post-MAST3 (β = 0.44, 95% CI: [0.23 to 0.65]). Network amplitude in the DMN was best explained by a model including stress only (model 1), with increases in activation at post-MAST1 (β = 0.08, 95% CI: [0.00 to 0.16]). There were no improvements on the null model for the Temp–Ins–Amyg network or any of the between-network connectivity pairs and the interactive model (model 3 above) did not fit any of the amplitude or connectivity data better than the stress + group model (model 2), suggesting either that neural response to stress did not vary systematically across groups or this study was underpowered to identify a significant interaction effect. The optimal base model parameters are shown in Table 1.

Fig. 3: Stress and diagnosis affect network amplitude.
figure 3

A Frontoparietal network amplitude decreased by stress. B Frontoparietal network amplitude decreased in MDD and rMDD. C Salience network amplitude decreased in MDD and rMDD. D vmPFC–Str–ACC network amplitude increased by stress. HC: N = 17; MDD: N = 18; rMDD: N = 14. The dark line inside the box represents the median. The top of the box is 75th percentile and the bottom of the box is the 25th percentile. The endpoints of the lines (i.e., whiskers) are at a distance of 1.5× interquartile range (the distance between 25th and 75th percentiles).

Full size image
Table 1 Effects of group and stress on network amplitude response to stress in MDD and rMDD.
Full size table

Association with affective ratings, cortisol, and GABA

To leverage the multimodal data, we compared extended mixed-effects models to examine the relationship between network amplitude/connectivity and the additional predictors rACC/dlPFC GABA, cortisol AUC, and affective ratings. Building on the initial group/stress analyses, the optimum base model (null model/model 1/model 2) was augmented with the additional predictors in a stepwise manner, and model AICs were compared. Models including AUCG and AUCI were directly compared using the same criteria. FPN network amplitude was best explained by a model that included rACC GABA, dlPFC GABA, and cortisol AUCG (R2 = 0.32), whereas a model including only rACC GABA best-explained vmPFC–Str–ACC activation (R2 = 0.36). For SN, the model included rACC GABA, dlPFC GABA, Cortisol AUCG, state anxiety change, and negative affect change (R2 = 0.36). Extended models for the Temp–Ins–Amyg and DMN networks did not improve on the optimum base model (Supplemental Table S4).

Dimensional analyses of GABA and network connectivity/amplitude

Given the significant decrease in rACC GABA in the MDD group compared to HC, further analyses examined relationships between rACC GABA and network connectivity or amplitude. We adopted a dimensional approach, regressing connectivity/amplitude values on fixed effects of stress (pre-MAST, post-MAST1, post-MAST2, and post-MAST3), rACC GABA, and their interaction, with a per-participant random adjustment to the fixed intercept. Rostral ACC GABA was positively associated with between-network connectivity in two network pairs across all conditions; FPN↔Temp–Ins–Amyg (β = 2.507, 95% CI: [0.256 to 4.757]) and the Temp–Ins–Amyg↔DMN (β = 2.226, 95% CI: [0.141 to 4.311]). There were no interactions between rACC GABA and stress. Crucially, the rACC MRS voxel overlapped the Temp–Ins–Amyg network (Fig. 4).

Fig. 4: Multimodal intersection of fMRI and MRS data.
figure 4

Overlap between Temp–Ins–Amyg network and approximate location of rostral anterior cingulate MRS voxel (35 × 20 × 25 mm3).

Full size image

Discussion

Acute negative stress significantly increased cortisol across all females, although females with MDD had lower cortisol AUCG than those with rMDD. In addition, females with MDD were characterized by decreased rACC GABA compared to healthy females. Three networks of interest were affected by acute stress across all groups. Specifically, the FPN network was characterized by decreased activation during the MIST task under acute stress, whereas the vmPFC–Str–ACC network and DMN networks showed increased stress-related activation. Irrespective of stress, females with MDD and rMDD were characterized by decreased FPN and SN activation compared to HCs. Critically, levels of GABA in the rACC were positively associated with activation in the SN and vmPFC-Str-ACC and negatively associated with activation in the FPN. Finally, dimensional analyses showed that rACC GABA was associated with increased connectivity between limbic regions (Temp–Ins–Amyg) and the FPN and DMN.

Together, these findings point to a key role for rACC GABA in downstream network activation and cortico-limbic between-network connectivity under negative stress. Rostral ACC GABA was reduced in MDD while females with MDD and rMDD showed decreased FPN and SN activation. This suggests that a deficit in rACC GABA may reflect impaired inhibition, which may be associated with deficient top-down control and aberrant salience mapping. Similar to prior findings in females with MDD [61], we observed a blunted cortisol response to stress in MDD. This finding may also be related to the number of depressive episodes or young age of the females, given that we previously reported differences in negative stress circuitry responses and hormonal physiology dependent on the age of females with MDD [2, 11]. Collectively, we speculate that these patterns highlight blunted stress response in MDD which may lead to decreased SN activation, thereby lowering demand on the FPN to exert top-down control. This was further supported by dimensional analyses showing that rACC GABA was associated with increases in connectivity between a limbic network (Temp–Ins–Amyg) and the FPN and DMN.

Healthy stress response circuitry involves salience mapping (SN) and top-down control over limbic regions (Temp–Ins–Amyg) by increasing connectivity with the FPN. Current findings suggest that the females with MDD and rMDD had reduced FPN and SN activation overall, highlighting impaired top-down control and compromised salience mapping. Findings showed increased stress-related amplitude across all groups in a network including vmPFC, striatum, and ACC. The vmPFC coordinates behavioral and physiological stress responses across multiple temporal and contextual domains [62] and is critical to active coping in acute stress response, which can become blunted with repeated life stress [63]. This is reinforced by prior findings in MDD and rMDD [30] showing reduced vmPFC activation to negative stress, although we did not see vmPFC group differences.

GABA MRS data add a new dimension to the role of rACC inhibition in stress, with the addition of rACC GABA improving predictive models of amplitude in the FPN, SN, and vmPFC–Str–ACC. Taking a dimensional approach, higher rACC GABA was associated with higher connectivity between the Temp–Ins–Amyg network (which overlaps the MRS voxel) and two networks (FPN and DMN) implicated in MDD. The Temp–Ins–Amyg network contains regions from the cingulo-opercular network, which is characterized as a critical hub with pervasive other network connectivity [15]. Decreased GABAergic inhibition in the rACC in MDD could result in impaired communication between regions involved in the stress response which manifests as reduced connectivity of this network including limbic regions (Temp–Ins–Amyg) with the FPN and DMN. We suggest that stress-induced activation in the Temp–Ins–Amyg is critically affected by its connectivity with other regions, namely the FPN and the DMN, and that rACC GABA may reduce the influence of the limbic system on this critical stress hub. A recent trial of the positive allosteric modulator of GABAA receptors has shown promising results in the treatment of MDD [64]. The additional complexity of the model that should be noted is that the measure of total tissue GABA does not necessarily reflect inhibitory function. Finally, the lack of findings for dlPFC GABA could be explained through its indirect role in regulating the stress response via top-down cognitive control, rather than the more direct role of the rACC in inhibiting arousal.

Limitations of the study include potential stress from IV insertion, a lack of behavioral measures or group differences in between-network connectivity, participants taking contraception, modest sample size, and shorter duration blocks used for this task-based fMRI study than are typically used in resting-state designs. However, our repeated measures design, the inclusion of only females, and controlling of menstrual cycle and time of day likely improved our statistical power. All our reported effects were the main effects of stress, group, or rACC GABA. There were no interactions of stress X group, suggesting that the clinical groups did not respond to stress differently to the HCs. There are also potential networks that we did not examine, due to the need to reduce multiple comparisons. Even so, the SN model, group differences in rACC GABA and their association with network connectivity did not survive correction for a number of regions/networks examined. A stepwise model comparison approach was taken with the aim of parsimony and inclusion of all data points. However, this approach is less hypothesis-driven than more traditional approaches. Nevertheless, brain activation, connectivity, and physiology differences between groups provide insights into understanding differential strategies to maintain similar behavioral responses. Trauma was not examined in this study and could be a key factor in stress sensitivity. Moreover, the young age of the sample, lack of comorbidities, and limited ecological validity of acute laboratory stressors limit generalizability. Finally, future studies could leverage advances in functional MRS [65] to examine how changes in GABA unfold over the course of a stressor.

In sum, we found that stress increased amplitude in a network including the vmPFC, striatum, and ACC, and in the DMN, and reduced amplitude in the FPN. In addition, compared to HC, females with MDD and rMDD showed decreased network amplitude in the SN and FPN, implicating impaired trait-like top-down control and salience mapping. Importantly, rostral ACC GABA was associated with network amplitude in the SN, FPN, and vmPFC-Str-ACC, and connectivity between an overlapping limbic network (Temp–Ins–Amyg) and two circuits (FPN and DMN) critical to MDD pathophysiology. Together, these novel findings suggest that inhibitory GABAergic mechanisms have downstream effects on activation of, and connectivity between, circuits implicated in negative stress and pathophysiology in MDD.