Bilinearity in Spatiotemporal Integration of Synaptic Inputs

Neurons process information via integration of synaptic inputs from dendrites. Many experimental results demonstrate dendritic integration could be highly nonlinear, yet few theoretical analyses have been performed to obtain a precise quantitative characterization analytically. Based on asymptotic analysis of a two-compartment passive cable model, given a pair of time-dependent synaptic conductance inputs, we derive a bilinear spatiotemporal dendritic integration rule. The summed somatic potential can be well approximated by the linear summation of the two postsynaptic potentials elicited separately, plus a third additional bilinear term proportional to their product with a proportionality coefficient . The rule is valid for a pair of synaptic inputs of all types, including excitation-inhibition, excitation-excitation, and inhibition-inhibition. In addition, the rule is valid during the whole dendritic integration process for a pair of synaptic inputs with arbitrary input time differences and input locations. The coefficient is demonstrated to be nearly independent of the input strengths but is dependent on input times and input locations. This rule is then verified through simulation of a realistic pyramidal neuron model and in electrophysiological experiments of rat hippocampal CA1 neurons. The rule is further generalized to describe the spatiotemporal dendritic integration of multiple excitatory and inhibitory synaptic inputs. The integration of multiple inputs can be decomposed into the sum of all possible pairwise integration, where each paired integration obeys the bilinear rule. This decomposition leads to a graph representation of dendritic integration, which can be viewed as functionally sparse.


Introduction
For information processing, a neuron receives and integrates thousands of synaptic inputs from its dendrites and then induces the change of its membrane potential at the soma. This process is usually known as dendritic integration [1][2][3]. The dendritic integration of synaptic inputs is crucial for neuronal computation [2][3][4]. For example, the integration of excitatory and inhibitory inputs has been found to enhance motion detection [5], regularize spiking patterns [6], and achieve optimal information coding [7] in many sensory systems. They have also been suggested to be able to fine tune information processing within the brain, such as the modulation of frequency [8] and the improvement of the robustness [9] of gamma oscillations. In order to understand how information is processed in neuronal networks in the brain, it is important to understand the computational rules that govern the dendritic integration of synaptic inputs.
Dendritic integration has been brought into focus with active experimental investigations (see reviews [1,10] and references therein). There have also been many theoretical developments based on physiologically realistic neuron models [11,12]. Among those works, only a few investigate quantitative dendritic integration rules for a pair of excitatory and inhibitory inputs [3,13] and there has yet to be an extensive investigation of the integration of a pair of excitatory inputs or a pair of inhibitory inputs. In this work, we propose a precise quantitative rule to characterize the dendritic integration for all types of synaptic inputs and validate this rule via realistic neuron modeling and electrophysiological experiments.
We first develop a theoretical approach to quantitatively characterize the spatiotemporal dendritic integration. Initially, we introduce an idealized two-compartment passive cable model to understand the mathematical structure of the dendritic integration rule. We then verify the rule by taking into account the complicated dendritic geometry and active ion channels. For time-dependent synaptic conductance inputs, we develop an asymptotic approach to analytically solve the cable model. In this approach, the membrane potential is represented by an asymptotic expansion with respect to the input strengths. Consequently, a hierarchy of cable-type equations with different orders can be derived from the cable model. These equations can be analytically solved order by order using the Green's function method. The asymptotic solution to the second order approximation is shown to be in excellent agreement with the numerical solutions of the original cable model with physiologically realistic parameters.
Based on our asymptotic approach, we obtain a new theoretical result, namely, a nonlinear spatiotemporal dendritic integration rule for a pair of synaptic inputs: the summed somatic potential (SSP) V S can be well approximated by the summation of the two postsynaptic potentials V 1 and V 2 elicited separately, plus an additional third nonlinear term proportional to their product, i.e., V S (t)~V 1 (t)zV 2 (t)zk(t)V 1 (t)V 2 (t): The proportionality coefficient k encodes the spatiotemporal information of the input signals, including the input locations and the input arrival times. In addition, we demonstrate that the coefficient k is nearly independent of the input strengths. Because the correction term kV 1 V 2 to the linear summation of V 1 and V 2 takes a bilinear form, we will refer to the rule (1) as the bilinear spatiotemporal dendritic integration rule. In the remainder of the article, unless otherwise specified, all the membrane potentials will be referred to those measured at the soma.
We note that our bilinear integration rule is consistent with recent experimental observations [3]. In the experiments [3], the rule was examined at the time when the excitatory postsynaptic potential (EPSP) measured at the soma reaches its peak for a pair of excitatory and inhibitory inputs elicited concurrently. We demonstrate that our bilinear integration rule is more general than that in Ref. [3]: (i) our rule holds for a pair of excitatory and inhibitory inputs that can arrive at different times; (ii) our rule is also valid at any time and is not limited to the peak time of the EPSP; (iii) our rule is general for all types of paired synaptic input integration, including excitatory-inhibitory, excitatory-excitatory and inhibitory-inhibitory inputs.
Our bilinear integration rule is derived from the two-compartment passive cable model. We then validate the rule in a biologically realistic pyramidal neuron model with active ion channels embedded. The simulation results from the realistic model are consistent with the rule derived from the passive cable model. We further validate the rule in electrophysiological experiments in rat hippocampal CA1 pyramidal neurons. All of our results suggest that the form of the bilinear integration rule is preserved in the presence of active dendrites.
As mentioned previously, there are thousands of synaptic inputs received by a neuron in the brain. We therefore further apply our analysis to describe the dendritic integration of multiple synaptic inputs. We demonstrate that the spatiotemporal dendritic integration of all synaptic inputs can be decomposed into the sum of all possible pairwise dendritic integration, and each pair obeys the bilinear integration rule (1), i.e., where V S denotes the SSP, V p E denotes the p th individual EPSP, V q I denotes the q th individual inhibitory postsynaptic potential (IPSP), k ij EI , k kl EE , and k mn II are the corresponding proportionality coefficients with superscripts denoting the index of the synaptic inputs. We then confirm the bilinear integration rule (2) numerically using realistic neuron modeling. The decomposition of multiple inputs integration in rule (2) leads to a graph representation of the dendritic integration. Each node in the graph corresponds to a synaptic input location, and each edge connecting two nodes represents the bilinear term for a pair of synaptic inputs given at the corresponding locations. This graph evolves with time, and is all-to-all connected when stimuli are given at all synaptic sites simultaneously. However, based on simulation results and experimental observations, we can estimate that there are only a small number of activated synaptic integration, or edges in the graph, within a short time interval. Therefore, the graph representing the dendritic integration can indeed be functionally sparse.
Finally, we comment that, in general, it is theoretically challenging to analytically describe the dynamical response of a neuron with dendritic structures under time-dependent synaptic conductance inputs. One simple approach to circumvent this difficulty is to analyze the steady state of neuronal input-output relationships by assuming that both the synaptic conductance and the membrane potential are constant [3,12]. Such analyses can be applied to study dendritic integration, but they usually oversimplify the description of the spatial integration, and fail to describe the temporal integration. Another approach to circumvent the difficulty is to study the cable model [14,15] analytically or numerically. For the subthreshold regime, in which voltage-gated channels are weakly activated, the dendrites can be considered as a passive cable. Along the cable, the membrane potential is linearly dependent on injected current input. This linearity enables one to use the Green's function method to analytically obtain the membrane potential with externally injected current. In contrast, the membrane potential depends nonlinearly on the synaptic conductance input [12]. This nonlinearity greatly complicates mathematical analyses. Therefore, in order to solve the cable model analytically, one usually makes the approximation of

