Central amygdala circuit dynamics underlying the benzodiazepine anxiolytic effect

Benzodiazepines (BZDs) have been a standard treatment for anxiety disorders for decades, but the neuronal circuit interactions mediating their anxiolytic effect remain largely unknown. Here, we find that systemic BZDs modulate central amygdala (CEA) microcircuit activity to gate amygdala output. Combining connectome data with immediate early gene (IEG) activation maps, we identified the CEA as a primary site for diazepam (DZP) anxiolytic action. Deep brain calcium imaging revealed that brain-wide DZP interactions shifted neuronal activity in CEA microcircuits. Chemogenetic silencing showed that PKCδ+/SST− neurons in the lateral CEA (CEAl) are necessary and sufficient to induce the DZP anxiolytic effect. We propose that BZDs block the relay of aversive signals through the CEA, in part by local binding to CEAl SST+/PKCδ− neurons and reshaping intra-CEA circuit dynamics. This work delineates a strategy to identify biomedically relevant circuit interactions of clinical drugs and highlights the critical role for CEA circuitry in the pathophysiology of anxiety.


Subject history
For the c-Fos screen, animals were exposed to the appropriate anxiogenic stimuli and euthanized. Their brains were prepared for histological assays. Home cage animals were taken from their home cage and euthanized directly. To test the involvement of PKC + /SSTneurons in the central lateral amygdala, animals underwent surgery for virus injection, later tested first in the EPM and then in a reward-fear conditioning paradigm. After the last test, animals were euthanized and virus expression their brains was histologically assessed. For testing the interaction of DZP and PVT-CEAl projection using fMRI, animals underwent surgery for virus injection and fiber implantation and were then subjected to the fMRI scans. After the last test, animals were euthanized and virus expression in their brains was histologically assessed. For deep brain calcium imaging experiments, mice underwent surgery for virus injection and microendoscopic fiber implantation and were later fear conditioned and re-exposed to the aversive conditioning context. For neuronal population sequencing and electrophysiological recordings, mice were taken from their home cages and euthanized directly.
When injections were finished, the needle was kept in place for an additional 5 min to minimize leakage when the needle was pulled out. The scalp was then closed with 2 to 3 stitches, and the animal placed into a clean cage by itself, where its body temperature was kept constant using a heating pad. Once the animal was awake, it was transferred back to its home cage, where it received enrofloxacin (100mg/ml, Baytril®, Bayer Austria) and carprofen (Rimadyl®, 50mg/ml; Pfizer Austria) via drinking water for 14 days. Animals were allowed to recuperate for at least 5 weeks before the first experiments to allow for optimal virus expression. For detailed virus information see Supplementary Table S1.

Pharmacology
Clozapine-n-oxide (CNO; Sigma-Aldrich Co. LLC) was dissolved in 1xPBS and used on the same day for intraperitoneal (i.p.) injection at 5 mg/kg. This dose was chosen to ensure effective silencing. CNO is metabolized to clozapine in vivo, potentially resulting in unwanted psychoactive effects mediated by clozapine.
However, the clozapine level back-converted from CNO is quite low (a CNO:clozapine ratio of about 20:1 (3)), and CNO doses of up to 10 mg/kg do not affect locomotion behavior in the timeframe of the experiment (3). In our study, CNO treatment did not alter the total distance travelled or the time in the open arms of the EPM when compared to controls -the most critical behavioral parameters in the context of this study (cf. Fig. 4c,f, saline GFP group and Supplementary Fig. S1a,d, saline group). However, to account for potential side effects, we chose an experimental design in which all groups received an equal dose of CNO. Diazepam (Gewacalm®, Nycomed Austria GmbH) was diluted in 1xPBS just before i.p. injection at 1 mg/kg. The injected volume was 5 ml/kg for each substance. Behavioral tests were performed 30 min after drug injection.

Anxiety tests
A standard elevated plus maze (EPM) assay was used. The apparatus was raised 40 cm from the ground and For the c-Fos screen, 15 min after injection animals were placed in a novel chamber for 10 min where only the appropriate cohorts received 10 foot-shocks, 0.5mA, of 1s at randomized intervals of 20-100s. Mice were then transferred back to their home cage where they either stayed for 5 min before they were exposed to the EPM, or for 95 min until euthanasia.

