- The drplot function now accepts custom xlim and ylim values

- Confidence intervals for the EC50 can now (again) be generated by 
  drfit by giving the argument conf=TRUE
- Update of the INDEX



git-svn-id: http://kriemhild.uft.uni-bremen.de/svn/drfit@62 d1b72e20-2ee0-0310-a1c4-ad5adbbefcdc

1 files changed, 22 insertions, 4 deletions

diff --git a/R/drfit.R b/R/drfit.R
index ff03d90..1bbffa8 100644
--- a/R/drfit.R
+++ b/R/drfit.R
@@ -1,6 +1,7 @@
 drfit <- function(data, startlogED50 = NA, chooseone=TRUE,
         probit = TRUE, logit = FALSE, weibull = FALSE,
-        linlogit = FALSE, linlogitWrong = NA, allWrong = NA,
+        linlogit = FALSE, conf = FALSE,
+        linlogitWrong = NA, allWrong = NA,
         s0 = 0.5, b0 = 2, f0 = 0)
 {
     if(!is.null(data$ok)) data <- subset(data,ok!="no fit") # Don't use data with 
@@ -20,6 +21,7 @@ drfit <- function(data, startlogED50 = NA, chooseone=TRUE,
     sigma <- array()                      # the standard deviation of the residuals
     logED50 <- vector()
     stderrlogED50 <- vector()
+    conflogED50 <- vector()
     a <- b <- c <- vector()
 
     splitted <- split(data,data$substance)
@@ -47,7 +49,7 @@ drfit <- function(data, startlogED50 = NA, chooseone=TRUE,
             n <- ndl <- 0
         } else {
             ndl <- length(levels(factor(tmp$dose)))
-            n <- round(length(tmp$response)/ndl)
+            n <- length(tmp$response)
             if (is.na(startlogED50[i])){
                 w <- 1/abs(tmp$response - 0.3)
                 startlogED50[[i]] <- sum(w * log10(tmp$dose))/sum(w)
@@ -86,12 +88,14 @@ drfit <- function(data, startlogED50 = NA, chooseone=TRUE,
                                 mtype[[ri]] <- "no fit"
                                 logED50[[ri]] <- NA
                                 stderrlogED50[[ri]] <- NA
+                                conflogED50[[ri]] <- NA
                                 a[[ri]] <- NA
                                 b[[ri]] <- NA
                                 c[[ri]] <- NA
                             } else {
                                 mtype[[ri]] <- "linlogit"
                                 stderrlogED50[[ri]] <- s$parameters["logED50","Std. Error"]
+                                conflogED50[[ri]] <- stderrlogED50[[ri]] * qt(0.975, n - 3)
                                 a[[ri]] <- coef(m)[["logED50"]]
                                 b[[ri]] <- coef(m)[["b"]]
                                 c[[ri]] <- coef(m)[["f"]]
@@ -122,11 +126,13 @@ drfit <- function(data, startlogED50 = NA, chooseone=TRUE,
                                     mtype[[ri]] <- "no fit"
                                     logED50[[ri]] <- NA
                                     stderrlogED50[[ri]] <- NA
+                                    conflogED50[[ri]] <- NA
                                     a[[ri]] <- NA
                                     b[[ri]] <- NA
                                 } else {
                                     mtype[[ri]] <- "probit"
                                     stderrlogED50[[ri]] <- s$parameters["logED50","Std. Error"]
+                                    conflogED50[[ri]] <- stderrlogED50[[ri]] * qt(0.975, n - 2)
                                     a[[ri]] <- coef(m)[["logED50"]]
                                     b[[ri]] <- coef(m)[["scale"]]
                                 }
@@ -158,11 +164,13 @@ drfit <- function(data, startlogED50 = NA, chooseone=TRUE,
                                     mtype[[ri]] <- "no fit"
                                     logED50[[ri]] <- NA
                                     stderrlogED50[[ri]] <- NA
+                                    conflogED50[[ri]] <- NA
                                     a[[ri]] <- NA
                                     b[[ri]] <- NA
                                 } else {
                                     mtype[[ri]] <- "logit"
                                     stderrlogED50[[ri]] <- s$parameters["logED50","Std. Error"]
+                                    conflogED50[[ri]] <- stderrlogED50[[ri]] * qt(0.975, n - 2)
                                 }
                             }
                         }
@@ -196,11 +204,13 @@ drfit <- function(data, startlogED50 = NA, chooseone=TRUE,
                                     mtype[[ri]] <- "no fit"
                                     logED50[[ri]] <- NA
                                     stderrlogED50[[ri]] <- NA
+                                    conflogED50[[ri]] <- NA
                                     a[[ri]] <- NA
                                     b[[ri]] <- NA
                                 } else {
                                     mtype[[ri]] <- "weibull"
                                     stderrlogED50[[ri]] <- NA
+                                    conflogED50[[ri]] <- stderrlogED50[[ri]] * qt(0.975, n - 2)
                                 }
                             }
                         }
@@ -238,13 +248,21 @@ drfit <- function(data, startlogED50 = NA, chooseone=TRUE,
             sigma[[ri]] <- NA
             logED50[[ri]] <- NA
             stderrlogED50[[ri]] <- NA
+            conflogED50[[ri]] <- NA
             a[[ri]] <- NA
             b[[ri]] <- NA
             c[[ri]] <- NA
         }
     }
-    results <- data.frame(rsubstance, rndl, rn, rlld, rlhd, mtype, logED50, stderrlogED50, runit, sigma, a, b)
-    names(results) <- c("Substance","ndl","n","lld","lhd","mtype","logED50","std","unit","sigma","a","b")
+    if (conf)
+    {
+        results <- data.frame(rsubstance, rndl, rn, rlld, rlhd, mtype, logED50, conflogED50, runit, sigma, a, b)
+        names(results) <- c("Substance","ndl","n","lld","lhd","mtype","logED50","conf","unit","sigma","a","b")
+    } else {
+        results <- data.frame(rsubstance, rndl, rn, rlld, rlhd, mtype, logED50, stderrlogED50, runit, sigma, a, b)
+        names(results) <- c("Substance","ndl","n","lld","lhd","mtype","logED50","std","unit","sigma","a","b")
+    } 
+
     if (linlogit) {
         results$c <- c
     }