INTRODUCTION

TOP ABSTRACT INTRODUCTION MATERIALS AND METHODS RESULTS DISCUSSION REFERENCES

THE MINIMAL MODEL of glucose disappearance (5) is widely used to estimate metabolic indexes of glucose effectiveness (S_G) and insulin sensitivity (S_I) in humans in both normal and physiopathological conditions. Since 1994, ~280 peer-reviewed scientific reports from the intermediary metabolism community have used the minimal model approach or have discussed related methodology, demonstrating substantial interest in the method; 65 of these reports have appeared in the past two years. Use of the minimal model for S_I determination is often preferred over the more invasive and labor-intensive glucose clamp technique, especially in large clinical trials and epidemiological studies (26). Minimal model-based insulin sensitivity (S_I) functions very well as a surrogate measure of cardiovascular disease risk in normal subjects (18) and in type 2 diabetics (17). Currently, several studies are using it to search for diabetes-relevant genes in the human genome (15, 16). Caumo et al. (7) have criticized this method because of its likely undermodeling of plasma glucose kinetics. However, its rising popularity among several circles warrants some effort to improve not only the data collection portion of the experiment but also the related procedures of data analysis, because physiologically significant conclusions will rest on both of these aspects.

Because of its simplicity, the standard intravenous glucose tolerance test (IVGTT) protocol (sIVGTT), based on the injection of a glucose bolus at time 0 and subsequent sampling for 3 or 4 h, is by far the most popular experimental protocol. However, the precision of S_G and S_I estimates is not always satisfactory, especially with reduced sampling schedules. The insulin-modified (mIVGTT) protocol (13, 37), which involves an infusion of insulin between 20 and 25 min after the glucose bolus, has allowed considerable improvement in the precision of the individual estimates of S_G and S_I. However, the modified protocol becomes substantially more complex than the standard procedure and may lead to hypoglycemic episodes in normal volunteers; also, it does not carry over easily to pediatric populations. In addition, the presence of exogenous insulin makes the determination of concomitant -cell activity difficult, and Pacini et al. (21) have recently suggested that, in studies comparing groups, use of the sIVGTT is more desirable. Improvement of the precision of the estimates in the sIVGTT is therefore paramount.

In the standard approach to minimal modeling, each subject is analyzed individually, and the additional information that all subjects belong to a homogeneous population is never exploited, even when available. With a single exception (10), minimal model analysis of large groups has so far neglected the fact that each individual belongs to a population of subjects who share certain quantitative traits. Nevertheless, the fact that all subjects in similarly performed studies have similar characteristics, like body weight, age, medications, etc., carries a certain amount of information. Thus the techniques of population analysis could bring a significant contribution to the minimal model method. Population analysis is the methodology used to quantify between-subject (also called intersubject in the literature) variability relative to a given population model. Its use is widespread in pharmacokinetic/pharmacodynamic studies, because it helps solve problems crucial in clinical studies such as sparse sampling and protocol deviation. Similar issues will become increasingly important in physiology as metabolic studies like the IVGTT move out of the investigative stage (small number of subjects) into the clinical arena (large number of subjects). The statistics literature has especially studied population analysis, developing the related methodologies of analysis of repeated measurement data, nonlinear mixed-effects modeling, and analysis of longitudinal data (9, 28).

In this work, we will describe the results of population analysis applied to glucose minimal modeling, and we will use the population results to provide Bayesian priors for the individual estimates of the metabolic indexes. The database consists of 16 sIVGTTs performed in normal subjects; we analyze these data with both a frequent (FSS) and a reduced sampling schedule (RSS). The use of an RSS is of interest in experimental design, because it potentially reduces experiment costs, laboratory worker biohazard, and patient discomfort. We will show that the use of a population analysis iterative algorithm, the iterative two-stage (ITS) approach, allows significant improvement in the precision of individual parameter estimates for each subject, even with RSS.

	MATERIALS AND METHODS

TOP ABSTRACT INTRODUCTION MATERIALS AND METHODS RESULTS DISCUSSION REFERENCES

Data. The data consist of a series of sIVGTTs performed in a population of 16 young adults. Experiments were performed at the Washington University School of Medicine, St. Louis, MO, and the Department of Metabolic Diseases, University of Padova, Padua, Italy. We published these data previously in Vicini et al. (35), to which we refer the reader for experiment details. Briefly, a glucose bolus (between 300 and 330 mg/kg) was administered at time 0. Samples were drawn with an FSS of 30 samples, from which an optimal RSS of only 14 samples was extracted (0, 2, 3, 4, 5, 8, 10, 12, 14, 16, 18, 20, 24, 28, 32, 40, 45, 50, 60, 70, 80, 90, 100, 110, 120, 140, 160, 180, 210, 240 min; the RSS in italics). The RSS is based on optimal glucose sampling schedule studies (8, 25). Because different laboratories gathered the data, the sampling schedule differed slightly among subjects; the ones we report are only indicative.

