Modeling of free fatty acid dynamics: insulin and nicotinic acid resistance under acute and chronic treatments

Nicotinic acid (NiAc) is a potent inhibitor of adipose tissue lipolysis. Acute administration results in a rapid reduction of plasma free fatty acid (FFA) concentrations. Sustained NiAc exposure is associated with tolerance development (drug resistance) and complete adaptation (FFA returning to pretreatment levels). We conducted a meta-analysis on a rich pre-clinical data set of the NiAc–FFA interaction to establish the acute and chronic exposure-response relations from a macro perspective. The data were analyzed using a nonlinear mixed-effects framework. We also developed a new turnover model that describes the adaptation seen in plasma FFA concentrations in lean Sprague–Dawley and obese Zucker rats following acute and chronic NiAc exposure. The adaptive mechanisms within the system were described using integral control systems and dynamic efficacies in the traditional \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$I_{\text{max}}$$\end{document}Imax model. Insulin was incorporated in parallel with NiAc as the main endogenous co-variate of FFA dynamics. The model captured profound insulin resistance and complete drug resistance in obese rats. The efficacy of NiAc as an inhibitor of FFA release went from 1 to approximately 0 during sustained exposure in obese rats. The potency of NiAc as an inhibitor of insulin and of FFA release was estimated to be 0.338 and 0.436 \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${\upmu {\text{M}}}$$\end{document}μM, respectively, in obese rats. A range of dosing regimens was analyzed and predictions made for optimizing NiAc delivery to minimize FFA exposure. Given the exposure levels of the experiments, the importance of washout periods in-between NiAc infusions was illustrated. The washout periods should be \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\sim$$\end{document}∼2 h longer than the infusions in order to optimize 24 h lowering of FFA in rats. However, the predicted concentration-response relationships suggests that higher AUC reductions might be attained at lower NiAc exposures.


