A computational model of postprandial adipose tissue lipid metabolism derived using human arteriovenous stable isotope tracer data

Given the association of disturbances in non-esterified fatty acid (NEFA) metabolism with the development of Type 2 Diabetes and Non-Alcoholic Fatty Liver Disease, computational models of glucose-insulin dynamics have been extended to account for the interplay with NEFA. In this study, we use arteriovenous measurement across the subcutaneous adipose tissue during a mixed meal challenge test to evaluate the performance and underlying assumptions of three existing models of adipose tissue metabolism and construct a new, refined model of adipose tissue metabolism. Our model introduces new terms, explicitly accounting for the conversion of glucose to glyceraldehye-3-phosphate, the postprandial influx of glycerol into the adipose tissue, and several physiologically relevant delays in insulin signalling in order to better describe the measured adipose tissues fluxes. We then applied our refined model to human adipose tissue flux data collected before and after a diet intervention as part of the Yoyo study, to quantify the effects of caloric restriction on postprandial adipose tissue metabolism. Significant increases were observed in the model parameters describing the rate of uptake and release of both glycerol and NEFA. Additionally, decreases in the model’s delay in insulin signalling parameters indicates there is an improvement in adipose tissue insulin sensitivity following caloric restriction.


Introduction
The adipose tissue plays a key role in the regulation of plasma triglyceride and non-esterified fatty acid (NEFA) concentrations [1,2]. In the fasting state, lipolysis of triglycerides stored in the adipose tissue delivers NEFA to the plasma, where it can be taken up by other tissues, including skeletal muscle, heart, and liver [3]. In the postprandial state dietary derived chylomicron-triglycerides are removed from the circulation by insulin stimulated hydrolysis by lipoprotein lipase (LPL) at the endothelial cell wall [4]. NEFA released by LPL lipolysis is primarily taken up by the adipose tissue and may subsequently be stored as triglyceride [5]. Disturbances in the regulation of lipid metabolism by the adipose tissue can lead to raised plasma NEFA concentrations resulting in ectopic fat deposition. The accumulation of fat in insulin sensitive tissues such as skeletal muscle and liver has been related to the development of insulin resistance and, consequently, increases the risk for developing Type 2 Diabetes and cardiovascular disease [6][7][8][9]. The contribution of disturbed NEFA metabolism to the development of insulin resistance highlights the need to include adipose tissue NEFA metabolism when investigating disturbances in the glucose-insulin regulatory system during disease development [10].
Several computational models exist which describe adipose tissue NEFA metabolism; from smaller models focussing on NEFA dynamics between the plasma and adipose tissue alone [11,12] to larger models describing several metabolites across multiple tissues [10,13,14]. Typically, these models have been parameterised using values reported in literature and their performance assessed by comparison with measured plasma metabolite concentrations during an oral glucose or mixed meal challenge test. However, with traditional venous sampling it is not possible to differentiate between the contributions of individual tissues (i.e. liver or skeletal muscle) to the whole body plasma metabolite concentration. Consequently, it is not possible to fully assess the structure of the model and assumptions upon which the model has been constructed. Therefore, more detailed measurements are needed, providing information about contributions of individual tissues to postprandial changes in plasma NEFA concentrations.
In recent years, new quantitative knowledge has been gained on human adipose tissue NEFA dynamics in both the fasting and postprandial state through the use of stable isotope tracers coupled with arteriovenous sampling [5,15,16]. The simultaneous sampling of blood from an artery (or arterialised hand vein) and a vein draining a specific tissue depot (e.g. the abdominal subcutaneous adipose tissue) allows for the calculation of metabolite fluxes across this tissue. In combination with the use of stable isotope tracers, it is now possible to further untangle postprandial lipid metabolism and to quantify rates of appearance of triglyceride in the plasma, rates of LPL lipolysis of lipoprotein derived triglyceride, and spill-over of LPL derived NEFA into the plasma, while also providing insights into the behaviour of secondary metabolites such as glycerol.
Arteriovenous measurements provide invaluable information in the evaluation and refinement existing models of adipose tissue fatty acid metabolism. One example where such data has been successfully used is the 2008 study of Kim et al. [12], who modelled adipose tissue metabolism with a focus on the regulation of lipolysis. However, their model focused primarily on the fasting state, as well as the effect of an epinephrine infusion, and did not consider the postprandial state nor the influence of insulin.
In the present study, we used human in vivo measurements of net triglyceride, NEFA, glucose, and glycerol fluxes across the abdominal subcutaneous adipose tissue along with arterial insulin collected in the fasting and postprandial state during a high fat mixed meal challenge test to evaluate three existing models of postprandial adipose tissue dynamics. The inclusion of a palmitate stable isotope tracer in the meal provides a first opportunity to evaluate model terms describing the postprandial spill-over of LPL derived NEFA on human data. Ultimately, as no existing model could sufficiently describe the calculated metabolite fluxes, we construct a new refined model of postprandial adipose tissue metabolism which could be parameterised by data; introducing novel terms explicitly accounting for the conversion of glucose to G-3-P for use in re-esterification and the postprandial uptake of glycerol into the adipose tissue. The resulting model was then used to quantitatively estimate the impact of caloric restriction induced weight loss on subcutaneous adipose tissue fatty acid dynamics in a population of sixteen overweight or obese individuals that participated in the Yoyo study [17]. The model also allows for the prediction of the dynamic response of several reactions, which were not directly measured, to the ingestion of a meal (e.g. uptake and release of NEFA and glycerol by the adipose compartment), providing additional insights into the alterations in adipose tissue metabolism following caloric restriction which contribute to the observed changes in the calculated adipose tissue fluxes.

