A stage-structured population model for activity-dependent dendritic spines

ABSTRACT Here we present a novel application of stage-structured population modelling to explore the properties of neuronal dendrites with spines. Dendritic spines are small protrusions that emanate from the dendritic shaft of several functionally important neurons in the cerebral cortex. They are the postsynaptic sites of over 90% of excitatory synapses in the mammalian brain. Here, we formulate a stage-structured population model of a passive dendrite with activity-dependent spines using a continuum approach. This computational study models three dynamic populations of activity-dependent spine types, corresponding to the anatomical categories of stubby, mushroom, and thin spines. In this stage-structured population model, transitions between spine type populations are driven by calcium levels that depend on local electrical activity. We explore the influence of the changing spine populations and spine types on the development of electrical propagation pathways in response to repetitive synaptic input, and which input frequencies are best for facilitating these pathways.


Introduction
Efforts to unravel the functional role of dendritic spines, especially in relation to learning and memory, have led to some of the most exciting research in neuroscience. First described in 1888 by Ramón y Cajal [6], dendritic spines are bacterial size membrane protrusions from the surface of the neuron, consisting of a bulbous or mushroom-shaped head connected to the dendritic shaft by a neck or stem. Over 90% of excitatory synapses in the cerebral cortex terminate on spines and there are over ten trillion spines (among the hundred billion neurons) in the human brain [16]. Their structures are dynamic, in that they can physically morph into different spine classifications in response to various physiological, behavioural or pathological states. This synaptic activity-dependent structural plasticity exists over a vast range of time scales, from minutes to days or weeks [8,15,16].
Although spines may have a continuum of anatomical configurations [1], studies have traditionally grouped spines based on their shapes [12,18]. The groups are generally : Magnified region from A shows synapses between presynaptic axon terminals and three kinds of postsynaptic dendritic processes: stubby s (type-I), mushroom m (type-II), and thin t (type-III) spines. Note that at a fixed location along the dendrite more than one spine type may protrude with synaptic connections to the same or different axons. classified as (see Figure 1) stubby (Type-I), mushroom (Type-II), and thin (Type-III). In early optical studies of cortical pyramidal cells it was observed that spine morphology correlates with distance from the cell body [12,18]; proximal to the soma, where the dendrite has the largest diameter, spines have short necks and appear stubby, whereas distally, where the dendrite is thinner, spines have thinner and longer stems. At intermediate distances, a variety of spine types were seen including the intermediate mushroom-shaped spines. It is important to note that the accuracy of these observations and the subsequent conclusions have been challenged by more recent imaging and statistical analyses. In 1982, motivated by these early optical studies, Rinzel [22] proposed that the observed anatomical configurations presented an apparent paradox: why should the membrane potential response initiated at a distant synaptic site, which is severely attenuated due to cable properties alone, suffer additional attenuation due to a long, thin, and therefore more electrically resistive, spine stem?
Rinzel, in conjunction with earlier work of Rall [19,21], provided a possible explanation based on analysis and computations with a passive cable model. They found that if the resistance to current flow in the stem is correlated with the input resistance to the cable at the spine base, then attenuation is maximally sensitive to changes in the stem resistance. Rinzel suggested that if such an operating range exists for adjusting synaptic potency, then it is not surprising to find thinner spine stems at more distal sites and stubbier spines near the soma. This analysis of passive spines, which predated what is now known about calcium-mediated restructuring of spines, did provide an explanation for how the spines arrived at their spatial configuration. Although Rinzel's mathematical analysis provided an explanation that added support to the results of the optical studies, today it is viewed that the spatial distribution of spines depicted in those early investigations were not accurate. However, the value of Rinzel's effort provided an important contribution that linked the electrical properties of passive spines to how spines may spatially configure along the dendrite.
In this paper we develop the framework for a stage-structured population model to illustrate how spine configurations observed in dendrites (e.g. cortical pyramidal cells) may form and how the morphology can evolve in response to patterns of synaptic activity. Conversely, for excitable spines, we show how an emerging pattern of spine morphology along the dendrite can influence the generation and propagation of an impulse (emerging signalling properties). The model is based on experimental findings suggesting that spines are compartments that regulate local changes in the calcium concentration independently of the parent dendrite [11,13], and that moderate calcium concentration levels cause spine elongation but higher concentrations lead to spine shrinkage [23]. We break from previous modelling efforts [7,25,26] by specifically tracking the population density of three types of spines: stubby, mushroom, and thin, and allow for the possibility that multiple spine types may emerge radially from the membrane at any point along the dendrite. The model is stage-structured with spine types representing reversible growth stages. Stage transitions between spine types depend on intraspine calcium levels that increase or decrease in response to patterns of electrical activity in the cell. For example, high calcium levels in the thin population reduce the spine density of that population as thin spines change to mushroom spines. The primary motivation for introducing a stage-structured approach is that it opens the door for direct comparison of computational results with morphological analysis based on experiments.
The paper is organized as follows. The model is formulated in Section 2. In Section 3 we derive a reduced system of ordinary differential equations for passive spines to investigate how dendritic diameter and the frequency of synaptic input influences spine morphology. In Section 4 we demonstrate, using the full PDE continuum model, that synaptic input to excitable spines can drive spine restructuring within the input region as well as far outside of the region. We show that this is accomplished by activity-dependent spine restructuring that increases the density of mushroom spines, forging a pathway for propagation. Section 5 is Discussion.

