Weighted Deming – Rocke-Lorenzato - known sigma, kappa

Description

This code fits the weighted Deming regression on predicate readings (X) and test readings (Y), with user-supplied Rocke-Lorenzato ("RL") parameters sigma (\sigma) and kappa (\kappa).

Usage

PWD_RL(X, Y, sigma, kappa, lambda=1, alpha=NA, beta=NA, mu=NA, epsilon=1e-8)

Arguments

Details

The Rocke-Lorenzato precision profile model assumes the following forms for the variances, with proportionality constant \lambda:

• predicate precision profile model: g_i = var(X_i) = \lambda\left(\sigma^2 + \left[\kappa\cdot \mu_i\right]^2\right) and

• test precision profile model: h_i = var(Y_i) = \sigma^2 + \left[\kappa\cdot (\alpha + \beta\mu_i)\right]^2.

The algorithm uses maximum likelihood estimation. Proportionality constant \lambda is assumed to be known or estimated externally.

Value

A list containing the following components:

Author(s)

Douglas M. Hawkins, Jessica J. Kraker krakerjj@uwec.edu

References

Hawkins DM and Kraker JJ (in press). Precision Profile Weighted Deming Regression for Methods Comparison. The Journal of Applied Laboratory Medicine. doi:10.1093/jalm/jfaf183

Hawkins DM (2014). A Model for Assay Precision. Statistics in Biopharmaceutical Research, 6, 263-269. doi:10.1080/19466315.2014.899511

Examples

# library
library(ppwdeming)

# parameter specifications
sigma <- 1
kappa <- 0.08
alpha <- 1
beta  <- 1.1
true  <- 8*10^((0:99)/99)
truey <- alpha+beta*true
# simulate single sample - set seed for reproducibility
set.seed(1039)
# specifications for predicate method
X     <- sigma*rnorm(100)+true *(1+kappa*rnorm(100))
# specifications for test method
Y     <- sigma*rnorm(100)+truey*(1+kappa*rnorm(100))

# fit RL with given sigma and kappa
RL_results <- PWD_RL(X,Y,sigma,kappa)
cat("\nWith given sigma and kappa, the estimated intercept is",
    signif(RL_results$alpha,4), "and the estimated slope is",
    signif(RL_results$beta,4), "\n")

Estimate of Variance Profile Functions (proportional)

Description

This code estimates the variance profiles, assumed proportional, of the Rocke-Lorenzato form; also provides the resulting weighted Deming fit and residuals.

Usage

PWD_get_gh(X, Y, lambda = 1, rho=NA, alpha=NA, beta=NA, mu=NA,
           epsilon = 1e-8, printem=FALSE)

Arguments

Details

This workhorse routine optimizes the likelihood in the unknown g, h setting over its n+4 parameters (the two Rocke-Lorenzato precision profile parameters \sigma and \kappa, the intercept \alpha and slope \beta, and the n latent true concentrations \mu_i).

That is, the assumed forms are:

• predicate precision profile model: g_i = var(X_i) = \lambda\left(\sigma^2 + \left[\kappa\cdot \mu_i\right]^2\right) and

• test precision profile model: h_i = var(Y_i) = \sigma^2 + \left[\kappa\cdot (\alpha + \beta\mu_i)\right]^2.

The search algorithm implements an efficient approach via reparameterization to the ratio \rho = \frac{\sigma}{\kappa}.

If initial estimates are not provided, the parameters are initialized as:

• alpha and beta are initially intercept and slope from simple linear regression;

• rho is initialized as the vector c(0.01, 1, 100); and

• mu is initialized as the values of X.

Value

A list containing the following components:

Author(s)

Douglas M. Hawkins, Jessica J. Kraker krakerjj@uwec.edu

References

Hawkins DM and Kraker JJ (in press). Precision Profile Weighted Deming Regression for Methods Comparison. The Journal of Applied Laboratory Medicine. doi:10.1093/jalm/jfaf183

