Modeling fructose-load-induced hepatic de-novo lipogenesis by model simplification

Abstract

Hepatic de-novo lipogenesis is a metabolic process implemented in the pathogenesis of type 2 diabetes. Clinically, the rate of this process can be ascertained by use of labeled acetate and stimulation by fructose administration. A systems pharmacology model of this process is desirable because it facilitates the description, analysis, and prediction of this experiment. Due to the multiple enzymes involved in de-novo lipogenesis, and the limited data, it is desirable to use single functional expressions to encapsulate the flux between multiple enzymes. To accomplish this we developed a novel simplification technique which uses the available information about the properties of the individual enzymes to bound the parameters of a single governing ‘transfer function’. This method should be applicable to any model with linear chains of enzymes that are well stimulated. We validated this approach with computational simulations and analytical justification in a limiting case. Using this technique we generated a simple model of hepatic de-novo lipogenesis in these experimental conditions that matched prior data. This model can be used to assess pharmacological intervention at specific points on this pathway. We have demonstrated this with prospective simulation of acetyl-CoA carboxylase inhibition. This simplification technique suggests how the constituent properties of an enzymatic chain of reactions gives rise to the sensitivity (to substrate) of the pathway as a whole.

Keywords

Linear chain de-novo lipogenesis model simplification

Introduction

De-novo lipogenesis (DNL) is a complicated hormonally regulated pathway that is active in humans in liver and adipose tissue,¹ and to a lesser extent in skeletal muscle.² This pathway is one mechanism for converting carbohydrates into fat, which is eventually exported, or stored, as triglycerides (TGs). Hepatic DNL has been implicated in the pathogenesis of diabetes and hepatic steatosis.^3,4 In type 2 diabetics, where homeostatic mechanisms fail to control the level of glucose and gluconeogeneic precursors,⁵ hepatic DNL may act as a way of diverting these substrates away from gluconeogenesis and hepatic glucose output.⁶ However, the increased lipid load is an independent risk factor for cardiovascular disease.⁷

A mechanistic, quantitative, understanding of the relationship between substrate load and lipid accumulation via DNL is fundamental to fully appreciating the role of DNL in disease. Fructose is a lipogenic dietary substrate that is up-taken into hepatocytes rapidly on the first pass through the liver.⁸ Hence, one technique to stimulate and assay DNL is repeated oral dosing of fructose in the presence of exogenous C14 labeled acetate. Then DNL is measured as the fraction of plasma TGs that have incorporated C14 acetate, using the MIDA method. Hudgins et al.⁹ reported that this technique increased de-novo lipids in the plasma by approximately 15%. The purpose of this paper is to develop and use a mathematical ‘systems pharmacology’ model of this protocol to gain a quantitative understanding of the process. This will allow us, in the future, to predict the effect of putative pharmacological agents in this context.

The development of systems pharmacology models is a time-consuming process. This is because they attempt to faithfully represent key physiological mechanisms mathematically. Hence, they are often large and complex, requiring data from multiple species and experiments. A crucial challenge within the development process is accommodating a balance between mechanistic detail and available experimental knowledge to inform that detail. Therefore, we believe that tools which resolve this tension by systematic simplification of models are particularly valuable. Often the reduction to simplified representations occurs in the development and design phase of the process. In this case it is usually difficult to retain the correspondence to physiology and experimental data.

In the protocol applied by Hudgins et al.,⁹ it is plausible to assume that metabolism of fructose to plasma TGs can be well approximated by a linear network (Figure 1). The most direct route from fructose to plasma TGs involves 16 intermediaries: see the supplementary file. However, even in a model where each reaction is described by Hill kinetics this will give rise to at least a 32 parameter model, whereas if we seek to simply retain the boxed species (Figure 1), and construct transfer functions between them, the size of the model will be reduced dramatically. As such, it is desirable to design a simplified representative model of this process incorporating key mechanistic features and known experimental data regarding enzyme properties of the individual steps, while not modeling each enzyme explicitly.

Figure 1.

The shortest route from fructose to plasma triglycerides (TGs) includes 16 steps (see Supplementary Table 2). In the experimental conditions of Hudgins et al.,⁹ where hepatic fructose exposure is high for a period of approximately 8 hours, we assume that the overwhelming majority of fructose-derived TGs has been synthesized by this pathway.

