"anova.drclist" <- function(object, ..., details = TRUE, test = NULL)
{
objects <- list(object, ...)
if (length(objects) > 2) {stop("Only two models can be compared")}
# if (inherits(object, "bindrc")) # the argument 'test="F"' is not used
# {
# obj1 <- objects[[1]]$loglik
# obj2 <- objects[[2]]$loglik
# rowNames <- c("Model 1", "Model 2")
#
#
# loglik <- c(obj1[3], obj2[3])
# dfModel <- c(obj1[4], obj2[4])
# testStat <- (2*abs(loglik[1] - loglik[2]))
# dfDiff <- c(NA, abs(dfModel[1] - dfModel[2]))
#
# pVal <- c(NA, 1 - pchisq(testStat, dfDiff[2]))
# testStat <- c(NA, testStat)
#
# headName <- "Analysis of deviance table\n"
# colNames <- c("ModelDf", "Loglik", "Df", "LR value", "p value")
#
# } else {
## Testing two models against each other
obj1 <- objects[[1]]
obj2 <- objects[[2]]
rowNames <- c("1st model", "2nd model")
# sumObj1 <- summary(obj1)
# sumObj2 <- summary(obj2)
if ( !(obj1$"type"==obj2$"type") ) {stop("The two models are based on different types on data")}
if (obj1$"type" == "binomial" && (is.null(test)) ) {test <- "Chisq"}
if (obj1$"type" == "continuous" && (is.null(test)) ) {test <- "F"}
if (!(test == "F")) # chis-square based test
{
loglik <- c(logLik(obj1), logLik(obj2)) # c(sumObj1[[4]][1], sumObj2[[4]][1])
dfModel <- c(attr(logLik(obj1), "df"), attr(logLik(obj2), "df")) # c(sumObj1[[4]][2], sumObj2[[4]][2])
testStat <- (2*abs(loglik[1] - loglik[2]))
dfDiff <- c(NA, abs(dfModel[1] - dfModel[2]))
pVal <- c(NA, 1 - pchisq(testStat, dfDiff[2]))
testStat <- c(NA, testStat)
headName <- "ANOVA-like table\n"
colNames <- c("ModelDf", "Loglik", "Df", "LR value", "p value")
} else { # F-test
# sumVec1 <- obj1$"summary"
# sumVec2 <- obj2$"summary"
df1 <- df.residual(obj1)
df2 <- df.residual(obj2)
if (df2 > df1)
{
objTemp <- obj1
obj1 <- obj2
obj2 <- objTemp
df1 <- df.residual(obj1)
df2 <- df.residual(obj2)
rowNames <- rowNames[c(2, 1)]
}
# if (sumVec2[6]>sumVec1[6])
# {sumTemp <- sumVec1; sumVec1 <- sumVec2; sumVec2 <- sumTemp; rowNames <- rowNames[c(2,1)]}
# loglik <- c(sumVec1[5],sumVec2[5])
loglik <- c(obj1$"summary"[4], obj2$"summary"[4]) # use an extractor "rss()" instead?
# dfModel <- c(sumVec1[6], sumVec2[6])
dfModel <- c(df1, df2)
# loglikTemp <- c(sumVec1[5], sumVec2[5])
# loglik <- c((loglikTemp[1]-loglikTemp[2])/(dfModel[1]-dfModel[2]), loglikTemp[2]/dfModel[2])
# dfModel <- c(sumVec1[6], sumVec2[6])
# dfDiff <- c((loglik[1]-loglik[2])/(dfModel[1]-dfModel[2]), loglik[2]/dfModel[2])
# dfDiff <- c(NA, dfModel[1] - dfModel[2])
dfDiff <- c(NA, df1 - df2)
# testStat <- dfDiff[1]/dfDiff[2]
testStat <- ((loglik[1] - loglik[2]) / dfDiff[2]) / (loglik[2] / df2)
pVal <- c(NA, 1 - pf(testStat, dfDiff[2], df2))
testStat <- c(NA, testStat)
headName <- "ANOVA table\n"
colNames <- c("ModelDf", "RSS", "Df", "F value", "p value")
}
# dataFra <- data.frame(dfModel, loglik, dfDiff, testStat, pVal)
if (details)
{
## Specifying the models
cat("\n")
collapse1 <- obj1[[8]]$collapse
if (is.null(collapse1)) {collapse1 <- obj1[[8]]$pmodels}
if (!is.null(obj1$"pmodelsText")) {collapse1 <- obj1$"pmodelsText"}
# print(collapse1)
# if (is.null(collapse1)) {collapse1 <- paste(deparse(obj1[[8]]$assayNo), "(for all parameters)")} else {collapse1 <- deparse(obj1[[8]]$collapse)}
if (is.null(collapse1))
{
if (is.null(obj1[[8]]$curve))
{
collapse1 <- "1 (for all parameters)"
} else {
collapse1 <- paste(deparse(obj1[[8]]$curve), "(for all parameters)")
}
} else {
# collapse1 <- paste(deparse(obj1[[8]]$collapse), collapse="")
if (!is.character(collapse1))
{
collapse1 <- paste(deparse(collapse1), collapse = "")
}
# collapse1 <- paste(deparse(collapse1), collapse = "")
collapse1 <- gsub(" ", "", collapse1, fixed = TRUE) # removing extra spaces
}
# pos <- 1
# if (is.data.frame(eval(collapse1)))
pos <- regexpr("data.frame(", collapse1, fixed = TRUE)
if (pos > 0)
{
# collapse1 <- deparse(obj1[[8]]$collapse)
# collapse1 <- substring(collapse1, pos+11, nchar(collapse1)-1)
collapse1 <- substring(collapse1, 12, nchar(collapse1)-1)
}
# if (is.list(eval(collapse1)))
pos <- regexpr("list(", collapse1, fixed=TRUE)
if (pos > 0)
{
# collapse1 <- deparse(obj1[[8]]$collapse)
# collapse1 <- substring(collapse1, pos+5, nchar(collapse1)-1)
collapse1 <- substring(collapse1, 6, nchar(collapse1)-1)
}
collapse2 <- obj2[[8]]$collapse
if (is.null(collapse2)) {collapse2 <- obj2[[8]]$pmodels}
if (!is.null(obj2$"pmodelsText")) {collapse2 <- obj2$"pmodelsText"}
# print(collapse2)
# if (is.null(collapse2)) {collapse2 <- paste(deparse(obj2[[8]]$assayNo), "(for all parameters)")} else {collapse2 <- deparse(obj2[[8]]$collapse)}
if (is.null(collapse2))
{
if (is.null(obj2[[8]]$curve))
{
collapse2 <- "1 (for all parameters)"
} else {
collapse2 <- paste(deparse(obj2[[8]]$curve), "(for all parameters)")
}
} else {
# collapse2 <- paste(deparse(obj2[[8]]$collapse), collapse = "")
if (!is.character(collapse2))
{
collapse2 <- paste(deparse(collapse2), collapse = "")
}
# collapse2 <- paste(deparse(collapse2), collapse = "")
collapse2 <- gsub(" ", "", collapse2, fixed = TRUE) # removing extra spaces
}
# if (is.data.frame(eval(collapse2)))
pos <- regexpr("data.frame(", collapse2, fixed = TRUE)
if (pos > 0)
{
# collapse2 <- deparse(obj2[[8]]$collapse)
# collapse2 <- substring(collapse2, pos+11, nchar(collapse2)-1)
collapse2 <- substring(collapse2, 12, nchar(collapse2) - 1)
}
# if (is.list(eval(collapse2)))
pos <- regexpr("list(", collapse2, fixed = TRUE)
if (pos > 0)
{
# collapse2 <- deparse(obj2[[8]]$collapse)
# collapse2 <- substring(collapse2, pos+5, nchar(collapse2)-1)
collapse2 <- substring(collapse2, 6, nchar(collapse2) - 1)
}
## Omitting collapse line if content is the same in both lines
if (identical(collapse1, collapse2)) {colLine <- FALSE} else {colLine <- TRUE}
if ( (!is.null(obj1[[8]]$pmodels)) || (!is.null(obj2[[8]]$pmodels)) )
{
fctStart <- " fct: "
colStart <- " pmodels: "
} else {
fctStart <- " fct: "
colStart <- " pmodels: "
# changed July 4 2011
}
cat("1st model\n")
fctInfo <- ifelse(is.null(obj1$"text"), deparse(obj1[[8]]$fct), obj1$"text")
cat(paste(fctStart, fctInfo, "\n", sep = ""))
# cat(paste(fctStart, deparse(obj1[[8]]$fct), "\n", sep = ""))
# if (colLine) {cat(paste(" collapse: ", collapse1, "\n", sep=""))}
if (colLine) {cat(paste(colStart, collapse1, "\n", sep = ""))}
cat("2nd model\n")
fctInfo <- ifelse(is.null(obj2$"text"), deparse(obj2[[8]]$fct), obj2$"text")
cat(paste(fctStart, fctInfo, "\n", sep = ""))
# cat(paste(fctStart, deparse(obj2[[8]]$fct), "\n", sep = ""))
# if (colLine) {cat(paste(" collapse: ", collapse2, "\n", sep=""))}
if (colLine) {cat(paste(colStart, collapse2, "\n", sep = ""))}
cat("\n")
}
# }
dataFra <- data.frame(dfModel, loglik, dfDiff, testStat, pVal)
dimnames(dataFra) <- list(rowNames, colNames)
structure(dataFra, heading = headName, class = c("anova", "data.frame"))
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.