Study design and work-flow
The study design and work-flow is visualised in Fig 1 and detailed below. In summary, three existing models of postprandial adipose tissue glucose, fatty acid, and triglyceride dynamics were evaluated using arteriovenous measurements across the abdominal subcutaneous adipose tissue. Subsequently, a refined model was constructed by either selecting the best fitting model term of the three existing models or, in cases where none of the existing models allowed for an adequate fit to the data, new or modified physiological mechanisms were introduced. Finally, model parameters were estimated before and after weight loss in the Yoyo study [17] and compared, taking their confidence intervals into account.
(BMI > 25kg/m 2 ) and obese (BMI > 30kg/m 2 ) individuals (BMI 28-35 kg/m 2 ) were randomly assigned to one of two diet interventions designed to achieve a weight loss of approximately 10%, a low calorie diet of 1250 kcal/day over twelve weeks (slow weight loss) or a very low calorie diet of 500 kcal/day for five weeks (rapid weight loss). Following the diet intervention all participants underwent a four week weight stabilisation period with a diet based on the energy requirements of each participant. As part of this study, sixteen individuals underwent a high fat mixed-meal challenge test (milk shake) containing 100 mg of [U-13 C] palmitate tracer at baseline and following the weight stabilisation period [17]. Fasting samples were collected at -30 and 0 minutes. The meal was consumed at 0 minutes, participants were asked to consume the shake within 10 minutes, and samples were collected at 60, 120, 180, 240, 300 minutes postprandially from an arterialised dorsal hand vein and the superficial epigastric vein draining the abdominal subcutaneous adipose tissue. Abdominal subcutaneous adipose tissue blood flow was measured using the 133 Xe washout technique [19,20]. study. Existing model terms are first evaluated based on their ability to describe the experimentally measured data. In addition, as the refined model is constructed with the aim that parameters could be estimated directly from measured data, existing model terms are also compared using identifiability analysis, using the Profile Likelihood, sensitivity analysis, and statistical evaluation using Akaike Information Criterion. The refined model of adipose tissue metabolism is constructed either by selecting the best fitting of the existing model terms, or by introducing modifications and, where necessary, novel physiological mechanisms. The resulting refined model is then parameterised using data from both the baseline and following weight stabilisation from the Yoyo study, and the values compared. https://doi.org/10.1371/journal.pcbi.1007400.g001

Ethics statement
All subjects gave their written informed consent before participation in the Yoyo Study. The Yoyo Study was performed according to the Declaration of Helsinki and was approved by the Medical Ethics Committee of Maastricht University Medical Centre (METC 11-3-066) with the approval number NL38099.068.11.

Calculations
Individual metabolite fluxes across the abdominal subcutaneous adipose tissue were computed for NEFA, triglyceride, glucose, and glycerol by multiplying the arterio-venous plasma concentration difference of each measured metabolite by the adipose tissue blood flow [17].
Plasma tracer concentrations in NEFA and TG fractions were computed as described in the original study [17]. The rate of fractional spill-over of NEFA derived from LPL lipolysis of chylomicron triglyceride was calculated as one minus the rate of fractional extraction as described by Bickerton et al. [5] Two indices of insulin resistance were calculated using the arterial metabolite measurements. HOMA-IR is a measure of whole body insulin resistance, and is calculated as fasting plasma glucose(mmol/l) times the fasting plasma insulin (μ U/ml) divided by 22.4 [21]. Several studied have evaluated appropriate HOMA-IR cutt-off values for determining insulin resistance, with identified cut-off calues ranging from 1.85 to 2.01 depending on the study population [22]. ADIPO-IR is a surrogate measure of adipose tissue specific insulin resistance and is calculated as fasting plasma NEFA concentration (mmol/l) times the fasting plasma insulin concentration (pmol/l) [23].

Existing models
Models of postprandial lipid metabolism were assessed to identify models that described insulin mediated adipose tissue specific dynamics of glucose, triglyceride, NEFA, and/or glycerol. Three published models describing complementary aspects of adipose tissue metabolism were included in the analysis, namely the Jelic [11], Pratt [13], and Sips [10] models. The Pratt Model is a large, multi-compartmental, computational model describing the dynamics of several metabolite species including glucose, NEFA, glycerol, pyruvate, and both endogenous and dietary triglyceride in the fasting and postprandial states in plasma, skeletal muscle, adipose tissue, and liver [13]. The Jelic model is a two compartment, physiology-based, mathematical model describing NEFA dynamics across the plasma and a lumped interstitial adipose space in the postprandial state [11]. The Jelic model also accounts for two physiologically relevant delays in insulin signalling not present in the Pratt model. The model parameters for both the Jelic and Pratt models have been taken from literature and both models have been validated using oral glucose tolerance test (OGTT) and mix meal challenge test data in lean and abdominally obese individuals. The Sips model [10] is an extension of the Dalla-Man model [24] of postprandial glucose-insulin interplay to include NEFA kinetics. Unlike the previous two models, all parameters in the Sips model have been estimated from data. Where possible parameter values were maintained at the values provided for the lean healthy individual in the Dalla-Man model, as these parameters had been validated using gold standard triple-glucose tracer data [20]. The remaining 21 parameters were estimated from experimental data consisting of plasma measurements of several metabolites from individuals under various clamp conditions [25,26] and also a frequently sampled oral glucose tolerance test (OGTT) and oral lipid tolerance test [10].