Model reduction will often be an esoteric case-dependent exercise. However, there are some general approaches that have been developed. For example, one way to simplify a model is to separate time-scales, so that interactions faster than the time-scale of interest can be assumed to be in a quasi steady-state.^10–12 Another approach is to reduce parameters based on Hankel singular values, and/or singular value decomposition.¹³ This method ranks the parameters based on their effect on the input–output response of the model. The approach we are proposing here is a novel contribution because it uses available experimental information to establish approximate values for lumped parameters that otherwise do not have a direct biochemical and experimental meaning.

An obvious way to simplify a model is to compact linear-chain motifs (synonymously, path graph or linear tree) of a network into the input and output of the motif. If we believe the intermediaries of this graph are not rate-limiting and we do not have time-course data for those nodes, then all we are really interested in, or can constrain, is the ‘transfer function’ from the first to the last node (Figure 2).

Figure 2.

Approximating linear chains with transfer functions. Often it is desirable to reduce the ‘full model’ motif, to a ‘simplified model’ motif. In this case, the input species $X_{1}$ and the output species $X_{n}$ are retained in the simplification, and the intermediary fluxes approximated by the total flux, captured by $T (X_{1})$ , which must approximate $J_{n - 1} (X_{n - 1})$ .

For enzymatic reaction networks a biologically relevant transfer function is a Hill function,¹⁴ which arises naturally in the case of Michaelis–Menten kinetics. This rate law is typically a reasonable assumption for metabolic networks, but for signal transduction networks, where interactions are more complex (for example, when proteins can be both enzymes and substrates) there are superior approximations.¹⁵ Here we assume that Hill functions (with Hill coefficient one) can capture the experimental dynamics and still be predictive, even if the assumptions required for Michaelis–Menten kinetics are not precisely met. Hill functions, similar to metabolic reactions, are saturable, and can give rise to rich dynamics,¹⁶ which we presume includes the requisite model behavior. Here we consider Hill functions with coefficient one, but our analysis shows that similar relationships to the one we develop will hold true at higher values.

Our approach to simplifying linear-chain motifs uses prior knowledge of the system under consideration. For metabolic pathways it is often the case that the maximum rate of a reaction ( $V_{\max}$ ) is not known, but enzyme–substrate affinity is known (related to the Michaelis constant $K_{m}$ in the case of Michaelis–Menten kinetics). This is because $K_{m}$ can be measured in vitro, whereas $V_{\max}$ is a function of in-vivo enzyme expression and regulation, and is less readily known in vivo. Furthermore, in dysregulated disease states $V_{\max}$ for given reactions could be, for example, greatly attenuated. In systems pharmacology models these parameters are often using experiments that capture the whole system behavior (i.e. the system the model is describing), rather than informed from direct measurements of the behavior of constitutive component of the model.

A linear chain of n reactions modeled with Michaelis–Menten kinetics will have 2n parameters. By introducing a transfer function from the first to last node we can simplify this motif to being determined by either two or four parameters, depending on if we assume conservation of mass (so that the flux leaving the first node is matched by flux into the last). In this case, mass is not conserved in the model because it is ‘lost’ to the intermediary substrates that we are not considering. However, due to the linear structure of the pathway, we gain a constraint that the maximum rate, $V_{T}$ , of the transfer function should be less than or equal to $V_{1}$ (the maximum flux out of the first node).

While this simplification is often necessary, it is still highly desirable to use all prior knowledge to constrain the transfer function. As mentioned, it is unusual for detailed information regarding the $V_{i}$ to be available (however, the maximum rate, $V_{T}$ , of the transfer function is more readily measured or constrained by experimental data). Although there are graph theoretic approaches to solving models with simpler flux terms,¹⁷ as far as we are aware there are no existing methods for using data about individual reactions for constraining $K_{T}$ , even though detailed measurements of these values are often available. How to use available information to estimate $K_{T}$ is the focus of this work.

To motivate our approach for simplification we first introduce a novel model of hepatic DNL which includes the aforementioned ‘transfer functions’. Following this we describe our novel approach to parameterizing these functions, which is justified by both computational and analytical rationale.

Modeling hepatic DNL in response to high fructose levels

Model equations

The model detailed in Figure 3 was described mathematically as follows. It was assumed that all reactions could be captured by Hill functions (i.e. can be approximated by Michaelis–Menten-type kinetics). This approximation will not always be valid, and certainly the Michaelis–Menten interpretation of their meaning is not necessarily valid unless specific assumptions are met. We assume that in conditions where the pathway is operating unidirectionally (i.e. highly stimulated at the initial node) that this functional form sufficiently captures the flux from node to node.

Figure 3.

Simplified model of DNL in fructose-stimulated conditions. This model structure was designed based on known biology, experimental protocol,⁹ and the simplification of a larger metabolic pathway (Figure 2).

