Monte Carlo Simulation

All Monte Carlo simulations were performed using C programming routines written in-house for the two-compartment model shown in Figure 1A, which is identical to that published previously (Shen et al, 1999). Metabolic fluxes were assumed as follows (in μmol/g/min): ⁿV_TCA = 0.71, ^aV_TCA = 0.06, V_cyc = 0.32, V_ana = 0.04, V_dil(Lac) = 0.05, V_dil(Gln) = 0.14 (see Figure 1A for their definitions), which are the mean metabolic fluxes of resting human brain (Shen et al, 1999). Metabolic concentrations were assumed as follows (in mmol/L): [glutamate] ([Glu]) = 9.1, [glutamine] ([Gln]) = 4.3, based on previous ¹³C MRS measurements of human brain (Gruetter et al, 1994; Mason et al, 2003, 2007; Shen et al, 1999). The metabolic model consisted of two coupled differential equations which were solved numerically as described previously (Shen et al, 1999), yielding turnover curves of glutamate C4 and glutamine C4 as labeled from [1-¹³C]glucose.

The Monte Carlo simulation procedure was performed as follows: synthetic curves for glutamate C4 and glutamine C4 were first generated by solving the two coupled differential equations. Then, noise of normal distribution with a predefined standard deviation (σ) was added to the synthetic curves. At the early portion of the curves, if a resultant turnover data point was less than the natural abundance concentration because of addition of random noise, it was forced to take the corresponding natural abundance value. A total of 160 mins of infusion time was simulated with 32 data points per curve (equivalent to a temporal resolution of 5 mins in MRS data acquisition). Two noise levels were investigated, σ = 0.1 and 0.2 μmol/g, which spanned the range of noise levels observed within the experimental data sets. A total of 100 synthetic data sets for each simulation (i.e., different noise realization but with the same σ and other conditions) were generated using the Monte Carlo method. Each data set was then fitted using the metabolic model shown in Figure 1A to obtain the best fit values for each of the free metabolic fluxes. All metabolic fluxes were constrained between 0 and 10 μmol/g/min unless specified otherwise based on previous experiments and constraints. Minimization was performed using the simulated annealing algorithm (Shen et al, 1999) which is well-suited to finding the global minimum in a multidimensional error space (Kirkpatrick et al, 1983). The best fit was defined to be the combination of the free metabolic fluxes which minimizes the sum of squares of the difference between the Monte Carlo data point and the corresponding value calculated using the two-compartment model.

Calculated Values of V_cyc, the Impact of Noise, and Constraining Other Fluxes

Table 1 summarizes the results of Monte Carlo simulations performed for each predetermined condition (noise level and constraints) using the model depicted in Figure 1A which includes astroglial dilution. The constraints examined were on the rate of the glial TCA cycle (^aV_TCA) and on the rate of glutamine efflux (V_efflux) with V_efflux equal to 20% of the flux of glutamine synthesis (V_Gln). With no constraints and a noise level of σ = 0.1 μmol/g a mean V_cyc of 0.27 μmol/g/min was calculated with a standard deviation of 37%. At the higher noise level of σ = 0.2 μmol/g a lower value of V_cyc (0.16 μmol/g/min) was determined, as the fitting tends to calculate a very high value of ^aV_TCA. When ^aV_TCA was constrained to be < 0.1 μmol/g, the resulting value of V_cyc was 0.32 μmol/g/min for both noise levels and with a reduced standard deviation. An additional reduction in the standard deviation was obtained by also constraining the value for the net rate of glutamine efflux. Constraining the values for the dilution of the glutamate and glutamine C4 pools further reduced the standard deviation of the calculations.

V_dil(Gln) and the Sensitivity of ¹³C Turnover Curves to V_cyc