Model evaluation and refinement
The three existing models (Jelic, Pratt, Sips) were decomposed into relevant subunits describing each measured metabolite's dynamics across and within the adipose tissue using a so-called divide and conquer approach [27] (Fig 1, S1 Table). These model terms were then evaluated by comparison to the calculated metabolite flux from the A-V data (average μ = (μ 1 ,. . ., μ 7 ) and standard deviation σ = (σ 1 ,. . ., σ 7 )) according to the following error function C(p): Here M(p, t i ) is the model prediction at time t i , t = (−30, 0, 60, 120, 180, 240, 300), for a parameter set p.
Optimal parameter sets for each term were obtained through non-linear regression, by finding the parameter set which best described the measured flux, minimising the error term C(p). Optimisation of parameters was performed using lsqnonlin (MATLAB 2014b, The MathWorks Inc., Natick, Massachusetts, United States) a local, gradient-based least square solver. The parameter values supplied in the original model studies were used as initial values for the parameter search. Parameterisation of the final, refined model for comparison of baseline and following weight stabilisation was performed using a combination of a global and local search algorithm. Controlled Random Search [28] with 250 randomly selected initial parameter sets was utilised to search the parameter space in order to provide a good initial value for lsqnonlin.
Terms describing the uptake and release of metabolites from tissues other than the adipose tissue, for example appearance from the gut and insulin secretion, are not described in this model. Available measured arterial concentrations of insulin, glucose, NEFA, triglyceride, and glycerol are supplied as dependent inputs to the model terms [29].
In the cases when the existing models were not capable of describing the experimentally estimated fluxes, underlying assumptions and model terms were evaluated based on existing biological knowledge and modified accordingly to provide an improved description of postprandial metabolite dynamics across the adipose tissue. Additionally, all model terms were examined statistically, using the Akaike Information Criterion (AIC c ) corrected for small sample sizes [30], to select the most parsimonious model that is both biologically sound and can adequately describe the measured flux data.
The above procedure was performed to evaluate terms describing the triglyceride, glucose, glycerol, and NEFA fluxes across the adipose tissue. Terms describing the fractional spill-over of NEFA by LPL lipolysis could also be assessed for the first time using the [U-13 C] palmitate tracer data.

Parameter identifiability
Model terms were also evaluated for the identifiability of parameters, as a primary objective of this study was to provide a model which could be parameterised by experimental data and could therefore be used to quantify adipose tissue metabolism from meal challenge test data. Identifiability of model parameters was evaluated using Profile Likelihood Analysis [31], whereby one parameter was varied iteratively from its optimal value and the remaining parameters were re-estimated. For an identifiable parameter the error measure C(p) would be expected to increase as the parameter value deviated from its optimum [31,32]. 95% confidence intervals for the estimated parameter values have been estimated using the parameter value covariance matrix approximated using the Jacobian matrix provided as output of the lsqnonlin algorithm used in the parameter estimation procedure.

Results
The three existing models of postprandial adipose tissue metabolism implemented in this study were not capable of describing the calculated postprandial adipose tissue metabolite fluxes. Consequently, a refined model of adipose tissue postprandial metabolism was constructed. The refined model is formulated as a two compartment model consisting of a blood plasma compartment and a lumped interstitial adipose compartment reflecting the sampling method used in the arteriovenous data (Fig 2). In the following, we show the results of the Structure of refined mathematical model of adipose tissue metabolism. The refined model of adipose tissue metabolism consists of a two compartmental model, describing dynamics between the plasma and a lumped interstitial adipose space. The insulin stimulated LPL lipolysis of circulating triglycerides releases glycerol and NEFA. The hydrolysed NEFA passes into the adipose space, with an insulin-dependent fraction spilling over into the plasma. The insulin inhibited lipolysis of the triglycerides stored within the adipose tissue releases NEFA and glycerol. It is assumed that this glycerol cannot be recycled within the adipose space and enters the plasma for transportation to the liver. Glucose passes into the adipose space at an insulin dependent and independent rate. Glucose is converted to G-3-P and provides the glycerol backbone necessary for re-esterification of NEFA for storage as triglyceride within the adipose space. Novel model terms, introduced in this analysis, are shown in red, existing models are shown in black. Reactions that are stimulated by insulin are depicted with a green insulin symbol. Reactions that are inhibited by insulin are shown with a red insulin symbol.
https://doi.org/10.1371/journal.pcbi.1007400.g002 model evaluation and refinement separately for each measured metabolite using the adipose tissue flux data before the weight loss intervention.

