Significance Statement

  • This work combined in vitro experiments and computer simulations to study the interplay between homeostatic structural plasticity and synaptic scaling.

  • We observed a non-linear relationship between spine numbers and neural activity, where partial or complete inhibition of synaptic transmissions led to increased and reduced spine density, respec-tively.

  • Partial inhibition increased spine sizes regardless of their initial sizes, while complete inhibition selectively increased relatively large spines.

  • A bi-phasic spine-number-based structural plasticity rule reconciled the divergent experimental results in activity-dependent spine density changes.

  • We used a framework of engineering systems and complex systems to analyze the roles of the bi-phasic structural plasticity in maintaining the robustness of the neural network — firing rate homeostasis.

  • The bi-phasic rule partially uses the negative feedback strategy and serves as a redundant and heterogeneous mechanism to the synaptic-weight-based homeostatic synaptic scaling rule.

  • Both rules are based on intracellular calcium concentration which is the integral signal of neural activity. We, therefore, proposed a critical role for integral feedback control in firing rate home-ostasis.

Introduction

To survive in a dynamic world, animals must properly respond to novel and familiar environmental cues. The former demands their brain networks to sustain activity at a critical state, where information can be transmitted into action potentials and conveyed as neural avalanches.16 The latter describes experience-based learning such that subtle memory cues could trigger excessive responses in the corresponding pathway during recall. These two mechanisms are reflected by firing rate homeostasis and associative learning, respectively.

Robust firing rate homeostasis has been reported in rodents’ visual cortices and hippocampi in vivo, where neural activities get restored at either individual neuron79 or population level10, 11 within a few days after perturbation. Long before that, however, the idea of firing rate homeostasis has been hover-ing in theories since long-term synaptic potentiation (LTP) was observed.12 LTP is a positive-feedback modification of synaptic transmission among excitatory neurons (“neurons that fire together, wire to-gether”13), which was postulated by Donald Hebb14 and thus termed Hebbian plasticity. Its associative properties allow synapses to preserve memory traces as increased synaptic weights but risk pulling the network dynamics from equilibrium status into overexcitation or silence.15 Restoring the neural network to its equilibrium dynamics without erasing memories in this condition is therefore critical and attests to the robustness of the complex system—the human brain network.1621

Robustness is a ubiquitous concept in many complex systems or engineering systems, meaning that when external perturbations forcefully change its dynamics, the system adapts to return to its original attractor or moves to a new attractor state that maintains its functions.19, 20 Several features have been proven essential for maintaining robustness, including negative feedback control, redundancy, and heterogeneity. Homeostatic synaptic scaling falls perfectly into the negative feedback strategy.22 Upon long-term manipulation of neural activity,22 membrane potential,23 or calcium concentration,24 neurons proportionally down-or upscale their incoming synaptic weights in a compensatory way. Redundancy and heterogeneity encompass various mechanisms that achieve the same goal in case one or two fail. While redundancy and heterogeneity have been identified in biological systems in the form of genetic buffering and convergent molecular circuits,19, 20, 25 their representation in activity-dependent plasticity has not been discussed explicitly. Induction of homosynaptic long-term depression (LTD)26 with the Bienenstock– Cooper–Munro (BCM) rule27 or STDP rule,2831 or induction of heterosynaptic LTD via synaptic tagging and capture32 are in agreement with the idea of stabilizing firing rate. Inhibitory plasticity also contributes to re-establishing a balance between excitation and inhibition.3335 All these mechanisms support the notion that in the functional plasticity regime, there are local or global mechanisms that are redundant and heterogeneous to the homeostatic synaptic scaling.

However, whether such redundancy is also reflected in the structural substrate of functional trans-mission remains elusive. Structural plasticity includes changes in sizes or numbers of dendritic spines and axonal boutons, numbers of synapses, and network connectivity, which are all pivotal to functional transmission and memory capacity.36 Meanwhile, previous work on activity-dependent changes in spine density showed both homeostatic and non-homeostatic regulations upon activity perturbation.37 In the present study, we used whole-cell patch-clamp recordings and time-lapse imaging experiments, which in-formed the computational implementation of a bi-phasic activity-dependent (i.e., calcium-based) rule for synapse-number-based structural plasticity that reconciles various alterations of spine numbers. Then, with computer simulations, we showed that this structural plasticity rule is competitive, redundant, and compensatory to the synaptic-weight-based homeostatic synaptic scaling rule, in different time scales and magnitude scenarios. We demonstrated that calcium concentration-based integral feedback control, in the formats of structural plasticity and synaptic scaling, enables robust adaption of network activity.21, 38, 39

Results

Integral feedback mechanisms in regulating firing rate homeostasis

Neural network dynamics are determined by the product of external inputs and intrinsic connectiv-ity. External stimulation and inhibition increase or decrease neural activity deterministically in a static network. Hebbian plasticity amplifies network activity by potentiating synaptic weights via positive feedback. In contrast, homeostatic synaptic scaling and homeostatic structural plasticity reduce input-triggered responses by modifying synaptic weights or numbers in the presence or absence of Hebbian plasticity (Fig. 1A). Mechanism-wise, homeostatic synaptic scaling and homeostatic structural plasticity seem redundant in regulating network dynamics. However, compared to homeostatic synaptic scaling, ex-perimental spine density analyses yielded ambiguous, even contradicting results, especially upon chronic activity inhibition37 (Fig. 1B). Highly heterogeneous experimental systems partially explain the inconclu-sive results. Nevertheless, in a complex biological neural network where homeostatic synaptic scaling is constantly subject to and meanwhile modulates neural activity, changes in synaptic weights are expected to interact with the activity-dependent structural remodeling of networks. From the perspective of a com-plex system, we proposed to view synaptic scaling and spine-number-based structural plasticity as two integral feedback mechanisms (Fig. 1C), which use intracellular calcium concentration (C(t)) to inform their regulatory behaviors. We approximated their interactions systematically using computer simula-tions. Specifically for structural plasticity, different working models are available where neural activity either monotonically or non-monotonically governs the neurites’ outgrowth and retraction, and hence synapse formation and loss. Additional information is therefore needed to inform the model selection.

Integral feedback mechanisms in regulating firing rate homeostasis.

(Ai) Neural net-work activity is determined by external inputs in a static network. (Aii) Hebbian plasticity amplifies the network responses to external inputs via a positive feedback mechanism. (Aiii) Homeostatic plasticity restores setpoint activity via a negative feedback mechanism, which is key to firing rate homeostasis. (Bi) Hebbian functional plasticity strengthens recurrent connectivity by potentiating the weights of certain groups of synapses. (Bii) Homeostatic synaptic scaling proportionally down-scale or up-scale all synaptic strengths upon chronic excitation and inhibition. (Biii) Homeostatic structural plasticity presents homeo-static spine loss upon chronic excitation, while divergent changes in spine density upon chronic inhibition have been observed. (C) Synaptic scaling and structural plasticity use intracellular calcium concentration ([Ca2+]i, C(t)) to track neural activity (AP, action potential, S(t)). Calcium concentration updates each time with calcium influx (βCa) upon the arrival of an action potential and decays with a time constant τCa. Homeostatic synaptic scaling is implemented as a weight-dependent rule which updates the synaptic weight w(t) with a scaling factor ρ. The discrepancy from the setpoint ɛ determines the direction of weight scaling. Structural plasticity is also calcium-dependent and governs the growth and retraction of axonal boutons and dendritic spines by the setpoint value. Two examples of the structural plasticity rule are presented here. They are either linearly or non-linearly dependent on the intracellular calcium concentration.

Heterogeneous modification of spine sizes and numbers upon chronic synaptic inhibition

