Y<-c(69.38,
69.72,
69.58,
69.50,
69.48,
69.56,
69.90,
69.60,
69.80,
69.70,
69.50,
69.40,
69.40,
70.02,
69.62,
69.78,
69.70,
69.46,
69.50,
69.68,
69.94,
69.52,
69.90,
69.92,
69.50,
69.42,
69.64,
69.88)
Piezas<-rep(1:7,4)
Inspector<-rep(1:2,each=14)
df<-data.frame(Piezas=as.factor(Piezas),Inspector=as.factor(Inspector),Y)
library(lme4)
## Loading required package: Matrix
modelo<-lmer(Y~(1|Piezas)+(1|Inspector)+(1|Piezas:Inspector),data=df)
#modelo<-lm(Y~Piezas*Inspector,data=df)
summary(modelo)
## Linear mixed model fit by REML ['lmerMod']
## Formula: Y ~ (1 | Piezas) + (1 | Inspector) + (1 | Piezas:Inspector)
##    Data: df
## 
## REML criterion at convergence: -33.9
## 
## Scaled residuals: 
##      Min       1Q   Median       3Q      Max 
## -1.69924 -0.55074  0.06051  0.44442  2.00487 
## 
## Random effects:
##  Groups           Name        Variance Std.Dev.
##  Piezas:Inspector (Intercept) 0.000331 0.01819 
##  Piezas           (Intercept) 0.030600 0.17493 
##  Inspector        (Intercept) 0.001598 0.03997 
##  Residual                     0.007200 0.08485 
## Number of obs: 28, groups:  Piezas:Inspector, 14; Piezas, 7; Inspector, 2
## 
## Fixed effects:
##             Estimate Std. Error t value
## (Intercept) 69.64286    0.07383   943.3
with(df, interaction.plot(x.factor = Piezas, trace.factor = Inspector, 
                              response = Y))

with(df, interaction.plot(x.factor = Inspector, trace.factor = Piezas, 
                              response = Y))

modelo1<-lm(Y~Piezas*Inspector,data=df)
fit<-aov(modelo1)
summary(fit)
##                  Df Sum Sq Mean Sq F value   Pr(>F)    
## Piezas            6 0.7816 0.13026  18.092 7.41e-06 ***
## Inspector         1 0.0302 0.03023   4.198   0.0597 .  
## Piezas:Inspector  6 0.0472 0.00786   1.092   0.4137    
## Residuals        14 0.1008 0.00720                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Ejemplo de la Web

Let us now consider the analysis about an inheritance study with beef animals. Five sires (male animals) were each mated to a separate group of dams (female animals). Birth weight of eight male calves (from different dams) in each of the five sire groups were recorded. We want to find out what part of the total variation is due to variation between different sires.

We first create the data set and visualize it.

weight <- c(61, 100,  56, 113,  99, 103,  75,  62,  ## sire 1
            75, 102,  95, 103,  98, 115,  98,  94,  ## sire 2
            58,  60,  60,  57,  57,  59,  54, 100,  ## sire 3
            57,  56,  67,  59,  58, 121, 101, 101,  ## sire 4
            59,  46, 120, 115, 115,  93, 105,  75)  ## sire 5
sire    <- factor(rep(1:5, each = 8))
animals <- data.frame(weight, sire)
str(animals)
## 'data.frame':    40 obs. of  2 variables:
##  $ weight: num  61 100 56 113 99 103 75 62 75 102 ...
##  $ sire  : Factor w/ 5 levels "1","2","3","4",..: 1 1 1 1 1 1 1 1 2 2 ...
library(lme4)
fit.sire <- lmer(weight ~ (1 | sire), data = animals)
summary(fit.sire)
## Linear mixed model fit by REML ['lmerMod']
## Formula: weight ~ (1 | sire)
##    Data: animals
## 
## REML criterion at convergence: 358.2
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -1.9593 -0.7459 -0.1581  0.8143  1.9421 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  sire     (Intercept) 116.7    10.81   
##  Residual             463.8    21.54   
## Number of obs: 40, groups:  sire, 5
## 
## Fixed effects:
##             Estimate Std. Error t value
## (Intercept)   82.550      5.911   13.96