Triglyceride flux
Circulating plasma triglyceride is hydrolysed by LPL at the endothelial wall releasing NEFA and glycerol. The Pratt model describes this reaction as being proportional to the concentration of triglyceride present in circulating lipoproteins (both chylomicron and VLDL) and occurring at a basal and insulin stimulated rate using linear terms [13], with the plasma insulin concentration having an instantaneous effect. The Jelic model accounts for the saturation of enzyme mediated LPL lipolysis using Michaelis-Menten kinetics, with the rate of hydrolysis of circulating triglyceride being dependent on the plasma triglyceride concentration and stimulated by a delayed insulin effect [11]. The transcription of LPL and subsequent secretion of the LPL protein to the endothelial wall are known to be stimulated by insulin, however, these processes take some time [4]. As a result, the Jelic model introduced a three-fold insulin delay term, with the time delay parameter set to 240 mins [11]. The Sips model makes use of the Jelic model term and parameters [10].
Using the parameter values as specified in the original publications, neither the Jelic nor Pratt models are capable of describing the mean measured triglyceride flux across the adipose tissue ( Fig 3A). Allowing the parameters to be estimated from the data, the Jelic model provides an improved fit. A more detailed investigation demonstrated that the insulin delay is the essential component missing in the Pratt model (S2 and S3 Figs). As a result, the optimised triglyceride flux term ([TG] flux ) shown below was derived, making use of the insulin dependent linear approximation from the Pratt model while introducing the insulin delay term of the Jelic model.
Where K ad is the rate parameter for LPL lipolysis, [TG art ] is the arterial triglyceride concentration and [I LPL ] is the delayed insulin signal, modelled with a three compartmental delay (Eq 12) with a time delay parameter τ LPL also estimated from data. (Fig 4) Fractional spill-over NEFA released by LPL lipolysis of circulating triglyceride is taken up by the adipose tissue. However, this process is inherently leaky with a proportion of the released NEFA spilling-over directly into the plasma [33,34]. Arteriovenous studies combined with NEFA stable isotope tracers have previously demonstrated that the fractional spill-over of NEFA from LPL lipolysis increases later in the postprandial period [5,33]. The Jelic model does not describe the spill-over of NEFA derived from LPL lipolysis [11]. The Pratt model describes fractional spill-over at a constant rate of 25% [13]. The Sips model describes the fractional spill-over of LPL derived NEFA as occurring at a basal and insulin inhibited rate, using the same delayed insulin signal as LPL lipolysis with a 240 min time delay parameter [10].
The Sips term could not describe the increasing postprandial fractional spill-over measured in the Yoyo Study data (Fig 3B). The 240 min insulin time delay served to over-damp the insulin signal. As a result, the insulin delay was removed. Further analysis of the Sips term with use of the Akaike Information Criterion led to the removal of the basal fractional spill-over rate, yielding the following optimised fractional spill-over term.
Here I B is the basal plasma insulin concentration, [I art ] the arterial insulin concentration and D spill an estimated constant ( Fig 4B).

Glucose flux
Glucose enters the adipose tissue along the concentration gradient facilitated by insulin dependent (GLUT4) and insulin independent (GLUT1) transporters [35]. Within the adipocyte, glucose is quickly converted to glucose-6-phosphate (G-6-P) which serves to trap the glucose within the adipocyte. Glucose is the primary source of the glyceraldehyde-3-phosphate (G-3-P) backbone needed for re-esterification of adipocyte NEFA [36]. Both the Pratt [13] and the Sips model [10] describes the uptake of glucose by the adipose tissue as occurring at an insulin dependent and independent rate. However, the Pratt model uses direct plasma insulin stimulation rather than accounting for a delay in the effect of insulin signalling which induces the translocation of GLUT4 transporters from the transport vesicles of the cell to the membrane. The Jelic model does not account for glucose dynamics. With the introduction of a time delay in insulin stimulation the linear term from the Pratt model was capable of describing the glucose flux. In addition, it was necessary to explicitly account for the conversion of glucose to G-6-P in order to remedy an erroneous prediction of a postprandial glucose efflux from the interstitial adipose space (S3C Fig).
Where [G art ] is the arterial glucose concentration, [I AT ] is the adipose tissue delayed insulin signal described using a threefold delay with a time delay parameter τ AT (Eq 12), GLUT1 and GLUT4 are constants estimated from the data (Fig 4C).

Glycerol flux
Glycerol is released into the plasma by LPL lipolysis of circulating triglyceride at the endothelial wall [35]. Glycerol is also released by lipolysis of triglyceride stored within the adipose tissue (denoted as ATL lipolysis [37]). It is commonly assumed tht due to the inactivity of glycerol kinase in the adipose tissue [38], all glycerol released by ATL lipolysis enters the plasma for transportation to the liver. Postprandial glycerol dynamics are not accounted for in the Jelic [11] nor Sips models [10]. The Pratt model accounts for the release of glycerol by ATL lipolysis within the adipose tissue which appears directly in the liver compartment [13]. Glycerol release into plasma by LPL lipolysis is not described in this model. Analysis of the postprandial glycerol and triglyceride fluxes from the YoYo Study indicate the uptake of glycerol into the interstitial adipose space in the postprandial state ( Fig 5). As a result, the Pratt glycerol term was extended to account for the concentration gradient diffusion of plasma glycerol into and from the interstitial adipose compartment.
Where ATL lipolysis is described as occurring at a basal (B ATL ) and saturable insulin inhibited rate using terms from the Jelic model (ATL max , K ATL ) [11]. Where [I AT ] is the delayed interstitial adipose compartment insulin signal as described for the glucose flux above. (Figs  3D and 4D)