A stage population model for dendritic spines
We formulate a continuum model for the stage transition between stubby, mushroom, and thin dendritic spines. New in this model is that multiple spine types may coexist simultaneously at points along the dendrite (see Figure 1), and the model allows for the study of population dynamics between spine types. A continuum formulation [3] is used for tracking the interaction between the many activity-dependent (active or passive) spines and for studying the impact of their individual and collective dynamics on the output properties of the dendrite.
To begin, consider a passive dendritic cable of electrotonic length l (μm), with both ends sealed, studded with three distinct populations of dendritic spines. The spine densities N i , i = 1, 2, 3 represent the number of spines of type i per unit physical length. Over a short segment x, the three spine types deliver current xN i I ssi to the dendrite, where I ssi represents the current flowing through the stem of type i. The stem current (I ssi ) is expressed as an I · R voltage drop across the spine stem resistance R ssi (M ), given by where V shi and V d (mV) are respectively the membrane potential of spine head i and the dendritic shaft. The spine stem is modelled as in previous studies [3,24] as a lumped Ohmic resistor, neglecting the stem's membrane and cable properties. If the potential in the spine head is larger than the potential in the dendritic shaft (V shi > V d ), then I ssi > 0 and the current is flowing from spine head i to the spine base. Conversely, if the potential in the base is larger than the potential in the head (V d > V shi ), then I ssi < 0 and the current flow is in the opposite direction. If I ssi = 0, then no current is passing through that spine stem. The electrical potential V d (x, t) in a passive dendrite studded with 3 i=1 N i spines per unit physical length satisfies the cable equation Here R i ( · cm) is the specific cytoplasmic resistivity (here subscript i refers to intracellular, not spine type); R m ( · cm 2 ) is the resistance across a unit area of passive membrane; C m (μF/cm 2 ) is the specific membrane capacitance and d (μm) is the diameter of the dendrite. Parameter values for the cable are in Table 1.
It is convenient to rewrite the cable equation in terms of dimensionless (electrotonic) length. After multiplying through by R m /π d, we substitute into Equation (2) the Here, n i represents the number of spines of type i over length λ (or number of spines per unit electrotonic length; denoted here by spines/e.l.). It is assumed throughout this paper that both ends of the dendrite are sealed, with a uniform resting potential of zero in the cable and the spine heads. Each spine head is modelled as an isopotential compartment with surface area A sh (μm 2 ) and specific membrane capacitance C m (μF/cm 2 ). The capacitance of an individual spine head is therefore C sh = A sh C m (μF). It is well known that different spine types have different spine head areas. However, we assume here that mushroom, thin, and stubby spines have the same spine head geometry, but have different capacities for sequestering and releasing intraspine calcium stores. An equation for the membrane potential in a single spine of type i is obtained from a current balance relation for the capacitive, ionic, spine stem and synaptic currents given by The term I ioni represents ionic currents passing through the head membrane, and I syni represents synaptic current. If the membrane is passive, then I ioni = V shi /R sh , where R sh is the membrane resistance of the head. When modelling excitable membrane in the spine heads, we used Hodgkin-Huxley kinetics [9] for voltage-dependent ion channel currents: Here, g Na , g K , g L and V Na , V K , V L are maximal conductances and reversal potentials, respectively, for sodium, potassium and leakage currents. We followed [24] and used increased channel densities (γ = 2.5) and a temperature of 22 • C.
We simulate the activation of synapses by applying to all spine types i in an activation region, X 0 ≤ X ≤ X 0 + X, where V syn is the synaptic reversal potential (assumed the same value in each spine type) and g syn is a brief synaptic conductance generated by the α-function which reaches a peak synaptic conductance g p , t p ms after activation. Equations (6) and (7) model synaptic currents observed in experiments and are similar to the equations for synaptic current used in other models [3,20]. Parameter values for spine heads, ionic and synaptic currents may be found in Table 1.
Next we combine experimental evidence for how calcium drives structural change in spine shape [13,16,17,23] with morphological evidence that spines tend to organize in defined groups by shapes, classified as stubby, mushroom, or thin [12,14,18]. Although research into the precise mechanism and associated detailed mathematical models for how spines change their shape is ongoing (for example, see [10]), we define the spine stem current (I ssi ) as the electrical influence, in spine type i, that drives changes in intraspine calcium levels. When there is no activity (I ssi = 0) we assume that calcium levels decay to a minimum value C min exponentially. The calcium equations are where constants κ i are important scale factors that control the steady-state calcium levels in response to synaptic input. These scale-factors quantify a prediction of how much calcium a particular spine type is able to release from intracellular stores, which may depend on the spine head area and specialized membrane apparatus associated with that spine type. Figure 2(A) is a modification of Figure 1 in [7] (adapted originally from the conceptual model by [23]). In Figure 2(A), when the function f (C a ) is positive, for intermediate values of calcium C 1 < C a < C 2 , spines become longer and thinner through a transition from stubby to mushroom to thin. For moderately low (C 0 < C a < C 1 ) or high (C a > C 2 ) calcium values, the transition is from thin to mushroom to stubby. All spine transitions are considered sequential; e.g. a stubby spine undergoing transition to a thin spine, must first go through the mushroom spine stage, and vice versa. In this formulation we assume that there is no spine pruning or formation of new spines, so the total density of spines is conserved (See Discussion). For very low calcium concentrations, it is assumed that spine transitions do not occur. The stoichiometric representation of the transition between spine C 2 ) f > 0 and spines grow through transitions from stubby to mushroom to thin (s→m→t). For moderately low or high values of calcium (C 0 < C a < C 1 or C a > C 2 ), f < 0 and spines shrink through transitions in the opposite direction (t→m→s). For very low calcium (C a < C 0 ), f = 0 and spine structure is static. (B): Stoichiometric representation of the transitions between spine types quantified by densities n 1 , n 2 , and n 3 . The calcium dependent rates K ± i are defined in Equation (12).
types is shown in Figure 2(B) and the corresponding equations are Here, the rates K ± i are dependent on calcium and defined by where δ i , i = 1, 2, 3, are constants. The above definition of the calcium-dependent rates are consistent with the dynamics described by Figure 2(A). For example, suppose that at some location x along the dendrite, mushroom spines have calcium levels in the intermediate range (C 1 < C a < C 2 ) but the other two spine types have high levels of calcium (C a > C 2 ). We expect that the mushroom spines would see a positive contribution to their density n 2 from thin transitions to mushroom (K − 3 = 0), yet we would not expect any stubby spines to transition to mushroom (K + 1 = 0). Since mushroom spines have intermediate calcium levels, they would be making a transition to thin spines but not to stubby (K + 2 = 0, but K − 2 = 0). Whether or not the mushroom spine density ultimately increases, decreases, or stays the same, depends on the balance between the terms K − 3 n 3 and K + 2 n 2 . Summarizing, the stage population model for three types of dendritic spines is where the K ± i 's are defined above in Equation (12) and the currents in Equation (14): I ssi , I ioni , and I syni , are defined in Equations (1), (5), and (6) respectively.
In this work, we choose the functional form of f to be (see Figure 2(A)) System (13)-(18) is nonlinear even if the spine head membrane is passive, since the rates K ± i are nonlinearly dependent on calcium. However, it is worth noting that for passive spine heads, Equations (13)-(14) constitute a linear subsystem which acts on a fast time scale (milliseconds) compared to subsystem (15)- (18), which behaves on a slow time scale (seconds). In our simulations, the synaptic input will be repeated at least every 10 ms, a period long enough to allow the synaptic inputs I syni and potentials V shi and V d to return to rest (see Figure 3(A)). However, intraspine calcium values do not return to rest (Figure 3(B)) due to their relatively slow decay rates. The steady-state calcium values for each spine type is governed by their respective scale-factors κ i . Figure 3(B) shows a case where the steady-state calcium level of the thin spines is largest, followed by mushroom and stubby. The mechanisms for handling calcium in a spine, along with spine morphology, ultimately determines the magnitude and time course of calcium (see [28] for a more thorough discussion and issues surrounding calcium handling in a spine and their parent dendrite). Here we assume, for simplicity, that the time course of the calcium concentration in each spine type is similar, but the magnitudes may differ. We examine all permutations of steady-state calcium levels. The model does not take into account the affect of calcium in the parent dendrite. Finally, it is important to point out that the biology does not preclude extremely large disparities in time scale between fast and slow subsystems. For practical computational reasons we chose time scales for the slow system, δ i and i , just small enough so that decreasing the order of magnitude of these parameters does not qualitatively affect the dynamical properties of the system (see [27]).
For spatially uniform synaptic input and spatially uniform distribution of spines, the partial differential equation system reduces to a system of ordinary differential equations which were integrated and cross-checked using standard solvers. Extensions to the model, such as how to incorporate the formation and loss of synapses, and how to add spine populations to handle more than one input frequency at points along the dendrite, are addressed in the discussion.