Introduction
Nicotinic acid (NiAc; or niacin) has long been used to treat dyslipidemia [1,2]. When given in large doses (1-3 g/day), NiAc improves the plasma lipid profile by reducing total cholesterol, triglycerides, low-density lipoprotein cholesterol, and very-low-density lipoprotein cholesterol, and increasing levels of high-density lipoprotein cholesterol [3]. Moreover, by binding to the G-protein coupled receptor GPR109A, NiAc potently inhibits lipolysis in adipose tissue, leading to decreased plasma free fatty acid (FFA) concentrations [4,5]. The mechanisms of NiAc-induced antilipolysis have been thoroughly analyzed in previous studies [6][7][8][9]. Chronically elevated plasma FFA concentrations are associated with several metabolic diseases, including insulin resistance [10][11][12]; NiAc-induced FFA lowering is a potential approach to ameliorating these conditions. However, the clinically applied dosing regimens have not been designed to lower FFA; rather, the goal has been to ameliorate dyslipidemia [13].
Although acute administration of NiAc results in rapid reduction in FFA concentrations [2,14], long-term infusions are associated with tolerance development (drug resistance) and plasma FFA concentrations returning to pre-treatment levels (complete adaptation) [15]. Furthermore, abrupt cessation of the NiAc infusions produces an FFA rebound that overshoots the pre-infusion levels [9,15]. Numerous studies have sought to quantitatively determine the acute concentration-response relationship between NiAc and FFA [14,[16][17][18][19][20][21][22][23]. The acute NiAc-induced FFA response has been successfully characterized using pharmacokinetic/pharmacodynamic (PK/PD) models, but these models fail to describe the complete return of FFA to pretreatment levels associated with chronic NiAc treatment. Thus, an improved model is required in order to predict optimal treatment regimens, aimed to achieve durable NiAc-induced FFA lowering.
In this study, we sought to further develop the concepts used in previous analyses to develop a more general NiAc-FFA interaction model-applicable to a large set of dosing regimens and NiAc exposure durations. The model was also aimed at quantitatively determining the impact of disease on the FFA-insulin system and to provide predictions for optimal drug delivery. We conducted a meta-analysis on a rich pre-clinical data set of the interaction between NiAc and FFA, as well as insulin, in a nonlinear mixed-effects (NLME) modeling framework. Using various routes and modes of NiAc provocations, we collected concentration-time course data of NiAc (drug kinetics), insulin and FFA (drug-induced dynamics). Experiments were done both in lean Sprague-Dawley and obese Zucker rats-allowing disease impact to be evaluated. Furthermore, by including insulin as a covariate of the FFA response, we could quantitatively analyze the endogenous antilipolytic effects of insulin [24] under NiAc provocations. Moreover, optimal dosing regimens, consisting of constant rate infusion periods followed by washout periods, were investigated.

Animals
Male Sprague Dawley (lean) and Zucker rats (fa/fa, obese) were purchased from (conscious groups) Harlan Laboratories B.V. (The Netherlands) or (anesthetized groups) Charles River Laboratories (USA). Experimental procedures were approved by the local Ethics Committee for Animal Experimentation (Gothenburg region, Sweden). Rats were housed in an Association for Assessment and Accreditation of Laboratory Animal Care accredited facility with environmental control: 20-22 C, relative humidity 40-60%, and 12 h light-dark cycle. During acclimatization ( ! 5 days), animals were housed in groups of 5 with free access to both water and standard rodent chow (R70, Laktamin AB, Stockholm, Sweden).

Surgical preparations
To prevent potential infections in conjunction with surgery, oral antibiotics were given 1 day before pump/catheter surgery and then once daily for 3 days (sulfamethoxazole and trimethoprim 40 ? 8 mg mL À1 ; Bactrim Ò , 0.2mL / animal, Roche Ltd, Basel, Switzerland). Surgery was performed under isoflurane (Forene Ò , Abbott Scandinavia AB, Solna, Sweden) anesthesia, with body temperature maintained at 37 C. For NiAc/saline administration, a programmable mini pump (iPrecio Ò SMP200 Micro Infusion Pump, Primetech Corporation, Tokyo, Japan) was implanted subcutaneously, via a dorsal skin incision. To allow blood sampling during the terminal experiment (conscious animals only), a polyurethane catheter (Instech Laboratories Inc, Plymouth Meeting, PA USA) was placed in the right jugular vein via an incision in the neck. In order to maintain its patency up to the acute experiment, the jugular catheter was filled with sterile 45.5% (wt/wt) PVP (polyvinylpyrrolidone, K30, MW $ 40,000 Fluka, Sigma-Aldrich, Sweden) dissolved in a sodium-citrate solution (20.6 mmol), sealed and exteriorized at the nape of the neck. Each animal received a post-operative, subcutaneous analgesic injection (buprenorphine, Temgesic Ò , 1.85 lg kg À1 , RB Pharmaceuticals Ltd, Berkshire, GB). Animals were then housed individually and allowed three days of recovery before the start of the pre-programmed pump infusion. Throughout the study, body weight and general health status were monitored and recorded daily.

Nicotinic acid exposure selection and formulation
A key aspect of the study design was to achieve plateau plasma nicotinic acid (NiAc) concentrations corresponding to therapeutically relevant levels in the rat ( $ 1 lM), based on the relationship between plasma NiAc levels and FFA lowering [16]. For intravenous infusions (i.v.), NiAc (pyridine-3-carboxylic acid, Sigma-Aldrich, St. Louis, MO, USA) was dissolved in sterile saline. For subcutaneous (s.c.) infusions, NiAc was dissolved in sterile water and adjusted to physiological pH using sodium hydroxide. Vehicle, for control animals, consisted of sodium chloride solutions at equimolar concentrations. Freshly prepared formulations were loaded into the infusion pump (see below) via a 0.2 lm sterile filter (Acrodisc Ò , Pall Corporation, Ann Arbor, MI, USA) just before pump implantation.

Experimental protocols
Conscious animals (NiAc naïve, Cont. NiAc and Inter. NiAc groups) Both lean and obese animals were divided into 3 dose groups and NiAc was given either acutely (NiAc naïve) or following 5 days of either continuous (Cont. NiAc) or intermittent (Inter. NiAc) administration. Each dose group was matched with corresponding saline infused controls. NiAc infusions were given subcutaneously at 170 nmol min À1 kg À1 . The intermittent infusion protocol was programmed as a 12 h onoff cycle (infusion on at 13:00). Following overnight fast, in the morning of the acute experimental day, the jugular catheter was connected to a swivel system to enable blood sampling in unrestrained animals. Jugular catheter patency was maintained by continuous infusion (5 lmol min À1 ) of sodium-citrate solution (20.6 mM). After a 3-4 h adaptation period, at $ 12:00, the basal phase of the acute experiment commenced with 2-3 blood samples drawn between -60 and -5 min, relative to start of NiAc/saline infusion (note that, in the Cont. NiAc groups, infusion pumps were on throughout this sampling period). Blood samples (16-17/ animal) were drawn under an 8 h experimental period. Samples, 30-150 ll (with total loss less than 5% of blood volume), were collected in potassium-EDTA tubes, centrifuged and plasma stored at -80 C pending analysis for NiAc, FFA and insulin.

Anesthetized animals (NiAc Off and NiAc Stp-Dwn 12 h infusion groups)
Before the infusions began, lean and obese rats were fasted for 8 h. On the day of the acute study, at 01:00 (corresponding to time = 0 h), the implanted pre-programmed pump began infusing NiAc at a constant rate of 170 nmol min À1 kg À1 for 12 h. At 8.5 h animals were anesthetized (Na-thiobutabarbitol, Inactin Ò , 180 mg kg À1 , i.p., RBI, Natick, MA, USA), underwent a tracheotomy with PE 240 tubing, and breathed spontaneously. One catheter (PE 50 tubing) was placed in the left carotid artery for blood sampling and for recording arterial blood pressure and heart rate. One catheter (PE 10 tubing) was placed in the right external jugular vein to infuse top-up doses of anesthetic. The arterial catheter patency was maintained by continuous infusion of sodium-citrate (20.6 mM in saline, 5 ll min À1 ) from shortly after carotid catheterization until the experiment ended. Body temperature was monitored using a rectal thermocouple and maintained at 37.5 C by means of servo controlled external heating. After surgery, animals were allowed a stabilization period of at least 1.5 h and blood sampling began at 11.0 h. At 12.0 h, NiAc infusion was either programmed to switch off (NiAc Off) or to decrease in a step-wise manner, with final switch-off at 15.5 h (NiAc Stp-Dwn). The step-down NiAc infusion rates were 88.9, 58.3, 43.7, 34.0, 24.3, 17.0, and 9.7 nmol min À1 kg À1 . All NiAc protocols were matched with saline-infused controls. Blood samples (18/animal) were drawn during a 6 h experimental period. Samples, 30-150 ll (with total loss less than 5% of blood volume), were collected in potassium-EDTA tubes, centrifuged, and plasma was stored at -80 C pending analysis for NiAc, FFA and insulin.

Anesthetized animals (NiAc Off and NiAc Stp-Dwn 1 h infusion groups)
After an overnight fast, lean and obese rats were anesthetized and surgically prepared, as described above. They were allowed a stabilization period after surgery of at least 1.5 h. Two basal blood samples were obtained, after which an i.v. NiAc infusion was given at a constant rate (170 nmol min À1 kg À1 ) for 1.0 h (the start of infusion was taken as time = 0 h). The NiAc infusion was then either switched off (NiAc-Off 1 h) or decreased in a step-wise manner, with final switch-off at 4.5 h (NiAc Stp-Dwn 1 h). The step-down NiAc infusion rates were: 31.1, 20.4, 15.3, 11.9, 8.50, 5.95 and 3.40 nmol min À1 kg À1 . All NiAc protocols were matched with saline infused controls. Blood samples (13-18/animal) were drawn during a 6 h experimental period. Samples, 30-150 ll (with total loss less than 5% of blood volume), were collected in potassium-EDTA tubes, centrifuged, and plasma was stored at -80 C pending analysis for NiAc, FFA, and insulin. All of the experimental groups are summarized in Table 1.

Analytical methods
Plasma FFA was analyzed using an enzymatic colorimetric method (Wako Chemicals GmbH, Neuss, Germany). Plasma insulin from obese rats was analyzed with a radioimmunoassay kit (rat insulin RIA kit, Millipore Corporation, St. Charles, Missouri, USA). Plasma insulin concentrations from lean rats were determined using a colorimetric ELISA kit (Ultra Sensitive Rat Insulin ELISA Kit, Crystal Chem INC, Downers Grove, IL, USA). The ELISA was used for lean rats to minimize blood sample volume (only 5 ll plasma required vs. $ 50 ll plasma for RIA). The RIA was used for the obese rats because their high lipid levels in plasma interfere with the ELISA but not the RIA measurement. Due to the hyperinsulinemia in the obese rats only 5 ll of plasma was required. For lean-rat plasma (with low lipid levels) the absolute insulin measurements are equivalent for the RIA and ELISA assays, according to an in-house comparison. Plasma NiAc concentrations were analyzed using LC-MS/MS with a hydrophilic interaction liquid chromatography (HILIC) approach, separated on a 50Â2.1 mm Biobasic AX column, with 5 lm particles (Thermo Hypersil-Keystone, Runcorn, Cheshire, UK) as previously described [16].

Model development
The exposure (PK) and biomarker (PD) models were developed sequentially because of the interaction between the model components; the kinetics of NiAc are assumed to be unaffected by insulin and FFA, whereas NiAc inhibits the release of both insulin and FFA. Furthermore, due to its antilipolytic effect, insulin affects FFA release. The interactions between the three models (NiAc, insulin, and FFA) are illustrated in Fig. 1b, and model interactions of previously published NiAc-FFA models [14,19,20,23] are illustrated in Fig. 1a for comparison. When a sub-model had been estimated, the random effects were fixed to the Empirical Bayes Estimates (EBE) and used as covariates in the subsequent sub-model.

Disease modeling and inter-study variability
The PK and PD were significantly different between lean (normal) and obese (diseased) rats and, consequently, these groups were modeled separately. Furthermore, the animal experiments were done under different conditions (separate time periods, anesthetized/conscious animals) which may have provoked different dynamic behaviors. To account for this, inter-study variability was included in the models in the form of fixed-study effects [25].

Notation conventions
To improve readability and enable the reader to differentiate between separate sub-model parameters, PD (insulin and FFA) model parameters are labeled with a subscript, indicating to which model they belong. For example, the turnover rate of FFA will be referred to as k inF and the turnover rate of insulin is k inI (i.e., F for FFA and I for insulin). Parameters that link NiAc, insulin, and FFA are Fig. 1 Schematic illustration of how the dependency between NiAc and FFA has been modeled in previous studies (a) and how the dependencies between NiAc, insulin, and FFA were modeled in this study (b). Solid lines represent fluxes while dashed lines represent control. NiAc inhibits the turnover of insulin (1). Insulin, in turn, has feedback mechanisms that inhibits its turnover (2) and stimulates its fractional turnover (3). Both NiAc (4) and insulin (5) inhibit the turnover of FFA. In this study, FFA has a single feedback mechanism which inhibits its turnover (6), while in previous studies, FFA was modeled using an additional feedback mechanism which stimulates its fractional turnover (7) labeled with both sub-model subscripts (e.g., potency of NiAc as an FFA inhibitor will be called IC 50NF , whereby the N is for NiAc and the F is for FFA).

NiAc exposure model
The pharmacokinetic properties of NiAc have been thoroughly characterized in previous studies [14, 16-19, 21-23, 26]. Ahlström et al. [16] introduced a two-compartment disposition model with parallel nonlinear (Michaelis-Menten) elimination for lean Sprague-Dawley rats, and a one-compartmental model with a single nonlinear elimination for obese Zucker rats (a schematic illustration of the PK models is given in Fig. 2).

Lean rats
In lean rats, the NiAc disposition is given by where C p ðtÞ is the observed NiAc concentration in the central plasma compartment and C t ðtÞ is the concentration in the peripheral tissue compartment (derivations of the initial conditions for these compartments are given in Appendix 2), and V p and V t are, respectively, the volumes of distribution of the plasma and tissue compartments. The parameters V max1 and K m1 are the maximal elimination rate and the Michaelis constant of the first pathway, and V max2 and K m2 are the maximal elimination rate and the Michaelis constant of the second pathway (low and high affinity pathway, respectively). Furthermore, Cl d is the inter-compartmental distribution, Synt the endogenous NiAc synthesis, and InputðtÞ is a time-dependent function determined by the route of administration according to where Inf. rate is the infusion rate, A sc ðtÞ is the amount of drug in the subcutaneous compartment, and k a is the absorption rate from the subcutaneous compartment to plasma. The rate of change of A sc ðtÞ is given by with initial condition A sc ð0Þ ¼ 0. Here, Pump rate represents the infusion rate from a subcutaneous mini-pump. The mini-pump was surgically implanted seven days before the final acute experiment. During this period, when the pump is not infusing, interstitial tissue fluid may diffuse into the tip of the catheter, diluting the NiAc dosing solution, whilst the solution is leaking into the tissue. Consequently, a concentration gradient may form, resulting in an apparently lower initial infusion rate compared to the pre-programmed setting (particularly pronounced in lean NiAc naïve rats, see Fig. 7a). To capture this, the pump infusion rate is modeled as where Inf. rate is the programmed infusion rate of the pump, d is a lumped diffusion parameter, and t 0 is the pump inactivation time (in this case 7 days). Here erf is the error function [27]. The derivation of the Pump rate is given in Appendix 1. Given NiAc's low molecular weight (123.11 g/mol), bioavailability from the subcutaneous compartment was assumed to be equal to unity.

Obese rats
For obese Zucker rats, the NiAc disposition is given by where C p ðtÞ is the NiAc concentration in the central plasma compartment, V c the volume of distribution, V max1 the maximal elimination rate, K m1 the Michaelis constant, and Synt the endogenous synthesis. The term InputðtÞ is the same as for the lean rats (the relations given in Eqs. 3, 4 and 5).

Between-subject and residual variability
The modeling was performed in an NLME framework to capture the between-subject variability seen in the exposure-time data. The parameters that varied within the population were k a , V max1 , and Synt, though Synt varied only in lean rats. These were assumed to be log-normally distributed in order to keep the parameter values positive. However, the five-day continuous infusion group of obese rats did not have exposure data. Consequently, these rats were assumed to behave like the estimated median individual. The residual variability was normally distributed and modeled using a proportional error model.

Estimated parameters
Because of sparse sampling, all parameter values could not be estimated from the data. By applying an a priori sensitivity analysis [28,29], we identified the parameters that had the greatest influence on the output. These were then estimated from the data and the remaining parameters were obtained from the literature [23]. The population parameters estimated from the data were k a , d, and V max1 .

Insulin turnover model
The primary aim of the insulin model was to establish smooth trajectories that would accurately describe the insulin-time courses under various provocations of NiAc, rather than describe all of the mechanistic aspects of insulin dynamics. To this end, the model structure was kept as simple as possible. The insulin model could subsequently be used to provide an input to the FFA model, enabling a quantitative analysis of the antilipolytic effects of insulin. Given this premise, a phenomenologically based modeling approach was applied. Under the assumption that NiAc perturbs insulin, the characteristics seen in the data were used to establish an insulin model with NiAc as input. The characteristic behavior of the data for acute and long-term NiAc provocations in lean and obese rats is illustrated in Fig. 3. Attributes seen include indirect action, tolerance (drug resistance), rebound, and complete adaptation (insulin levels returning to pre-treatment levels). Data with similar properties as those seen in the acute experiments (Fig. 3a, c) were modeled using turnover equations with moderator feedback control [14,30]. Furthermore, to capture the different long-term adaptive behaviors with (Fig. 3b), and without ( Fig. 3d) rebound, a 'NiAc action compartment' was included, as well as an integral feedback control. The insulin dynamics are given by with initial conditions where I(t) denotes the observed insulin level, and M 1I ðtÞ and M 2I ðtÞ the first and second moderator compartments, respectively. The parameters k inI and k outI are the turnover rate and fractional turnover rate of insulin, respectively, and k tolI is the fractional turnover rate of the moderators. The regulator compartment R I ðtÞ is given by where k inRI is the turnover rate, k outRI the fractional turnover rate, and I(t) the insulin concentration. The regulator compartment is initially at steady-state with By integrating Eq. 12, the dynamics of R I ðtÞ can be expressed as Hence, by construction, R I ðtÞ represents the output of an insulin-driven integral feedback controller [31] with I 0 as the set-point and k outRI as the integral gain parameter (k outRI will from here on be referred to as the integral gain parameter). The integral feedback controller will ensure that insulin levels return to the baseline I 0 , despite persistent external effects on insulin turnover and fractional turnover. The inhibitory NiAc function on insulin is given by where IC 50NI is the potency of NiAc on insulin and n the Hill coefficient of the inhibitory function. The term E NI ðN I ðtÞÞ represents the drug efficacy, which is fixed for lean rats and dependent on the concentration in a hypothetical NiAc action compartment, N I ðtÞ, for obese rats, according to where I maxNI is the initial efficacy of NiAc on insulin, N 50I the potency of the NiAc action compartment, S NI the longterm NiAc efficacy loss, and c the corresponding Hill coefficient of the efficacy relation. The dynamics of N I are in turn given by with N I ð0Þ ¼ C p ð0Þ. Here k NI is the turnover rate of the NiAc action concentration.
The NiAc action compartment is initially at steady-state with the plasma NiAc compartment C p . As infusions begin, and the plasma compartment concentration increases, N I ðtÞ increases until it reaches the steady-state NiAc concentration N ss ðtÞ ¼ C pss . With increasing levels in the NiAc action compartment, EðN I ðtÞÞ decreases to a minimum of 1 À S NI and, consequently, the efficacy of NiAc as an insulin inhibitor is down-regulated. In other words, the system has developed tolerance to the drug. The turnover rate k NI determines the rate at which tolerance develops. A schematic illustration of the insulin model is given in Fig. 4.

Between-subject, inter-study, and residual variability
Individual variations seen in the insulin data were incorporated in the model by allowing the parameters I 0 , k tolI , and IC 50NI to vary in the population. As in the PK model, these parameters were assumed to be log-normally distributed. The choice of these parameters was guided by an a priori sensitivity analysis. Moreover, the parameters I 0 and k tolI varied over study groups according to fixed-study effects on both the mean and individual parameter distributions [25]. In other words, for S, the number of groups, the parameter I 0 for an individual j was modeled as where Study k ¼ 1 if individual j is in group k and 0 otherwise. The residual variability was modeled using an additive model (with normally distributed errors). Mechanistic FFA model The model suggested in this study (schematically illustrated in Fig. 5) is founded on preceding approaches [14,19,20,23]; however, insulin has been included as the main endogenous regulator of FFA as insulin provides a homeostatic force on the system-thereby keeping FFA levels in the vicinity of its baseline concentration. Furthermore, the NiAc efficacy is dynamic in that it is decreasing during long-term infusions, which allows for complete systemic adaptation -a feature apparent in the data [32,33]. The characteristic behavior of the data, for acute and chronic NiAc provocations in lean and obese rats, is illustrated in Fig. 6. Attributes observed include indirect response, tolerance (drug resistance), rebound, and complete adaptation (FFA concentrations returning to pretreatment levels). The behavior observed in the acute experiment (Fig. 6a, c) has been described by turnover equations with moderator feedback (as described for the insulin system). The long-term behavior, and in particular the adaptations with, and without, rebound, is captured by dynamic NiAc efficacy and an insulin-controlled regulator. The FFA model is given by with initial condition Here, F(t) denotes the observed FFA level, k inF the turnover rate, and k outF the fractional turnover rate. The moderator compartment M F is given by with initial condition where the parameter k tolF represents the turnover rate of the moderator compartment. The moderator compartment provides a feedback mechanism for the turnover of FFA, that strives to dampen deviations from the baseline response. The regulator compartment R F ðtÞ, that links insulin dynamics to FFA release, is similar to that of the insulin model (Eq. 12) and is given by where k inRF is the turnover rate, k outRF the fractional turnover rate, and I(t) the insulin concentration. As for the insulin regulator, R F ðtÞ represents the output of an insulindriven integral controller with I 0 as the set-point and k outRF as the integral gain parameter. The contribution of this integral controller during acute and chronic NiAc treatments in lean and obese rats is illustrated in Fig. 9b. The inhibitory NiAc function on FFA (similar to that for the insulin model, Eq. 15), is given by Tolerance and rebound are captured by the moderator compartment M F , which acts on the turnover rate of FFA. The regulator compartment R acts on the turnover rate of FFA and the fractional turnover rate of R is affected by insulin Fig. 4 Mechanisms of insulin dynamics. The parameters k inI and k outI represent the turnover rate and fractional turnover rate, respectively. The turnover of insulin is inhibited by the NiAc action function H NI ðC p Þ. Tolerance and rebound is captured by the moderator compartments M 1I and M 2I , which act on the turnover rate and fractional turnover rate of insulin, respectively. The regulator R I , representing an integral feedback controller, acts on the turnover rate of insulin, in that it strives to maintain insulin baseline, I 0 , despite persistent external effects on the turnover where IC 50NF is the potency of NiAc as an inhibitor of FFA release and m is the Hill coefficient. The drug efficacy is dynamic and changes (down-regulates) during long-term infusions of NiAc. The efficacy is given by where I maxNF is the initial efficacy of NiAc on FFA, N 50F the potency of the NiAc action compartment, S NF the longterm NiAc efficacy loss, / the Hill coefficient, and N F ðtÞ the concentration in the NiAc action compartment. The dynamics of the NiAc action compartment are in turn described by with initial condition N F ð0Þ ¼ C p ð0Þ. Here, the parameter k NF is the turnover rate of the NiAc action state.

Between-subject, inter-study, and residual variability
Random effects were again selected using an a priori sensitivity analysis. The parameters that varied in the population were F 0 , k tolF , and IC 50NF (according to a lognormal distribution). Moreover, inter-study variability was incorporated in the model according to a fixed-study effect (as described for the insulin model). The parameters that varied between experimental groups were F 0 and k tolF . The residual variability was modeled using an additive model (with normally distributed errors).

Numerical analysis
The NLME modelling and simulations and the identifiability analysis were performed using Wolfram Mathematica (Wolfram Research, Inc., Mathematica, Version 10.3, Champaign, IL (2014).

Identifiability analysis
All population model structures analyzed in this study were proven to be structurally locally identifiable in a fixed effects setting (identifiability of the population model (fixed effects) implies identifiability of the statistical model (random effects) [34]). The identifiability analysis was performed using the Exact Arithmetic Rank (EAR) approach [35][36][37]-implemented in the Identifiabil-ityAnalysis Wolfram Mathematica package, developed by the Fraunhofer-Chalmers Centre. The EAR algorithm requires that all states and system parameters are rational functions of their arguments. This requirement is not fulfilled in the insulin and FFA systems (for example, the state space variable C p is raised to the power of n in Eq. 15). An illustrative example of how this requirement can be achieved is provided in the Appendix 3.

Selection of random effect parameters
An a priori sensitivity analysis was used to guide selection of the random parameters [28,29]. The system output sensitivity, with respect to the parameters, was analyzed and the parameters were ranked accordingly. The parameters with the highest sensitivity, given by the absolute value of the partial derivative of the system output with respect to a specific parameter evaluated at a given point in the parameter space, were considered random in the model.

Parameter estimation
Parameter estimates for the NLME models were computed by maximizing the first-order conditional estimation (FOCE) approximation of the population likelihood. This was done using a method developed and implemented in Mathematica 10 (Wolfram Research) at the Fraunhofer-Chalmers Research Centre for Industrial Mathematics (Gothenburg, Sweden) [38], which combines exact gradients of the FOCE likelihood based on the so-called sensitivity equations with the Boyden-Fletcher-Goldfarb-Shanno optimization algorithm [39]. Parameter standard errors were derived using the Hessian of the approximate population likelihood with respect to the parameters, evaluated at the point estimate. The Hessian was computed using finite differences of the exact gradients. From the steady-state relations in the insulin and FFA models, dependencies were derived which enabled the parameters k inI , k inF , k inRI , and k inRF to be expressed in terms of other model parameters (derivation given in appendix 2). Consequently, these parameters were redundant and could be replaced in the parameter estimation. Furthermore, some parameters were initially estimated to be very close to their physiological limit (e.g. I maxNI ¼ 0:9999 % 1 for obese rats) and were consequently fixed for numerical stability. Finally, to simplify the parameter estimation, some parameters were fixed (e.g., S NI ¼ 1 for obese rats). This is motivated by the complete systemic adaptation apparent in the long-term insulin-time data (obese rats), implying that S NI must be 1 (The fixed parameters are given in Table 2)). The long-term NiAc efficacy loss for lean rats was initially estimated to be % 0, whereby this part was omitted in the final model.

Results
The parameter estimation for the three sub-models (NiAc, insulin, and FFA) was performed sequentially, as described in the Model development section. The estimates and between-subject variabilities (expressed in CV%), both with corresponding relative standard errors (RSE%), for normal Sprague-Dawley rats and obese Zucker rats are given in Table 2. Weighted summaries [25] are presented for the parameters that varied between studies. The resulting models were qualitatively evaluated using visual predictive check (VPC) plots [40]; illustrating the data, the model predicted median individual, and 90% Monte Carlo prediction intervals generated from the models [40,41]. The VPC's are shown in Fig. 7 for lean Sprague-Dawley rats and in Fig. 8 for obese Zucker rats. The VPC's are generated from the PK, insulin, and FFA models for all provocations of NiAc.

Pharmacokinetic model
The pharmacokinetic system reached a steady-state concentration of about 1lM for all protocols both in lean and obese rats (first column in Fig. 7 and first column in Fig. 8). The steady-state was attained faster with intravenous than with subcutaneous administration. When infusions were terminated, the drug was cleared from the system within minutes and the NiAc concentration approached the endogenous level.
The absorption from the subcutaneous compartment had a half-lives of 0.16 and 0.13 h for lean and obese rats, respectively. At steady-state, the elimination of NiAc from the plasma compartment in lean rats was approximately three times faster for the high affinity pathway than the low affinity one. Moreover, the drug elimination rate from the plasma at steady-state was $ 20 and $ 25lmolkg À1 h À1 for lean and obese rats, respectively. The lumped diffusion coefficient was estimated to be 77 and 62 h À1=2 for lean and obese rats, respectively, implying that the NiAc dosing solution was diluted during the first $ 1.5 h.

Insulin model
The insulin concentration was suppressed below its baseline value at all provocations of NiAc. The suppression was more pronounced at an early stage of the infusions; at later stages, the insulin concentrations drifted back towards their baselines (cf. Figs. 7b-h or 8b-h). After the infusions were terminated, the insulin concentration rebounded before reaching its baseline value. Rebound was highest in the rats receiving the 12 h Off protocols (Figs. 7k, 8k) and was less pronounced in those receiving step-down protocols (Figs. 7n, t, 8n, t). In obese rats, the insulin concentrations returned to their baselines after long-term infusions of NiAc and did not rebound after the extended infusions were terminated (Fig. 8h).
The median baseline concentrations across groups were 0.233 and 3.51 nM for lean and obese rats, respectively. The estimates of the individual groups ranged between 0.151 and 0.264 nM (inter-study variability of 33% 1 ) for lean rats and 2.69-4.82 nM (inter-study variability of 33% 1 ) for obese rats. The inter-study variability was not correlated to the anesthetic condition of the rats. The median turnover half-lives of insulin, for the moderator, the integral controller, and the NiAc action level for lean and obese rats are given in Table 3. For lean rats, the efficacy of NiAc on insulin inhibition, I maxNI , was estimated to be 0.793; consequently, NiAc cannot completely inhibit insulin release. The established NiAc exposure was about 1lM which is approximately three times the NiAc potency related to inhibition of insulin (IC 50NI ¼0.338 lM). This implies that the inhibitory function was saturated at steadystate. The estimated Hill coefficient n indicates a steep NiAc concentration-insulin response relationship at steadystate. Furthermore, for obese rats, the efficacy was fixed to 1 (described in the parameter estimation section) and the corresponding potency was high since the IC 50NI was low in comparison to the NiAc steady-state exposure. However, since the estimated Hill coefficient was 0:84\1, indicating a gentle NiAc-concentration insulin-response relationship at steady state, the NiAc concentrations never reached levels high enough to saturate the inhibitory function. The estimated N 50I of the NiAc action compartment was lower than the steady-state NiAc concentration (0.897\1) and the Hill coefficient of the dynamic efficacy was estimated to be 18.9 (suggesting an all or non-response). This implies that the efficacy was completely down-regulated at the end of the long-term experiments in obese rats, implying no NiAc inhibition on insulin release.

FFA model
The FFA concentration was suppressed below its baseline value for all provocations of NiAc. Suppression was more pronounced initially during NiAc infusions. FFA concentrations then drifted back towards their baseline values (cf. Figs. 7c-f or 8c-i). After the infusions were terminated, the FFA concentrations rebounded before reaching their baseline values. The step-down protocols resulted in less rebound. The FFA concentrations returned to their baselines during extended exposure of NiAc in lean and obese rats (Figs. 7f, 8i, respectively). However, as the long-term exposure was terminated, rebound occured in lean, but not in obese rats (Figs. 7f, 8i). The median baseline concentrations across groups were 0.707 and 1.14 mM with corresponding ranges of 0.652-0.801 mM (inter-study variability of 11% 2 ) in lean and 0.789-1.22 mM (interstudy variability of 20% 2 ) obese rats, respectively. The inter-study variability was not correlated to the anesthetic condition of the rats. The median turnover half-lives of the FFA, the moderator, the integral controller, and the NiAc action level are given in Table 3.
Model simulations of the acute and chronic action of NiAc (H NF ðC p ðtÞÞ), the insulin-driven integral controller (R F ðtÞalso referred to as the insulin-controlled regulator), and the moderator feedback (M 0F =M F ðtÞ) on the FFA turnover are illustrated in Fig. 9 with the corresponding acute and chronic FFA responses. The simulated NiAc concentrations were set to be around 1lM (corresponding to the experimental exposures). The turnover rate of FFA was initially inhibited about 80% (lean) and 70% (obese) by the NiAc infusion (Fig. 9a). Upon the extended exposure to NiAc (120 h) the inhibitory action on the turnover rate was decreased by approximately 13% due to intrinsic tolerance mechanisms (Eq. 26). In obese rats, the NiAc action vanished completely (Fig. 9a). The insulin-driven controller (Eq. 14) provides a stimulatory action on the FFA turnover rate as insulin concentrations fall below the baseline (Fig.9b). The positive (stimulatory) action increases from 100% (at baseline) to about 200% after extended (120 h) exposure to NiAc in lean rats. The insulin action is totally abolished at equilibrium (120 h) in obese rats (Fig. 9b). The positive impact of the moderator is seen acutely both in lean and obese rats, whereas the moderator action has receded in obese rats at equilibrium (Fig. 9c). In lean rats, the combined inhibitory (NiAc) and stimulatory (insulin, moderator) action on the FFA turnover rate causes a rebound in FFA response when the NiAc infusion is stopped (Fig. 9d). This is not seen in obese rats since all of the NiAc, insulin, and moderator actions are back to baseline at equilibrium (120 h).
The potency of the NiAc action compartment was low for lean rats, since N 50F [ [ NðtÞ at equilibrium (120 h), and high for obese rats, since N 50F \NðtÞ at equilibrium (120 h). Consequently, the loss of efficacy was low in lean rats and high in obese ones. Furthermore, the estimate of Hill coefficient / ¼ 8:37 for obese rats suggests an all-ornone efficacy loss in obese rats.

Model predictions
The resulting model was used to predict 24 h FFA lowering (AUC 24 ) for a range of protocols at steady-state (i.e., after multiple dosing). Here, the lowering over a period T is given by b Fig. 7 Visual predictive checks for lean Sprague-Dawley rats. The first column shows the PK fit, the second column the insulin, and the third column the FFA. The rows represent the different protocols of NiAc (as described in the Experimental protocols section The protocol design consisted of a 0.25-12 h NiAc exposure period followed by a 0-12 h washout period. The NiAc infusions were designed to generate concentrations around the therapeutically relevant level ( $ 1lM) that was used in the experiments [16]. The predicted AUC 24 and proportional reduction, in comparison to baseline levels, on a median obese rat are given in Fig. 10. The model predicted an optimal dosing strategy of $ 2 h longer washout period than the exposure period and the maximal AUC reduction is 5.60 mM h. These results were consistent when AUC 24 was predicted for outliers with high/low potencies (IC 50NF ), baseline responses (F 0 ), and/or moderator turnover rates (k tolF ). Furthermore, the NiAc-concentration FFA-response relationship at steady-state, for obese rats, is illustrated in Fig. 11, with the corresponding NiAc-, insulin-, and moderator actions. The largest reductions in FFA exposure occurs within a 'window of opportunity', with NiAc concentrations between the IC 50NF and the N 50F . The predicted AUC 24 at steady-state, for the optimal NiAc exposure of 0.500lM, was 7.40 mM h.

Discussion
In this study, we applied a population modeling approach to a unique pre-clinical data set containing FFA-, insulin-, and NiAc-time courses obtained from acute and chronic provocations of NiAc in lean and obese rats. The aim was to identify a general model, from a macro perspective, to be used to predict optimal NiAc exposure profiles and to generate durable chronic dosing regimens.
Recent experimental data from long-term NiAc protocols have illustrated adaptation, with FFA concentration returning to its baseline value [15]. Those findings challenge previous models [18,23] for long-term dosing predictions.
The PK/PD model applied in this study was derived based on previous models [18,20,23], but we added, crucial mechanistic components in order to describe two different kinds of complete adaptation, one seen in lean rats and one in obese. In particular, insulin was included as an endogenous regulator of the turnover of FFA [9]. To this end, a separate insulin model was developed to describe the insulin dynamics for all provocations of NiAc.

Model characteristics
A key observation of the FFA-time course data in our study is that lean and obese rats both acquire complete adaptation, with FFA levels returning towards baseline at equilibrium (Fig. 6). The post-infusion FFA rebound in lean rats implies a NiAc-sensitive system. However, the rebound was less pronounced than that for the 1-12 h experiments. Consequently, the inhibiting effect of NiAc has been down-regulated during the long-term infusions. For obese rats, the rebound completely vanished, which implies a NiAc-insensitive system. This suggests two separate mechanisms that result in complete adaptation. To describe adaptation with drug effect (seen in lean rats), an insulin-driven integral feedback control was incorporated into the model [39]. Control system techniques have previously been applied in glucose-insulin models in form of PID (proportional-integral-derivative) controllers [42]. The set-point of the insulin-driven integral controller is the insulin baseline and, hence, its action reflects the deviations from baseline insulin concentrations. For the FFA model, the controller represents the antilipolytic effect of insulin [24] via a regulator; as the insulin level increases, the antilipolysis will be more pronounced since the elimination rate of R F will increase and, consequently, R F 's stimulation of FFA release will be lowered (and vice versa). The traditional I max -equation was modified with a dynamic efficacy function in order to capture the phenomena of NiAc resistance. The dynamic efficacy may represent the down-regulation of the PDE3B gene expression [7].
The impact of the dynamic efficacy and the integral controller ( Fig. 9) show that both adaptive actions push FFA concentrations back towards baseline (at equilibrium) b Fig. 8 Visual predictive checks for obese Zucker rats. The first column shows the PK fit, the second column the insulin, and the third column the FFA. The rows represent the different protocols of NiAc (as described in the Experimental protocols section). The dots represent the data, with colors indicating separate individuals, the black line the estimated median individual, and the grey area the 90% population prediction interval. No exposure data were available from the Cont. protocol (g). despite ongoing NiAc exposure. The insulin-driven controller has less effect in obese as compared to lean rats in spite of 10-fold higher insulin concentrations in the former group. This reflects the insulin resistance of the obese rats [7].

Model evaluation
The x-axis represents the infusion time and the y-axis represents the washout time. The model predicts an optimal infusion regimen of $ 2-h longer washout period than the infusion. The maximal AUC reduction is 5.60 mM h turnover rate is operating on a significantly shorter timescale than the feedback and adaptive mechanisms. This resulted in low precision of k outF (RSE% of 140 and 120, respectively, in lean and obese rats). To achieve higher precision a denser sampling of the FFA time-course is needed during the initial infusion phases. On the other hand, the drug resistance in lean rats was a slower process (half-life of the NiAc action of more than 100 h). This resulted in low precision of the parameters linked to efficacy loss (RSE% of 160 for k NF and 190 for S NF ). Some estimated median responses of the insulin and FFA models were under-or over-predicted (cf. Fig. 7l or Fig. 7k). This is most likely due to the low number of individuals per trial (5-10), implying that every 4-9th population median will be estimated below or above the trends seen in the individual data [43]. Furthermore, some predicted 90% population spans also under-or over-predicted the response for the insulin and FFA model (cf. Fig. 7l or Fig. 7k). This is most likely due to correlations of the between-subject parameter variabilities, which were not captured because diagonal covariance matrices were used in the FFA model [44]. When sampling from the resulting distributions to generate the VPC's, non-feasible parameter combinations may occur which render a skewed population [44]. Diagonal covariance matrices were chosen in order to simplify the parameter estimation.
The turnover processes in the insulin and FFA models operate over completely different time scales (Table 3). The general trend in both models is that the moderator and integral controller processes are slower in obese rats. For example, the insulin-driven integral controller has a halflife that is 50 times longer in obese as compared to lean rats. These control processes are tightly controlled in lean animals, and are probably elongated by the tremendous hyperinsulinemia (nearly 10-fold greater pathological concentrations), and corresponding insulin resistance, in obese rats. In the FFA model, turnover of FFA is more than 100 times faster than moderator feedback and integral control. Furthermore, turnover of the NiAc action compartment had a half-life of more than 100 h for lean rats, thus spanning the entire duration of the experiment. The corresponding half-life in obese rats was 18 h. Hence, obese rats reach complete adaptation much faster than lean rats.
Due to the nature of the FFA dynamics, with tolerance development and rebound post-infusions, constructing an optimal dosing protocol is challenging. By selecting an inappropriate dosing regimen, the NiAc provocation can yield an increased FFA exposure in comparison to controls (Fig. 10-the negative AUC area). Given a NiAc exposure of $ 1lM (the exposure that was used in the experiments), there is an optimal strategy whereby washout periods are 2 h longer than infusion periods; this is illustrated by Fig. 10. These strategies were consistent in that they were tested on the median individual and on 90% quantiles (i.e., individuals that had the 5 and 95% quantiles of the parameters that varied in the population: IC 50NF , F 0 , and k tolF ). However, rather surprisingly, higher reductions in AUC were attained at constant NiAc exposure (without washout periods) at lower concentrations, where the maximal reduction was 7.40 mM h (in comparison to the maximal reduction of 5.60 mM attained from the exposure/ washout protocols).

Conclusions
Our study presents a novel NiAc-insulin-FFA model that could accurately describe the concentration-response relations seen during acute and chronic NiAc treatment in lean Sprague-Dawley and obese Zucker rats. In particular the model could describe two different types of adaptive changes. This was done by applying an insulin-driven integral controller and a dynamic efficacy in the traditional Fig. 11 Predicted steady-state concentration-response (lefthand y-axis) and concentrationaction (right-hand y-axis) relationships for obese rats. The black line represents the FFA response, the blue line the inhibitory action of NiAc on FFA turnover, the red line the insulin action on FFA turnover, and the purple line the moderator action on the FFA turnover. The IC 50F and the N 50F concentrations are given by the black vertical lines (Color figure online) I max model. The dynamics of these methods make them suitable for a range of tolerance scenarios. Finally, the model was used to simulate infusion-washout regimens, with a NiAc exposure of 1lM, in order to estimate the 24 h lowering of FFA. Given the targeted exposure, the importance of incorporating washout periods in-between infusions was illustrated. However, the predicted concentration-response relationship suggests that higher reductions in AUC could be obtained by using lower NiAc concentrations. These findings should be experimentally verified in future studies.
NiAc is close to zero, i.e. C p ð0Þ % 0, the inhibitory NiAc function on insulin is H NI ðC p ð0ÞÞ % 1 and we obtain dI dt ¼ k inI À k outI Á I 0 ¼ 0 () ð38Þ Thus, we can eliminate k inI from the estimation procedure and obtain an estimate for this parameter from Eq. 39. Assuming that the system in Eq. 20 is at an initial steady-state with dF dt ¼ 0 and that the initial plasma concentration of NiAc is close to zero, i.e. C p ð0Þ % 0, then the inhibitory NiAc function on insulin is H NF ðC p ð0ÞÞ % 1 and we obtain dF dt Thus, we can eliminate k inF from the estimation procedure and obtain an estimate for this parameter from Eq. 41. [45]