In the results of the Monte Carlo simulations listed in Table 1, V_dil(Gln) was shown to have a large influence on the sensitivity of the glutamine C4 turnover curve to the value of V_cyc. Incorporating both the V_dil(Gln) flux into the model and constraining its value significantly affected the relative standard deviation of V_cyc. The relationship between V_dil(Gln) and the sensitivity of the glutamine C4 turnover curve to the value of V_cyc may seem surprising at first. Intuitively, this relationship can be easily explained using a simplified version of the two-compartment model (Figure 1B), in which all ¹³C labels flow from pyruvate C3 to neuronal glutamate C4, which is in exchange with astroglial glutamine C4. This simplified model represents the dominant metabolic pathway for the trafficking of ¹³C labels in brain when [1-¹³C] or [1,6-¹³C₂]glucose infusion is used. For this simple model, the kinetics of ¹³C label incorporation into glutamine C4 is described by the following differential equation describing the bidirectional label dilution of glutamine C4:

where fe represents the fractional enrichment. The incoming astroglial dilution flux brings in ¹²C-glutamine C4 and does not explicitly appear in Equation (1). At isotopic steady state,

With V_dil(Gln) = 0, the fractional enrichment of the downstream glutamine C4 approaches that of the upstream glutamate C4 after the initial infusion period, thereby losing sensitivity to the label-feeding flux (V_cyc) with time and leading to an overall reduced sensitivity to V_cyc. At isotopic steady state, the fractional enrichments of glutamine C4 and glutamate C4 are equal for V_dil(Gln) = 0; thus, the isotopic steady state glutamine C4 curve bears no relationship with V_cyc. However, with the V_dil(Gln) term included, Equation (2) of the simplified model clearly shows that the entire glutamine C4 turnover curve (fe(Gln C4) versus time) is sensitive to V_cyc, even at isotopic steady state. The relationship fe(Gln C4) < fe(Glu C4) is maintained at isotopic steady state. Unless V_cycV_dil(Gln), the extent of the label dilution at glutamate C4 will not be the same as that at glutamine C4. As illustrated by Equations (1 and 2), V_dil(Gln) affects the turnover of glutamine C4 quite differently from V_cyc, allowing the two fluxes to be separately determined.

The sensitivity of the glutamine C4 turnover curve to V_cyc, as illustrated by the simplified model shown in Figure 1B and the values given in Table 1, was confirmed using a direct sensitivity analysis of the two-compartment model shown in Figure 1A. Figures 3A and 3B show the effect of changing the value of V_cyc by ±50% (from 0.32 to 0.16 and 0.48 μmol/g/min) while maintaining all other parameters unchanged. In the case of V_dil(Gln) = 0, the end point glutamine C4 intensities converge as the isotopic steady state is approached, and where the fractional enrichments of glutamine C4 and glutamate C4 are equal regardless of V_cyc. Much larger changes in the simulated glutamine C4 turnover curve was clearly seen in Figure 3B for V_dil(Gln) = 0.14 μmol/g/min, indicating that the glutamine C4 enrichment curve is very sensitive to changes in V_cyc, especially so when V_cyc is not much larger than 0.32 μmol/g/min. Note that, in Figure 3B, both the glutamate C4 and glutamine C4 intensities are lower than in Figure 3A because of the additional dilution and due to the glutamine-to-glutamate half of the glutamate–glutamine cycling which passes label dilution from glutamine C4 to glutamate C4. The glutamate C4 turnover curve is also responsive (in the opposite direction) to changes in V_cyc. In fact, the relative constancy of the glutamine C4-to-glutamate C4 ratio after the initial infusion period and its significant sensitivity to V_cyc indicate that this ratio could be used as a surrogate marker for V_cyc in clinical studies where relatively long scan times to capture the full labeling kinetics is difficult or impossible to implement.

Open in a separate window

Figure 3

The effect of changing the value of V_cyc by ±50% (from 0.32 to 0.16 and 0.48 μmol/g/min) on the glutamate C4 and glutamine C4 turnover curves while maintaining all other parameters unchanged. (A) V_dil(Gln) = 0. (B) V_dil(Gln) = 0.14 μmol/g/min. Red line: V_cyc = 0.32 μmol/g/min; black line: V_cyc = 0.48 μmol/g/min; blue line: V_cyc = 0.16 μmol/g/min.

