aboutsummaryrefslogtreecommitdiff
path: root/R
diff options
context:
space:
mode:
authorranke <ranke@d1b72e20-2ee0-0310-a1c4-ad5adbbefcdc>2014-03-18 11:21:20 +0000
committerranke <ranke@d1b72e20-2ee0-0310-a1c4-ad5adbbefcdc>2014-03-18 11:21:20 +0000
commitec7e9cce71d3b46a84be7268da1de4b0953d54c6 (patch)
treeccb734b11acaf2d3cbfc0ac5a117360e25c714bf /R
parent189aad5e39df513075fc7c7a64e646b1a354bc4c (diff)
Improvements on the EDx calculations in drcfit. However, in cases that there is no
hormesis in the data, drfit often fails to fit the linlogit model, and drcfit often fails to find the EDx values. git-svn-id: http://kriemhild.uft.uni-bremen.de/svn/drfit@100 d1b72e20-2ee0-0310-a1c4-ad5adbbefcdc
Diffstat (limited to 'R')
-rw-r--r--R/drcfit.R34
-rw-r--r--R/drfit.R3
-rw-r--r--R/drplot.R5
3 files changed, 25 insertions, 17 deletions
diff --git a/R/drcfit.R b/R/drcfit.R
index ff9e6bb..8c52f6d 100644
--- a/R/drcfit.R
+++ b/R/drcfit.R
@@ -37,6 +37,9 @@ drcfit <- function(data, chooseone=TRUE,
unit <- ""
message("\n",i,": No data\n")
}
+ if (length(unit) == 0) {
+ unit <- NA
+ }
if (length(unit) > 1) {
message("More than one unit for substance ",i,", halting\n\n")
break
@@ -81,23 +84,24 @@ drcfit <- function(data, chooseone=TRUE,
runit[[ri]] <- unit
rlld[[ri]] <- log10(lowestdose)
rlhd[[ri]] <- log10(highestdose)
- ED50 <- try(ED(m, 50, interval = "delta", display = FALSE))
+ logED50[[ri]] <- NA
+ logED50low[[ri]] <- NA
+ logED50high[[ri]] <- NA
+ a[[ri]] <- coef(m)[[2]]
+ b[[ri]] <- coef(m)[[1]]
+ c[[ri]] <- coef(m)[[3]]
+ ED50 <- try(ED(m, 50, interval = "delta",
+ lower = lowestdose / 10,
+ upper = highestdose * 10,
+ display = FALSE))
if (!inherits(ED50, "try-error")) {
logED50[[ri]] <- log10(ED50["1:50", "Estimate"])
logED50low[[ri]] <- log10(ED50["1:50", "Lower"])
logED50high[[ri]] <- log10(ED50["1:50", "Upper"])
if (logED50[[ri]] > rlhd[[ri]]) {
mtype[[ri]] <- "no fit"
- logED50[[ri]] <- NA
- logED50low[[ri]] <- NA
- logED50high[[ri]] <- NA
}
- } else {
- mtype[[ri]] <- "no ED50"
- }
- a[[ri]] <- coef(m)[[2]]
- b[[ri]] <- coef(m)[[1]]
- c[[ri]] <- coef(m)[[3]]
+ }
}
}
}
@@ -119,21 +123,19 @@ drcfit <- function(data, chooseone=TRUE,
rlld[[ri]] <- log10(lowestdose)
rlhd[[ri]] <- log10(highestdose)
logED50[[ri]] <- log10(coef(m)[[2]])
+ a[[ri]] <- coef(m)[[2]]
+ b[[ri]] <- coef(m)[[1]]
c[[ri]] <- NA
if (logED50[[ri]] > rlhd[[ri]]) {
mtype[[ri]] <- "no fit"
logED50[[ri]] <- NA
logED50low[[ri]] <- NA
logED50high[[ri]] <- NA
- a[[ri]] <- NA
- b[[ri]] <- NA
} else {
mtype[[ri]] <- "probit"
ED50 <- ED(m, 50, interval = "delta", display = FALSE)
logED50low[[ri]] <- log10(ED50["1:50", "Lower"])
logED50high[[ri]] <- log10(ED50["1:50", "Upper"])
- a[[ri]] <- coef(m)[[2]]
- b[[ri]] <- coef(m)[[1]]
}
}
}
@@ -155,6 +157,8 @@ drcfit <- function(data, chooseone=TRUE,
rlld[[ri]] <- log10(lowestdose)
rlhd[[ri]] <- log10(highestdose)
logED50[[ri]] <- log10(coef(m)[[2]])
+ a[[ri]] <- coef(m)[[2]]
+ b[[ri]] <- coef(m)[[1]]
c[[ri]] <- NA
if (logED50[[ri]] > rlhd[[ri]]) {
mtype[[ri]] <- "no fit"
@@ -168,8 +172,6 @@ drcfit <- function(data, chooseone=TRUE,
ED50 <- ED(m, 50, interval = "delta", display = FALSE)
logED50low[[ri]] <- log10(ED50["1:50", "Lower"])
logED50high[[ri]] <- log10(ED50["1:50", "Upper"])
- a[[ri]] <- coef(m)[[2]]
- b[[ri]] <- coef(m)[[1]]
}
}
}
diff --git a/R/drfit.R b/R/drfit.R
index ad7d16c..7fa53bf 100644
--- a/R/drfit.R
+++ b/R/drfit.R
@@ -41,6 +41,9 @@ drfit <- function(data, startlogED50 = NA, chooseone=TRUE,
unit <- ""
message("\n",i,": No data\n")
}
+ if (length(unit) == 0) {
+ unit <- NA
+ }
if (length(unit) > 1) {
message("More than one unit for substance ",i,", halting\n\n")
break
diff --git a/R/drplot.R b/R/drplot.R
index ded7aca..11c4d42 100644
--- a/R/drplot.R
+++ b/R/drplot.R
@@ -204,7 +204,10 @@ drplot <- function(drresults, data,
fit.row = fit.rows[j]
if (overlay) nl <- nl + 1 else nl = j
lty <- ltys[nl]
- if (drresults[[fit.row, "mtype"]] %in% c("probit", "logit", "weibull", "linlogit"))
+ if (drresults[[fit.row, "mtype"]] %in% c("probit",
+ "logit",
+ "weibull",
+ "linlogit"))
{
m <- attr(drresults, "models")[[as.numeric(fit.row)]]
of <- function(x) {

Contact - Imprint