Skip to content
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 1 addition & 1 deletion .travis.yml
Original file line number Diff line number Diff line change
Expand Up @@ -9,7 +9,7 @@ r:
- devel

r_github_packages:
- DoseResponse/drcData
# - DoseResponse/drcData
- r-lib/covr

after_success:
Expand Down
8 changes: 5 additions & 3 deletions DESCRIPTION
Original file line number Diff line number Diff line change
@@ -1,11 +1,12 @@
Package: drc
Version: 3.2-0
Date: 2021-01-13
Date: 2023-04-25
Title: Analysis of Dose-Response Data
Authors@R: c(person("Christian" , "Ritz", role = c("aut", "cre"), email="ritz@bioassay.dk"),
person("Jens", "Streibig", middle="C.", email="streibig@bioassay.dk", role = "aut"))
person("Jens", "Streibig", middle="C.", email="streibig@bioassay.dk", role = "aut"),
person("Ryan", "Ward", role = "ctb"))
Maintainer: Christian Ritz <ritz@bioassay.dk>
Depends: R (>= 2.0.0), MASS, stats, drcData
Depends: R (>= 2.0.0), MASS, stats
Imports: car, gtools, multcomp, plotrix, sandwich, scales
LazyLoad: yes
LazyData: yes
Expand All @@ -14,3 +15,4 @@ License: GPL-2
URL: http://www.r-project.org, http://www.bioassay.dk
RoxygenNote: 6.1.1
BugReports: https://github.com/DoseResponse/drc/issues/

2 changes: 1 addition & 1 deletion NAMESPACE
Original file line number Diff line number Diff line change
Expand Up @@ -39,7 +39,7 @@ importFrom(scales, alpha)
#importFrom(stats, cooks.distance, hatvalues)
#exportMethods(BIC)

import(drcData)
# import(drcData)
importFrom("graphics", "abline", "axTicks", "axis", "lines", "par",
"plot", "points", "polygon", "segments", "text")

Expand Down
13 changes: 11 additions & 2 deletions R/drmOpt.R
Original file line number Diff line number Diff line change
@@ -1,3 +1,12 @@
generateOptimControl <- function(maxIt, relTol, parscale, traceVal, method) {
if (method == "L-BFGS-B") {
control <- list(maxit = maxIt, factr = relTol, parscale = parscale, trace = traceVal)
} else {
control <- list(maxit = maxIt, reltol = relTol, parscale = parscale, trace = traceVal)
}
return(control)
}