To examine whether spine density is linearly dependent on neural activity or not, we used 200 nM and 50 µM NBQX (2,3-dioxo-6-nitro-7-sulfamoyl-benzo[f]quinoxaline), a competitive AMPA receptor antago-nist, and studied CA1 pyramidal neurons of entorhinal-hippocampal tissue cultures (Fig. 2A). Whole-cell patch-clamp recording of AMPA receptor-mediated spontaneous excitatory postsynaptic currents (sEP-SCs) demonstrated a significant reduction in mean sEPSC amplitude and frequency with 200 nM NBQX and a near-to-complete blockade of excitatory synaptic transmission with 50 µM NBQX (Fig. 2Bi-Bii). Then the two concentrations were applied to Thy1-eGFP cultures for three days, and individual dendritic segments were followed before and after chronic activity inhibition with the help of time-lapse imaging (Fig. 2C). Analysis showed that 200 nM NBQX treatment increased spine density (p = 0.003, Wilcoxon test), while 50 µM NBQX treatment reduced spine density (p = 0.008, Wilcoxon test), in comparison to their baseline densities (Fig. 2D). No significant changes were observed in the control cultures (p = 0.06, Wilcoxon test). Our results indicated a non-linear activity dependency for spine numbers.

Heterogeneous modification of spine numbers and sizes upon three-day synaptic inhibition via NBQX.

(A) Example CA1 pyramidal neuron recorded in the entorhinal-hippocampal tissue culture for probing the effects of different NBQX concentrations. Scale bar: 500 µm. (Bi-Bii) Frequency and averaged amplitudes of sEPSC events for each recorded neuron in three groups (N = 18 for the control group, N = 18 for the 200 nM NBQX-treated group, N = 12 for the 50 µM NBQX-treated group). (C) Example Thy1-eGPF culture and the example dendritic segments from the radiatum layer (rad.) before and after three-day treatment. Scale bar: 200 µm and 5 µm. (D) Spine density at baseline and after the three-day treatment. All values were normalized by the corresponding baseline values. Lines with light shades are raw data (solid and dashed lines represent increased or reduced spine density, respectively). Dark-shaded lines with error bars are each group’s means and the standard error of the means (s.e.m.s). (N = 19 for the control group, N = 24 for the 200 nM-treated group, N = 33 for the 50 µM-treated group) (E) Cumulative distribution function of spine sizes before and after the three-day treatment (N = 489 for the control group, N = 736 for the 200 nM-treated group, N = 675 for 50 µM-treated group). Inset shows the corresponding averages of normalized spine sizes. (F) Normalized changes in spine sizes grouped by their initial spine sizes. Values on the x-axis are the upper limits of each group. (G) Each segment’s average change in spine sizes against its initial spine density. The marker size labels the net change in spine density over a three-day course. (H) Each segment’s average changes in spine sizes against its change in spine density. (I) The table summarises the experimental data, while the graph displays the extrapolated relationship between neural activity and spine densities or sizes.

We also tracked the size of individual spines before and after the treatment. The spine sizes in all three groups followed a long tail distribution in which the majority were small (Supplementary Fig. 1A) while their alterations after the three-day treatment (Δ spine size) displayed a normal distribution (Supplemen-tary Fig. 1B). Population-wise (Fig. 2E), we observed a significant increase and a significant reduction in spine sizes, respectively, after 200 nM (p < 0.001, Wilcoxon test; 99.9%CI = [9.09, 25.8], LLM) and 50 µM NBQX treatment (p < 0.001, Wilcoxon test; 99.9%CI = [9.09, 25.8], LLM). However, we noticed that after normalizing the size alterations by their initial sizes, all three groups showed an increased ten-dency (Fig. 2E inset and Supplementary Fig. 1C), which we suspected was inflated by the enlargement of large numbers of small spines. Therefore, we sorted the spine size alterations based on their initial sizes in Fig. 2F. Indeed, control segments (black curve) showed a natural fluctuation of spine sizes that small spines tend to grow over the three-day period while large spines shrink. When partial inhibition via 200 nM NBQX (orange curve) was applied, the overall dynamics shifted toward increase regardless of their initial sizes. However, complete inhibition via 50 µM NBQX (light blue curve) only tended to increase the sizes of a small number of large spines, but counterbalanced the natural increase of small spines and caused shrinkage of middle-sized spines. These data confirmed a general homeostatic adjust-ment of spine sizes upon partial and complete inhibition. However, their initial sizes seem to mediate the fate of individual spines in response to complete inhibition.

To provide insights into the interaction between spine density and size, which may hint at the rela-tionship between structural plasticity and synaptic scaling, we further quantified the average changes of all spine sizes for individual dendritic segments (spine size). We observed a positive correlation be-tween the average changes in spine sizes and their baseline spine densities under 200 nM NBQX treatment (Fig. 2G, r(24) = 0.65, p = 0.00064, Pearson’s correlation; Supplementary Fig. 2A) but not in control or 50 µM NBQX treated groups. When plotting the average changes of spine sizes against the corre-sponding alterations of spine densities for all dendritic segments, we found that weak inhibition (200 nM NBQX) triggered both spine density increase and spine size enlargement in most of the segments (see orange dots clustered in the upper right quadrant in Fig. 2H). These data implied a potentially critical role for segments with dense spines and a synergistic effect between structural and functional plasticity upon weak inhibition. In the case of 50 µM NBQX treatment, however, the size-density data points were widely scattered in four quadrants (light blue dots in Fig. 2H). This result supported our observation of heterogeneous modification of spine sizes upon strong inhibition. In summary, our imaging results suggested a non-linear relationship between network activity and spine density or spine sizes, which may result from the complex interaction between synaptic scaling and structural plasticity upon weak and strong activity deprivation (Fig. 2I).

Stabilizing and characterizing a biphasic structural plasticity rule

Our results regarding spine density alterations demonstrated a non-linear (biphasic) relationship between neural activity and spine numbers, where partial inhibition resulted in an increase in spine numbers and complete inhibition resulted in reduced spine numbers. Because changes in spine numbers or density are functionally correlated to synapse formation, loss, and rewiring, we decided to implement structural plasticity rules in point neurons to study synaptic rewiring in a large homogeneous spiking neuron network without considering the heterogeneous neural morphology (Fig. 3A). The vanilla linear growth rule was implemented for comparison reasons (Fig. 3Bi). We adopted a Gaussian-shaped growth rule of synaptic element numbers with two setpoints to represent the observed biphasic dependency for simplification. Particularly, by setting the first setpoint at η = 0 or η > 0, we achieved two slightly different variants of the Gaussian rule (Fig. 3Bii-Biii). We then grew three neural networks where the synapses among excitatory neurons were respectively subject to the three mentioned rules (Fig. 1C). As reported before,11, 4042 the linear-rule-guided network developed smoothly to a homogeneous sparsely connected network (10% of connection probability, equivalent to 1 000 excitatory synapses for each neuron) with an average firing rate around the target rate, i.e., the setpoint value ɛ = 7.9 Hz (Fig. 1Di). The Gaussian rule with a zero setpoint (η = 0) also developed the network to an equilibrium state similar to the linear rule (Fig. 1Dii). However, the Gaussian rule with two non-zero setpoints (η = 0.7 and ɛ = 7.9) did not permit proper network development (Fig. 1Diii).

Growing a neural network with three distinct structural plasticity rules.

