diff --git a/ChangeLog b/ChangeLog index 983484d..940808f 100644 --- a/ChangeLog +++ b/ChangeLog @@ -1,11 +1,20 @@ +2016-11-21 vs. 1.2-6 - Ivan Bezerra Allaman + * Added the adjusted.pvalue argument in the TukeyC function. + * Removed dispersion argument of the TukeyC function and added in the plot function. + * Added three news methods in the package: print.TukeyC, TukeyC.formula and TukeyC.lm. + * Substitute the TukeyC.aov function by TukeyC.lm. + * Substitute the m.inf.1a, m.inf.1b, m.inf.2a, m.inf.2b, m.inf.3a and m.inf.3b by m.infos.lm, m.infos.aovlist, m.infos.nest.lm and m.infos.nest.aovlist. + * The TukeyC.nest.lm and TukeyC.nest.aovlist are internal. + * The TukeyC.nest functions were otimized. + * Added the dependency of the doBy package because the LSmeans function. + * Added the official url development of the package in description file. + 2014-08-16 vs. 1.2-5 - José Cláudio Faria - * A new function 'cv' was added. It allows to perform the coefficient - of the experiment variation for lm, aov and aovlist objetcs. + * A new function 'cv' was added. It allows to perform the coefficient of the experiment variation for lm, aov and aovlist objetcs. * Small changes to meet the requirements of CRAN. 2013-11-21 vs. 1.2-4 - José Cláudio Faria - * All calls to the old method model.frame.aovlist now call the generic - to meet the requirements of CRAN. + * All calls to the old method model.frame.aovlist now call the generic to meet the requirements of CRAN. * Released to CRAN. 2013-04-16 vs. 1.2-3 - José Cláudio Faria @@ -18,8 +27,7 @@ * Restrict to testers. 2012-12-10 vs. 1.1-1 - José Cláudio Faria - * A small bug was fixed in the function TukeyC.nest.aovlist.R. - The lines that were with "[:punct:]" (wrongly) were changed to "[[:punct:]]" (correctly). + * A small bug was fixed in the function TukeyC.nest.aovlist.R. The lines that were with "[:punct:]" (wrongly) were changed to "[[:punct:]]" (correctly). * New data: SPET (Split-plot in time). * Restrict to testers. @@ -58,3 +66,4 @@ 2010-05-27 vs. 1.0-0 - José Cláudio Faria * Restrict to testers. +GVIM 168 43 0 27 diff --git a/DESCRIPTION b/DESCRIPTION index cdc308d..d05c196 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,15 +1,17 @@ Package: TukeyC Type: Package Title: Conventional Tukey Test -Version: 1.1-5 -Date: 2014-08-16 +Version: 1.2-6 +Date: 2016-11-21 Author: José Cláudio Faria Enio G. Jelihovschi Ivan Bezerra Allaman Maintainer: José Cláudio Faria -Depends: R (>= 2.6.0), base -Description: Perform the conventional Tukey test from aov and aovlist - objects +Depends: R (>= 2.6.0), doBy +Description: Perform the conventional Tukey test from formula, aov, aovlist and lm objects. License: GPL (>= 2) +URL: https://github.com/jcfaria/TukeyC Encoding: latin1 LazyLoad: yes +NeedsCompilation: no +Packaged: 2016-11-23 00:59:54 UTC; ivan diff --git a/NAMESPACE b/NAMESPACE index 2328f3c..dba6ce7 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -1,5 +1,19 @@ -# Default NAMESPACE created by R -# Remove the previous line if you edit this file +export(TukeyC) +export(cv) -# Export all names -exportPattern(".") +import( + doBy, + utils, + stats, + graphics +) + +# S3 methods +S3method(TukeyC, formula) +S3method(TukeyC, lm) +S3method(TukeyC, aovlist) +S3method(TukeyC, nest.lm) +S3method(TukeyC, nest.aovlist) +S3method(summary, TukeyC) +S3method(print, TukeyC) +S3method(plot, TukeyC) diff --git a/R/TukeyC.R b/R/TukeyC.R index 28fdece..445a605 100644 --- a/R/TukeyC.R +++ b/R/TukeyC.R @@ -1,6 +1,4 @@ ## ## Main function S3 based ## - -TukeyC <- - function(x, ...) UseMethod('TukeyC') +TukeyC <- function(x, ...) UseMethod('TukeyC') diff --git a/R/TukeyC.aov.R b/R/TukeyC.aov.R deleted file mode 100644 index cd1dbe8..0000000 --- a/R/TukeyC.aov.R +++ /dev/null @@ -1,60 +0,0 @@ -## -## S3 method to 'aov' object -## - -TukeyC.aov <- function(x, - which=NULL, - sig.level=.05, - round=2, - dispersion=c('mm', 's', 'se'), ...) -{ - if(is.null(which)) - which <- names(x$model)[2] - - mt <- model.tables(x, - "means") # summary tables for model fits - - if(is.null(mt$n)) - stop("No factors in the fitted model!") - - tabs <- mt$tables[-1][which] # specified group means - - r <- mt$n[names(tabs)][[which]] # groups and its number of replicates - - bal <- ifelse(length(r) == 1, - TRUE, - FALSE) # is (or not) balanced - - MSE <- sum(resid(x)^2) / x$df.residual - - nms <- names(tabs[[which]]) - - ord <- order(as.vector(tabs[[which]]), - decreasing=TRUE) - - m.inf <- m.inf.1a(x, - which, - dispersion) - - rownames(m.inf) <- nms - - m.inf <- m.inf[order(m.inf[,1], - decreasing=TRUE),] - - dfr <- x$df.residual # residual degrees of freedom - - out <- make.TukeyC.test(r=r, - MSE=MSE, - m.inf=m.inf, - ord=ord, - sig.level=sig.level, - dfr=dfr, - bal=bal, - mt=mt, - round) - - class(out) <- c('TukeyC', - 'list') - - return(out) -} diff --git a/R/TukeyC.aovlist.R b/R/TukeyC.aovlist.R index d3f4ffd..3ed7939 100644 --- a/R/TukeyC.aovlist.R +++ b/R/TukeyC.aovlist.R @@ -3,56 +3,206 @@ ## TukeyC.aovlist <- function(x, - which, - error, - sig.level=.05, - round=2, - dispersion=c('mm', 's', 'se'), ...) + which = NULL, + fl1 = NULL, + fl2 = NULL, + error = NULL, + sig.level = .05, + round = 2, + adjusted.pvalue = 'none', ...) { - mt <- model.tables(x, - "means") # summary tables for model fits - if(is.null(mt$n)) - stop("No factors in the fitted model!") - tabs <- mt$tables[-1][which] # specified group means + # Interações com erro experimental + if(!is.null(fl1) & is.null(error)){ - r <- mt$n[names(tabs)][[which]] # groups and its number of replicates + pos_error <- length(names(x)) + SSE <- deviance(x[[pos_error]]) # experimental error + dfr <- df.residual(x[[pos_error]])# experimental error + MSE <- SSE/dfr - bal <- ifelse(length(r) == 1, - TRUE, - FALSE) # is (or not) balanced + class(x) <- c('nest.aovlist',class(x)) - MSE <- sum(resid(x[[error]])^2) / x[[error]][[8]] + res <- TukeyC(x = x, + which = which, + fl1 = fl1, + fl2 = fl2, + MSE = MSE, + dfr = dfr, + sig.level = sig.level, + round = round, + adjusted.pvalue = adjusted.pvalue, + ...) - nms <- names(tabs[[which]]) + class(res) <- c('TukeyC', + 'list') - ord <- order(as.vector(tabs[[which]]), - decreasing=TRUE) + return(res) - m.inf <- m.inf.1b(x, - which, - dispersion) - - rownames(m.inf) <- nms + } - m.inf <- m.inf[order(m.inf[,1], - decreasing=TRUE),] + # Interações com outros erros + if(!is.null(fl1) & !is.null(error)){ - dfr <- x[[error]][[8]] # residual degrees of freedom + many_errors <- unlist(strsplit(error, + '\\/')) - out <- make.TukeyC.test(r=r, - MSE=MSE, - m.inf=m.inf, - ord=ord, - sig.level=sig.level, - dfr=dfr, - bal=bal, - mt=mt, - round) + n_errors <- length(many_errors) - class(out) <- c('TukeyC', + if(n_errors > 1){# combinação de erros!!! + + aux_SSE <- NULL + aux_dfr <- NULL + + for(i in 1:n_errors){ + + aux_SSE[i] <- deviance(x[[many_errors[i]]]) + aux_dfr[i] <- df.residual(x[[many_errors[i]]]) + } + + aux_MSE <- aux_SSE/aux_dfr + + factors <- unlist(strsplit(which, + '[[:punct:]]')) + + aux_levels <- attr(x,'xlevel') + + aux_levels1 <- lapply(aux_levels, + length) + + levelss <- unlist(aux_levels1[factors]) + + if(length(levelss) == 2){ + + cp <- c(levelss[1]-1, + 1) + + MSE <- (cp%*%aux_MSE)/levelss[1] + + numer <- (cp%*%aux_MSE)^2 + denom <- (cp[1]*aux_MSE[1])^2/aux_dfr[1] + aux_MSE[2]^2/aux_dfr[2] + dfr <- numer/denom + + } else { + + cp <- c(levelss[2]*(levelss[1]-1), + levelss[2]-1, + 1) + + MSE <- (cp%*%aux_MSE)/prod(levelss[1:2]) + + numer <- (cp%*%aux_MSE)^2 + denom <- (cp[1]*aux_MSE[1])^2/aux_dfr[1] + (cp[2]*aux_MSE[2])^2/aux_dfr[2] + aux_MSE[3]^2/aux_dfr[3] + dfr <- numer/denom + + } + }else{# não há combinação de erros!!! + + SSE <- deviance(x[[error]]) # experimental error + dfr <- df.residual(x[[error]])# experimental error + MSE <- SSE/dfr + + } + class(x) <- c('nest.aovlist',class(x)) + + res <- TukeyC(x = x, + which = which, + fl1 = fl1, + fl2 = fl2, + MSE = MSE, + dfr = dfr, + sig.level = sig.level, + round = round, + adjusted.pvalue = adjusted.pvalue, + ...) + + class(res) <- c('TukeyC', + 'list') + + return(res) + + } + + # Aqui não há interesse em interações!!! + if(is.null(fl1) & !is.null(error)){ + + SSE <- deviance(x[[error]]) # experimental error + dfr <- df.residual(x[[error]])# experimental error + MSE <- SSE/dfr + + } else {# Erro experimental + + pos_error <- length(names(x)) + SSE <- deviance(x[[pos_error]]) # experimental error + dfr <- df.residual(x[[pos_error]])# experimental error + MSE <- SSE/dfr + + } + + my <- as.character(attr(x,'terms')[[2]]) + + #m1 <- gsub('\\:','\\+', which) + + #forminter <- as.formula(paste(my, '~', m1)) + forminter <- as.formula(paste(my, + '~', + which)) + + dat <- model.frame(x) + + aux_mt <- aggregate(forminter, + data = dat, + function(x) c(means = mean(x), + r = length(x))) + # + # aux_r <- aggregate(forminter, + # data = dat, + # function(x) r = length(x)) + # + # reps <- aux_r[[my]] + # + # aux_mt <- LSmeans(x, + # effect = which) + # + # aux_mt1 <- aux_mt$coef[,1] + # + # aux_mt2 <- data.frame(means = aux_mt1, + # reps = reps) + # + # row.names(aux_mt2) <- aux_r[,1] + # + # mt <- aux_mt2[order(aux_mt2[,1], + # decreasing = TRUE),] + # + aux_mt1 <- aux_mt[order(aux_mt[[my]][,1], + decreasing = TRUE),] + + mt <- data.frame(which = aux_mt1[,1], + means = aux_mt1[[my]][,1], + reps = aux_mt1[[my]][,2]) + + row.names(mt) <- aux_mt1[,1] + + out <- make.TukeyC.test(obj = mt, + MSE = MSE, + sig.level = sig.level, + dfr = dfr, + round = round, + adjusted.pvalue = adjusted.pvalue) + + m.inf <- m.infos.aovlist(x = x, + my = my, + forminter = forminter, + which = which, + sig.level = sig.level, + aux_mt = aux_mt, + MSE = MSE) + + res <- list(out = out, + info = m.inf) + + class(res) <- c('TukeyC', 'list') - return(out) + return(res) } diff --git a/R/TukeyC.default.R b/R/TukeyC.default.R index 19cd4a1..2ff9511 100644 --- a/R/TukeyC.default.R +++ b/R/TukeyC.default.R @@ -3,45 +3,11 @@ ## TukeyC.default <- function(x, - y=NULL, - model, - which, - error, - sig.level=.05, - round=2, - dispersion=c('mm', 's', 'se'), ...) + ...) { - if (is.data.frame(y)) - y <- as.matrix(y[, 1]) # manova is not contemplated - else - stopifnot(is.atomic(y)) - if (is.matrix(x) || is.atomic(x)) - x <- as.data.frame(x) + stop(paste("class", + class(x), + "objects are not valid for TukeyC" )) - if(!is.null(y)) - dat <- as.data.frame(cbind(x, - y)) - else - dat <- x - - av <- eval(substitute(aov(fo, - dat), - list(fo=formula(model)))) - - if(class(av)[1] == 'aov') - res <- TukeyC.aov(x=av, - which=which, - sig.level=sig.level, - round, - dispersion=dispersion) - else - res <- TukeyC.aovlist(x=av, - which=which, - error=error, - sig.level=sig.level, - round, - dispersion=dispersion) - - invisible(res) } diff --git a/R/TukeyC.formula.R b/R/TukeyC.formula.R new file mode 100644 index 0000000..62e8055 --- /dev/null +++ b/R/TukeyC.formula.R @@ -0,0 +1,73 @@ +## +## S3 method to design matrix and response variable or data.frame objects +## + +TukeyC.formula <- function(formula, + data = NULL, + which = NULL, + fl1 = NULL, + fl2 = NULL, + error = NULL, + sig.level = .05, + round = 2, + adjusted.pvalue = 'none', ...) +{ + + aux <- regexpr("Error", + as.character(formula), + perl=TRUE) + aux_err <- regmatches(as.character(formula), + aux) + + if(length(aux_err) == 0){ + + model <- lm(formula, + data = data) + + res <- TukeyC(x = model, + which = which, + fl1 = fl1, + fl2 = fl2, + error = error, + sig.level = sig.level, + round = round, + adjusted.pvalue = adjusted.pvalue, + ...) + + invisible(res) + + } else { + + basee <- aov(formula, + data = data) + + oc <- attr(basee, + "call") + Terms <- attr(basee, + "terms") + indError <- attr(Terms, + "specials")$Error + errorterm <- attr(Terms, + "variables")[[1 + indError]] + form <- update.formula(Terms, + paste(". ~ .-", + deparse(errorterm, + width.cutoff = 500L, backtick = TRUE), "+", deparse(errorterm[[2L]], + width.cutoff = 500L, backtick = TRUE))) + model <- lm(form, + data = data) + + res <- TukeyC(x = model, + which = which, + fl1 = fl1, + fl2 = fl2, + error = error, + sig.level = sig.level, + round = round, + adjusted.pvalue = adjusted.pvalue, + ...) + + invisible(res) + + } +} diff --git a/R/TukeyC.lm.R b/R/TukeyC.lm.R new file mode 100644 index 0000000..68ebee6 --- /dev/null +++ b/R/TukeyC.lm.R @@ -0,0 +1,209 @@ +TukeyC.lm <- function(x, + which = NULL, + fl1 = NULL, + fl2 = NULL, + error = NULL, + sig.level = .05, + round = 2, + adjusted.pvalue = 'none', ...) +{ + + if(is.null(which)){ + + which <- names(x$model)[2] + + } + + # Interações com erro experimental + if(!is.null(fl1) & is.null(error)){ + + SSE <- deviance(x) # sum square error + dfr <- df.residual(x) # residual degrees of freedom + MSE <- SSE/dfr + + class(x) <- c('nest.lm',class(x)) + + res <- TukeyC(x = x, + which = which, + fl1 = fl1, + fl2 = fl2, + MSE = MSE, + dfr = dfr, + sig.level = sig.level, + round = round, + adjusted.pvalue = adjusted.pvalue, + ...) + + class(res) <- c('TukeyC', + 'list') + + return(res) + } + + # Interações com outros erros + if(!is.null(fl1) & !is.null(error)){ + + many_errors <- unlist(strsplit(error, + '\\/')) + + ifelse(any(many_errors == 'Within'), + many_errors <- gsub('Within', + 'Residuals', + many_errors), + many_errors <- many_errors) + + n_errors <- length(many_errors) + + if(n_errors > 1){# combinação de erros!!! + + aux_MSE <- NULL + aux_dfr <- NULL + anov <- anova(x) + + for(i in 1:n_errors){ + + aux_MSE[i] <- anov[rownames(anov) == many_errors[i],][1,3] + aux_dfr[i] <- anov[rownames(anov) == many_errors[i],][1,1] + } + + factors <- unlist(strsplit(which, + '[[:punct:]]')) + + aux_levels <- x$xlevel + + aux_levels1 <- lapply(aux_levels, + length) + + levelss <- unlist(aux_levels1[factors]) + + if(length(levelss) == 2){ + + cp <- c(levelss[1]-1, + 1) + + MSE <- (cp%*%aux_MSE)/levelss[1] + + numer <- (cp%*%aux_MSE)^2 + denom <- (cp[1]*aux_MSE[1])^2/aux_dfr[1] + aux_MSE[2]^2/aux_dfr[2] + dfr <- numer/denom + + } else { + + cp <- c(levelss[2]*(levelss[1]-1), + levelss[2]-1, + 1) + + MSE <- (cp%*%aux_MSE)/prod(levelss[1:2]) + + numer <- (cp%*%aux_MSE)^2 + denom <- (cp[1]*aux_MSE[1])^2/aux_dfr[1] + (cp[2]*aux_MSE[2])^2/aux_dfr[2] + aux_MSE[3]^2/aux_dfr[3] + dfr <- numer/denom + + } + + } else {# não há combinação de erros!!! + + anov <- anova(x) + SSE <- anov[rownames(anov) == error,][1,2] # sum square error + dfr <- anov[rownames(anov) == error,][1,1] # residual degrees of freedom + MSE <- SSE/dfr + + } + + class(x) <- c('nest.lm',class(x)) + + res <- TukeyC(x = x, + which = which, + fl1 = fl1, + fl2 = fl2, + MSE = MSE, + dfr = dfr, + sig.level = sig.level, + round = round, + adjusted.pvalue = adjusted.pvalue, + ...) + + class(res) <- c('TukeyC', + 'list') + + return(res) + + } + + # Aqui não há interesse em interações!!!! + if(is.null(fl1) & !is.null(error)){#Um erro de interesse do usuário + + anov <- anova(x) + SSE <- anov[rownames(anov) == error,][1,2] # sum square error + dfr <- anov[rownames(anov) == error,][1,1] # residual degrees of freedom + MSE <- SSE/dfr + + } else { #Erro experimental + + SSE <- deviance(x) # sum square error + dfr <- df.residual(x) # residual degrees of freedom + MSE <- SSE/dfr + + } + + my <- as.character(formula(x)[[2]]) + #m1 <- gsub('\\:','\\+', which) + + #forminter <- as.formula(paste(my, '~', m1)) + forminter <- as.formula(paste(my, + '~', + which)) + + # aux_mt <- aggregate(forminter, + # data = x$model, + # function(x) c(means = mean(x), + # r = length(x))) + # + # mt <- aux_mt[order(aux_mt[[my]][,1], + # decreasing = TRUE),] + # + # row.names(mt) <- mt[,1] + # + + aux_r <- aggregate(forminter, + data = x$model, + function(x) r = length(x)) + reps <- aux_r[[my]] + + aux_mt <- suppressWarnings(LSmeans(x, + effect = which, + level = 1 - sig.level)) + + aux_mt1 <- aux_mt$coef[,1] + + aux_mt2 <- data.frame(means = aux_mt1, + reps = reps) + + row.names(aux_mt2) <- aux_r[,1] + + mt <- aux_mt2[order(aux_mt2[,1], + decreasing = TRUE),] + + out <- make.TukeyC.test(obj = mt, + MSE = MSE, + sig.level = sig.level, + dfr = dfr, + round = round, + adjusted.pvalue = adjusted.pvalue) + + m.inf <- m.infos.lm(x = x, + my = my, + forminter = forminter, + which = which, + sig.level = sig.level, + aux_mt = aux_mt) + + res <- list(out = out, + info = m.inf) + + class(res) <- c('TukeyC', + 'list') + + return(res) + +} diff --git a/R/TukeyC.lm.R_bug b/R/TukeyC.lm.R_bug new file mode 100644 index 0000000..e6136a8 --- /dev/null +++ b/R/TukeyC.lm.R_bug @@ -0,0 +1,203 @@ +TukeyC.lm <- function(x, + which = NULL, + fl1 = NULL, + fl2 = NULL, + error = NULL, + sig.level = .05, + round = 2, + adjusted.pvalue = 'none', ...) +{ + + if(is.null(which)){ + + which <- names(x$model)[2] + + } + + # Interações com erro experimental + if(!is.null(fl1) & is.null(error)){ + + SSE <- deviance(x) # sum square error + dfr <- df.residual(x) # residual degrees of freedom + MSE <- SSE/dfr + + class(x) <- c('nest.lm',class(x)) + + res <- TukeyC(x = x, + which = which, + fl1 = fl1, + fl2 = fl2, + MSE = MSE, + dfr = dfr, + sig.level = sig.level, + round = round, + adjusted.pvalue = adjusted.pvalue, + ...) + + class(res) <- c('TukeyC', + 'list') + + return(res) + } + + # Interações com outros erros + if(!is.null(fl1) & !is.null(error)){ + + many_errors <- unlist(strsplit(error, + '\\/')) + + n_errors <- length(many_errors) + + if(n_errors > 1){# combinação de erros!!! + + aux_MSE <- NULL + aux_dfr <- NULL + anov <- anova(x) + + for(i in 1:n_errors){ + + aux_MSE[i] <- anov[rownames(anov) == many_errors[i],][1,3] + aux_dfr[i] <- anov[rownames(anov) == many_errors[i],][1,1] + } + + factors <- unlist(strsplit(which, + '[[:punct:]]')) + + aux_levels <- x$xlevel + + aux_levels1 <- lapply(aux_levels, + length) + + levelss <- unlist(aux_levels1[factors]) + + if(length(levelss) == 2){ + + cp <- c(levelss[1]-1, + 1) + + MSE <- (cp%*%aux_MSE)/levelss[1] + + numer <- (cp%*%aux_MSE)^2 + denom <- (cp[1]*aux_MSE[1])^2/aux_dfr[1] + aux_MSE[2]^2/aux_dfr[2] + dfr <- numer/denom + + } else { + + cp <- c(levelss[2]*(levelss[1]-1), + levelss[2]-1, + 1) + + MSE <- (cp%*%aux_MSE)/prod(levelss[1:2]) + + numer <- (cp%*%aux_MSE)^2 + denom <- (cp[1]*aux_MSE[1])^2/aux_dfr[1] + (cp[2]*aux_MSE[2])^2/aux_dfr[2] + aux_MSE[3]^2/aux_dfr[3] + dfr <- numer/denom + + } + + } else {# não há combinação de erros!!! + + anov <- anova(x) + SSE <- anov[rownames(anov) == error,][1,2] # sum square error + dfr <- anov[rownames(anov) == error,][1,1] # residual degrees of freedom + MSE <- SSE/dfr + + } + + class(x) <- c('nest.lm',class(x)) + + res <- TukeyC(x = x, + which = which, + fl1 = fl1, + fl2 = fl2, + MSE = MSE, + dfr = dfr, + sig.level = sig.level, + round = round, + adjusted.pvalue = adjusted.pvalue, + ...) + + class(res) <- c('TukeyC', + 'list') + + return(res) + + } + + # Aqui não há interesse em interações!!!! + if(is.null(fl1) & !is.null(error)){#Um erro de interesse do usuário + + anov <- anova(x) + SSE <- anov[rownames(anov) == error,][1,2] # sum square error + dfr <- anov[rownames(anov) == error,][1,1] # residual degrees of freedom + MSE <- SSE/dfr + + } else { #Erro experimental + + SSE <- deviance(x) # sum square error + dfr <- df.residual(x) # residual degrees of freedom + MSE <- SSE/dfr + + } + + my <- as.character(formula(x)[[2]]) + #m1 <- gsub('\\:','\\+', which) + + #forminter <- as.formula(paste(my, '~', m1)) + forminter <- as.formula(paste(my, + '~', + which)) + + # aux_mt <- aggregate(forminter, + # data = x$model, + # function(x) c(means = mean(x), + # r = length(x))) + # + # mt <- aux_mt[order(aux_mt[[my]][,1], + # decreasing = TRUE),] + # + # row.names(mt) <- mt[,1] + # + + aux_r <- aggregate(forminter, + data = x$model, + function(x) r = length(x)) + reps <- aux_r[[my]] + + aux_mt <- suppressWarnings(LSmeans(x, + effect = which, + level = 1 - sig.level)) + + aux_mt1 <- aux_mt$coef[,1] + + aux_mt2 <- data.frame(means = aux_mt1, + reps = reps) + + row.names(aux_mt2) <- aux_r[,1] + + mt <- aux_mt2[order(aux_mt2[,1], + decreasing = TRUE),] + + out <- make.TukeyC.test(obj = mt, + MSE = MSE, + sig.level = sig.level, + dfr = dfr, + round = round, + adjusted.pvalue = adjusted.pvalue) + + m.inf <- m.infos.lm(x = x, + my = my, + forminter = forminter, + which = which, + sig.level = sig.level, + aux_mt = aux_mt) + + res <- list(out = out, + info = m.inf) + + class(res) <- c('TukeyC', + 'list') + + return(res) + +} diff --git a/R/TukeyC.nest.R b/R/TukeyC.nest.R deleted file mode 100644 index edf91b3..0000000 --- a/R/TukeyC.nest.R +++ /dev/null @@ -1,6 +0,0 @@ -## -## Main function S3 based -## - -TukeyC.nest <- - function(x, ...) UseMethod('TukeyC.nest') diff --git a/R/TukeyC.nest.aov.R b/R/TukeyC.nest.aov.R deleted file mode 100644 index 4a02452..0000000 --- a/R/TukeyC.nest.aov.R +++ /dev/null @@ -1,521 +0,0 @@ -## -## S3 method to 'aov' object -## - -TukeyC.nest.aov <- function (x, - which, - fl1, - fl2=0, - sig.level=.05, - round=2, - dispersion=c('mm', 's', 'se'), ...) -{ - mt <- model.tables(x, - "means") # summary tables for model fits - - if(is.null(mt$n)) - stop("No factors in the fitted model!") - - nfa <- names(mt$tables)[-1] # nomes dos fatores - - nf1 <- unlist(strsplit(which, - split=':'))[1] # nome do primeiro fator do which - - nf2 <- unlist(strsplit(which, - split=':'))[2] # nome do segundo fator do which - - nf3 <- unlist(strsplit(which, - split=':'))[3] # nome do terceiro fator do which - - MSE <- deviance(x)/df.residual(x) - - group <- NULL - - group2 <- NULL - - if(fl2 == 0){ - # MODELO SOMENTE COM DOIS FATORES - if(length(nfa[grep('[[:punct:]]', - nfa)]) == 1 && # condição necessária para certificar-se que no modelo há somente uma interação! - which != nfa[grep('[[:punct:]]', - nfa)]) { - whichn <- paste(nf2, - nf1, - sep=':') # necessário pois no modelo original estamos em uma ordem inversa! - - r <- mt$n[names(mt$tables)][[whichn]] # groups and its number of replicates - - bal <- ifelse(length(r) == 1, - TRUE, - FALSE) # is (or not) balanced - - m <- as.vector(mt$tables[whichn][[whichn]][,fl1]) # pegando as médias de interesse - - which1 <- names(dimnames(mt$tables[whichn][[whichn]]))[2] # corresponde ao primeiro fator do seu 'which' - - which2 <- names(dimnames(mt$tables[whichn][[whichn]]))[1] # corresponde ao segundo fator do seu 'which' - - m.inf <- m.inf.2a(x, - which1, - which2, - dispersion) - - f1 <- levels(x$model[,which2]) # correspondem aos fatores que se quer comparar! - - f2 <- levels(x$model[,which1])[fl1] # corresponde ao fator onde se está fazendo o desdobramento! - - m.inf <- subset(m.inf, - group == f2)[,2] - - rownames(m.inf) <- paste(f2, - f1, - sep='/') - - ord <- order(m, - decreasing=TRUE) - - m.inf <- cbind(m.inf[,1][ord], - m.inf[,2][ord], - m.inf[,3][ord]) - - tab <- mt$tables[whichn][[whichn]] - } else if(length(nfa[grep('[[:punct:]]', - nfa)]) == 1 && - which == nfa[grep('[[:punct:]]', - nfa)]) { - - r <- mt$n[names(mt$tables)][[which]] # groups and its number of replicates - - bal <- ifelse(length(r) == 1, - TRUE, - FALSE) # is (or not) balanced - - m <- as.vector(mt$tables[which][[which]][fl1,]) # pegando as médias de interesse - - which1 <- names(dimnames(mt$tables[which][[which]]))[1] # corresponde ao primeiro fator do seu 'which' - - which2 <- names(dimnames(mt$tables[which][[which]]))[2] # corresponde ao segundo fator do seu 'which' - - m.inf <- m.inf.2a(x, - which1, - which2, - dispersion) - - f1 <- levels(x$model[,which2]) # correspondem aos fatores que se quer comparar! - - f2 <- levels(x$model[,which1])[fl1] # corresponde ao fator onde se está fazendo o desdobramento! - - m.inf <- subset(m.inf, - group == f2)[,2] - - rownames(m.inf) <- paste(f2, - f1, - sep='/') - - ord <- order(m, - decreasing=TRUE) - - m.inf <- cbind(m.inf[,1][ord], - m.inf[,2][ord], - m.inf[,3][ord]) - - tab <- mt$tables[which][[which]] - } else - # MODELO COM TRES FATORES - if(length(nfa[grep('[[:punct:]]', - nfa)]) != 1 && - which == nfa[grep('[[:punct:]]', - nfa)][1] | - which == nfa[grep('[[:punct:]]', - nfa)][2] | - which == nfa[grep('[[:punct:]]', - nfa)][3] ) { - - r <- mt$n[names(mt$tables)][[which]] # groups and its number of replicates - - bal <- ifelse(length(r) == 1, - TRUE, - FALSE) # is (or not) balanced - - m <- as.vector(mt$tables[which][[which]][fl1,]) - - which1 <- names(dimnames(mt$tables[which][[which]]))[1] - - which2 <- names(dimnames(mt$tables[which][[which]]))[2] - - m.inf <- m.inf.2a(x, - which1, - which2, - dispersion) - - f1 <- levels(x$model[,which2]) - - f2 <- levels(x$model[,which1])[fl1] - - m.inf <- subset(m.inf, - group == f2)[,2] - - rownames(m.inf) <- paste(f2, - f1, - sep='/') - - ord <- order(m, decreasing=TRUE) - - m.inf <- cbind(m.inf[,1][ord], - m.inf[,2][ord], - m.inf[,3][ord]) - - tab <- mt$tables[which][[which]] - } else if(length(nfa[grep('[[:punct:]]', - nfa)]) != 1 && - which != nfa[grep('[[:punct:]]', - nfa)][1]) { - - t1 <- unlist(strsplit(which, - split=':'))[1] - - t2 <- unlist(strsplit(which, - split=':'))[2] - - whichn <- paste(t2, - t1, - sep=':') - - r <- mt$n[names(mt$tables)][[whichn]] # groups and its number of replicates - - bal <- ifelse(length(r) == 1, - TRUE, - FALSE) # is (or not) balanced - - m <- as.vector(mt$tables[whichn][[whichn]][,fl1]) - - which1 <- names(dimnames(mt$tables[whichn][[whichn]]))[2] - - which2 <- names(dimnames(mt$tables[whichn][[whichn]]))[1] - - m.inf <- m.inf.2a(x, - which1, - which2, - dispersion) - - f1 <- levels(x$model[,which2]) - - f2 <- levels(x$model[,which1])[fl1] - - m.inf <- subset(m.inf, - group == f2)[,2] - - rownames(m.inf) <- paste(f2, - f1, - sep='/') - - ord <- order(m, - decreasing=TRUE) - - m.inf <- cbind(m.inf[,1][ord], - m.inf[,2][ord], - m.inf[,3][ord]) - - tab <- mt$tables[whichn][[whichn]] - } - } else if(fl2 != 0) { - natri <- nfa[grep('[[:punct:]].{,100}[[:punct:]]', - nfa)] - - nt1 <- unlist(strsplit(natri, - split=':'))[1] - - nt2 <- unlist(strsplit(natri, - split=':'))[2] - - nt3 <- unlist(strsplit(natri, - split=':'))[3] - if(which == natri){ - r <- mt$n[names(mt$tables)][[which]] # groups and its number of replicates - - bal <- ifelse(length(r) == 1, - TRUE, - FALSE) # is (or not) balanced - - which1 <- names(dimnames(mt$tables[which][[which]]))[1] - - which2 <- names(dimnames(mt$tables[which][[which]]))[2] - - which3 <- names(dimnames(mt$tables[which][[which]]))[3] - - m.inf <- m.inf.3a(x, - which1=which3, - which2=which2, - which3=which1, - dispersion) - - f1 <- levels(x$model[,which3]) - - f2 <- levels(x$model[,which2])[fl2] - - f3 <- levels(x$model[,which1])[fl1] - - m.inf <- subset(m.inf, - group == f2 & - group2 == f3)[,3] - - rownames(m.inf) <- paste(f3, - f2, - f1, - sep='/') - - ord <- order(as.vector(m.inf[,1]), - decreasing=TRUE) - - m.inf <- cbind(m.inf[,1][ord], - m.inf[,2][ord], - m.inf[,3][ord]) - - tab <- mt$tables[which][[which]] - } else if(which != natri){ - if(nf1 == nt2 && nf2 == nt3) { - - r <- mt$n[names(mt$tables)][[natri]] # groups and its number of replicates - - bal <- ifelse(length(r) == 1, - TRUE, - FALSE) # is (or not) balanced - - m <- as.vector(mt$tables[natri][[natri]][,fl1,fl2]) - - which1 <- names(dimnames(mt$tables[natri][[natri]]))[1] - - which2 <- names(dimnames(mt$tables[natri][[natri]]))[2] - - which3 <- names(dimnames(mt$tables[natri][[natri]]))[3] - - m.inf <- m.inf.3a(x, - which1=which1, - which2=which3, - which3=which2, - dispersion) - - f1 <- levels(x$model[,which1]) - - f2 <- levels(x$model[,which3])[fl2] - - f3 <- levels(x$model[,which2])[fl1] - - m.inf <- subset(m.inf, - group == f2 & - group2 == f3)[,3] - - rownames(m.inf) <- paste(f3, - f2, - f1, - sep='/') - - ord <- order(m, - decreasing=TRUE) - - m.inf <- cbind(m.inf[,1][ord], - m.inf[,2][ord], - m.inf[,3][ord]) - - tab <- mt$tables[natri][[natri]] - } else if(nf1 == nt2 && nf2 == nt1) { - r <- mt$n[names(mt$tables)][[natri]] # groups and its number of replicates - - bal <- ifelse(length(r) == 1, - TRUE, - FALSE) # is (or not) balanced - - which1 <- names(dimnames(mt$tables[natri][[natri]]))[1] - - which2 <- names(dimnames(mt$tables[natri][[natri]]))[2] - - which3 <- names(dimnames(mt$tables[natri][[natri]]))[3] - - m.inf <- m.inf.3a(x, - which1=which3, - which2=which1, - which3=which2, - dispersion) - - f1 <- levels(x$model[,which3]) - - f2 <- levels(x$model[,which1])[fl2] - - f3 <- levels(x$model[,which2])[fl1] - - m.inf <- subset(m.inf, - group == f2 & - group2 == f3)[,3] - - rownames(m.inf) <- paste(f3, - f2, - f1, - sep='/') - - ord <- order(as.vector(m.inf[,1]), - decreasing=TRUE) - - m.inf <- cbind(m.inf[,1][ord], - m.inf[,2][ord], - m.inf[,3][ord]) - - tab <- mt$tables[natri][[natri]] - } else if(nf1 == nt1 && nf2 == nt3) { - r <- mt$n[names(mt$tables)][[natri]] # groups and its number of replicates - - bal <- ifelse(length(r) == 1, - TRUE, - FALSE) # is (or not) balanced - - m <- as.vector(mt$tables[natri][[natri]][fl1,,fl2]) - - which1 <- names(dimnames(mt$tables[natri][[natri]]))[1] - - which2 <- names(dimnames(mt$tables[natri][[natri]]))[2] - - which3 <- names(dimnames(mt$tables[natri][[natri]]))[3] - - m.inf <- m.inf.3a(x, - which1=which2, - which2=which3, - which3=which1, - dispersion) - - f1 <- levels(x$model[,which2]) - - f2 <- levels(x$model[,which3])[fl2] - - f3 <- levels(x$model[,which1])[fl1] - - m.inf <- subset(m.inf, - group == f2 & - group2 == f3)[,3] - - rownames(m.inf) <- paste(f3, - f2, - f1, - sep='/') - - ord <- order(m, - decreasing=TRUE) - - m.inf <- cbind(m.inf[,1][ord], - m.inf[,2][ord], - m.inf[,3][ord]) - - tab <- mt$tables[natri][[natri]] - } else if(nf1 == nt3 && nf2 == nt2) { - r <- mt$n[names(mt$tables)][[natri]] # groups and its number of replicates - - bal <- ifelse(length(r) == 1, - TRUE, - FALSE) # is (or not) balanced - - m <- as.vector(mt$tables[natri][[natri]][,fl2,fl1]) - - which1 <- names(dimnames(mt$tables[natri][[natri]]))[1] - - which2 <- names(dimnames(mt$tables[natri][[natri]]))[2] - - which3 <- names(dimnames(mt$tables[natri][[natri]]))[3] - - m.inf <- m.inf.3a(x, - which1=which1, - which2=which2, - which3=which3, - dispersion) - - f1 <- levels(x$model[,which1]) - - f2 <- levels(x$model[,which2])[fl2] - - f3 <- levels(x$model[,which3])[fl1] - - m.inf <- subset(m.inf, - group == f2 & - group2 == f3)[,3] - - rownames(m.inf) <- paste(f3, - f2, - f1, - sep='/') - - ord <- order(m, - decreasing=TRUE) - - m.inf <- cbind(m.inf[,1][ord], - m.inf[,2][ord], - m.inf[,3][ord]) - - tab <- mt$tables[natri][[natri]] - } else { - r <- mt$n[names(mt$tables)][[natri]] # groups and its number of replicates - - bal <- ifelse(length(r) == 1, - TRUE, - FALSE) # is (or not) balanced - - m <- as.vector(mt$tables[natri][[natri]][fl2,,fl1]) - - which1 <- names(dimnames(mt$tables[natri][[natri]]))[1] - - which2 <- names(dimnames(mt$tables[natri][[natri]]))[2] - - which3 <- names(dimnames(mt$tables[natri][[natri]]))[3] - - m.inf <- m.inf.3a(x, - which1=which2, - which2=which1, - which3=which3, - dispersion) - - f1 <- levels(x$model[,which2]) - - f2 <- levels(x$model[,which1])[fl2] - - f3 <- levels(x$model[,which3])[fl1] - - m.inf <- subset(m.inf, - group == f2 & - group2 == f3)[,3] - - rownames(m.inf) <- paste(f3, - f2, - f1, - sep='/') - - ord <- order(m, - decreasing=TRUE) - - m.inf <- cbind(m.inf[,1][ord], - m.inf[,2][ord], - m.inf[,3][ord]) - - tab <- mt$tables[natri][[natri]] - } - } - } - colnames(m.inf) <- c('mean', - 'min', - 'max') - - mMSE <- MSE/r - - dfr <- df.residual(x) # residual degrees of freedom - - out <- make.TukeyC.test(r=r, - MSE=MSE, - m.inf=m.inf, - ord=ord, - sig.level=sig.level, - dfr=dfr, - bal=bal, - mt=mt, - round) - - class(out) <- c('TukeyC.nest', - 'TukeyC', - 'list') - - return(out) -} diff --git a/R/TukeyC.nest.aovlist.R b/R/TukeyC.nest.aovlist.R index e329ea2..2568608 100644 --- a/R/TukeyC.nest.aovlist.R +++ b/R/TukeyC.nest.aovlist.R @@ -1,514 +1,119 @@ ## ## S3 method to 'aovlist' object ## - TukeyC.nest.aovlist <- function(x, which, - error, fl1, - fl2=0, - sig.level=.05, - round=2, - dispersion=c('mm', 's', 'se'), ...) + fl2, + MSE, + dfr, + sig.level, + round, + adjusted.pvalue,...) { - mt <- model.tables(x, - "means") # summary tables for model fits - - if(is.null(mt$n)) - stop("No factors in the fitted model!") - - MSE <- sum(resid(x[[error]])^2) / x[[error]][[8]] - - nfa <- names(mt$tables)[-1] # nomes dos fatores + my <- as.character(attr(x,'terms')[[2]]) + m1 <- gsub('\\:', + '\\+', + which) + # m2 <- unlist(strsplit(which, + # '[[:punct:]]')) + # + forminter <- as.formula(paste(my, '~', m1)) + + dat <- model.frame(x) + + aux_mt1 <- aggregate(forminter, + data = dat, + function(x) c(means = mean(x), + r = length(x))) + + aux_mt2 <- aux_mt1[order(aux_mt1[[my]][,1], + decreasing = TRUE),] + + aux_mt3 <- data.frame(aux_mt2[1:length(names(aux_mt2))-1], + means = aux_mt2[[my]][,1], + reps = aux_mt2[[my]][,2]) + + #row.names(mt) <- aux_mt1[,1] + + #names(aux_mt2) <- gsub(my,'x',names(aux_mt2)) + # + # aux_r <- aggregate(forminter, + # data = x$model, + # function(x) r = length(x)) + # reps <- aux_r[[my]] + # + # aux_mt <- LSmeans(x, + # effect = m2) + # + # aux_mt1 <- aux_mt$coef[,1] + # + # aux_mt2 <- data.frame(aux_r[1:length(names(aux_r))-1], + # means = aux_mt1, + # reps = reps) + # + # aux_mt3 <- aux_mt2[order(aux_mt2[['means']], + # decreasing = TRUE),] + # nf1 <- unlist(strsplit(which, - split=':'))[1] # nome do primeiro fator do which + split = ':'))[1] # nome do primeiro fator do which nf2 <- unlist(strsplit(which, - split=':'))[2] # nome do segundo fator do which + split = ':'))[2] # nome do segundo fator do which nf3 <- unlist(strsplit(which, - split=':'))[3] # nome do terceiro fator do which - group <- NULL - - group2 <- NULL - - # SPLIT-PLOT - if (fl2 == 0){ - # MODELO SOMENTE COM DOIS FATORES - if(length(nfa[grep('[[:punct:]]', - nfa)]) == 1 && # condição necessária para certificar-se que no modelo há somente uma interação! - which != nfa[grep('[[:punct:]]', - nfa)]) { - whichn <- paste(nf2, - nf1, - sep=':') # necessário pois no modelo original estamos em uma ordem inversa! - - r <- mt$n[names(mt$tables)][[whichn]] # groups and its number of replicates - - bal <- ifelse(length(r) == 1, - TRUE, - FALSE) # is (or not) balanced - - m <- as.vector(mt$tables[whichn][[whichn]][,fl1]) # pegando as médias de interesse - - which1 <- names(dimnames(mt$tables[whichn][[whichn]]))[2] # corresponde ao primeiro fator do seu 'which' - - which2 <- names(dimnames(mt$tables[whichn][[whichn]]))[1] # corresponde ao segundo fator do seu 'which' - m.inf <- m.inf.2b(x, - which1, - which2, - dispersion) - - f1 <- levels(model.frame(x)[,which2]) # correspondem aos fatores que se quer comparar! - - f2 <- levels(model.frame(x)[,which1])[fl1] # corresponde ao fator onde se está fazendo o desdobramento! - - m.inf <- subset(m.inf, - group == f2)[,2] - - rownames(m.inf) <- paste(f2, - f1, - sep='/') - - ord <- order(m, - decreasing=TRUE) - - m.inf <- cbind(m.inf[,1][ord], - m.inf[,2][ord], - m.inf[,3][ord]) - - tab <- mt$tables[whichn][[whichn]] - } else if(length(nfa[grep('[[:punct:]]', - nfa)]) == 1 && - which == nfa[grep('[[:punct:]]', - nfa)]) { - r <- mt$n[names(mt$tables)][[which]] # groups and its number of replicates - - bal <- ifelse(length(r) == 1, - TRUE, - FALSE) # is (or not) balanced - - m <- as.vector(mt$tables[which][[which]][fl1,]) # pegando as médias de interesse - - which1 <- names(dimnames(mt$tables[which][[which]]))[1] # corresponde ao primeiro fator do seu 'which' - - which2 <- names(dimnames(mt$tables[which][[which]]))[2] # corresponde ao segundo fator do seu 'which' - - m.inf <- m.inf.2b(x, - which1, - which2, - dispersion) - - f1 <- levels(model.frame(x)[,which2]) # correspondem aos fatores que se quer comparar! - - f2 <- levels(model.frame(x)[,which1])[fl1] # corresponde ao fator onde se está fazendo o desdobramento! - - m.inf <- subset(m.inf, - group == f2)[,2] - - rownames(m.inf) <- paste(f2, - f1, - sep='/') - ord <- order(m, - decreasing=TRUE) - - m.inf <- cbind(m.inf[,1][ord], - m.inf[,2][ord], - m.inf[,3][ord]) - - tab <- mt$tables[which][[which]] - } else - # MODELO COM TRES FATORES - if(length(nfa[grep('[[:punct:]]', - nfa)]) != 1 && - which == nfa[grep('[[:punct:]]', - nfa)][1] | - which == nfa[grep('[[:punct:]]', - nfa)][2] | - which == nfa[grep('[[:punct:]]', - nfa)][3]) { - r <- mt$n[names(mt$tables)][[which]] # groups and its number of replicates - - bal <- ifelse(length(r) == 1, - TRUE, - FALSE) # is (or not) balanced - - m <- as.vector(mt$tables[which][[which]][fl1,]) - - which1 <- names(dimnames(mt$tables[which][[which]]))[1] - - which2 <- names(dimnames(mt$tables[which][[which]]))[2] - - m.inf <- m.inf.2b(x, - which1, - which2, - dispersion) - - f1 <- levels(model.frame(x)[,which2]) - - f2 <- levels(model.frame(x)[,which1])[fl1] - - m.inf <- subset(m.inf, - group == f2)[,2] - - rownames(m.inf) <- paste(f2, - f1, - sep='/') - - ord <- order(m, - decreasing=TRUE) - - m.inf <- cbind(m.inf[,1][ord], - m.inf[,2][ord], - m.inf[,3][ord]) - - tab <- mt$tables[which][[which]] - } else if(length(nfa[grep('[[:punct:]]', - nfa)]) != 1 && - which != nfa[grep('[[:punct:]]', - nfa)][1]) { - t1 <- unlist(strsplit(which, - split=':'))[1] - - t2 <- unlist(strsplit(which, - split=':'))[2] - - whichn <- paste(t2, - t1, - sep=':') - - r <- mt$n[names(mt$tables)][[whichn]] # groups and its number of replicates - - bal <- ifelse(length(r) == 1, - TRUE, - FALSE) # is (or not) balanced - - m <- as.vector(mt$tables[whichn][[whichn]][,fl1]) - - which1 <- names(dimnames(mt$tables[whichn][[whichn]]))[2] - - which2 <- names(dimnames(mt$tables[whichn][[whichn]]))[1] - - m.inf <- m.inf.2b(x, - which1, - which2, - dispersion) - - f1 <- levels(model.frame(x)[,which2]) - - f2 <- levels(model.frame(x)[,which1])[fl1] - - m.inf <- subset(m.inf, - group == f2)[,2] - - rownames(m.inf) <- paste(f2, - f1, - sep='/') - - ord <- order(m, - decreasing=TRUE) - - m.inf <- cbind(m.inf[,1][ord], - m.inf[,2][ord], - m.inf[,3][ord]) - - tab <- mt$tables[whichn][[whichn]] - } - } else if(fl2 != 0) { - natri <- nfa[grep('[[:punct:]].{,100}[[:punct:]]', - nfa)] - nt1 <- unlist(strsplit(natri, - split=':'))[1] - nt2 <- unlist(strsplit(natri, - split=':'))[2] - nt3 <- unlist(strsplit(natri, - split=':'))[3] - if(which == natri){ - r <- mt$n[names(mt$tables)][[which]] # groups and its number of replicates - - bal <- ifelse(length(r) == 1, - TRUE, - FALSE) # is (or not) balanced - - which1 <- names(dimnames(mt$tables[which][[which]]))[1] - - which2 <- names(dimnames(mt$tables[which][[which]]))[2] - - which3 <- names(dimnames(mt$tables[which][[which]]))[3] - - m.inf <- m.inf.3b(x, - which1=which3, - which2=which2, - which3=which1, - dispersion) - - f1 <- levels(model.frame(x)[,which3]) - - f2 <- levels(model.frame(x)[,which2])[fl2] - - f3 <- levels(model.frame(x)[,which1])[fl1] - - m.inf <- subset(m.inf, - group == f2 & - group2 == f3)[,3] - - rownames(m.inf) <- paste(f3, - f2, - f1, - sep='/') - - ord <- order(as.vector(m.inf[,1]), - decreasing=TRUE) - - m.inf <- cbind(m.inf[,1][ord], - m.inf[,2][ord], - m.inf[,3][ord]) - - tab <- mt$tables[which][[which]] - } else if(which != natri){ - if(nf1 == nt2 && nf2 == nt3){ - r <- mt$n[names(mt$tables)][[natri]]# groups and its number of replicates - - bal <- ifelse(length(r) == 1, - TRUE, - FALSE) # is (or not) balanced - - m <- as.vector(mt$tables[natri][[natri]][,fl1,fl2]) - - which1 <- names(dimnames(mt$tables[natri][[natri]]))[1] - - which2 <- names(dimnames(mt$tables[natri][[natri]]))[2] - - which3 <- names(dimnames(mt$tables[natri][[natri]]))[3] - - m.inf <- m.inf.3b(x, - which1=which1, - which2=which3, - which3=which2, - dispersion) - - f1 <- levels(model.frame(x)[,which1]) - - f2 <- levels(model.frame(x)[,which3])[fl2] - - f3 <- levels(model.frame(x)[,which2])[fl1] - - m.inf <- subset(m.inf, - group == f2 & - group2 == f3)[,3] - - rownames(m.inf) <- paste(f3, - f2, - f1, - sep='/') - - ord <- order(m, - decreasing=TRUE) - - m.inf <- cbind(m.inf[,1][ord], - m.inf[,2][ord], - m.inf[,3][ord]) - - tab <- mt$tables[natri][[natri]] - } else if(nf1 == nt2 && nf2 == nt1) { - - r <- mt$n[names(mt$tables)][[natri]] # groups and its number of replicates - - bal <- ifelse(length(r) == 1, - TRUE, - FALSE) # is (or not) balanced - - which1 <- names(dimnames(mt$tables[natri][[natri]]))[1] - - which2 <- names(dimnames(mt$tables[natri][[natri]]))[2] - - which3 <- names(dimnames(mt$tables[natri][[natri]]))[3] - - m.inf <- m.inf.3b(x, - which1=which3, - which2=which1, - which3=which2, - dispersion) - - f1 <- levels(model.frame(x)[,which3]) - - f2 <- levels(model.frame(x)[,which1])[fl2] - - f3 <- levels(model.frame(x)[,which2])[fl1] - - m.inf <- subset(m.inf, - group == f2 & - group2 == f3)[,3] - - rownames(m.inf) <- paste(f3, - f2, - f1, - sep='/') - - ord <- order(as.vector(m.inf[,1]), - decreasing=TRUE) - - m.inf <- cbind(m.inf[,1][ord], - m.inf[,2][ord], - m.inf[,3][ord]) - - tab <- mt$tables[natri][[natri]] - } else if(nf1 == nt1 && nf2 == nt3) { - r <- mt$n[names(mt$tables)][[natri]] # groups and its number of replicates - - bal <- ifelse(length(r) == 1, - TRUE, - FALSE) # is (or not) balanced - - m <- as.vector(mt$tables[natri][[natri]][fl1,,fl2]) - - which1 <- names(dimnames(mt$tables[natri][[natri]]))[1] - - which2 <- names(dimnames(mt$tables[natri][[natri]]))[2] - - which3 <- names(dimnames(mt$tables[natri][[natri]]))[3] - - m.inf <- m.inf.3b(x, - which1=which2, - which2=which3, - which3=which1, - dispersion) - - f1 <- levels(model.frame(x)[,which2]) - - f2 <- levels(model.frame(x)[,which3])[fl2] - - f3 <- levels(model.frame(x)[,which1])[fl1] - - m.inf <- subset(m.inf, - group == f2 & - group2 == f3)[,3] - - rownames(m.inf) <- paste(f3, - f2, - f1, - sep='/') - - ord <- order(m, - decreasing=TRUE) - - m.inf <-cbind(m.inf[,1][ord], - m.inf[,2][ord], - m.inf[,3][ord]) - - tab <- mt$tables[natri][[natri]] - } else if(nf1 == nt3 && nf2 == nt2) { - r <- mt$n[names(mt$tables)][[natri]] # groups and its number of replicates - - bal <- ifelse(length(r) == 1, - TRUE, - FALSE) # is (or not) balanced - - m <- as.vector(mt$tables[natri][[natri]][,fl2,fl1]) - - which1 <- names(dimnames(mt$tables[natri][[natri]]))[1] - - which2 <- names(dimnames(mt$tables[natri][[natri]]))[2] - - which3 <- names(dimnames(mt$tables[natri][[natri]]))[3] - - m.inf <- m.inf.3b(x, - which1=which1, - which2=which2, - which3=which3, - dispersion) - - f1 <- levels(model.frame(x)[,which1]) - - f2 <- levels(model.frame(x)[,which2])[fl2] - - f3 <- levels(model.frame(x)[,which3])[fl1] - - m.inf <- subset(m.inf, - group == f2 & - group2 == f3)[,3] - - rownames(m.inf) <- paste(f3, - f2, - f1, - sep='/') - - ord <- order(m, - decreasing=TRUE) - - m.inf <- cbind(m.inf[,1][ord], - m.inf[,2][ord], - m.inf[,3][ord]) - - tab <- mt$tables[natri][[natri]] - } else { - r <- mt$n[names(mt$tables)][[natri]] # groups and its number of replicates - - bal <- ifelse(length(r) == 1, - TRUE, - FALSE) # is (or not) balanced - - m <- as.vector(mt$tables[natri][[natri]][fl2,,fl1]) - - which1 <- names(dimnames(mt$tables[natri][[natri]]))[1] - - which2 <- names(dimnames(mt$tables[natri][[natri]]))[2] + split = ':'))[3] # nome do terceiro fator do which - which3 <- names(dimnames(mt$tables[natri][[natri]]))[3] + if(is.null(fl2)){ + # Interesse apenas na interação dupla + f1 <- levels(model.frame(x)[,nf2]) # correspondem aos fatores que se quer comparar! - m.inf <- m.inf.3b(x, - which1=which2, - which2=which1, - which3=which3, - dispersion) + f2 <- levels(model.frame(x)[,nf1])[fl1] # corresponde ao fator onde se está fazendo o desdobramento! - f1 <- levels(model.frame(x)[,which2]) + mt <- subset(aux_mt3, + eval(parse(text = nf1)) == f2) # pegando as médias de interesse - f2 <- levels(model.frame(x)[,which1])[fl2] + row.names(mt) <- paste(f2, + f1, + sep='/') + } # Interesse na interação tripla + else { - f3 <- levels(model.frame(x)[,which3])[fl1] + f1 <- levels(model.frame(x)[,nf3]) - m.inf <- subset(m.inf, - group == f2 & - group2 == f3)[,3] + f2 <- levels(model.frame(x)[,nf2])[fl2] - rownames(m.inf) <- paste(f3, - f2, - f1, - sep='/') + f3 <- levels(model.frame(x)[,nf1])[fl1] - ord <- order(m, - decreasing=TRUE) + mt <- subset(aux_mt3, + eval(parse(text = nf1)) == f3 & eval(parse(text=nf2)) == f2) # pegando as médias de interesse - m.inf <-cbind(m.inf[,1][ord], - m.inf[,2][ord], - m.inf[,3][ord]) + row.names(mt) <- paste(f3, + f2, + f1, + sep='/') - tab <- mt$tables[natri][[natri]] - } - } - } - colnames(m.inf) <- c('mean', - 'min', - 'max') - mMSE <- MSE / r + } - dfr <- x[[error]]$df.residual # residual degrees of freedom + out <- make.TukeyC.test(obj = mt, + MSE = MSE, + sig.level = sig.level, + dfr = dfr, + round = round, + adjusted.pvalue = adjusted.pvalue) - out <- make.TukeyC.test(r=r, - MSE=MSE, - m.inf=m.inf, - ord=ord, - sig.level=sig.level, - dfr=dfr, - bal=bal, - mt=mt, - round) + m.inf <- m.infos.nest.aovlist(x = x, + my = my, + forminter = forminter, + which = which, + fl1 = fl1, + fl2 = fl2, + sig.level = sig.level, + aux_mt = aux_mt1, + MSE = MSE) - class(out) <- c('TukeyC.nest', - 'TukeyC', - 'list') + res <- list(out = out, + info = m.inf) - invisible(out) } diff --git a/R/TukeyC.nest.default.R b/R/TukeyC.nest.default.R deleted file mode 100644 index 2e4a50a..0000000 --- a/R/TukeyC.nest.default.R +++ /dev/null @@ -1,51 +0,0 @@ -## -## S3 method to design matrix and response variable or data.frame objects -## - -TukeyC.nest.default <- function(x, - y=NULL, - model, - which, - error, - fl1, - fl2=0, - sig.level=.05, - round=2, - dispersion=c('mm', 's', 'se'), ...) -{ - - if (is.data.frame(y)) - y <- as.matrix(y[, 1]) # manova is not contemplated - else - stopifnot(is.atomic(y)) - - if (is.matrix(x) || is.atomic(x)) - x <- as.data.frame(x) - - if(!is.null(y)) - dat <- as.data.frame(cbind(x, - y)) - else - dat <- x - - av <- eval(substitute(aov(fo, - x), - list(fo=formula(model)))) - - if(class(av)[1] == 'aov') - res <- TukeyC.nest.aov(x=av, - which=which, - fl1=fl1, - fl2=fl2, - sig.level=sig.level, - round) - else - res <- TukeyC.nest.aovlist(x=av, - which=which, - error=error, - fl1=fl1, - fl2=fl2, - sig.level=sig.level, - round) - invisible(res) -} diff --git a/R/TukeyC.nest.lm.R b/R/TukeyC.nest.lm.R new file mode 100644 index 0000000..9382e23 --- /dev/null +++ b/R/TukeyC.nest.lm.R @@ -0,0 +1,110 @@ +TukeyC.nest.lm <- function(x, + which, + fl1, + fl2, + MSE, + dfr, + sig.level, + round, + adjusted.pvalue, ...) +{ + + my <- as.character(formula(x)[[2]]) + #my <- as.character(attr(x,'terms')[[2]]) + m1 <- gsub('\\:','\\+', which) + m2 <- unlist(strsplit(which, + '[[:punct:]]')) + + forminter <- as.formula(paste(my, '~', m1)) + + # aux_mt1 <- aggregate(forminter, + # data = x$model, + # function(x) c(means = mean(x), + # r = length(x))) + # + # aux_mt2 <- aux_mt1[order(aux_mt1[[my]][,1], + # decreasing = TRUE),] + # + # names(aux_mt2) <- gsub(my,'x',names(aux_mt2)) + # + + aux_r <- aggregate(forminter, + data = x$model, + function(x) r = length(x)) + reps <- aux_r[[my]] + + aux_mt <- suppressWarnings(LSmeans(x, + effect = m2)) + + aux_mt1 <- aux_mt$coef[,1] + + aux_mt2 <- data.frame(aux_r[1:length(names(aux_r))-1], + means = aux_mt1, + reps = reps) + + aux_mt3 <- aux_mt2[order(aux_mt2[['means']], + decreasing = TRUE),] + + nf1 <- unlist(strsplit(which, + split = ':'))[1] # nome do primeiro fator do which + + nf2 <- unlist(strsplit(which, + split = ':'))[2] # nome do segundo fator do which + + nf3 <- unlist(strsplit(which, + split = ':'))[3] # nome do terceiro fator do which + + if(is.null(fl2)){ + # Interesse apenas na interação dupla + + f1 <- levels(x$model[,nf2]) # correspondem aos fatores que se quer comparar! + + f2 <- levels(x$model[,nf1])[fl1] # corresponde ao fator onde se está fazendo o desdobramento! + + # mt <- subset(aux_mt2, + # eval(parse(text = nf1)) == f2) # pegando as médias de interesse + mt <- subset(aux_mt3, + eval(parse(text = nf1)) == f2) # pegando as médias de interesse + + row.names(mt) <- paste(mt[,1], + mt[,2], + sep='/') + } # Interesse na interação tripla + else { + + f1 <- levels(x$model[,nf3]) + + f2 <- levels(x$model[,nf2])[fl2] + + f3 <- levels(x$model[,nf1])[fl1] + + mt <- subset(aux_mt3, + eval(parse(text = nf1)) == f3 & eval(parse(text=nf2)) == f2) # pegando as médias de interesse + + row.names(mt) <- paste(mt[,1], + mt[,2], + mt[,3], + sep='/') + + } + + out <- make.TukeyC.test(obj = mt, + MSE = MSE, + sig.level = sig.level, + dfr = dfr, + round = round, + adjusted.pvalue = adjusted.pvalue) + + m.inf <- m.infos.nest.lm(x = x, + my = my, + forminter = forminter, + which = which, + fl1 = fl1, + fl2 = fl2, + sig.level = sig.level, + aux_mt = aux_mt) + + res <- list(out = out, + info = m.inf) + +} diff --git a/R/m.inf.1a.R b/R/m.inf.1a.R deleted file mode 100644 index 492c56b..0000000 --- a/R/m.inf.1a.R +++ /dev/null @@ -1,25 +0,0 @@ -m.inf.1a <- function(x, - which, - dispersion=c('mm', 's', 'se')) -{ - switch(match.arg(dispersion), - mm = { - m.inf <- aggregate(x$model[,1], - by=list(x$model[[which]]), - function(x) c(mean=mean(x), - min=min(x), - max=max(x)))[,2] - }, s = { - m.inf <- aggregate(x$model[,1], - by=list(x$model[[which]]), - function(x) c(mean=mean(x), - 'm - s'=mean(x) - sd(x), - 'm + s'=mean(x) + sd(x)))[,2] - }, se= { - m.inf <- aggregate(x$model[,1], - by=list(x$model[[which]]), - function(x) c(mean=mean(x), - 'm - se'=mean(x) - (sd(x) / sqrt(length(x))), - 'm + se'=mean(x) + (sd(x) / sqrt(length(x)))))[,2] - }) -} diff --git a/R/m.inf.1b.R b/R/m.inf.1b.R deleted file mode 100644 index 1422b56..0000000 --- a/R/m.inf.1b.R +++ /dev/null @@ -1,25 +0,0 @@ -m.inf.1b <- function(x, - which, - dispersion=c('mm', 's', 'se')) -{ - switch(match.arg(dispersion), - mm = { - m.inf <- aggregate(model.frame(x)[,1], - by=list(model.frame(x)[[which]]), - function(x) c(mean=mean(x), - min=min(x), - max=max(x)))[,2] - }, s = { - m.inf <- aggregate(model.frame(x)[,1], - by=list(model.frame(x)[[which]]), - function(x) c(mean=mean(x), - 'm - s'=mean(x) - sd(x), - 'm + s'=mean(x) + sd(x)))[,2] - }, se= { - m.inf <- aggregate(model.frame(x)[,1], - by=list(model.frame(x)[[which]]), - function(x) c(mean=mean(x), - 'm - se'=mean(x) - (sd(x) / sqrt(length(x))), - 'm + se'=mean(x) + (sd(x) / sqrt(length(x)))))[,2] - }) -} diff --git a/R/m.inf.2a.R b/R/m.inf.2a.R deleted file mode 100644 index 47cf7b3..0000000 --- a/R/m.inf.2a.R +++ /dev/null @@ -1,29 +0,0 @@ -m.inf.2a <- function(x, - which1, - which2, - dispersion=c('mm', 's', 'se')) -{ - switch(match.arg(dispersion), - mm = { - m.inf <- aggregate(x$model[,1], - by=list(x$model[[which2]], - group=x$model[[which1]]), - function(x) c(mean=mean(x), - min=min(x), - max=max(x)))[,2:3] - }, s = { - m.inf <- aggregate(x$model[,1], - by=list(x$model[[which2]], - group=x$model[[which1]]), - function(x) c(mean=mean(x), - 'm - s'=mean(x) - sd(x), - 'm + s'=mean(x) + sd(x)))[,2:3] - }, se= { - m.inf <- aggregate(x$model[,1], - by=list(x$model[[which2]], - group=x$model[[which1]]), - function(x) c(mean=mean(x), - se.min=mean(x) - (sd(x) / sqrt(length(x))), - se.max=mean(x) + (sd(x) / sqrt(length(x)))))[,2:3] - }) -} diff --git a/R/m.inf.2b.R b/R/m.inf.2b.R deleted file mode 100644 index 000981d..0000000 --- a/R/m.inf.2b.R +++ /dev/null @@ -1,29 +0,0 @@ -m.inf.2b <- function(x, - which1, - which2, - dispersion=c('mm', 's', 'se')) -{ - switch(match.arg(dispersion), - mm = { - m.inf <- aggregate(model.frame(x)[,1], - by=list(model.frame(x)[[which2]], - group=model.frame(x)[[which1]]), - function(x) c(mean=mean(x), - min=min(x), - max=max(x)))[,2:3] - }, s = { - m.inf <- aggregate(model.frame(x)[,1], - by=list(model.frame(x)[[which2]], - group=model.frame(x)[[which1]]), - function(x) c(mean=mean(x), - 'm - s'=mean(x) - sd(x), - 'm + s'=mean(x) + sd(x)))[,2:3] - }, se= { - m.inf <- aggregate(model.frame(x)[,1], - by=list(model.frame(x)[[which2]], - group=model.frame(x)[[which1]]), - function(x) c(mean=mean(x), - se.min=mean(x) - (sd(x) / sqrt(length(x))), - se.max=mean(x) + (sd(x) / sqrt(length(x)))))[,2:3] - }) -} diff --git a/R/m.inf.3a.R b/R/m.inf.3a.R deleted file mode 100644 index 194622a..0000000 --- a/R/m.inf.3a.R +++ /dev/null @@ -1,33 +0,0 @@ -m.inf.3a <- function(x, - which1, - which2, - which3, - dispersion=c('mm', 's', 'se')) -{ - switch(match.arg(dispersion), - mm = { - m.inf <- aggregate(x$model[,1], - by=list(x$model[[which1]], - group=x$model[[which2]], - group2=x$model[[which3]]), - function(x) c(mean=mean(x), - min=min(x), - max=max(x)))[,2:4] - }, s = { - m.inf <- aggregate(x$model[,1], - by=list(x$model[[which1]], - group=x$model[[which2]], - group2=x$model[[which3]]), - function(x) c(mean=mean(x), - 'm - s'=mean(x) - sd(x), - 'm + s'=mean(x) + sd(x)))[,2:4] - }, se= { - m.inf <- aggregate(x$model[,1], - by=list(x$model[[which1]], - group=x$model[[which2]], - group2=x$model[[which3]]), - function(x) c(mean=mean(x), - se.min=mean(x) - (sd(x) / sqrt(length(x))), - se.max=mean(x) + (sd(x) / sqrt(length(x)))))[,2:4] - }) -} diff --git a/R/m.inf.3b.R b/R/m.inf.3b.R deleted file mode 100644 index 6ffe1d8..0000000 --- a/R/m.inf.3b.R +++ /dev/null @@ -1,33 +0,0 @@ -m.inf.3b <- function(x, - which1, - which2, - which3, - dispersion=c('mm', 's', 'se')) -{ - switch(match.arg(dispersion), - mm = { - m.inf <- aggregate(model.frame(x)[,1], - by=list(model.frame(x)[[which1]], - group=model.frame(x)[[which2]], - group2=model.frame(x)[[which3]]), - function(x) c(mean=mean(x), - min=min(x), - max=max(x)))[,2:4] - }, s = { - m.inf <- aggregate(model.frame(x)[,1], - by=list(model.frame(x)[[which1]], - group=model.frame(x)[[which2]], - group2=model.frame(x)[[which3]]), - function(x) c(mean=mean(x), - 'm - s'=mean(x) - sd(x), - 'm + s'=mean(x) + sd(x)))[,2:4] - }, se= { - m.inf <- aggregate(model.frame(x)[,1], - by=list(model.frame(x)[[which1]], - group=model.frame(x)[[which2]], - group2=model.frame(x)[[which3]]), - function(x) c(mean=mean(x), - se.min=mean(x) - (sd(x) / sqrt(length(x))), - se.max=mean(x) + (sd(x) / sqrt(length(x)))))[,2:4] - }) -} diff --git a/R/m.infos.aovlist.R b/R/m.infos.aovlist.R new file mode 100644 index 0000000..0fc58a6 --- /dev/null +++ b/R/m.infos.aovlist.R @@ -0,0 +1,33 @@ +m.infos.aovlist <- function(x, + my, + forminter, + which, + sig.level, + aux_mt, + MSE, + ...) +{ + + aux_m.inf <- aggregate(forminter, + data = model.frame(x), + function(x) c(min = min(x), + max = max(x), + sd = sd(x), + se = sqrt(MSE/length(x)))) + + aux_m.inf1 <- data.frame(aux_m.inf[names(aux_m.inf)!=my], + means = aux_mt[[my]][,1], + aux_m.inf[[my]][,1:2], + 'linf_sd' = aux_mt[[my]][,1] - aux_m.inf[[my]][,3], + 'lsup_sd' = aux_mt[[my]][,1] + aux_m.inf[[my]][,3], + 'linf_se' = aux_mt[[my]][,1] - abs(qt(sig.level,aux_mt[[my]][,2]))*aux_m.inf[[my]][,4], + 'lsup_se' = aux_mt[[my]][,1] + abs(qt(sig.level,aux_mt[[my]][,2]))*aux_m.inf[[my]][,4]) + + aux_m.inf2 <- aux_m.inf1[order(aux_m.inf1[['means']], + decreasing = TRUE),] + + m.inf <- list(Means = aux_m.inf2[,c(1:2)], + mm = aux_m.inf2[,c(1,3:4)], + sd = aux_m.inf2[,c(1,5:6)], + ic = aux_m.inf2[,c(1,7:8)]) +} diff --git a/R/m.infos.lm.R b/R/m.infos.lm.R new file mode 100644 index 0000000..6861957 --- /dev/null +++ b/R/m.infos.lm.R @@ -0,0 +1,31 @@ +m.infos.lm <- function(x, + my, + forminter, + which, + sig.level, + aux_mt, + ...) +{ + + aux_m.inf <- aggregate(forminter, + data = x$model, + function(x) c(min = min(x), + max = max(x), + sd = sd(x))) + + aux_m.inf1 <- data.frame(groups = aux_m.inf[names(aux_m.inf)!=my], + means = aux_mt$coef[,1], + aux_m.inf[[my]][,1:2], + 'linf_sd' = aux_mt$coef[,1] - aux_m.inf[[my]][,3], + 'lsup_sd' = aux_mt$coef[,1] + aux_m.inf[[my]][,3], + 'linf_se' = aux_mt$coef[,1] - abs(qt(sig.level,aux_mt$coef[,3]))*aux_mt$coef[,2], + 'lsup_se' = aux_mt$coef[,1] + abs(qt(sig.level,aux_mt$coef[,3]))*aux_mt$coef[,2]) + + aux_m.inf2 <- aux_m.inf1[order(aux_m.inf1[['means']], + decreasing = TRUE),] + + m.inf <- list(Means = aux_m.inf2[,c(1:2)], + mm = aux_m.inf2[,c(1,3:4)], + sd = aux_m.inf2[,c(1,5:6)], + ic = aux_m.inf2[,c(1,7:8)]) +} diff --git a/R/m.infos.nest.aovlist.R b/R/m.infos.nest.aovlist.R new file mode 100644 index 0000000..53f5692 --- /dev/null +++ b/R/m.infos.nest.aovlist.R @@ -0,0 +1,67 @@ +m.infos.nest.aovlist <- function(x, + my, + forminter, + which, + fl1, + fl2, + sig.level, + aux_mt, + MSE, + ...) +{ + + aux_m.inf <- aggregate(forminter, + data = model.frame(x), + function(x) c(min = min(x), + max = max(x), + sd = sd(x), + se = sqrt(MSE/length(x)))) + + aux_m.inf1 <- data.frame(aux_m.inf[names(aux_m.inf)!=my], + means = aux_mt[[my]][,1], + aux_m.inf[[my]][,1:2], + 'linf_sd' = aux_mt[[my]][,1] - aux_m.inf[[my]][,3], + 'lsup_sd' = aux_mt[[my]][,1] + aux_m.inf[[my]][,3], + 'linf_se' = aux_mt[[my]][,1] - abs(qt(sig.level,aux_mt[[my]][,2]))*aux_m.inf[[my]][,4], + 'lsup_se' = aux_mt[[my]][,1] + abs(qt(sig.level,aux_mt[[my]][,2]))*aux_m.inf[[my]][,4]) + + aux_m.inf2 <- aux_m.inf1[order(aux_m.inf1[['means']], + decreasing = TRUE),] + + nf1 <- unlist(strsplit(which, + split = ':'))[1] # nome do primeiro fator do which + + nf2 <- unlist(strsplit(which, + split = ':'))[2] # nome do segundo fator do which + + nf3 <- unlist(strsplit(which, + split = ':'))[3] # nome do terceiro fator do which + + if(is.null(fl2)){ + + + f2 <- levels(model.frame(x)[,nf1])[fl1] # corresponde ao fator onde se está fazendo o desdobramento! + + aux_m.inf21 <- subset(aux_m.inf2, + eval(parse(text = nf1)) == f2) # pegando as médias de interesse + + m.inf <- list(Means = aux_m.inf21[,c(1:3)], + mm = aux_m.inf21[,c(1:2,4:5)], + sd = aux_m.inf21[,c(1:2,6:7)], + ic = aux_m.inf21[,c(1:2,8:9)]) + + } else { + + f2 <- levels(model.frame(x)[,nf2])[fl2] + + f3 <- levels(model.frame(x)[,nf1])[fl1] + + aux_m.inf21 <- subset(aux_m.inf2, + eval(parse(text = nf1)) == f3 & eval(parse(text = nf2)) == f2) # pegando as médias de interesse + + m.inf <- list(Means = aux_m.inf21[,c(1:4)], + mm = aux_m.inf21[,c(1:3,5:6)], + sd = aux_m.inf21[,c(1:3,7:8)], + ic = aux_m.inf21[,c(1:3,9:10)]) + } +} diff --git a/R/m.infos.nest.lm.R b/R/m.infos.nest.lm.R new file mode 100644 index 0000000..b9ad6e3 --- /dev/null +++ b/R/m.infos.nest.lm.R @@ -0,0 +1,65 @@ +m.infos.nest.lm <- function(x, + my, + forminter, + which, + fl1, + fl2, + sig.level, + aux_mt, + ...) +{ + + aux_m.inf <- aggregate(forminter, + data = x$model, + function(x) c(min=min(x), + max=max(x), + sd=sd(x))) + + aux_m.inf1 <- data.frame(aux_m.inf[names(aux_m.inf) != my], + means = aux_mt$coef[,1], + aux_m.inf[[my]][,1:2], + 'linf_sd' = aux_mt$coef[,1] - aux_m.inf[[my]][,3], + 'lsup_sd' = aux_mt$coef[,1] + aux_m.inf[[my]][,3], + 'linf_se' = aux_mt$coef[,1] - abs(qt(sig.level,aux_mt$coef[,3]))*aux_mt$coef[,2], + 'lsup_se' = aux_mt$coef[,1] + abs(qt(sig.level,aux_mt$coef[,3]))*aux_mt$coef[,2]) + + aux_m.inf2 <- aux_m.inf1[order(aux_m.inf1[['means']], + decreasing = TRUE),] + + nf1 <- unlist(strsplit(which, + split = ':'))[1] # nome do primeiro fator do which + + nf2 <- unlist(strsplit(which, + split = ':'))[2] # nome do segundo fator do which + + nf3 <- unlist(strsplit(which, + split = ':'))[3] # nome do terceiro fator do which + + if(is.null(fl2)){ + + f2 <- levels(x$model[,nf1])[fl1] # corresponde ao fator onde se está fazendo o desdobramento! + + aux_m.inf21 <- subset(aux_m.inf2, + eval(parse(text = nf1)) == f2) # pegando as médias de interesse + + m.inf <- list(Means = aux_m.inf21[,c(1:3)], + mm = aux_m.inf21[,c(1:2,4:5)], + sd = aux_m.inf21[,c(1:2,6:7)], + ic = aux_m.inf21[,c(1:2,8:9)]) + + } else { + + f2 <- levels(x$model[,nf2])[fl2] + + f3 <- levels(x$model[,nf1])[fl1] + + aux_m.inf21 <- subset(aux_m.inf2, + eval(parse(text = nf1)) == f3 & eval(parse(text = nf2)) == f2) # pegando as médias de interesse + + m.inf <- list(Means = aux_m.inf21[,c(1:4)], + mm = aux_m.inf21[,c(1:3,5:6)], + sd = aux_m.inf21[,c(1:3,7:8)], + ic = aux_m.inf21[,c(1:3,9:10)]) + + } +} diff --git a/R/make.TukeyC.test.R b/R/make.TukeyC.test.R index 321f3d0..404934d 100644 --- a/R/make.TukeyC.test.R +++ b/R/make.TukeyC.test.R @@ -1,88 +1,91 @@ ## ## Function to perform Tukey test ## - -make.TukeyC.test <- function(r=r, - MSE=MSE, - m.inf=m.inf, - ord=ord, - sig.level=sig.level, - dfr=dfr, - bal=bal, - mt=mt, - round=round) +make.TukeyC.test <- function(obj, + MSE, + sig.level, + dfr, + round, + adjusted.pvalue) { - if (length(r) < nrow(m.inf)) # expand n to the correct length if necessary - r <- rep.int(r, - nrow(m.inf)) - - r <- r[ord] - m.tmp <- m.inf[, 1] + #vari <- names(obj)[length(obj)] + #means <- obj[[vari]][,1] + means <- obj[['means']] + #names(means) <- obj[[vari]][,2] # necessary to vece object + names(means) <- obj[['reps']] - names(m.tmp) <- r + vece <- outer(X = means, + Y = means, # (v)ariance (e)stimation of (c)ontrast (e)stimation + function(X, Y) MSE * (1/as.numeric(names(X)) + 1/(as.numeric(names(Y))))) - vece <- outer(X=m.tmp, - Y=m.tmp, # (v)ariance (e)stimation of (c)ontrast (e)stimation - function(X, Y) MSE * (1/as.numeric(names(X)) + (1/as.numeric(names(Y))))) + qTukey <- qtukey(p = sig.level, + nmeans = dim(obj)[1], + df = dfr, + lower.tail = FALSE) - qTukey <- qtukey(p=sig.level, - nmeans=nrow(m.inf), - df=dfr, - lower.tail=FALSE) + msd <- qTukey * sqrt(1/2 * vece) # minimum significative difference - if (!bal) { - msd <- qTukey * sqrt(1/2 * vece) # minimum significative difference + diag(msd) <- 0 - diag(msd) <- 0 + dimnames(msd) <- list(row.names(obj), + row.names(obj)) - dimnames(msd) <- list(rownames(m.inf), - rownames(m.inf)) - } - else - msd <- qTukey * sqrt(1/2 * vece)[1,1] + names(means) <- row.names(obj) # necessary to difm object - m <- m.inf[,1] - - difm <- abs(outer(m, - m, + difm <- abs(outer(means, + means, "-")) # means difference + dif <- difm - msd - dif <- ifelse(dif <= 0, - FALSE, - TRUE) + # dif <- ifelse(dif <= 0, + # FALSE, + # TRUE) + # + + # The below estimates the probability of observed difference betweeen means be significative + # Matrix of the difference of means above diagonal and respective p-values of the Tukey test + # below diagonal + pval_tukey <- ptukey(q = difm[lower.tri(difm)] / sqrt(1/2 * vece[lower.tri(vece)]), + nmeans = dim(obj)[1], + df = dfr, + lower.tail = FALSE) + + adjusted_pvalue <- p.adjust(pval_tukey, + method = adjusted.pvalue) + + new_dif <- difm + new_dif[lower.tri(new_dif)] <- adjusted_pvalue + new_dif1 <- t(new_dif) + new_dif1[lower.tri(new_dif1)] <- adjusted_pvalue + diag(new_dif1) <- 1 + + new_dif2 <- ifelse(new_dif1 >= sig.level, + FALSE, + TRUE) - res <- make.TukeyC.groups(dif) + out_groups <- make.TukeyC.groups(new_dif2) - res <- cbind(format(round(m, + res <- cbind(format(round(means, round), - nsmall=2), - res) + nsmall = 2), + out_groups) colnames(res) <- c('Means', paste('G', 1:(ncol(res) - 1), sep='')) - if (bal) r <- r[1] + difm[lower.tri(difm)] <- adjusted_pvalue - # The below estimates the probability of observed difference betweeen means be significative - # Matrix of the difference of means above diagonal and respective p-values of the Tukey test - # below diagonal - difm[lower.tri(difm)] <- ptukey(q=difm[lower.tri(difm)] / sqrt(1/2 * vece[lower.tri(vece)]), - nmeans=nrow(m.inf), - df=dfr, - lower.tail=FALSE) diag(difm) <- 0 # To be sure! - out <- list(Table = mt, - Means = m.inf, - Result = as.data.frame(res), + out <- list(Result = as.data.frame(res), Sig.level = sig.level, Diff_Prob = round(difm, 3), MSD = round(msd, 3), - Replicates = r) + Replicates = obj[['reps']]) return(out) } diff --git a/R/plot.TukeyC.R b/R/plot.TukeyC.R index bf78cd1..6ba123b 100644 --- a/R/plot.TukeyC.R +++ b/R/plot.TukeyC.R @@ -3,22 +3,25 @@ ## plot.TukeyC <- function(x, - result=TRUE, - replicates=TRUE, - pch=19, - col=NULL, - xlab=NULL, - ylab=NULL, - xlim=NULL, - ylim=NULL, - id.lab=NULL, - id.las=1, - rl=TRUE, - rl.lty=3, - rl.col='gray', - mm=TRUE, - mm.lty=1, - title='', ...) + result = TRUE, + replicates = TRUE, + pch = 19, + col = NULL, + xlab = NULL, + ylab = NULL, + xlim = NULL, + ylim = NULL, + id.lab = NULL, + id.las = 1, + rl = TRUE, + rl.lty = 3, + rl.col = 'gray', + dispersion = c('none', + 'mm', + 'sd', + 'ci'), + dispersion.lty = 1, + title = '', ...) { fun <- function(m) { a <- rep('\n', length(m)) @@ -34,14 +37,14 @@ plot.TukeyC <- function(x, if(is.null(ylab)) ylab <- 'Means' - means <- x$Means[, 1] + means <- x$info$Means[['means']] if(replicates) - r <- x$Replicates + r <- x$out$Replicates groups <- 1:length(means) - m.res <- t(x$Result[, 2:ncol(x$Result)]) + m.res <- t(x$out$Result[, 2:ncol(x$out$Result)]) if(dim(m.res)[1] != 1) { m.res <- apply(m.res, 2, fun) @@ -53,7 +56,21 @@ plot.TukeyC <- function(x, else id.groups <- m.res - minmax <- x$Means[, 2:3] + #minmax <- x$Means[, 2:3] + minmax1 <- x$info$mm[['min']] + minmax2 <- x$info$mm[['max']] + minmax <- data.frame(minmax1, + minmax2) + + sd1 <- x$info$sd[['linf_sd']] + sd2 <- x$info$sd[['lsup_sd']] + sdd <- data.frame(sd1, + sd2) + + ic1 <- x$info$ic[['linf_se']] + ic2 <- x$info$ic[['lsup_se']] + ic <- data.frame(ic1, + ic2) if(is.null(col)) col <- 'black' @@ -70,17 +87,24 @@ plot.TukeyC <- function(x, ylim <- c(min(minmax[, 1]), max(minmax[, 2])) - old.par <- par(mai=c(1, 1, 1.25, 1)) + # By J.C.Faria + ngroups <- dim(x$out$Result)[2] - 1 + if(ngroups > 3){ + op <- par('mar') # Original par('mar') + np <- op # A copy + np[3] <- ngroups + 1 # Changing top to show all letters + par(mar=np) # Setting new par('mar') + } plot(groups, means, - pch=pch, - col=col, - xlab=xlab, - ylab=ylab, - xlim=xlim, - ylim=ylim, - axes=FALSE, ...) + pch = pch, + col = col, + xlab = xlab, + ylab = ylab, + xlim = xlim, + ylim = ylim, + axes = FALSE, ...) if(rl == TRUE) segments(rep(-0.5, @@ -88,51 +112,79 @@ plot.TukeyC <- function(x, means, groups, means, - lty=rl.lty, - col=rl.col, ...) - - if(mm == TRUE) - segments(groups, - minmax[, 2], - groups, - minmax[, 1], - lty=mm.lty, - col=col, ...) + lty = rl.lty, + col = rl.col, ...) + + # if(mm == TRUE) + # segments(groups, + # minmax[, 2], + # groups, + # minmax[, 1], + # lty=mm.lty, + # col=col, ...) + switch(match.arg(dispersion), + mm = { + segments(groups, + minmax[,1], + groups, + minmax[,2], + lty = dispersion.lty, + col = col, ...) + }, + sd = { + segments(groups, + sdd[,1], + groups, + sdd[,2], + lty = dispersion.lty, + col = col, ...) + }, + ci = { + segments(groups, + ic[,1], + groups, + ic[,2], + lty = dispersion.lty, + col = col, ...) + }, + none = NULL) axis(2, - at=round(seq(ylim[1], - ylim[2], - length.out=5), - 1)) + at = round(seq(ylim[1], + ylim[2], + length.out = 5), + 1)) if(is.null(id.lab)) - id.lab <- rownames(x$Means) + id.lab <- names(x$out$Result[,1]) axis(1, - at=1:length(means), - labels=id.lab, - las=id.las, - col.axis=FALSE, ...) + at = 1:length(means), + labels = id.lab, + las = id.las, + col.axis = FALSE, ...) if(result) axis(3, - at=1:length(means), - labels=id.groups, ...) + at = 1:length(means), + labels = id.groups, ...) if(replicates) - text(x=1:length(means), - y=min(ylim), - labels=r, - pos=3, ...) + text(x = 1:length(means), + y = min(ylim), + labels = r, + pos = 3, ...) - mtext(text=id.lab, - side=1, - line=1, - at=1:length(means), - las=id.las, ...) + mtext(text = id.lab, + side = 1, + line = 1, + at = 1:length(means), + las = id.las, ...) title(title, ...) - par(old.par) + # By J.C.Faria + if(ngroups > 3) + par(mar=op) # Restoring the original par('mar') } diff --git a/R/print.TukeyC.R b/R/print.TukeyC.R new file mode 100644 index 0000000..709be62 --- /dev/null +++ b/R/print.TukeyC.R @@ -0,0 +1,23 @@ +## +## S3 method to print 'TukeyC' object +## + +print.TukeyC <- function(x, + ...) +{ + cat('Results\n') + print(x$out$Result, + ...) + + cat('\nSig.level\n', + x$out$Sig.level) + + cat('\n\nDiff_Prob\n') + print(x$out$Diff_Prob, + ...) + + cat('\nMSD\n') + print(x$out$MSD, + ...) + +} diff --git a/R/summary.TukeyC.R b/R/summary.TukeyC.R index 4ed1038..5806a0b 100644 --- a/R/summary.TukeyC.R +++ b/R/summary.TukeyC.R @@ -11,19 +11,19 @@ summary.TukeyC <- function(object, if(complete){ cat('Goups of means at sig.level =', - object$Sig.level, + object$out$Sig.level, '\n') - print(object$Result) + print(object$out$Result) cat('\nMatrix of the difference of means above diagonal and') cat('\nrespective p-values of the Tukey test below diagonal values\n') - print(object$Diff_Prob) + print(object$out$Diff_Prob) } else { - res <- object$Result + res <- object$out$Result res } diff --git a/demo/CRD.R b/demo/CRD.R index 49ea58f..40d3871 100644 --- a/demo/CRD.R +++ b/demo/CRD.R @@ -3,67 +3,26 @@ ## ## The parameters can be: vectors, design matrix and the response variable, -## data.frame or aov +## data.frame, aov or lm. -## Example 1 +## Example 1: an small experiment library(TukeyC) data(CRD1) -## From: vectors x and y - balanced +## From: formula - balanced tk1 <- with(CRD1, - TukeyC(x=x, - y=y, - model='y ~ x', - which='x')) + TukeyC(y ~ x, + dfm)) tk1 summary(tk1) plot(tk1) -## From: vectors x and y - unbalanced -tk1u <- with(CRD1, - TukeyC(x=x[-1], - y=y[-1], - model='y ~ x', - which='x', - dispersion='s')) -tk1u -summary(tk1u) -plot(tk1u) - -## From: design matrix (dm) and response variable (y) - balanced -tk2 <- with(CRD1, - TukeyC(x=dm, - y=y, - model='y ~ x', - which='x', - dispersion='se')) -tk2 -summary(tk2) - -## From: design matrix (dm) and response variable (y) - unbalanced -tk2u <- with(CRD1, - TukeyC(x=dm[-1,], - y=y[-1], - model='y ~ x', - which='x')) -tk2u -summary(tk2u) - -## From: data.frame (dfm) - balanced -tk3 <- with(CRD1, - TukeyC(x=dfm, - model='y ~ x', - which='x')) -tk3 -summary(tk3) - -## From: data.frame (dfm) - unbalanced -tk3u <- with(CRD1, - TukeyC(x=dfm[-1,], - model='y ~ x', - which='x')) -tk3u -summary(tk3u) +## From: formula - unbalanced +utk1 <- with(CRD1, + TukeyC(y ~ x, + dfm[-1,])) +utk1 +summary(utk1) ## From: aov - balanced av1 <- with(CRD1, @@ -71,86 +30,36 @@ av1 <- with(CRD1, data=dfm)) summary(av1) -tk4 <- TukeyC(x=av1, - which='x') -tk4 -summary(tk4) - -## From: aov - unbalanced -av1u <- with(CRD1, - aov(y ~ x, - data=dfm[-1,])) -summary(av1u) +tk2 <- TukeyC(av1) +tk2 +summary(tk2) -tk4u <- TukeyC(x=av1u, - which='x') -tk4u -summary(tk4u) +## From: lm - unbalanced +ulm1 <- with(CRD1, + lm(y ~ x, + data=dfm[-1,])) +summary(ulm1) +utk2 <- TukeyC(ulm1) +utk2 +summary(utk2) -## Example 2 -library(TukeyC) +## Example 2: a lot of groups data(CRD2) -## From: vectors x and y - balanced -tk5 <- with(CRD2, - TukeyC(x=x, - y=y, - model='y ~ x', - which='x', round=3)) -tk5 -summary(tk5) -plot(tk5, - id.las=2, - rl=FALSE) - -## From: vectors x and y - unbalanced -tk5u <- with(CRD2, - TukeyC(x=x[-1], - y=y[-1], - model='y ~ x', - which='x', - round=3)) -tk5u -summary(tk5u) - -## From: design matrix (dm) and response variable (y) - balanced -tk6 <- with(CRD2, - TukeyC(x=dm, - y=y, - model='y ~ x', - which='x')) -tk6 -summary(tk6) - -## From: design matrix (dm) and response variable (y) - unbalanced -tk6u <- with(CRD2, - TukeyC(x=dm[-1,], - y=y[-1], - model='y ~ x', - which='x')) -tk6u -summary(tk6u) -plot(tk6u, +## From: data.frame (dfm) - balanced +tk3 <- with(CRD2, + TukeyC(y ~ x, + dfm)) +plot(tk3, id.las=2, rl=FALSE) -## From: data.frame (dfm) - balanced -tk7 <- with(CRD2, - TukeyC(x=dfm, - model='y ~ x', - which='x')) -tk7 -summary(tk7) - ## From: data.frame (dfm) - unbalanced -tk7u <- with(CRD2, - TukeyC(x=dfm[-1,], - model='y ~ x', - which='x')) -tk7u -summary(tk7u) -plot(tk7u, +utk3 <- with(CRD2, + TukeyC(y ~ x, + dfm[-1,])) +plot(utk3, id.las=2, rl=FALSE) @@ -160,21 +69,19 @@ av2 <- with(CRD2, data=dfm)) summary(av2) -tk8 <- TukeyC(x=av2, - which='x') -tk8 -summary(tk8) - -## From: aov - unbalanced -av2u <- with(CRD2, - aov(y ~ x, - data=dfm[-1,])) -summary(av2u) - -tk8u <- TukeyC(x=av2u, - which='x') -tk8u -summary(tk8u) -plot(tk8u, +tk4 <- TukeyC(av2) +plot(tk4, + id.las=2, + rl=FALSE) + +## From: lm - unbalanced +ulm2 <- with(CRD2, + lm(y ~ x, + data=dfm[-1,])) +summary(ulm2) + +utk8 <- TukeyC(ulm2) + +plot(utk8, id.las=2, rl=FALSE) diff --git a/demo/FE.R b/demo/FE.R index a005d37..7be75ab 100644 --- a/demo/FE.R +++ b/demo/FE.R @@ -2,229 +2,67 @@ ## Example: Factorial Experiment (FE) ## -## The parameters can be: design matrix and the response variable, -## data.frame or aov +## The parameters can be: vectors, design matrix and the response variable, +## data.frame, aov or lm. ## Note: Upper case for factors and lowercase for levels library(TukeyC) data(FE) -## From: design matrix (dm) and response variable (y) -## Main factor: N +## From: formula +## Nested: k1/p2/N +## The indices (1, 2, ...) are used to set the level of the factor tk1 <- with(FE, - TukeyC(x=dm, - y=y, - model='y ~ blk + N*P*K', - which='N')) + TukeyC(y ~ blk + N*P*K, + dfm, + which='K:P:N', + fl1=1, + fl2=2)) summary(tk1) -plot(tk1) -## Main factor: P +## Nested: k2/p1/N tk2 <- with(FE, - TukeyC(x=dm, - y=y, - model='y ~ blk + N*P*K', - which='P', - dispersion='s')) + TukeyC(y ~ blk + N*P*K, + dfm, + which='K:P:N', + fl1=2, + fl2=1)) summary(tk2) -plot(tk2) - -## Main factor: K -tk3 <- with(FE, - TukeyC(x=dm, - y=y, - model='y ~ blk + N*P*K', - which='K', - dispersion='se')) -summary(tk3) -plot(tk3) - -## Nested: p1/N -## Testing N inside of level one of P -ntk1 <- with(FE, - TukeyC.nest(x=dm, - y=y, - model='y ~ blk + N*P*K', - which='P:N', - fl1=1)) -summary(ntk1) -plot(ntk1) - -## Nested: p2/N -ntk2 <- with(FE, - TukeyC.nest(x=dm, - y=y, - model='y ~ blk + N*P*K', - which='P:N', - fl1=2)) -summary(ntk2) -plot(ntk2) - -## Nested: k1/N -ntk3 <- with(FE, - TukeyC.nest(x=dm, - y=y, - model='y ~ blk + N*P*K', - which='K:N', - fl1=1)) -summary(ntk3) - -## Nested: k2/N -ntk4 <- with(FE, - TukeyC.nest(x=dm, - y=y, - model='y ~ blk + N*P*K', - which='K:N', - fl1=2)) -summary(ntk4) - -## Nested: k1/P -ntk5 <- with(FE, - TukeyC.nest(x=dm, - y=y, - model='y ~ blk + N*P*K', - which='K:P', - fl1=1)) -summary(ntk5) - -## Nested: k2/P -ntk6 <- with(FE, - TukeyC.nest(x=dm, - y=y, - model='y ~ blk + N*P*K', - which='K:P', - fl1=2)) -summary(ntk6) - -## Nested: k1/p1/N -## Testing N inside of level one of K and level one of P -ntk7 <- with(FE, - TukeyC.nest(x=dm, - y=y, - model='y ~ blk + N*P*K', - which='K:P:N', - fl1=1, - fl2=1)) -summary(ntk7) # p-value=.053 to be different -plot(ntk7) - -ntk7 <- with(FE, - TukeyC.nest(x=dm, - y=y, - model='y ~ blk + N*P*K', - which='K:P:N', - fl1=1, - fl2=1, - sig.level=.054)) -summary(ntk7) # The p-value is very informative! -plot(ntk7) - -## Nested: k2/p2/N -ntk8 <- with(FE, - TukeyC.nest(x=dm, - y=y, - model='y ~ blk + N*P*K', - which='K:P:N', - fl1=2, - fl2=2)) -summary(ntk8) -plot(ntk8) - -## Nested: k1/n1/P -ntk9 <- with(FE, - TukeyC.nest(x=dm, - y=y, - model='y ~ blk + P*N*K', - which='K:N:P', - fl1=1, - fl2=1)) -summary(ntk9) - -## Nested: k2/n2/P -ntk8 <- with(FE, - TukeyC.nest(x=dm, - y=y, - model='y ~ blk + P*N*K', - which='K:N:P', - fl1=2, - fl2=2)) -summary(ntk8) - -## Nested: p1/n1/K -ntk10 <- with(FE, - TukeyC.nest(x=dm, - y=y, - model='y ~ blk + K*N*P', - which='P:N:K', - fl1=1, - fl2=1)) -summary(ntk10) - -## Nested: p2/n2/K -ntk11 <- with(FE, - TukeyC.nest(x=dm, - y=y, - model='y ~ blk + K*N*P', - which='P:N:K', - fl1=2, - fl2=2)) -summary(ntk11) - - -## From: data.frame -## Nested: k2/p1/N -ntk12 <- with(FE, - TukeyC.nest(x=dfm, - model='y ~ blk + N*P*K', - which='K:P:N', - fl1=1, - fl2=2)) -summary(ntk12) - -## Nested: k1/p2/N -ntk13 <- with(FE, - TukeyC.nest(x=dfm, - model='y ~ blk + N*P*K', - which='K:P:N', - fl1=2, - fl2=1)) -summary(ntk13) - ## From aov -nav1 <- with(FE, - aov(y ~ blk + N*P*K , - data=dfm)) -summary(nav1) +av1 <- with(FE, + aov(y ~ blk + N*P*K , + data=dfm)) +summary(av1) ## Main factor: N -ntk14 <- TukeyC(nav1, - which='N') -summary(ntk14) +tk3 <- TukeyC(av1, + which='N') +summary(tk3) ## Nested: k1/P -ntk15 <- TukeyC.nest(nav1, - which='K:P', - fl1=1) -summary(ntk15) +tk4 <- TukeyC(av1, + which='K:P', + fl1=1) +summary(tk4) -## Nested: k2/p1/N -ntk16 <- TukeyC.nest(nav1, - which='K:P:N', - fl1=1, - fl2=2) -summary(ntk16) - -# Changing the order of factors (test) -nav2 <- with(FE, - aov(y ~ blk + K*N* P, - data=dfm)) -summary(nav2) +## Nested: k1/p2/N +tk4 <- TukeyC(av1, + which='K:P:N', + fl1=1, + fl2=2) +summary(tk4) + +# Changing the order of factors (for test only) +av2 <- with(FE, + aov(y ~ blk + K*N*P, + data=dfm)) +summary(av2) ## Nested: p1/n1/K -ntk17 <- TukeyC.nest(nav2, - which='P:N:K', - fl1=1, - fl2=1) -summary(ntk17) +tk5 <- TukeyC(av2, + which='P:N:K', + fl1=1, + fl2=1) +summary(tk5) diff --git a/demo/LSD.R b/demo/LSD.R index 6380169..faca069 100644 --- a/demo/LSD.R +++ b/demo/LSD.R @@ -2,36 +2,27 @@ ## Example: Latin Squares Design (LSD) ## -## The parameters can be: design matrix and the response variable, -## data.frame or aov +## The parameters can be: vectors, design matrix and the response variable, +## data.frame, aov or lm. library(TukeyC) data(LSD) -## From: design matrix (dm) and response variable (y) +## From: formula +## Testing tra tk1 <- with(LSD, - TukeyC(x=dm, - y=y, - model='y ~ rows + cols + tra', - which='tra')) -summary(tk1) # p-value E-A = .051 -plot(tk1) - -tk1 <- with(LSD, - TukeyC(x=dm, - y=y, - model='y ~ rows + cols + tra', + TukeyC(y ~ rows + cols + tra, + dfm, which='tra', sig.level=.052)) summary(tk1) -plot(tk1) -## From: data.frame +## From: formula +## Testing rows tk2 <- with(LSD, - TukeyC(x=dfm, - model='y ~ rows + cols + tra', - which='tra', - sig.level=.052)) + TukeyC(y ~ rows + cols + tra, + dfm, + which='rows')) summary(tk2) ## From: aov @@ -40,7 +31,15 @@ av1 <- with(LSD, data=dfm)) summary(av1) +## From: aov +## Testing tra tk3 <- TukeyC(av1, which='tra', sig.level=.052) summary(tk3) + +## From: aov +## Testing cols +tk4 <- TukeyC(av1, + which='cols') +summary(tk4) diff --git a/demo/RCBD.R b/demo/RCBD.R index c4d2ab0..a8dcacb 100644 --- a/demo/RCBD.R +++ b/demo/RCBD.R @@ -2,52 +2,43 @@ ## Example: Randomized Complete Block Design (RCBD) ## -## The parameters can be: design matrix and the response variable, -## data.frame or aov +## The parameters can be: vectors, design matrix and the response variable, +## data.frame, aov or lm. library(TukeyC) data(RCBD) -## Design matrix (dm) and response variable (y) +## From: data.frame (dfm), which='tra' tk1 <- with(RCBD, - TukeyC(x=dm, - y=y, - model='y ~ blk + tra', + TukeyC(y ~ blk + tra, + dfm, which='tra')) summary(tk1) plot(tk1) -## From: data.frame (dfm), which='tra' +## From: formula, which='blk' implicit (due to be the first arg of the model) tk2 <- with(RCBD, - TukeyC(x=dfm, - model='y ~ blk + tra', - which='tra')) + TukeyC(y ~ blk + tra, + dfm)) summary(tk2) +plot(tk2) -## From: data.frame (dfm), which='blk' -tk3 <- with(RCBD, - TukeyC(x=dfm, - model='y ~ blk + tra', - which='blk')) -summary(tk3) -plot(tk3) -## From: aov av1 <- with(RCBD, aov(y ~ blk + tra, data=dfm)) summary(av1) -## From: aov, which='blk' implicit -tk4 <- TukeyC(x=av1) -summary(tk4) +## From: aov, which='blk' implicit (due to be the first arg of the model) +tk3 <- TukeyC(av1) +summary(tk3) ## From: aov, which='blk' explicit -tk5 <- TukeyC(x=av1, +tk4 <- TukeyC(x=av1, which='blk') -summary(tk5) +summary(tk4) ## From: aov, which='tra' explicit -tk6 <- TukeyC(x=av1, +tk5 <- TukeyC(x=av1, which='tra') -summary(tk6) +summary(tk5) diff --git a/demo/SPE.R b/demo/SPE.R new file mode 100644 index 0000000..ba8390d --- /dev/null +++ b/demo/SPE.R @@ -0,0 +1,101 @@ +## +## Example: Split-plot Experiment (SPE) +## + +## The parameters can be: vectors, design matrix and the response variable, +## data.frame, aov or lm. + +## Note: Upper case for factors and lowercase for levels + +library(TukeyC) +data(SPE) + +## From: formula +## Main factor: SP +tk1 <- with(SPE, + TukeyC(y ~ blk + P*SP + Error(blk/P), + data=dfm, + which='SP')) +summary(tk1) + +## Nested: p1/SP +tk2 <- with(SPE, + TukeyC(y ~ blk + P*SP + Error(blk/P), + data=dfm, + which='P:SP', + fl1=1 )) +summary(tk2) +plot(tk2) + +## Nested: sp1/P - it is necessary to inform how to combinate the errors +tk3 <- with(SPE, + TukeyC(y ~ blk + P*SP + Error(blk/P), + data=dfm, + which='SP:P', + error='Within/blk:P', + fl1=1)) +summary(tk3) + +## From: lm +lm1 <- with(SPE, + lm(y ~ blk*P + P*SP, + dfm)) + +tk4 <- TukeyC(lm1, + which='P:SP', + fl1=1) + +summary(tk4) + +tk5 <- TukeyC(lm1, + which='SP:P', + error='Within/blk:P', + fl1=1) +summary(tk5) + + +## From: aov +av1 <- with(SPE, + aov(y ~ blk + P*SP + Error(blk/P), + data=dfm)) +summary(av1) + +## Main factor: SP +tk6 <- TukeyC(av1, + which='SP', + sig.level=0.1) +summary(tk6) + +## Main factor: P +## It is necessary to inform the appropriate error for the test +tk7 <- TukeyC(av1, + which='P', + error='blk:P') + +summary(tk7) + +## Nested: p1/SP +## Testing SP inside of level one of P +tk8 <- TukeyC(av1, + which='P:SP', + fl1=1) +summary(tk8) + +## Nested: p2/SP +tk9 <- TukeyC(av1, + which='P:SP', + fl1=2) +summary(tk9) + +## Nested: p3/SP +tk10 <- TukeyC(av1, + which='P:SP', + fl1=3) +summary(tk10) + +## Nested: sp1/P - it is necessary to inform how to combinate the errors +tk11 <- TukeyC(av1, + which='SP:P', + error='Within/blk:P', + fl1=1) +summary(tk11) diff --git a/demo/SPE.r b/demo/SPE.r deleted file mode 100644 index 416280e..0000000 --- a/demo/SPE.r +++ /dev/null @@ -1,106 +0,0 @@ -## -## Example: Split-plot Experiment (SPE) -## - -## The parameters can be: design matrix and the response variable, -## data.frame or aov - -## Note: Upper case for factors and lowercase for levels - -library(TukeyC) -data(SPE) - -## From: design matrix (dm) and response variable (y) -## Main factor: P -tk1 <- with(SPE, - TukeyC(x=dm, - y=y, - model='y ~ blk + P*SP + Error(blk/P)', - which='P', - error='blk:P')) -summary(tk1) -plot(tk1) - -## Main factor: SP -tk2 <- with(SPE, - TukeyC(x=dm, - y=y, - model='y ~ blk + P*SP + Error(blk/P)', - which='SP', - error ='Within')) -summary(tk2) -plot(tk2) - -## Nested: p1/SP -## Testing SP inside of level one of P -tkn1 <- with(SPE, - TukeyC.nest(x=dm, - y=y, - model='y ~ blk + P*SP + Error(blk/P)', - which='P:SP', - error='Within', - fl1=1)) -summary(tkn1) - - -## From: data.frame -## Main factor: SP -tk3 <- with(SPE, - TukeyC(dfm, - model='y ~ blk + P*SP + Error(blk/P)', - which='SP', - error='Within', - sig.level=0.025)) -summary(tk3) - -## Nested: p1/SP -tkn2 <- with(SPE, - TukeyC.nest(dfm, - model='y ~ blk + P*SP + Error(blk/P)', - which='P:SP', - error='Within', - fl1=1 )) -summary(tkn2) -plot(tkn2) - -## From: aovlist -av1 <- with(SPE, - aov(y ~ blk + P*SP + Error(blk/P), - data=dfm)) -summary(av1) - -## Main factor: SP -tk4 <- TukeyC(av1, - which='SP', - error='Within', - sig.level=0.1) -summary(tk4) - -## Main factor: P -tk5 <- TukeyC(av1, - which='P', - error='blk:P') - -summary(tk5) - -## Nested: p1/SP -## Testing SP inside of level one of P -tkn3 <- TukeyC.nest(av1, - which='P:SP', - error='Within', - fl1=1) -summary(tkn3) - -## Nested: p2/SP -tkn4 <- TukeyC.nest(av1, - which='P:SP', - error='Within', - fl1=2) -summary(tkn4) - -## Nested: p3/SP -tkn5 <- TukeyC.nest(av1, - which='P:SP', - error='Within', - fl1=3) -summary(tkn5) diff --git a/demo/SPET.R b/demo/SPET.R new file mode 100644 index 0000000..bff5a46 --- /dev/null +++ b/demo/SPET.R @@ -0,0 +1,124 @@ +## +## Example: Split-plot Experiment in time (SPET) +## + +## The parameters can be: vectors, design matrix and the response variable, +## data.frame, aov or lm. + +## Note: The factors are in uppercase and its levels in lowercase! + +library(TukeyC) +data(SPET) + +## From: formula +## Main factor: year +tk1 <- with(SPET, + TukeyC(y ~ blk + tra*year + Error(blk/tra), + dfm, + which='year')) +summary(tk1) + +## Nested: crotjuncea/year +tk2 <- with(SPET, + TukeyC(y ~ blk + tra*year + Error(blk/tra), + dfm, + which='tra:year', + fl1=2)) +summary(tk2) + +## Nested: year_1/tra +## It is necessary to inform how to combinate the errors +tk3 <- with(SPET, + TukeyC(y ~ blk + tra*year + Error(blk/tra), + dfm, + which='year:tra', + error='Within/blk:tra', + fl1=1)) +summary(tk3) + + +## From: lm +lm1 <- with(SPET, + lm(y ~ blk*tra + tra*year, + data=dfm)) + +## Nested: tra1/year +tk4 <- TukeyC(lm1, + which='tra:year', + fl1=1) + +summary(tk4) + +## Nested: year1/tra +## It is necessary to inform how to combinate the errors +tk5 <- TukeyC(lm1, + which='year:tra', + error='Within/blk:tra', + fl1=1) +summary(tk5, + complete=FALSE) + +## Nested: year2/tra +## It is necessary to inform how to combinate the errors +tk6 <- TukeyC(lm1, + which='year:tra', + error='Within/blk:tra', + fl1=2) +summary(tk6, + complete=FALSE) + +## From: aov +av1 <- with(SPET, + aov(y ~ blk + tra*year + Error(blk/tra), + data=dfm)) +summary(av1) + +## Main factor: year +tk7 <- TukeyC(av1, + which='year') +summary(tk7) + +## Main factor: tra +## It is necessary to inform the appropriate error for the test +tk8 <- TukeyC(av1, + which='tra', + error='blk:tra') +summary(tk8, + complete=FALSE) + +## Nested: crotjuncea/year +tk9 <- TukeyC(av1, + which='tra:year', + fl1=2) +summary(tk9) + +## Nested: guandu/year +tk10 <- TukeyC(av1, + which='tra:year', + fl1=4) +summary(tk10) + +## Nested: year_1/tra - it is necessary to inform how to combinate the errors +tk11 <- TukeyC(av1, + which='year:tra', + error='Within/blk:tra', + fl1=1) +summary(tk11, + complete=FALSE) + +op <- par(mar=c(6, 3, 3, 2)) +plot(tk10, + id.las=2, + xlab='') + +## Nested: year_2/tra - it is necessary to inform how to combinate the errors +tk12 <- TukeyC(av1, + which='year:tra', + error='Within/blk:tra', + fl1=2) +summary(tk12) +op <- par(mar=c(7, 3, 3, 2)) +plot(tk12, + id.las=2, + xlab='') +par(op) diff --git a/demo/SPET.r b/demo/SPET.r deleted file mode 100644 index d73c2d3..0000000 --- a/demo/SPET.r +++ /dev/null @@ -1,141 +0,0 @@ -## -## Example: Split-plot Experiment in time (SPET) -## - -## The parameters can be: design matrix and the response variable, -## data.frame or aov - -## Note: The factors are in uppercase and its levels in lowercase! - -library(TukeyC) -data(SPET) - -## From: design matrix (dm) and response variable (y) -## Main factor: tra -tk1 <- with(SPET, - TukeyC(x=dm, - y=y, - model='y ~ blk + tra*year + Error(blk/tra)', - which='tra', - error='blk:tra')) -summary(tk1) - -plot(tk1) - -## Main factor: year -tk2 <- with(SPET, - TukeyC(x=dm, - y=y, - model='y ~ blk + tra*year + Error(blk/tra)', - which='year', - error='Within')) -summary(tk2) -plot(tk2, - title='Main effect: year') - -## Nested: crotjuncea/year -tkn1 <- with(SPET, - TukeyC.nest(x=dm, - y=y, - model='y ~ blk + tra*year + Error(blk/tra)', - which='tra:year', - error='Within', - fl1=2)) -summary(tkn1) -plot(tkn1, - title='Effect: crotjuncea/year') - - -## From: data.frame -## Main factor: year -tk3 <- with(SPET, - TukeyC(dfm, - model='y ~ blk + tra*year + Error(blk/tra)', - which='year', - error='Within')) -summary(tk3) -plot(tk3, - title='Main effect: year') - -## Nested: crotjuncea/year -tkn2 <- with(SPET, - TukeyC.nest(dfm, - model='y ~ blk + tra*year + Error(blk/tra)', - which='tra:year', - error ='Within', - fl1=2)) -summary(tkn2) -plot(tkn2, - title='Effect: crotjuncea/year') - -## Nested: year_1/tra -tkn3 <- with(SPET, - TukeyC.nest(dfm, - model='y ~ blk + tra*year + Error(blk/tra)', - which='year:tra', - error ='Within', - fl1=1)) -summary(tkn3) -plot(tkn3, - title='Effect: year_1/tra') - - - -## From: aovlist -av1 <- with(SPET, - aov(y ~ blk + tra*year + Error(blk/tra), - data=dfm)) -summary(av1) - -## Main factor: year -tk4 <- TukeyC(av1, - which='year', - error='Within') -summary(tk4) -plot(tk4, title='Main effect: year') - -## Main factor: tra -tk5 <- TukeyC(av1, - which='tra', - error='blk:tra') -summary(tk5) -plot(tk5, - title='Main effect: tra') - -## Nested: crotjuncea/year -tkn4 <- TukeyC.nest(av1, - which='tra:year', - error='Within', - fl1=2) -summary(tkn4) -plot(tkn4, - title='Effect: crotjuncea/year') - -## Nested: guandu/year -tkn5 <- TukeyC.nest(av1, - which='tra:year', - error='Within', - fl1=4) -summary(tkn5) -plot(tkn5, - title='Effect: guandu/year') - - -## Nested: year_1/tra -tkn6 <- TukeyC.nest(av1, - which='year:tra', - error='Within', - fl1=1) -summary(tkn6) -plot(tkn6, - title='Effect: year_1/tra') - -## Nested: year_2/tra -tkn7 <- TukeyC.nest(av1, - which='year:tra', - error='Within', - fl1=2) -summary(tkn7) -plot(tkn7, - title='Effect: year_2/tra') - diff --git a/demo/SSPE.R b/demo/SSPE.R new file mode 100644 index 0000000..628a39d --- /dev/null +++ b/demo/SSPE.R @@ -0,0 +1,132 @@ +## +## Example: Split-split-plot Experiment (SSPE) +## + +## The parameters can be: vectors, design matrix and the response variable, +## data.frame, aov or lm. + +## Note: Upper case for factors and lowercase for levels + +library(TukeyC) +data(SSPE) + +## From: formula +## Main factor: P +## It is necessary to inform the appropriate error for the test +tk1 <- with(SSPE, + TukeyC(y ~ blk + P*SP*SSP + Error(blk/P/SP), + dfm, + which='P', + error='blk:P')) +summary(tk1) + +## Nested: p2/SP +## It is necessary to inform the appropriate error for the test +tk2 <- with(SSPE, + TukeyC(y ~ blk + P*SP*SSP + Error(blk/P/SP), + dfm, + which='P:SP', + error='blk:P:SP', + fl1=2)) +summary(tk2) + +## Nested: p2/SSP +tk3 <- with(SSPE, + TukeyC(y ~ blk + P*SP*SSP + Error(blk/P/SP), + dfm, + which='P:SSP', + fl1=2)) +summary(tk3) +plot(tk3) + + +## From: lm +lm1 <- with(SSPE, + lm(y ~ blk*P + blk*P*SP + P*SP*SSP, + data=dfm)) +summary(lm1) + +## Main factor: P +## It is necessary to inform the appropriate error for the test +tk4 <- TukeyC(lm1, + which='P', + error='blk:P') +summary(tk4) + +## Main factor: SP +tk5 <- TukeyC(lm1, + which='SP', + error='blk:P:SP') +summary(tk5) + +## Main factor: SSP +tk6 <- TukeyC(lm1, + which='SSP') +summary(tk6) + +## Nested: p1/SP +## It is necessary to inform the appropriate error for the test +tk7 <- TukeyC(lm1, + which='P:SP', + error='blk:P:SP', + fl1=1) +summary(tk7) + + +## From: aov +av1 <- with(SSPE, + aov(y ~ blk + P*SP*SSP + Error(blk/P/SP), + data=dfm)) +summary(av1) + +## Main factor: P +## It is necessary to inform the appropriate error for the test +tk8 <- TukeyC(av1, + which='P', + error='blk:P') +summary(tk8) + +## Main factor: SSP +tk9 <- TukeyC(av1, + which='SSP') +summary(tk9) + +## Nested: p1/SP +## It is necessary to inform the appropriate error for the test +tk10 <- TukeyC(av1, + which='P:SP', + error='blk:P:SP', + fl1=1) +summary(tk10) + +## Nested: p2/SP +tk11 <- TukeyC(av1, + which='P:SP', + error='blk:P:SP', + fl1=2) +summary(tk11) + +## Nested: Pi/SPi/SSP (at various levels of P and SP) +tk12 <- TukeyC(av1, + which='P:SP:SSP', + fl1=1, + fl2=1) +summary(tk12) + +tk13 <- TukeyC(av1, + which='P:SP:SSP', + fl1=2, + fl2=1) +summary(tk13) + +tk14 <- TukeyC(av1, + which='P:SP:SSP', + fl1=3, + fl2=3) +summary(tk14) + +tk15 <- TukeyC(av1, + which='P:SP:SSP', + fl1=2, + fl2=3) +summary(tk15) diff --git a/demo/SSPE.r b/demo/SSPE.r deleted file mode 100644 index 1ecb38d..0000000 --- a/demo/SSPE.r +++ /dev/null @@ -1,141 +0,0 @@ -## -## Example: Split-split-plot Experiment (SSPE) -## - -## Note: Upper case for factors and lowercase for levels - -library(TukeyC) -data(SSPE) - -## From: design matrix (dm) and response variable (y) -## Main factor: P -tk1 <- with(SSPE, - TukeyC(dm, - y, - model='y ~ blk + P*SP*SSP + Error(blk/P/SP)', - which='P', - error='blk:P')) -summary(tk1) -plot(tk1) - -# Main factor: SP -tk2 <- with(SSPE, - TukeyC(dm, - y, - model='y ~ blk + P*SP*SSP + Error(blk/P/SP)', - which='SP', - error='blk:P:SP')) -summary(tk2) -plot(tk2) - -# Main factor: SSP -tk3 <- with(SSPE, - TukeyC(dm, - y, - model='y ~ blk + P*SP*SSP + Error(blk/P/SP)', - which='SSP', - error='Within')) -summary(tk3) -plot(tk3) - -## Nested: p1/SP -tkn1 <- with(SSPE, - TukeyC.nest(dm, - y, - model='y ~ blk + P*SP*SSP + Error(blk/P/SP)', - which='P:SP', - error='blk:P:SP', - fl1=1)) -summary(tkn1) - - -## From: data.frame -## Main factor: P -tk4 <- with(SSPE, - TukeyC(dfm, - model='y ~ blk + P*SP*SSP + Error(blk/P/SP)', - which='P', - error='blk:P')) -summary(tk4) - -## Nested: p2/SP -tkn2 <- with(SSPE, - TukeyC.nest(dfm, - model='y ~ blk + P*SP*SSP + Error(blk/P/SP)', - which='P:SP', - error='blk:P:SP', - fl1=2)) -summary(tkn2) - -## Nested: p2/SSP -tkn3 <- with(SSPE, - TukeyC.nest(dfm, - model='y ~ blk + P*SP*SSP + Error(blk/P/SP)', - which='P:SSP', - error='Within', - fl1=2)) -summary(tkn3) -plot(tkn3) - - -## From: aovlist -av <- with(SSPE, - aov(y ~ blk + P*SP*SSP + Error(blk/P/SP), - data=dfm)) -summary(av) - -## Main factor: P -tk5 <- TukeyC(av, - which='P', - error='blk:P') -summary(tk5) - -## Main factor: SSP -tk6 <- TukeyC(av, - which='SSP', - error='Within') -summary(tk6) - -## Nested: p1/SP -tkn4 <- TukeyC.nest(av, - which='P:SP', - error='blk:P:SP', - fl1=1) -summary(tkn4) - -## Nested: p2/SP -tkn5 <- TukeyC.nest(av, - which='P:SP', - error='blk:P:SP', - fl1=2) -summary(tkn5) - -## Nested: Pi/SPi/SSP (at various levels of P and SP) -tkn6 <- TukeyC.nest(av, - which='P:SP:SSP', - error='Within', - fl1=1, - fl2=1) -summary(tkn6) -plot(tkn6) - -tkn7 <- TukeyC.nest(av, - which='P:SP:SSP', - error='Within', - fl1=2, - fl2=1) -summary(tkn7) - -tkn8 <- TukeyC.nest(av, - which='P:SP:SSP', - error='Within', - fl1=3, - fl2=3) -summary(tkn8) - -tkn9 <- TukeyC.nest(av, - which='P:SP:SSP', - error='Within', - fl1=2, - fl2=3) -summary(tkn9) diff --git a/man/TukeyC-internal.Rd b/man/TukeyC-internal.Rd new file mode 100644 index 0000000..9c0458e --- /dev/null +++ b/man/TukeyC-internal.Rd @@ -0,0 +1,9 @@ +\name{TukeyC-internal} +\title{Internal TukeyC functions} +\alias{m.infos.aovlist} +\alias{m.infos.lm} +\alias{m.infos.nest.aovlist} +\alias{m.infos.nest.lm} +\description{Internal TukeyC functions.} +\details{These are not to be called by the user and are undocumented.} +\keyword{internal} diff --git a/man/TukeyC-package.Rd b/man/TukeyC-package.Rd index cd0b904..037f009 100644 --- a/man/TukeyC-package.Rd +++ b/man/TukeyC-package.Rd @@ -1,5 +1,6 @@ -\name{TukeyC-package} +\name{Tukey test} \alias{TukeyC-package} + \docType{package} \title{ @@ -64,313 +65,3 @@ \keyword{package} -\examples{ - ## - ## Examples: Completely Randomized Design (CRD) - ## More details: demo(package='TukeyC') - ## - - ## The parameters can be: vectors, design matrix and the response variable, - ## data.frame or aov - data(CRD2) - - ## From: design matrix (dm) and response variable (y) - balanced - tk1 <- with(CRD2, - TukeyC(x=dm, - y=y, - model='y ~ x', - which='x')) - summary(tk1) - plot(tk1, - id.las=2, - rl=FALSE) - - ## From: design matrix (dm) and response variable (y) - unbalanced - tk1u <- with(CRD2, - TukeyC(x=dm[-1,], - y=y[-1], - model='y ~ x', - which='x', - dispersion='s')) - summary(tk1u) - plot(tk1u) - - ## From: data.frame (dfm) - balanced - tk2 <- with(CRD2, - TukeyC(x=dfm, - model='y ~ x', - which='x', - dispersion='se')) - summary(tk2) - plot(tk2) - - ## From: data.frame (dfm) - balanced - tk2u <- with(CRD2, - TukeyC(x=dfm[-1,], - model='y ~ x', - which='x')) - summary(tk2u) - - ## From: aov - balanced - av <- with(CRD2, - aov(y ~ x, - data=dfm)) - summary(av) - - tk3 <- TukeyC(x=av, - which='x') - summary(tk3) - - ## From: aov - unbalanced - avu <- with(CRD2, - aov(y ~ x, - data=dfm[-1,])) - summary(avu) - - tk3u <- TukeyC(x=avu, - which='x') - summary(tk3u) - - ## - ## Example: Randomized Complete Block Design (RCBD) - ## More details: demo(package='TukeyC') - ## - - ## The parameters can be: design matrix and the response variable, - ## data.frame or aov - - data(RCBD) - - ## Design matrix (dm) and response variable (y) - tk1 <- with(RCBD, - TukeyC(x=dm, - y=y, - model='y ~ blk + tra', - which='tra')) - summary(tk1) - plot(tk1) - - ## From: data.frame (dfm), which='tra' - tk2 <- with(RCBD, - TukeyC(x=dfm, - model='y ~ blk + tra', - which='tra')) - summary(tk2) - - ## - ## Example: Latin Squares Design (LSD) - ## More details: demo(package='TukeyC') - ## - - ## The parameters can be: design matrix and the response variable, - ## data.frame or aov - - data(LSD) - - ## From: design matrix (dm) and response variable (y) - tk1 <- with(LSD, - TukeyC(x=dm, - y=y, - model='y ~ rows + cols + tra', - which='tra')) - summary(tk1) - plot(tk1) - - ## From: data.frame - tk2 <- with(LSD, - TukeyC(x=dfm, - model='y ~ rows + cols + tra', - which='tra')) - summary(tk2) - - ## From: aov - av <- with(LSD, - aov(y ~ rows + cols + tra, - data=dfm)) - summary(av) - - tk3 <- TukeyC(av, - which='tra') - summary(tk3) - - ## - ## Example: Factorial Experiment (FE) - ## More details: demo(package='TukeyC') - ## - - ## The parameters can be: design matrix and the response variable, - ## data.frame or aov - - data(FE) - - ## From: design matrix (dm) and response variable (y) - ## Main factor: N - tk1 <- with(FE, - TukeyC(x=dm, - y=y, - model='y ~ blk + N*P*K', - which='N')) - summary(tk1) - plot(tk1) - - ## Nested: p1/N - ## Testing N inside of level one of P - ntk1 <- with(FE, - TukeyC.nest(x=dm, - y=y, - model='y ~ blk + N*P*K', - which='P:N', - fl1=1)) - summary(ntk1) - - ## Nested: k1/P - ntk2 <- with(FE, - TukeyC.nest(x=dm, - y=y, - model='y ~ blk + N*P*K', - which='K:P', - fl1=1)) - summary(ntk2) - - ## Nested: k1/p1/N - ## Testing N inside of level one of K and level one of P - ntk3 <- with(FE, - TukeyC.nest(x=dm, - y=y, - model='y ~ blk + N*P*K', - which='K:P:N', - fl1=1, - fl2=1)) - summary(ntk3) - - ## Nested: k2/n2/P - ntk4 <- with(FE, - TukeyC.nest(x=dm, - y=y, - model='y ~ blk + N*P*K', - which='K:N:P', - fl1=2, - fl2=2)) - summary(ntk4) - - ## Nested: p1/n1/K - ntk5 <- with(FE, - TukeyC.nest(x=dm, - y=y, - model='y ~ blk + N*P*K', - which='P:N:K', - fl1=1, - fl2=1)) - summary(ntk5) - - ## - ## Example: Split-plot Experiment (SPE) - ## More details: demo(package='TukeyC') - ## - - data(SPE) - - ## The parameters can be: design matrix and the response variable, - ## data.frame or aov - - ## From: design matrix (dm) and response variable (y) - ## Main factor: P - tk1 <- with(SPE, - TukeyC(x=dm, - y=y, - model='y ~ blk + P*SP + Error(blk/P)', - which='P', - error='blk:P')) - summary(tk1) - - ## Main factor: SP - tk2 <- with(SPE, - TukeyC(x=dm, - y=y, - model='y ~ blk + P*SP + Error(blk/P)', - which='SP', - error='Within')) - summary(tk2) - plot(tk2) - - ## Nested: p=1/sp - tkn1 <- with(SPE, - TukeyC.nest(x=dm, - y=y, - model='y ~ blk + P*SP + Error(blk/P)', - which='P:SP', - error='Within', - fl1=1 )) - summary(tkn1) - - ## - ## Example: Split-split-plot Experiment (SSPE) - ## More details: demo(package='TukeyC') - ## - - data(SSPE) - - ## From: design matrix (dm) and response variable (y) - ## Main factor: P - tk1 <- with(SSPE, - TukeyC(dm, - y, - model='y ~ blk + P*SP*SSP + Error(blk/P/SP)', - which='P', - error='blk:P')) - summary(tk1) - - # Main factor: SP - tk2 <- with(SSPE, - TukeyC(dm, - y, - model='y ~ blk + P*SP*SSP + Error(blk/P/SP)', - which='SP', - error='blk:P:SP')) - summary(tk2) - - # Main factor: SSP - tk3 <- with(SSPE, - TukeyC(dm, - y, - model='y ~ blk + P*SP*SSP + Error(blk/P/SP)', - which='SSP', - error='Within')) - summary(tk3) - plot(tk3) - - ## Nested: p1/SP - tkn1 <- with(SSPE, - TukeyC.nest(dm, - y, - model='y ~ blk + P*SP*SSP + Error(blk/P/SP)', - which='P:SP', - error='blk:P:SP', - fl1=1)) - summary(tkn1) - - ## From: aovlist - av <- with(SSPE, - aov(y ~ blk + P*SP*SSP + Error(blk/P/SP), - data=dfm)) - summary(av) - - ## Nested: p1/sp1/SSP - ## Testing SSP inside of level one of P and level one of SP - tkn6 <- TukeyC.nest(av, - which='P:SP:SSP', - error='Within', - fl1=1, - fl2=1) - summary(tkn6) - plot(tkn6) - - ## Nested: p2/sp1/SSP - tkn7 <- TukeyC.nest(av, - which='P:SP:SSP', - error='Within', - fl1=2, - fl2=1) - summary(tkn7) - plot(tkn7) -} diff --git a/man/TukeyC.Rd b/man/TukeyC.Rd index d77f42b..cc96c77 100644 --- a/man/TukeyC.Rd +++ b/man/TukeyC.Rd @@ -1,7 +1,8 @@ \name{TukeyC} \alias{TukeyC} \alias{TukeyC.default} -\alias{TukeyC.aov} +\alias{TukeyC.lm} +\alias{TukeyC.formula} \alias{TukeyC.aovlist} \title{ @@ -9,71 +10,77 @@ } \description{ - These are methods for objects of class \code{vector}, \code{matrix} or - \code{data.frame} joined as default, \code{aov} and \code{aovlist} for - single experiments. + These are methods for objects of class \code{formula}, \code{lm}, \code{aov} and \code{aovlist} for single, factorial, split-plot and split-split-plot experiments. } \usage{ - \method{TukeyC}{default}(x, - y=NULL, - model, - which, - error, - sig.level=.05, - round=2, - dispersion=c('mm', 's', 'se'), \dots) - - \method{TukeyC}{aov}(x, - which=NULL, - sig.level=.05, - round=2, - dispersion=c('mm', 's', 'se'), \dots) - - \method{TukeyC}{aovlist}(x, - which, - error, - sig.level=.05, - round=2, - dispersion=c('mm', 's', 'se'), \dots) +TukeyC(x,...) + +\method{TukeyC}{formula}(formula, + data = NULL, + which = NULL, + fl1 = NULL, + fl2 = NULL, + error = NULL, + sig.level = .05, + round = 2, + adjusted.pvalue = 'none', + \dots) + +\method{TukeyC}{lm}(x, + which = NULL, + fl1 = NULL, + fl2 = NULL, + error = NULL, + sig.level = .05, + round = 2, + adjusted.pvalue = 'none', + \dots) + +\method{TukeyC}{aovlist}(x, + which = NULL, + fl1 = NULL, + fl2 = NULL, + error = NULL, + sig.level = .05, + round = 2, + adjusted.pvalue = 'none', + \dots) + } \arguments{ - \item{x}{A design matrix, \code{data.frame} or an \code{aov} object.} - \item{y}{A vector of response variable. It is necessary to inform this - parameter only if \samp{x} represent the design matrix.} + \item{x,formula}{A \code{formula}, \code{lm}, \code{aov} and \code{aovlist} object. Objects of the \code{formula} class follow ``response variable ~ predicted variable.} + \item{data}{A object of the \code{data.frame} class. Use only objects of \code{formula} class.} \item{which}{The name of the treatment to be used in the comparison. The name must be inside quoting marks.} - \item{model}{If \samp{x} is a \code{data.frame} object, the model to be used in the - aov must be specified.} - \item{error}{The error to be considered.} + \item{fl1}{A vector of length 1 giving the level of the first factor in nesting order tested.} + \item{fl2}{A vector of length 1 giving the level of the second factor in nesting order tested.} + \item{error}{The error to be considered. If from experiment at split plot or split-split plot pay attention! See details!} \item{sig.level}{Level of Significance used in the TukeyC algorithm to create the groups of means. The default value is 0.05.} \item{round}{Integer indicating the number of decimal places.} - \item{dispersion}{The dispersion to be considered to the means. \cr - The possible vaues are: \samp{'mm'} = \samp{minimum and maximum}, - \samp{'s'} = \samp{standart deviation}, - \samp{'se'} = \samp{standart deviation of the mean}.} + \item{adjusted.pvalue}{Method for adjusting p values (see \code{p.adjust} to more details). The possible values are: \code{"holm"}, \code{"hochberg"}, \code{"hommel"}, \code{"bonferroni"}, \code{"BH"}, \code{"BY"}, \code{"fdr"} and \code{"none"}. The default is \code{"none"}. } \item{\dots}{Potential further arguments (required by generic).} } \details{ The function \code{TukeyC} returns an object of class \code{TukeyC} - respectivally containing the groups of means plus other + containing the groups of means plus other necessary variables for summary and plot. - + The generic functions \code{summary} and \code{plot} are used to obtain and print a summary and a plot of the results. + + The error arguments may be used whenever the user want a specific error other than the experimental error. At the split plot and split-split plot experiment, combination of error may be specified with "/" in the sequence of the \code{which} argument. For example, a object of \code{aovlist} class, a possible combination would be \code{error = 'Within/blk:plot'} at case block split plot experiment with \code{which = 'subplot:plot'} argument. } \value{ The function \code{TukeyC} returns a list of the class \code{TukeyC} with the slots: - \item{av}{A \code{list} storing the result of \code{aov}.} - \item{groups}{A vector of length equal the number of factor levels marking the groups generated.} - \item{nms}{A vector of the labels of the factor levels.} - \item{ord}{A vector which keeps the position of the means of the factor levels in decreasing order.} - \item{m.inf}{A matrix which keeps the means and dispersion of the factor levels in decreasing order.} - \item{sig.level}{A vector of length 1 giving the level of significance of the test.} + \item{Result}{A \code{data.frame} storing the result of Tukey test.} + \item{Sig.level}{A scalar giving the level of significance of the test.} + \item{Diff_Prob}{A \code{matrix} at the lower diagonal with p-values and upper diagonal with means differences.} + \item{MSD}{A \code{matrix} with minimum significance differences by Tukey methodology. If balanced data, then all values are equal.} } \author{ @@ -86,140 +93,280 @@ Miller, R.G. (1981) \emph{Simultaneous Statistical Inference}. Springer. Ramalho M.A.P, Ferreira D.F & Oliveira A.C. (2000) \emph{Experimentação em Genética - e Melhoramento de Plantas}. Editora UFLA. + e Melhoramento de Plantas}. Editora UFLA. Steel, R.G., Torrie, J.H & Dickey D.A. (1997) \emph{Principles and procedures of statistics: - a biometrical approach}. Third Edition. - + a biometrical approach}. Third Edition. + Yandell, B.S. (1997) \emph{Practical Data Analysis for Designed Experiments}. Chapman & Hall. } \examples{ - ## - ## Examples: Completely Randomized Design (CRD) - ## More details: demo(package='TukeyC') - ## - - ## The parameters can be: vectors, design matrix and the response variable, - ## data.frame or aov - data(CRD2) - - ## From: design matrix (dm) and response variable (y) - tk1 <- with(CRD2, - TukeyC(x=dm, - y=y, - model='y ~ x', - which='x', - id.trim=5)) - summary(tk1) - - ## From: data.frame (dfm) - tk2 <- with(CRD2, - TukeyC(x=dfm, - model='y ~ x', - which='x', - id.trim=5)) - summary(tk2) - - ## From: aov - av <- with(CRD2, - aov(y ~ x, - data=dfm)) - summary(av) - - tk3 <- with(CRD2, - TukeyC(x=av, - which='x', - id.trim=5)) - summary(tk3) - - ## - ## Example: Randomized Complete Block Design (RCBD) - ## More details: demo(package='TukeyC') - ## - - ## The parameters can be: design matrix and the response variable, - ## data.frame or aov - - data(RCBD) - - ## Design matrix (dm) and response variable (y) - tk1 <- with(RCBD, - TukeyC(x=dm, - y=y, - model='y ~ blk + tra', - which='tra')) - summary(tk1) - - ## From: data.frame (dfm), which='tra' - tk2 <- with(RCBD, - TukeyC(x=dfm, - model='y ~ blk + tra', - which='tra')) - summary(tk2) - - ## - ## Example: Latin Squares Design (LSD) - ## More details: demo(package='TukeyC') - ## - - ## The parameters can be: design matrix and the response variable, - ## data.frame or aov - - data(LSD) - - ## From: design matrix (dm) and response variable (y) - tk1 <- with(LSD, - TukeyC(x=dm, - y=y, - model='y ~ rows + cols + tra', - which='tra')) - summary(tk1) - - ## From: data.frame - tk2 <- with(LSD, - TukeyC(x=dfm, - model='y ~ rows + cols + tra', - which='tra')) - summary(tk2) - - ## From: aov - av <- with(LSD, - aov(y ~ rows + cols + tra, - data=dfm)) - summary(av) - - tk3 <- TukeyC(av, - which='tra') - summary(tk3) - - ## - ## Example: Factorial Experiment (FE) - ## More details: demo(package='TukeyC') - ## - - ## The parameters can be: design matrix and the response variable, - ## data.frame or aov - - data(FE) - ## From: design matrix (dm) and response variable (y) - ## Main factor: N - tk1 <- with(FE, - TukeyC(x=dm, - y=y, - model='y ~ blk + N*P*K', - which='N')) - summary(tk1) - - ## Nested: p1/N - ntk1 <- with(FE, - TukeyC.nest(x=dm, - y=y, - model='y ~ blk + N*P*K', - which='P:N', - fl1=1)) - summary(ntk1) +## +## Examples: Completely Randomized Design (CRD) +## More details: demo(package='TukeyC') +## + +## The parameters can be: formula, aov, lm and aovlist. + +data(RCBD) + +## From: formula +tk1 <- with(RCBD, + TukeyC(y ~ blk + tra, + data=dfm, + which='tra')) +summary(tk1) + +## +## Example: Latin Squares Design (LSD) +## More details: demo(package='TukeyC') +## + +data(LSD) + +## From: formula +tk2 <- with(LSD, + TukeyC(y ~ rows + cols + tra, + data=dfm, + which='tra')) +summary(tk2) + +## From: aov +av1 <- with(LSD, + aov(y ~ rows + cols + tra, + data=dfm)) + +tk3 <- TukeyC(av1, + which='tra') +summary(tk3) + +## From: lm +lm1 <- with(LSD, + lm(y ~ rows + cols + tra, + data=dfm)) + +tk4 <- TukeyC(lm1, + which='tra') +summary(tk4) + +## +## Example: Factorial Experiment (FE) +## More details: demo(package='TukeyC') +## + +data(FE) +## From: formula +## Main factor: N +tk5 <- with(FE, + TukeyC(y ~ blk + N*P*K, + data=dfm, + which='N')) +summary(tk5) + +## Nested: p1/N +# From: formula +ntk1 <- with(FE, + TukeyC(y ~ blk + N*P*K, + data=dfm, + which='P:N', + fl1=1)) +summary(ntk1) + +## Nested: p2/N +# From: lm +lm2 <- with(FE, + lm(y ~ blk + N*P*K, + dfm)) + +ntk2 <- with(FE, + TukeyC(lm2, + which='P:N', + fl1=2)) +summary(ntk2) + +## Nested: n1/P +# From: aov +av2 <- with(FE, + aov(y ~ blk + N*P*K, + dfm)) + +ntk3 <- with(FE, + TukeyC(av2, + which='N:P', + fl1=1)) +summary(ntk3) + +## +## Example: Split-plot Experiment (SPET) +## More details: demo(package='TukeyC') +## +data(SPET) + +lm2 <- with(SPET, + lm(y ~ blk*tra + tra*year, + dfm)) + +# crotgrantiana/year +sptk1 <- TukeyC(lm2, + which='tra:year', + fl1=1) +summary(sptk1) + +# year1/tra +# It is necessary to set year error with trat error in the order of the "which" argument. +# It is necessary to inform how to combinate the errors +tk6 <- with(FE, + TukeyC(lm2, + which='year:tra', + error='Residuals/blk:tra', + fl1=1)) +summary(tk6) + + +## Example: Split-split-plot Experiment (SSPE) +## More details: demo(package='TukeyC') +## + +data(SSPE) +## From: formula +## Main factor: P +## It is necessary to inform the appropriate error for the test +tk7 <- with(SSPE, + TukeyC(y ~ blk + P*SP*SSP + Error(blk/P/SP), + data=dfm, + which='P', + error='blk:P')) +summary(tk7) + +## Main factor: SP +## It is necessary to inform the appropriate error for the test +tk8 <- with(SSPE, + TukeyC(y ~ blk + P*SP*SSP + Error(blk/P/SP), + data=dfm, + which='SP', + error='blk:P:SP')) +summary(tk8) + +## Main factor: SSP +tk9 <- with(SSPE, + TukeyC(y ~ blk + P*SP*SSP + Error(blk/P/SP), + data=dfm, + which='SSP')) +summary(tk9) + +## From: aov +## Main factor: SSP +av1 <- with(SSPE, + aov(y ~ blk + P*SP*SSP + Error(blk/P/SP), + data=dfm)) + +tk10 <- TukeyC(av1, + which='SSP') +summary(tk10) + +## Nested: p1/SP +## It is necessary to inform the appropriate error for the test +tk11 <- TukeyC(av1, + which='P:SP', + error='blk:P:SP', + fl1=1) +summary(tk11) + +## Nested: p1/SSP +tk12 <- TukeyC(av1, + which='P:SSP', + fl1=1) +summary(tk12) + +## Nested: p1/sp1/SSP +## Testing SSP inside of level one of P and level one of SP +tk13 <- TukeyC(av1, + which='P:SP:SSP', + fl1=1, + fl2=1) +summary(tk13) + +## Nested: p2/sp1/SSP +tk14 <- TukeyC(av1, + which='P:SP:SSP', + fl1=2, + fl2=1) +summary(tk14) + +## Nested: sp1/P +## It is necessary to inform the appropriate error for the test +tk15 <- TukeyC(av1, + which='SP:P', + error='blk:P:SP/blk:P', + fl1=1) + +summary(tk15) + +## Nested: ssp1/SP +tk16 <- TukeyC(av1, + which='SSP:SP', + error='Within/blk:P:SP', + fl1=1) +summary(tk16) + +## Nested: ssp1/sp1/P +## It is necessary to inform the appropriate error for the test +tk17 <- TukeyC(av1, + which='SSP:SP:P', + error='Within/blk:P:SP/blk:P', + fl1=1, + fl2=1) +summary(tk17) + +## UNBALANCED DATA +## The average are adjusted by "Least-Square-Means" methodology. +## From: formula +data(CRD2) + +uCRD2 <- CRD2$dfm +uCRD2[c(3, 5, 10, 44, 45), 3] <- NA + +utk1 <- TukeyC(y ~ x, + data=uCRD2, + which='x') +summary(utk1) + +## From: lm +ulm1 <- lm(y ~ x, + data=uCRD2) + +utk2 <- TukeyC(ulm1, + which='x') +summary(utk2) + + +## Factorial Experiments +## Nested: p1/N +# From: lm + +uFE <- FE$dfm +uFE[c(3, 6, 7, 20, 31, 32), 5] <- NA + +ulm2 <- lm(y ~ blk + N*P*K, + uFE) + +## Nested: p1/N +utk3 <- TukeyC(ulm2, + data=uFE, + which='P:N', + fl1=1) +summary(utk3) + +## Nested: p2/n2/K +utk4 <- TukeyC(ulm2, + data=uFE, + which='P:N:K', + fl1=2, + fl2=2) +summary(utk4) + } \keyword{package} diff --git a/man/TukeyC.nest.Rd b/man/TukeyC.nest.Rd deleted file mode 100644 index ff6c10f..0000000 --- a/man/TukeyC.nest.Rd +++ /dev/null @@ -1,186 +0,0 @@ -\name{TukeyC.nest} -\alias{TukeyC.nest} -\alias{TukeyC.nest.default} -\alias{TukeyC.nest.aov} -\alias{TukeyC.nest.aovlist} - -\title{ - The TukeyC test for Factorial, Split-plot and Split-SPlit plot Experiments -} - -\description{ - These are methods for objects of class \code{vector}, \code{matrix} or - \code{data.frame} joined as default, \code{aov} and \code{aovlist} for - factorial, split-plot and split-split-plot experiments. -} - -\usage{ - \method{TukeyC.nest}{default}(x, - y=NULL, - model, - which, - error, - fl1, - fl2=0, - sig.level=.05, - round=2, - dispersion=c('mm', 's', 'se'), \dots) - - \method{TukeyC.nest}{aov}(x, - which, - fl1, - fl2=0, - sig.level=.05, - round=2, - dispersion=c('mm', 's', 'se'), \dots) - - \method{TukeyC.nest}{aovlist}(x, - which, - error, - fl1, - fl2=0, - sig.level=.05, - round=2, - dispersion=c('mm', 's', 'se'), \dots) -} - -\arguments{ - \item{x}{A design matrix, \code{data.frame} or an \code{aov} object.} - \item{y}{A vector of response variable. It is necessary to inform this - parameter only if \code{x} represent the design matrix.} - \item{which}{The name of the treatment to be used in the comparison. - The name must be inside quoting marks.} - \item{model}{If x is a \code{data.frame} object, the model to be used in the - aov must be specified.} - \item{fl1}{A vector of length 1 giving the level of the first factor in nesting order tested.} - \item{fl2}{A vector of length 1 giving the level of the second factor in nesting order tested.} - \item{error}{The error to be considered, only in case of split-plots experiments.} - \item{sig.level}{Level of Significance used in the TukeyC algorithm to create - the groups of means. The default value is 0.05.} - \item{round}{Integer indicating the number of decimal places.} - \item{dispersion}{The dispersion to be considered to the means. \cr - The possible vaues are: \samp{'mm'} = \samp{minimum and maximum}, - \samp{'s'} = \samp{standart deviation}, - \samp{'se'} = \samp{standart deviation of the mean}.} - \item{\dots}{Potential further arguments (required by generic).} -} - -\details{ - The function \code{TukeyC.nest} returns an object of class - \code{TukeyC.nest} containing the groups of means plus other - necessary variables for summary and plot. - - The generic functions \code{summary} and \code{plot} are used to obtain and - print a summary and a plot of the results. -} - -\value{ - The function \code{TukeyC.nest} returns a list of the class \code{TukeyC.nest} with the slots: - \item{av}{A \code{list} storing the result of \code{aov}.} - \item{groups}{A vector of length equal the number of factor levels marking the groups generated.} - \item{nms}{A vector of the labels of the factor levels.} - \item{ord}{A vector which keeps the position of the means of the factor levels in decreasing order.} - \item{m.inf}{A matrix which keeps the means and dispersion of the factor levels in decreasing order.} - \item{sig.level}{A vector of length 1 giving the level of significance of the test.} - \item{r}{A vector of length 1 giving the number of replicates.} - \item{which}{The name of the factor whose levels were tested.} - \item{tab}{An array keeping the names of the factors and factor levels and also the mean value of the repetitions for every combination of factor levels.} - \item{fl1}{A vector of length 1 giving the level of the first factor in nesting order tested.} - \item{fl2}{A vector of length 1 giving the level of the second factor in nesting order tested.} -} - -\author{ - José Cláudio Faria (\email{joseclaudio.faria@gmail.com})\cr - Enio Jelihovschi (\email{eniojelihovs@gmail.com})\cr - Ivan Bezerra Allaman (\email{ivanalaman@gmail.com}) -} - -\references{ - Miller, R.G. (1981) \emph{Simultaneous Statistical Inference}. Springer. - - Ramalho M.A.P, Ferreira D.F & Oliveira A.C. (2000) \emph{Experimentação em Genética - e Melhoramento de Plantas}. Editora UFLA. - - Steel, R.G., Torrie, J.H & Dickey D.A. (1997) \emph{Principles and procedures of statistics: - a biometrical approach}. Third Edition. - - Yandell, B.S. (1997) \emph{Practical Data Analysis for Designed Experiments}. - Chapman & Hall. -} - -\examples{ - ## - ## Example: Split-split-plot Experiment (SSPE) - ## More details: demo(package='TukeyC') - ## - - data(SSPE) - ## From: design matrix (dm) and response variable (y) - ## Main factor: P - tk1 <- with(SSPE, - TukeyC(dm, - y, - model='y ~ blk + P*SP*SSP + Error(blk/P/SP)', - which='P', - error='blk:P')) - summary(tk1) - plot(tk1) - - # Main factor: SP - tk2 <- with(SSPE, - TukeyC(dm, - y, - model='y ~ blk + P*SP*SSP + Error(blk/P/SP)', - which='SP', - error='blk:P:SP', - dispersion='s')) - summary(tk2) - plot(tk2) - - # Main factor: SSP - tk3 <- with(SSPE, - TukeyC(dm, - y, - model='y ~ blk + P*SP*SSP + Error(blk/P/SP)', - which='SSP', - error='Within', - dispersion='se')) - summary(tk3) - plot(tk3) - - ## Nested: p1/SP - tkn1 <- with(SSPE, - TukeyC.nest(dm, - y, - model='y ~ blk + P*SP*SSP + Error(blk/P/SP)', - which='P:SP', - error='blk:P:SP', - fl1=1)) - summary(tkn1) - - ## From: aovlist - av <- with(SSPE, - aov(y ~ blk + P*SP*SSP + Error(blk/P/SP), - data=dfm)) - summary(av) - - ## Nested: p1/sp1/SSP - ## Testing SSP inside of level one of P and level one of SP - tkn2 <- TukeyC.nest(av, - which='P:SP:SSP', - error='Within', - fl1=1, - fl2=1) - summary(tkn2) - - ## Nested: p2/sp1/SSP - tkn3 <- TukeyC.nest(av, - which='P:SP:SSP', - error='Within', - fl1=2, - fl2=1) - summary(tkn3) -} - -\keyword{package} - diff --git a/man/m.inf.Rd b/man/m.inf.Rd deleted file mode 100644 index 5ce1b2a..0000000 --- a/man/m.inf.Rd +++ /dev/null @@ -1,72 +0,0 @@ -\name{m.inf} -\alias{m.inf} -\alias{m.inf.1a} -\alias{m.inf.1b} -\alias{m.inf.2a} -\alias{m.inf.2b} -\alias{m.inf.3a} -\alias{m.inf.3b} - -\title{ - Calculates the Means and Dispersion -} - -\description{ - Calculates the means and dispersion for one, two and three factor(s) by \code{model} and \code{model.frame} functions. -} - -\usage{ - m.inf.1a(x, - which, - dispersion=c('mm', 's', 'se')) - m.inf.1b(x, - which, - dispersion=c('mm', 's', 'se')) - m.inf.2a(x, - which1, - which2, - dispersion=c('mm', 's', 'se')) - m.inf.2b(x, - which1, - which2, - dispersion=c('mm', 's', 'se')) - m.inf.3a(x, - which1, - which2, - which3, - dispersion=c('mm', 's', 'se')) - m.inf.3b(x, - which1, - which2, - which3, - dispersion=c('mm', 's', 'se')) -} - -\arguments{ - \item{x}{A \code{SK} object.} - \item{which}{The name of the treatment to be used in the comparison. - For all the value is determined internally by the package.} - \item{which1}{The name of the treatment to be used in the comparison. - For all the value is determined internally by the package.} - \item{which2}{The name of the treatment to be used in the comparison. - For all the value is determined internally by the package.} - \item{which3}{The name of the treatment to be used in the comparison. - For all the value is determined internally by the package.} - \item{dispersion}{The dispersion to be considered to the means. \cr - The possible vaues are: \samp{'mm'} = \samp{minimum and maximum}, - \samp{'s'} = \samp{standart deviation}, - \samp{'se'} = \samp{standart deviation of the mean}.} -} - -\author{ - Enio Jelihovschi (\email{eniojelihovs@gmail.com})\cr - José Cláudio Faria (\email{joseclaudio.faria@gmail.com})\cr - Ivan Bezerra Allaman (\email{ivanalaman@gmail.com})\cr -} - -\note{ - This function is mainly for internal use in the \pkg{TukeyC} package. -} - -\keyword{m.inf} - diff --git a/man/make.TukeyC.test.Rd b/man/make.TukeyC.test.Rd index 3955109..8615401 100644 --- a/man/make.TukeyC.test.Rd +++ b/man/make.TukeyC.test.Rd @@ -10,33 +10,25 @@ } \usage{ - make.TukeyC.test(r=r, - MSE=MSE, - m.inf=m.inf, - ord=ord, - sig.level=sig.level, - dfr=dfr, - bal=bal, - mt=mt, - round=round) + make.TukeyC.test(obj, + MSE, + sig.level, + dfr, + round, + adjusted.pvalue) } \arguments{ - \item{r}{A vector of the number of repicates of each level of the factor being tested.} + \item{obj}{A data.frame with the means and replicate of the factors.} \item{MSE}{A vector of length 1 giving the mean squared error.} - \item{m.inf}{A \code{matrix} of the levels of the factor being tested in decreasing order.} - \item{ord}{A vector of ordered levels of the factor being tested.} \item{sig.level}{A vector of length 1 giving the level of significance of the test.} \item{dfr}{A vector of length 1 giving the degrees of freedom of \samp{MSE}.} - \item{bal}{A vector of length 1 giving the information wheter the experiment is or not balanced.} - \item{mt}{The model table.} \item{round}{Integer indicating the number of decimal places.} + \item{adjusted.pvalue}{A vector of pvalues obtained of the ptukey function.} } \value{ - A list with 7 slots containing the most important results of the test performed: - \item{Table}{Table of means of the factors} - \item{Means}{The means of the factors} + A list with 5 slots containing the most important results of the test performed: \item{Result}{The result of the Tukey test} \item{Sig.Level}{The significance of the test} \item{Diff_Prob}{A \code{matrix} with the observed means differences (\code{upper.tri}) and respective probalities (\code{lower.tri})} diff --git a/man/plot.TukeyC.Rd b/man/plot.TukeyC.Rd index 36d0867..f49f5fc 100644 --- a/man/plot.TukeyC.Rd +++ b/man/plot.TukeyC.Rd @@ -10,23 +10,23 @@ } \usage{ - \method{plot}{TukeyC}(x, - result=TRUE, - replicates=TRUE, - pch=19, - col=NULL, - xlab=NULL, - ylab=NULL, - xlim=NULL, - ylim=NULL, - id.lab=NULL, - id.las=1, - rl=TRUE, - rl.lty=3, - rl.col='gray', - mm=TRUE, - mm.lty=1, - title='', \dots) +\method{plot}{TukeyC}(x, + result = TRUE, + replicates = TRUE, + pch = 19, + col = NULL, + xlab = NULL, + ylab = NULL, + xlim = NULL, + ylim = NULL, + id.lab = NULL, + id.las = 1, + rl = TRUE, + rl.lty = 3, + rl.col = 'gray', + dispersion = c('none','mm','sd','ci'), + dispersion.lty = 1, + title = '', \dots) } \arguments{ @@ -44,8 +44,8 @@ \item{rl}{Horizontal line connecting the circle to the \samp{y} axis.} \item{rl.lty}{Line type of \samp{rl}.} \item{rl.col}{Line color of \samp{rl}.} - \item{mm}{Vertical line through the circle (mean value) linking the minimum to the maximum of the factor level values corresponding to that mean value.} - \item{mm.lty}{Line type of mm.} + \item{dispersion}{Vertical line through the circle (mean value) linking the minimum to the maximum of the factor level values corresponding to that mean value. Other options are: sd (standard deviation), ci (confidence interval) and none.} + \item{dispersion.lty}{Line type of dispersion.} \item{title}{A title for the plot.} \item{\dots}{Optional plotting parameters.} } @@ -53,8 +53,8 @@ \details{ The \code{plot.TukeyC} function is a S3 method to plot \samp{Tukey} and \code{TukeyC.nest} objetcs. It generates a serie of points (the means) and a - vertical line showing the minimum e maximum of the values corresponding to - each group mean. + vertical line showing the dispersion of the values corresponding to + each group mean. } \author{ @@ -72,56 +72,55 @@ } \examples{ - ## - ## Examples: Completely Randomized Design (CRD) - ## More details: demo(package='TukeyC') - ## - - library(TukeyC) - data(CRD2) - - ## From: vectors x and y - tk1 <- with(CRD2, - TukeyC(x=x, - y=y, - model='y ~ x', - which='x')) - plot(tk1, - id.las=2, - rl=FALSE) - - ## From: design matrix (dm) and response variable (y) - tk2 <- with(CRD2, - TukeyC(x=dm, - y=y, - model='y ~ x', - which='x')) - plot(tk2, - mm.lty=3, - id.las=2, - rl=FALSE) - - ## From: data.frame (dfm) - tk3 <- with(CRD2, - TukeyC(x=dfm, - model='y ~ x', - which='x')) - plot(tk3, - id.las=2, - rl=FALSE) - - ## From: aov - av <- with(CRD2, - aov(y ~ x, - data=dfm)) - summary(av) - - tk4 <- with(CRD2, - TukeyC(x=av, - which='x')) - plot(tk4, - rl=FALSE, - id.las=2) +## +## Examples: Completely Randomized Design (CRD) +## More details: demo(package='TukeyC') +## + +library(TukeyC) +data(CRD2) + +## From: formula +tk1 <- with(CRD2, + TukeyC(y ~ x, + data=dfm, + which='x')) + +old.par <- par(mar=c(6, 3, 6, 2)) +plot(tk1, + id.las=2) + +plot(tk1, + dispersion='sd', + mm.lty=3, + id.las=2, + rl=FALSE) + +## From: aov +av <- with(CRD2, + aov(y ~ x, + data=dfm)) +summary(av) + +tk2 <- TukeyC(x=av, + which='x') +plot(tk2, + rl=FALSE, + id.las=2) + +# From: lm +av_lm <- with(CRD2, + lm(y ~ x, + data=dfm)) + +tk3 <- TukeyC(x=av_lm, + which='x') +plot(tk3, + dispersion='ci', + id.las=2, + rl=FALSE) + +par(old.par) } \keyword{package} diff --git a/man/print.TukeyC.Rd b/man/print.TukeyC.Rd new file mode 100644 index 0000000..334bd9d --- /dev/null +++ b/man/print.TukeyC.Rd @@ -0,0 +1,43 @@ +\name{print.TukeyC} +\alias{print.TukeyC} + +\title{ + Print Method for \code{TukeyC} objects. +} + +\description{ + Returns (and prints) a list for objects of class \code{TukeyC}. +} + +\usage{ + +\method{print}{TukeyC}(x, ...) + +} + +\arguments{ + \item{x}{A given object of the class \code{TukeyC}.} + \item{\dots }{Further arguments (require by generic).} +} + +\author{ + Enio G. Jelihovschi (\email{eniojelihovs@gmail.com})\cr + Ivan Bezerra Allaman (\email{ivanalaman@gmail.com})\cr +} + +\seealso{\code{\link{TukeyC}} +} + +\examples{ +data(RCBD) + +tk <- with(RCBD, + TukeyC(y ~ blk + tra, + data=dfm, + which='tra')) +tk + +} + +\keyword{package} +\keyword{TukeyC} diff --git a/man/summary.Rd b/man/summary.Rd index d578099..218d9c5 100644 --- a/man/summary.Rd +++ b/man/summary.Rd @@ -41,15 +41,13 @@ ## More details: demo(package='TukeyC') ## - ## The parameters can be: vectors, design matrix and the response variable, - ## data.frame or aov + ## The parameters can be: formula, aov, lm and aovlist data(CRD2) - ## From: design matrix (dm) and response variable (y) + ## From: formula tk1 <- with(CRD2, - TukeyC(x=dm, - y=y, - model='y ~ x', + TukeyC(y ~ x, + data=dfm, which='x', id.trim=5)) summary(tk1) @@ -59,18 +57,15 @@ ## More details: demo(package='TukeyC') ## - ## The parameters can be: design matrix and the response variable, - ## data.frame or aov - + ## The parameters can be: formula, aov, lm and aovlist data(RCBD) - ## Design matrix (dm) and response variable (y) - tk1 <- with(RCBD, - TukeyC(x=dm, - y=y, - model='y ~ blk + tra', + ## From: formula + tk2 <- with(RCBD, + TukeyC(y ~ blk + tra, + data=dfm, which='tra')) - summary(tk1) + summary(tk2) ## ## Example: Latin Squares Design (LSD) @@ -82,13 +77,12 @@ data(LSD) - ## From: design matrix (dm) and response variable (y) - tk1 <- with(LSD, - TukeyC(x=dm, - y=y, - model='y ~ rows + cols + tra', + ## From: formula + tk3 <- with(LSD, + TukeyC(y ~ rows + cols + tra, + data=dfm, which='tra')) - summary(tk1) + summary(tk3) ## ## Example: Factorial Experiment (FE) @@ -101,52 +95,47 @@ data(FE) ## From: design matrix (dm) and response variable (y) ## Main factor: N - tk1 <- with(FE, - TukeyC(x=dm, - y=y, - model='y ~ blk + N*P*K', + tk4 <- with(FE, + TukeyC(y ~ blk + N*P*K, + data=dfm, which='N')) - summary(tk1) + summary(tk4) ## Nested: p1/N ## Testing N inside of level one of P ntk1 <- with(FE, - TukeyC.nest(x=dm, - y=y, - model='y ~ blk + N*P*K', - which='P:N', - fl1=1)) + TukeyC(y ~ blk + N*P*K, + data=dfm, + which='P:N', + fl1=1)) summary(ntk1) ## Nested: k1/p1/N ## Testing N inside of level one of K and level one of P ntk2 <- with(FE, - TukeyC.nest(x=dm, - y=y, - model='y ~ blk + N*P*K', - which='K:P:N', - fl1=1, - fl2=1)) + TukeyC(y ~ blk + N*P*K, + data=dfm, + which='K:P:N', + fl1=1, + fl2=1)) summary(ntk2) ## Nested: k2/n2/P ntk3 <- with(FE, - TukeyC.nest(x=dm, - y=y, - model='y ~ blk + N*P*K', - which='K:N:P', - fl1=2, - fl2=2)) + TukeyC(y ~ blk + N*P*K, + data=dfm, + which='K:N:P', + fl1=2, + fl2=2)) summary(ntk3) ## Nested: p1/n1/K ntk4 <- with(FE, - TukeyC.nest(x=dm, - y=y, - model='y ~ blk + N*P*K', - which='P:N:K', - fl1=1, - fl2=1)) + TukeyC(y ~ blk + N*P*K, + data=dfm, + which='P:N:K', + fl1=1, + fl2=1)) summary(ntk4) ## @@ -155,91 +144,88 @@ ## data(SPE) - - ## The parameters can be: design matrix and the response variable, - ## data.frame or aov - - ## From: design matrix (dm) and response variable (y) + ## From: formula ## Main factor: P + ## It is necessary to inform the appropriate error for the test tk1 <- with(SPE, - TukeyC(x=dm, - y=y, - model='y ~ blk + P*SP + Error(blk/P)', + TukeyC(y ~ blk + P*SP + Error(blk/P), + data=dfm, which='P', error='blk:P')) summary(tk1) ## Nested: p1/SP tkn1 <- with(SPE, - TukeyC.nest(x=dm, - y=y, - model='y ~ blk + P*SP + Error(blk/P)', - which='P:SP', - error='Within', - fl1=1 )) + TukeyC(y ~ blk + P*SP + Error(blk/P), + data=dfm, + which='P:SP', + fl1=1 )) summary(tkn1) - data(SSPE) - - ## From: design matrix (dm) and response variable (y) + ## From: formula ## Main factor: P + ## It is necessary to inform the appropriate error for the test + data(SSPE) + tk1 <- with(SSPE, - TukeyC(dm, - y, - model='y ~ blk + P*SP*SSP + Error(blk/P/SP)', + TukeyC(y ~ blk + P*SP*SSP + Error(blk/P/SP), + data=dfm, which='P', error='blk:P')) summary(tk1) - - # Main factor: SP + + ## Main factor: SP + ## It is necessary to inform the appropriate error for the test tk2 <- with(SSPE, - TukeyC(dm, - y, - model='y ~ blk + P*SP*SSP + Error(blk/P/SP)', + TukeyC(y ~ blk + P*SP*SSP + Error(blk/P/SP), + data=dfm, which='SP', error='blk:P:SP')) summary(tk2) - - # Main factor: SSP + + ## Main factor: SSP tk3 <- with(SSPE, - TukeyC(dm, - y, - model='y ~ blk + P*SP*SSP + Error(blk/P/SP)', - which='SSP', - error='Within')) + TukeyC(y ~ blk + P*SP*SSP + Error(blk/P/SP), + data=dfm, + which='SSP')) summary(tk3) ## Nested: p1/SSP tkn1 <- with(SSPE, - TukeyC.nest(dm, - y, - model='y ~ blk + P*SP*SSP + Error(blk/P/SP)', - which='P:SSP', - error='blk:P:SP', - fl1=1)) + TukeyC(y ~ blk + P*SP*SSP + Error(blk/P/SP), + data=dfm, + which='P:SSP', + fl1=1)) summary(tkn1) ## From: aovlist av <- with(SSPE, aov(y ~ blk + P*SP*SSP + Error(blk/P/SP), - data=dfm)) + data=dfm)) summary(av) ## Nested: P1/SP1/SSP - tkn2 <- TukeyC.nest(av, - which='P:SP:SSP', - error='Within', - fl1=1, - fl2=1) + tkn2 <- TukeyC(av, + which='P:SP:SSP', + fl1=1, + fl2=1) summary(tkn2) ## Nested: P2/SP1/SSP - tkn3 <- TukeyC.nest(av, - which='P:SP:SSP', - error='Within', - fl1=2, - fl2=1) + tkn3 <- TukeyC(av, + which='P:SP:SSP', + fl1=2, + fl2=1) summary(tkn3) + + ## Nested: SSP2/P1/SP - it is necessary to inform how to combinate the errors + tkn4 <- TukeyC(av, + which='SSP:P:SP', + fl1=2, + fl2=1, + error='Within/blk:P/blk:P:SP') + summary(tkn4) + } \keyword{package}