aboutsummaryrefslogtreecommitdiff
path: root/R
diff options
context:
space:
mode:
Diffstat (limited to 'R')
-rw-r--r--R/drfit.R227
1 files changed, 118 insertions, 109 deletions
diff --git a/R/drfit.R b/R/drfit.R
index 9eee87a..1488cd4 100644
--- a/R/drfit.R
+++ b/R/drfit.R
@@ -101,50 +101,19 @@ drfit <- function(data, startlogED50 = NA, chooseone=TRUE,
lhd <- log10(highestdose)
lld <- log10(lowestdose)
responseathighestdose <- mean(subset(tmp,dose==highestdose)$response)
+ responseatlowestdose <- mean(subset(tmp,dose==lowestdose)$response)
if (responseathighestdose < 0.5) {
inactive <- FALSE
-
- if (linlogit &&
- length(subset(linlogitWrong,linlogitWrong == i))==0 &&
- length(subset(allWrong,allWrong == i))==0) {
- m <- try(nls(response ~ linlogitf(dose,1,f,logED50,b),
- data=tmp,
- start=list(f=f0,logED50=startlogED50[[i]],b=b0)))
- if (!inherits(m, "try-error")) {
- fit <- TRUE
- ri <- ri + 1
- s <- summary(m)
- sigma[[ri]] <- s$sigma
- rsubstance[[ri]] <- i
- rndl[[ri]] <- ndl
- rn[[ri]] <- n
- runit[[ri]] <- unit
- rlld[[ri]] <- log10(lowestdose)
- rlhd[[ri]] <- log10(highestdose)
- logED50[[ri]] <- coef(m)[["logED50"]]
- if (logED50[[ri]] > rlhd[[ri]]) {
- mtype[[ri]] <- "no fit"
- logED50[[ri]] <- NA
- stderrlogED50[[ri]] <- NA
- a[[ri]] <- NA
- b[[ri]] <- NA
- c[[ri]] <- NA
- } else {
- mtype[[ri]] <- "linlogit"
- stderrlogED50[[ri]] <- s$parameters["logED50","Std. Error"]
- a[[ri]] <- coef(m)[["logED50"]]
- b[[ri]] <- coef(m)[["b"]]
- c[[ri]] <- coef(m)[["f"]]
- }
- }
- }
-
- if (probit &&
- length(subset(allWrong,allWrong == i))==0) {
- m <- try(nls(response ~ pnorm(-log10(dose),-logED50,scale),
+ if (responseatlowestdose < 0.5) {
+ active <- TRUE
+ } else {
+ active <- FALSE
+ if (linlogit &&
+ length(subset(linlogitWrong,linlogitWrong == i))==0 &&
+ length(subset(allWrong,allWrong == i))==0) {
+ m <- try(nls(response ~ linlogitf(dose,1,f,logED50,b),
data=tmp,
- start=list(logED50=startlogED50[[i]],scale=1)))
- if (chooseone==FALSE || fit==FALSE) {
+ start=list(f=f0,logED50=startlogED50[[i]],b=b0)))
if (!inherits(m, "try-error")) {
fit <- TRUE
ri <- ri + 1
@@ -157,95 +126,131 @@ drfit <- function(data, startlogED50 = NA, chooseone=TRUE,
rlld[[ri]] <- log10(lowestdose)
rlhd[[ri]] <- log10(highestdose)
logED50[[ri]] <- coef(m)[["logED50"]]
- c[[ri]] <- NA
if (logED50[[ri]] > rlhd[[ri]]) {
mtype[[ri]] <- "no fit"
logED50[[ri]] <- NA
stderrlogED50[[ri]] <- NA
a[[ri]] <- NA
b[[ri]] <- NA
+ c[[ri]] <- NA
} else {
- mtype[[ri]] <- "probit"
+ mtype[[ri]] <- "linlogit"
stderrlogED50[[ri]] <- s$parameters["logED50","Std. Error"]
a[[ri]] <- coef(m)[["logED50"]]
- b[[ri]] <- coef(m)[["scale"]]
+ b[[ri]] <- coef(m)[["b"]]
+ c[[ri]] <- coef(m)[["f"]]
}
}
}
- }
- if (logit &&
- length(subset(allWrong,allWrong == i))==0) {
- m <- try(nls(response ~ plogis(-log10(dose),-logED50,scale),
- data=tmp,
- start=list(logED50=startlogED50[[i]],scale=1)))
- if (chooseone==FALSE || fit==FALSE) {
- if (!inherits(m, "try-error")) {
- fit <- TRUE
- ri <- ri + 1
- s <- summary(m)
- sigma[[ri]] <- s$sigma
- rsubstance[[ri]] <- i
- rndl[[ri]] <- ndl
- rn[[ri]] <- n
- runit[[ri]] <- unit
- rlld[[ri]] <- log10(lowestdose)
- rlhd[[ri]] <- log10(highestdose)
- logED50[[ri]] <- a[[ri]] <- coef(m)[["logED50"]]
- b[[ri]] <- coef(m)[["scale"]]
- c[[ri]] <- NA
- if (logED50[[ri]] > rlhd[[ri]]) {
- mtype[[ri]] <- "no fit"
- logED50[[ri]] <- NA
- stderrlogED50[[ri]] <- NA
- a[[ri]] <- NA
- b[[ri]] <- NA
- } else {
- mtype[[ri]] <- "logit"
- stderrlogED50[[ri]] <- s$parameters["logED50","Std. Error"]
+ if (probit &&
+ length(subset(allWrong,allWrong == i))==0) {
+ m <- try(nls(response ~ pnorm(-log10(dose),-logED50,scale),
+ data=tmp,
+ start=list(logED50=startlogED50[[i]],scale=1)))
+ if (chooseone==FALSE || fit==FALSE) {
+ if (!inherits(m, "try-error")) {
+ fit <- TRUE
+ ri <- ri + 1
+ s <- summary(m)
+ sigma[[ri]] <- s$sigma
+ rsubstance[[ri]] <- i
+ rndl[[ri]] <- ndl
+ rn[[ri]] <- n
+ runit[[ri]] <- unit
+ rlld[[ri]] <- log10(lowestdose)
+ rlhd[[ri]] <- log10(highestdose)
+ logED50[[ri]] <- coef(m)[["logED50"]]
+ c[[ri]] <- NA
+ if (logED50[[ri]] > rlhd[[ri]]) {
+ mtype[[ri]] <- "no fit"
+ logED50[[ri]] <- NA
+ stderrlogED50[[ri]] <- NA
+ a[[ri]] <- NA
+ b[[ri]] <- NA
+ } else {
+ mtype[[ri]] <- "probit"
+ stderrlogED50[[ri]] <- s$parameters["logED50","Std. Error"]
+ a[[ri]] <- coef(m)[["logED50"]]
+ b[[ri]] <- coef(m)[["scale"]]
+ }
}
}
}
- }
- if (weibull &&
- length(subset(allWrong,allWrong == i))==0) {
- m <- try(nls(response ~ pweibull(-log10(dose)+location,shape),
- data=tmp,
- start=list(location=startlogED50[[i]],shape=s0)))
- if (chooseone==FALSE || fit==FALSE) {
- if (!inherits(m, "try-error")) {
- fit <- TRUE
- ri <- ri + 1
- s <- summary(m)
- sigma[[ri]] <- s$sigma
- rsubstance[[ri]] <- i
- rndl[[ri]] <- ndl
- rn[[ri]] <- n
- runit[[ri]] <- unit
- rlld[[ri]] <- log10(lowestdose)
- rlhd[[ri]] <- log10(highestdose)
- a[[ri]] <- coef(m)[["location"]]
- b[[ri]] <- coef(m)[["shape"]]
- sqrdev <- function(logdose) {
- (0.5 - pweibull( - logdose + a[[ri]], b[[ri]]))^2
+ if (logit &&
+ length(subset(allWrong,allWrong == i))==0) {
+ m <- try(nls(response ~ plogis(-log10(dose),-logED50,scale),
+ data=tmp,
+ start=list(logED50=startlogED50[[i]],scale=1)))
+ if (chooseone==FALSE || fit==FALSE) {
+ if (!inherits(m, "try-error")) {
+ fit <- TRUE
+ ri <- ri + 1
+ s <- summary(m)
+ sigma[[ri]] <- s$sigma
+ rsubstance[[ri]] <- i
+ rndl[[ri]] <- ndl
+ rn[[ri]] <- n
+ runit[[ri]] <- unit
+ rlld[[ri]] <- log10(lowestdose)
+ rlhd[[ri]] <- log10(highestdose)
+ logED50[[ri]] <- a[[ri]] <- coef(m)[["logED50"]]
+ b[[ri]] <- coef(m)[["scale"]]
+ c[[ri]] <- NA
+ if (logED50[[ri]] > rlhd[[ri]]) {
+ mtype[[ri]] <- "no fit"
+ logED50[[ri]] <- NA
+ stderrlogED50[[ri]] <- NA
+ a[[ri]] <- NA
+ b[[ri]] <- NA
+ } else {
+ mtype[[ri]] <- "logit"
+ stderrlogED50[[ri]] <- s$parameters["logED50","Std. Error"]
+ }
}
- logED50[[ri]] <- nlm(sqrdev,startlogED50[[i]])$estimate
- c[[ri]] <- NA
- if (logED50[[ri]] > rlhd[[ri]]) {
- mtype[[ri]] <- "no fit"
- logED50[[ri]] <- NA
- stderrlogED50[[ri]] <- NA
- a[[ri]] <- NA
- b[[ri]] <- NA
- } else {
- mtype[[ri]] <- "weibull"
- stderrlogED50[[ri]] <- NA
+ }
+ }
+
+ if (weibull &&
+ length(subset(allWrong,allWrong == i))==0) {
+ m <- try(nls(response ~ pweibull(-log10(dose)+location,shape),
+ data=tmp,
+ start=list(location=startlogED50[[i]],shape=s0)))
+ if (chooseone==FALSE || fit==FALSE) {
+ if (!inherits(m, "try-error")) {
+ fit <- TRUE
+ ri <- ri + 1
+ s <- summary(m)
+ sigma[[ri]] <- s$sigma
+ rsubstance[[ri]] <- i
+ rndl[[ri]] <- ndl
+ rn[[ri]] <- n
+ runit[[ri]] <- unit
+ rlld[[ri]] <- log10(lowestdose)
+ rlhd[[ri]] <- log10(highestdose)
+ a[[ri]] <- coef(m)[["location"]]
+ b[[ri]] <- coef(m)[["shape"]]
+ sqrdev <- function(logdose) {
+ (0.5 - pweibull( - logdose + a[[ri]], b[[ri]]))^2
+ }
+ logED50[[ri]] <- nlm(sqrdev,startlogED50[[i]])$estimate
+ c[[ri]] <- NA
+ if (logED50[[ri]] > rlhd[[ri]]) {
+ mtype[[ri]] <- "no fit"
+ logED50[[ri]] <- NA
+ stderrlogED50[[ri]] <- NA
+ a[[ri]] <- NA
+ b[[ri]] <- NA
+ } else {
+ mtype[[ri]] <- "weibull"
+ stderrlogED50[[ri]] <- NA
+ }
}
}
}
- }
+ }
} else {
inactive <- TRUE
@@ -267,7 +272,11 @@ drfit <- function(data, startlogED50 = NA, chooseone=TRUE,
if (inactive) {
mtype[[ri]] <- "inactive"
} else {
- mtype[[ri]] <- "no fit"
+ if (active) {
+ mtype[[ri]] <- "active"
+ } else {
+ mtype[[ri]] <- "no fit"
+ }
}
}
sigma[[ri]] <- NA

Contact - Imprint