(A) The point neuron model was used to study structural plasticity, where dendritic morphology is reduced. Dendritic spines are represented by pink sticks on the soma; axonal boutons are represented by empty or solid half circles. An empty circle with a dashed line labels an axon during retraction. Calcium concentration is linearly correlated with neural firing rate in our implementation, so neural activity would be used in the rest of the manuscript to reflect the hidden calcium dynamics. (Bi-Biii) Three activity-dependent growth rules of structural plasticity regulate the change of synaptic element numbers. (Bi) Linear rule with one setpoint (ɛ = 7.9). (Bii) Gaussian rule with two setpoints that one is zero (ɛ = 7.9 and η = 0). (Biii) Gaussian rule with two non-zero setpoints (ɛ = 7.9 and η = 0.7). Three shades indicate 100%, 50%, or 10% of the original growth rate (ν). Positive and negative values indicate, respectively, the speed of outgrowth and retraction of synaptic elements. (C) The neural network architecture of the Brunel network. 10 000 excitatory (blue) and 2 500 inhibitory neurons (red) are stimulated by external Poissonian inputs. All I-I, I-E, and E-I synapses are hard-wired with 10% probability. E-E synapses are subject to structural plasticity rules. (Di-Diii) Temporal dynamics of neural activity and network connectivity (Γ) during growth, respectively, guided by three rules. If not otherwise stated, the curve and the shaded area in activity plots represent the mean and standard deviation of the neural activity for the inhibitory population (I) and excitatory population (E). The network developed to an equilibrium state (Γ = 10%) in Di and Dii but not in Diii. The firing rates distribution and network connectivity matrices of the chosen time points, indicated by solid triangles, are included in Supplementary Figure 3 for Di-Dii, and in Figure 4 for Diii.

A closer scrutinization of the neural firing rates and network connectivity at a late growth stage showed that half of the excitatory neurons were silent (Fig. 4A), and all excitatory neurons were connected inhomogeneously where almost half of them were isolated (Fig. 4Bi-Bii). Correlating the firing rates and synapse numbers of individual excitatory neurons confirmed that the silent neurons were isolated while neurons that reached the target rate possessed a regular synapse number close to that observed in the linear rule (1 000 synapses, Fig. 4C). These observations are consistent with the properties of the Gaussian growth rules, indicating that one setpoint is stable (ɛ = 7.9), while the other one is either stable if the synapse numbers remain unchanged after silencing (η = 0) or unstable if disconnection occurs after silencing (η > 0). We thus refer to the Gaussian rule with two stable setpoints as stable Gaussian rule and the other Gaussian as biphasic Gaussian rule in the manuscript.

Silent neurons remained isolated in the network regulated by the biphasic Gaussian rule.

(A) Histogram of firing rates of excitatory and inhibitory neurons sampled at the time point indicated in Figure 3Diii. Almost half excitatory population was silent. The blue vertical line labels the mean firing rate of non-silent neurons. The orange dashed vertical line indicates the target firing rate (ɛ = 7.9 Hz). (Bi-Bii) Network connectivity matrix and the distribution of synapse numbers that individual excitatory neurons have. (C) Correlation heatmap between neural activity and synapse number of individual excitatory neurons. Neurons that were silent did not form synapses either. Neurons that fired around the target rate formed around 1 000 synapses from other active excitatory neurons. (D) Network architecture when facilitating current (Ifacilitating) was injected to boost the network development. (E) Temporal dynamics of neural activity and network connectivity when damping facilitating current were injected. The small inset shows the firing rate distributions of both excitatory and inhibitory neurons at a chosen time point (solid triangle). The facilitating current started at 750 pA and decayed linearly to zero at 4 000 s. (F) Different starting values of facilitating currents ended with different network connectivities. We used 750 pA throughout the manuscript.

Inspired by the observation that neural excitability facilitates circuit development,43 we increased the excitability of excitatory neurons during the growth period by injecting facilitating current (Ifacilitating) to lift the average membrane potential closer to the threshold potential (Fig. 4D). With the help of facilitating current, the neural network with biphasic Gaussian rule grew to the equilibrium state and retained the dynamics after the facilitating current decayed to zero (Fig. 4E). Given that different inten-sities of facilitating current resulted in different equilibrium connectivity (Fig. 4F), we used the intensity that achieved the same network connectivity (10%) and firing dynamics as the linear rule for further exploration.

Biphasic Gaussian rule reconciles divergent structural changes triggered by activity perturbation

As a proof-of-concept experiment, we systematically tuned the input strength to a subgroup of excitatory neurons (10% of the total excitatory population) from 0% to 200% fold of the original intensity (FOI, Fig. 5A) after the growth period to mimic input deprivation and stimulation. Example traces in Fig. 5Bi-Biii,C showed the protocols, neural activity, and network connectivity when the network was subject to the biphasic Gaussian rule. Intuitively, stimulation increased the neural activity of the subpopulation (S) (left panel in Bii) and triggered a reduction in synaptic connectivity among the stimulated neurons (S-S) and between the stimulated and non-stimulated excitatory neurons (S-E; left panel in Biii). The connectivity matrices in panel C display the final network connectivity. The homeostatic disconnection eventually restored the activity level of the stimulated neurons to the setpoint value (orange line in Bii). In contrast, the magnitude of activity deprivation determined whether the connectivity underwent homeostatic increase or silencing-induced degeneration. Both weak deprivation and silencing decreased the neural activity of the subpopulation (middle and right panels in Bii), while the network only managed to restore its network activity to the setpoint value by homeostatically increasing the network connectivity under weak deprivation. Silencing protocol induced disconnection of the silent neurons (right panels in Biii and C).

Divergent regulation of network connectivity upon stimulation and deprivation by three structural plasticity rules.

(A) A subpopulation (10%) of excitatory neurons (S) was subject to activity perturbation. All E-E, S-S, and E-S synapses are subject to the biphasic Gaussian rule. (Bi) Activity perturbation protocol. Three different folds of the original intensity (FOI) of the Poisson generator were used as examples to represent stimulation (110% FOI), weak deprivation (95% FOI), and silencing (0% FOI). (Bii) Temporal dynamics of neural activity of the subpopulation (S), the rest of the excitatory neurons (E), and inhibitory neurons (I) under corresponding protocols. (Biii) Temporal evolution of the overall network connectivity and connectivity of different subgroups under corresponding protocols. Synaptic connection probability from E to S (E-S) is identical to that from S to E (S-E) here so only S-E traces are shown. (C) Network connectivity matrices at the end of three protocols. (D) Average incoming synapse numbers of S neurons under different FOIs. Empty green circles are data from networks under extreme stimulation or inhibition, where neural activity and network connectivity dynamics were unstable. (E) Examples of two neurons that received the same external inputs but have different calcium decay time constants (τCa). The upper panel shows the membrane potential; the middle panel shows the spike trains of the two neurons; the lower panel displays the integrated calcium concentration over time. (F) Connectivity traces of subnetwork upon silencing under three different conditions.

Systematic analysis of incoming synapse numbers confirmed the biphasic dependency of the biphasic Gaussian rule (dark yellow curve in Fig. 5D) and the monotonic property of the linear rule (green curve). The stable Gaussian rule, which has two stable setpoints (η = 0), however, showed an intermediate tendency that external stimulation and weak deprivation triggered a homeostatic reduction or increase in synapse numbers, whereas strong deprivation left the network silent and intact (light yellow curve). Among the three rules, the biphasic Gaussian captured the homeostatic properties of the linear rule and allowed for silencing-induced spine loss, as reported in our experiments and many other previous experiments.4448

