## Abstract

The oral glucose tolerance test and meal tolerance test are common clinical tests of the glucose-insulin system. Several mathematical models have been suggested as means to extract information about β-cell function from data from oral tolerance tests. Any such model needs to be fairly simple but should at the same time be linked to the underlying biology of the insulin-secreting β-cells. The scope of the present work is to present a way to make such a connection using a recent model describing intracellular mechanisms. We show how the three main components of oral minimal secretion models, derivative control, proportional control, and delay, are related to subcellular events, thus providing mechanistic underpinning of the assumptions of the minimal models.

- secretory granules
- β-cell function
- C-peptide secretion
- threshold distribution
- mathematical models

pancreatic β-cells are responsible for secretion of insulin in response to elevated plasma glucose levels, and the malfunctioning of these cells is a major contributor to the development of diabetes. The oral glucose tolerance test (OGTT) and the meal tolerance test are common clinical tests, where plasma glucose, insulin, and C-peptide concentrations are measured after ingestion of the sugar or nutrients. Several mathematical models have been suggested as means to extract information about β-cell function from these kinds of data (3, 8, 16, 21, 32). Any such model needs to be fairly simple to ensure identification of its parameters, but at the same time one wishes to link the minimal mathematical description to the underlying biology of the insulin-secreting β-cells. The scope of the present work is to present a way to make such a connection using a recent model (25) describing intracellular events.

Oral minimal models differ in their formulation. The plasma glucose concentration is the main controller of secretion, either instantaneously (16, 21) or with a delay (3, 8, 32). This proportional control has been denoted as static secretion (3, 32). In the case of no delay, a correcting factor is needed to modify the degree of glucose control (21). Some models include additional control by the time derivative of the glucose concentration (3, 21, 32), whereas other authors do not include such a dynamic secretion term (8, 16). It has been shown that derivative control via the dynamic secretion term is necessary for fitting data from in vivo studies (3, 4, 32).

Derivative control can arise from the threshold hypothesis, as explained by Grodsky (14) and in greater detail by Licko (19). Grodsky (14) proposed that insulin is located in “packets,” plausibly the insulin-containing granules but possibly entire β-cells. Some of the insulin is stored in a reserve pool, whereas other insulin packets are located in a labile pool, ready for release in response to glucose. The labile pool is responsible for the first phase of insulin secretion following a step in the glucose concentration, whereas the reserve pool is responsible for creating a sustained second phase (14). This basic distinction has been, at least partly, confirmed when the packets are identified with granules (9, 23). Moreover, Grodsky (14) assumed that the labile pool is heterogeneous in the sense that the packets in the pool have different thresholds with respect to glucose, beyond which they release their content. This assumption was necessary for explaining the so-called staircase experiment, where the glucose concentration was stepped up and each step gave rise to a peak of insulin.

There is no supporting evidence for granules having different thresholds (22), but Grodsky (14) mentioned that cells apparently have different thresholds based on electrophysiological measurements. Later, Jonkers and Henquin (17) showed that the number of active β-cells is a sigmoidal function of the glucose concentration, as assumed by Grodsky (14) for the threshold distribution. Recently, we showed how to unify the threshold distribution for cells with the pool description for granules (25), thus providing an updated version of Grodsky's model (14), which takes into account more of the recent knowledge of β-cell biology.

The present work will show how derivative control arises in our recent mechanistic model (25) and will investigate how the three main components of oral minimal models, derivative control, proportional control, and delay, are related to the subcellular events. This will provide mechanistic underpinning of the assumptions of the minimal models. We will investigate the relative contribution of dynamic and static insulin secretion during a meal test using mathematical analysis and computer simulations.

## METHODS

Toffolo et al. (32) and Breda et al. (3) proposed the following minimal model of insulin secretion:
_{b} is the basal rate of secretion. The static rate of secretion, SR_{s} = Y, due to mobilization with a delay with time constant *T*, is controlled by plasma glucose concentration G as

The dynamic secretion rate is controlled by the glucose rate of change dG/d*t* as
*k* are measurements of static and dynamic β-cell sensitivity, respectively. According to *Eq. 3*, the dynamic control is maximum when glucose increases above its basal value G_{b}; then it decreases linearly with glucose concentration and vanishes when G exceeds the threshold level G_{t}. In combination with a model of C-peptide kinetics (see below), the minimal secretion model allows estimation of its parameters from OGTT and meal test C-peptide data and yields indices of β-cell function: the delay *T*, the basal sensitivity index SR_{b}/G_{b}, the static sensitivity index β, and the dynamic sensitivity index related to *k* normalized appropriately to the range of glucose concentrations (3, 4, 32). Because of identification problems, it is often assumed that *k*(G) = *k*_{d} is glucose independent, which corresponds formally to setting G_{t} equal to infinity in *Eq. 3*.