Minimal model of glucose disappearance. We describe the minimal model of glucose disappearance (5) by

(1)

where D is the glucose dose (mg/kg body wt), Q(t) (mg/kg) is glucose mass in plasma, G(t) (mg/dl) is glucose concentration, I(t) (µU/ml) is insulin concentration, G_b (= Q_b/V) and I_b are their basal values, and X(t) is insulin action (min¹). Insulin concentration acts as a known input (without error) in the second equation, and we estimate model parameters by fitting the model response, G(t), to glucose concentration data. The model has four uniquely identifiable parameters: S_G (min¹), S_I (min¹ · µU¹ · ml), p₂ (min¹), the insulin action parameter, and V (dl kg¹), the glucose distribution volume. S_G and S_I are the minimal model indexes of glucose effectiveness and insulin sensitivity, respectively, and reflect the effect of glucose and insulin on both glucose disposal and production. In particular, S_G measures the ability of glucose per se, at basal insulin, to stimulate glucose disposal and to inhibit glucose production. Similarly, S_I measures the ability of insulin to enhance the glucose per se stimulation of glucose disposal and the glucose per se inhibition of glucose production. The minimal model parameter vector is thus p = [S_G, S_I, p₂, V].

Standard two-stage approach. The typical method of estimating the mean (first-order moment) and the variance (second-order moment) of the population distribution of S_G and S_I consists of estimating S_G and S_I separately in each subject via nonlinear regression and then calculating the sample mean (µ) and covariance (²) of all the S_G and S_I estimates. The results are grouped according to the population (e.g., a population of normal controls and a population of diabetic patients should be analyzed separately). This method is often called the standard two-stage (STS) approach. The values of µ and ² represent the mean and variance (or rather, the first- and second-order moments, if we cannot assume normality) of the population distribution. For the STS method, we used weighted nonlinear least squares, as implemented in the kinetic analysis software SAAM II (1), to determine each subject's parameters. We minimized the weighted residuals sum of squares with respect to the vector of model parameters for a subject j, p_j

(2)

where N_j is the number of data points available for the j^th subject, t_i,j and G_i,j^OBS are the i^th time point and data point, respectively, of the j^th subject, ²_i,j is the variance of the measurement error of the i^th data point, and G(p_j,t_i,j) is the minimal model prediction of glucose concentration for a given p_j.

(3)

View larger version (23K): [in this window] [in a new window]   Fig. 1.   Plot of individual data points () and individual minimal model fits () for the frequent sampling schedule (FSS) and standard two-stage (STS) analysis.

(4)

(5)

ITS approach. The ITS is a parametric, iterative population analysis method based on the concepts of population prior knowledge and maximum a posteriori (MAP) probability empirical Bayes estimation. Steimer et al. (34) proposed it as a computationally attractive alternative to the nonlinear mixed-effects model approach (9). The steps of the ITS follow.

(6)

(7)

(8)

(9)

Statistical methods. Because of both the nonnormality of the sample and the heterogeneity of variances between the FSS and RSS groups (see Ref. 30, p. 296), we compared the S_I and S_G values obtained from different estimators and different sampling schedules by means of the Wilcoxon signed-ranks test for matched pairs (Ref. 30, p. 128). We set the significance level at P = 0.05, and we used the nonparametric Spearman rank correlation to measure agreement.

	RESULTS

TOP ABSTRACT INTRODUCTION MATERIALS AND METHODS RESULTS DISCUSSION REFERENCES

We report here the STS and the ITS results for both sampling schedules (FSS and RSS). Table 1 (FSS by use of STS and ITS) and Table 2 (RSS by use of STS and ITS) show individual results for all of the subjects for all of the minimal model parameters, whereas Fig. 2 displays a visual summary of our findings on the percent precisions of parameter estimates.

View larger version (20K): [in this window] [in a new window]   Fig. 2.   Precision of glucose effectiveness (SG) and insulin sensitivity (SI) estimates with different sampling schedules [FSS or reduced sampling schedule (RSS)] and population analysis approaches [STS and iterative two stage (ITS)]. Mean values are shown.

We first applied the STS method to both FSS and RSS data. We detected a significant difference (P < 0.05) between S_I determined with STS analysis of the FSS compared with the RSS data. To investigate the effect of the number of samples on population analysis, we used linear regression to compare the individual estimates across sampling schedules. From Fig. 3 we can infer that S_I is fairly stable (except when its value is very low, as in subject 11), whereas S_G is slightly less stable. With the FSS data, the average precision of the estimates, calculated as the mean of all individual precisions, was 51, 33, 73, and 7%, respectively. However, the average precision of the estimates worsened considerably with the RSS data (79, 91, 227, and 11%, respectively). Therefore, the RSS considerably hampers individual parameter reliability.