Author Summary
A neuron, as a fundamental unit of brain computation, exhibits extraordinary computational power in processing input signals from neighboring neurons. It usually integrates thousands of synaptic inputs from its dendrites to achieve information processing. This process is known as dendritic integration. To elucidate information coding, it is important to investigate quantitative spatiotemporal dendritic integration rules. However, there has yet to be extensive experimental investigations to quantitatively describe dendritic integration. Meanwhile, most theoretical neuron models considering time-dependent synaptic inputs are difficult to solve analytically, thus impossible to be used to quantify dendritic integration. In this work, we develop a mathematical method to analytically solve a two-compartment neuron model with time-dependent synaptic inputs. Using these solutions, we derive a quantitative rule to capture the dendritic integration of all types, including excitation-inhibition, excitation-excitation, inhibition-inhibition, and multiple excitatory and inhibitory inputs. We then validate our dendritic integration rule through both realistic neuron modeling and electrophysiological experiments. We conclude that the general spatiotemporal dendritic integration structure can be well characterized by our dendritic integration rule. We finally demonstrate that the rule leads to a graph representation of dendritic integration that exhibits functionally sparse properties.
constant synaptic conductance [16,17]. The approximation can help investigate some aspects of dendritic integration, however, the approximation in such a case is not sufficiently realistic because the synaptic conductances in vivo are generally time-dependent. On the other hand, one can study the dendritic integration in the cable model numerically. The compartmental modeling approach [14] enables one to solve the cable model with time-dependent synaptic inputs. This approach has been used to investigate many aspects of dendritic integration. For instance, it was discovered computationally that dendritic integration of excitatory inputs obeys a certain qualitative rule, i.e., EPSPs are first integrated nonlinearly at individual branches before summed linearly at the soma [18,19], which was verified later in experiments [20,21]. Clearly, the computational approach can help gain insights into various phenomena of spatiotemporal dynamics observed at the dendrites, however, a deep, comprehensive understanding often requires analytical approaches. Note that this point has also been emphasized in Ref. [22]. Here, our analytical asymptotic method can solve the cable model with time-dependent synaptic inputs analytically and reveal a precise quantitative spatiotemporal dendritic integration rule, as will be further illustrated below.

Results
We first study the dendritic integration of a pair of excitatory and inhibitory inputs (E-I), a pair of excitatory inputs (E-E), and a pair of inhibitory inputs (I-I) case by case. In each case, we first analytically derive the bilinear integration rule from the twocompartment passive cable model, and then validate the bilinear integration rule using the realistic model of a pyramidal neuron with both active channels and dendritic branches; we further validate the bilinear integration rule in electrophysiological experiments in rat CA1 pyramidal neurons. We then derive the bilinear integration rule for multiple excitatory and inhibitory inputs, and validate this rule in the simulation. Based on our bilinear integration rule for multi-inputs, we finally propose a graph representation of dendritic integration.

