From f35bd986bed9d46adefff898e71d4716258334bc Mon Sep 17 00:00:00 2001
From: Walmes Zeviani <walmeszeviani@gmail.com>
Date: Tue, 25 Aug 2015 20:41:40 -0300
Subject: [PATCH] Add biasBox function, Box bias measure for nonlinear
regression.
---
NAMESPACE | 1 +
R/biasBox.R | 95 ++++++++++++++++++++++++++++++++++++++++++++++++++
man/biasBox.Rd | 74 +++++++++++++++++++++++++++++++++++++++
3 files changed, 170 insertions(+)
create mode 100644 R/biasBox.R
create mode 100644 man/biasBox.Rd
diff --git a/NAMESPACE b/NAMESPACE
index 07e3e94..ced143c 100644
--- a/NAMESPACE
+++ b/NAMESPACE
@@ -1,5 +1,6 @@
# Generated by roxygen2 (4.1.1): do not edit by hand
+export(biasBox)
export(polyGui)
export(subsetDropAll)
export(twoStripCombined)
diff --git a/R/biasBox.R b/R/biasBox.R
new file mode 100644
index 0000000..3d31353
--- /dev/null
+++ b/R/biasBox.R
@@ -0,0 +1,95 @@
+#' @title Box bias measaure for nonlinear regression models
+#'
+#' @name biasBox
+#'
+#' @description This function calculates the asymptotic Box bias measure for
+#' nonlinear regression models defined by Box (1971). See the references section.
+#'
+#' @param nls.obj An object of class \code{nls} that has gradient and
+#' hessian attributes. See \link[stats]{deriv3}.
+#'
+#' @return A list with three named elements:
+#' \itemize{
+#' \item \code{Absolute_bias} is the absolute bias.
+#' \item \code{Relative_theta} is the bias in relation to the
+#' pontual parameter estimates, in percentage.
+#' \item \code{Relative_std.error} is the bias in relation to the
+#' precision of the potual parameter estimates, in this case, its
+#' standard error, in percentage.
+#' }
+#'
+#' @author Walmes Zeviani, \email{walmes@@ufpr.br}
+#'
+#' @references Box, M. J. (1971). Bias in nonlinear estimation. Journal
+#' of Royal Statistical Society. Serie B. Methodological, $33 (2), 171–201.
+#'
+#' @export
+#' @examples
+#'
+#' library(lattice)
+#'
+#' data(Puromycin)
+#'
+#' xyplot(rate~conc, groups=state, data=Puromycin)
+#'
+#' da <- subset(Puromycin, state=="treated")
+#'
+#' ##-------------------------------------------
+#' ## Model 1: Michaelis-Menten.
+#'
+#' model1 <- deriv3(expr=~A*conc/(B+conc),
+#' namevec=c("A", "B"),
+#' function.arg=function(conc, A, B){ NULL })
+#'
+#' m1 <- nls(rate~model1(conc, A, B), data=da,
+#' start=list(A=200, B=0.05))
+#' coef(m1)
+#' bb1 <- biasBox(m1)
+#'
+#' ##-------------------------------------------
+#' ## Model 2: monomolecular.
+#'
+#' model2 <- deriv3(expr=~A*(1-exp(-log(2)*conc/B)),
+#' namevec=c("A", "B"),
+#' function.arg=function(conc, A, B){ NULL })
+#'
+#' m2 <- nls(rate~model2(conc, A, B), data=da,
+#' start=list(A=200, B=0.05))
+#' coef(m2)
+#' bb2 <- biasBox(m2)
+#'
+#' ##-------------------------------------------
+#' ## Bias side by side.
+#'
+#' cbind(do.call(rbind, bb1),
+#' do.call(rbind, bb2))
+#'
+#'
+biasBox <- function(nls.obj){
+ smm.obj <- summary(nls.obj)
+ theta <- smm.obj$coef[,1]
+ sd.theta <- smm.obj$coef[,2]
+ sig <- smm.obj$sigma
+ F <- attr(nls.obj$m$fitted(), "gradient")
+ H <- attr(nls.obj$m$fitted(), "hessian")
+ if (is.null(F) | is.null(H)){
+ stop("Models doesn't have a gradient/hessian attributes.")
+ }
+ n <- nrow(F)
+ FtF <- crossprod(F)
+ iFtF <- solve(FtF)
+ d <- -(sig^2/2)*
+ sapply(1:n,
+ function(x){
+ sum(diag(iFtF%*%H[x, , ]))
+ })
+ bias <- as.vector(iFtF%*%t(F)%*%d)
+ names(bias) <- names(coef(nls.obj))
+ bias.sd <- 100*bias/sd.theta
+ bias.th <- 100*bias/theta
+ L <- list("Absolute_bias"=bias,
+ "Relative_theta"=bias.th,
+ "Relative_std.error"=bias.sd)
+ class(L) <- "biasBox"
+ return(L)
+}
diff --git a/man/biasBox.Rd b/man/biasBox.Rd
new file mode 100644
index 0000000..c7ecc49
--- /dev/null
+++ b/man/biasBox.Rd
@@ -0,0 +1,74 @@
+% Generated by roxygen2 (4.1.1): do not edit by hand
+% Please edit documentation in R/biasBox.R
+\name{biasBox}
+\alias{biasBox}
+\title{Box bias measaure for nonlinear regression models}
+\usage{
+biasBox(nls.obj)
+}
+\arguments{
+\item{nls.obj}{An object of class \code{nls} that has gradient and
+hessian attributes. See \link[stats]{deriv3}.}
+}
+\value{
+A list with three named elements:
+\itemize{
+ \item \code{Absolute_bias} is the absolute bias.
+ \item \code{Relative_theta} is the bias in relation to the
+pontual parameter estimates, in percentage.
+ \item \code{Relative_std.error} is the bias in relation to the
+precision of the potual parameter estimates, in this case, its
+standard error, in percentage.
+}
+}
+\description{
+This function calculates the asymptotic Box bias measure for
+nonlinear regression models defined by Box (1971). See the references section.
+}
+\examples{
+library(lattice)
+
+data(Puromycin)
+
+xyplot(rate~conc, groups=state, data=Puromycin)
+
+da <- subset(Puromycin, state=="treated")
+
+##-------------------------------------------
+## Model 1: Michaelis-Menten.
+
+model1 <- deriv3(expr=~A*conc/(B+conc),
+ namevec=c("A", "B"),
+ function.arg=function(conc, A, B){ NULL })
+
+m1 <- nls(rate~model1(conc, A, B), data=da,
+ start=list(A=200, B=0.05))
+coef(m1)
+bb1 <- biasBox(m1)
+
+##-------------------------------------------
+## Model 2: monomolecular.
+
+model2 <- deriv3(expr=~A*(1-exp(-log(2)*conc/B)),
+ namevec=c("A", "B"),
+ function.arg=function(conc, A, B){ NULL })
+
+m2 <- nls(rate~model2(conc, A, B), data=da,
+ start=list(A=200, B=0.05))
+coef(m2)
+bb2 <- biasBox(m2)
+
+##-------------------------------------------
+## Bias side by side.
+
+cbind(do.call(rbind, bb1),
+ do.call(rbind, bb2))
+}
+\author{
+Walmes Zeviani, \email{walmes@ufpr.br}
+}
+\references{
+Box, M. J. (1971). Bias in nonlinear estimation. Journal
+of Royal Statistical Society. Serie B. Methodological, $33 (2), 171–201.
+}
+
--
GitLab