View larger version (15K): [in this window] [in a new window]   Fig. 3.   Impact of the RSS on the SG and SI estimates from the standard intravenous glucose tolerance test (sIVGTT). Estimates are compared with the corresponding estimates obtained under FSS. R2 values reported are Spearman rank correlation results.

We then analyzed FSS and RSS data with the ITS method. Results for FSS data were comparable to the STS results. The population spread was slightly smaller than the STS estimate, as expected from theory, i.e., the STS overestimates the population covariance. Interestingly, we found no statistically significant difference between STS-FSS and ITS-RSS insulin sensitivity estimates, showing that the RSS is not robust unless the data are analyzed with a population analysis method. We can speculate that, in a data-rich situation, STS and ITS would perform similarly; however, the average precision of the estimates improved dramatically with ITS (24, 18, 35, and 4%, respectively). This approach maximizes the information available from both the rich data set of the FSS and the knowledge of the population mean. The ITS, however, reveals its power in a (relatively) data-poor situation like the RSS. The average precision of the estimates in this case was very good: 23, 27, 37, and 4%, respectively. These numbers are similar to those obtained with ITS-FSS, but the process employed fewer data points. The gain of the ITS with respect to the STS, with both sampling schedules, is quite evident from the summary of the estimated parameter precision in all four cases shown in Fig. 2.

Last, we investigated the performance of the estimation methods applied to the same sampling schedule. To give an idea of how individual estimates from the STS map into individual estimates from the ITS, we compared the individual STS and ITS parameter estimates via linear regression; rank correlation of S_G and S_I was always between 0.97 and 1.00.

	DISCUSSION

TOP ABSTRACT INTRODUCTION MATERIALS AND METHODS RESULTS DISCUSSION REFERENCES

Population analysis methods find their natural application in the analysis of data-poor clinical studies, as when the number of samples available for each individual is rather small (e.g., 3 or 4). The natural arena is that of pharmacokinetics and pharmacodynamics, but there have been some recent contributions in the metabolism arena (24). We have found here that both rough (like the STS) and somewhat more sophisticated (like the ITS) methods give similar results for the population mean and variance. However, the population-based approach obtained the same results for population mean and variance with one-half of the total number of blood samples (an attractive feature). Moreover, population analysis allows a significant gain in the individual precisions of S_G and S_I (Fig. 2), thus giving more confidence in the overall results. A population method like the ITS allows the calculation of much more precise individual estimates for all of the subjects, even for those where the STS parameter precisions become unacceptable. For instance, our results demonstrate that an RSS is feasible for the sIVGTT, provided that a population-based method is used for data analysis; in fact, we found a statistically significant difference between the FSS and RSS S_I estimates with STS analysis. Note that the overestimation error associated with the STS estimate of between-individual variability increases as within-individual variability (related to the precision of parameter estimates and also referred to as intraindividual variability) increases (which is particularly evident when comparing the estimate of  for p₂ in Table 2, STS vs. ITS). This is the reason for the failure of STS when applied to the RSS case. Therefore, in an experimental protocol with very precisely estimated individual parameters, as in tracer studies (8), or with the mIVGTT (36), the STS would give quite reliable estimates of the population mean and variance, because within-individual variability would have only a small confounding effect. In contrast, in a situation where within-individual variability is as large as or larger than between-individual variability (i.e., the parameters are poorly estimated), the two effects would confound each other. Thus the sIVGTT protocol with RSS benefits considerably from population analysis.

The minimal model community widely uses measured values instead of predicted values when calculating ²_i,j in Eq. 2. This practice, however, introduces measurement error in the weights. Inspection of weighted residuals from the individual fits indicates that model misspecification is small (data not shown), thus suggesting that the difference would not be significant. Nevertheless, an extended least squares (ELS) or maximum likelihood criterion, where a term in the objective function with the logarithm of the data variance balances the use of predictions in the weighting, would better account for the variance in the data (2, 3)

(10)

Equation 10 gives twice the negative logarithm of the likelihood function of the data. If ²(p_j,t_ij) did not depend on the parameters and were constant across time (constant standard deviation), then the two methods (weighted and extended least squares) would be entirely equivalent.