Lesion-or denervation-induced plasticity is of considerable interest in activity-dependent structural plasticity and clinical-relevant scenarios.4951 The biphasic Gaussian rule thus provides a model to study denervation-induced degeneration and subsequent regeneration. Nevertheless, in the temporal evolution of network connectivity, we noticed a small increment right after the silencing, which seemed biologically unrealistic. This might have resulted from the fast growth rate or the residual effects of previous firing activity on calcium concentration. In Fig. 5E, we showed that two neurons receiving the same external inputs presented identical membrane potential dynamics and spike trains, while the difference in calcium time constant (τCa) led to a gap in accumulated calcium concentration. Because the synaptic element number depends on the calcium concentration, which reflects the firing rate, the shorter time constant (τCa = 1 s, orange curve in 5F) carried less activity history and resulted in a smoother reduction of synaptic connectivity. Similar effects were achieved by using a slower growth rate (10% of the original rate) for structural plasticity (light blue dashed curve, 5F). Here we showed with computer simulations that the time scales of calcium dynamics and structural plasticity interfered. In reality, the decay time constant of calcium concentration is difficult to estimate due to multiple sources of calcium signals across time and space, while structural plasticity is widely accepted as a rather slow process occurring within hours to days. Therefore, we decided to use the slow growth rate and the original calcium time constant to further study the activity-dependent structural plasticity. The fast growth rate was only used in Fig. 6D-F to accelerate spine turnover in order to save simulation time when we need to display the dynamics over a longer period.

Activity perturbation and recurrent connectivity shaped the evolution of network connectivity.

(A) In the default network, we used the same growth rule for axonal boutons (A, light brown curve) and dendritic spines (D, pink curve). Alternatively, different η values could be used for axonal and dendritic elements. A silencing protocol was applied. (B) Neural activity of S, E, and I neurons under two conditions. (C) The network connectivity matrices at the end of the silencing protocol under two conditions. (D) Protocols used to examine the effects of recurrent connectivity and external stimulation. (E-F) Time courses of neural activity and connectivity upon silencing and external stimulation. Under Protocol 1, the growth rate in the left panel was ten times faster than that in the middle panel. In the right panel, external stimulation intensity was doubled.

Activity perturbation and recurrent connectivity shape the evolution of spine number

Although silencing-induced degeneration captured the feature of deprivation-induced spine loss, it was rarely observed that input deprivation led to complete neural isolation. Rather, spine recovery has been observed after spine loss.44, 46, 48 Early computer simulations with the Gaussian rule also reported such “physiological” recovery phenomena by using distinct rules for axonal bouton and dendritic spines52, 53 (Fig. 6A). We thus used different η values for bouton and spine growth and applied the same silencing protocol. The network grew properly to the equilibrium state with altered rules but failed to recover synapse numbers (Fig. 6B), only asymmetric connectivity was observed in the input and output synapses of the deprived subpopulation (Fig. 6C). However, a closer look at their network structure revealed the critical role of distance-dependent connectivity patterns in the early work. When distance-dependent connectivity was applied, neurons close to the border of the inhibited area would still receive active inputs from the non-inhibited neighbors. Therefore, these neurons were less deprived than neurons close to the center, even though their external inputs were equally removed.

Accordingly, we were inspired to examine the recurrent inputs’ role in connectivity recovery. We applied the silencing protocol and delivered external stimulation to the deprived subnetwork. In Fig. 6E (Protocol 1), the network in the left panel was simulated at a ten times faster speed than in the middle panels, such that the recurrent connectivity at the timing of stimulation was much lower in the left than in the middle panel. As stated, network dynamics are determined by the product of external inputs and internal connectivity, such that the triggered activity in the deprived subnetwork remained silent in the left panel but was reactivated in the middle panel, despite the same stimulation strength. The connectivity trace in Fig. 6F also confirmed that the residual recurrent connectivity amplified the external stimulation and initiated synapse regeneration and rewiring. However, if the external stimulation were too high (protocol 2), the over-amplification effects of the recurrent connectivity would push the system from silencing-induced degeneration to over-excitation-induced homeostatic degeneration (Fig. 6E,F, right panels). Our simulation results suggested that using recurrent connectivity, external stimulation could modulate neural activity and shape synapse number and network connectivity.

Synaptic scaling is redundant and complementary to biphasic Gaussian rule in maintaining firing rate homeostasis

Following the idea of activity manipulation, homeostatic synaptic scaling could, in theory, sculpt recur-rent connectivity by continuously adjusting functional synaptic transmission during input deprivation. Particularly, homeostatic strengthening of excitatory synapses could increase neural activity and thus may push the activity-dependent structural plasticity from silencing-induced synapse loss to rewiring, as in the case of external stimulation.

To examine this hypothesis, we implemented a monotonic rule for calcium-based synaptic scaling. We grew the network with the biphasic Gaussian structural rule to the equilibrium state and turned on the calcium-based synaptic scaling at the same time when the silencing protocol was initiated (Fig. 7A). To specify the contribution of structural plasticity and synaptic scaling, two types of connectivity were analyzed: (1) structural connectivity (Γstruc.), which was calculated by synapse numbers only; and (2) effective connectivity (Γeffec.), which was calculated by multiplying synapse numbers with the weight of individual synapses. In the scenario where synaptic scaling was absent, the synaptic weights were uniform. Thus, the two types of connectivity were identical (Fig. 7C, left panels). When allowing for weak scaling (Fig. 7B,C, ρ = 0.01), the effective connectivity to the deprived neurons (E-S) increased over time, while its impact on firing rate was negligible and therefore failed to restore synapse numbers and network connectivity. By increasing the scaling strength in the right panel (ρ = 0.02), we observed a higher temporal increment of effective connectivity in E-S synapses which managed to reactive the deprived neurons from silent status and therefore re-initiated the synapse regeneration and rewiring.

Homeostatic synaptic scaling sculpted effective connectivity to interfere with struc-tural plasticity.

(A) Protocol of silencing and synaptic scaling enabling. Three different scaling strengths were applied, that ρ = 0 is for w/o scaling, ρ = 0.01 and ρ = 0.02 represent weak and strong scaling, respectively. (B-C) Time courses of network activity and connectivity. Γstruc. is synapse-number-based structural connectivity. Γeffec. denotes the effective connectivity which multiplies synapse numbers and synaptic weights. (Di-Dii) Structural and effective connectivity matrices of the whole network at t2. (E) Time courses of firing activity and average synaptic weights for an example active excitatory neuron (blue) and an example silent excitatory neuron (grey). In the weight plot, synapses from other silent neurons (solid grey line, S) and from other active excitatory neurons (dashed grey line, E) are specified. (E) Raster plots of 100 selected inhibitory neurons (red), excitatory neurons (blue), and silent neurons (grey) from the network at t1, t2, and t3 labeled in panel B.

We further plotted the structural connectivity matrix (Fig. 7Di) and the effective connectivity ma-trix (Fig. 7Dii) at the time point t2. We showed that although the synapse numbers to the deprived subgroup decreased, the remaining synapses which received inputs from the non-deprived neurons gained weight. This observation was confirmed by the temporal evolution of the mean synaptic weights of two example neurons being deprived (gray) or non-deprived (blue) in Fig. 7E. Our simulation results con-firmed that structural plasticity and synaptic scaling react oppositely but complement each other under extreme activity deprivation. With proper scaling strength, synaptic scaling could initiate synapse regen-eration by carving recurrent connectivity to make the maximum use of recurrent inputs from neighbor neurons. On the other hand, the two rules behaved redundantly under stimulation as expected, that homeostatic disconnection occurred less when homeostatic synaptic down-scaling took place at the same time (Supplementary Figure 4).

While proper scaling strength could restore the synapse numbers and neural activity after silencing, we wondered whether the temporal network dynamics were rescued after synapse regeneration. We sampled three distinct time points (before silencing, during silencing, and after regeneration) indicated in Fig. 7B and plotted the raster plots of inhibitory neurons (red, I), non-deprived excitatory neurons (blue, E), and deprived subpopulation (grey, S). As shown in Fig. 7F, the neurons’ spiking activity was asynchronous and irregular before deprivation (t1), and the deprived neurons became silent as expected during silencing (t2). Intriguingly, synaptic scaling rescued the average firing rate of the deprived neurons, but their activity pattern became highly synchronized at t3. The highly synchronized activity observed may be attributed to the potentiation of the residual active synapses, which dominated the inputs afterward. Our results implied that although the interaction between structural plasticity and synaptic scaling triggered natural synapse rewiring after silencing or denervation, non-homogeneous connectivity among the deprived neurons may introduce unwanted dynamics, which could be of interest in the context of post-traumatic, i.e., injury-triggered-epilepsy.

