Skip Navigation


Biostatistics Advance Access originally published online on July 14, 2005
Biostatistics 2006 7(1):115-129; doi:10.1093/biostatistics/kxi044
This Article
Right arrow Abstract Freely available
Right arrow FREE Full Text (PDF) Freely available
Right arrow Supplementary data
Right arrow All Versions of this Article:
7/1/115    most recent
kxi044v1
Right arrow Alert me when this article is cited
Right arrow Alert me if a correction is posted
Services
Right arrow Email this article to a friend
Right arrow Similar articles in this journal
Right arrow Similar articles in PubMed
Right arrow Alert me to new issues of the journal
Right arrow Add to My Personal Archive
Right arrow Download to citation manager
Right arrowRequest Permissions
Right arrow Disclaimer
Google Scholar
Right arrow Articles by Ogden, R. T.
Right arrow Articles by Tarpey, T.
Right arrow Search for Related Content
PubMed
Right arrow PubMed Citation
Right arrow Articles by Ogden, R. T.
Right arrow Articles by Tarpey, T.
Social Bookmarking
Add to CiteULike   Add to Connotea   Add to Del.icio.us  
What's this?

© The Author 2005. Published by Oxford University Press. All rights reserved. For permissions, please e-mail: journals.permissions@oupjournals.org.

Estimation in regression models with externally estimated parameters

R. Todd Ogden*

Department of Biostatistics, Columbia University, 6th floor, 722 West 168th Street, New York, NY 10032, USA to166{at}columbia.edu

Thaddeus Tarpey

Department of Mathematics and Statistics, Wright State University, Dayton, OH 45435, USA

* To whom correspondence should be addressed.


   SUMMARY
TOP
 SUMMARY
1. INTRODUCTION
 2. THE BASIC MODEL
 3. ASYMPTOTIC RESULTS
 4. BOOTSTRAP STANDARD ERRORS
 5. APPLICATIONS
 6. DISCUSSION
 REFERENCES
 
In many regression applications, some of the model parameters are estimated from separate data sources. Typically, these estimates are plugged into the regression model and the remainder of the parameters is estimated from the primary data source. This situation arises frequently in compartment modeling when there is an external input function to the system. This paper provides asymptotic and bootstrap-based approaches for accounting for all sources of variability when computing standard errors for estimated regression model parameters. Examples and simulations are provided to motivate and illustrate the ideas.

Keywords: Bootstrap; Compartment models; PET imaging


   1. INTRODUCTION
TOP
 SUMMARY
 1. INTRODUCTION
2. THE BASIC MODEL
 3. ASYMPTOTIC RESULTS
 4. BOOTSTRAP STANDARD ERRORS
 5. APPLICATIONS
 6. DISCUSSION
 REFERENCES
 
It is a general statistical practice to report standard errors associated with estimated model parameters. Typically, model parameters are estimated by a single data source, but there are many examples in which two or more sources are used. A common example occurs in chemistry and physics applications where nonlinear models are defined in terms of well-known constants (e.g. Planck's constant and Faraday's constant) whose estimates have been previously established. The uncertainties associated with such physical constants are generally very small. However, in other applications, some of the model parameters must be estimated by external data sources and the uncertainties associated with these estimates can be quite large. We shall call the parameters estimated by external data sets ‘input’ parameters. In routine applications, the variability associated with the estimation of input parameters is often ignored.

The following illustrate examples of regression models requiring input parameters that are estimated by separate data sets.

EXAMPLE 1 (POSITRON EMISSION TOMOGRAPHY IMAGING FOR NEURORECEPTOR MAPPING) Positron emission tomography (PET) imaging of the brain often involves kinetic modeling (see, e.g. Gunn et al., 2001Go) or one of its variants (Ogden, 2003Go; Gunn et al., 2002Go). The primary data source is measurements made by the PET camera of the concentration of a ligand throughout the brain. In such an application, the objective is to estimate the density of the target receptor throughout the brain. This is accomplished by modeling the kinetics of the ligand and estimating the relevant rate constants. The density of the receptor may be estimated by computing a function of the estimated rate parameters. These estimated densities may then be compared among different patient populations or treatment groups. Since there is substantial noise involved with imaging and modeling, in addition to large intrasubject variability, these comparisons often lack power. With standard error estimates, power can be dramatically increased.

In general kinetic modeling, the ligand concentration at a given brain location at time t, denoted C(t), is given in terms of a convolution:

(1.1)

where g(t;ß) represents the kinetics of the brain compartmental system and involves unknown rate parameters ß. The input function Cp(t;{alpha}) represents the concentration of the ligand in the plasma at time t and the input parameters {alpha} are estimated separately using data obtained directly from blood samples drawn during the imaging experiment.