Bilinear Rule for E-I Integration
We begin to study the spatiotemporal dendritic integration of a pair of excitatory and inhibitory inputs. An analytical derivation of the bilinear integration rule is described in the section of Derivation of the Rule. The details of the cable model used in the derivation can be found in the section of Materials and Methods. The validation of the bilinear integration rule using the realistic neuron modeling and electrophysiological experiments is described in the section of Validation of the Rule. The spatial dependence of the coefficient k in the rule is described in the section of Spatial Dependence of k EI .
Derivation of the rule. Most neurons possess complicated dendritic morphology, however, for simplicity, we start to investigate the spatiotemporal dendritic integration rule with an idealized two-compartment passive cable model, in which a spherical soma is connected to an unbranched cylindrical dendrite with finite length l and diameter d. The distance between a dendritic location and the soma is denoted by x. Given an excitatory input at location x~x E and at time t~t E , and an inhibitory input at location x~x I and at time t~t I , the membrane potential dynamics is governed by where u is the membrane potential with respect to the resting potential on the dendrite, c is the membrane capacitance per unit area, r a is the axial resistivity, and L is the leak conductance per unit area. The excitatory and inhibitory input strengths f E and f I control the amplitude of EPSP and IPSP, respectively. The variables E and I denote the unitary excitatory and inhibitory conductances per unit area with their peak value normalized to unity [described by Equation (29)]. e E and e I are the corresponding reversal potentials, respectively. The values of the parameters in the model can be found in the section of Materials and Methods.
By assuming that one dendritic end is sealed and the other dendritic end connects to the spherical soma, the boundary conditions are given by where S is the surface area of the soma. The initial condition is simply set as for a neuron at its resting state. For the physiological regime, the corresponding synaptic input strengths, f E and f I , are relatively small in the model [Equations (3)- (6)]. To be specific, for the amplitude of an EPSP less than 6mV and the amplitude of an IPSP less than {3mV , the corresponding strengths f E and f I are considered to be small. Therefore, we can represent the membrane potential as an asymptotic series in powers of f E and f I as follows: Substituting Equation (7) into Equations (3)-(6), order by order, we can obtain its asymptotic solutions u mn (x,t). The solutions to the second order (mznƒ2) are described below (see the section of Materials and Methods for a detailed calculation). For the zeroth order, we have The solution corresponds to the fact that the membrane potential response remains at its resting state when there is no stimuli presented (f E~0 , f I~0 ). The first order excitation O(f E ) is where 'Ã' denotes convolution in time, and G(x,y,t) is the Green's function of the cable equation given a d impulse input (the analytical expression is presented in the section of Materials and Methods). Note that the input e E E at x E can be viewed as the synaptic current when the local membrane potential is maintained at the resting state. The second order of excitation where the expression of u 10 (x E ,t) is given in Equation (9). Note that the input { E u 10 at x E can be viewed as the synaptic current, which is the product of the local conductance E with the first order local membrane potential u 10 . Similarly, we have the first and second order inhibition O(f I ) and O(f 2 I ): For the order of O(f E f I ), we have where u 10 (x,t) and u 01 (x,t) are given by Equations (9) and (11), respectively. On account of the fact that e E~7 0mV (relative to the resting potential) is nearly an order of magnitude larger than De I D~10mV , u 11 in Equation (13) can be further simplified as which indicates that the nonlinear integration effect mainly originates from the outward synaptic current, i.e., I u 10 , induced by the first order EPSP measured at the inhibitory input site x I , i.e., u 10 (x I ,t).
Numerical simulation of the cable model indicates that the second order asymptotic approximation is sufficiently accurate in capturing the model's solution of physiologically realistic membrane potentials, as demonstrated in Fig. 1. Therefore, if only an individual excitatory input is given [f I~0 in Equation (3)], the corresponding EPSP measured at the soma, denoted by V E , can be approximated by Similarly, if only an individual inhibitory input is given [f E~0 in Equation (3)], the corresponding IPSP measured at the soma, denoted by V I , can be approximated by If both the excitatory and inhibitory inputs are given at x E and x I , the corresponding SSP measured at the soma, denoted by V S , can be approximated by We define the difference between the SSP, and the linear summation of EPSP and IPSP as From Equations (15), (16) and (17), we have If e I is set to be at the resting potential, we can show that the value of V I indicating hyperpolarization vanishes, while the value of V SC stays nearly the same. Therefore, V SC is mainly caused by the shunting effect and thus is referred to as the shunting component (SC). In our analysis, the SC is the leading order of the nonlinear integration between excitation and inhibition. If we further define the shunting coefficient k EI as then k EI is nearly independent of the amplitude of EPSP and IPSP, because the input strengths f E and f I , which determine the amplitudes of EPSP and IPSP, cancel each other in both denominator and numerator in Equation (20b).
From Equation (20a), we have the following spatiotemporal dendritic integration rule The shunting coefficient k EI is nearly independent of the amplitude of EPSP and IPSP. In addition, k EI depends on the location of excitatory and inhibitory inputs x E and x I . For a fixed pair of excitatory and inhibitory input locations, k EI is a function of both time t and the arrival time difference between the excitatory and inhibitory input t~t E {t I , as illustrated in Fig. 2.
Validation of the Rule. As the bilinear integration rule (21) is derived from the idealized passive neuron model, we need to investigate its validity for a realistic neuron, which has active ion channels embedded in its tree-like dendrites.
We first perform the simulation of a biologically realistic pyramidal neuron with active channels (morphology shown in Fig. 2). The details of the model and the related computational method can be found in the section of Materials and Methods. The simulation results are summarized below.
For the case of concurrent inputs (Here, we use ''concurrent'' to denote the case when t E~tI ), as shown in Fig. 3A, when the excitatory and inhibitory inputs are elicited concurrently at different locations on the dendritic trunk, the SSP is found to be always smaller than the linear sum of the EPSP and the IPSP when elicited separately. In this case, the bilinear integration rule (21) holds at the time t Ã when the EPSP reaches its peak value. We can vary f E to control the amplitude of EPSP (less than 6mV ) and vary f I to control the amplitude of IPSP (less than {3mV ). For fixed input strengths f E and f I , we obtain the set of time courses of the EPSP, the IPSP, and the corresponding SSP. Using 9 sets of such data with different input strengths, we find that the SC amplitude V SC (t Ã ) depends linearly on the product of the EPSP and IPSP amplitudes, i.e., V E (t Ã )V I (t Ã ). The excellent linear fitting in Fig. 4A shows that the slope k EI (t Ã ) is independent of the amplitude of EPSP and IPSP. This result is consistent with the experimental observation [3].  Here t Ã denotes the time when EPSP reaches its peak value. (B) the IPSP is elicited 20ms before the EPSP. The results are obtained in the realistic pyramidal neuron model simulation which is described in detail in the section of Materials and Methods. The excitatory input is given at the location x E~2 83mm and the inhibitory input is given at the location For the same case of concurrent inputs, we then calculate k EI at any time t, instead of the peak time t Ã . We compute k EI (t) within the interval t Ã {svtvt Ã zs, where s~10ms. We choose this 20ms interval because the amplitude of the EPSP is relatively large within that interval. As a consequence, we can avoid a small denominator in calculating k EI (t) and improve the numerical accuracy. At each fixed time t, k EI (t) is estimated by linear regression using the same 9 sets of data used to validate the bilinear integration rule at t Ã . As demonstrated in Fig. 4B, at each time t, there is a small error bar for the slope k EI estimation and the R 2 value is very close to 1. Both facts indicate an excellent linear fitting of V E (t)V I (t) vs. V SC (t). Therefore, the shunting coefficient k EI is nearly independent of the amplitude of EPSP (V E ) and IPSP (V I ) at any time. As expected, the error bar for the slope estimation increases dramatically far away from the peak time due to the fact that the numerical accuracy is low when EPSP and IPSP are small, in particular, when they approach zero. However, the SC amplitude in this case is sufficiently small, and can thus be neglected. Therefore, the bilinear integration rule can naturally be considered valid with k EI~0 mV {1 .
For the case of nonconcurrent inputs, when the onset of the inhibitory input is 20ms earlier than that of the excitatory input ( Fig. 3B), our numerical results show that the bilinear integration rule (21) still holds, as shown in Fig. 4C-D. The rule is also confirmed for any excitatory and inhibitory input locations arbitrarily distributed on the dendritic trees.
We next perform electrophysiological experiments to validate our bilinear integration rule (21). The details of the experimental procedure can be found in the section of Materials and Methods.
In experiments, the excitatory input is given at *100mm with the EPSP amplitude less than 8mV and the inhibitory input is given at *50mm with the IPSP amplitude less than {3mV. For the case of concurrent inputs, we found that V SC depends linearly on V E V I at the time when EPSP reaches its peak, as shown in Fig. 4E, and at a non-peak time, as shown in Fig. 4F. Therefore, the slope k EI is nearly independent of EPSP and IPSP amplitudes. For the case of nonconcurrent inputs, when the IPSP is elicited 20ms earlier than the EPSP, the linear relationship between V SC and V E V I still holds, as shown in Fig. 4G-H, except that the value R 2 of the regression is smaller (0.77 to 0.81) than those in the concurrent case (0.90 to 0.99). Therefore, it can be seen from the above that, the bilinear integration rule (21) is confirmed in rat hippocampal CA1 pyramidal neurons.
Note that the bilinear integration rule is derived from an idealized passive neuron case. Interestingly, our results from the simulation and experiments demonstrate that the structure of the rule is preserved in the presence of both active channels and dendritic branches.
Spatial dependence of k EI . Although we have obtained the form of the bilinear integration rule (21), how the value of the shunting coefficient k EI depends on the input location is difficult to analyze directly from its explicit form [Equation (20b)]. Here we investigate the spatial dependence of input location for k EI , which may partially reveal the way in which k EI encodes spatial integration information.
In our previous study [23], a spatial rule for k EI as a function of input locations has been proposed based on a theoretical analysis of a multi-compartment cable model. Under the situation when the excitatory and inhibitory inputs are given concurrently, the spatial profile of k EI at the time when EPSP reaches its peak is characterized by the following spatial rule: For a fixed inhibitory input location, the I path (marked by green on the dendritic trees in Fig. 5) is defined as the path between the soma and the inhibitory input; Along the I path, k EI is predicted to increase as the distance between excitatory input location and the soma increases; for any branch (including the trunk) connecting to the I path, k EI is predicted to be constant for all excitatory input sites on the branch.
The prediction of this spatial k EI rule is consistent with our electrophysiological experimental results as shown in Fig. 5 (also see Ref. [3]). We next use the spatial k EI rule to explain these experimental results.
In Fig. 5A and 5B, an inhibitory input is given on an oblique branch and two excitatory inputs are given at two locations on the same oblique branch. The corresponding shunting coefficients k EI were estimated based on Equation (20a). In our experiment, no significant difference was found between the values of k EI when the two excitatory locations were both distal, as shown in Fig. 5A. This experimental observation can be easily understood by our spatial rule for k EI because the two excitatory inputs are on the same branch connected to the I path. By our rule above, the shunting coefficient k EI should be the same on such a branch.
In contrast, as shown in Fig. 5B, for two excitatory inputs at proximal locations, k EI was experimentally found significantly smaller for the excitatory input closer to the soma. This is consistent with our rule because k EI is predicted to be an increasing function of the distance between the excitatory input location and the soma for this case.
In Fig. 5C, an inhibitory input is given on the apical trunk and two distal excitatory inputs are given at either the trunk or a branch. For this case, the k EI values were found to be nearly constant in our experiment. This is the case in which the two branches where the two excitatory inputs are located connect to the I path with the same branching point. Therefore, this experimental observation can be understood through our rule that all k EI on the two branches are the same as the one at the branching point.
In Fig. 5D, an inhibitory input is given at an oblique branch and two distal excitatory inputs are given at different branches. For this case, the value of k EI for inhibitory and excitatory inputs located at the same branch was found in our experiment to be significantly larger than the case in which inhibitory and excitatory inputs are located at different branches. For the case when inhibitory and excitatory inputs located at different branches, by our rule, k EI equals to the value at the branching point on the I path. Because this branching point is closer to the soma than the other excitatory input location on the I path, the increase of k EI along the I path predicted by our rule explains the experimental observation.
Our theory [23] further predicts that, in principle, it is possible for an inhibitory input located at a branch to shunt the excitatory input on other branches with a large value of k EI . This generalizes the conclusion in Ref. [3], in which the shunting inhibition is mainly confined within the branch where the inhibitory input is located. In Ref. [3], the simulation investigation is focused on the case when inhibitory input is elicited on a branch close to soma. In this case, k EI for excitatory inputs on different branches is found to be small in simulation. This observation is consistent with our rule. According to our rule, the value of k EI for an excitatory input on the branch is equal to that for excitatory input on the branching point connected to the I path. Due to the branching point being very close to the soma in their simulation, k EI can be small. In general, when the inhibitory input is elicited at a location far from the soma, for the excitatory input on some other distal dendritic branch, the distance between the branching point and soma can be large, which leads to a large value of k EI .

Bilinear Rules for E-E & I-I Integration
So far we have addressed the dendritic integration for a pair of excitatory and inhibitory inputs. A natural question arises: how does a neuron integrate a pair of time-dependent synaptic conductance inputs with identical type?
The dendritic integration of excitatory inputs has been extensively investigated in experiments (reviewed in Ref. [1]), yet a precise quantitative characterization is still lacking. According to our idealized cable model, given a pair of excitatory inputs with input strengths f E1 and f E2 at locations x~x E1 and x~x E2 and at times t~t E1 and t~t E2 , the dynamics of the membrane potential on the dendrite is governed by the following equation: with the initial and boundary conditions the same as given in Equations (4)- (6). Similarly, we can represent its solution as an The data marked by grey squares were collected from 7 neurons in our experiments, and lines connect data from the same neuron. The data marked by black squares are the average of the data marked by grey squares. The error bar indicates one standard deviation. In all figure panels, the locations of the inhibitory input (I) and excitatory inputs (E1 and E2) are marked by a blue dot and red dots, respectively. The I path is marked by green. (A) The inhibitory input I at an oblique branch: k EI is nearly constant for two distal E1 and E2 on the same branch. (B) As in (A) except that E1 and E2 are more proximal than I. k EI is significantly different at E1 and E2 sites. (C) The inhibitory input I at the trunk: k EI is nearly constant between E1 at the trunk and E2 at the oblique branch. (D) The inhibitory input I at an oblique branch: k EI is significantly different between E1 and E2, where E1 is on the same branch as I and E2 is on a different branch. doi:10.1371/journal.pcbi.1004014.g005 Spatiotemporal Integration Rule of Synaptic Inputs PLOS Computational Biology | www.ploscompbiol.org asymptotic series and solve it order by order to obtain the following bilinear integration rule: where V E1 and V E2 are EPSPs induced by two individual excitatory inputs, and V S is the SSP when the two excitatory inputs are present. Similar to the case of a pair of excitatory and inhibitory inputs, the shunting coefficient k EE (t) only depends on the excitatory input locations and the input time difference. It does not depend on the EPSPs' amplitudes. Here k EE will still be referred to as a shunting coefficient because the origin of the nonlinear integration for the paired excitatory inputs is exactly the same as that for the paired excitatory and inhibitory inputs from the passive cable model.
The bilinear integration rule (23) is found to be consistent with the numerical results obtained using the same realistic pyramidal neuron model as the one used in the section of Bilinear Rule for E-I Integration. For a pair of excitatory inputs with their locations fixed on the dendritic trunk, the rule holds when the amplitude of each EPSP is less than 2mV . For the case of concurrent inputs, at the time t Ã when one of the EPSPs reaches its peak value is found to be linearly dependent of V E1 (t Ã )V E2 (t Ã ), as shown in Fig. 6A. This linear relationship indicates k EE is independent of the amplitudes of the two EPSPs. In addition, as shown in Fig. 6B, the bilinear integration rule is numerically verified in the time interval t Ã {svtvt Ã zs, for s~10ms, within which the amplitude of EPSPs are relatively large. For the case of nonconcurrent inputs, the bilinear integration rule is also numerically verified in the same way, as shown in Fig. 6C-D.
In addition, we find that when the input strengths become sufficiently strong so as to make the depolarized membrane potential too large, i.e. V E1 V E2 w5mV 2 , there is a deviation from the bilinear integration rule (23). This deviation can be ascribed to the voltage-gated ionic channel activities in our realistic pyramidal neuron model. After blocking the active channels, the rule becomes valid with a different value of k EE for large EPSPs amplitudes, as shown in Fig. 7. However, we note that, regardless of input strengths, the amplitude of SC is always two orders of magnitude smaller than the amplitude of SSP. Therefore, the integration of two excitatory inputs can be naturally approximated by the linear summation of two individual EPSPs, i.e. V SC~0 mV . We then perform electrophysiological experiments with a pair of excitatory synaptic inputs to confirm the linear summation. As expected, this linear summation is also observed in our experiments for both concurrent and nonconcurrent input cases, as shown in Fig. 6E and 6F, respectively. Note that, the linear summation is also consistent with experimental observations as reported in Ref. [24].
Similarly, for a pair of inhibitory inputs, we can arrive at the following bilinear integration rule from the cable model: where V I1 and V I2 are IPSPs induced by two individual inhibitory inputs, and V S is the SSP when the two inhibitory inputs are present. Here, k II (t) is the shunting coefficient that is independent of the IPSPs amplitudes but is dependent on the input time difference and input locations. The above bilinear integration rule (24) is consistent with our numerical results using the realistic pyramidal neuron model, as shown in Fig. 8A-D. Our electrophysiological experimental observations further confirm this rule, as shown in Fig. 8E-H.

Bilinear Rule for Multi-input Integration
In the previous sections, we have discussed the integration of a pair of synaptic inputs. In vivo, a neuron receives thousands of excitatory and inhibitory inputs from dendrites [2]. Therefore, we now address the question of whether the integration rule derived for a pair of synaptic inputs can be generalized to the case of multiple inputs.
Our theoretical analysis shows that, for multiple inputs, the SSP can be approximated by the linear sum of all individual EPSPs and IPSPs, plus the bilinear interactions between all the paired inputs with shunting coefficients k ij EI , k kl EE , and k mn II respectively (the superscript labels the synaptic inputs), i.e., We next validate the rule (25) using the realistic pyramidal neuron model. It has been reported that, for a CA1 neuron, inhibitory inputs are locally concentrated on the proximal dendrites while excitatory inputs are broadly distributed on the entire dendrites [25]. Based on this observation, we randomly choose 15 excitatory input locations and 5 inhibitory input locations on the model neuron's dendrites (Fig. 9A). In the simulation, all inputs are elicited starting randomly from 0ms to 100ms. In order to compare Equation (25) with the SSP simulated in the realistic neuron model, we first measure k EI , k EE , and k II pair by pair for all possible pairs. We then record all membrane potential traces V 1 E ,V 2 E ,:::,V 15 E and V 1 I ,V 2 I ,:::,V 5 I induced by the corresponding individual synaptic inputs. Our results show that the SSP measured from our simulation is indeed given by the bilinear integration rule (25), as shown in Fig. 9B and 9C. In contrast, the SSP in our numerical simulation deviates significantly from the linear summation of all individual EPSPs and IPSPs.

Graph Representation of Dendritic Integration
According to our bilinear integration rule (25), the dendritic integration of multiple synaptic inputs can be decomposed into the summation of all possible pairwise dendritic integration. Therefore, we can map dendritic computation in a dendritic tree onto a graph. Each dendritic site corresponds to a node in the graph and the corresponding shunting component is mapped to the weight of the edge connecting the two nodes. We refer to such a graph as a dendritic graph. The dendritic graph is an all-to-all connected graph if all stimuli are given concurrently (Fig. 10A). However, the dendritic integration for all possible pairs of synaptic inputs is usually not activated concurrently in realistic situations. For instance, if the arrival time difference between two inputs is sufficiently large, there is no interaction between them. The activated level of the nonlinear dendritic integration for a pair of synaptic inputs can be quantified by the SC amplitude-the weight of the edge in the graph. The simulation result shows that the number of activated edges at any time is relatively small on the dendritic graph (Fig. 10B-D), compared with the total number of edges on the all-to-all connected graph (Fig. 10A). Therefore, for the case of a hippocampal pyramidal neuron, the dendritic graph could be functionally sparse in time. The functional sparsity of a dendritic graph may also exist in neocortical pyramidal neurons. In vivo, a cortical pyramidal neuron receives about *10 4 synaptic inputs [26]. Most of them are from other cortical neurons [27,28], which typically fire about 10 spikes per second in awake animals [29,30]. Thus, the neuron can be expected to receive *10 5 synaptic inputs per second. The average number of synaptic inputs within 10ms (membrane potential time constants in vivo) is *10 3 . The number of activated dendritic integration pairs within the 10ms interval is *10 6 , which is relatively small compared with the total possible synaptic integration pairs *10 8 . Therefore, the activated integrations or edges in the dendritic graph within a short time window can be indeed functionally sparse (*10 {2 ). In general, the neuronal firing rates vary across different cell types, cortical regions, brain states and so on. Therefore, based on the above estimate, in an average sense, the graph of dendritic integration is functionally sparse.

Discussion
Our bilinear dendritic integration rule (21) is consistent with the rule previously reported [3], but is more general in the following aspects: (i) Our dendritic integration rule holds at any time and is not limited to the time when the EPSP reaches its peak value. (ii) The rule holds when the two inputs are even nonconcurrent. This situation often occurs because the excitatory and inhibitory inputs may not always arrive at precisely the same time. (iii) The form of the rule can be extended to describe the integration between a pair of excitatory inputs, a pair of inhibitory inputs, and even multiple inputs of mixed-types. The spatiotemporal information of synaptic inputs interaction is coded in the shunting coefficient, which is a function of the input locations and input arrival time difference.
Our bilinear integration rule holds in the subthreshold regime for a large range of membrane potential. When we derive the bilinear rule from the passive cable model, we assume that the   (25) are plotted against the SSP from the simulation of the realistic neuron model. Here, the data are points on the corresponding curves from ten trials sampled uniformly from 0ms to 100ms. For comparison, the slope of the grey line is unity. It can be observed that the red dots fall on the grey line. This indicates that the predicted SSP is equal to the simulated SSP at any time. doi:10.1371/journal.pcbi.1004014.g009 input strengths or the amplitudes of membrane potentials require to be small. This assumption forms the basis of the asymptotic analysis, because the second order asymptotic solutions of EPSP, IPSP and SSP converge to their exact solutions as the asymptotic parameters f E and f I (denoting the excitatory and inhibitory input strengths) approach zero. In general, in the passive cable model, the bilinear rule will be more accurate for small amplitudes of EPSPs and IPSPs than large amplitudes. Importantly, the assumption holds naturally that in the physiological regime when EPSP amplitude is less than 6mV and IPSP amplitude is less than -3mV, f E and f I are small O(10 {3 ). However, even for EPSP amplitude close to the threshold, i.e., 10mV, which is unusually large physiologically, we can show that the second order asymptotic solution can still well approximate the EPSP with a relative error less than 5%. Thus the bilinear rule is still valid for large depolarizations near the threshold. The validity of the bilinear rule for large membrane potentials is also confirmed in both simulations and experiments. In particular, in the analysis of our experimental data, to validate the bilinear rule, we have already included all the data when the EPSP amplitude is below and close to the threshold because we have only excluded those data corresponding to the case when a neuron fires.
Our bilinear dendritic integration rule (21) is derived from the passive cable model. However, the simulation results and the experimental observations demonstrate that the form of dendritic integration is preserved for active dendrites. Additional simulation results show that for the same input locations, the shunting coefficients are generally larger on the active dendrites than those on the passive dendrites with all active channels blocked. We also note that the value of k in simulation is different from the value measured in experiments. This difference may arise from the fact that some parameters of the passive membrane properties, such as the membrane leak conductance, may not be exactly the same as those in the biological neuron, and we have only used a limited set of ion channels in simulation compared with those in the biological neuron. In addition, the input locations in the simulation and the experiments are different, which may also contribute to this derivation. However, the bilinear form is a universal feature in both simulation and experiment.
By fixing excitatory input location while varying inhibitory input location, our model exhibits that there exists a region in the distal dendritic trunk within which the shunting inhibition can be more powerful, i.e, a larger k EI , than in proximal dendrites. This result is consistent with what is reported in Ref. [31]. Compared with Ref. [31], our work provides a different perspective of dendritic computation. In their work, the multiple inhibitory inputs can induce a global shunting effect on the dendrites. However, if we focus on the shunting effect only at the soma instead of the dendrites, our theory shows that all the interactions among multiple inputs can then be decomposed into pairwise interactions, as described by the bilinear integration rule (25). In addition, in this work, we focus on the somatic membrane potential that is directly related to the generation of an action potential. However, it is also important to investigate the local integration of membrane potentials measured at a dendritic site instead of that measured at the soma. Asymptotic analysis of the cable model can show that our bilinear integration rule is still valid for the description of the integration on the dendrites. On the dendrites, the broadly distributed dendritic spines with high neck resistances [32,33] will filter a postsynaptic potential to a few millivolts on a branch [34,35]. Within this regime our bilinear integration rule is valid. Note that our rule may fail to capture the supralinear integration of synaptic inputs measured on the dendrites during the generation of a dendritic spike [36]. However, if the integration is measured at the soma, our rule remains valid even when there is a dendritic spike induced by a strong excitatory input and an inhibitory synaptic input on different branches [3].
The bilinear integration rule (25) can help improve the computational efficiency in a simulation of neuronal network with dendritic structures. By our results, once the shunting coefficients for all pairs of input locations are measured, we can predict the neuronal response at the soma by the bilinear integration rule (25). By taking advantage of this, one can establish library-based algorithms to simulate the membrane potential dynamics of a biologically realistic neuron. An example of a library-based algorithm can be found in Ref. [37]. To be specific, based on the full simulation of a realistic neuron model, we can measure the timedependent shunting coefficient as a function of the arrival time difference and input locations for all possible pairs of synaptic inputs and record them in a library in advance. For a particular simulation task, given the specific synaptic inputs on the dendrites, we can then search the library for the corresponding shunting coefficients to compute the neuronal response according to the bilinear integration rule (25) directly. In such a computational framework, one can avoid directly solving partial differential equations that govern the spatiotemporal dynamics of dendrites and greatly reduces the computational cost for large-scale simulations of networks of neurons incorporating dendritic integration.

Ethics Statement
The animal-use protocol was approved by the Animal Management Committee of the State Key Laboratory of

The Cable Model
We consider an idealized passive neuron whose isotropic spherical soma is attached to an unbranched cylindric dendrite with finite length l and diameter d. Each small segment in the neuron can be viewed as an RC circuit with a constant capacitance and leak conductance density [11,38]. The current conservation within a segment ½x,xzDx on the dendrite leads to where u is the membrane potential with respect to the resting potential on the dendrite, c is the membrane capacitance per unit area, and L is the leak conductance per unit area. Here, I syn is the synaptic current given by: where G E and G I are excitatory and inhibitory synaptic conductance per unit area and e E and e I are their reversal potentials, respectively. When excitatory inputs are elicited at M E dendritic sites and inhibitory inputs are elicited at M I dendritic sites, we have where q~E,I. For a synaptic input of type q, f ij q is the input strength of the j th input at the i th location, t ij q is the arrival time of the j th input at the i th location, x i q is the i th input location. The unitary conductance is often modeled as with the peak value normalized to unity by the normalization factor N q , and with s qr and s qd as rise and decay time constants, respectively [38]. Here H(t) is a Heaviside function. The axial current I(x) can be derived based on the Ohm's law, where r a is the axial resistivity. Taking the limit Dx?0, Equation (26) becomes our unbranched dendritic cable model, In particular, for a pair of excitatory and inhibitory inputs with strength f E and f I received at x E and x I , and at time t E and t I , respectively, we have Similarly, for a pair of excitatory or inhibitory inputs with strengths f q1 and f q2 received at x q1 and x q2 , and at time t q1 and t q2 (q~E,I), respectively, we have For the boundary condition of the cable model [Equation (31)], we assume one end of the dendrite is sealed: For the other end connecting to the soma, which can also be modeled as an RC circuit, by the law of current conservation, we have where S is the somatic membrane area, and u s is the somatic membrane potential. The dendritic current flowing to the soma, I dend , takes the form of Equation (30) at x~0. Because the membrane potential is continuous at the connection point we arrive at the other boundary condition at x~0: For a resting neuron, the initial condition is simply set as Green's Function In the absence of synaptic inputs, Equation (31) is a linear system. Using a d impulse input, its Green's function G(x,y,t) can be obtained from with the following boundary conditions and initial condition, c LG(0,y,t) Lt~{ L G(0,y,t)z pd 2 4Sr a LG(x,y,t) Combining the two boundary conditions (B(g,s) is thus eliminated), we have where whose denominator is denoted as f(s) for later discussions. For the inverse Laplace transform, we need to deal with singular points that are given by the roots of f(s)~0. It can be easily verified that these singularities are simple poles and LH(j,g,s) is analytic at infinity. Then LH(j,g,s) can be written as where H n (j,g) is a constant coefficient in the complex s domain, and s~{k n are the singular points. Then taking the inverse Laplace transform of Equation (44), we obtain Now we only need to solve k n and H n (j,g) in Equation (45) to obtain the Green's function of Equation (40). We solve the singular points s~{k n first. Defining w n~{ i ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ffi {k n z L p l, f(s)~0 yields whose roots can be determined numerically. There are solutions for w n with (n{1=2)pvw n v(nz1=2)p for n §1 and w 0~0 : Next, to determine the factors H n (j,g), we use the residue theorem for integrals. For a contour C n that winds in the counter-clockwise direction around the pole s~{k n , and that does not include any other singular points, the integral of LH(j,g,s) on this contour is given by Using Equations (42)(43)(44) and (47), we obtain H n (j,g)~cD n cos w n (1{j=l) where D n~2 clzclw {1 n sin (w n ) cos (w n )z2 cos 2 (w n ) for n §0. The solution of the original Green's function for Equation (39) can now be expressed as

Asymptotic Analysis
We first consider the case when a pair of excitatory and inhibitory inputs are received by a neuron. Similar results can be obtained for a pair of excitatory inputs and a pair of inhibitory inputs. For the physiological regime (the amplitude of an EPSP being less than 6mV and the amplitude of an IPSP being less than {3mV ), the corresponding required input strengths f E and f I are relatively small. Therefore, given an excitatory input at location x~x E and time t~t E , and an inhibitory input at location x~x I and time t~t I , we represent u(x,t) as an asymptotic series in the powers of f E and f I , Substituting Equation (51) into the cable equation (31), order by order, we obtain a set of differential equations. For the zerothorder, we have Using the boundary and initial conditions [Equations (34), (37), and (38)], the solution is simply u 00~0 : Because u 10 is given by Equation (55), the solution of Equation (56) is Similarly, we can have the first and second order inhibitory solutions,

Hippocampal Slice Preparation and Electrophysiology
The experimental measurements of summation of EPSPs or IPSPs in single hippocampal CA1 pyramidal cells in the acute brain slice followed a method described in Ref. [3], with some modifications. A brief description of modified experimental procedure is as follows. Acute hippocampal slices (350mm thick) were prepared from Sprague Dawley rats (postnatal day 14-16), using a vibratome (VT1200, Leica). The slices were incubated at 34uC for 30 min before transferring to a recording chamber perfused with the aCSF solution (2ml/min; 30-32uC). The aCSF contained (in mM) 125 NaCl, 3 KCl, 2 CaCl 2 , 2 MgSO 4 , 1.25 NaH 2 PO 4 , 1.3 sodium ascorbate, 0.6 sodium pyruvate, 26 NaHCO 3 , and 11 D-glucose, and was saturated with gas containing 95% O 2 and 5% CO 2 (pH 7.4). Whole-cell recording was made from the soma of CA1 pyramidal cells using glass micropipettes under an upright microscope (BX51WI, Olympus) equipped with the DIC optics and an infrared camera (IR-1000E, DAGE-MTI). The intra-micropipette solution contained (in mM) 145 K-gluconate, 5 KCl, 10 HEPES, 10 disodium phosphocreatine, 4 Mg 2 ATP, 0.3 Na 2 GTP, and 0.2 EGTA (pH 7.3), together with fluorescent dye Alexa Fluor 488 (20mM, Invitrogen) to visualize the dendritic trees. Pipette resistance was about 3-4 MV, and the access resistance during the whole-cell recording was normally less than 20 MV. The same method for microiontophoretic application of extracelluar glutamate or GABA at the apical dendrite of CA1 pyramidal cells was used to elicit rapid membrane depolarizations (EPSPs) and hyperpolarizations (IPSPs). For all three experimental configurations (EPSP-IPSP, EPSP-EPSP and IPSP-IPSP summation), two micro-iontophoretic pipettes were placed at dendritic locations 100mm and 50mm from the soma, respectively, in particular for the EPSP-IPSP summation GABA iontophoretic pipette was always placed at the more proximal location than glutamate iontophoretic pipette was placed. For each recorded cell, an electrode was placed at the soma to set the resting membrane potential to about {60mV in order to obtain a driving force of 15{20mV for inhibitory GABA inputs. Electrical signals of individual and summed iontophoretic responses were amplified and filtered at 3 kHz (low pass) by a patch clamp amplifier (MultiClamp 700B, Molecular Devices), digitalized (100 kHz) by an AD-DA converter (Digidata 1440A, Molecular Devices), and acquired by a pClamp 10.3 (Molecular Devices) into a computer for further analysis.

Data Processing
In order to study the dendritic integration of a pair of excitatory and inhibitory inputs, for fixed input locations and strengths, a moving average technique with time lag 1 ms was first applied to smooth each individual trace of the EPSP, IPSP, and SSP recorded in our experiments. After smoothing, we measured the amplitudes of EPSP, IPSP, and SSP at different times, including those when EPSP reached its peak value, and denoted them by V E , V I , and V S , respectively. By varying the excitatory and inhibitory input strengths, we measured values of V E , V I and V S . We then constructed a scatter plot of {V E V I vs. {V SC~{ (V S {V E {V I ) at different times. We divided the range of {V E V I into approximately 10 bins and averaged all the data points ({V E V I ,{V SC ) within each bin. The number of bins was chosen to ensure at least 8 data points were used for averaging. However, the qualitative results were not sensitive to the number of bins (e.g, from 6 bins to 16 bins). Using the Curve Fitting Toolbox in Matlab Version 7.14, we finally fitted the averaged data points by a linear function V SC~kEI V E V I , from which the slope k EI was estimated together with its 95% confidence interval. For the dendritic integration of a pair of identical type, the same data processing procedure was followed.