Notation

We use F, P, $[F A]$ , and $[T G]$ to denote fructose, pyruvate, fatty acids, and TGs, respectively. Square brackets denote species concentration. Subscripts refer to different pools for each species (bolus, plasma, or liver). A superscript D denotes the metabolites which have derived from DNL (identified experimentally with C14 acetate) a superscript T denotes total metabolite concentration (unlabeled + labeled):

\frac{d F_{B}}{d t} = D . δ (T, t) - β F_{B}

(1)

\frac{d [F_{P}]}{d t} = \frac{(1 - ρ_{1}) β F_{B}}{V_{P}^{F}} - γ [F_{P}] + α_{1}

(2)

\frac{d [F_{L}]}{d t} = \frac{ρ_{1} β F_{B}}{V_{L}^{F}} - J_{1}

(3)

\frac{d [P_{L}]}{d t} = T_{1} - J_{2} + α_{2}

(4)

\frac{d {[F A]}_{L}}{d t} = (1 - R) T_{2} - r_{1} J_{3} - σ_{1} {[F A]}_{L} + α_{3}

(5)

\frac{d {[T G]}_{L}}{d t} = r_{1} T_{3} - r_{2} J_{4}

(6)

\frac{d {[T G]}_{P}}{d t} = V_{r} r_{2} J_{4} - σ_{2} {[T G]}_{P} + α_{4}

(7)

\frac{d {[F A]}_{L}^{D}}{d t} = R T_{2} - (1 - r_{1}) J_{3} - σ_{1} {[F A]}_{L}^{D}

(8)

\frac{d {[T G]}_{L}^{D}}{d t} = (1 - r_{1}) T_{1} - (1 - r_{2}) J_{4}

(9)

\frac{d {[T G]}_{P}^{D}}{d t} = V_{r} (1 - r_{2}) J_{4} - σ_{2} {[T G]}_{P}^{D}

(10)

Rationale for equations (1)–(3)

These equations predict plasma and liver fructose following oral administration of fructose at discrete times given by the set $T$ . Here D is the dose of fructose administered at each time $T_{i} \in T$ , and $δ$ is the dirac delta function. We use $F_{B}$ to denote the amount of fructose bolus that has not yet appeared in plasma in units of mmol and $[F_{P}]$ is the plasma fructose concentration in units of mmol/L. Here $V_{p}^{F}$ is the volume of distribution of fructose in plasma. Plasma fructose is cleared at a rate $γ$ in units of $\min^{- 1}$ . In this model, the fructose bolus is partitioned into either liver or plasma $ρ_{1}$ is the fraction absorbed into the liver, with volume of distribution $V_{L}^{F}$ . Here $V_{r}$ is the ratio of liver and plasma volume of distributions for TGs.

Rationale for equations (4)–(10)

The flux terms in these expression are defined as outlined in Figure 2, more specifically:

J_{1} = \frac{V_{1} [F_{L}]}{K_{1} + [F_{L}]}

(11)

T_{1} = \frac{V_{2} [F_{L}]}{K_{2} + [F_{L}]}

(12)

J_{2} = \frac{V_{3} [P_{L}]}{K_{3} + [P_{L}]}

(13)

T_{2} = \frac{V_{4} [P_{L}]}{K_{4} + [P_{L}]}

(14)

J_{3} = \frac{V_{5} {[F A]}_{L}^{T}}{K_{5} + {[F A]}_{L}^{T}}

(15)

T_{3} = \frac{V_{6} {[F A]}_{L}^{T}}{K_{6} + {[F A]}_{L}^{T}}

(16)

J_{4} = \frac{V_{7} {[T G]}_{L}^{T}}{K_{7} + {[T G]}_{L}^{T}}

(17)

When applied, $J_{3}$ , $T_{3}$ , and $J_{4}$ are scaled by $r_{1}$ and $r_{2}$ which are defined as the fraction unlabeled metabolite in each case (see equations (6), (7), (9), and (10)). Due to the constraints $T_{i} ⩽ J_{i}$ (see equation (22)), we can rewrite $V_{2} = ρ_{2} V_{1}$ , $V_{4} = ρ_{3} V_{3}$ and $V_{6} = ρ_{4} J_{5}$ , where $ρ_{i} \in [0, 1]$ .