EXAMPLE 2 (ZINC METABOLISM) Foster et al. (1979)Go used a kinetic model to study the early phases of zinc metabolism in humans using data from Aamodt et al. (1979)Go. In this study, a subject ingests a zinc isotope and the concentration of zinc is modeled over time in various compartments. The objective is to construct a model for the movement of zinc throughout the body. Results from such a study are often displayed as a table of estimated rate parameters. Standard errors would help in assessing the identifiability of each parameter and provide a guide in determining an appropriate model.

Assuming that the rate at which zinc enters the urine is proportional to the amount of zinc in the plasma, the amount of administered zinc in the urine compartment is modeled as

(1.2)

The plasma over time is modeled as a sum of four exponentials:

(1.3)

In applications, the plasma concentration curve is fitted using standard nonlinear least-squares methodology on combined data from several subjects. Subsequently, estimates for the input parameters {alpha} are plugged into (1.2) in order to estimate subject-specific parameters for the urine compartment.

EXAMPLE 3 (ENVIRONMENTAL STUDIES OF TOXICANTS) The concentration of a toxicant in an organism measured over time C(t) is often modeled as

(1.4)

where the parameters ßu and ßd are the uptake and elimination rates of the toxicant, respectively, in the organism and Ce is an input parameter to the model representing the toxicant level in the environment (e.g. Bailer et al. (2000Go; Wheeler and Bailer (2003Go; Lotufo et al. (2000Go; Rand et al. (1995Go). The primary objective in such studies is often to determine the bioconcentration factor (BCF), defined as the ratio of ßu to ßd. The BCF is a key concern of environmental scientists as it quantifies how much chemical accumulates in an organism relative to the environmental exposures. The constant Ce is estimated separately from water or sediment samples and the uncertainty in estimating Ce is generally ignored when ascertaining the uncertainty in the estimated BCF. However, the validity of statistical inferences in environmental studies concerning changes of chemical accumulations in organisms requires accurate standard errors for estimated BCFs, which in turn requires the incorporation of the variability in estimating Ce.

All these examples fall under the general heading of compartment models (e.g. Matis and Wehrly, 1994Go; Jacquez, 1985Go) which include input parameters that are estimated from separate data sources.

The problem of estimating input parameters from a separate data set is closely related to models where separate data are available for covariates measured with error. For example, Higgins et al. (1997)Go consider a two-step approach using maximum likelihood estimation where time-varying covariates are estimated in the first step and are then used to obtain estimates of the parameters of primary importance. A parametric bootstrap approach is used to obtain standard errors that account for the uncertainty of estimating the covariates. The covariates in their setup are analogous to our input parameters. Wu (2002)Go considers a similar model with the added complication of censoring, proposing a joint modeling approach that simultaneously estimates the covariates and model parameters using a Monte Carlo expectation-maximization (EM) algorithm. The estimated standard errors are then obtained by approximating the information matrix.

In Section 2, a general regression model is provided for situations involving input parameters estimated from separate data sources. Section 3 gives asymptotic results for the sampling distribution of the estimated parameters when the input parameters are estimated. The resulting asymptotic covariance matrix can be used to provide approximate standard errors of model parameters. Section 4 describes a bootstrap alternative to asymptotic standard error estimation which is particularly useful for small- or moderate-sized samples. The proposed methodologies are explored for two situations via simulation in Section 5, and some general discussion is given in Section 6.


   2. THE BASIC MODEL
TOP
 SUMMARY
 1. INTRODUCTION
 2. THE BASIC MODEL
3. ASYMPTOTIC RESULTS
 4. BOOTSTRAP STANDARD ERRORS
 5. APPLICATIONS
 6. DISCUSSION
 REFERENCES
 
The three examples mentioned in Section 1 are special cases of the following general regression model: Let Y1, ..., Yn denote a sample from a known regression model depending on unknown parameter vectors {alpha} and ß:

(2.1)

In (2.1), the xi are the fixed, known values of predictors and the {epsilon}i are a sequence of i.i.d. mean zero random variables with Var({epsilon}i) = {sigma}2. The parameters {alpha} are input parameters which are to be estimated from a separate data source. Primary interest is in estimating ß and computing reliable standard errors for .

It is common practice to estimate ß from the Y data in (2.1) using standard regression techniques with {alpha} replaced by its estimated value . This is appropriate for the purposes of point estimation, but without taking into account the variability in , estimated standard errors of will tend to be underestimated, perhaps severely so. If the final analysis consists of comparing estimated ß values for a sample of individual subjects across, for instance, treatment groups, then reliable standard errors can dramatically enhance the power of such comparisons.

Let (x1, Y1), ..., (xn, Yn) represent data observed from (2.1). If {alpha} were known, then the covariance matrix of can be estimated using standard least-squares theory with () = 2(X'X)–1, where X is the matrix of partial derivatives of f with respect to ß evaluated at . Because {alpha} is not known in practice, this covariance matrix will tend to produce variances that are too small.

The scenario of interest here is that an estimate of {alpha}, as well as an estimated variance–covariance matrix for , is available from another source. For this situation, the standard least-squares approach is to determine the value of ß minimizing Defining

the estimator of ß in practice is = H(, Y) whose value therefore depends on realizations of two other random variables. If the estimator of ß in the case of known {alpha} is approximately unbiased for ß and the function H(, Y) is well approximated by a first-order Taylor series about {alpha} and {alpha} must be estimated, then

showing that is approximately unbiased provided that is approximately unbiased for {alpha} and independent of Y, and that the expectations exist.

However, the variability in the is amplified when an estimator of {alpha} is used in place of its true value. If is used in the estimation of ß, then the well-known variance decomposition formula can be used to relate Var() to Var(|):

(2.2)

Note that Var()–E[Var(|)] is a positive semidefinite matrix, indicating that the conventional approach based on Var(|) underestimates the true variance–covariance matrix of on average.


   3. ASYMPTOTIC RESULTS
TOP
 SUMMARY
 1. INTRODUCTION
 2. THE BASIC MODEL
 3. ASYMPTOTIC RESULTS
4. BOOTSTRAP STANDARD ERRORS
 5. APPLICATIONS
 6. DISCUSSION
 REFERENCES
 
In this section, the asymptotic distribution for the least-squares estimator of ß in (2.1) is derived when the {alpha} parameters are estimated from a separate data set. The asymptotic results derived here are related to previous work on pseudo maximum likelihood by Gong and Samaniego (1981)Go and on estimating equations by Carroll et al. (1995)Go. Gong and Samaniego (1981)Go derive asymptotic normality results similar to our result below except that they consider maximum likelihood estimation with a single data source where a set of nuisance parameters is estimated using a method different from maximum likelihood. Carroll et al. (1995)Go ((1995)Go, pp. 267–269) consider the estimation of ß based on a random sample of n individuals, from which estimators of {alpha} are also computed.

We shall assume that the function f in (2.1) is differentiable with respect to the parameters. Let {alpha} = ({alpha}1, ..., {alpha}p)' and ß = (ß1, ..., ßq)'. Denote by {alpha}* and ß* the ‘true’ values of the parameters. Suppose that is a strongly consistent estimator of {alpha} that is based on another data set of sample size n1 and that

where {Rightarrow} represents convergence in distribution. We shall assume that the sample sizes n1 and n go to infinity at the same rate with where 0 < c < {infty}. Furthermore, we assume that is independent of {epsilon}1, ..., {epsilon}n. Some necessary additional notation follows:

and

where {partial}f(xi; {alpha}*, ß*)/{partial}ß represents the column vector of length q with the j th element given by {partial}f(xi; {alpha}, ß)/{partial}ßj evaluated at the point {alpha} = {alpha}* and ß = ß*. (The derivative with respect to ß' is the transpose of this.) The matrices {Omega}1 and {Omega}2 have dimensions q x q and p x q, respectively. We assume that the above limits exist and that {Omega}1 is nonsingular. We also assume appropriate regularity conditions ensuring that is consistent for ß in the case that {alpha} is known [for additional details, see Chapter 12 of Seber and Wild (1989)Go]. The following result demonstrates the asymptotic normality of . The proof of the theorem, given in the supplementary material (www.biostatistics.oxfordjournals.org), follows the general approach taken by Amemiya (1983)Go.

THEOREM 3.1 Assuming the appropriate regularity conditions, is asymptotically normal with mean 0 and variance–covariance matrix

The square root of the diagonal of the asymptotic covariance matrix can be used to obtain standard errors for the rate parameters. In practice, this asymptotic covariance matrix must be estimated. In particular, the estimated covariance matrix of can be used in place of {Psi}. Similarly, {Omega}1 and {Omega}2 can be estimated by taking the appropriate partial derivatives and evaluating them using the parameter estimates and averaging over all observations.

It is worth noting that the variance–covariance matrix of based on the asymptotic formula above has two components just as in (2.2). The first component is ordinarily computed as the estimated variance–covariance matrix in any nonlinear least-squares setting—the ‘naive’ estimator of the standard error in which the variability of the estimated input function is neglected. The second component arises from the variability in estimating {alpha}. This component depends on the factor c = limn,n1->{infty}n/n1 which can be estimated by the ratio of the actual sample sizes. Thus, if the input parameters are estimated with a large sample size n1 relative to n, then the constant c will be small and the variability of will play a diminished role in affecting the variability of .

Theorem 3.1 can be generalized to handle situations where the i.i.d. assumption may not hold (e.g. nonconstant variance or longitudinal studies) or where the input data set is not independent of the primary data set. Consider the following (possibly nonlinear) model for the input data set

where {alpha} represents the input parameters as before, the {xi}i are mean zero errors, and the function h is twice differentiable with respect to {alpha}. Stack the input errors {xi} = ({xi}1, ..., {xi}n1)' on the errors {epsilon} from (2.1) and let the covariance matrix for the joint distribution of ({xi}', {epsilon}')' be given by

Let Gn = [{partial}f/{partial}ß] denote the n x q matrix of partial derivatives of the nonlinear regression function f evaluated at x1, ..., xn. Similarly, let Hn1 = [{partial}h/{partial}{alpha}] denote the n1 x p matrix of partial derivatives of h with respect to {alpha} evaluated at w1, ..., wn1. Implementing the asymptotic results of Amemiya (1983)Go ((1983)Go, p. 354) for nonlinear regression with correlated errors and incorporating the covariance between the input data and primary data, it follows that as n1 and n go to infinity, is asymptotically normal with mean 0 and covariance matrix

(3.1)

where

The proof of (3.1) is outlined in the supplementary material (www.biostatistics.oxfordjournals.org).


   4. BOOTSTRAP STANDARD ERRORS
TOP
 SUMMARY
 1. INTRODUCTION
 2. THE BASIC MODEL
 3. ASYMPTOTIC RESULTS
 4. BOOTSTRAP STANDARD ERRORS
5. APPLICATIONS
 6. DISCUSSION
 REFERENCES
 
An alternative to the asymptotic approach from the previous section is to use the bootstrap (Efron and Tibshirani, 1993Go) for estimating the standard errors. In some applications, the asymptotic results may be difficult to obtain if the computation of the partial derivative matrices needed for the asymptotic covariance matrix is complex. Also, standard errors based on asymptotic computations may not provide good estimates for data sets of small to moderate size. The implementation of the bootstrap approach must incorporate the variability of the estimation of the input parameters into the calculation of standard errors for the rate parameters.

If the points xi are fixed by design, then it is more appropriate to compute bootstrap samples of the Y data using the residuals from the fit rather than resampling the raw data (see Efron and Tibshirani, 1993Go, see, Section 9.5). The proper resampling scheme to generate resampled input parameter estimates depends on the nature of the input data. If raw input data are not available, a parametric bootstrap can be applied, provided both an estimate of {alpha} and its estimated variance–covariance matrix are available. The choice of parametric bootstrap algorithm should depend upon the nature of the input data, but if it can be reasonably assumed that is normally distributed with mean {alpha}* and variance–covariance matrix then each resampled set of estimated parameters could be generated from an distribution (replacing {Psi} with its estimate where appropriate).

There are a number of variants of the general bootstrap algorithm that would be appropriate for this situation, depending on whether original data or residuals are to be resampled, or if a parametric bootstrap is to be utilized. For instance, if data are thought to be serially correlated, a moving blocks bootstrap algorithm (Künsch, 1989Go; Liu and Singh, 1992Go) may be appropriate. This would involve resampling blocks of adjacent residuals rather than sampling from the residuals without regard to their order as in the usual bootstrap algorithm. This helps to preserve the correlation structure inherent in the original data.

For any selected set of bootstrap schemes, the general algorithm is:

Bootstrap Algorithm 1 (one-stage bootstrap)

  1. For each k = 1, ..., B,
    (a) obtain a bootstrap sample of the {alpha} data and calculate the resulting k value;
    (b) obtain a bootstrap sample of the Y data and compute k from these data by plugging in the value of k in the expression for f(xi; {alpha}, ß).

  2. The bootstrap standard error is the standard deviation of components of 1, ..., B.

An alternative approach would be to conduct a two-stage (nested) bootstrap. This would allow for the estimation of each component of the variance decomposition formula (2.2), which may be useful for purposes of design as well as for ascertaining the individual sources of variability of the parameter estimate. In the case that the resampling for the Y data is based on residuals from a fitted curve involving , this will also allow for a richer set of residuals for the Y data (since a new set of residuals would be computed for each bootstrapped value of ). This algorithm is given by:

Bootstrap Algorithm 2 (two-stage bootstrap)

  1. The outer loop: For each k = 1, ..., B1,
    (a) obtain a bootstrap sample of the {alpha} data and calculate the resulting k value.
    (b) The inner loop: For each {ell} = 1, ..., B2,
    (i) using the computed value of k, obtain a bootstrap sample of the Y data and compute {ell}(k).

    (c) Compute k, the average of the {ell}(k) values.
    (d) Compute k, the sample variance–covariance matrix of the {ell}(k) values.

  2. Estimate Var(E[|]) by computing the sample variance–covariance matrix of the k values.
  3. Estimate E[Var(|)] by computing the matrix average of the k s.
  4. The bootstrap estimate of the (unconditional) variance–covariance matrix of is the sum of these two components.

Note that in the two-stage algorithm, each k is an estimate of the corresponding E[|k] and k is an estimate of Var(|k).


   5. APPLICATIONS
TOP
 SUMMARY
 1. INTRODUCTION
 2. THE BASIC MODEL
 3. ASYMPTOTIC RESULTS
 4. BOOTSTRAP STANDARD ERRORS
 5. APPLICATIONS
6. DISCUSSION
 REFERENCES
 
We return to Examples 1 and 2 from Section 1 to illustrate the computation of standard errors when input parameters are estimated via an external data set. Simulations were performed to examine the behavior of estimated standard errors using the asymptotic formula as given in the theorem in Section 3, both bootstrap approaches as described in Section 4, and the ‘naive’ approach (in which the variability in estimating the input parameters is neglected). We begin with the relatively simple model of zinc metabolism in urine as described in Example 2 of Section 1.

5.1 Zinc metabolism

To provide a simple illustration of the methods described in this paper, data on zinc concentration in urine (1.2) and in plasma (1.3) were simulated using parameter values based on the results reported in Foster et al. (1979)Go. As is common in these situations, we assumed a multiplicative error and performed the fitting using a log-transformed model:

In this simple case, the estimate of ß is just the arithmetic mean of the values The input parameters in (1.3) were estimated using nonlinear least squares. Gaussian errors were generated throughout, and all simulations were run using Matlab software. A plasma input function was estimated based on averaged data from 17 subjects who ingested a zinc isotope. Concentration of the isotope in other compartments (liver, urine, etc.) was measured for the same 17 subjects at several time points over 5 days and averaged across subjects. For the simulation, parameter values were chosen to conform with the data reported in Aamodt et al. (1979)Go. In particular, input parameter values for the plasma were set at {alpha} = (0.790, 176, 0.175, 73.4, 0.022, 5.87, 0.013, 0.053)', and the value of the rate parameter ß was set to –1.46. The standard deviation for the plasma data was set to 0.0015 and that for the urine data was set to 0.02.

The sampling distribution of was simulated in two cases and is shown in Figure 1. The dashed curve in Figure 1 shows the distribution of the estimated rate parameter in the case when the input parameters are known without error while the solid curve shows the distribution of when the input parameters are estimated. In both cases, is approximately unbiased, but its variability is much greater when the input function is estimated.



View larger version (11K):
[in this window]
[in a new window]
  		Fig. 1. Distributions of simulated values of for the case in which the input function is known (dashed curve) and estimated (solid curve). The vertical line indicates the true value of ß.

 
Next, simulations were run in order to compare the standard errors of the naive, asymptotic, and bootstrap approaches. Two separate simulations were run for the zinc model differing only in the number of time points at which plasma samples were to be taken. In the first run, the time points were set to the original time points as in Aamodt et al. (1979)Go; in the second (to allow for better comparison of the asymptotic expression for the standard error), the plasma sampling rate was increased.

The first simulation run used the original time points from the Foster et al. paper (for both the plasma and the urine data), giving n1 = 16 time points for the plasma data and n = 12 urine samples over the course of about 5 days. A total of 1000 simulated data sets were generated. For the one-stage bootstrap algorithm, we used B = 200. For the two-stage bootstrap, we set B1 = B2 = 200. The results of this simulation are displayed in Figure 2. Corresponding summary statistics are provided in Table 1. The five smoothed histograms in Figure 2 were generated from the 1000 simulated data sets. Note that these are all on the standard deviation scale. The solid curve represents the naive standard error estimates (neglecting the variability of the input data). The short-length dashed curve represents the distribution of the average of the values as described in Section 4, the bootstrap estimate of E[Var(|)]. The dotted curve shows the distribution of the estimates of Var(E[|]), the sample variance–covariance matrix of the values. Combining these two gives the medium-length dashed curve representing the total bootstrap estimate of the standard errors from the two-stage algorithm. Finally, the long-length dashed curve shows the estimated distribution of standard errors computed using the one-stage bootstrap algorithm. The vertical line represents the ‘true’ standard error, computed as the standard deviation of the values from a separate simulation of 10 000 data sets.



View larger version (17K):
[in this window]
[in a new window]
  		Fig. 2. Smoothed histograms of the results from the simulation study based on the original data as reported in Section 5.1. The vertical line represents the ‘true’ standard error for this situation, as computed by a separate simulation study of 10 000 simulated data sets.

 

View this table:
[in this window]
[in a new window]
  		Table 1. Summary statistics of simulation results presented in Figure 2

 
It is clear in Figure 2 that the naive approach to standard error estimation badly underestimates the truth. The two-stage bootstrap algorithm demonstrates that most of the overall variability of is due to the variability in the estimates of {alpha}, so the Var(E[|]) term tends to dominate the overall estimate of the standard error. There is little difference between the one-stage and the two-stage bootstrap standard errors in this example. Figure 2 shows that both bootstrap approaches show some negative bias. This may be due in part to the small sample size, in particular, the necessity to estimate eight plasma parameters from only 16 data points using nonlinear regression.

The asymptotic formula given in Section 3 when applied to the simulation data using original time points gave very large and unreliable standard error estimates for about 1–2% of the simulated data sets. This is due to numerical instability of estimating the derivatives for the plasma function based on a small sample size, as the denominator of the estimated derivatives can take values near 0. If these values are excluded, the distribution of the asymptotic standard error estimates closely resembles the bootstrap results.

In order to compare the performance of the asymptotic formula with the bootstrap approach, a second simulation study was conducted using n = 12 urine samples [again, as in the original Aamodt et al. (1979)Go paper] but n1 = 80 plasma samples (computed by equally subdividing each interval in the original plasma data). In this new setting, 1000 simulated data sets were generated and standard error estimates were computed using both bootstrap algorithms as well as the asymptotic formulas. The results of this simulation are shown in Figure 3 and summary statistics are given in Table 2.



View larger version (19K):
[in this window]
[in a new window]
  		Fig. 3. Smoothed histograms of the results from the simulation study with increased plasma sampling as reported in Section 5.1. The vertical line represents the ‘true’ standard error for this situation, as computed by a separate simulation study of 10 000 simulated data sets.

 

View this table:
[in this window]
[in a new window]
  		Table 2. Summary statistics of simulation results presented in Figure 3

 
Figure 3 displays the same information as Figure 2, but for the new, more frequently sampled plasma data. In this situation, the bootstrap standard error is more balanced between its two components and the two components taken together provide a nearly unbiased estimate of the ‘true’ standard error. The one-stage bootstrap algorithm is close but exhibits some negative bias. Again, as may be expected, the naive estimator conforms closely to the E[Var(|)] term.

For the setting with increased plasma sampling (n1 = 80), the distribution of the bootstrap estimates nearly coincides with that of the asymptotic estimates, and both are approximately unbiased for the true standard error.

5.2 PET imaging

In the PET imaging application described in Example 1 of Section 1, the simplest model in common use is a two-compartment model [(one compartment for the plasma and one for a brain region of interest (ROI)]. For model (1.1), the g(t; ß) function is a simple exponential, so the concentration of the ligand in a given brain location is given by for unknown rate parameters k1 and k2. The plasma concentration function is commonly modeled by a piecewise expression:

(5.1)

Plasma data (consisting of measured concentration based on total radioactivity counts corrected for rate of decay and metabolism) are taken over time from arterial samples while the subject is being scanned:

The parameters involved in the expression of Cp(t) in (5.1) are estimated using nonlinear least squares.

The convolution defining C(t) can be worked out analytically as

Through the PET imaging modality, concentration is measured for each of several time-ordered ‘frames’ at each location on a 3D grid. In some applications, aggregate measures are computed for several anatomically defined ROI by combining all data from voxels contained within the ROI. Assuming the two-compartment model, the brain data are given by

One outcome measure of interest is the ‘total volume of distribution,’ defined as VT = K1/k2.

Simulation parameters are set to mimic real data from PET imaging experiments as in Parsey et al. (2000)Go. The ‘true’ plasma and brain regression functions are displayed in Figure 4.



View larger version (10K):
[in this window]
[in a new window]
  		Fig. 4. ‘True’ plasma (left panel) and ‘true’ brain (right panel) concentrations, taken from real PET data.

 
This PET example is more complicated than the zinc example and requires a number of modifications to the bootstrap procedures. First, in practice, the plasma data points are measured together with an estimate of their standard errors. These errors are used to determine weights in a weighted least-squares fitting algorithm that substantially improves the quality of the fit and the stability of estimates. The weights, taken from an actual study, are used to generate the plasma data in the simulation. Additional improvement of the bootstrap procedure was achieved by centering the residuals at zero and dividing the residuals by where hii represents the nonlinear regression analogue to the ‘hat’ matrix; see Davison and Hinkley (1997)Go (Chapter 6 (1997)Go, Chapter 6) for a discussion of this procedure in the linear regression case.

Figure 5 and Table 3 summarize the results of the simulation. As in the zinc example, the standard error of the naive estimator badly underestimates the true standard error. On the other hand, both the bootstrap and asymptotic results are reasonably close to the true standard error.



View larger version (15K):
[in this window]
[in a new window]
  		Fig. 5. Smoothed histograms of one-stage bootstrap and asymptotic standard errors for the PET simulation study as described in Section 5.2. The vertical line represents the ‘true’ standard error for this situation, as computed by a separate simulation study of 10 000 simulated data sets.

 

View this table:
[in this window]
[in a new window]
  		Table 3. Summary statistics of simulation results presented in Figure 5

 
For the purposes of illustration, these methods were applied to data from the ventral striatum of a PET imaging subject. For this study with an estimated VT of 43.8, the bootstrap estimate of the standard error was 11.1 and the asymptotic estimate was 12.6. The naive estimates of the standard error neglecting the plasma effect are 9.01 for the bootstrap and 8.60 for the asymptotic estimator. This example illustrates that in current PET imaging practice, the naively estimated standard errors for parameter estimates such as VT are underestimating the true standard error.


   6. DISCUSSION
TOP
 SUMMARY
 1. INTRODUCTION
 2. THE BASIC MODEL
 3. ASYMPTOTIC RESULTS
 4. BOOTSTRAP STANDARD ERRORS
 5. APPLICATIONS
 6. DISCUSSION
REFERENCES
 
This paper proposes two general methods for the estimates of standard errors in regression models defined in terms of input parameters that are estimated separately. The usual method of estimation neglects the variability in estimating the input parameters. In contrast, the methods described here take the variability of the estimated input parameters into account and consequently yield much better estimates of the true variability in the estimated model parameters.

The asymptotic formulas gave reasonably good standard error estimates even for relatively small sample sizes. Computationally, the asymptotic standard errors can be computed very fast. The quick computation makes the asymptotic standard error approach very attractive in applications like PET imaging (Section 5.2) that involve model fittings at hundreds of thousands of voxels. However, a drawback of the asymptotic approach is that the coding routine can require painstaking calculations and is prone to error.

On the other hand, the bootstrap approaches tend to work quite well for a variety of settings. Coding the bootstrap is straightforward, but, of course, they require considerably more computational expense than the asymptotic approach.


   ACKNOWLEDGMENTS
 
The authors are grateful to the editor, associate editor, and referee for helpful and insightful comments that have led to a substantial improvement of this paper.


   REFERENCES
TOP
 SUMMARY
 1. INTRODUCTION
 2. THE BASIC MODEL
 3. ASYMPTOTIC RESULTS
 4. BOOTSTRAP STANDARD ERRORS
 5. APPLICATIONS
 6. DISCUSSION
 REFERENCES
 

    AAMODT, R. L., RUMBLE, W. F., JOHNSTON, G. S., FOSTER, D. AND HENKIN, R. I. (1979). Zinc metabolism in humans after oral and intravenous administration of Zn69-m. The American Journal of Clinical Nutrition 32, 559–569.[Abstract/Free Full Text]

    AMEMIYA, T. (1983). Non-linear regression models. In Griliches, Z. and Intriligator, M. D. (eds), Handbook of Econometrics, Volume 1. Amsterdam: North Holland, pp. 333–389.

    BAILER, A. J., WALKER, S. E. AND VENUS, K. J. (2000). Estimating and testing bioconcentration factors. Environmental Toxicology and Chemistry 19, 2338–2340.[CrossRef]

    CARROLL, R. J., RUPPERT, D. AND STEFANSKI, L. A. (1995). Measurement Error in Nonlinear Models. London: Chapman and Hall.

    DAVISON, A. C. AND HINKLEY, D. V. (1997). Bootstrap Methods and Their Application. Cambridge: Cambridge University Press.

    EFRON, B. AND TIBSHIRANI, R. J. (1993). An Introduction to the Bootstrap. London: Chapman and Hall.

    FOSTER, D. M., AAMODT, R. L., HENKIN, R. I. AND BERMAN, M. (1979). Zinc metabolism in humans: a kinetic model. American Journal of Physiology 237, R340–R349.

    GONG, G. AND SAMANIEGO, F. J. (1981). Pseudo maximum likelihood estimation: theory and applications. Annals of Statistics 9, 861–869.

    GUNN, R. N., GUNN, S. R. AND CUNNINGHAM, V. J. (2001). Positron emission tomography compartmental models. Journal of Cerebral Blood Flow and Metabolism 21, 635–652.[CrossRef][Medline]

    GUNN, R. N., GUNN, S. R., TURKHEIMER, F. E., ASTON, J. A. D. AND CUNNINGHAM, V. J. (2002). Positron emission tomography compartmental models: a basis pursuit strategy for kinetic modeling. Journal of Cerebral Blood Flow and Metabolism 22, 1425–1439.[Medline]

    HIGGINS, K. M., DAVIDIAN, M. AND GILTINAN, D. M. (1997). A two-step approach to measurement error in time-dependent covariates in nonlinear mixed-effects models, with applications to IGF-I pharmacokinetics. Journal of the American Statistical Association 92, 436–448.

    JACQUEZ, J. A. (1985). Compartmental Analysis in Biology and Medicine, 2nd edition. Ann Arbor, MI: University of Michigan Press.

    KÜNSCH, H. R. (1989). The jackknife and the bootstrap for general stationary observations. Annals of Statistics 17, 1217–1241.

    LIU, R. Y. AND SINGH, K. (1992). Moving blocks jackknife and bootstrap capture weak dependence. In LePage, R. and Billard, L. (eds), Exploring the Limits of Bootstrap. New York: Wiley, pp. 225–248.

    LOTUFO, G. R., LANDRUM, P. F., GEDEON, M. L., TIGUE, E. A. AND HERCHE, L. R. (2000). Comparative toxicity and toxicokinetics of DDT and its major metabolites in freshwater amphipods. Environmental Toxicology and Chemistry 19, 369–379.

    MATIS, J. H. AND WEHRLY, T. E. (1994). Compartmental models of ecological and environmental systems. In Patil, G. P. and Rao, C. R. (eds), Handbook of Statistics, Volume 12. Amsterdam: Elsevier, pp. 583–613.

    OGDEN, R. T. (2003). Estimation of kinetic parameters in graphical analysis of PET imaging data. Statistics in Medicine 22, 3557–3568.[CrossRef][ISI][Medline]

    PARSEY, R. V., KEGELES, L. S., HWANG, D. R., SIMPSON, N., ABI-DARGHAM, A., MAWLAWI, O., SLIFSTEIN, M., VAN HEERTUM, R. L., MANN, J. J. AND LARUELLE, M. (2000). In vivo quantification of brain serotonin transporters in humans using [11C]McN 5652. Journal of Nuclear Medicine 41, 1465–1477.[Abstract/Free Full Text]

    RAND, G. M., WELLS, P. G. AND MCCARTY, L. S. (1995). Introduction to aquatic toxicology. In Rand, G. M. (ed.), Fundamentals of Aquatic Toxicology: Effects, Environmental Fate, and Risk Assessment. Washington DC: Taylor and Francis, pp. 3–67.

    SEBER, G. A. F. AND WILD, C. J. (1989). Nonlinear Regression. New York: Wiley.

    WHEELER, M. W. AND BAILER, A. J. (2003). A simulation study of methods for constructing confidence intervals for bioaccumulation factors. Environmental Toxicology and Chemistry 22, 921–927.[Medline]

    WU, L. (2002). A joint model for nonlinear mixed-effects models with censoring and covariates measured with error, with application to AIDS studies. Journal of the American Statistical Association 97, 955–964.[CrossRef]

    Received May 19, 2004; revised October 20, 2004; revised April 8, 2005; revised June 14, 2005; revised June 22, 2005; accepted for publication July 13, 2005.


    Add to CiteULike CiteULike   Add to Connotea Connotea   Add to Del.icio.us Del.icio.us    What's this?


    This article has been cited by other articles:


    		Home page
    		 JNMHome page
    M. S. Milak, A. J. Severance, R. T. Ogden, J. Prabhakaran, J.S. D. Kumar, V. J. Majo, J. J. Mann, and R. V. Parsey
    Modeling Considerations for 11C-CUMI-101, an Agonist Radiotracer for Imaging Serotonin 1A Receptor In Vivo with PET
    J. Nucl. Med., April 1, 2008; 49(4): 587 - 596.
    [Abstract] [Full Text] [PDF]

    This Article
    Right arrow Abstract Freely available
    Right arrow FREE Full Text (PDF) Freely available
    Right arrow Supplementary data
    Right arrow All Versions of this Article:
    7/1/115    most recent
    kxi044v1
    Right arrow Alert me when this article is cited
    Right arrow Alert me if a correction is posted
    Services
    Right arrow Email this article to a friend
    Right arrow Similar articles in this journal
    Right arrow Similar articles in PubMed
    Right arrow Alert me to new issues of the journal
    Right arrow Add to My Personal Archive
    Right arrow Download to citation manager
    Right arrowRequest Permissions
    Right arrow Disclaimer
    Google Scholar
    Right arrow Articles by Ogden, R. T.
    Right arrow Articles by Tarpey, T.
    Right arrow Search for Related Content
    PubMed
    Right arrow PubMed Citation
    Right arrow Articles by Ogden, R. T.
    Right arrow Articles by Tarpey, T.
    Social Bookmarking
    Add to CiteULike   Add to Connotea   Add to Del.icio.us  
    What's this?