#Power Curve ########################################################
# Independently from b, the curve is 0 for x = 0
# The second form adds a displacement on Y axis,
# so that y != 0 when x = 0
powerCurve.fun <- function(predictor, a, b) {
a * ( predictor ^ b )
}
powerCurveNO.fun <- function(predictor, a, b, c) {
a * ( predictor ^ b ) + c
}
powerCurve.Init <- function(mCall, LHS, data, ...) {
xy <- sortedXyData(mCall[["predictor"]], LHS, data)
lmFit <- lm(log(xy[, "y"]) ~ log(xy[, "x"]))
coefs <- coef(lmFit)
a <- exp(coefs[1])
b <- coefs[2]
value <- c(a, b)
names(value) <- mCall[c("a", "b")]
value
}
NLS.powerCurve <- selfStart(powerCurve.fun, powerCurve.Init, parameters=c("a", "b"))
# powerCurveNO.Init <- function(mCall, LHS, data, ...) {
# xy <- sortedXyData(mCall[["predictor"]], LHS, data)
# pseud
# pseudoY <- log(xy[, "y"])
# pseudoX <- log(xy[, "x"])
# lmFit <- lm(pseudoY ~ pseudoX)
# coefs <- coef(lmFit)
# a <- exp(coefs[1])
# b <- coefs[2]
# value <- c(a, b)
# names(value) <- mCall[c("a", "b")]
# value
# }
#
# NLS.powerCurveNO <- selfStart(powerCurve.fun, powerCurve.Init, parameters=c("a", "b"))
DRC.powerCurve <- function(fixed = c(NA, NA), names = c("a", "b"))
{
## Checking arguments
numParm <- 2
if (!is.character(names) | !(length(names) == numParm)) {stop("Not correct 'names' argument")}
if (!(length(fixed) == numParm)) {stop("Not correct 'fixed' argument")}
## Fixing parameters (using argument 'fixed')
notFixed <- is.na(fixed)
parmVec <- rep(0, numParm)
parmVec[!notFixed] <- fixed[!notFixed]
## Defining the non-linear function
fct <- function(x, parm)
{
parmMat <- matrix(parmVec, nrow(parm), numParm, byrow = TRUE)
parmMat[, notFixed] <- parm
a <- parmMat[, 1]; b <- parmMat[, 2]
a * x ^(b)
}
## Defining self starter function
ssfct <- function(dataf)
{
x <- dataf[, 1]
y <- dataf[, 2]
#regression on pseudo y values
pseudoY <- log( y + 0.00001)
pseudoX <- log(x)
coefs <- coef( lm(pseudoY ~ pseudoX) )
a <- exp(coefs[1])
b <- coefs[2]
return(c(a, b)[notFixed])
}
## Defining names
pnames <- names[notFixed]
## Defining derivatives
deriv1 <- function(x, parms){
parmMat <- matrix(parmVec, nrow(parms),
numParm, byrow = TRUE)
parmMat[, notFixed] <- parms
# Approximation by using finite differences
a <- as.numeric(parmMat[,1])
b <- as.numeric(parmMat[,2])
d1.1 <- expoDecay.fun(x, a, b)
d1.2 <- expoDecay.fun(x, (a + 10e-7), b)
d1 <- (d1.2 - d1.1)/10e-7
d2.1 <- expoDecay.fun(x, a, b)
d2.2 <- expoDecay.fun(x, a, (b + 10e-7) )
d2 <- (d2.2 - d2.1)/10e-7
cbind(d1, d2)[notFixed]
}
## Defining the first derivative (in x=dose)
derivx <- function(x, parm)
{
parmMat <- matrix(parmVec, nrow(parm), numParm, byrow = TRUE)
parmMat[, notFixed] <- parm
a <- as.numeric(parmMat[,1])
b <- as.numeric(parmMat[,2])
d1.1 <- expoGrowth.fun(x, a, b)
d1.2 <- expoGrowth.fun((x + 10e-7), a, b)
d1 <- (d1.2 - d1.1)/10e-7
d1
}
## Defining the ED function
## Defining the inverse function
## Defining descriptive text
text <- "Power curve (Freundlich equation)"
## Returning the function with self starter and names
returnList <- list(fct = fct, ssfct = ssfct, names = pnames, text = text,
noParm = sum(is.na(fixed)),
deriv1 = deriv1, derivx = derivx)
class(returnList) <- "drcMean"
invisible(returnList)
}
"DRC.powerCurveNO" <- function(fixed = c(NA, NA, NA), names = c("a", "b", "c"))
{
## Checking arguments
numParm <- 3
if (!is.character(names) | !(length(names) == numParm)) {stop("Not correct 'names' argument")}
if (!(length(fixed) == numParm)) {stop("Not correct 'fixed' argument")}
## Fixing parameters (using argument 'fixed')
notFixed <- is.na(fixed)
parmVec <- rep(0, numParm)
parmVec[!notFixed] <- fixed[!notFixed]
## Defining the non-linear function
fct <- function(x, parm)
{
parmMat <- matrix(parmVec, nrow(parm), numParm, byrow = TRUE)
parmMat[, notFixed] <- parm
a <- parmMat[, 1]; b <- parmMat[, 2]; c <- parmMat[, 3]
a * x ^(b) + c
}
## Defining self starter function
# ssfct <- function(dataf)
# {
# x <- dataf[, 1]
# y <- dataf[, 2]
#
# #regression on pseudo y values
# pseudoY <- log( y + 0.00001)
# pseudoX <- log(x)
# coefs <- coef( lm(pseudoY ~ pseudoX) )
# a <- exp(coefs[1])
#
# b <- coefs[2]
#
# return(c(a, b)[notFixed])
# }
## Defining names
pnames <- names[notFixed]
## Defining derivatives
## Defining the ED function
## Defining the inverse function
## Defining descriptive text
text <- "Power curve not passing for origin"
## Returning the function with self starter and names
returnList <- list(fct = fct, names = pnames, text = text, noParm = sum(is.na(fixed)))
class(returnList) <- "drcMean"
invisible(returnList)
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.