NEFA flux
In the fasting state NEFA is released by the lipolysis of stored triglyceride within the adipocyte and passes into the plasma for delivery to other tissues. The lipolysis of stored triglyceride is inhibited by increasing insulin concentration in the postprandial state [39,40] which also stimulates the rate of re-esterification of free NEFA within the adipocyte for storage as triglyceride [41]. In addition, insulin stimulates the LPL lipolysis of circulating triglyceride, removing excess dietary triglyceride from the plasma. While the majority of NEFA released by LPL lipolysis passes into the adipocyte for re-esterification a fraction spills-over directly into the plasma [33,34]. Terms describing the fractional spill-over, LPL and ATL lipolysis have previously been derived. Terms describing the rate of re-esterification were compared using the Yoyo study NEFA flux. The Pratt model describes the rate of re-esterification as a linear term, directly stimulated by plasma insulin and adipose glucose concentration [13]. The Jelic model uses saturable kinetics coupled with a delayed insulin stimulation [11]. Both models used passive diffusion of NEFA between the plasma and interstitial adipose compartment. The Sips model does not explicitly describe the adipose compartment, instead using a linear term which accounts for net lipolysis and re-esterification of NEFA within the adipose tissue [10].
While the Jelic term could produce a good describtion of the NEFA flux when fit to the data (S2 Fig), it assumed an unlimited supply of G-3-P and required the estimation of five parameters from the data (S1 Table). The linear re-esterification term of the Pratt model was modified to include the delayed interstitial adipose compartment insulin signal. In addition, terms explicitly accounting for the production of the G-3-P backbone necessary for re-esterification from glucose taken into the adipose tissue were introduced, providing a sufficient description of the measured postprandial NEFA flux (Figs 3E and 4E).

d½G-3-P AT
Re-esterification occurs at an insulin stimulated linear rate K reester dependent on the available concentrations of G-3-P and NEFA. A portion of the glucose taken up by the adipose space (frac use ) is converted to G-6-P and then to adipose tissue G-3-P ([G-3-P AT ]) which will be available for use in re-esterification, the remaining glucose leaves the system for use by other cellular functions which are not explicitly described here. The production of adipose tissue G-3-P is described by the term [G-3-P pro ] with a two compartmental delay governed by the delay constant τ G-3-P .

Insulin delays
Insulin stimulates the translocation, secretion, and in some cases transcription of several enzymes and transport proteins involved in adipose tissue metabolism. These processes occur over time spans of several minutes to several hours. The effectiveness of the measured plasma insulin signal may be further damped by the presence of insulin resistance, a condition where the responsiveness of tissues to plasma insulin is reduced. Rather than describing the full sequence of reactions involved the dampening of the plasma insulin signal is approximated using a three compartmental delay, as in the Jelic model [11].
Where [I art ] is the measured arterial insulin concentration and τ is the respective time delay parameter. Two insulin time delays are simulated, a long time delay for LPL lipolysis (τ = τ LPL ) and a shorter interstitial adipose compartment time delay (τ = τ AT ). The structure of the insulin delay terms are not altered from the original model [11]. However, the time delay parameters are not fixed, as in the Jelic and Sips models, and are estimated from the data.