Population dynamics of passive spines
In this section we reduce system (13)- (18), with passive spines, to a system of ordinary differential equations to investigate how (1) dendritic diameter and (2) synaptic input frequency influence spine morphology along the dendrite. The reduction is motivated by computational studies which show that spatially uniform repetitive synaptic input to passive spines, within an input region, drive spine restructuring primarily within that input region, with little or no effect on membrane potential and structural change outside the region [26]. This is in stark contrast to excitable spines, where theoretical studies suggest that synaptic activity can drive restructuring and membrane potential far outside the input region [7,26]. For the passive spine case, we use the above insights to model an isolated segment of cable with sealed ends corresponding to the endpoints of the input region. If the passive cable segment receives synchronous synaptic input, and initially each spine type is uniformly distributed, model (13)- (18) reduces to the following system of ordinary differential equations: To investigate how dendritic diameter influences spine morphology in a passive cable we use system (20)- (25) to simulate calcium dynamics and changes in spine densities. Figure 4 shows the calcium dynamics and spine densities for stubby, mushroom, and thin spine types, for three different dendritic diameters: d = 0.5, 1.0, 1.5 μm. A synaptic input frequency of 50 Hz (synaptic input ever 20 ms) is applied to all three spine types. Initially the spine densities are equal, n i = 20 spines per electrotonic length (e.l.). In this simulation parameters κ i are set so that thin spines accumulate the largest intraspine calcium, followed by mushroom and stubby. For d = 0.5 the calcium for all three spine types does not exceed C 2 = 300 nM in response to repetitive synaptic input. Calcium levels remain in the range C 1 < C a < C 2 . According to Figure 2(A) the direction of spine transition is s → m → t. Indeed, as Figure 4(A) (right) shows, the thin spines dominate while mushroom and stubby spine densities decay toward zero. It is important to note that the total spine density remains constant at all times, since we are modelling transitions between spine types. As the diameter increases, the calcium levels for all three spine types shift upward. Figure 4(B) (left) shows the case where the thin spine calcium concentration rises above C 2 , mushroom calcium is partially above C 2 , while stubby calcium remains below this critical value. Figure 2(A) indicates that both thin and stubby spines should transition to mushroom spines. The mushroom spine calcium concentration oscillates above and below C 2 , spending most of the time above C 2 . From Figure 2(A), a larger fraction of mushroom spines should become stubby, while a smaller fraction becomes thin. This scenario creates the possibility for coexistence between spine types. Figure 4(B) (right) demonstrates that this is the case; for large time, the mushrooms dominate stubbys, while the thin spines decay to low but nonzero values.
To better understand how coexistence comes about and obtain an analytic estimate of the spine densities for this case, consider subsystem (23)- (25) for large time in conjunction with Figure 4(B) (left). To simplify the analysis assume that the function f (C a ) in Figure 4 (Equation (19)) is replaced by a piecewise continuous function that transitions between 1 and −1 at C 2 . Since the large time calcium concentration for stubby spines falls below C 2 , from Equation (12) the rate K + 1 = δ 1 . The thin spine concentration lies above C 2 , so K − 3 = δ 3 . The mushroom spine concentration straddles C 2 spending fraction p of its time above and (1 − p) below the critical value; therefore K − 2 = pδ 2 and K + 2 = (1 − p)δ 2 . In the simulation δ i = δ, i = 1, 2, 3 (δ = 2 × 10 −3 ms −1 ), and the approximating linear subsystem is with initial conditions n i (0) = 20, i = 1, 2, 3. The solution to this linear system is ⎛ For large time ⎛ demonstrating coexistence. In Figure 4(B) (left), p ≈ 0.9, and from Equation (28) the long time spine density values are estimated as n 2 ≈ 30 (mushroom type), n 1 ≈ 27 (stubby type) and n 3 ≈ 3 (thin type), which agrees well with the numerical results in Figure 4(B) (right). Note from Equation (28)) that the density of stubby and thin spines depend on the fraction of time the calcium values for mushroom spines stay above or below C 2 , whereas the mushroom spines approach a final density (half the total density) that is independent of its straddling configuration. For diameter d = 1.5 in Figure 4(C) (left), both thin and mushroom calcium concentrations eventually exceed C 2 , but this time the calcium for the stubby spines straddle C 2 . In this case, K + 1 = (1 − p)δ, K − 2 = δ, K + 2 = 0, and K − 3 = δ, and for large time In Figure 4(C) (left), p ≈ 0.8, and from Equation (29), the long time spine density values for this case are estimated as n 1 ≈ 50 (stubby type), n 2 ≈ 10 (mushroom type) and n 3 ≈ 0 (thin type), which agree well with Figure 4(C) (right). Figure 5 depicts bifurcation diagrams that plot the steady-state values of spine densities for t, m, and s over the dendritic cable diameter range 0.5 ≤ d ≤ 2, for a synaptic input frequency of 50 Hz. We refer to each panel in this figure as bifurcation diagrams since steady-states are computed for each value of d, as opposed to using the full PDE model to compute the steady-state solution over a cable that has a continuously increasing diameter (a flared cable). In the figure, scale factors κ i are chosen to represent six permutations of the steady-state calcium levels for the three spine types. To clarify, consider the steady-state C a levels for spine types t, m, and s in Figure 4  In all panels of Figure 5, thin spines dominate for d < 0.7 whereas stubby spines dominate for d > 1.4. This result, for passive spines, suggests that for a cable that tapers from the soma, stubby spines dominate near the soma whereas thin spines dominate distally. This result for passive spines is consistent with the theoretical results of Rinzel [22]. Finally, note that for medium range diameters multiple spine types can coexist.
In Figure 6 the bifurcation parameter is now the synaptic input frequency rather than the cable diameter. In this figure the cable diameter is fixed (d = 1.5). The bifurcation structures of Figure 6 are remarkably similar to the corresponding structures in Figure 5. Figure 6 shows for higher frequencies of synaptic input the stubby spines dominate the other spine types due to excessive levels of calcium, whereas in the low frequency range thin spines dominate. In the middle range of frequency we see a coexistence between spine types. Figure 4 and 5 suggest the possibility of an interesting interplay between a hard-wired geometric parameter (cable diameter) and a dynamic variable (the frequency of synaptic input). A possible implication of this interplay is addressed in the Discussion.