The model by Hovorka et al. (16) is given by SR = Y_{∞} (G), whereas the model by Cretti et al. (8) has SR = Y, as in *Eq. 2*; i.e., it is identical to *Eq. 1*, with the dynamic term equal to zero. The model by Mari et al. (21) is given by SR = SR_{d} + P(*t*) f(G), with SR_{d} as in *Eq. 3* and f(G) a nonlinear dose-response function that is corrected by a time-varying function, P(*t*).

Our subcellular model (25) includes mobilization of secretory granules from a reserve pool to the cell periphery, where they attach to the plasma membrane (docking). The granules can mature further (priming), thus entering the “readily releasable pool” (RRP). Calcium influx triggers membrane fusion and subsequent insulin release. We included the possibility of so-called kiss-and-run exocytosis, where the fusion pore reseals before the granule cargo is released. The glucose-dependent increase in the number of cells showing a calcium signal (17) was included by distinguishing between readily releasable granules in silent and active cells. Therefore, the RRP is heterogeneous in the sense that only granules residing in cells with a threshold for calcium activity below the ambient glucose concentration are allowed to fuse. An overview of the model is given in Fig. 1.

Mobilization is assumed to depend on glucose but with a delay, τ, similar to Grodsky (14):

The docked pool develops according to the mass-balance equation

The RRP is described by a time-varying density function h(g,*t*), indicating the amount of insulin in the RRP in β-cells with a threshold between g and g + dg. Granules are primed with rate *p ^{+}* and are assumed to lose the capacity of rapid exocytosis with rate

*p*

^{−}. Moreover, if the granule is in a triggered β-cell, i.e., a cell with a threshold below the glucose concentration, it will fuse with rate

*f*. This leads to the equation

^{+}Here, θ(G − g) is the Heaviside unit step function, which is 1 for G > g and zero otherwise, indicating that fusion occurs only when the threshold is reached. The priming flux *p ^{+}*D distributes among cells according to the fraction of cells with threshold g described by the time-independent function φ(g). Thus, priming is assumed to occur with the same rate in all cells, but we take into account the fraction of cells with the corresponding threshold.

The secretion rate can be expressed as
_{b} is basal secretion, *m* is the rate constant of release, and F is the size of the fused pool, which follows
*f ^{+}* is the rate constant of fusion and

*κ*is the kiss-and-run rate. The integral represents the amount of insulin in the RRP in cells with a threshold below G. For further details of the model, we refer to the original article (25).

For the simulation studies, parameters of φ were chosen to get a decent fit of Φ(G) = ∫_{0}^{G} φ(g) dg to the fraction of cell clusters showing calcium activity at different glucose concentrations taken from Jonkers and Henquin (Fig. 2 in Ref. 17), since β-cell clusters are more physiological than isolated cells due to cell coupling. Hence, in the model description above, “active cells” should be interpreted as “cells in active clusters,” and φ describes accordingly the fraction of cells located in active clusters. This question is treated further in the conclusions. Parameters of M_{∞}(G) were chosen to reproduce data on the glucose dependence of mobilization (Fig. 10 in Ref. 14). Both Φ and M_{∞} were assumed to be Hill functions of G. However, to have a more direct control of basal secretion, we assumed that Φ vanished below G = 5 mM, which is approximately the mean basal glucose concentration in the subjects used for parameter adjustments (see below). This was obtained by shifting the function 5 mM to the right by setting

We use the following expression for steady-state mobilization:

Other parameters were adjusted to give a satisfactory simulation of mean C-peptide concentration in healthy individuals following a meal (2). For this purpose, it was necessary to complement the above model with the classical two-compartment model of C-peptide kinetics (11) given by

