From 51d591316daab73ca0f7ec3cf0b0c5925005d29b Mon Sep 17 00:00:00 2001 From: ranke Date: Fri, 25 Feb 2005 21:43:31 +0000 Subject: Some bugfixes and the addition of the linearlogisWrong argument to the drfit function. git-svn-id: http://kriemhild.uft.uni-bremen.de/svn/drfit@16 d1b72e20-2ee0-0310-a1c4-ad5adbbefcdc --- DESCRIPTION | 4 +- R/drfit.R | 177 +++++++++++++++++++++++++++------------------------- man/drfit.Rd | 10 ++- man/drplot.Rd | 32 +++++----- man/linearlogisf.Rd | 6 +- 5 files changed, 118 insertions(+), 111 deletions(-) diff --git a/DESCRIPTION b/DESCRIPTION index 15ffa95..69f3d1d 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,6 +1,6 @@ Package: drfit -Version: 0.03-13 -Date: 2005-02-24 +Version: 0.03-16 +Date: 2005-02-25 Title: Dose-response data evaluation Author: Johannes Ranke Maintainer: Johannes Ranke diff --git a/R/drfit.R b/R/drfit.R index b328171..9516180 100644 --- a/R/drfit.R +++ b/R/drfit.R @@ -38,7 +38,8 @@ linearlogisf <- function(x,k,f,mu,b) drfit <- function(data, startlogEC50 = NA, chooseone=TRUE, lognorm = TRUE, logis = FALSE, - linearlogis = FALSE, b0 = 2, f0 = 0) + linearlogis = FALSE, linearlogisWrong = NA, + b0 = 2, f0 = 0) { substances <- levels(data$substance) unit <- levels(as.factor(data$unit)) @@ -67,6 +68,7 @@ drfit <- function(data, startlogEC50 = NA, chooseone=TRUE, } else { nodata = FALSE } + rix <- ri if (!nodata) { if (is.na(startlogEC50[i])){ w <- 1/abs(tmp$response - 0.3) @@ -77,11 +79,11 @@ drfit <- function(data, startlogEC50 = NA, chooseone=TRUE, lhd <- log10(highestdose) lld <- log10(lowestdose) responseathighestdose <- mean(subset(tmp,dose==highestdose)$response) - rix <- ri if (responseathighestdose < 0.5) { inactive <- FALSE - if (linearlogis) { + if (linearlogis && + length(subset(linearlogisWrong,linearlogisWrong == i))==0) { m <- try(nls(response ~ linearlogisf(dose,1,f,logEC50,b), data=tmp, start=list(f=f0,logEC50=startlogEC50[[i]],b=b0))) @@ -223,8 +225,7 @@ drplot <- function(drresults, data, dtype = "std", alpha = 0.95, } # Get the plot limits on the x-axis (log of the dose) - if(is.data.frame(data)) - { + if(is.data.frame(data)) { if (min(data$dose == 0)) { cat("At least one of the dose levels is 0 - this is not a valid dose.") } else { @@ -277,94 +278,98 @@ drplot <- function(drresults, data, dtype = "std", alpha = 0.95, for (i in dsubstances) { n <- n + 1 tmp <- splitted[[i]] - color <- colors[[n]] - # Prepare the single graphs if an overlay is not requested - if (!overlay) - { - if (postscript) { - filename = paste(path,fileprefix,sub(" ","_",gsub("([\(\) ])", "", i)),".eps",sep="") - postscript(file=filename, - paper="special",width=7,height=7,horizontal=FALSE,pointsize=12) - cat("Created File: ",filename,"\n") - } - if (png) { - filename = paste(path,fileprefix,sub(" ","_",gsub("([\(\) ])", "", i)),".png",sep="") - png(filename=filename, - width=500, height=500,pointsize=12) - cat("Created File: ",filename,"\n") - } - if (!postscript && !png) { - x11(width=7,height=7,pointsize=12) - } - - plot(0,type="n", - xlim=c(lld - 0.5, lhd + 2), - ylim= c(-0.1, hr + 0.2), - xlab=paste("Decadic Logarithm of the dose in ", unit), - ylab="Normalized response") - } - if (!overlay) legend(lhd - 1, hr + 0.1, i,lty = 1, col = color) - tmp$dosefactor <- factor(tmp$dose) # necessary because the old - # factor has all levels, not - # only the ones tested with - # this substance - # Plot the data, if requested - if (dtype != "none") { - if (dtype == "raw") { - points(log10(tmp$dose),tmp$response,col=color) - } else { - splitresponses <- split(tmp$response,tmp$dosefactor) - means <- sapply(splitresponses,mean) - lengths <- sapply(splitresponses,length) - vars <- sapply(splitresponses,var) - standarddeviations <- sqrt(vars) - } - if (dtype == "std") - { - tops <- means + standarddeviations - bottoms <- means - standarddeviations - } - if (dtype == "conf") + if (length(tmp$response) != 0) { + color <- colors[[n]] + # Prepare the single graphs if an overlay is not requested + if (!overlay) { - confidencedeltas <- qt((1 + alpha)/2, lengths - 1) * sqrt(vars) - tops <- means + confidencedeltas - bottoms <- means - confidencedeltas - } - if (dtype != "raw") - { - x <- log10(as.numeric(levels(tmp$dosefactor))) - segments(x,bottoms,x,tops,col=color) - points(x,means,col=color) - smidge <- 0.05 - segments(x - smidge,bottoms,x + smidge,bottoms,col=color) - segments(x - smidge,tops,x + smidge,tops,col=color) + if (postscript) { + filename = paste(path,fileprefix,sub(" ","_",gsub("([\(\) ])", "", i)),".eps",sep="") + postscript(file=filename, + paper="special",width=7,height=7,horizontal=FALSE,pointsize=12) + cat("Created File: ",filename,"\n") + } + if (png) { + filename = paste(path,fileprefix,sub(" ","_",gsub("([\(\) ])", "", i)),".png",sep="") + png(filename=filename, + width=500, height=500,pointsize=12) + cat("Created File: ",filename,"\n") + } + if (!postscript && !png) { + x11(width=7,height=7,pointsize=12) + } + + plot(0,type="n", + xlim=c(lld - 0.5, lhd + 2), + ylim= c(-0.1, hr + 0.2), + xlab=paste("Decadic Logarithm of the dose in ", unit), + ylab="Normalized response") } - } - - # Plot the fits, if there are any - fits <- subset(drresults,Substance == i) - nf <- length(fits$Substance) # number of fits to plot - if (nf > 0) { - for (j in 1:nf) - { - logEC50 <- fits[j,"logEC50"] - mtype <- as.character(fits[j, "mtype"]) - if (mtype == "lognorm") { - slope <- fits[j,"slope"] - plot(function(x) pnorm(-x,-logEC50,slope),lld - 0.5, lhd + 2, add=TRUE,col=color) + if (!overlay) legend(lhd - 1, hr + 0.1, i,lty = 1, col = color) + tmp$dosefactor <- factor(tmp$dose) # necessary because the old + # factor has all levels, not + # only the ones tested with + # this substance + # Plot the data, if requested + if (dtype != "none") { + if (dtype == "raw") { + points(log10(tmp$dose),tmp$response,col=color) + } else { + splitresponses <- split(tmp$response,tmp$dosefactor) + means <- sapply(splitresponses,mean) + lengths <- sapply(splitresponses,length) + vars <- sapply(splitresponses,var) + standarddeviations <- sqrt(vars) + } + if (dtype == "std") + { + tops <- means + standarddeviations + bottoms <- means - standarddeviations + } + if (dtype == "conf") + { + confidencedeltas <- qt((1 + alpha)/2, lengths - 1) * sqrt(vars) + tops <- means + confidencedeltas + bottoms <- means - confidencedeltas } - if (mtype == "logis") { - slope <- fits[j,"slope"] - plot(function(x) plogis(-x,-logEC50,slope),lld - 0.5, lhd + 2, add=TRUE,col=color) + if (dtype != "raw") + { + x <- log10(as.numeric(levels(tmp$dosefactor))) + segments(x,bottoms,x,tops,col=color) + points(x,means,col=color) + smidge <- 0.05 + segments(x - smidge,bottoms,x + smidge,bottoms,col=color) + segments(x - smidge,tops,x + smidge,tops,col=color) } - if (mtype == "linearlogis") { - plot(function(x) linearlogisf(10^x,1,fits[j,"f"],fits[j,"logEC50"],fits[j,"b"]), - lld - 0.5, lhd + 2, - add=TRUE,col=color) + } + + # Plot the fits, if there are any + fits <- subset(drresults,Substance == i) + nf <- length(fits$Substance) # number of fits to plot + if (nf > 0) { + for (j in 1:nf) + { + logEC50 <- fits[j,"logEC50"] + mtype <- as.character(fits[j, "mtype"]) + if (mtype == "lognorm") { + slope <- fits[j,"slope"] + plot(function(x) pnorm(-x,-logEC50,slope),lld - 0.5, lhd + 2, add=TRUE,col=color) + } + if (mtype == "logis") { + slope <- fits[j,"slope"] + plot(function(x) plogis(-x,-logEC50,slope),lld - 0.5, lhd + 2, add=TRUE,col=color) + } + if (mtype == "linearlogis") { + plot(function(x) linearlogisf(10^x,1,fits[j,"f"],fits[j,"logEC50"],fits[j,"b"]), + lld - 0.5, lhd + 2, + add=TRUE,col=color) + } } } + if (!overlay && (postscript || png)) dev.off() + } else { + cat("No data for ",i,"\n") } - if (!overlay && (postscript || png)) dev.off() } } if (overlay) legend(lhd - 1, hr + 0.1, dsubstances,lty = 1, col = colors) diff --git a/man/drfit.Rd b/man/drfit.Rd index 30f3fd4..ddcf80b 100644 --- a/man/drfit.Rd +++ b/man/drfit.Rd @@ -7,7 +7,7 @@ } \usage{ drfit(data, startlogEC50 = NA, chooseone = TRUE, lognorm = TRUE, logis = FALSE, - linearlogis = FALSE, b0 = 2, f0 = 0) + linearlogis = FALSE, linearlogisWrong = NA, b0 = 2, f0 = 0) } \arguments{ \item{data}{ @@ -30,8 +30,12 @@ A boolean defining if cumulative densitiy curves of logistic distributions are fitted to the data. Default is FALSE.} \item{linearlogis}{ - A boolean defining if the linear-logistic function as defined by van Ewijk and Hoekstra - 1993 is fitted to the data. Default is FALSE.} + A boolean defining if the linear-logistic function + \code{\link{linearlogisf}} as defined by van Ewijk and Hoekstra 1993 is + fitted to the data. Default is FALSE.} + \item{linearlogisWrong}{ + An optional vector containing the names of the substances for which the + linearlogis function produces a wrong fit.} \item{chooseone}{ If TRUE (default), the models are tried in the order linearlogis, logis and lognorm, and the first model that produces a valid fit is used. Usually this will be the one diff --git a/man/drplot.Rd b/man/drplot.Rd index 418de15..784b65b 100644 --- a/man/drplot.Rd +++ b/man/drplot.Rd @@ -7,7 +7,7 @@ } \usage{ drplot(drresults, data, dtype, alpha, path, fileprefix, overlay, - postscript, color, datacolors, fitcolors) + postscript, png, bw, colors) } \arguments{ \item{drresults}{ @@ -46,23 +46,21 @@ will all be put into the last graph in this case. } \item{postscript}{ - If TRUE, (a) postscript graph(s) will be created. Otherwise, graphics will be + If TRUE, (a) postscript graph(s) will be created. Otherwise, and if + the png argument is also FALSE, graphics will be displayed with a screen graphics device. } - \item{color}{ - If TRUE, a sensible selection of colors will be attempted. If false, everything - will be drawn in black + \item{png}{ + If TRUE, (a) png graph(s) will be created. Otherwise, and if the + postscript argument is also FALSE, graphics will be displayed with a + screen graphics device. } - \item{datacolors}{ - This is a vector of colors, defaulting to 1:8, used for plotting the data. + \item{bw}{ + A boolean deciding if the plots will be black and white or not. Default + is TRUE. } - \item{fitcolors}{ - Here you can specify a palette for the colors of the dose-response fits. The - default value is "default", which produces the default palette, if the - number of fits to be plotted is 8 or less. Otherwise, rainbow colors - will be plotted. Unless there is more than one fit per substance to be plotted, - or the number of fits is larger than 8, the fitcolors will match the - datacolors. + \item{colors}{ + This is a vector of colors, defaulting to 1:8, used for plotting the data. } } \value{ @@ -73,9 +71,9 @@ } \note{ - Turn off the colors if you don't like them and don't want to fiddle with - them. Treatment of legends is quite bad. Be sure all devices are closed - (e.g. by calling \code{dev.off()}) before calling \code{drplot} again. + Treatment of legends is not really well solved. Be sure all devices are + closed (e.g. by calling \code{dev.off()}) before calling \code{drplot} + again after a failure. } \examples{ data(antifoul) diff --git a/man/linearlogisf.Rd b/man/linearlogisf.Rd index d49dc9a..1ebe8be 100644 --- a/man/linearlogisf.Rd +++ b/man/linearlogisf.Rd @@ -23,9 +23,9 @@ \value{ The response at dose x. } -\examples{ - -} +\references{ + van Ewijk, P. H. and Hoekstra, J. A. (1993) \emph{Ecotox Environ Safety} + \bold{25} 25-32} \author{ Johannes Ranke \email{jranke@uni-bremen.de} -- cgit v1.2.1