Hudgins et al.⁹ measured DNL by using C14 acetate, which identifies the de-novo lipids using the MIDA method. In the model we can consider that $R (t)$ is the percentage of fatty acids from DNL which are ‘labeled’ at time $t$ , and that $R (t)$ is dependent solely on how much of the acetate pool is labeled. Note that in actuality the MIDA method is a calculation to correct for the combinatorics of combining endogenous and exogenously derived acetyl-CoA to form lipids, however in the model we can identify the metabolites originating from de-novo lipids directly. To allow for the tracer kinetics, we calculate $R (t)$ , which is at steady state until the acetate infusion is stopped, as

\begin{array}{l} R = \frac{A_{L} (t)}{A_{E} (t) + A_{L} (t)} \\ R = \frac{1}{1 + \frac{k_{E}}{k_{I}} e^{λ (t - T)}} \end{array}

where $A_{L} (t) = (k_{I} / λ) e^{- λ (t - T)}$ and $A_{E} (t) = k_{E} / λ$ are the labeled and endogenous acetate, respectively, and T is the time that infusion is stopped. Prior to time T we assume that R is at steady state $= \frac{K_{I}}{K_{I} + K_{E}}$ . The endogenous rate of appearance of acetate ( $k_{E}$ ) and its clearance ( $λ$ ) are taken from literature values.

Parameter estimates

Parameter ranges, units, and final model fits are reported in the supplementary material. Here we detail our strategy for generating estimates for the parameter ranges, and fitting the model.

For equations (1)–(3) parameters $β$ , $γ$ , and $α_{1}$ were independently fit using data from an oral fructose test.¹⁸ Fractional absorption of fructose, $ρ$ , was estimated as 0.72, based on a reported 72% of oral fructose being up taken into the liver of fasted rats, and in humans at least half of intravenously delivered fructose being metabolized hepatically.¹⁹ Here $α_{1}$ was estimated assuming a low basal concentration of plasma fructose of 51 µmol/L.¹⁸

Feasible ranges for the zeroth-order source terms in the model ( $α_{i}$ ) were made as follows. Here $α_{1}$ was assumed to be known with high confidence (see above) and $α_{2}$ , a source term for non-fructose derived pyruvate, was estimated by noting that rate of hepatic glucose output is approximately 0.8 mmol/min,²⁰ a significant portion of this will be from glycogen. Assuming that $α_{2}$ would be within a reasonable range of this value, we constrained $α_{2}$ to be within ${[10}^{- 6}, 5 \times 0.8 \times V_{L}^{F}]$ mmol/L/min. The uptake rate of fatty acids, $α_{3}$ , was estimated to be a maximum of 0.6 mmol/min,²¹ similarly a wide including this value was confined: ${[10}^{- 6}, 2 \times 0.6 / V_{L}^{F}]$ mmol/L/min. The basal source term for plasma TGs, $α_{4}$ was estimated from the estimate for clearance $σ_{2}$ (estimated from the fructose DNL experiment⁹), and reported values for plasma TG concentrations (1.2 mmol/L).²² Total whole body fatty acid oxidation was estimated to be between $[0, 995]$ mg/kg/day as reported by Aarsland et al.,²³ depending on substrate availability. Given the unique nature of the protocol being simulated here, and the large variability in the reported values, this range was used to calculate the first-order rate constant ( $σ_{2}$ ), by dividing by reported typical hepatic fatty acid levels.²⁴

Apart from the Michaelis constants in the transfer functions ( $K_{2}$ , $K_{4}$ , $K_{6}$ ), for which we describe a novel approach in Section 3, the remaining parameters of the DNL model were estimated using a bounded simulated annealing routine (simulannealingbnd²⁵) which is an algorithm for bounded global optimization. For each parameter a range was estimated to pass to the routine.

If a range was not estimable from experimental data, a reasonably wide-range around the reported value was passed to the optimization routine: this range spanned one, two, or three orders of magnitude, the latter being assigned to the $V_{i}$ parameters. This is representative of the large degree of uncertainty in their value.

Computational approach to estimate parameters of the transfer functions

In this section we describe a novel computational approach to estimate a value of $K_{2}$ , $K_{4}$ , and $K_{6}$ in the transfer functions (equations (12), (14), and (16)). These parameters have no direct biological interpretation because the transfer functions (12), (16), and (16) are lumped representations of multiple metabolic steps. Nevertheless, if there are data regarding the individual reactions then it is desirable to apply that data to inform the possible values of these parameters.

Below we show that with limited assumptions about the nature of the pathway and experimental conditions, a rule can be constructed to estimate the expected value, $E (K_{T})$ , of the transfer-function parameter $K_{T}$ , based solely on the values of the $K_{i}$ (the Michaelis–Menten constant of each enzymatic/metabolic step, see Figure 1). This rule is analogous to the well-known result that $V_{T} = m i n (V_{i})$ , i.e. the maximum rate of the pathway is controlled by the rate-limiting step.