Here, *compartment 1*, accessible to measurement, represents plasma and rapidly equilibrating tissues, whereas *compartment 2* represents tissues in slow exchange with plasma. The factor 43 changes units from micrograms of insulin per pancreas per minute to picomolar concentrations of C-peptide per minute, assuming a volume of *compartment 1* of 4 liters (32) and equimolar secretion of insulin and C-peptide (36). Parameters of C-peptide kinetics were set to mean values from (32) *k*_{01} = 0.06, *k*_{21} = 0.0559, and *k*_{12} = 0.0492. We used the mean glucose profile from the data set (2) as input to the combined secretion-C-peptide model with SR_{b} = *k*_{01} CP_{1}(0)/43 and did a prerun to start from a steady state. All parameters of the secretion model are given in Table 1. Equations were solved with a fifth-order Dormand-Prince integrator, as implemented in XPPAUT (12).

## RESULTS

#### Theory.

Due to the threshold on cells, the mechanistic model (25) possesses derivative control, which can be seen as follows. Rewriting *Eq. 8*, the fused pool develops according to
*t*) = ∫_{0}^{G} h(g, *t*) dg follows (using *Eq. 6*)

Assuming that D, G, and dG/d*t* change much more slowly than the time scale of the rates *f ^{+}*,

*κ*, and

*m*, it is reasonable to expect quasi-equilibrium in

*Eqs. 9*and

*10*, and we obtain from

*Eq. 7*that the secretion rate is given by

Thus, SR is composed of the three main ingredients of the OGTT models: *1*) a static term, which includes *2*) a delay, τ, due to the delayed refilling of the docked pool, D, and *3*) a dynamic term proportional to dG/d*t*. Compare with *Eqs. 1–3* and Refs. 3, 21, and 32.

We note that for dG/d*t* < 0, cells with threshold G are already active and will therefore have near-empty RRP, so h(G,*t*) will be small. This can also be seen from an approximation to h(G), obtained by assuming quasi-equilibrium in *Eq. 6*. Then,
*f ^{+}* is much larger than

*p*and

^{+}*p*

^{−}, we observe that a significant derivative control is obtained only when dG/d

*t*< 0, in agreement with

*Eq. 3*. We observe also that the secretion rate in

*Eq. 11*is insensitive to changes in the fusion rate

*f*, as long as fusion is faster than the other processes, since the fraction in square brackets in C is approximately near 1, because

^{+}*f*is much larger than

^{+}*p*

^{−}. This fact was observed numerically by Pedersen et al.(25).

To obtain a model with no derivative control (8, 16), we would need to have a very small RRP, and thus low values for h, e.g., due to the forward priming rate *p ^{+}* being much smaller than the backward priming rate

*p*

^{−}, as seen from

*Eq. 12*. This might be the situation in diabetics, which show a reduced first phase of secretion (10, 13).

Direct control by G with no delay would be obtained if the docked pool D (or a similar pool refilling the RRP) does not change its size notably during the experiment. The time-varying correction function (21), P(*t*), would likely correspond to a time-varying priming rate, which is not controlled directly by glucose.

#### Simulations.

To test the validity of our simplifying assumptions, we compare the secretion rate in response to the mean plasma glucose concentration G following a meal (2). We calculate dG/d*t* numerically and simulate the secretion rate from the full mechanistic model as well as from the approximating expression given by *Eq. 11* (Fig. 2). We see that the approximation to the full model (Fig. 2, dotted curve) does an excellent job in predicting the secretion rate from the full model (Fig. 2, solid curve). This allows us to separate the secretion rate in its dynamic (Fig. 2, dashed curve), static (Fig. 2, dashed-dotted curve), and basal (data not shown) parts.

The dynamic term given in *Eq. 3* of the OGTT secretion model depends on G through the function *k*, whereas derivative control in the approximation (*Eq. 11*) to the model is modified mainly through h, which is nearly proportional to the function φ, as seen in *Eq. 12*. Figure 3 shows φ and *k*, scaled to have area under the curve equal to 1 as φ, assuming G_{t} = 12 mM for the value where *k* intersects the *x*-axis. In in vivo studies, G_{t} can often not be estimated and has in many studies been set to infinity, meaning no glucose dependence of *k*. In the cases where it was possible to estimate G_{t} (3), values of 9.8 and 11.9 mM were found, whereas an up-and-down glucose infusion test gave G_{t} values of 11.2, 12.4, and 28.4 mM when the parameter was identifiable (32). Moreover, Toschi et al. (33) observed plateauing of the first-phase insulin secretion response to steps in plasma glucose concentration in humans at ∼12 mM, and Henquin et al. (15) similarly found that insulin secretion from human islets reached maximum rates between 10 and 15 mM glucose. These results indicate that the RRP in humans is emptied by glucose concentrations of ∼12 mM, as is the case for rats (14). This was also suggested from an intravenous glucose tolerance test study, where amounts of either 300 or 500 mg/kg glucose were injected, but the total releasable pool (∼RRP) was estimated to be similar in the two scenarios, suggesting that glucose concentrations of <15 mM are sufficient to empty the RRP (31). We see from Fig. 3 that, under this assumption, φ and *k* and their integrated functions are indeed similar. We suggest fixing G_{t} = 12 mM, or a similar value, in studies where its value cannot be estimated reliably rather than assuming no glucose dependence of *k*.