Context fear conditioning
Fear conditioning occurred in 4 identical experimental chambers (16.5 cm wide x 16.5 cm deep x 30.5 cm high, H10-11M-TC, Coulbourn Instruments, Whitehall, PA, USA) encased in a sound-attenuating shell. Above the chamber, a custom-made house light provided illumination (around 10 lux), an infrared spotlight (Kemo Electronic, Geestland, Germany) improved animal detection, a speaker (Audiocomm, Vienna, Austria) provided sounds designed using Audacity software (http://www.audacityteam.org) with a maximal sampling frequency of 192 kHz and played from a Terratec sound card (Alsdorf / Germany). A video camera (Basler, Ahrensburg, Germany) monitored the animal's behavior. Shocks were delivered via a stainless-steel shock floor (H10-11M-TC-SF, Coulbourn Instruments), and chambers were cleaned with diluted, lemon-scented cleaning solution before each mouse. The house lights were turned off. The fear conditioning session started with a 2 min baseline, followed by conditioned stimulus (CS) presentation (3 kHz tone, total duration of 10 s, 75 dB) which was again immediately followed by a foot shock (1 s, 0.5 mA) delivered to the floor via an external shocker (H13-15, Coulbourn Instruments). The session lasted 10.5 min and consisted of 5 trials with a randomized ITI of 100 ± 30 s.
Freezing was scored with Ethovision XT 8 (Noldus Information Technology, Wageningen, the Netherlands) and defined as a lack of all movement except for respiratory-related movements (4). To control the success of fear conditioning, freezing to each aversive-CS presentation was compared to a 2-min baseline period. Videos were recorded at a rate of 20 frames/s and analyzed with Ethovision using an immobility threshold of 2.5% pixel change per sample and sampling of 0.5 s (determined by eye, blind to experimental conditions). For calcium imaging, mice were re-exposed to the conditioning context on two separate days. Before each session, mice received i.p. injections of either DZP (1 mg/kg, day 8) or saline (day 10). When comparing baseline freezing levels from day 1 (before conditioning) to freezing in the fear context (day 10), mice showed increased freezing to the context (mean ± SEM: baseline = 3.17±2.06%, context = 23.48±5.97%; n = 11 mice; paired ttest: t = 3.448, df = 10, P = 0.03).

Histology and Immunohistochemistry
To assess c-Fos or virus expression, mice were deeply anesthetized and perfused transcardially with 10 ml of 10U Heparin/ml PBS and 30 ml of 4% paraformaldehyde. In the case of c-Fos expression, perfusions were timed to exactly 90 min after completion of behavioral tests or 120 min after i.p. injection. For all animals, brains were extracted and left overnight in 30% sucrose (wt/vol) at 4°C. The following morning, they were transferred to a 1:1 mixture of 30% sucrose and Tissue Tech ® O.C.T., followed by a transfer to plain O.C.T. These were all in the 'no injection' (homecage) group and only used for illustration in Fig. 1b, but not for any further analysis. Please see Table S4 for details. To assess virus spread, mean fluorescence was quantified using Fiji ImageJ (6). A Zeiss LSM 710 Spectral confocal microscope was used to produce exemplary scans for figures. See Table S3

Cell Counting
For analyzing images, a semi-automated, machine learning based approach was used. After brain regions were first marked in Pannoramic Viewer, they were exported and analyzed in Definiens Developer XD. First, the nuclei were segmented on DAPI using a LoG-Filter and watershed algorithms. On a set of training images, nuclei objects were manually classified as positive and negative and used as samples for training a machine learning algorithm. Mean intensity, standard deviation and local contrast parameters were extracted from these as input for a decision tree based classifier. This classifier was then applied to all nuclei objects on a larger set of images to define positive cells. To assess nuclear and cytoplasmic co-labeling, images were hand-counted by a blinded investigator. Raw c-Fos numbers obtained with this staining and quantification protocol (CEA ~ 4 -60 / mm 2 , BLA ~ 25 -100 / mm 2 , IL ~ 0 -250 / mm 2 ) were in the typical range for c-Fos IHC studies (7)(8)(9)(10). To compare c-Fos expression levels, the number of c-Fos positive nuclei was then normalized to the number of DAPI positive nuclei in any given region.

c-Fos Network Analysis
The functional network analysis was performed using an in-house application implemented in Python 3 and a series of open source libraries available at http://dx.doi.org/10.5281/zenodo.60880 under the GPLv3 License.
Input to the application consists of connectome data describing the structural connectivity of the studied regions, c-Fos expression data under different treatment conditions, a consistent name mapping scheme for region names in the connectome and the c-Fos data, and a configuration file that guides the data processing and algorithm execution.
First, the c-Fos expression data (containing multiple samples per region and mouse) was loaded and aggregated by first calculating the mean expression over each region per mouse, followed by a mean over all mice for each region. This procedure produced one average expression value per region and treatment group, which were then converted to z-scores for each region over all treatment groups: one for all groups including the control group (i.e. home cage mice) and one without the control group. The z-scored data containing the no-injection (i.e. home-cage) control group was used solely to calculate hierarchical clustering using Euclidian distance and the average linkage clustering agglomeration method. The data without the no-injection group were used for further analysis.
Second, the c-Fos data were used to compute correlation matrices describing the correlation between expression rates of two brain regions under with-drug and without-drug treatment conditions (Fig. S3a calculating the difference of  between drug and no drug measures (i.e.  = DZP -saline). Fourth, the connectome was imported from the Allen Mouse Brain Atlas, whereby recursive projections were filtered out, region names were mapped according to user specification, and multiple projections between the same two regions were aggregated by taking their mean. We used the left-brain projections from the Allen Mouse Brain Connectivity Atlas. From the connectome, a subset of brain regions was extracted that contains only regions present in the c-Fos expression data. This connectome subset was then normalized to the outgoing edges of each region (out = 1 , out being the connection strength of outgoing connections).
Next, con was calculated by taking the c-Fos correlation values calculated in the second step and scaling them by the normalized connectome generated in the third step. In order to obtain a single measure for each node in the network, the sum of all scaled correlation pairs was calculated from that node to its first order neighbors con = . This then served as a proxy, combining functional and structural information as well as reflecting the change in correlation between two regions. The effect of the drug (DZP) is analyzed by calculating the difference of con between drug and no drug measures (i.e. con = con(DZP) -con(saline)). As con is a network measure on a directed graph, it is possible to perform a node analysis only based on outgoing edges (i.e. con
The cell was allowed to reestablish constant activity post break-in for 5 minutes and a baseline was recorded for at least 5 minutes in whole-cell voltage-clamp configuration (-70mV), followed by application of Diazepam (10M) to the bath. Frequency and shape of spontaneous inhibitory post-synaptic currents (sIPSCs) were analyzed using Clampfit software (Molecular devices). Two min bins of sIPSCs from before and after addition of diazepam were analyzed for each cell. Internal solution contained the following: 130mM KCl, 10mM HEPES, 10mM EGTA, 5mM CaCl2, 2mM MgCl2, 2mM Na2ATP, 0.4mM NaGTP, 5mM Na2Phosphocreatine.

Statistical analysis.
Sample size. The sample size for all experiments was based on previous experiments and published literature (13,14).
Exclusion of data points. Virus-injected mice that did not show virus expression in post-mortem histological classification were excluded from all data analyses (see Supplementary Table S3). Further, one mouse (saline M4 cohort) had to be excluded from the analysis of the DREADD EPM experiment because it jumped off the apparatus. These routine exclusion criteria had been established prior to testing.
Randomization. Within each experiment, animals were assigned to treatment groups at random, and assignment was balanced over cages.
Blinding. For the c-Fos analysis, cell counting was either performed by an experimenter blinded to the experimental conditions or computationally. For slice electrophysiology, the experimenter was not blinded to drug treatment. For deep brain calcium imaging experiments, the experimenters were not blinded to drug treatment, but analysis was performed computationally. For behavioral experiments, the experimenter was not blinded, but analysis was performed blinded or computationally. For fMRI experiments, the experimenter was not blinded to drug treatment, but was blinded to viral expression.
Normality. For c-Fos, deep brain calcium imaging and gene expression data, normality of the underlying populations was assumed. Ex-vivo electrophysiology, mouse behavior, and virus expression data were tested for normality using D'Agostino & Pearson normality test. If normality was not detected, non-parametric statistical tests were used. For fMRI, standard statistical workflows were used, and data were treated as normally distributed.
Variation. Variance is represented by individual data points in each graph. Corrections for unequal variances were applied whenever possible.
Unless indicated otherwise, all statistical analyses were performed in Graph Pad Prism ® (Version 7). The statistical tests used are indicated in the figures legends, two-tailed, and were corrected for multiple comparisons whenever applicable. Wherever significance is not made explicit, the test did not reach statistical significance.

fMRI literature search
To find brain regions involved in the effect of BZD on fear and anxiety, PubMed (https://www.ncbi.nlm.nih.gov/pubmed/) and the Cochrane Library (http://www.cochranelibrary.com) were searched (October 2016, search terms below). Results not relevant for fear or anxiety were excluded. From the remaining literature (see Table S2), a meta score for each region was calculated as the sum of all region-specific BZD effects found in this literature set. Regions were rank ordered based on that score.

Search terms
PubMed: (   Fig. S5e). In order to compare the number of connections across stimulation condition and experimental groups the average matrices were binarized using FDR threshold (≤ 0.05). These binarized matrices were used to create networks (i.e. graphs). These graphs were visualized in a custom-made module in Amira based on the 3D coordinates of each brain structure -node -taken from our 3D mouse brain atlas derived from Paxinos et al. (5). Nodes without edges (orphan nodes) were omitted. The nodesize represents the sum of edges, i.e. significant connections with all other brain structures. For a better appreciation of the 3D node distribution, a transparent mouse brain surface of an average anatomical standard mouse brain was rendered and used as anatomical reference. Only the nodes of the limbic subnet, guided by the c-Fos data (amygdala hip area, amygdala piriform, anterior amygdala, basolateral amygdaloid nucleus, basomedial amygdaloid nucleus, bed nucleus of stria terminalis, central nucleus of the amygdala, cingulate cortex, cortical amygdala, infralimbic cortex, insular cortex, medial amygdaloid nucleus, paraventricular hypothalamic nucleus, periaqueductal gray, prelimbic cortex, septum, sublenticular extended amygdala) were colour coded (Supplementary Fig. S6b).
To emphasize the differences of the experimental conditions to the control condition "GFP saline", the average limbic subnet matrix of GFP saline was subtracted from the average limbic subnet matrices of the other groups resulting in networks with differences in functional connectivity strength as edge weight between the nodes.
These difference subnets were displayed using a forced based projection algorithm (Kamadai-Kawa in NWB, NWB Team, 2006). Thus, nodes sharing edges, here differences in functional connectivity, with each other are clustered together (Supplementary Fig. S6e). The node-size represents the sum of the absolute differences of all edges of that node. The thickness of each edge codes its weight, the absolute difference in connectivity strength. Red edges indicate an increase of connectivity strength compared to GFP saline and blue edges a decrease (Fig. S6e). The colour of each node was chosen anatomically (see Supplementary Fig. S6b).
Additionally, we took advantage of the non-invasive and high 3D resolution of the MRI modality for spatially referencing the placement of the optical fiber. Therefore, we Supplementary Fig. S8f,g).

Deep brain calcium imaging
Deep brain calcium imaging was performed using a Vista HD 2.0 in vivo Rodent Brain Imaging System sessions. Traces from these units were low-pass filtered at 0.5 Hz (Fig. 3b, right, bottom). Calcium events were detected at a peak threshold > 6 S.D and decay time  > 0.5s. (Fig. 3b, right, top). Units were projected onto mean movies (Fig. 3b, left). Neuronal activity was expressed as events score, computed as the cumulative amplitude of above threshold calcium events for each neuron in each session (Fig. 3c, top). With this workflow, we analyzed neurons from a total of 3 animals (CEAl SST + /PKC -), 4 animals (CEAl PKC + /SST -) and 4 animals (CEm), yielding a total 24 units from 2 (CEAl SST + /PKC -), 14 units from 3 animals (CEAl PKC + /SST -) and 33 units from 4 animals (CEm) with significant activity (passing PCA/ICA detection). Population activity was expressed as average of the cell-wise event scores for each session. These units were classified based on d activity BZD-activity saline session in cells with increasing (event score BZD-saline >0), decreasing (event score BZD-saline <0) or unaffected (event score BZD-saline =0) activity, respectively (Fig. 3c, bottom).

Code availability
The code used for the c-Fos network analysis is available at http://dx.doi.org/10.5281/zenodo.60880 under the GPLv3 License.   Table S4 for exact number of samples.  Supplementary Fig. S6 Optogenetic activation of PVT-CEAl modulates aversive brain state globally and is modulated by DZP. a BOLD activation patterns for the GFP and ChR2 with saline and DZP treatments at 3 representative positions (Bregma 1.94, -1.58 and -3.08) as a result of the initial GLM group analysis. b Brain-wide BOLD fMRI functional networks recruited by optogenetic activation of PVT-CEAl projections and systemic DZP application together with corresponding GFP controls. Animals virally expressing GFP (n = 4) or ChR2 (n = 3) in the PVT were laser-stimulated in the CEAl. The limbic functional subnetwork is colorized. Node size codes for sum of all edges, i.e. significant connections of that node. A 3D surface of an average mouse brain was rendered transparently. Dashed lines depict the overall network similarities of d. c Cluster analysis of brain-wide functional networks in c by column-normalized node strength. d Similarity assessment of the networks in c across the different experimental groups. Values are Pearson correlation coefficients of the respective connectivity matrices. e Kamadai-Kawai plots of differential networks, i.e. subtracting the GFP saline as control, representing DZP effects (left), effects of PVT-CEl stimulation (right) and the interaction of the two (middle). Here, the node-size represents the sum of the absolute differences of all edges of that node. The thickness of an edge codes its weight or the absolute difference in connectivity strength. Red edges indicate an increase of connectivity strength compared to GFP saline and blue edges a decrease. Note that DZP treatment rearranges functional coupling of the amygdala network and increases ILA, amygdala and PAG interaction (left). This is partially reverted by PVT-CEAl optogenetic activation under DZP treatment (middle), which reduces intra-amygdala, ILA, and PAG dominance (right). Increased interaction between amygdala and PAG may arise through DZP suppression, while PVT-CEAl activation activates the CEAm inhibitory output to the brainstem. Color code in b, e represents anatomical groups. Size represents node strength or correlation b and absolute differences of functional connectivity e. R indicates right hemisphere. AA -anterior amygdala; ACA -anterior cingulate area; Acb -nucleus accumbens; adHC -anteriodorsal hippocampus; AHtA -amygdala-hippocampus transition area; AI -agranular insular cortex; Amy -amygdala; APtA -amygdala-piriform transition area; Au -auditory cortex; BLA -basolateral amygdala; BMA -basomedial amygdala; BSTbed nucleus of the stria terminalis; CEA -central amygdala, CEAl -central amygdala lateral subdivision; CEAm -central amygdala, medial subdivision; COA -cortical amygdala; CoM -mammillary bodies; DG -dentate gyrus; ILA -infralimbic area; LS -lateral and medial septum; M2 -secondary motor cortex; MEA -medial amygdala; Orb -orbitofrontal cortex; OT -olfactory tubercle; PAG -periaqueductal gray; Pir -piriform cortex; PL -prelimbic area; Pt -parietal association cortex;PVH -paraventricular hypothalamus; Rh -entorhinal cortex; RS -retrosplenial cortex; S1 -primary sensory cortex; S1BF -primary sensory cortex barrel field; S2 -secondary somatosensory cortex; SLEA -sublenticular extended amygdala; SN -substantia nigra; vHC -ventral hippocampus; Vis -visual cortex; VTA -ventral tegmental area. Note that in all cohorts virus marker fluorescence is higher in CEAl than in BLA, and that virus fluorescence in the BLA is never above PKCbackground, indicating that the virus did not spread outside the CEAl. Bars are means ± s.e.m. Significance levels between groups at * Q < 0.05, ** Q < 0.01, *** Q < 0.001, **** Q < 0.0001. Supplementary Fig. S8 Virus expression and optical fiber placement in fMRI cohorts. a AAV::EGFP injections into the PVT of C57BL/6 mice. Data taken from Allen Mouse Brain Connectivity Atlas. Experiment number: 120875111 (http://connectivity.brain-map.org). b,c Confocal scans showing typical expression (left), as well as infection rate (right) of AAV::GFP (b) and AAV::ChR2 (c). d,e Confocal scans showing typical PVT-CEAl projections in AAV::GFP (d) and AAV::ChR2 (e) injected animals. Note that terminal fields of virally infected PVT neurons within the temporal lobe are highly specific for the CEAl (a,d,e). f,g Optical fiber placement. f 3D reconstruction of average fiber placement in the experimental animals. Red area indicates the 3D standard deviation of fiber positions. Of note, the representation of the optical fiber in MR images is considerably overestimated due to susceptibility artefacts of the fiber itself. The fiber positions were determined from the anatomical reference scans for each animal. g Location of the fiber tip in each animal. R indicates right hemisphere.

Supplementary
Supplementary Fig. S9 A circuit framework for the BZD anxiolytic effect. Our proposed model of CEA circuitry and BZD anxiolytic effect, representing the major neuronal types in CEAl. We find that DZP modulates CEA signaling both by global (center-left) and local (center-right) mechanisms, which synergize (right) to promote an anxiolytic state in CEA circuitry. Inhibitory interactions within CEA may underlie these dynamic changes.
Supplementary Figure S10 Translational evaluation of the BZD anxiolytic mechanism. Evaluation of amygdala inhibitory gating for the BZD anxiolytic effect in mice and humans. Note the strong amygdala interaction with DZP (|con|) in mice, high expression of the (GABAA-2) subunits in mice (cf. Fig. 2k) and humans (data from Allen Mouse , http://mouse.brain-map.org and Allen Human Brain Databases, http://human.brain-map.org), and evidence from the literature (see methods) for the amygdala mediating BZD anxiolytic effect from human fMRI (gray). Regions were ranked based on the absolute values in each category. Grey value intensity indicates rank above the median for functional connectivity, gene expression data and literature evidence (see methods). Boxes indicate where regions were pooled due to lower anatomical resolution in databases or literature. Human GABAA-2 expression includes the medial septal nuclei. ACAd -anterior cingulate area, dorsal part; ACAvanterior cingulate area, ventral part; AI -agranular insular area; BLA -basolateral amygdalar nucleus; BST -bed nuclei of the stria terminalis; CEA -central amygdalar nucleus; ILA -infralimbic area; LA -lateral amygdalar nucleus; LSc -lateral septal nucleus, caudodorsal part; LSr -lateral septal nucleus, rostroventral part; LSv -lateral septal nucleus, ventral part; PAG -periaqueductal gray; PL -prelimbic area; PVH -paraventricular hypothalamic nucleus; PVT -paraventricular nucleus of the thalamus. Gabra2 -GABAA receptor subunit 2. Sample number for each region and condition in the c-Fos screen. A sample number of n indicates samples from n different animals, each from at least 2 histological slices. Please note that the 'home cage' group was not used for any analysis, but only for illustration in Fig. 1b. ACAd -anterior cingulate area, dorsal part; ACAv -anterior cingulate area, ventral part; AI -agranular insular area; BLA -basolateral amygdalar nucleus; BST -bed nuclei of the stria terminalis; CEA -central amygdalar nucleus; ILA -infralimbic area; LA -lateral amygdalar nucleus; LSc -lateral septal nucleus, caudodorsal part; LSr -lateral septal nucleus, rostroventral part; LSvlateral septal nucleus, ventral part; PAG -periaqueductal gray; PL -prelimbic area; PVH -paraventricular hypothalamic nucleus; PVT -paraventricular nucleus of the thalamus