It is imperative to point out that the rule developed in this section only provides an estimate, or bounds (depending on prior knowledge), for $E (K_{T})$ that is required to be refined based on fitting of experimental data (specifically data where input and output to the system are measured). Nevertheless, this rule provides a rationale for incorporating prior knowledge about the motif into the lumped parameter estimation. For the model introduced in Section 2 this method was used to estimate a plausible range for the fitted $K_{T}$ . This range was then used in the fitting procedure as a bound on $K_{T}$ . To the best of the authors’ knowledge there was no other way to construct a relatively narrow, but physiologically relevant, range for $K_{T}$ without this analysis.

Our strategy for approaching this problem was to simulate the full model in Figure 2, and attempt to fit the flux into the final node as function of the concentration of the first node (Figure 2, inset). In this manner we could attempt to establish a relationship between the set ${K_{i}}$ and $K_{T}$ . Computational simulation (see Appendix A.2) illustrated a reproducible relationship between $Σ K_{i}$ and $K_{T}$ . Specifically this relationship is that $\log (E (K_{T}))$ is a linear function of $\log (Σ K_{i})$ , where $E (K_{T})$ is the expected value of the transfer function (assuming that the $V_{i}$ are unknown). Figure 4 is an example of this relationship for a chain of length 8. Note that this relationship provides an order-of-magnitude estimate for $K_{T}$ which, although broad, is an improvement on the a priori range for human enzyme $K_{m}$ (which can be estimated as ≈10 order of magnitudes, as assessed by the BRENDA tool²⁶). This estimate provides a rationale for establishing a suitable, and computationally feasible, range to pass to a constrained optimization routine (or to evaluate model performance using unconstrained optimization). While this estimate may not lead to improved model fits, it will lead to a model where the parameter choices have improved experimental justification.

Figure 4.

Predicting transfer function kinetic parameters from individual flux parameters. For a high stimulus level and given chain length ( $n = 8$ ) $\log (Σ K_{i})$ is linearly correlated ( $r^{2} = 0.52$ ) with $\log (K_{T})$ .

To investigate how robust this relationship is, we repeated this procedure for different lengths of chains (n = 3 to 24) and stimulus strengths (over five orders of magnitude). We noticed that the correlation between $\log (E (K_{T}))$ and $\log (Σ K_{i})$ was obscured at lower stimulus levels (Figure 5A), but was largely independent of the length of the chain (Figure 5B and C). The y-intercept of the regression line is determined by the stimulus level (Figure 5C).

Figure 5.

The effect of chain length and stimulus level on the relationship between $K_{T}$ and ∑ $K_{i}$ ,see Figure 4. A relationship and correlation between transfer function kinetic parameters from individual flux parameters is maintained over a wide range of chain-lengths and stimulus level. For weak stimuli the correlation is lost.

In the appendix we provide an analytical analysis of this problem, which provides a direct interpretation of why this relationship holds in a limiting case.

Simplification parameters for the hepatic DNL model

An unknown in applying this technique to the DNL model is the stimulus strength relative to the $E (V_{i})$ . Given the high levels of fructose loading in the clinical experiment,⁹ and based on our analytical analysis, we assumed a range of [0,1] for the y-intercept (equation (21)). This corresponds to a range of stimulus of $S \approx E (V_{i})$ (the stimulus is likely to be saturating some reactions on the pathway) to $S \approx 10 E (V_{i})$ (the stimulus is very likely to be saturating reactions on the pathway).

Application to and simulation of a model of hepatic DNL