Rocke DM, Lorenzato S (2012). A Two-Component Model for Measurement Error in Analytical Chemistry. Technometrics, 37:2:176-184.

Examples

# library
library(ppwdeming)

# parameter specifications
sigma <- 1
kappa <- 0.08
alpha <- 1
beta  <- 1.1
true  <- 8*10^((0:99)/99)
truey <- alpha+beta*true
# simulate single sample - set seed for reproducibility
set.seed(1039)
# specifications for predicate method
X     <- sigma*rnorm(100)+true *(1+kappa*rnorm(100))
# specifications for test method
Y     <- sigma*rnorm(100)+truey*(1+kappa*rnorm(100))

# fit with RL precision profile to estimate parameters
RL_gh_fit  <- PWD_get_gh(X,Y,printem=TRUE)
# RL precision profile estimated parameters
cat("\nsigmahat=", signif(RL_gh_fit$sigma,6),
    "and kappahat=", signif(RL_gh_fit$kappa,6), "\n")
# with estimated linear coefficients
cat("\nalphahat=", signif(RL_gh_fit$alpha,6),
    "and betahat=", signif(RL_gh_fit$beta,6), "\n")

Weighted Deming Regression – Inference

Description

This routine fits the regression, uses the jackknife to get its precision, and optionally prints it out. Implements Rocke-Lorenzato as the variance profile model.

Usage

PWD_inference(X, Y, lambda=1, rho=NA, alpha=NA, beta=NA, mu=NA, MDL=NA,
              epsilon=1e-8, printem=FALSE)

Arguments

Details

For the linear model relating the predicate and test readings, the standard errors of the estimators \hat{\alpha}, \hat{\beta}, and their covariance are estimated by the jackknife. The point estimates of the intercept and slope are output, along with their standard errors and covariance.

These estimates are further used to estimate the predictions at the input MDL.

Value

A list containing the following components:

Author(s)

Douglas M. Hawkins, Jessica J. Kraker krakerjj@uwec.edu

References

Hawkins DM and Kraker JJ (in press). Precision Profile Weighted Deming Regression for Methods Comparison. The Journal of Applied Laboratory Medicine. doi:10.1093/jalm/jfaf183

Efron, B (1982). The jackknife, the bootstrap and other resampling plans. Society for Industrial and Applied Mathematics.

Examples

# library
library(ppwdeming)

# parameter specifications
sigma <- 1
kappa <- 0.08
alpha <- 1
beta  <- 1.1
true  <- 8*10^((0:99)/99)
truey <- alpha+beta*true
# simulate single sample - set seed for reproducibility
set.seed(1039)
# specifications for predicate method
X     <- sigma*rnorm(100)+true *(1+kappa*rnorm(100))
# specifications for test method
Y     <- sigma*rnorm(100)+truey*(1+kappa*rnorm(100))

# fit with RL precision profile to estimate parameters and variability
RL_inf <- PWD_inference(X,Y,MDL=12,printem=TRUE)

Weighted Deming Regression – general weights

Description

This code is used for the setting of known precision profiles implemented in user-provided R functions called gfun and hfun.

Usage

PWD_known(X, Y, gfun, hfun, gparms, hparms, epsilon=1e-8,
          MDL=NA, getCI=TRUE, printem=FALSE)

Arguments

Details

The functions gfun and hfun are allowed as inputs, to support flexibility in specification of the forms of these variance functions. The known precision profiles specified by the functions gfun and hfun, when provided with estimated vectors of \mu and \alpha + \beta\mu respectively and with any required parameters, will produce the vectors g and h. These vectors are then integrated into the iterative estimation of the slope and intercept of the linear relationship between predicate and test readings.

Value

A list containing the following components:

Author(s)

Douglas M. Hawkins, Jessica J. Kraker krakerjj@uwec.edu

Examples

# library
library(ppwdeming)

