Tres sucedaneos com diferente concentracao de proteina (tratamentos) foram testados em bezerros da raca holandes:

O experimento consistiu de um delineamento em blocos (9) completos casualizados.

Session

pkg <- c("agricolae", "MASS", "ggplot2", "Rmisc", "plyr")
sapply(pkg, library, character.only=TRUE, logical.return=TRUE)
## agricolae      MASS   ggplot2     Rmisc      plyr 
##      TRUE      TRUE      TRUE      TRUE      TRUE
#setwd("/media/DATA/Dropbox/Varios-BRA/Vachi/vachi_R/exp2") #PC Juanchi
#source('/media/DATA/Dropbox/MyR/Sources/boxplot.with.outlier.label.r', encoding='UTF-8') # pc juanchi

setwd("~/Dropbox/MyR/Análise otros/vachi_R/exp2")#lab juanchi
source('~/Dropbox/MyR/Sources/boxplot.with.outlier.label.r', encoding='UTF-8') # lab

#setwd("C:/Users/evangelina/Dropbox/Varios-BRA/Vachi/vachi_R/exp2") #pc vachi
#source('C:/Users/evangelina/Dropbox/Varios-BRA/Vachi/vachi_R/exp2/boxplot.with.outlier.label.r', encoding='UTF-8')

Dataset

#list.files()
dat = read.csv("d_diarreia e vida.csv", h=T, s=";", d=",")

Exploracao

##   ani trat bloc dd dv
## 1   2    B    2 37 56
## 2   3    C    1 23 56
## 3   4    C    3  4 17
## 4   5    B    3 40 56
## 5   6    A    1 35 56
## 6   7    A    3 24 56
## 'data.frame':    33 obs. of  5 variables:
##  $ ani : int  2 3 4 5 6 7 8 9 10 11 ...
##  $ trat: Factor w/ 3 levels "A","B","C": 2 3 3 2 1 1 2 3 1 2 ...
##  $ bloc: int  2 1 3 3 1 3 1 2 2 6 ...
##  $ dd  : int  37 23 4 40 35 24 39 48 41 39 ...
##  $ dv  : int  56 56 17 56 56 56 56 56 56 56 ...
##     bloc
## trat 1 2 3 4 5 6 7 8 9
##    A 1 1 1 1 2 1 1 2 1
##    B 1 1 1 1 2 1 1 1 1
##    C 1 1 1 1 3 1 1 2 1
## Loading required package: TeachingDemos

Anova

# Ajuste de modelo
mod = aov(dd~ trat + bloc, data=dat[-12,])

# Diagnostico
par(mfrow=c(2,2))
plot(mod, which=1:2); boxcox(mod)
car::influencePlot(mod); layout(1)

##       StudRes       Hat     CookD
## 14  2.2221839 0.3936189 0.4953667
## 15 -2.4034049 0.1871486 0.3138464
## 17  0.8999226 0.5437751 0.2975802
anova(mod) # cuadro anova
## Analysis of Variance Table
## 
## Response: dd
##           Df Sum Sq Mean Sq F value   Pr(>F)   
## trat       2 1828.9  914.46  7.0200 0.004628 **
## bloc       8 1121.7  140.22  1.0764 0.415926   
## Residuals 21 2735.6  130.26                    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# comparaciones multiples
TukeyHSD(x=mod, 'trat', conf.level=0.95)
##   Tukey multiple comparisons of means
##     95% family-wise confidence level
## 
## Fit: aov(formula = dd ~ trat + bloc, data = dat[-12, ])
## 
## $trat
##           diff        lwr        upr     p adj
## B-A   7.484848  -5.445494 20.4151912 0.3301994
## C-A -11.015152 -23.023674  0.9933713 0.0760201
## C-B -18.500000 -31.185582 -5.8144178 0.0038523
(tuk = HSD.test(mod, 'trat', group=TRUE))
## $statistics
##       Mean       CV MSerror      HSD r.harmonic
##   35.15625 32.46469 130.265 12.54754   10.51327
## 
## $parameters
##   Df ntr StudentizedRange
##   21   3         3.564625
## 
## $means
##         dd       std  r Min Max
## A 37.18182  8.244006 11  24  53
## B 44.66667  6.204837  9  37  54
## C 26.16667 16.151743 12   3  50
## 
## $comparison
## NULL
## 
## $groups
##   trt    means  M
## 1   B 44.66667  a
## 2   A 37.18182 ab
## 3   C 26.16667  b

Plots

(sumdat <- summarySEwithin(dat[-12,], measurevar="dd", withinvars="trat",
                           na.rm=FALSE, conf.interval=.95))
##   trat  N       dd        sd       se        ci
## 1    A 11 37.18182 10.096804 3.044301  6.783125
## 2    B  9 44.66667  7.599342 2.533114  5.841371
## 3    C 12 26.16667 19.781764 5.710503 12.568733
df = data.frame(
let = tuk$groups[with(tuk$groups, order(trt)), ],
me = tuk$means$dd,
se = tuk$means$std/sqrt(tuk$means$r))
df
##   let.trt let.means let.M       me       se
## 2       A  37.18182    ab 37.18182 2.485661
## 1       B  44.66667     a 44.66667 2.068279
## 3       C  26.16667     b 26.16667 4.662606
(p1<-ggplot(df, aes(let.trt, me)) + theme_bw() +
    geom_point(size=3) + 
    geom_errorbar(aes(x = let.trt, ymin = me-se, 
                      ymax = me+se),colour="black", width=.1) +  
     geom_line() +
     xlab("Sucedáneos") + ylab("Dias com diarreia")  +
     annotate("text", x = (1:3)-0.1, y = sumdat$dd,
              label = df$let.M)
 
)
## geom_path: Each group consist of only one observation. Do you need to adjust the group aesthetic?