To demonstrate this technique we applied it to the model of hepatic DNL presented in Figure 5. Note that this model represents a significant simplification of the 16 enzyme path (Figure 1). The individual reactions and estimates for their $K_{m}$ are given in Supplementary Table 2. In this case, we have assumed that we can daisy chain this approximation together by applying this simplification three times. This allowed us to take account of the major sources and sinks in this pathway, and retain a suitable level of detail for future mechanistic questions. Applying the method above to this data gave a range of values for the $K_{T}$ for each of the transfer functions in this model (Supplementary Table 3). The range for each $K_{i}$ is generated by both experimental uncertainty and uncertainty estimating the stimulus level relative to ${V_{i}}$ in the pathway. We could have considered the additional uncertainty due to $Var (K_{T})$ , given that equation (21), only provides an estimate for the expected value of $K_{T}$ , but this was not necessary in this example. In this case we took the maximum value of the y-intercept in the relationship as being in the range [0,1], which in our prior simulations (Figure 5) corresponded to a range where the stimulus level ≈[0,1] times the average $V_{i}$ in the pathway. This range represents mildly saturating to highly saturating conditions, which seemed reasonable based on observed experimental data. Since we only established a search range for the $K_{i}$ parameters, and since there were still unconstrained parameters in the model (for example, the $V m a x$ parameters), it was necessary to fit the model against further experimental, dynamic, data. To do this we applied this model to previously published human data, where subjects were repeatedly fed fructose over an 8 hour time period and the appearance of de-novo lipids in blood was measured.⁹ The model could be fit to the available data but not uniquely (Figure 6). This is typical for a model of this size applied to limited human data. However, this did test that estimating the $K_{i}$ parameters of the transfer functions did not lead to constraints that were discordant with available dynamic data. Given the remaining uncertainty we captured multiple parameter sets (≈3000), so that we can propagate the uncertainty forward for future simulations (Supplementary Table 1).

Figure 6.

Summary of model performance compared with clinical data from Hudgins et al.9 (black points). Multiple fits (≈1800) were found that could describe the data. The upper-quartile of fits were found to cluster into four distinct groups (cluster centroid, solid red lines). The shaded region summarizes the upper and lower bounds of all fits with outliers removed (defined as the top and bottom percentiles).

Because novel inhibitors of acetyl-CoA carboxylase (ACC) have been developed,^27,28 we wanted to investigate the effect of inhibiting ACC in this model. This was simulated by re-running the model with a reduction of parameter $V_{4}$ by 90% and 50%. Recalling the simplification paradigm, this parameter is the $V_{\max}$ of the rate limiting step of the reactions between pyruvate and fatty acids (Figure 3). In this case, the rate-limiting step of this part of the pathway, fat synthesis, is ACC.²⁹ Hence, this modification represents inhibition of ACC.

Figure 7 shows the effect of ACC inhibition in the model. Inhibition of ACC is predicted to have a potent effect on fractional DNL as assayed by this clinical methodology. The model also predicts that in a clinical trial the variability in response will be higher for doses that induce 50% inhibition than placebo or higher doses (as shown by wider prediction intervals for the 50% case). This prediction illustrates the utility in generating multiple parameter sets to explore variability in response, and also has implications for how such a study should be conducted. However, this effect could be masked by other sources of variability.

Figure 7.

Evaluation of the model response to 50% and 90% inhibition of acetyl-CoA carboxylase (ACC) in comparison with the baseline response (shaded area is the 10–90% prediction interval). The model predicts that the fractional DNL response to fructose stimulation can be effectively modulated by inhibition of ACC. Variability in predicted response at 50% inhibition is significantly larger than the baseline and 90% case.

Since the expression of ACC can be tightly controlled,³⁰ it is interesting to interpret these results as the effect of modulating ACC expression to 50% or 10% of baseline (mathematically, these are the same simulations as Figure 7). In this interpretation the results supports the hypothesis that regulation of ACC by SREBP1c³⁰ is a crucial component of the lipogenic regulatory program controlled by SREBPC1 and regulation of ACC expression directly controls the contribution of DNL to the liver and plasma TG pools.

Discussion

Here we have developed a novel systems pharmacology model of hepatic DNL. In particular, this model is designed to support, and analyze, clinical experiments where hepatic DNL is measured. As such, this model could be a powerful tool to assess novel pharmacological agents that may affect DNL. For example it could be used to infer pharmacologically meaningful parameters such as IC50, or $I_{\max}$ from clinical data, as indicated in Figure 7.

As a focused simulation of a specific experimental protocol we have high confidence in the model when applied in this context, but in other contexts the model would require further validation and refinement before application. However, we believe that within this context the model is a good representation of the underlying biology because of the novel simplification approach we developed. This approach links underlying biochemical constants with a parameter of a lumped representation that hitherto had no biologically interpretable meaning.

The model we have presented here is an acute view of fructose-induced DNL and tracks a major metabolic route of fructose in these conditions. However, fructose metabolism also has a major signaling role which can alter downstream metabolic pathways. These signaling effects include activation of SREBP1c which has a profound effect on the entire lipogenic pathway.³⁰ Our assumption is that the acute nature of the study allows us to ignore these effects during the course of the simulated trial. However, future work could include analysis and modification of the model to describe patients from a different population, such as those with non-alcoholic fatty liver disease. This population consumes large amounts of dietary fructose and have a higher fasting DNL level. Such an analysis like could facilitate insight to the precise biochemical parameters that are inducing increased DNL in these patients, and how they might respond in this clinical protocol.