The importance of V_dil(Gln) in the sensitivity of the ¹³C turnover curves to the value of V_cyc was further confirmed by the expected finding of a relatively strong positive correlation between V_cyc and V_dil(Gln) using the metabolic model in Figure 1A in cases where no constraints were applied (Table 1, the first two rows). The corresponding Pearson's product–moment correlation coefficient for V_cyc~V_dil(Gln) was 0.87 for σ = 0.1 μmol/g and 0.74 for σ = 0.2 μmol/g, respectively. Note that in Figure 3, the incorporation of V_dil(Gln) into the model also causes a decrease in labeled glutamate C4 concentration because of the glutamine-to-glutamate half of the glutamate–glutamine cycle. Because ¹³C label dilution at the level of pyruvate/lactate (V_dil(Lac)) competes with V_dil(Gln) in fitting to the glutamate and glutamine C4 enrichment levels, a negative correlation is expected between V_dil(Lac) and V_dil(Gln) and, by inference, between V_cyc and V_dil(Lac) for the first two cases given in Table 1 (the first two rows). The corresponding Pearson's correlation coefficient for V_cyc~V_dil(Lac) is −0.89 for σ = 0.1 μmol/g and −0.42 for σ = 0.2 μmol/g, respectively.

Effects of Constraining other Metabolic Fluxes

A limitation of obtaining metabolic fluxes from glutamate and glutamine C4 labeling curves for [1-¹³C]glucose as precursor is that the rate information on anaplerosis and the astroglial TCA cycle resides primarily in the early part of the glutamate and glutamine labeling curves, placing a premium on the signal-to-noise ratio of the acquired data. An error in the calculation of these fluxes will propagate into the V_cyc calculation, as shown in Table 1. With no constraints the V_cyc values obtained were 0.27±0.12 and 0.16±0.11 μmol/g/min for noise levels of 0.1 and 0.2 μmol/g, respectively. The corresponding ^aV_TCA values were 0.19±0.27 and 1.11±1.95 μmol/g/min, respectively. The downward trend in the calculated value of V_cyc was due to the high value of ^aV_TCA. At the extreme limit where ^aV_TCA is infinitely fast, it can be shown that the value found for V_cyc will approach asymptotically a low value of ~0.16 μmol/g/min because of the very rapid labeling of the astroglial glutamate precursors. Table 1 shows that the bias in the mean value of V_cyc was removed when ^aV_TCA was singularly constrained (with an upper limit of 0.1 μmol/g/min). Fitting the simulated glutamate and glutamine C4 time courses with σ = 0.1 μmol/g and the constraint ^aV_TCA≤0.15 or ≤0.20 μmol/g/min gave values of V_cyc of 0.31±0.10 or 0.29±0.10 μmol/g/min, respectively. We note that these constraints reflect values of ^aV_TCA above those found by direct measurement in vivo (Blüml et al, 2002; Gruetter et al, 1998, 2001; Lebon et al, 2002). Furthermore, the Lebon et al (2002) and Blüml et al (2002) studies used glial-specific substrates, which enhanced the accuracy of the measurement.

Another limitation of extracting metabolic fluxes from glutamate and glutamine C4 data and [1-¹³C] or [1,6-¹³C₂]glucose infusion is that there is no information on the rate of ¹³C label exchange between the mitochondria and cytosol (V_x). In the model used by Shen et al (1999), V_x was assumed to be much greater than ⁿV_TCA. The assumption V_xⁿV_TCA is supported by the studies of Patel et al (2004) and de Graaf et al (2004). Nevertheless, it should be pointed out that V_x has little effect on the astroglial isotopic dilution flux V_dil(Gln). Regardless of V_x, the extent of astroglial dilution of glutamine C4 is determined by the relative magnitude of label influx (mainly V_cyc × fe(Glu C4)) versus label efflux (V_dil(Gln) × fe(Gln C4) + V_cyc × fe(Gln C4)). Equation (1) (and therefore its solution) does not change upon addition of V_x into the model shown in Figure 1B.