Application of adipose model-Quantification of caloric restriction
The resulting refined model was then used to quantify the effect of caloric restriction using meal challenge test data collected before and after a diet intervention as part of the Yoyo study. Parameter values before and after caloric restriction. Parameter sets were estimated from data collected at baseline and after a weight stabilisation period following caloric restriction (Table 1, the complete parameter sets are provided in S2 Table). The rate parameters for both glycerol and NEFA concentration gradient based diffusion from the plasma to the adipose space (P GLY and P NEFA ) increase significantly following caloric restriction, indicating an increase in the rate of glycerol and NEFA transport in and out of the adipose tissue. A large, although not significant, decrease can be seen in the insulin delay parameter for LPL lipolysis of circulating triglyceride, τ LPL , (157 mins at baseline to 113 mins following weight stabilisation) which is not accompanied by a strong change in the value for K ad , the rate parameter for that term. In addition, there is a four minute decrease in the adipose tissue insulin delay parameter τ AT . The reduction in the insulin delay parameters indicates a less damped response to insulin signalling following caloric restriction, indicative of improved insulin sensitivity.
Decomposition of glycerol and NEFA model predictions. The estimated values for P GLY and P NEFA increased significantly following caloric restriction, but it remained unclear if the observed differences in the glycerol and NEFA fluxes following caloric restriction (Fig 4) could be entirely attributed to these increases in the rates of uptake and release of both metabolites. In Fig 6, our refined model's predictions are used to decompose the calculated adipose glycerol and NEFA fluxes over the duration of the mixed meal challenge test into their constituent reactions. This allows for visualisation of how the rates of specific, unmeasured reactions change following caloric restriction. The model predicts an increase in the rate of ATL lipolysis within the adipose tissue following caloric restriction (Fig 6B and 6E, red lines), with a predicted 56% increase in the rate of lipolysis of stored triglyceride (0.229 μmol/100 ml tissue/min at baseline to 0.359 μmol/100 ml tissue/min following caloric restriction), resulting in the observed increase in efflux of both NEFA and glycerol from the adipose tissue in the fasting state ( Fig  6A and 6D, purple lines). An increase in ATL lipolysis would be expected following such a period of caloric restriction as the triglyceride stores in the adipose tissue are hydrolysed to supply NEFA for use in other tissues. In the later postprandial phase, the lower concentrations of circulating triglyceride follow caloric restriction (total AUC 642.2 μmol a baseline to just 463 μmol following the diet intervention), results in the reduced efflux of glycerol and NEFA from LPL lipolysis after 120 minutes (Fig 6A and 6D, yellow lines).
Interestingly, the peak in the rate of LPL occurs approximately 40 mins earlier (230 mins at baseline versus 190 mins following caloric restriction) reflecting the estimated decrease in τ LPL (Fig 6A and 6E, yellow lines). Note, there is no decrease in the peak time of measured arterial insulin (S4E Fig), however there is a decrease in both the fasting and postprandial plasma insulin concentrations (total AUC 6821.5 μU of insulin at baseline to 5291.8 μU following caloric restriction). Despite the decrease in circulating concentrations of insulin the postprandial inhibition of ATL lipolysis and stimulation of re-esterification by insulin occur at faster rates. The magnitude of slope describing the rate of postprandial inhibition of ATL lipolysis increases from 0.019 at baseline to 0.0023 following caloric restriction. Concurrently, the rate of stimulation of re-esterification increases from 0.0055 at baseline to 0.0078 following the diet intervention. These increases are indicative of improved responsiveness to circulating insulin, as reflected with the decrease in the estimated values for both insulin delay parameters.
Alternative measures of insulin resistance. Both whole body and adipose tissue specific insulin resistance were assessed using the HOMA-IR and ADIPO-IR indices respectively. HOMA-IR values decreased significantly from a mean values of 2.27 at baseline, which would be considered insulin resistant, to 1.51 following the diet intervention (p = 0.036). ADIPO-IR values also decreased significantly from 47.27 at baseline to 28.90 following caloric restriction (p = 0.006). Supporting the model prediction of an improvement in insulin sesntivity.