Systems pharmacology models are a useful tool to understand biological systems, disease pathophysiology, and the potential effects of pharmacological intervention. As more biological detail is added to a model it becomes less constrained and the balance between the number of parameters and available data becomes less favorable. Therefore it becomes necessary to simplify the model. However, there is a balance between redressing this while retaining enough mechanistic detail to satisfactorily address the biological questions and hypotheses of interest. In this case it is important that all available data is used to establish parameter values.

The model reduction technique applied here simplifies a model while utilizing available data about individual enzyme kinetics. This technique supplements existing approaches of constraining, or reducing, a model based on dynamic data. By bounding the search for lumped parameter values into a biologically and mathematically justified space, physiologically reasonable parameter sets can be identified. Our analysis and simulations show that the Michaelis constant of the transfer function is related to the sum of the individual Michaelis constants and the stimulus level relative to the distribution of $V_{i}$ in the pathway. We have demonstrated the validity of this approach by analysis of a limiting discontinuous case, and simulation of many parameterizations of an example motif in a continuous case. This motif was based upon a chain of enzymes that can be described by Michaelis–Menten kinetics. Typically, and in the DNL example we give, the regulation of metabolic enzymes and their kinetics are far more complicated than the relatively simple Michaelis–Menten representation (for example, product inhibition, multiple-substrate reactions, co-operative binding). Given our analytical results we hypothesize that these results should be valid as long as the pathway is saturable, and unidirectional (or is so in the conditions of interest). Application of this approach to enzymatic chains with mixed kinetics is a future avenue of research.

We have focused on this method in order to develop rules for fitting parameters. However, there is a biological interpretation of our result which generates an interesting hypothesis for future testing. We can consider the transfer function as a summary of the kinetics of the pathway. As such, $K_{T}$ is a measure of the sensitivity of the pathway to stimulus. If $K_{T}$ is small, the pathway is very sensitive and fully activated for small stimuli. Conversely, if $K_{T}$ is large the pathway is not sensitive. Equation (21) is then interpreted as follows: the sensitivity of a pathway decreases (larger $K_{T}$ ) in both long reaction chains or in high stimulus conditions, which seems intuitively obvious. The inverse of this is that, according to this analysis, increased enzyme expression sensitizes a pathway. This analysis also indicates that therapeutic intervention designed to improve sensitivity to substrate should be focused on the ‘sensitivity limiting step’ (i.e. the reaction with the highest $K_{m}$ ). In the future we hope to test these concepts experimentally.

Footnotes

Appendix

Acknowledgements

We would like to acknowledge Gianluca Nucci (Pfizer Inc.) for his support and motivation of this work. We would also like to thank Theodore Rieger (Pfizer Inc.), and William C. Thompson (formerly Pfizer Inc., currently SAS Institute Inc.) for valuable input and feedback on this work.

Author Contributions

RJA performed the analysis, wrote the manuscript, and approved the final manuscript. CJM reviewed and approved the final manuscript.

Peer review:

Five peer reviewers contributed to the peer review report. Reviewers’ reports totaled 673 words, excluding any confidential comments to the academic editor.

Funding:

The author(s) disclosed receipt of the following financial support for the research, authorship, and/or publication of this article: RJA, and CJM were employees of Pfizer Inc. during the completion and/or analyses of these studies. Pfizer Inc. has provided all funding for these studies and was responsible for the study design, data analysis, decision to publish, and preparation of the manuscript.

References

Schutz

Concept of fat balance in human obesity revisited with particular reference to de novo lipogenesis

Int J Obesity 2004;28:S3–S11.

Katsuhiko

Haowei

. Muscle lipogenesis balances insulin sensitivity and strength through calcium signaling J Clin Investig 2013;123:1229–1240.

Fatima

Lisa

Shahida

Hubert

Nousheen

De novo lipogenesis in health and disease

Metabolism 2014;63:895–902.

Ferré

Foufelle

Hepatic steatosis: a role for de novo lipogenesis and the transcription factor SREBP-1c. Diabetes Obesity Metabolism 2010;12(Suppl. 2):83–92.

Leahy

JL.

Pathogenesis of type 2 diabetes mellitus. Arch Med Res 2005;36:197–209.

Hellerstein

MK.

De novo lipogenesis in humans: metabolic and regulatory aspects

Eur J Clin Nutrition 1999;1:S53–65.

Jackson

Poppitt

Minihane

AM.

Postprandial lipemia and cardiovascular disease risk: Interrelationships between dietary, physiological and genetic determinants