An additional reduction in the standard deviation was obtained by constraining glutamine efflux from the brain based on experimental measurements of glutamine efflux, ammonia fixation, and isotopic labeling of glutamine through the anaplerotic pathway, which mirrors the approach taken in previous work (Mason et al, 2003, 2007; Shen et al, 1998; Sibson et al, 2001). To explore the effect of this parameter we constrained V_efflux as 0.2V_Gln based on the average of several measurements of the human brain (Gruetter et al, 2001; Mason et al, 2007; Shen et al, 1999), where V_Gln = V_cyc + V_efflux. We note that this is a high value because the efflux of glutamine is considerably lower than the rate of anaplerosis from which the value of V_efflux is derived in human MRS studies (Lebon et al, 2002; Mason et al, 2007). This constraint reduced the standard deviation in V_cyc by preventing the assignment of most of V_Gln to the anaplerotic flux (V_ana) due to the relatively low sensitivity of the glutamate and glutamine C4 curves to ¹³C labeling through this pathway. A secondary effect of this constraint was revealed as a downward bias in the value of V_cyc due to the assumed value of V_efflux being higher than what was used to generate the simulated curves. However, a direct measurement may be obtained by use of a glial-specific substrate to measure the rate of anaplerosis, providing a maximum estimate of V_efflux (Sibson et al, 2001; Mason et al, 2007) or direct NH₃ incorporation into glutamine (Shen et al, 1998).

It is significant to note that in previously published studies where fitted glutamate and glutamine C4 time courses were used in the determination of ⁿV_TCA and/or V_cyc, the fluxes V_x, V_efflux/V_ana, and ^aV_TCA were constrained either explicitly or implicitly (Gruetter et al, 1998, 2001; Shen et al, 1999; Sibson et al, 1997, 1998). Our simulations indicate that the improvements obtained with these constraints are significant. As discussed below the uncertainties we derive are highly consistent with the reported standard deviations in all of the human studies published to date (Table 3).

Table 3

Experimental measurement of V_cyc in the resting human brain from literature^a

Reference	Vcyc	nVTCA	Vcyc/nVTCA
Gruetter et al (1998)	0.32 (0.09)	0.86 (0.06)	0.37b
Shen et al (1999)	0.32 (0.05)	0.71 (0.08)	0.45b
Gruetter et al (2001)	0.17 (0.05)	0.57 (0.06)	0.30b
Chhina et al (2001)c	0.25 (0.08)d	0.75 (0.06)d	0.33b
Lebon et al (2002)	0.32 (0.07)e	0.80f	0.39 (0.09)
Mason et al (2003)g	0.21 (0.04)	0.52 (0.12)	0.40b
Mason et al (2007)	0.34 (0.16)	0.72 (0.25)	0.47b

Open in a separate window

V_cyc, cerebral rate of the glutamate–glutamine cycle; V_TCA, trichloroacetic acid cycle flux.

^aAll flux rates are in the units of μmol/g/min. Data in () are standard deviation in the same unit.

^bCalculated by the ratio of the mean V_cyc and ⁿV_TCA.

^cResting brain data only.

^dCalculated from the means of the reported V_cyc values in Chhina et al (2001).

^eCalculated by assuming ⁿV_TCA = 0.80 μmol/g/min.

^fAssumed ⁿV_TCA.

^gFrom V_Gln (i.v.) and V_cyc = 0.9 V_Gln as in Mason et al, 2003.

Comparison with Experimentally Derived Metabolic Rates and Uncertainty

Table 3 shows the published mean and standard deviation of ⁿV_TCA and V_cyc for all of the human ¹³C MRS studies to date that analyzed their data using a two-compartment model. In addition the ratio V_cyc/ⁿV_TCA is given for each study. The studies of Shen et al (1999) and later ones (Mason et al, 2003, 2007) included the astroglial dilution flux. We note that the standard deviations in our previously published reports are consistent with the results of the simulations presented here. Lebon et al (2002) used the steady state enrichments from the metabolism of [2-¹³C]acetate, a glial-specific tracer, to determine the ratio of V_cyc/ⁿV_TCA. The use of the steady state enrichment data to compute the ratio (relative fluxes) removes the dependence on the early portions of the labeling curves. As shown in Table 3, the value of the ratio, V_cyc/ⁿV_TCA, measured with [1-¹³C]glucose as tracer is consistent with that measured using the glial-specific [2-¹³C]acetate tracer.

This work was supported by the Intramural Program of NIH, NIMH, NS37527 (NINDS), and DK27121 (NIDDK).