Lactate supply overtakes glucose when neural computational and cognitive loads scale up

Significance Brain activity and performance are constrained by energy metabolism. Glucose and lactate have been proposed as energy substrates for neurons allocated to various forms of memory. We demonstrate that glucose and lactate metabolism are differentially engaged in neuronal fueling depending on the activity-dependent plasticity and behavioral complexity. These results reconcile a decades-long debate regarding the nature of the neuroenergetic sources used by synaptic activity with the potential of inspiring future lines of research regarding neuroenergetic rules. The brain has high energy demands, and alterations in neuroenergetics are hallmarks of several neuropathologies. A better knowledge of the cellular and molecular mechanisms of neuroenergetics, as reported here, may be instructive in targeting energy metabolism deficits as a therapeutic approach for neurodegenerative diseases.

was placed back into the arena for object discrimination and were exposed to the exchanged position of two of the presented objects, C-B-A-D for half of the animals and A-D-C-B for another half, for 10 min. NOR and OiP performance analysis. Rats were videotaped during familiarization and test sessions. Videos were analyzed and time spent on exploration each object (sniffing or licking) was measured during familiarization and test sessions. For assessing NOR and OiP memory performance, the calculation of object discrimination, the exploration time of the novel object was expressed as percentage of the total exploration time of familiar and novel objects during familiarization and test sessions. Then, we calculated the relative time of exploration per object and the preference index (%) in control (saline-injected rats) and in oxamate (Na-oxamate-injected rats).
In vivo electrophysiology in behaving rats during NOR and OiP task e-fEPSP were measured in the left CA1 over the 3-day behavioral assessment in rats subjected to NOR (familiarization: A-A and test: A-B) or OiP (familiarization: A-B-C-D and test: C-B-A-D) tasks. A recording wire (stainless steel, 0.005" diameter) was stereotaxically implanted under pentobarbital-ketamine anesthesia (pentobarbital: 30 mg/kg i.p., Ceva Santé Animale, Libourne, France; ketamine: 27.5 mg/kg, i.m., Imalgène, Mérial, Lyon, France) in CA1 stratum radiatum along the guide cannula (anteroposterior: -4.1mm, mediolateral: 2.5, dorsoventral : 3.3), and a bipolar stimulating electrode (2 twisted wires, same as the recording electrode) in the ipsilateral Schaffer collaterals (anteroposterior: -4.4mm, mediolateral: 4.4, dorsoventral : 3.4). The rat body temperature was maintained using at 37°C during surgery. Rats were allowed at least a week of recovery before electrophysiological recordings and behavioral testing. Recordings were performed during 10-15 min in a 24x44cm Plexiglas box, 60 min before familiarization (before the oxamate or saline injection) and, 2 and 24 hours after familiarization. Schaffer e-fEPSPs were amplified, and acquired at 20 kHz using a KJE-1001 system (Amplipex, Szeged, Hungary). Test pulses (150-200 μs, 20-450 µA) were evoked every 30 seconds using a square pulse stimulator and stimulus isolator (Model 2100, AM-Systems, Sequim, WA, USA). Responses were analyzed offline using custom Matlab codes (2019b, The Mathworks, Natick, MA, USA) scripts. Traces was smoothed (linear averaging) over 0.15, 0.25 or 0.5 ms depending on the noise level. Detection of local extrema was performed to define the start and peak of the response, and to measure the e-fEPSP amplitude. Statistical analyses (two-tailed Student t-test) and calculation of plasticity ratio were made on 25 trials. The average baseline of e-fEPSP amplitude obtained at t=-60 min was used to normalize every individual e-fEPSP amplitude at +2 and +24 hours, which were then averaged to obtain the plasticity ratio. The F-statistic and p-value associated with the average vectors were obtained from MANCOVA test (two dependent variables: preference index and plasticity ratio; one factor: belonging to the saline-or oxamate-injected rat group). Pearson's correlation coefficients were computed for paired plasticity ratios measured at 2 or 24 hours after familiarization phase.

Histology
After the completion of NOR and OiP experiments, cannulas and recording and stimulation electrodes positions were examined. Rats were anaesthetized (sodium pentobarbital, 150mg⁄kg i.p.) and transcardially perfused with saline followed by 4% ice-cold paraformaldehyde in 0.1 M phosphate buffer. Brains were extracted, post-fixed in PFA (4% in PBS) for 2 days at 4°C and then cryoprotected in 30% sucrose solution for one week. Brains were cut in coronal 50 µm sections using a cryotome (HM400, Microm Microtech, Francheville, France) and slices were maintained in 0.1 M potassium-PBS (pH=7.4). The sections were Nisslstained with thionin (Neurotrace 500/525, Thermofischer, Waltham, MA, USA), and images were acquired using a stereozoom microscope (Axiozoom, Zeiss, Oberkochen, Germany) and processed in ImageJ. Only the rats with cannula tips located bilaterally within the dorsal hippocampi were included in the data analysis.

Mathematical model
In the mathematical model, though, parameter calibration imposes to specify which cell type one considers. The vast majority of published models and data on the subject is specific of astrocytes; we opted for astrocytes as glial cells in the model. Our model therefore simulates the network of signaling and metabolic reactions occurring in a postsynaptic neuronal terminal and an interacting astrocyte shown in Figure 2A.
Stimulations. The membrane potential of the postsynaptic compartment " is given by: with the indicator function 1 [6,I] ( ) = 1 if ∈ [ , ], 0 otherwise , 2 H is the time of the beginning of i th depolarization of the postsynaptic neuron, max is the amplitude of the depolarization current in the postsynaptic neuron and dur its duration. step is the time step of the step current, is the attenuation of the bAP at the postsynaptic compartment and the models for the ionic currents 0123 , 3143 , 56+ ( Fig. 2A)  where ⌊ ⌋ denotes the integer part of , mod( , ) = − ⌊ / ⌋ is its modulo j and = mod( − 1,40). On the postsynapse, experiments show that during a presynaptic TBS, the probability for a presynaptic spike to trigger a postsynaptic bAP is circa 0.5, therefore we set max = 0 (no postsynaptic current injected in TBS) and set max = 0 for odd s. STDP. The variation of postsynaptic cytosolic calcium concentration is computed as: where " is the concentration of ATP in the postsynaptic compartment, max , max and max are, respectively, the amplitudes of the calcium-gated LTD and LTP and of the ATP-gated depotentiation and start , start and Thr their respective thresholds. * sets the value of the unstable (intermediate) steadystate of the bistable and € the time scale to reach the two stables steady-states = 0 and = 1. Finally, the synaptic weight is taken an affine function of the signaling state variable: Astrocyte-neuron lactate shuttle. To derive the time course of " in eq. (6) as a function of the pre-and post-synaptic stimulations, we used the model developed by Jolivet et al (2015) (4). In both the postsynaptic terminal and a nearby interacting astrocyte ( Fig. 2A), the model accounts for glycolysis, i.e. the production of pyruvate from glucose with glyceraldehyde-3 phosphate (GAP) and phosphoenolpyruvate (PEP) as intermediates as well as pyruvate formation from lactate by lactate dehydrogenase (LDH). The model also incorporates the production of NADH from pyruvate by the TCA cycle in the mitochondria, resulting ATP production by the electron transport chain as well as NADH shuttling from the cytosol to the mitochondria. Glucose and lactate are exchanged between the cytosol of the two compartments and the periplasmic extracellular medium via transporters, including MCT2 for lactate transport from the periplasmic volume to the postsynaptic compartment. Diffusive transport also occurs between the periplasmic volume and reservoir solutions with fixed concentrations (bath solution or blood capillaries of the slices). In addition, the model accounts for the activity-dependent dynamics of sodium ions in both compartments. Presynaptic stimulations trigger Na influx via EEAT2 channels in the astrocyte and voltage-gated sodium channels (VGSC) and (partly) AMPA receptors in the postsynaptic compartment. Sodium is then pumped back by Na,K-ATPases (ATPase), consuming ATP in the process.
We have implemented the model of Jolivet et al (2015) (4) in extenso. Our only change has been to substitute the Hodgkin-Huxley equation used for the postsynaptic membrane voltage in (30) by our eq. (1). This made us slightly adapt the maximum conductance for the postsynaptic VGSC, the strength of the effect of the presynaptic stimulation on the astrocyte and the leak conductance. All the remaining of the model, i.e. the 28 ODEs and the 80+ remaining parameters were taken unchanged from (4).

Parameter estimation.
Coupling the two models (STDP + astrocyte-neuron lactate shuttle) left us with 27 parameters to estimate (SI Appendix Table S2). We used a subset of our experimental data (training set), whereas validation was carried out by checking the accuracy of the model output for experimental conditions that were not used in parameter estimation (pharmacology perturbation experiments, changes of extracellular glucose concentrations, different stimulation protocols). See Supplementary Information for a complete description of the model as well as parameter estimation strategy and values.

Code availability:
Computer code for the model is publicly available at https://gitlab.inria.fr/hberry/anls_stdp

Mathematical model
Ionic currents of the postsynaptic compartment membrane. The AMPA current was modeled as a double exponential conductance (5): where cde H is the time of the i th presynaptic stimulation, 3143max the maximal AMPA conductance and 3143[ and 3143a the two-time scales. The NMDAR was also modelled with a double exponential (6): where the magnesium block = Z1 + -Mg aE š/3.57exp(−0.062 " )\^[. The VGCC current was modelled after (7) as a high-voltage L-type current: with activation and inactivation functions given by: The voltage-gated sodium current was taken from (4): where the steady-state activation was … = $ ( $ + $ ) ⁄ with $ = −0.1 ( " + 33) -exp]−0.1( " + 33)_ − 1š ⁄ and $ = 4exp(− ( " + 58) 12 ⁄ ) and the inactivation was given by: Note that 06 does not contribute to the evolution of the membrane potential in eq.(1) since eq.(1) does not include the spike currents but emulates the bAP using a simple increment of the potential at the spiking time followed by an exponential decay.
Astrocyte-neuron lactate shuttle. We implemented the model of Jolivet et al (2015) (4) where ® denotes the concentration of chemical species S in compartment = { , } for astrocyte or postsynaptic neuronal, respectively. We give below a complete depiction of the equations and parameters of the model, but interested readers should refer to the original paper for further specifics.
The cytosolic sodium concentration is given by with the Nernst potential 06® = ⁄ log( · / ® ), 06® the leak conductance and pump® the maximal rate of the Na-K-ATPase pump in compartment x. stim® is the sum of sodium influxes triggered by the stimulations. In the postsynaptic compartment: stim" = where N is the total number of presynaptic stimulations and f their frequency (4). Cytosolic glucose concentrations are given by where · is the extracellular concentration of glucose, i.e. in the pericellular volume, while Â is itsconstant -concentration in the large (reservoir) volume of the bath solution and/or the blood vessels. ®Å is the constant of glucose transport between compartments x and y, ½´4{´® is the maximal rate for the part of glycolysis from hexokinase to phosphofructokinase (lumped as a single equivalent reaction).
Glyceraldehyde-3 phosphate (GAP) concentrations are obtained with: where 4¹´® is the maximal rate for the part of glycolysis from GAP dehydrogenase to enolase, ® Â is the concentration of NADH in the cytosol of compartment x and N is the total concentration of NADH, i.e. the sum of NADH and NAD+ concentration in the cytosol. Phosphoenolpyruvate (PEP) concentration obeys: with 4´® the pyruvate kinase rate in x. Next, the dynamics of pyruvate concentration (PYR) is given by: Here +2½on® and +2½off® are the forward and reverse rates, respectively, of lactate dehydrogenase (LDH) in compartment x, ® $ is the concentration of NADH in the mitochondrion of x and mitoin® is the maximal rate of the TCA cycle. Lactate (LAC) dynamics is obtained via: where ®Å is the constant of lactate transport between compartments x and y and · is the extracellular concentration of lactate, i.e. in the pericellular volume. Likewise, in the astrocyte: where ®Å is the constant of lactate transport between compartments x and y and Â is the extracellular concentration of lactate in the reservoir volume. Likewise, the concentration of lactose in the pericellular extracellular medium is given by: and that of glucose: with ®Å the ratio between the volume of compartment x and that of compartment y. The concentration of cytosolic NADH ( ® Â ) is obtained by integration of: where 032½® is the maximal rate for NADH shuttling from cytosol to the mitochondria sub-compartment, is the relative mitochondria volume, ^® = where mitoout® is the maximal rate of the electron transport chain in x and a ® the oxygen concentration in this compartment. ATP concentration in the postsynaptic compartment is given by: with 3"46Ú·Ú" a parameter accounting for ATPase activities outside Na-K-ATPases, CKon" and CKoff" the forward and backward rates of creatine kinase, respectively, and C the total concentration of creatine plus phosphocreatine.
" " ⁄ , the ratio between deoxyAMP and deoxyATP is computed with where 3´ is the adenylate kinase equilibrium constant, = ® + ® + ® is the total adenine nucleotide concentration and ® = 3´a + 4 3´( ® ⁄ − 1). Similarly, the concentration of ADP is computed from that of ATP using ® = 0.5 ® ]− 3´+ Ý ® _. Now, in the astrocyte, the concentration of ATP is given by: where phosphocreatine concentrations are given by Finally, oxygen concentrations are obtained through: where Â6Þ ® ⁄ is the oxygen transport constant from the reservoir to the cells and 0 a Â oxygen concentration in the reservoir.
Numerical integration. The system of ODEs consisting of eq.(1)-(7) and eq.(S1)-(S25) was integrated using a variable-order adaptive-stepsize stiff solver (ode15s in Matlab®), interrupting integration at every stimulation events ( if relevant) to account for the discontinuities of eq.(2-4). The system was first integrated for 2x10 5 seconds in the absence of any electrical stimulation to guaranty that the initial state (before electrical stimulation) is the stable steady-state. After electrical simulation (TBS or STDP), the system was integrated for 45 min. The state of the system 45 min after the stimulation sets the value of the synaptic weight (eq.7) after the stimulation, thus the potential expression of LTP.
Parameter estimation. To calibrate the model, we used the following strategy. All the parameters of equations (S3) to (S25) kept the values set by (4) (SI Appendix Table S2) except for 06max (eq. S6), L (eq.1) and the prefactor of eq.(S10), that were modified as explained below. Therefore, parameter estimation was almost completely restricted to the parameters of eq. (1-7). Together, this represented 27 parameters to estimate (SI Appendix Table S2). Parameter values were estimated on a subset of our available experimental data (see below) whereas validation was carried out by checking the accuracy of the model output for experimental conditions that were not used in parameter estimations.