Discussion
We present a refined mathematical model of postprandial adipose tissue insulin mediated dynamics of glucose, NEFA, triglyceride, and glycerol. Our refined model elucidates the strong effect of insulin on adipose tissue lipid metabolism, building on the work of Jelic [11], Pratt [13], and Sips [10], with all measured adipose tissue metabolite fluxes being stimulated or inhibited by insulin. Our model also introduces several novel terms. We explicitly account for the conversion of glucose to G-3-P, which serves to trap glucose within the adipose space. Thereby, counteracting the erroneous predictions of a postprandial efflux of glucose from the adipose compartment in the existing models due to glucose accumulation. We also introduce a term accounting for the concentration gradient dependent uptake and release of glycerol. This mechanism had been postulated by Coppack et al. in their 2005 model [42], but has not been accounted for in the three subsequently published models that we have compared. Moreover, inclusion of a [U-13 C] palmitate stable isotope tracer in the ingested meal allowed us to estimate the postprandial fractional spill-over of LPL derived NEFA into the plasma. The corresponding model term has, to the best of our knowledge, not previously been validated using experimental data.
Arteriovenous measurements in the fasting and postprandial conditions across the abdominal subcutaneous adipose tissue were employed to evaluate three existing models of postprandial adipose tissue metabolism and the assumptions upon which they have been constructed. While the existing models included in this analysis perform well in describing dynamics between triglycerides, NEFA, and glucose at the whole-body level, none of the existing models could describe the adipose tissue specific fluxes well, leading to the construction of our refined model. Methods such as area under the curve (AUC) are most often employed to compare time series of metabolite concentrations during challenge tests, with each metabolite evaluated independently of the others. However, AUC fails to capture the dynamic properties of the postprandial metabolite curves. It is possible for several different response curves to have the same AUC value. With computational models the dynamic responses of all measured metabolites are used to parameterise a model in a physiologically meaningful way. Consequently, computational models could also prove to be a powerful tool for the interpretation and quantification of dynamic time series of data, as in the case of meal challenge test data. Our refined model was applied to time series data collected before and after a weight loss intervention study [17] to investigate the effects of caloric restriction on adipose tissue metabolism. Comparison of parameter values estimated from data collected before and after caloric restriction indicated a significant increase in the rate coefficient for the concentration gradient based transport of both glycerol and NEFA in and out of the adipose compartment. Use of the model to decompose the glycerol and NEFA flux predictions into their constituent reactions predicted an increase in the rate of ATL lipolysis within the adipose tissue following caloric restriction, resulting in the measured increased efflux of NEFA from the adipose tissue in the fasting state. This would be expected following a period of caloric restriction as triglyceride stored in the adipose tissue is hydrolysed and the resulting NEFA is transported for use in other tissues. A large, but non-significant, decrease was observed in the insulin time delay parameter for LPL lipolysis of circulating triglyceride, τ LPL , from 157 minutes at baseline to 113 minutes following caloric restriction, the effect of which can be observed with the peak in the rate of LPL lipolysis occurring approximately 40 minutes earlier following weight stabilisation. This is accompanied by a non-significant four minute decrease in the general adipose tissue insulin delay τ AT . These reductions indicate a less damped response to insulin stimulation following the diet intervention, which is in line with previous studies reporting a reduction in whole body insulin resistance following prolonged caloric restriction [43]. Whole body insulin sensitivity for the sixteen individuals, measured using HOMA-IR, decreased significantly following caloric restriction, from a mean of 2.27 at baseline to 1.51 following the diet intervention (p = 0.036). In addition, adipose tissue specific insulin resistance, assessed using ADIPO-IR, decreased significantly following the diet intervention, further supporting our model's prediction of improved insulin sensitivity following caloric restriction. It is also of note that the LPL time delay estimates are more than one, and close to two hours, shorter than the fixed 240 minutes time delay proposed by Jelic et al. [11] and subsequently used in the Sips model [10]. While changes could be observed in several other parameter values before and after caloric restriction, these differences were not significant. Previous analysis of the Yoyo study data using more traditional techniques (incremental area under the curve) found only a significant decrease in the fasting and postprandial triglyceride flux, however this was accompanied by reduced arterial triglyceride concentration [17]. Thus, we see the additional insights which can be gained through the fitting of computational models in the quantification of dynamic time series data [44].
Multiple arteriovenous studies have used the glycerol flux to infer information regarding rate of ATL lipolysis under the assumption that glycerol release by LPL lipolysis, equivalent to the instantaneous triglyceride flux, remains within the plasma, and due to the negligible activity of glycerol kinase within the adipose tissue [38], all glycerol relased by ATL lipolysis enters the plasma for transport to the liver. Thereby, the rate of ATL lipolysis can be estimted as the difference between the measured glycerol and triglyceride fluxes [45,46]. However, analysis of the postprandial glycerol and triglyceride fluxes from the Yoyo Study indicates that there is an influx of glycerol into the interstitial adipose space in the postprandial state, in contradiction with this assumption (Fig 5). The introduction of a term accounting for the concentration gradient dependent uptake and release of glycerol by the adipose tissue, which was first proposed by Coppack et al. [41], allows the model to described the measured postprandial influx of glycerol. In addition, in order to achieve the best fit of the flux data, in particular the measured postprandial influx of glycerol from the plasma to the adipose tissue, it was necessary to introduce a sink term in adipose tissue glycerol which cannot be accounted for with current biological understanding. In future work it may be possible with the integration of adipose tissue specific measurements of other omics data to determine the cause of the glycerol disappearance.
In combination with arterio-venous measurements, the use of stable isotope tracers allow for the quantification of reactions that are not directly measurable, such as rates of appearance and turn over. The meal administered in the Yoyo study includes 100mg of [U-13 C] palmitate stable isotope tracer which allows for the estimation of the rate of fractional spill-over of NEFA from LPL lipolysis of dietary triglycerides using the method as described by Bickerton et al. [5]. As the Yoyo study utilizes a single tracer we do not have the data to evaluate chylomicron and VLDL triglyceride separately. Therefore, we have grouped the terms and used the available measurements of fractional spill-over of chylomicron derived NEFA to be representative of fractional spill-over from total circulating triglyceride pool. Moreover, we assume the dynamics of palmitate is representative for the generic NEFA pool. Incorporation of additional fatty acid tracers under different challenge conditions would allow for the extension of our model to describe the complex systemic interplay between different lipoprotein classes and NEFA species. Inclusion of a palmitate stable isotope tracer in a meal has been shown to label dietary derived triglycerides and NEFA for approximately the first 120 to 180 minutes [5,47]. Due to recycling of NEFA both within the adipose tissue and eventual incorporation of labelled NEFA into VLDL, estimates of fractional spill-over become less reliable beyond this time. This is reflected with the increases in the measured standard error of mean for later time points in part B of Fig 3. To this end, error measurements have been weighted by the standard deviation within the 16 individuals such that a time point with a larger standard deviation will have a lower priority when estimating parameter values.
In order to produce a model which could be parameterised from data, certain model terms were simplified from those presented in the original models, e.g. use of a linear term in place of Michaelis-Menten dynamics to describe LPL lipolysis of circulating triglycerides. A linear approximation is sufficient in this situation to describe this reaction, as saturation is not expected to occur given the measured arterial triglyceride concentrations within the Yoyo study. The use of Profile Likelihood Analysis indicates that six of the fourteen model parameters are identifiable given the data, with a further six having upper or lower bounds (S3 Fig). For future implementation of the model, fixation of parameter values for which reliable experimental estimates are available would improve the estimation of parameters of interest for the biological questions being asked. As with previous computational models of adipose tissue metabolism, we have considered whole body adipose tissue as one homogenous unit. While the arteriovenous sampling of the abdominal subcutaneous tissue technically restricts our model evaluations to a single adipose tissue depot, we believe that our model sufficiently captures postprandial adipose tissue dynamics given that the measured depot, the abdominal subcutaneous adipose tissue, is considered the primary site for the storage and release of NEFA into systemic circulation.
Arteriovenous measurements across other tissue depots, such as skeletal muscle, presents an opportunity to further evaluate and, where necessary, refine existing model terms describing contributions of other tissues to glucose and NEFA homeostasis, as in the case of adipose tissue metabolism in this study. In this study we present a work-flow where, through the use of dependent inputs in combination with a divide and conquer approach, the model evaluation procedure can be reduced to a non-linear regression problem, greatly reducing the computational time.
Given the association of ectopic fat deposition with the development of several metabolic disorders, such as Non-Alcoholic Fatty Liver Disease and Type 2 Diabetes, the study of the dysregulation of adipose tissue metabolism has garnered much attention in recent years. Our refined model is capable of simulating in silico the dynamics of human in vivo adipose tissue metabolism in both the fasting and postprandial state. Through the variation of model parameters or modulation of dependent inputs, such as the arterial insulin concentration, it is possible to simulate the dynamics in adipose tissue metabolism in response to different external stimuli and gain insight into potential sources of dysregulation in adipose metabolism. Furthermore, embedding our refined adipose tissue model into a larger whole body model, with terms describing the contributions of other tissues, would allow for the investigation how the model simulated dysregulation in adipose tissue metabolism would impact on the glucose and NEFA dynamics in other tissues and whole body metabolism.
In conclusion, we present a refined physiology-based computational model of adipose tissue metabolism which has been shown, using arteriovenous measurements across the abdominal subcutaneous adipose tissue, to outperform several existing models of adipose tissue metabolism. Our model elucidates the strong influence of insulin signalling on adipose tissue dynamics, particularly the cycling between storage and release of NEFA in the fasting and postprandial states. Application of our model to data collected before and after a diet intervention allows for quantification of the effect of caloric restriction on adipose tissue metabolism. Estimated parameter values indicate that the delays in insulin effectiveness in the system are not fixed, with the estimated time delays in insulin signalling not only differing from the values used in previous models, but also decreasing following the diet intervention suggesting an improvement in adipose tissue insulin sensitivity following caloric restriction.
Supporting information S1 Fig. Results of Profile Likelihood Analysis for each parameter of the refined model. The parameter values estimated from the baseline data are shown with a red cross and the value of C(p) resulting from iteratively adjusting the parameter value and re-estimating the parameter values indicated by the blue line. A parabola with the parameter estimate at it's base (as in the case of K ad , τ LPL , D spill , τ AT , P GLY , and P NEFA ) indicates an identifiable parameter. Profile likelihood for several other parameters (GLUT1, GLUT4, B ATL , ATL max and τ G3P ) indicate the existence of an upper bound, these parameters have been bound below by zero in the parameter estimation procedure for physiologically relevant reasons. Similarly the parameter describing the fractional usage of glucose to for G-3-P production in re-esterification is bound above by one. Finally two parameters (K ATL , and K reester ) appear to be practically non-identifiable. Given the product of K reester and the model predicted concentration of G-3-P in the adipose tissue are equal to the maximum rate of G3P production it is unsurprising that the parameter K reester in non-identifiable. Any change in the value of K reest in compensated for by a corresponding change in the model predicted concentration of adipose G-3-P. Jelic model is capable of producing quite a good fit of the triglyceride and NEFA fluxes. However, the Jelic model does not include terms to describe the frational spill over, glucose or glycerol fluxes. While the Pratt model can produce a qualitatively good fit of the glucose, NEFA, and glycerol fluxes it cannot produce a good fit of the triglyceride flux (A). The Jelic, Sips, and refined model descibed LPL lipolysis as being dependent on insulin, with a delayed insulin signal. The Pratt model also includes a term describing the stimulation of LPL lipolysis by insulin. However, it does not account for any delay in insulin signalling, with LPL lipolysis being directly stimulated by plasma insulin. The Pratt model also assumes that the contribution of insulin dependent lipolysis to the overall lipolysis of circulating triglyceride is neglible, consequently, the rate of LPL lipolysis is primarily determined by the circulating triglyceride concentration. (TIF) S3 Fig. Increasing the influence of insulin dependent LPL lipolysis in the Pratt model. In the above figure we increase the weight of the contribution of insulin dependent LPL lipolysis, such that insulin dependent LPL lipolysis accounts for 1% (blue line) and 10% (yellow line) of the total adipose tissue triglyceride flux while maintaining the other parameters at the values provided in the original publication. As the Pratt model uses direct plasma insulin stimulation rather than accounting for delays in insulin signalling, as in the Jelic, Sips, and refined models, the model simulated triglyceride flux begins to peak too early under the influence of plasma insulin. The refined model makes use of the LPL lipolysis term from the Pratt model, but introduces the three compartmental delay from the Jelic and Sips models. With this delay in the insulin signalling, the refined model can produce a good fit to the triglyceride flux data at baseline and following weight stabilsation (red line).  Table. Comparison of terms from the Jelic, Pratt, Sips, and refined adipose tissue models. Terms describing individual metabolite fluxes across the adipose tissue were extracted from the Jelic, Pratt, Sips, and refined model and compared. (PDF) S2 Table. Parameter values estimated for data collected at baseline and following caloric restriction. Complete set of parameter values estimated by fitting the refined model to the calculated adipose tissue metabolite fluxes at baseline and following weight stabilisation. 95% confidence intervals for parameter estimates are displayed in parentheses below the estimated value. The coloured boxes indicate the model term, or terms, in which each parameter appears. � All parameters were bound below by zero during parameter estimation, however the method for calculating the confidence intervals assumes the confidence interval is symmertric about the estimated parameter values. † terms that decribe a fractional value were also bound above by 1 during parameter estimation. (PDF) S1 File. Sensitivity analysis for refined model.