Hybrid combinations of synaptic scaling and structural plasticity resolve in-consistent results of spine density

Simulations of point neuron networks allowed us to treat individual synapses homogeneously regardless of the neural morphology. However, in a neuron with highly compartmentalized dendritic segments, we may encounter different local combinations of synaptic scaling strengths and structural plasticity time scales in different neural types, which we hypothesize contribute to the inconsistent results of observed spine density in experiments (Fig. 8A).

Systematic study of the interaction between synaptic scaling and structural plas-ticity in response to activity silencing.

(A) We hypothesize that different combinations of synaptic scaling (SS) and homeostatic structural plasticity (HSP) may apply to different neuron types or dendritic segments within the same neuron, or the same type of neuron under different experimental conditions, such that the empirically observed structural plasticity was highly heterogeneous. (Sample neurons were reconstructed based on CA1 pyramidal neuron and dentate gyrus granule cells previously recorded in our lab.) (B) Simulation protocol. We systematically changed the values for the growth rate of the HSP rule (ν) and the scaling strength of the SS rule. (C-D) Example traces of neural activity and struc-tural plasticity of the deprived subpopulation (S neurons) under different parameter combinations. (E-F) Connectivity discrepancies and firing rate discrepancies when different growth rates and scaling strengths were combined. All the discrepancies were calculated by estimating the area between the actual time course and the equilibrium connection probability (10%) or the target rate (ɛ) from the time of silencing until the end of the simulation.

To address this issue, we combined various structural plasticity growth rates (ν) with various synaptic scaling strengths (ρ) in different simulations (Fig. 8B). Example traces in Fig. 8C and D display different time courses of neural activity and structural connectivity of the deprived subpopulation. It varies from condition to condition. To obtain an overview of the impact of these parameter combinations on network rewiring, we first summated the discrepancy between the actual structural connectivity among the subnetwork and its equilibrium state value (10%) over all time points upon silencing until the end of the simulation in Fig. 8E. Cold colors indicate failure in connectivity restoration, and warm colors suggest success in rewiring or over-reconnection. Slower growth rates result in slower spine loss upon silencing and are equivalent to more residual synapses. Clearly, we saw that for slow growth rates, even weak synaptic scaling was sufficient to rescue the synapse loss; for fast growth rates, the small number of residual synapses would require much stronger synaptic scaling to achieve a similar level of synapse restoration.

Similarly, we summated the discrepancy between the actual neural activity and the target rate (7.9 Hz) for the deprived population over all time points since silencing to overview the recovery of neural activity in Fig. 8F. We saw that the combinations that failed to restore network connectivity in Fig. 7E presented negative discrepancies from the target rate (dark blue in Fig. 7F), suggesting a failure in reactivating the silent neurons. However, it is worth noticing that the combinations that succeeded in synapse regeneration in Fig. 7E (warm colors) did not necessarily achieve smooth firing rate restoration. Instead, only those parameters that achieved mild synapse regeneration (light pinkish colors panel E) presented smooth firing rate restoration. Fast growth rates of structural plasticity in combination with strong synaptic scaling strength may lead to very high neural activity (bright colors) which resembles epileptiform activity or result in a drop in the transient firing rate (dark blue colors). Our systematic study confirmed that calcium-based synaptic scaling and homeostatic structural plasticity interact with each other at various time scales and strengths, which may explain the inconsistent experimental results reported.37 Our results also suggested that fast-time-scaled homeostatic structural plasticity and strong synaptic scaling would harm rather than facilitate firing rate homeostasis.

Discussion

The mammalian brain is a complex system consisting of billions of neurons and non-neuron cells packed into layers, regions, and circuits.54 This system is highly robust such that ample external information does not cause catastrophic forgetting,55 nor will input deprivation destroy the network dynamics.79 The most extensively observed synaptic-weight-based plasticity in experiments, Hebbian plasticity and homeostatic synaptic scaling, are regarded as major components in balancing flexibility and stability.56, 57 Spine-number-based structural plasticity, however, remains less understood due to divergent experimental observations.58 Nevertheless, structural plasticity, which allows for synapse formation and rewiring, has been proven to greatly increase the memory capacity of the network36 and maintain its criticality,59, 60 suggesting its discernible role in network function and potential interaction with functional plasticity. This manuscript elaborated on the synapse-number-based homeostatic structural plasticity rule and its interaction with the synaptic-weight-based homeostatic synaptic scaling rule. Our time-lapse-imaging experiments revealed a non-linear dependency between spine density and neural activity, which informed a biphasic structural plasticity rule for the systematic study. Computer simulations showed that this biphasic rule and homeostatic synaptic scaling are redundant and interchangeable upon stimulation. However, they behave oppositely under silencing: neurons would reduce spine numbers while potentiating the remaining input synapses. Synaptic scaling may rescue silencing-induced spine loss by restoring neural activity. Our study reconciled varied experimental results and suggested that the redundancy and heterogeneity between homeostatic structural plasticity and synaptic scaling achieve robustness of the complex network of neurons—firing rate homeostasis.

Our work generally extended the scope of structural plasticity by exploring its functional interaction with synaptic scaling. Structural plasticity is less extensively modeled than functional plasticity, due to its high-dimensional morphological data. Many morphological changes may occur in neurons, ranging from dendritic tree branching pattern, axon length, and bouton size, to spine density and head volume, in response to environment and experience, disease and injury, or activity perturbation. These changes would influence network connectivity. Therefore, structural plasticity was frequently proposed to account for functional cortical rearrangement in humans after injury, amputation, or disease-related cognition decay, but only a few dimensions of neural morphology were modeled mathematically. Except for a few models about morphological development,6163 most structural plasticity models focus on synapse formation and rewiring, which was regarded as the key difference from functional plasticity models. However, early models inherited the synapse-specific view of Hebbian plasticity and modeled the growth and removal of individual synapses. In these models, less-frequently-used or small-weighted synapses were pruned from a highly-connected network, or from a network with regular random synapse growth.6468 They are variants of the correlation-based rule and also share Hebbian plasticity’s properties and limitations. Unlike synapse-specific correlation-based rules, Butz and van Ooyen proposed a cell-autonomous homeostatic structural plasticity rule that follows a setpoint firing rate or calcium concentration to create or delete spines and boutons randomly.52, 69 This spine-number-based rule is a structural analog to synaptic scaling in maintaining firing rate homeostasis. Apart from that, it preserves the correlation-based associative properties11, 41, 70 and allows for self-organized criticality in developing neural networks,59 which is beyond synaptic scaling and suggests an alternative solution to integrate Hebbian and homeostatic plasticity. Our work used the same model, but added biological features inspired by experiments and investigated its function in more biological scenarios where other plasticity rules also apply.