"drmOpt" <-
function(opfct, opdfct1, startVec, optMethod, constrained, warnVal,
upperLimits, lowerLimits, errorMessage, maxIt, relTol, opdfct2 = NULL, parmVec, traceVal, silentVal = TRUE,
Expand All @@ -21,7 +30,7 @@ matchCall)
{
nlsObj <- try(optim(startVec, opfct, opdfct1, hessian = hes, method = "L-BFGS-B",
lower = lowerLimits, upper = upperLimits,
control = list(maxit = maxIt, reltol = relTol, parscale = psVec)), silent = silentVal)
control = generateOptimControl(maxIt, relTol, psVec, traceVal, "L-BFGS-B")), silent = silentVal)
} else {
nlsObj <- try(optim(startVec, opfct, opdfct1, hessian = hes, method = optMethod,
control = list(maxit = maxIt, reltol = relTol, parscale = psVec)), silent = silentVal)
Expand Down Expand Up @@ -54,7 +63,7 @@ matchCall)
# print(opfct(startVec))
nlsObj <- try(optim(startVec, opfct, hessian = TRUE, method = "L-BFGS-B",
lower = lowerLimits, upper = upperLimits,
control = list(maxit = maxIt, parscale = psVec, reltol = relTol, trace = traceVal)), silent = silentVal)
control = generateOptimControl(maxIt, relTol, psVec, traceVal, "L-BFGS-B")), silent = silentVal)
# parscale is needed for the example in methionine.Rd
} else {
# psVec <- abs(startVec)
Expand Down
148 changes: 73 additions & 75 deletions R/predict.drc.R
Original file line number Diff line number Diff line change
@@ -1,18 +1,18 @@
"predict.drc" <- function(object, newdata, se.fit = FALSE,
interval = c("none", "confidence", "prediction", "ssd"),
level = 0.95, na.action = na.pass, od = FALSE, vcov. = vcov,
"predict.drc" <- function(object, newdata, se.fit = FALSE,
interval = c("none", "confidence", "prediction", "ssd"),
level = 0.95, na.action = na.pass, od = FALSE, vcov. = vcov,
ssdSEfct = NULL, constrain = TRUE, checkND = TRUE, ...)
{
## Checking arguments
interval <- match.arg(interval)
respType <- object[["type"]]

dataList <- object[["dataList"]]
dataList <- object[["dataList"]]
doseDim <- ncol(dataList[["dose"]])
if (is.null(doseDim)) {doseDim <- 1}
if (is.null(doseDim)) {doseDim <- 1}

## Assigning dataset from object if no data frame is provided
if (missing(newdata))
if (missing(newdata))
{
# predValues <- fitted(object) # not used
# newdata <- data.frame(object$data[, 1], object$data[, 3])
Expand All @@ -25,56 +25,56 @@
groupLevels <- as.character(dataList[["plotid"]])
} else {
groupLevels <- as.character(dataList[["curveid"]])
}
#
}
#
# if (identical(respType, "event"))
# {
# newdata <- data.frame(dataList[["dose"]], dataList[["plotid"]])
# } else {
# newdata <- data.frame(dataList[["dose"]], dataList[["curveid"]])
# }
} else {

if (checkND)
{
{
dName <- dataList[["names"]][["dName"]]
if (any(names(newdata) %in% dName))
{
doseVec <- newdata[, dName]
doseVec <- newdata[, dName]
} else {
doseVec <- newdata[, 1]
# warning("Dose variable not in 'newdata'")
}
} else {
doseVec <- newdata
}

cName <- dataList[["names"]][["cNames"]]
if (any(names(newdata) %in% cName))
{
groupLevels <- as.character(newdata[, cName])
# as.character() removes factor encoding
groupLevels <- as.character(newdata[, cName])
# as.character() removes factor encoding

} else {
groupLevels <- rep(1, nrow(newdata))
}
#
#
#
#
# if (ncol(newdata) < (doseDim + 1)) {newdata <- data.frame(newdata, rep(1, nrow(newdata)))}
# # ndncol <- ncol(newdata)
# # doseVec <- newdata[, 1:(ndncol-1)]
# doseVec <- newdata[, 1:doseDim]
# # groupLevels <- as.character(newdata[, ndncol]) # 'as.character()' used to suppress factor levels
# groupLevels <- as.character(newdata[, doseDim + 1]) # 'as.character()' used to suppress factor levels
# # groupLevels <- as.character(newdata[, ndncol]) # 'as.character()' used to suppress factor levels
# groupLevels <- as.character(newdata[, doseDim + 1]) # 'as.character()' used to suppress factor levels
}
noNewData <- length(groupLevels)

# if (ncol(newdata) < 2) {newdata <- data.frame(newdata, rep(1, nrow(newdata)))}
# if (ncol(newdata) > 2) {stop("More than 2 variables in 'newdata' argument")}
## Defining dose values -- dose in the first column!

## Defining dose values -- dose in the first column!
# doseVec <- newdata[, 1]
# groupLevels <- as.character(newdata[, 2]) # 'as.character()' used to suppress factor levels
# groupLevels <- as.character(newdata[, 2]) # 'as.character()' used to suppress factor levels
# noNewData <- length(doseVec)


Expand All @@ -86,56 +86,56 @@
}

## Retrieving matrix of parameter estimates
parmMat <- object[["parmMat"]]
parmMat <- object[["parmMat"]]
pm <- t(parmMat[, groupLevels, drop = FALSE])

# parmNames <- colnames(parmMat)
# lenCN <- length(parmNames)
# indVec <- 1:lenCN
# names(indVec) <- parmNames
# if (lenCN > 1)
# {
# indVec <- indVec[as.character(newdata[, 2])]
#
#
## groupLevels <- newdata[, 2]
# if (!all(is.numeric(groupLevels)))
# {
## pm <- parmMat[, as.character(groupLevels)] # 'as.character()' used to suppress factor levels
# pm <- parmMat[, groupLevels]
## pm <- parmMat[, as.character(groupLevels)] # 'as.character()' used to suppress factor levels
# pm <- parmMat[, groupLevels]
# } else {
# pm <- parmMat[, groupLevels]
# }
# pm <- parmMat[, groupLevels]
#
#
# } else {
# lenDV <- length(doseVec)
## indVec <- rep(1, lenDV)
# pm <- matrix(parmMat[, 1], length(parmMat[, 1]), lenDV)
# }
# }

# ## Checking for NAs in matrix of parameter estimates
# naVec <- rep(FALSE, lenCN)
# for (i in 1:lenCN)
# {
# naVec[i] <- any(is.na(parmMat[, i]))
# }
# parmMat <- parmMat[, !naVec, drop = FALSE]
# parmMat <- parmMat[, !naVec, drop = FALSE]


## Retrieving variance-covariance matrix
sumObj <- summary(object, od = od)
# varMat <- sumObj[["varMat"]]
vcovMat <- vcov.(object)
# varMat <- sumObj[["varMat"]]
vcovMat <- vcov.(object)

## Defining index matrix for parameter estimates
indexMat <- object[["indexMat"]]
## Calculating predicted values
# indexVec <- as.vector(indVec)
# print(indexVec)
# lenIV <- length(indexVec)

## Calculating predicted values
# indexVec <- as.vector(indVec)
# print(indexVec)
# lenIV <- length(indexVec)


# retMat <- matrix(0, lenIV, 4)
retMat <- matrix(0, noNewData, 4)
colnames(retMat) <- c("Prediction", "SE", "Lower", "Upper")
Expand All @@ -144,25 +144,25 @@
# print(doseVec)
retMat[, 1] <- objFct$"fct"(doseVec, pm)
# print(pm)

## Checking if derivatives are available
deriv1 <- objFct$"deriv1"
if (is.null(deriv1))
{
return(retMat[, 1])
}
return(retMat[, 1])
}