Forging a pathway for propagation
Similar to previous results [3], the stage-structured continuum spine model suggests that the propagation of a pulse along the cable is correlated with the density of spines with a medium range of stem resistance. Here, we have categorized these spines as the mushroom spine type. Also, higher densities of stubby and thin spines along with mushroom spines can enhance pulse propagation, but stubby and thin types without mushroom type spines do not lead to propagation.
In this section, we present results of numerical simulations for the model described by Equations (13)- (18), and Hodgkin-Huxley kinetics (Equation (5)) for ionic currents. We apply a repetitive synaptic input with a frequency of 50 Hz to the spine heads at one end (0 ≤ X ≤ 0.2) of a dimensionless cable of length L = 3 studded uniformly with 180 spines (60 of each type). Here, parameters that reflect the scaling of calcium dynamics are κ 1 = 1.7, κ 2 = 1.92, and κ 3 = 2.4 nM/(mA·ms). Biological evidence suggests that persistent synaptic activity leads to mushroom type spines that are more stable and resistant to morphological changes [4,5,29]. To reflect this, we reduce the rate of change for mushroom spines δ 2 from 0.001 to 0.0001 ms −1 . The model implementation employs a semi-discretization with a first-order numerical differentiation approximation for the second derivative and then solving the system of resulting ordinary differential equations using MATLAB's ODE15 solver. The results for this numerical approximation are represented in Figures 7-9. Figure 7(A) shows the dynamics of spines and membrane potential during the first 6 ms of initial synaptic input. The left panel shows the dynamics of spines along the cable, the  For active spines with activity-dependent densities, the peak head and dendritic potentials increase slowly as the mushroom spine density increases. Here we employ the same parameter values as Figure 7. The time evolution of the peak spine head potential and the dendritic membrane potential, along with the mushroom spine density, are plotted for three spatial locations: at the input site X = 0.12, at X = 1.2, and distant from the input site at X = 2. (A): Rise in peak potential in the mushroom spine head at three spatial points along the dendrite. (B): Peak dendritic potential rises with the same shape as the peak head potentials at the same spatial location, but with smaller peak values due to the spine stem resistance. (C): After 50 inputs, the mushroom type spine density increases at all locations, with the highest increase at the input site. middle panel is the membrane potential at mushroom spine's head and the right panel is the dynamics of dendrite membrane potential over time for three different locations X = 0.12, X = 1.2, and X = 2. Figure 7(B) left shows how the density of three spine types changes toward the domination of mushroom spine type along the cable. The middle and right panels demonstrate how the potential at mushroom spine's head and dendrite will change at three different locations of the cable. At the Figure 7(C) left the complete dominant of mushroom spine type is visible after 80 synaptic input. This domination of mushroom spines lead to prorogation of pulse along the dendrite. Figure 8 shows the change in the peak membrane potential at three different locations along the cable, X = 0.12, X = 1.2, and X = 2, as a result of the changes in the mushroom spine densities at these locations. Parameter values are the same as those used in simulations shown in Figure 7. When repetitive stimulation is applied to the spines heads in the initial portion of the cable, 0 ≤ X ≤ 0.2, the synaptic current term I syn in Equation (4) and Hodgkin-Huxley kinetics for ionic current trigger firing in the spine head. This increase in the spine head potential V sh increases the spine stem current I ss , resulting in a change in the calcium level at the spine head (Equation 13). Changes in the levels of calcium cause the dominance of the mushroom spine type for the chosen values of calcium scale-factors (κ 1 = 1.7, κ 2 = 1.92, κ 3 = 2.4 nM/(mA·ms)). This increase in I ss also increases the dendrite potential V d in Equation (3). Due to the voltage gradient, this increment propagates along the dendritic cable, changing the value of I ss outside of the stimulated area. This increase in turn generates spikes on the spine heads, and this chain reaction increases the number of mushroom spines throughout the cable, forging a pathway for propagation of a spike along the dendrite.
In Figure 9, we compare the propagation of impulses along the dendritic cable for different frequencies of synaptic stimulation. This figure depicts the propagation of the peak potential for mushroom type spine heads for different frequencies following 0.8 ms of simulation time. We increased the frequency from 20 to 100 Hz by 20 Hz increment, while the other parameters are the same as those in the Figure 7. For synaptic stimulation at a frequency of 20 Hz, the spike in potential for mushroom type spines reaches the location X = 1 after 0.8 ms. The maximum distance of the spike propagation obtained after 0.8 ms increases with increases in frequency, e.g for a frequency of 40 Hz, propagation reaches X = 1.5, and for a frequency of 60 Hz, propagation almost reaches X = 1.8. The spike reaches the other end of the cable at a input frequency of 100 Hz, and the visibly higher peak potential seen near the end of the cable is due to the sealed end boundary condition.