Particularly, we suggested a biphasic activity dependency for homeostatic structural plasticity with experi-ments. As reported in a recent review,37 spine density alterations observed in experiments seem not as homeostatic as synaptic scaling, especially under prolonged input deprivation. Heterogeneous experimental paradigms may cause inconsistency. Some used lesion methods, for instance, which could increase the secretion of neuroin-flammatory cytokines besides input deprivation.46, 47, 71, 72 Neuroinflammatory cytokines, such as TNF-α, have been shown to influence plasticity induction.7375 Non-traumatic deprivation methods, such as dark rearing76 or whisker trimming,77 avoided the confounding effects of injury-related cytokine release. But the resulting input deprivation was non-homogeneous on the whole dendritic tree, as apical and basal dendrites usually receive dis-tinct inputs. The same principle may also apply to pharmacological treatments targeting pre-synaptic parts78 or glutamate uncaging in the neighborhood of a few spines.79 However, in studies where drugs homogeneously treated postsynaptic neurons, changes in spine density remain inconsistent.8085 A closer look at their experi-mental systems suggested a potential role of deprivation strength. Using non-competitive or competitive drugs (with different concentrations), conducting in vivo or in vitro, and embedding within intact circuits or not all lead to the actual magnitudes of inhibition on the target neuron varying from experiment to experiment. Our recent work has demonstrated that depriving synaptic inputs to CA1 pyramidal neurons by lesioning entorhinal cortex or CA3 Schaffer collateral triggered distinct pathway-dependent cellular, structural, and molecular alterations, where the differentiated dominance of the two pathways may play a role among other factors.86 Given that, we hypothesized a biphasic activity-dependent rule to capture the rich scenarios of spine density alterations, that partial activity suppression shall promote spine growth while complete inhibition induces spine loss. Indeed, our time-lapse imaging results confirmed that partial inhibition (200 nM) increased spine density while complete inhibition (50 µM) reduced spine density, and informed our selection of a biphasic rule. Nevertheless, our results also suggested complex behaviors of spine size upon inhibition, which may influence our interpretation of spine density, and we will discuss this later.

In the shape of a Gaussian growth curve, the biphasic structural plasticity model utterly reconciled the rich scenarios of activity-dependent spine density alterations observed in experiments. As indicated by the growth rule, the unstable setpoint (η > 0) leads to varied results in synapse numbers (equivalent to spine numbers) after activity perturbation: increase, no change, or reduction. The silencing-induced synapse loss is especially beyond the reach of any linear (monophasic) rules. As a derived property, this rule enables phase transition by changing neural activity via manipulating external inputs or intrinsic connectivity. By applying external stimulation at proper strengths and timings, for instance, it is possible to reverse the silencing-induced synapse loss into homeostatic rewiring, which may provide insights for using non-invasive brain stimulation modalities to prevent stroke-induced neural lesions.87 However, Butz and van Ooyen have suggested that the Gaussian rule could achieve a natural transition from silencing-induced spine loss to rewiring if different rules for axonal boutons and dendritic spines are used accordingly,52, 53, 88 with which our simulation in Brunel network disagreed. Instead, we showed that the “recovery” observed in their study was the combined effect of distance-dependent connectivity and distinct axonal/dendritic rules. The distance-dependent connectivity created an inhibition gradient, that the neurons at the border were less inhibited than those close to the center, and triggered a gradual connectivity recovery with the Gaussian rule. Our systematic study demonstrated the fundamental effects of neural activity and resulting calcium concentration in regulating retracting and rewiring, regardless of whether it is driven by intrinsic connectivity or external inputs.

Homeostatic structural plasticity usually happens in tandem with functional synaptic scaling, such that they may influence each other’s behavior and make it hard to interpret the observed experimental results. In our experiments, we observed increased spine size and density after three days of 200 nM NBQX treatment. Spine size enlargement may raise the chance of small spines being detected as a spine. Therefore, the spine density increase observed in imaging experiments is a compound effect of spinogenesis, spine stability, and small spines’ potentiation. The same confounding effects happen to interpret synaptic scaling results. Turrigiano et al.22 observed changes in the amplitude but not the frequency of mEPSC after chronic inhibition and excitation, leading to the conclusion of weight scaling, while subsequent experiments reported either frequency and/or amplitude changes. The conversion of silent synapses into active ones89 may increase the frequency but reduce the average amplitude of mEPSCs or sEPSCs due to their small weights. Such changes, by definition, are synapse-number-based structural plasticity, but may not be reflected by imaging techniques, depending on the laser intensity or exposure time. Furthermore, we observed divergent size changes of small and large spines upon complete inhibition. Assuming that shrunk spines may disappear and enlarged big spines will contribute to increased sEPSC amplitudes, we would again obtain different and even opposite results from imaging experiments. Therefore, spine numbers and synaptic weights might be the two sides of the coin named “effective transmission”, and judging it by either side will hide the whole picture of how structural and functional plasticity interfere. Our systematic study tackled this question and demonstrated that the interaction between biphasic structural plasticity and monotonic synaptic scaling achieves an economical and robust control of firing rate homeostasis.

Synaptic scaling originally describes the compensatory changes of synaptic weights in response to activity perturbation, and several variants of mathematical rules were extrapolated and implemented. These models used target firing rate,9093 normalized synaptic weights,94 or regulated presynapses95, 96 to realize synaptic scaling, except one suggested an integral feedback control algorithm,97 where the authors used a variable with an expo-nential kernel to represent the internal status of neural spiking activity through time. This variable regulates synaptic weights accordingly to reach a target value. However, by the time they came up with the model, it was yet unknown that calcium concentration is the variable that tracks the internal trace of spiking activity.98 Taking advantage of this model, we implemented the calcium-based synaptic scaling rule and studied its interaction with the calcium-based structural plasticity model. Our results first suggested their redundancy. Assuming both rules are monophasic activity-dependent, they are fully redundant, and the one with a faster time scale would dominate. In other words, if synaptic scaling acts faster than structural plasticity, we shall easily observe mEPSC amplitude (or frequency) changes rather than altered spine density. Structural plasticity will join the game only if the counteracting effects of synaptic scaling are insufficient. However, supposing structural plasticity follows a biphasic rule while synaptic scaling follows a monotonic rule, they are only redundant upon stimulation and weak deprivation but bifurcate upon strong inactivity, as confirmed by our simulation. As a result, we may observe neurons with reduced spine density and enlarged spine head sizes in our experiments and others.82, 83 From an economical point of view, spine growth and synapse formation are constrained by the availability of proteins; the transcription and expression of proteins require the energy currency ATPs. Therefore, we postulate that it is more effective for neurons to potentiate the existing synapses than searching for new inputs with new spines. Our sys-tematic study revealed an interdependency between synaptic scaling and structural plasticity, which may explain divergent results obtained in different brain regions, neurons, dendritic segments, and experimental settings.

In summary, how homeostatic plasticity sculpts network connectivity, including synaptic weights and synapse numbers, resembles an integral feedback controller (Supplementary Figure 5), which is ubiquitous in engineering to achieve robust performance regardless of external perturbation. In our setting, both structural plasticity and synaptic scaling rules use calcium concentration as the control signal to drive connectivity kernel modification via negative feedback. Calcium concentration integrates neural firing rates over time and returns a filtered signal responding to its activity. Using this filtered signal rather than the raw firing rate to control the connectivity modification could address the cumulative temporal effect of error in estimating the firing rate. However, the detailed mechanisms of the two rules are slightly different; the former is closer to a Proportional, Integral, Differential (PID) controller, while the latter is a PI controller. Except for the shared integral (I) component, the proportional (P) component for structural plasticity is the growth rate of the growth rule,40 as it provides a direct modification on the number of synaptic elements created or deleted per time unit based on a magnitude relative to the difference between the expected and real values of calcium concentration. For synaptic scaling, the scaling factor (β) regulates the proportional effects. Finally, the unique derivative (D) component for structural plasticity is given by the steepness of the growth curve, which decides how fast or slow the changes should reach the neuron’s desired activity. The derivative component is critical to ensure a smooth and robust response to stimuli by avoiding overshooting and abrupt oscillations that are not biologically realistic, which is especially relevant in the Gaussian case, as the distance to the target firing rate has different accelerations depending on the position on the firing rate (or calcium concentration) axis. Therefore, although the exact growth rule shape remains unknown for homeostatic structural plasticity, we propose that a biphasic curve with changing slope provides a redundant and heterogeneous structural backup for synaptic scaling to maintain robust firing rate homeostasis. This view also prepared us to anticipate their interaction with other plasticity rules where sub-or supra-threshold calcium dynamics are involved.99, 100