Atherosclerosis 2012;220:22–33.

Stanhope

Griffen

Bair

Swarbrick

Keim

Havel

PJ.

Twenty-four-hour endocrine and metabolic profiles following consumption of high-fructose corn syrup-, sucrose-, fructose-, and glucose-sweetened beverages with meals. Amer J Clin Nutrition 2008;87:1194–1203.

Hudgins

Parker

Levine

Hellerstein

MK.

A dual sugar challenge test for lipogenic sensitivity to dietary fructose

J Clin Endocrinol Metabolism 2011;96:861–868.

10.

Krüger

Heinrich

. Model reduction and analysis of robustness for the Wnt/. BETA.-catenin signal transduction pathway Genome Informatics 2004;1:138–148.

11.

Segel

Slemrod

The quasi-steady-state assumption: a case study in perturbation

SIAM Rev 2006;31:446–477.

12.

Alon

Simplicity in biology

Nature 2007;446:497–497.

13.

Sun

Hahn

Parameter reduction for stable dynamical systems based on Hankel singular values and sensitivity analysis

Chem Eng Sci 2006;61:5393–5403.

14.

Aldridge

Burke

Lauffenburger

Sorger

PK.

Physicochemical modelling of cell signalling pathways. Nature Cell Biol 2006;8:1195–1203.

15.

Ciliberto

Capuani

Tyson

JJ.

Modeling networks of coupled enzymatic reactions using the total quasi-steady state approximation. PLoS Computat Biol 2007;3:e45.

16.

Allen

Elston

TC.

From physics to pharmacology?

Rep Prog Phys 2010;74:016601.

17.

Karmakar

Mandal

Graph theoretical analysis on the kinetic rate equations of linear chain and cyclic reaction networks. J Phys Chem A 2014;118:7672–7682.

18.

Chong

MF-F

Fielding

Frayn

. Mechanisms for the acute effect of fructose on postprandial lipemia. Amer J Clin Nutrition 2007;85:1511–1520.

19.

Mayes

PA.

Intermediary metabolism of fructose. Amer J Clin Nutrition 1993;58:754S–765S.

20.

Glauber

Wallace

Brechtel

Effects of fasting on plasma glucose and prolonged tracer measurement of hepatic glucose output in NIDDM

Diabetes 1987;10:1187–1194.

21.

Hallgreen

CE.

The interplay between glucose and fat metabolism: A biosimulation approach. PhD thesis, Technical University of Denmark, Copenhagen, 2009.

22.

Laws

Reaven

GM.

Evidence for an independent relationship between insulin resistance and fasting plasma HDL-cholesterol, triglyceride and insulin concentrations. Acta Medica Scandinavica 1992;231:25–30.

23.

Aarsland

Chinkes

Wolfe

RR.

Hepatic and whole-body fat synthesis in humans during carbohydrate overfeeding. Amer J Clin Nutrition 1997;65:1774–1782.

24.

Cairns

S R

Kark

A E

Peters

T J.

Raised hepatic free fatty acids in a patient with acute fatty liver after gastric surgery for morbid obesity. J Clin Pathol 1986;39:647–649.

25.

MATLAB. MATLAB and Global Optimization Toolbox Release 2014b, 2014.

26.

Schomburg

Chang

Schomburg

BRENDA, enzyme data and metabolic information. Nucleic Acids Res 2002;30:47–49.

27.

Griffith

Kung

Esler

. Decreasing the rate of metabolic ketone reduction in the discovery of a clinical acetyl-CoA carboxylase inhibitor for the treatment of diabetes. J Med Chem 2014;57:10512–10526.

28.

Harriman

Greenwood

Bhat

. Acetyl-CoA carboxylase inhibition by ND-630 reduces hepatic steatosis, improves insulin sensitivity, and modulates dyslipidemia in rats. Proc Natl Acad Sci 2016;113:E1796–E1805.

29.

Donaldson

. Regulation of fatty acid synthesis. Federation Proc 1979;38:2617–2621.

30.

Softic

Cohen

Kahn

CR.

Role of dietary fructose and hepatic de novo lipogenesis in fatty liver disease

Digest Dis Sci 2016;61:1282–1293.

Supplementary Material

Please find the following supplemental material available below.

For Open Access articles published under a Creative Commons License, all supplemental material carries the same license as the article it is associated with.

For non-Open Access articles published, all supplemental material carries a non-exclusive license, and permission requests for re-use of supplemental material or any part of supplemental material shall be sent directly to the copyright owner as specified in the copyright notice associated with the article.

0.17 MB

0.00 MB