> chem.mod0 <- aov(serum ~ method + day + method:day)
> summary(chem.mod0)
            Df Sum Sq Mean Sq F value   Pr(>F)   
method       1 329.42  329.42 22.8488 0.001391 **
day          3 431.44  143.81  9.9750 0.004442 **
method:day   3 184.83   61.61  4.2733 0.044602 * 
Residuals    8 115.34   14.42                    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 
> 
> chem.mod1 <- aov(serum ~ method + day + Error(method:day))
Warning message:
In aov(serum ~ method + day + Error(method:day)) :
  Error() model is singular
> summary(chem.mod1)

Error: method:day
          Df Sum Sq Mean Sq F value Pr(>F)
method     1 329.42  329.42  5.3468 0.1038
day        3 431.44  143.81  2.3342 0.2522
Residuals  3 184.83   61.61               

Error: Within
          Df Sum Sq Mean Sq F value Pr(>F)
Residuals  8 115.34  14.418               
> 
> library(nlme)
> chem.mod2 <- lme(fixed=serum~method,random=~1|day/method)
> 
> summary(chem.mod2)
Linear mixed-effects model fit by REML
 Data: NULL 
       AIC      BIC    logLik
  102.5048 105.7001 -46.25239

Random effects:
 Formula: ~1 | day
        (Intercept)
StdDev:    4.533303

 Formula: ~1 | method %in% day
        (Intercept) Residual
StdDev:     4.85764 3.797038

Fixed effects: serum ~ method 
              Value Std.Error DF  t-value p-value
(Intercept) 147.000  3.583164  8 41.02520  0.0000
method2      -9.075  3.924628  3 -2.31232  0.1038
 Correlation: 
        (Intr)
method2 -0.548

Standardized Within-Group Residuals:
        Min          Q1         Med          Q3         Max 
-1.58604086 -0.53477901 -0.01443907  0.23216813  1.41629914 

Number of Observations: 16
Number of Groups: 
            day method %in% day 
              4               8 
> intervals(chem.mod2)
Approximate 95% confidence intervals

 Fixed effects:
                lower    est.      upper
(Intercept) 138.73721 147.000 155.262790
method2     -21.56492  -9.075   3.414919
attr(,"label")
[1] "Fixed effects:"

 Random Effects:
  Level: day 
                    lower     est.    upper
sd((Intercept)) 0.9886142 4.533303 20.78752
  Level: method 
                   lower    est.    upper
sd((Intercept)) 1.690910 4.85764 13.95501

 Within-group standard error:
   lower     est.    upper 
2.326189 3.797038 6.197904 
> anova(chem.mod2)
            numDF denDF   F-value p-value
(Intercept)     1     8 2257.9766  <.0001
method          1     3    5.3468  0.1038
> 
> library(lme4)
> chem.mod3 <- lmer(serum~method+(1|day)+(1|method:day))
> summary(chem.mod3)
Linear mixed model fit by REML 
Formula: serum ~ method + (1 | day) + (1 | method:day) 
   AIC   BIC logLik deviance REMLdev
 102.5 106.4 -46.25    101.0    92.5
Random effects:
 Groups     Name        Variance Std.Dev.
 method:day (Intercept) 23.597   4.8576  
 day        (Intercept) 20.551   4.5333  
 Residual               14.417   3.7970  
Number of obs: 16, groups: method:day, 8; day, 4

Fixed effects:
            Estimate Std. Error t value
(Intercept)  147.000      3.583   41.03
method2       -9.075      3.925   -2.31

Correlation of Fixed Effects:
        (Intr)
method2 -0.548
> anova(chem.mod3)
Analysis of Variance Table
       Df Sum Sq Mean Sq F value
method  1 77.088  77.088  5.3468