Materials and Methods

Ethics statement

We used mouse pups postnatal at 3 to 5 days (P3-P5) from C57BL/6J (wild-type) and Thy1-eGFP mouse lines to prepare entorhinal-hippocampal tissue cultures in the current study. All animals were kept under a 12 h 12 h light-dark cycle with food and water provided ad-libitum. One male and one or two female(s) were kept within the same cage for breeding. All animal experiments were approved by the appropriate animal welfare committee and the animal welfare officer of Albert-Ludwigs-University Freiburg, Faculty of Medicine under X-21/01B and X-18/02C. All effort was made to reduce the pain or distress of animals.

Preparation of tissue cultures

Entorhinal-hippocampal tissue cultures were prepared as published before.101, 102 All tissue cultures were cul-tivated for at least 18 days inside the incubator with a humidified atmosphere (5% CO2 at 35 °C to reach an equilibrium status. The incubation medium consists of 50% (v/v) 1× minimum essential media (#21575 022, Thermo Fisher, USA), 25% (v/v) 1× basal medium eagle (#41010 026, Thermo Fisher, USA), 25% (v/v) heat-inactivated normal horse serum, 25 mM 1 M HEPES buffer solution (#15630 056, Gibco), 0.15% (w/v) sodium bicarbonate (#25080 060, Gibco), 0.65% (w/v) glucose (#RNBK3082, Sigma), 0.1 mg/ml streptomycin, 100 U/ml penicillin, and 2 mM glutamax (#35050 061, Gibco). The incubation medium was renewed 3 times per week. If not stated otherwise, the incubation medium applied to the tissue cultures was always pre-warmed to 35 °C and adjusted around pH = 7.38.

Experimental design

The main objective of the experiments in the current study was to probe whether there is a linear or non-linear dose-dependent regulation of dendritic spine density such that we can implement a spine-number-based structural plasticity rule close to biological reality for further systematic studies. Competitive AMPA receptor antagonist NBQX (2,3-dioxo-6-nitro-7-sulfamoyl-benzo[f]quinoxaline) at different concentrations was bath applied to wild-type cultures while recording synaptic transmission at CA1 pyramidal neurons. Two magnitudes of inhibition (“partial” and “complete”) were thus determined by the reduction of amplitudes and frequencies of spontaneous excitatory postsynaptic currents (sEPSCs). Then the same concentration ladders were used to treat Thy1-eGFP cultures for three days, where we tracked individual dendritic segments of CA1 pyramidal neurons before and after the three-day treatment to investigate whether two different magnitudes of AMPA-receptor inhibition resulted in divergent alterations of spine densities and spine sizes.

Whole-cell patch-clamp recordings

To probe the effects of NBQX administration at different concentrations on synaptic transmission, whole-cell patch-clamp recordings were conducted in CA1 pyramidal neurons. Recordings were performed at 35C. The bath solution contained (in mM) 126 NaCl, 2.5 KCl, 26 NaHCO3, 1.25 NaH2PO4, 2 CaCl2, 2 MgCl2, and 10 glucose (aCSF) and was continuously oxygenated with carbogen (5% CO2/95% O2). Glass patch pipettes had a tip resistance of 4 6 MΩ, filled with the internal solution which contained (in mM) 126 K-gluconate, 10 HEPES, 4 KCl, 4 ATP-Mg, 0.3 GTP − Na2, 10 PO-Creatine, 0.3% (w/v) biocytin. The internal solution was adjusted to pH = 7.25 with KOH and reached 290 mOsm with sucrose). We patched 6 neurons per culture to record the spontaneous excitatory postsynaptic currents (sEPSCs) of CA1 pyramidal neurons in the voltage-clamp mode at a holding potential of 70 mV. Series resistance was monitored before and after each recording and the neuron data was excluded if the series resistance went above 30 MΩ. Each neuron was recorded for 2 min.

NBQX treatment

NBQX (2,3-dioxo-6-nitro-7-sulfamoyl-benzo[f]quinoxaline) is a competitive antagonist of AMPA receptors.103 We chose two concentrations 200 nM and 50 µM and delivered by bath treatment to achieve a partial or complete inhibition of AMPA receptor currents. Wild-type cultures were either recorded in normal ACSF or ACSF that contained two different NBQX concentrations. In time-lapse imaging experiments, Thy1-eGFP cultures were treated with 200 nM and 50 µM by adding NBQX (Cat. No. 1044, Tocris Bioscience, Germany) in the incubation medium for three days. Whenever we changed the new medium, fresh NBQX was administrated accordingly.

Tissue fixation and immunohistochemical staining

Recorded cultures were fixed and stained for post hoc inspection. Cultures were fixed by immersing into 4% (w/v) paraformaldehyde (PFA) in 1× phosphate-buffered saline (PBS, 0.1 M, pH = 7.38) for 1 h and transferred into 1× PBS for storage at 4 °C after being washed in 1× PBS. Before staining, all fixed cultures were again washed three times with 1× PBS (3 × 10 min) to remove residual PFA. We incubated the cultures with Streptavidin 488 (1 : 1000, #S32354, Invitrogen, Thermo Fisher, USA) in 1× PBS with 10% (v/v) in normal goat serum and 0.05% (v/v) Triton X-100 at 4 °C overnight. In the next morning, cultures were rinsed with 1× PBS (3 × 10 min) and incubated with DAPI (1 : 2000) in 1× PBS for 20 min. After another 4 washes with 1× PBS (4 × 10 min), we mounted the cultures on glass slides with DAKO anti-fading mounting medium (#S302380 2, Agilent) for confocal microscope imaging.

Time-lapse imaging

To inspect whether neural spine densities and spine sizes were altered by the three-day administration of NBQX, we employed time-lapse imaging to follow the same apical dendritic segments of CA1 pyramidal neurons before and after treatment. Live cell imaging was performed at a Zeiss LSM800 microscope with 10× water-immersion (W N-Achroplan 10×/0.3 M27; 420947-9900-000, Carl Zeiss) and 63× water-immersion objectives (W Plan-Apochromat 63×/1,0 M27; 421480-9900-000, Carl Zeiss). Thy1-eGFP tissue cultures where clear CA1 pyramidal neurons could be identified were used in this experiment. Dendritic segments from the radiatum layer were imaged. We imaged individual dendritic segments prior to the NBQX treatment and again after the three-day treatment. During the imaging session, the membrane insert with 4 cultures was placed into a 35 mm petri dish filled with 5 ml incubation medium on a platform which was constantly maintained at 35 °C. We used the same pre-warmed and pH-adjusted incubation medium for imaging procedures.

Experimental data quantification

Spontaneous excitatory postsynaptic currents (sEPSCs) were analyzed using the automated event detection tool from the pClamp11 software package as previously described.104

z-stacked fluorescent images of Thy1-eGFP cultures were projected to create a 2D representation of individual dendritic segments. ImageJ plugin Spine Density Counter105 were used to count spine numbers and measure segment length, which estimates spine density. For the same dendritic segments imaged at different time points, special attention was paid to ensure the same starting and ending points were used. Post hoc visual inspection was applied to ensure the spine detection results were not strongly biased. Both raw spine density and normalized spine density by baseline were used in the analysis.

The same z-projected fluorescent images were used to track individual spines for spine size analysis. To eliminate the bias from drawing and automatic reconstruction, we drew circles manually around the spine to cut it from the dendrite, the spine size was estimated by measuring the signal intensity with an arbitrary unit of the drawn circle. The drawing and measurements were performed with FIJI ImageJ. Both the raw values and normalized values by baseline spine size were used in the analysis. Statistical methods were specified in the individual results section.

Statistical analysis

Dunn’s multiple comparison test was applied for statistical analysis regarding the sEPSC events among the three groups. For spine density analysis, the Wilcoxon test was applied to compare the values of each segment before and after the three-day treatment. If not otherwise stated, “ns” means no significant, “*” means p < 0.05, “**” means p < 0.01, “***” means p < 0.001. For spine size analysis, we first applied the Wilcoxon test for individual spines. To account for data clustering within each segment, linear mixed model (LMM) with spine_size timing was applied to compare the values before and after treatment with the segment ID as group factor as used in.42 The significance of LLM results was judged by whether the confidence interval (CI) crossed zero.

Neuron model and network model

We used the same spiking neuron model and network architecture as described before in41 and in.11, 42 Current-based leaky integrated-and-fire point neuron was used for both excitatory and inhibitory neurons. We build an inhibition-dominated network with 10 000 excitatory neurons and 2 500 inhibitory neurons.106 To simplify the scenario, we only grow the connections within the excitatory population (E-E) with the activity-based structural plasticity rule (see the Structural plasticity rule section below). Each inhibitory neuron was beforehand hard-wired randomly to receive synapses from 10% of the excitatory and inhibitory population. All details and parameters concerning neural and network models can be found in the Supplementary Materials. We performed the network simulations with NEST simulator 2.20.2 and NEST 3.0107 and MPI-based parallel computation. All the model parameters and protocols can be found in Supplementary Tables 1-7.

Structural plasticity rule

We enabled the growth, retraction, and rewiring of synapses among excitatory neurons with the help of structural plasticity rules. By definition, each excitatory neuron has multiple dendritic spines and axonal boutons, which are called synaptic elements. Synapses were formed by randomly matching free compatible synaptic elements. The growth and retraction of synaptic elements, or in other words, the number of synaptic elements, is governed by a growth rule. Three structural plasticity rules were explored in the current study: (i) linear growth rule; (ii) Gaussian growth rule with a zero setpoint and a non-zero setpoint; (iii) Gaussian growth rule with two non-zero setpoints. All three rules are determined by a function of calcium concentration that reflects neural activity,

where C(t) is the time evolution of calcium concentration. Calcium concentration decays with a time constant τCa and increases with calcium influx (βCa) upon the emission of an action potential S(t) of the postsynaptic neuron. This operation is performed internally by the NEST simulator and works as a low-pass filtered signal of the spiking activity of the neuron. The growth of synaptic elements is regulated differently depending on the calcium concentration and the shape of the growth rule.

Linear growth rule

The linear rule was first introduced in40 and systematically studied in an inhibitory-dominant neural network.11, 41, 42, 70

The number of synaptic elements (z(t)) is linearly dependent on the calcium concentration,

where ν is the growth rate, and ɛ is the target level of calcium concentration. Since calcium concentration reflects the neural activity loyally, the target level also suggests a setpoint of firing rate in the context of a certain neural network. As discussed before, when neurons fire below their target rate, they grow new synaptic elements and form new synapses. On the other hand, they break existing synapses and retract synaptic elements when they fire above the target rate (setpoint).

Gaussian growth rule

The Gaussian rule has a more complex dependency on the calcium concentration when neural activity is too low, as shown in Equation 3. This rule was suggested and explored initially in.52

where , and . In this rule, η and ɛ are two setpoints: ɛ is the stable setpoint as used in the linear rule. When the neuron fires above ɛ, it retracts synaptic elements as in the linear rule. η is another setpoint introduced specifically for the Gaussian rule, determining the regulation manner when the neuron activity drops below ɛ. When the neuron is firing below ɛ but above η, the number of synaptic elements present will undergo homeostatic outgrowth, but when the neuron is firing below η, neurons will break synapses and retract elements. In the case where η = 0, neurons cease to change synaptic elements when their firing rate drops to zero.

Homeostatic synaptic scaling

In order to achieve homeostatic synaptic scaling, we make use of a new synaptic model in NEST called scal-ing_synapse. In this synapse model, the weight of the synapse is regulated by the difference between a homeostatic setpoint and the calcium trace of the postsynaptic neuron,

where ρ is the scaling factor, and ɛ is the same target value as used in the homeostatic structural plasticity rule.

Activity perturbation

To examine different activity-induced scenarios on neural network connectivity, we performed systematic activity manipulation to a subnetwork of excitatory neurons (Nsub = 1 000), by changing its Poissonian input from 0% to 200% fold of the original intensity (FOI). All manipulations were performed at 6 000 s when the network has grown to the equilibrium state respectively with three structural plasticity rules. For the Gaussian rule with two non-zero setpoints, we applied damping current injection to the soma to facilitate its growth within the first 4 000 s and the network was simulated for another 2 000 s without any facilitating current.

Quantifying firing rate, network connectivity, and synapse number

Firing rate

Firing rate was calculated by the average spike count over a recording period for individual neurons. We have long interval (5 s) and short interval (1 s). Short intervals were only used within the short time window after activity perturbation to reveal its transient dynamics; otherwise, long intervals were used.

Network connectivity

Two types of connectivity were used in the present study. We used a N × N connectivity matrix (Aij) to represent the recurrent excitatory connections of our network, where columns and rows correspond to pre-and postsynaptic neurons. For structural connectivity, the entry Aij of the matrix represents the total number of synaptic connections from neuron j to neuron i. For effective connectivity, we integrated the synapse number with individual weights for each pair of neurons by summating the total weights. So the entry of the connectivity matrix is the equivalent number of unit synapses, by dividing the sum with a uniform weight 0.1 mV. To average the mean connectivity of the whole network or a subnetwork at any given time t, corresponding columns and rows of the connectivity matrix were selected and averaged by . Synapse numbers were calculated by the sum of the entry in the structural connectivity matrix.

Synapse number

Input and output synapse numbers, also called indegree and outdegree in another context, were calculated by summating the input and output synapse numbers of individual excitatory neurons based on the structural connectivity matrix.

Quantifying the discrepancies in firing rate and connectivity from the target values

To apply a systematic comparison among different combinations of the synaptic scaling strengths and structural plasticity growth rates, we summated the discrepancies in firing rate and connectivity from the target values for the subpopulation over time as an index of activity and connectivity recovery. All the discrepancies were calculated by estimating the area between the actual time course and the equilibrium connection probability (10%) or the target rate (ɛ) from the time of silencing until the end of the simulation. The method was explained in detail in Supplementary Fig. 6.

Supportive Information

The authors declare that they have no conflict of interest.

Abbreviated title

Redundancy and heterogeneity between homeostatic functional and structural plasticity

Acknowledgements

The work was supported by Deutsche Forschungsgemeinschaft (DFG; Project-ID 259373024 B14–CRC/TRR 167 to AV). The research leading to these results has received funding from the European Union’s Horizon 2020 Framework Programme for Research and Innovation under the Specific Grant Agreement No. 945539 (Human Brain Project SGA3). This research has also been partially funded by the Helmholtz Association through the Helmholtz Portfolio Theme Supercomputing and Modeling for the Human Brain. We acknowledge the use of Fenix Infrastructure resources, which are partially funded by the European Union’s Horizon 2020 research and innovation programme through the ICEI project under grant agreement No. 800858.

Author contributions

H.L. and A.V. conceived the project. M.L., H.L., and A.V. designed the experiments. M.L. performed the whole-cell patch-clamp recordings and analyzed the recorded data. H.L. performed the time-lapse imaging experiments and analyzed the spine data. H.L. established the bi-phasic HSP model in the network simulation and performed systematic explorations. S.D. established the HSS model and integrated it into the HSP model. H.L. and S.D. performed numerical simulations for the interaction between the two models. H.L. made the figures and drafted the manuscript. All authors contributed to the revision.