set.seed(123)
# Variable: Compuestos Fenolicos
cf = c(rnorm(40, 20, 1),
       rnorm(40, 12, 1.2),
       rnorm(40, 16, 1.15))

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

# Factor 2: Vierdades
varie = gl(2, 20, 120, c('var1', 'var2'))

#Covariable : Pesos
pes = runif(120,100,200)
pes = round(sort.int(pes, 6), 2)

df = data.frame(varie, tiemp, cf,pes)
head(df)
##   varie tiemp       cf    pes
## 1  var1     0 19.43952 100.63
## 2  var1     0 19.76982 102.00
## 3  var1     0 21.55871 100.94
## 4  var1     0 20.07051 100.82
## 5  var1     0 20.12929 103.37
## 6  var1     0 21.71506 103.49
mod1 = aov(cf ~ tiemp+varie+pes:pes, df)
summary(mod1)
##              Df Sum Sq Mean Sq F value Pr(>F)    
## tiemp         2 1297.1   648.5 663.336 <2e-16 ***
## varie         1    5.7     5.7   5.844 0.0172 *  
## pes           1    1.8     1.8   1.799 0.1824    
## Residuals   115  112.4     1.0                   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

#La covariable segun la tabla de covarianzas si influye en la concentracion de fenoles

mod1 = aov(cf ~ tiemp*varie*pes, df)
summary(mod1)
##                  Df Sum Sq Mean Sq F value Pr(>F)    
## tiemp             2 1297.1   648.5 643.586 <2e-16 ***
## varie             1    5.7     5.7   5.670  0.019 *  
## pes               1    1.8     1.8   1.746  0.189    
## tiemp:varie       2    2.8     1.4   1.396  0.252    
## tiemp:pes         2    0.4     0.2   0.176  0.839    
## varie:pes         1    0.2     0.2   0.203  0.653    
## tiemp:varie:pes   2    0.2     0.1   0.114  0.892    
## Residuals       108  108.8     1.0                   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1