Theory indicates that estimation of parameters with empirical Bayes estimation will result in greater precision (as calculated from the Fisher Information Matrix in Eqs. 3 and 9) than with the STS approach. However, the purpose of this study was to ascertain the magnitude and possible significance of such an improvement a priori unknown. Indeed, with extremely large variabilities of the metabolic indexes, we could expect negligible improvement. Also, we should note the (historical) fact that, in the clinic and even in large epidemiological studies (19), the STS is invariably used for S_G and S_I estimation, with both reduced and intensive sampling schedule (18-30 blood samples). We hope that this work can begin to introduce this large community to the concept that population-based methods and/or empirical Bayes approaches improve individual/group parameter quantification in the presence of large epidemiological data, and that it is feasible to use priors (when available and appropriate) to estimate individual patients' parameters (a practice sometimes referred to in pharmacokinetics as "Bayesian forecasting").

The improved precision of parameter estimates in this subject sample is due, at least in part, to the "shrinkage" of the individual estimates around the population mean (see Eq. 6). The not-so-hidden assumption is that a population mean indeed exists and that it is appropriate for all of the subjects involved in the study. Several levels of assumptions are embedded in kinetic modeling of population data: 1) that the mathematical model of the concentration time course is true; 2) that a certain measurement-error statistical structure exists (which allows the application of meaningful nonlinear regression); and 3) that, in the case of our version of the ITS, the shape of the population distribution of the parameters is multivariate normal. Examination of the consequences of all of these assumptions does not concern us here, since the present report is largely a "proof of concept" of the usefulness of population analysis in a reduced-data situation, and the ITS is only one possible algorithm for this purpose. We are thus not concerned here with either model accuracy or methods other than the ITS. However, some of the aforementioned assumptions, such as normality or unimodality of the population distribution of the parameters, can be relaxed, e.g., via mixture modeling, explicitly including covariants in the model, or performing nonparametric population analysis (see, e.g., Chapter 7 in Ref. 9). Davidian and Giltinan (9) comprehensively reviewed several different proposed estimation methods for population analysis of nonlinear models, all based on different degrees of approximation to the maximum-likelihood objective function. The attractiveness of the ITS derives from the need for a simple computational machinery based on successive iterations of the single-subject identification algorithm via a recursion formula. The ITS method is essentially an expectation-maximization (EM) algorithm (27). A related EM-type algorithm is the global two stage (GTS), which applies an iterative algorithm directly to estimated parameters and their within-individual covariances. Recently, Patron-Bizet et al. (22) compared the STS and the GTS used with a pharmacodynamic model and found that the GTS was superior. Steimer et al. (34) and Racine-Poon and Smith (23) compared the GTS and ITS, finding that they perform quite similarly. Although the ITS is computationally more intensive than the GTS, we chose to use the ITS here, because it allows a natural introduction of the concept of empirical Bayes individual estimation, which is at the root of the improved (asymptotic) precision. We have found the ITS method quite robust with respect to different starting values, with the main effect being the speed of convergence. After a few iterations, the ITS algorithm gets reasonably near its limit point. The literature suggests checking for convergence by calculation of the maximum-likelihood objective function, or an approximation thereof, at the point estimates (14).

Other available algorithms for the calculation of population parameters belong to the family of mixed-effects modeling, where the mean and the covariance of the population are fitted parameters in a maximum-likelihood context, and the computationally demanding maximum-likelihood integral is calculated via various approximations to the model function (29). Indeed, given that we have chosen a homogeneous population, mixed-effects modeling might be the best choice. Some of these approximations are currently available in the software NONMEM (28). Indeed, these algorithms are often preferred to the ITS, because they can be used even when individual estimates are not available in some subjects. They are based on the minimization of a well-defined objective function and return the precision of the estimates of population mean and covariance (and, in principle, also of each individual estimate). In contrast, the ITS algorithm does not readily allow calculation of the precision of population mean and variance. Because we were primarily interested in gauging the quantitative improvement of the individual posterior estimates, the ITS adequately met our purposes as a population analysis approach. Interestingly, De Gaetano et al. (10) used the first-order approximation to the population problem implemented in NONMEM to estimate minimal model population parameters in a group of 20 normal subjects by use of sIVGTT data. In that article, the precision of the population mean and covariance estimates improved with respect to the STS method. The mean and covariance estimates, however, did not change between NONMEM and STS (again, because the frequently sampled sIVGTT is a data-rich experiment). However, the authors did not report the values and precision of the individual parameter estimates (possibly because they are not reported by the NONMEM software).

Allowing for precise estimates of the parameters at the individual level in a reduced-data situation, and thus improving the a posteriori identifiability of the model parameters in a population, is only one of the possible applications of empirical Bayes estimation and population modeling in glucose metabolism and in intermediary metabolism in general. For example, there is ample pharmacokinetic/pharmacodynamic literature about the incorporation of demographic covariants into the mathematical model (12). Thus it is conceivable that population modeling of metabolic studies would, for instance, permit understanding the relationships of important features like visceral adiposity, weight, and age with metabolic indexes of insulin sensitivity and glucose effectiveness.