aboutsummaryrefslogtreecommitdiff
path: root/R
diff options
context:
space:
mode:
Diffstat (limited to 'R')
-rw-r--r--R/drfit.R88
1 files changed, 44 insertions, 44 deletions
diff --git a/R/drfit.R b/R/drfit.R
index d8631e6..e74667d 100644
--- a/R/drfit.R
+++ b/R/drfit.R
@@ -44,7 +44,7 @@ linlogitf <- function(x,k,f,mu,b)
k*(1 + f*x) / (1 + ((2*f*(10^mu) + 1) * ((x/(10^mu))^b)))
}
-drfit <- function(data, startlogEC50 = NA, chooseone=TRUE,
+drfit <- function(data, startlogED50 = NA, chooseone=TRUE,
probit = TRUE, logit = FALSE, weibull = FALSE,
linlogit = FALSE, linlogitWrong = NA, allWrong = NA,
s0 = 0.5, b0 = 2, f0 = 0)
@@ -62,8 +62,8 @@ drfit <- function(data, startlogEC50 = NA, chooseone=TRUE,
rlhd <- rlld <- vector() # highest and lowest doses tested
mtype <- array() # the modeltypes
sigma <- array() # the standard deviation of the residuals
- logEC50 <- vector()
- stderrlogEC50 <- vector()
+ logED50 <- vector()
+ stderrlogED50 <- vector()
a <- b <- c <- vector()
splitted <- split(data,data$substance)
@@ -89,9 +89,9 @@ drfit <- function(data, startlogEC50 = NA, chooseone=TRUE,
}
rix <- ri
if (!nodata) {
- if (is.na(startlogEC50[i])){
+ if (is.na(startlogED50[i])){
w <- 1/abs(tmp$response - 0.3)
- startlogEC50[[i]] <- sum(w * log10(tmp$dose))/sum(w)
+ startlogED50[[i]] <- sum(w * log10(tmp$dose))/sum(w)
}
highestdose <- max(tmp$dose)
lowestdose <- min(tmp$dose)
@@ -104,9 +104,9 @@ drfit <- function(data, startlogEC50 = NA, chooseone=TRUE,
if (linlogit &&
length(subset(linlogitWrong,linlogitWrong == i))==0 &&
length(subset(allWrong,allWrong == i))==0) {
- m <- try(nls(response ~ linlogitf(dose,1,f,logEC50,b),
+ m <- try(nls(response ~ linlogitf(dose,1,f,logED50,b),
data=tmp,
- start=list(f=f0,logEC50=startlogEC50[[i]],b=b0)))
+ start=list(f=f0,logED50=startlogED50[[i]],b=b0)))
if (!inherits(m, "try-error")) {
fit <- TRUE
ri <- ri + 1
@@ -117,18 +117,18 @@ drfit <- function(data, startlogEC50 = NA, chooseone=TRUE,
runit[[ri]] <- unit
rlld[[ri]] <- log10(lowestdose)
rlhd[[ri]] <- log10(highestdose)
- logEC50[[ri]] <- coef(m)[["logEC50"]]
- if (logEC50[[ri]] > rlhd[[ri]]) {
+ logED50[[ri]] <- coef(m)[["logED50"]]
+ if (logED50[[ri]] > rlhd[[ri]]) {
mtype[[ri]] <- "no fit"
- logEC50[[ri]] <- NA
- stderrlogEC50[[ri]] <- NA
+ logED50[[ri]] <- NA
+ stderrlogED50[[ri]] <- NA
a[[ri]] <- NA
b[[ri]] <- NA
c[[ri]] <- NA
} else {
mtype[[ri]] <- "linlogit"
- stderrlogEC50[[ri]] <- s$parameters["logEC50","Std. Error"]
- a[[ri]] <- coef(m)[["logEC50"]]
+ stderrlogED50[[ri]] <- s$parameters["logED50","Std. Error"]
+ a[[ri]] <- coef(m)[["logED50"]]
b[[ri]] <- coef(m)[["b"]]
c[[ri]] <- coef(m)[["f"]]
}
@@ -137,9 +137,9 @@ drfit <- function(data, startlogEC50 = NA, chooseone=TRUE,
if (probit &&
length(subset(allWrong,allWrong == i))==0) {
- m <- try(nls(response ~ pnorm(-log10(dose),-logEC50,scale),
+ m <- try(nls(response ~ pnorm(-log10(dose),-logED50,scale),
data=tmp,
- start=list(logEC50=startlogEC50[[i]],scale=1)))
+ start=list(logED50=startlogED50[[i]],scale=1)))
if (chooseone==FALSE || fit==FALSE) {
if (!inherits(m, "try-error")) {
fit <- TRUE
@@ -151,18 +151,18 @@ drfit <- function(data, startlogEC50 = NA, chooseone=TRUE,
runit[[ri]] <- unit
rlld[[ri]] <- log10(lowestdose)
rlhd[[ri]] <- log10(highestdose)
- logEC50[[ri]] <- coef(m)[["logEC50"]]
+ logED50[[ri]] <- coef(m)[["logED50"]]
c[[ri]] <- NA
- if (logEC50[[ri]] > rlhd[[ri]]) {
+ if (logED50[[ri]] > rlhd[[ri]]) {
mtype[[ri]] <- "no fit"
- logEC50[[ri]] <- NA
- stderrlogEC50[[ri]] <- NA
+ logED50[[ri]] <- NA
+ stderrlogED50[[ri]] <- NA
a[[ri]] <- NA
b[[ri]] <- NA
} else {
mtype[[ri]] <- "probit"
- stderrlogEC50[[ri]] <- s$parameters["logEC50","Std. Error"]
- a[[ri]] <- coef(m)[["logEC50"]]
+ stderrlogED50[[ri]] <- s$parameters["logED50","Std. Error"]
+ a[[ri]] <- coef(m)[["logED50"]]
b[[ri]] <- coef(m)[["scale"]]
}
}
@@ -171,9 +171,9 @@ drfit <- function(data, startlogEC50 = NA, chooseone=TRUE,
if (logit &&
length(subset(allWrong,allWrong == i))==0) {
- m <- try(nls(response ~ plogis(-log10(dose),-logEC50,scale),
+ m <- try(nls(response ~ plogis(-log10(dose),-logED50,scale),
data=tmp,
- start=list(logEC50=startlogEC50[[i]],scale=1)))
+ start=list(logED50=startlogED50[[i]],scale=1)))
if (chooseone==FALSE || fit==FALSE) {
if (!inherits(m, "try-error")) {
fit <- TRUE
@@ -185,18 +185,18 @@ drfit <- function(data, startlogEC50 = NA, chooseone=TRUE,
runit[[ri]] <- unit
rlld[[ri]] <- log10(lowestdose)
rlhd[[ri]] <- log10(highestdose)
- logEC50[[ri]] <- a[[ri]] <- coef(m)[["logEC50"]]
+ logED50[[ri]] <- a[[ri]] <- coef(m)[["logED50"]]
b[[ri]] <- coef(m)[["scale"]]
c[[ri]] <- NA
- if (logEC50[[ri]] > rlhd[[ri]]) {
+ if (logED50[[ri]] > rlhd[[ri]]) {
mtype[[ri]] <- "no fit"
- logEC50[[ri]] <- NA
- stderrlogEC50[[ri]] <- NA
+ logED50[[ri]] <- NA
+ stderrlogED50[[ri]] <- NA
a[[ri]] <- NA
b[[ri]] <- NA
} else {
mtype[[ri]] <- "logit"
- stderrlogEC50[[ri]] <- s$parameters["logEC50","Std. Error"]
+ stderrlogED50[[ri]] <- s$parameters["logED50","Std. Error"]
}
}
}
@@ -206,7 +206,7 @@ drfit <- function(data, startlogEC50 = NA, chooseone=TRUE,
length(subset(allWrong,allWrong == i))==0) {
m <- try(nls(response ~ pweibull(-log10(dose)+location,shape),
data=tmp,
- start=list(location=startlogEC50[[i]],shape=s0)))
+ start=list(location=startlogED50[[i]],shape=s0)))
if (chooseone==FALSE || fit==FALSE) {
if (!inherits(m, "try-error")) {
fit <- TRUE
@@ -223,17 +223,17 @@ drfit <- function(data, startlogEC50 = NA, chooseone=TRUE,
sqrdev <- function(logdose) {
(0.5 - pweibull( - logdose + a[[ri]], b[[ri]]))^2
}
- logEC50[[ri]] <- nlm(sqrdev,startlogEC50[[i]])$estimate
+ logED50[[ri]] <- nlm(sqrdev,startlogED50[[i]])$estimate
c[[ri]] <- NA
- if (logEC50[[ri]] > rlhd[[ri]]) {
+ if (logED50[[ri]] > rlhd[[ri]]) {
mtype[[ri]] <- "no fit"
- logEC50[[ri]] <- NA
- stderrlogEC50[[ri]] <- NA
+ logED50[[ri]] <- NA
+ stderrlogED50[[ri]] <- NA
a[[ri]] <- NA
b[[ri]] <- NA
} else {
mtype[[ri]] <- "weibull"
- stderrlogEC50[[ri]] <- NA
+ stderrlogED50[[ri]] <- NA
}
}
}
@@ -263,15 +263,15 @@ drfit <- function(data, startlogEC50 = NA, chooseone=TRUE,
}
}
sigma[[ri]] <- NA
- logEC50[[ri]] <- NA
- stderrlogEC50[[ri]] <- NA
+ logED50[[ri]] <- NA
+ stderrlogED50[[ri]] <- NA
a[[ri]] <- NA
b[[ri]] <- NA
c[[ri]] <- NA
}
}
- results <- data.frame(rsubstance, rn, rlld, rlhd, mtype, logEC50, stderrlogEC50, runit, sigma, a, b)
- names(results) <- c("Substance","n","lld","lhd","mtype","logEC50","std","unit","sigma","a","b")
+ results <- data.frame(rsubstance, rn, rlld, rlhd, mtype, logED50, stderrlogED50, runit, sigma, a, b)
+ names(results) <- c("Substance","n","lld","lhd","mtype","logED50","std","unit","sigma","a","b")
if (linlogit) {
results$c <- c
}
@@ -303,8 +303,8 @@ drplot <- function(drresults, data, dtype = "std", alpha = 0.95,
hr <- max(data$response)
dsubstances <- levels(data$substance)
} else {
- lld <- min(drresults[["logEC50"]],na.rm=TRUE) - 2
- lhd <- max(drresults[["logEC50"]],na.rm=TRUE) + 2
+ lld <- min(drresults[["logED50"]],na.rm=TRUE) - 2
+ lhd <- max(drresults[["logED50"]],na.rm=TRUE) + 2
if (length(subset(drresults,mtype=="linlogit")$Substance) != 0) {
hr <- 1.8
} else {
@@ -419,15 +419,15 @@ drplot <- function(drresults, data, dtype = "std", alpha = 0.95,
if (nf > 0) {
for (j in 1:nf)
{
- logEC50 <- fits[j,"logEC50"]
+ logED50 <- fits[j,"logED50"]
mtype <- as.character(fits[j, "mtype"])
if (mtype == "probit") {
scale <- fits[j,"b"]
- plot(function(x) pnorm(-x,-logEC50,scale),lld - 0.5, lhd + 2, add=TRUE,col=color)
+ plot(function(x) pnorm(-x,-logED50,scale),lld - 0.5, lhd + 2, add=TRUE,col=color)
}
if (mtype == "logit") {
scale <- fits[j,"b"]
- plot(function(x) plogis(-x,-logEC50,scale),lld - 0.5, lhd + 2, add=TRUE,col=color)
+ plot(function(x) plogis(-x,-logED50,scale),lld - 0.5, lhd + 2, add=TRUE,col=color)
}
if (mtype == "weibull") {
location <- fits[j,"a"]
@@ -435,7 +435,7 @@ drplot <- function(drresults, data, dtype = "std", alpha = 0.95,
plot(function(x) pweibull(-x+location,shape),lld - 0.5, lhd + 2, add=TRUE,col=color)
}
if (mtype == "linlogit") {
- plot(function(x) linlogitf(10^x,1,fits[j,"c"],fits[j,"logEC50"],fits[j,"b"]),
+ plot(function(x) linlogitf(10^x,1,fits[j,"c"],fits[j,"logED50"],fits[j,"b"]),
lld - 0.5, lhd + 2,
add=TRUE,col=color)
}

Contact - Imprint