1.
The parameters of the postsynaptic stimulations max , step , amp and as well as the delays [ and a and the leak conductance L were fitted to the experimental traces of the postsynaptic membrane potential as measured in the soma (SI Appendix Fig. S2B). The attenuation factor of the bAP between the soma and the postsynaptic compartment, AT, was set to a value that roughly corresponds to a synapse located at middistance between the soma and the dendrite end according to (8) (see their Fig 4D). 2. 3143max , 3143[ and 3143a were set to yield EPSPs of 2 mV amplitude, with short onset and a decay time scale around 10 ms, as measured in dendrites by (9). 0123[ and 0123a were set to yield a NMDA-component for the calcium influx that rises fast and takes roughly 200 ms to get back to zero and 0123max was fixed so that each presynaptic spike increases the cytosolic calcium level by 0.17 mM at -70 mV (10). 56+max was set so that one bAP (on top of the depolarizing current) triggers a calcium influx of 400 nM amplitude at the synapse in the absence of presynaptic stimulation (11). We also checked that with these parameter values, the amplitude of the calcium peak triggered by a postsynaptic stimulation comprising two bAPs is indeed roughly twice the amplitude obtained with a single bAP (SI Appendix Fig. S2B) as measured in our experimental setup both in spines and shafts (SI Appendix Fig. S6).

3.
In the original model (30), the membrane voltage of the presynaptic compartment is modelled as a Hodgkin-Huxley equation with variable Nernst potentials for Na. Because of the variable Nernst potential, this model exhibits strong spike-frequency adaptation so the postsynaptic neuron quickly ceases to emit spikes after the onset of the stimulation. Since our membrane voltage for the postsynaptic compartment does not exhibit spike-frequency adaptation, we had to adapt a pair of parameters to guaranty that the electrical stimulations employed in the experiments used to calibrate the model by (4) (their Fig. 4A, with experimental data taken from (11) will still yield the correct time course in our model. We therefore adapted the values of 06max and the prefactor of eq.(S10) so that the stimulation used in (11) (presynaptic stimulations for 20 sec) yielded in our model the same time-courses for 6 Â and " $ as their experimental measurements (their Fig. 4D). In particular (SI Appendix Fig. S2C) this stimulation leads to 1) a dip of " $ of around -10%, peaking around the end of the effective stimulation, and converging back to baseline after 15-20 min, and 2) a delayed overshoot of 6 Â that peaks at approx. +8% roughly 20 min after the end of the stimulation. Those dynamics reproduce previous experimental measurements (33), in both quantitative and qualitative terms.

4.
The parameters related to the dynamics of the synaptic weight were set as follows. The calcium thresholds start and start were set so that a "standard" STDP protocol (1 bAP at 1 Hz) triggers LTP for positive spike timings that are not larger than 30 ms and for more than approx. 20 pairings (SI Appendix Fig.  S2D). The values of max , max , * and € were fitted on the experimental measurements of the time course of the synaptic weight change for a STDP protocol with 1 bAP, 50 x at 0.5 Hz (Fig. 2B). Note that with 5-TBS, the experimental measurement of the final amplitude of the synaptic weight change was larger than with STDP (Fig. 2B), we thus adapted the value of € for 5-TBS (but kept the values of max , max and * to those obtained by fitting on STDP 0.5 Hz 50x).

5.
Finally, we fixed the ATP threshold Thr to a value that allows discriminating between STDP 0.5 Hz 50 pairings (LTP) stimulations and 5-TBS with oxamate (no LTP). In particular, the astrocyte-neuron lactate shuttle model predicts that neuronal ATP, " , falls below 2 mM for 5-TBS in the presence of oxamate (Fig.  2C), but remains well above 2 mM for 5-TBS in control conditions and STDP 50 pairings at 0.5 Hz (in control and with oxamate). In the absence of experimental data to set the value of max , we used a value large enough to cancel LTP for 5-TBS in the presence of oxamate.
Simulation of pharmacological experiments. The action of pharmacological agents was emulated by changing the corresponding parameters to the following values:    Figure 2). The model is initiated at time t=0 in control conditions. In control (t< 50s), ATP makes up roughly the totality of the 4 mM total adenosine phosphate concentration (top), the stationary level of pyruvate (Pyrn) and cytosolic NADH (NADHcn) in the neuron is large (middle), and neuronal glycolysis is low, as witnessed by the low activity of pyruvate kinase, PK (bottom). Oxamate addition in the neuron cytosol at t=50 sec rapidly switches the neuron to an oxidized redox state with a close to total depletion of cytosolic NADH. As a result, the ATP level drops well below 4 mM. Oxamate addition also results in neuronal glycolysis, as illustrated by the increase of pyruvate kinase activity. This restores significant levels of neuronal pyruvate and stabilizes ATP to roughly 50% of the total adenosine phosphate, i.e. slightly above 2 mM. Figure 2). Averaged time-course of the synaptic weight and EPSC amplitudes 45-55 min after STDP. When intracellular solution contained 2 mM of ATP and 5 mM phosphocreatine, LTP was not observed after STDP protocol in control or i-oxamate conditions. Only one significant LTP could be induced in control and i-oxamate conditions out of 6 and 7 cells, respectively. All data: mean±SEM. ns: not significant by two tailed t-test. See SI Appendix Table S3 for detailed data and statistics.

Fig. S5. LDH inhibition prevents 5-TBS-LTP and left unaffected STDP-LTP in adult mice.
Averaged time-course of the synaptic weight and EPSC amplitudes 50-60 min after 5-TBS or STDP in adult mice. Intracellular inhibition of LDH with i-oxamate show distinct effects on 5-TBS and STDP expression since 5-TBS did not induce plasticity whereas STDP triggered a potent LTP. Representative traces: 15 EPSCs averaged during baseline (grey) and 60 min (red) after protocols (arrows).All data: mean±SEM. *p<0.05; ns: not significant by two tailed t-test. See SI Appendix Table S3 for detailed data and statistics.  No difference for e-fEPSP amplitudes 24 hours after familiarization were observed between saline or oxamate-injected rats, with on average no plasticity for the NOR task (D) and a majority of LTD for the OiP task (C) (see average vectors in D and E). (F and G) Correlation between e-LFPs recorded 2 and 24 hours after familiarization in the NOR and OiP tasks, in relation with the preference index. NOR: significant positive correlation with the plasticity ratios decreasing from +2 to +24 hours, with larger LTPs remaining LTP at +24 hours, while smaller LTPs shifting into no plasticity or LTD (similar for both salineand oxamate-injected rats). OiP: no correlation between plasticity at +2 and +24 hours with LTP at +2 hours turned either into no plasticity or LTD while the majority of LTD at +2 hours persisted in LTD at +24 hours. All data are presented as mean±SD. *: p<0.05; **: p<0.01; ***: p<0.001; by two tailed t-test. See SI Appendix Tables S4C and S4D for detailed data and statistics.

5-TBS
STDP (50 pairings @ 0.5 Hz)  Basal Na-K-ATPase activity 0 (n) 0. 0687 (a) mM/s 1 : after compensation for the slight K + permeability of CaL channels, see 41 . 2 : n is for the postsynaptic neuronal compartment, a for the astrocytic one, e for the extracellular pericellular volume and c for constant reservoir concentrations.