## Calculating the quantile to be used in the confidence intervals
if (!identical(interval, "none"))
{
{
if (identical(respType, "continuous"))
{
tquan <- qt(1 - (1 - level)/2, df.residual(object))
tquan <- as.vector(qt(1 - (1 - level)/2, df.residual(object)))
} else {
tquan <- qnorm(1 - (1 - level)/2)
tquan <- as.vector(qnorm(1 - (1 - level)/2))
}
}
}

## Calculating standard errors and/or confidence intervals
if (se.fit || (!identical(interval, "none")))
{
Expand All @@ -171,8 +171,8 @@
{
estVec <- object[["dataList"]][["dose"]]
seVec <- object[["dataList"]][["weights"]]
if (is.null(ssdSEfct))

if (is.null(ssdSEfct))
{
# lmObj <- lm(seVec ~ estVec) # linear, not great
lmObj <- lm(log(seVec) ~ log(estVec))
Expand All @@ -186,24 +186,24 @@
# print(derivxRes)
# if (is.finite(derivxRes))
# {
# sumObjRV <- (derivxRes * sePred)^2
# sumObjRV <- (derivxRes * sePred)^2
# } else {
# sumObjRV <- 0 # setting Inf * 0 = 0 (would be NaN otherwise)
# }
# sumObjRV <- 0 # setting Inf * 0 = 0 (would be NaN otherwise)
# }
sumObjRV <- rep(0, length(derivxRes))
isFinDR <- is.finite(derivxRes)
isFinDR <- is.finite(derivxRes)
sumObjRV[isFinDR] <- ((derivxRes * sePred)^2)[isFinDR]

# print(sumObjRV)
}
}
if (identical(interval, "prediction"))
{
sumObjRV <- rep(sumObj$"resVar", noNewData)
}
}
#else {
# sumObjRV <- 0
# }
# rowIndex <- 1
# rowIndex <- 1
# for (i in indexVec)
# for (i in 1:ncol(indexMat))