# parameter specifications
alpha <- 1
beta  <- 1.1
true  <- 8*10^((0:99)/99)
truey <- alpha+beta*true
# forms of precision profiles
gfun    <- function(true, gparms) {
  gvals = gparms[1]+gparms[2]*true^gparms[3]
  gvals
}
hfun    <- function(true, hparms) {
  hvals = hparms[1]+hparms[2]*true^hparms[3]
  hvals
}

# Loosely motivated by Vitamin D data set
g     <- 4e-16+0.07*true^1.27
h     <- 6e-2+7e-5*truey^2.2
# simulate single sample - set seed for reproducibility
set.seed(1039)
# specifications for predicate method
X     <- true +sqrt(g)*rnorm(100)
# specifications for test method
Y     <- truey+sqrt(h)*rnorm(100)

# fit with to estimate linear parameters
pwd_known_fit <- PWD_known(X, Y, gfun, hfun,
                           gparms=c(4e-16, 0.07, 1.27),
                           hparms=c(6e-2, 7e-5, 2.2), MDL=12,
                           printem=TRUE)

Weighted Deming Regression – Outlier scanning

Description

This function tests for outliers from the fitted regression, and refits on a sanitized data set (with outliers removed).

Usage

PWD_outlier(X, Y, K, lambda=1, Pcut=0.01, rho=NA, alpha=NA, beta=NA, mu=NA,
            printem=FALSE)

Arguments

Details

The method is modeled on the Rosner sequential ESD outlier procedure and assumes the sample is large enough to ignore the effect of random variability in the parameter estimates on the distribution of the residuals.

Value

A list containing the following components:

Author(s)

Douglas M. Hawkins, Jessica J. Kraker krakerjj@uwec.edu

References

Hawkins DM and Kraker JJ (in press). Precision Profile Weighted Deming Regression for Methods Comparison. The Journal of Applied Laboratory Medicine. doi:10.1093/jalm/jfaf183

Hawkins DM (2008). Outliers in Wiley Encyclopedia of Clinical Trials, eds R. D’Agostino, L. Sullivan, and J. Massaro. Wiley, New York.

Examples

# library
library(ppwdeming)

# parameter specifications
sigma <- 1
kappa <- 0.08
alpha <- 1
beta  <- 1.1
true  <- 8*10^((0:99)/99)
truey <- alpha+beta*true
# simulate single sample - set seed for reproducibility
set.seed(1039)
# specifications for predicate method
X     <- sigma*rnorm(100)+true *(1+kappa*rnorm(100))
# specifications for test method
Y     <- sigma*rnorm(100)+truey*(1+kappa*rnorm(100))
# add some outliers
Y[c(1,2,100)] <- Y[c(1,2,100)] + c(7,4,-45)

# check for outliers, re-fit, and store output
outliers_assess <- PWD_outlier(X, Y, K=5, printem=TRUE)

Fit Rocke-Lorenzato profile model to residuals

Description

This routine fits the Rocke-Lorenzato precision profile model to the residuals from the fit.

Usage

PWD_resi(true, resi, epsilon=1e-8, printem=FALSE)

Arguments

Details

The Rocke-Lorenzato precision profile model is

SD^2 = \sigma_r^2 + (\kappa_r\cdot true)^2

for the residuals from a precision-profile model fit.

Under this model, the approach for reviewing residuals is to fit a variance profile model to the residuals r_i themselves. The output of this function includes a maximum-likelihood estimate of the remaining parameter in the special cases of:

• constant variance (\kappa_r = 0); and

• constant coefficient of variation (\sigma_r = 0).

Value

A list containing the following components:

Author(s)

Douglas M. Hawkins, Jessica J. Kraker krakerjj@uwec.edu

References

Hawkins DM and Kraker JJ (in press). Precision Profile Weighted Deming Regression for Methods Comparison. The Journal of Applied Laboratory Medicine. doi:10.1093/jalm/jfaf183