## CONCLUSIONS

Minimal models are necessary as “models to measure” for estimating parameters from in vivo data. Detailed mechanistic models, on the other hand, are “models to simulate” and provide insight into the underlying biological mechanisms. In the present work, we have shown how to connect these two types of models by using reasonable simplifications of the mechanistic model (25) to obtain an approximation to the model that is readily compared with oral minimal models (3, 8, 16, 21, 32). Such a comparison should allow the identification of candidate cellular mechanisms that underlie changes in indices of β-cell function. For example, we found that the threshold distribution φ corresponds approximately to the function *k* of the minimal oral model (3, 32) (*Eqs. 11* and *12* and Fig. 3). The importance of making *k* glucose independent by formally setting the parameter G_{t} equal to infinity in cases where G_{t} cannot be estimated should be investigated by comparing with results obtained with G_{t} fixed at 12 mM, which appears more correct from a mechanistic perspective, as discussed in the results.

Importantly, our analysis revealed that the mechanistic model contains the three main components of oral C-peptide minimal model (3, 32): *1*) proportional control, i.e., a static secretion term controlled by the glucose concentration; *2*) a delay; and *3*) derivative control, i.e., a dynamic secretion term controlled by the rate of change of the glucose concentration. The latter arises from the cellular activation thresholds (17) (see also Ref. 19).

The use of models with no dynamic term might be justified in subjects with a small RRP, as is often seen in diabetic patients. However, applying a model with a dynamic term to such patients will directly reveal the small RRP as a very low estimate of the size of derivative control expressed by the parameter *k*_{d} without assuming this a priori.

Assuming for a moment complete confidence in our mechanistic model, we conclude that a static term with no delay corresponds to a situation with a large pool refilling the RRP directly with a very short or no delay. In addition, since a static secretion term with no delay must be corrected by the time-varying correction function P(*t*) (21), such a glucose-controlled, immediate refilling should be complemented by, e.g., a time-varying priming rate that contains a glucose-insensitive part. It has been shown that this time-varying effect is not due to the incretin effect (6), but there exist a series of possible second messengers that might be responsible for such a modulation of the priming rate (29). However, restructuring of the cortical actin web permits mobilization of granules from the interior of the cell and is believed to be necessary for late insulin secretion (24, 35). Rac1, a member of Rho family GTPases regulating cytoskeletal organization, has been suggested to be involved in this process, and this molecule is activated 10–15 min after a rise in glucose stimulus (18, 34). Hence, a delay of this order might be inherent to mobilization and subsequent refilling of the RRP.

A fundamental assumption of our model is the fact that the β-cells activate at different glucose levels (17). In situ, the β-cells are coupled by gap junctions in the islets of Langerhans, and this coupling has been reported to reduce cell heterogeneity and promote interislet synchronization (27, 30). However, asynchronous Ca^{2+} oscillations have been seen in islets from rats (20) and humans (5). Moreover, by imposing a glucose gradient across single islets, Rocheleau et al. (28) showed that active cells were not able to activate otherwise silent cells; in other words, a cell was active if and only if its expected single-cell threshold was reached. In addition, even if coupled cells should show similar calcium levels, it is not a given that the cells will secrete insulin with the same rate (17), for example, due to different levels of cAMP and ATP, which are important for sensitizing the exocytotic machinery to the calcium signal (29). Such heterogeneity could provide the cell-to-cell variations used in the model.

More detailed models of granule dynamics than the one presented by Pedersen et al. (25) and analyzed here exist in the literature (1, 7, 26). It might be useful to do analyses of these models along with the ideas presented here to understand how processes such as calcium handling, metabolites, and granules with different properties relate to more macroscopical parameters.

## GRANTS

This work was supported by the European Union through the Network of Excellence BioSim (LSHB-CT-2004-005137). M. G. Pedersen received support from the Lundbeck Foundation.

## DISCLOSURES

No conflicts of interest are declared by the author(s).

- Copyright © 2010 the American Physiological Society