Articulo usado: http://www.scielo.org.co/pdf/rion/v31n1/0120-100X-rion-31-01-43.pdf

Hacer un analisis de covarianza

set.seed(123)

# covariable: peso inicial de remolacha 
per<-sort.int(runif(120,90,170))

# Bloqueo 
labo<- gl(2,10,120,c("biol","alim"))

# Factor 1: Tiempo de coccion
tiemp = gl(3,40,120,c(0,10,12))

# Factor 2: Variedades
varie = gl(2,20,120,c("var1","var2"))

# Variable respuesta: compuestos fenolicos 
cf = c(rnorm(40,20,1),
       rnorm(40,12,1.2),
       rnorm(40,16,1.15))

df2<- data.frame(per,labo,tiemp,varie,cf); head(df2)
##        per labo tiemp varie       cf
## 1 90.04998 biol     0  var1 20.37964
## 2 91.96909 biol     0  var1 19.49768
## 3 93.36476 biol     0  var1 19.66679
## 4 93.64452 biol     0  var1 18.98142
## 5 93.66649 biol     0  var1 18.92821
## 6 94.85765 biol     0  var1 20.30353
library(collapsibleTree)
## Warning: package 'collapsibleTree' was built under R version 4.0.4
collapsibleTreeSummary(
  df2,
  hierarchy = c('labo','varie','tiemp', "cf"))

Ancova modelo planteado

mod1<- aov(cf ~ labo + varie + tiemp + per); summary(mod1)
##              Df Sum Sq Mean Sq F value Pr(>F)    
## labo          1    0.1     0.1   0.068 0.7950    
## varie         1    5.4     5.4   4.686 0.0325 *  
## tiemp         2 1399.9   700.0 608.659 <2e-16 ***
## per           1    1.3     1.3   1.145 0.2869    
## Residuals   114  131.1     1.2                   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
res1<- mod1$residuals
# prueba de normalidad 
shapiro.test(res1)
## 
##  Shapiro-Wilk normality test
## 
## data:  res1
## W = 0.98861, p-value = 0.4184
plot(res1)

plot(mod1$fitted.values, res1)

Ancova con interaccion

mod2<- aov(cf ~ labo * varie * tiemp * per, df2); summary(mod2)
##                      Df Sum Sq Mean Sq F value Pr(>F)    
## labo                  1    0.1     0.1   0.062 0.8043    
## varie                 1    5.4     5.4   4.264 0.0416 *  
## tiemp                 2 1399.9   700.0 553.836 <2e-16 ***
## per                   1    1.3     1.3   1.042 0.3100    
## labo:varie            1    0.1     0.1   0.102 0.7502    
## labo:tiemp            2    1.0     0.5   0.386 0.6809    
## varie:tiemp           2    0.0     0.0   0.018 0.9826    
## labo:per              1    2.9     2.9   2.265 0.1356    
## varie:per             1    0.4     0.4   0.332 0.5656    
## tiemp:per             2    1.4     0.7   0.568 0.5683    
## labo:varie:tiemp      2    2.6     1.3   1.047 0.3549    
## labo:varie:per        1    0.0     0.0   0.000 0.9835    
## labo:tiemp:per        2    0.0     0.0   0.006 0.9940    
## varie:tiemp:per       2    0.7     0.3   0.272 0.7628    
## labo:varie:tiemp:per  2    0.6     0.3   0.219 0.8036    
## Residuals            96  121.3     1.3                   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
res2<- mod2$residuals
# prueba de normalidad 
shapiro.test(res2)
## 
##  Shapiro-Wilk normality test
## 
## data:  res2
## W = 0.9795, p-value = 0.06397

#Anova

mod3<- aov(cf ~ tiemp*varie,df2); summary(mod3)
##              Df Sum Sq Mean Sq F value Pr(>F)    
## tiemp         2 1399.9   700.0 602.599 <2e-16 ***
## varie         1    5.4     5.4   4.639 0.0334 *  
## tiemp:varie   2    0.1     0.0   0.033 0.9678    
## Residuals   114  132.4     1.2                   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
res3<- mod3$residuals
# prueba de normalidad 
shapiro.test(res3)
## 
##  Shapiro-Wilk normality test
## 
## data:  res3
## W = 0.98475, p-value = 0.1949
mod4<- aov(cf ~ labo*varie + tiemp, df2); summary(mod4)
##              Df Sum Sq Mean Sq F value Pr(>F)    
## labo          1    0.1     0.1   0.067 0.7959    
## varie         1    5.4     5.4   4.642 0.0333 *  
## tiemp         2 1399.9   700.0 602.906 <2e-16 ***
## labo:varie    1    0.1     0.1   0.056 0.8127    
## Residuals   114  132.4     1.2                   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
res4<- mod4$residuals
# prueba de normalidad 
shapiro.test(res4)
## 
##  Shapiro-Wilk normality test
## 
## data:  res4
## W = 0.98468, p-value = 0.1921
# Eficiencia bloqueo 
H<-0.1/132.4;H
## [1] 0.000755287

Con los datos anteriones se puede observar que el peso de los cubos de remolacha, osea nuestra covariable si influye en la concentracion de fenoles