Expand All @@ -214,17 +214,17 @@
for (rowIndex in 1:noNewData)
{
# parmInd <- indexMat[, i]
# print(indexVec)
# print(indexVec)
# print(varMat)
# print(parmInd)
# print(parmInd)

# varCov <- varMat[parmInd, parmInd]
# print(varCov)
# groupLevels <- newdata[, 2]
# parmInd <- indexMat[, groupLevels[rowIndex]]
# varCov <- varMat[parmInd, parmInd]

parmInd <- piMat[, rowIndex]
parmInd <- piMat[, rowIndex]
varCov <- vcovMat[parmInd, parmInd]

# parmChosen <- t(parmMat[, i, drop = FALSE])
Expand All @@ -233,28 +233,28 @@

dfEval <- deriv1(doseVec[rowIndex], pm[rowIndex, , drop = FALSE])
varVal <- dfEval %*% varCov %*% dfEval
retMat[rowIndex, 2] <- sqrt(varVal)
# retMat[rowIndex, 2] <- sqrt(dfEval %*% varCov %*% dfEval)
retMat[rowIndex, 2] <- sqrt(varVal)
# retMat[rowIndex, 2] <- sqrt(dfEval %*% varCov %*% dfEval)

if (!se.fit)
{
#retMat[rowIndex, 3:4] <- rep(retMat[rowIndex, 1], 2) +
#retMat[rowIndex, 3:4] <- rep(retMat[rowIndex, 1], 2) +
# (tquan * sqrt(varVal + sumObjRV[rowIndex])) * c(-1, 1)
retMat[rowIndex, 3] <- retMat[rowIndex, 1] - tquan * sqrt(varVal + sumObjRV[rowIndex])
retMat[rowIndex, 4] <- retMat[rowIndex, 1] + tquan * sqrt(varVal + sumObjRV[rowIndex])
}
retMat[rowIndex, 4] <- retMat[rowIndex, 1] + tquan * sqrt(varVal + sumObjRV[rowIndex])
}
# if (identical(interval, "confidence"))
# {
# retMat[rowIndex, 3] <- retMat[rowIndex, 1] - tquan * sqrt(varVal)
# retMat[rowIndex, 4] <- retMat[rowIndex, 1] + tquan * sqrt(varVal)
# retMat[rowIndex, 4] <- retMat[rowIndex, 1] + tquan * sqrt(varVal)
# }
# if (identical(interval, "prediction"))
# {
# sumObjRV <- sumObj$"resVar"
# retMat[rowIndex, 3] <- retMat[rowIndex, 1] - tquan * sqrt(varVal + sumObjRV)
# retMat[rowIndex, 4] <- retMat[rowIndex, 1] + tquan * sqrt(varVal + sumObjRV)
# }
# rowIndex <- rowIndex + 1
# retMat[rowIndex, 4] <- retMat[rowIndex, 1] + tquan * sqrt(varVal + sumObjRV)
# }
# rowIndex <- rowIndex + 1
}
}
## Imposing constraints on predicted values
Expand All @@ -270,13 +270,11 @@
retMat[, 3] <- pmax(retMat[, 3], 0)
}
}

## Keeping relevant indices
keepInd <- 1
if (se.fit) {keepInd <- c(keepInd, 2)}
if (!identical(interval, "none")) {keepInd <- c(keepInd, 3, 4)}

return(retMat[, keepInd]) # , drop = FALSE])
}


Loading