Discussion
The goal of this study is to provide a framework for modelling populations of multiple spine types and transitions between them. Novel to this study is the introduction of stage-structured dynamics applied to distinct populations of dendritic spines. The stage transition between spine types is formulated to depend on intraspine calcium levels that modulate in response to patterns of electrical activity. In this framework, the electrical activity is generated by repetitive synaptic input to spine heads. The magnitude of the spine stem current is used as the electrical influence that drives the changes in intracellular calcium levels [7,26].
In response to repetitive synaptic input, we find for passive stubby, mushroom and thin spine-types that thin spines emerge as dominant for small dendritic cable diameters, stubby spines dominate for large diameters, and multiple spine types including mushroom spines coexist at intermediate diameters ( Figure 5). A remarkably similar pattern or configuration of passive spines emerges when the dendritic diameter is fixed, but the frequency is varied. Here, thin spines dominate in response to low frequency input, stubby spines dominate at high frequencies, and mushroom spines coexist in the mid-frequency range ( Figure 6). This interplay between a hard-wired parameter (cable diameter) and a dynamic variable (the frequency of synaptic input) leads to an interesting interpretation. Figures 5 and 6 suggest the following: Consider a uniform constant diameter cable where all spines are stimulated by a repetitive synaptic input. However, suppose the frequency is spatially dependent; specifically, the frequency alternates periodically between high and low values along the cable. Figures 5 and 6 infer that the steady-state structure of emerging spine types would be equivalent to a dendritic cable with varicosities with all spines stimulated synchronously at the same frequency. The implication is that a steady subthreshold background pattern of synaptic input could have the virtual effect of morphing the cable geometry! When spines are excitable, the analysis and equations that determine the bifurcation structure in Figures 5 and 6 are no longer applicable. Using the full model (Equations (13)-(18)), we find that excitable spines drive the spread of membrane potential away from the input site, influencing the spine distribution pattern far from the input region. In this way, the model forges a pathway for impulse propagation due to the distal formation of mushroom spines. In our model, the mushroom spines have intermediate stem resistance values that sufficiently isolate the spine heads for generating local action potentials, but spine heads are not so isolated as to preclude action potential propagation. This result is consistent with early studies of spine impulse propagation dynamics [2,3]. Our computations studies reveal spatial distributions of mushroom spine densities along the cable that facilitate impulse propagation (Figures 7-9). It is important to add that currently it is extremely difficult to test the key predictions of this model experimentally. It is difficult to count spines and it is even more difficult to obtain accurate data identifying and counting spine types. With this said, we expect that as imaging and the precise control of synaptic input is perfected it will be possible to test the key predictions of this model.
In this stage-structured population model, all spine transitions are considered sequential; that is, a stubby spine undergoing transition to a thin spine, must first go through the mushroom spine stage, and vice versa. In this formulation we assume that there is no spine pruning or formation of new spines, so the total density of spines is conserved. However, previous models of activity-dependent changes in spine structure allowed for overall changes in spine density [7,26]. For example, in line with biological observations [23], the model described by Crook et al. [7] allows for the creation of new spines for intermediate levels of calcium, like those associated with spine elongation. Note that in the brain, newly formed spines begin as long, thin filopodia, which probe the environment prior to forming a functional synapse and becoming a thin spine. In contrast, for very low or very high levels of calcium like those typically associated with spine stem shortening, spines may completely disappear resulting in a decrease in spine density. These mechanisms could be added easily to the stage-structured formulation described here by allowing for increases in the density of thin spines to reflect the creation of new spines and decreases in the density of stubby spines to reflect the loss of spines, known as spine pruning, based on calcium levels.
In the simulation studies considered here, repetitive synaptic input at a single frequency is applied to spine heads in specific areas of the dendritic cable, similar to controlled experimental studies that aim to carefully manipulate spine structure in specific locations. In further computational studies, more biologically realistic synaptic inputs with multiple frequencies or complex spatio-temporal patterns of synaptic inputs could be applied to the model. In addition, the stage-structured population model formulation provides the means to administer inputs with different frequencies or patterns to different, specific spine types.