Hawkins DM (2014). A Model for Assay Precision. Statistics in Biopharmaceutical Research, 6, 263-269. http://dx.doi.org/10.1080/19466315.2014.899511

Examples

# library
library(ppwdeming)

# parameter specifications
sigma <- 1
kappa <- 0.08
alpha <- 1
beta  <- 1.1
true  <- 8*10^((0:99)/99)
truey <- alpha+beta*true
# simulate single sample - set seed for reproducibility
set.seed(1039)
# specifications for predicate method
X     <- sigma*rnorm(100)+true *(1+kappa*rnorm(100))
# specifications for test method
Y     <- sigma*rnorm(100)+truey*(1+kappa*rnorm(100))

# fit the model and store output
RL_gh_fit  <- PWD_get_gh(X,Y)
# run the residual analysis from the model output
post  <- PWD_resi(X, RL_gh_fit$resi, printem=TRUE)

Weighted Deming Regression

Description

This code fits the weighted Deming regression on predicate readings (X) and test readings (Y).

Usage

WD_General(X, Y, g, h, epsilon=1e-8)

Arguments

Details

This function is used when the variances of X and Y are already known, as can be the case from an adequate number of replicate X and Y readings of each sample. For input vectors g and h containing the variances of predicate readings X and test readings Y, respectively, iteratively fits weighted Deming regression.

Value

A list containing the following components:

Author(s)

Douglas M. Hawkins, Jessica J. Kraker krakerjj@uwec.edu

References

Ripley BD and Thompson M (1987). Regression techniques for the detection of analytical bias. Analyst, 112, 377-383.

Examples

# library
library(ppwdeming)

# parameter specifications
alpha <- 1
beta  <- 1.1
true  <- 8*10^((0:99)/99)
truey <- alpha+beta*true
# Loosely motivated by Vitamin D data set
g     <- 4e-16+0.07*true^1.27
h     <- 6e-2+7e-5*truey^2.2

# simulate single sample - set seed for reproducibility
set.seed(1039)
# specifications for predicate method
X     <- true +sqrt(g)*rnorm(100)
# specifications for test method
Y     <- truey+sqrt(h)*rnorm(100)

# fit with to estimate linear parameters
wd_fit <- WD_General(X,Y,g,h)
cat("\nWith given g and h, the estimated intercept is",
    signif(wd_fit$alpha,4), "and the estimated slope is",
    signif(wd_fit$beta,4), "\n")

Linnet proportional CV weighted Deming

Description

This routine, provided for convenience, makes Linnet’s constant CV fit.

Usage

WD_Linnet(X, Y, lambda=1, MDL=NA, getCI=TRUE, epsilon=1e-10, printem=FALSE)

Arguments

Details

Note that in cases where sigma happens to come out zero, Linnet’s constant CV fit differs from the precision-profile fit since the underlying precision profile models are not the same.

Value

A list containing the following components:

Author(s)

Douglas M. Hawkins, Jessica J. Kraker krakerjj@uwec.edu

References

Linnet K (1993). Evaluation of regression procedures for methods comparison studies. Clinical Chemistry, 39: 424-432.

Examples

# library
library(ppwdeming)

# parameter specifications
alpha <- 1
beta  <- 1.1
true  <- 8*10^((0:99)/99)
truey <- alpha+beta*true
kappa <- 0.1

# simulate single sample - set seed for reproducibility
set.seed(1039)
# specifications for predicate method
X     <- true *(1+kappa*rnorm(100))
# specifications for test method
Y     <- truey *(1+kappa*rnorm(100))

# fit with to estimate linear parameters
wd_fit <- WD_Linnet(X, Y, MDL=12, printem=TRUE)
cat("\nThe Linnet constant-CV estimated intercept is",
    signif(wd_fit$alpha,4), "and the estimated slope is",
    signif(wd_fit$beta,4), "\n")