Skip to content

ANCOVA with defined error terms

14 messages · hpdutra, Richard M. Heiberger

Original
hpdutra
# 
I am trying to run an ANCOVA with defined error terms.  Thus I have to use
AOV and not lm. 
my response variable is proportion of mice paw prints on track plates. These
plates were placed on plots that had vegetation and fruit manipulated to two
levels each (present or absent), and were sampled monthly for 14 months
(repeated measures). The fully crossed factor design was also blocked. My
sample size is just 3 and I cannot run the full model because I would ran
out of degrees of freedom. 
I measured raccoon activity on these plots as a covariate. So I am running
the following model.


summary(model<-(aov(mice~veget*fruit*time+(time*block)+(time*block*veget)+(time*block*fruit)+coon+Error(block/plot,
data = track))))

here is a sample of the data 
plot	veget	        fruit	     time	block 	coon     	        mice 
p1	Vremoved	Fintact	1	b1      	8.605276544	26.67179738
p2	Vintact	Fintact	1	b1 	     	16.64929378	45
p3	Vremoved	Fremoved	1	b1 	     	22.12023855	26.67179738
p4	Vintact	Fremoved	1	b1 	     	16.64929378	41.57117579
p5	Vintact	Fremoved	1	b2 	     	30.73522506	38.09196076
p6	Vremoved	Fintact	1	b2 	     	22.12023855	26.67179738
p7	Vintact	Fintact	1	b2 	     	8.605276544	8.605276544
p8	Vremoved	Fremoved	1	b2 	     	8.605276544	8.605276544
p9	Vintact	Fintact	1	b3 	     	8.605276544	16.64929378
p10	Vintact	Fremoved	1	b3 	        16.64929378	8.605276544
p11	Vremoved	Fintact	1	b3     	8.605276544	30.73522506
p12	Vremoved	Fremoved	1	b3  	     	8.605276544	16.64929378
p1	Vremoved	Fintact	2	b1       	22.12023855	41.57117579
p2	Vintact	Fintact	2	b1  	     	8.605276544	48.42882421

My question is. Is there a more simple way to run this mess using lm or not?
Richard M. Heiberger
# 
track.aov <- aov(mice ~ coon
                  + block*veget*fruit*time - block:veget:fruit:time
                  + Error(block/plot), data = track)
anova(track.aov)


I think this is what you are looking for.  This model in words says,
What is the effect of the four-way crossing after adjusting for
the covariate coon?  We don't have enough degrees of freedom for
the full crossing, so remove the four-way interaction.

anova() gives the sequential sums of squares of (multi-degree of
freedom) effects.  Be sure that the five variables (block, veget, fruit,
time, plot) are factors.  I am assuming time is a factor here.
hpdutra
# 
Thanks Richard,

I tried running the analysis the way you suggested but here is the error
that I get
+                 + Error(block/plot), data = track)
Warning message:
In aov(kotz.mice ~ kotz.coon + block * veget * fruit * time -
block:veget:fruit:time +  :
  Error() model is singular
Error in UseMethod("anova") : no applicable method for "anova"

I'm not sure you can use the function anova for an aov object. I tried
summary and it also didn't work. Did you mean something else?
Thanks
Richard M. Heiberger wrote:

  
    

      
      
    

  


  


      

    

    

      
      

        


  

        

      Richard M. Heiberger
    

    

      
      #
    

  


  

      
Yes, I meant summary().  anova() isn't defined for aovlist objects and 
summary() is.

            

      
      
      
    

        
You will need to investigate the singular Error() model.  You might want
to use the simpler Error model

Error(block)

My guess is the plot term is redundant.

What are kotz.mice?  That variable isn't in the model you showed.

For more detailed help you will need to
"provide commented, minimal, self-contained, reproducible code."


Rich

  
  


  


      

    

    
2 days later

    

      
      

        


  

        

      hpdutra
    

    

      
      #
    

  


  

      
Thank you Richard, this works. But, the model you suggested me lacks some
between subjects  interactions, namely:
 veget:time:block
 fruit:time:block. 
 According to Sokal and Rohlf I need to report those as well. 
Also, any ideas why the sum of squares on my model are different from yours

summary(mymodel<-(aov(mice~veget*fruit*time+(time*block)+(time*block*veget)+(time*block*fruit)+coon+Error(block/plot,
data = track))))
reports veg:fruit with a sum squares of 80.5

while yours 
summary <- (ricsmodel<-aov(mice ~ coon + block*veget*fruit*time -
block:veget:fruit:time+ Error(block/plot), data = track))
gives me a veget:fruit sum of squares of 253.2
the residuals degrees of freedom also differ my model gives 4 DF for the
between subjects while yourr does not report one for the between subjects. 
I know I'm asking too much, but I'm just trying to understand.
Thanks
Beto 
s

        
Richard M. Heiberger wrote:

            

      
      
      
    

        

      
      
    

        

  
    

      
      
    

  


  


      

    

    

      
      

        


  

        

      Richard M. Heiberger
    

    

      
      #
    

  


  

      
The three-way interactions you mention are included in the model formula
I suggested.  If they didn't appear in the expansion, it suggests
that you have some aliasing due to empty cells.

I can't do any more without your dataset.
You can post your dataset with random response values.
The exact data.frame for the predictors is needed.
Anything that is a factor must be a factor in what you send.
You can use dput to ensure accurate copying.


 > tmp <- data.frame(a=factor(rep(letters[1:3], 2)), y=rnorm(6))
 > tmp
   a          y
1 a -0.4313252
2 b -0.9065241
3 c -0.7285257
4 a  0.0368019
5 b  0.7982373
6 c -2.4712612
 > dput(tmp)
structure(list(a = structure(c(1L, 2L, 3L, 1L, 2L, 3L), .Label = c("a",
"b", "c"), class = "factor"), y = c(-0.431325155041393, -0.906524086347679,
-0.728525691910586, 0.0368018971284965, 0.798237317168781, -2.47126116272324
)), .Names = c("a", "y"), row.names = c(NA, -6L), class = "data.frame")
 >
 >
 > structure(list(a = structure(c(1L, 2L, 3L, 1L, 2L, 3L), .Label = c("a",
+ "b", "c"), class = "factor"), y = c(-0.431325155041393, 
-0.906524086347679,
+ -0.728525691910586, 0.0368018971284965, 0.798237317168781, 
-2.47126116272324
+ )), .Names = c("a", "y"), row.names = c(NA, -6L), class = "data.frame")
   a          y
1 a -0.4313252
2 b -0.9065241
3 c -0.7285257
4 a  0.0368019
5 b  0.7982373
6 c -2.4712612
 >

  
  


  


      

    

    

      
      

        


  

        

      hpdutra
    

    

      
      #
    

  


  

      
Hi Richard, there are no empty cells. 
I transform everything into factor, except the co-variate coon. 
Here is the full analysis with dput of the data.
I'm afraid I have not enough DF for the thre-way interaction using your
model as well. 12 plots divided in 3 blocks, each plot assigned to 2 crossed
treatments (veget and fruit two levels each), which means that my n=3. Plots
were sampled monthly for 14 months. 
I also attached the raw data, which is proportions of plates with paw
prints, but I had to use Kotz transformation to improve normality.
Thank you very much 

 track<-read.table(file.choose(), h=T)

            

      
      
      
    

        
Error: block
      Df Sum Sq Mean Sq
block  2 4581.5  2290.7

Error: block:plot
            Df  Sum Sq Mean Sq  F value  Pr(>F)  
veget        1 1602.91 1602.91 379.1653 0.03267 *
fruit        1   80.07   80.07  18.9411 0.14378  
coon         1  193.53  193.53  45.7791 0.09341 .
veget:fruit  1   40.26   40.26   9.5245 0.19948  
veget:block  2  355.53  177.76  42.0495 0.10840  
fruit:block  2   15.39    7.69   1.8200 0.46423  
Residuals    1    4.23    4.23                   
---
Signif. codes:  0 ?***? 0.001 ?**? 0.01 ?*? 0.05 ?.? 0.1 ? ? 1 

Error: Within
                 Df Sum Sq Mean Sq F value    Pr(>F)    
time             13 8451.5   650.1  8.7462 2.384e-06 ***
coon              1 1070.7  1070.7 14.4047 0.0008364 ***
veget:time       13 1214.4    93.4  1.2568 0.3006216    
fruit:time       13  474.7    36.5  0.4912 0.9090773    
time:block       26 2943.4   113.2  1.5230 0.1482787    
veget:fruit:time 13  726.1    55.9  0.7514 0.6992585    
veget:time:block 26 2194.9    84.4  1.1357 0.3762584    
fruit:time:block 26 2672.5   102.8  1.3828 0.2103929    
Residuals        25 1858.3    74.3                      
---
Signif. codes:  0 ?***? 0.001 ?**? 0.01 ?*? 0.05 ?.? 0.1 ? ? 1 
Warning message:
In aov(mice ~ veget * fruit * time + (time * block) + (time * block *  :
  Error() model is singular

            

      
      
      
    

        
structure(list(plot = structure(c(1L, 5L, 6L, 7L, 8L, 9L, 10L, 
11L, 12L, 2L, 3L, 4L, 1L, 5L, 6L, 7L, 8L, 9L, 10L, 11L, 12L, 
2L, 3L, 4L, 1L, 5L, 6L, 7L, 8L, 9L, 10L, 11L, 12L, 2L, 3L, 4L, 
1L, 5L, 6L, 7L, 8L, 9L, 10L, 11L, 12L, 2L, 3L, 4L, 1L, 5L, 6L, 
7L, 8L, 9L, 10L, 11L, 12L, 2L, 3L, 4L, 1L, 5L, 6L, 7L, 8L, 9L, 
10L, 11L, 12L, 2L, 3L, 4L, 1L, 5L, 6L, 7L, 8L, 9L, 10L, 11L, 
12L, 2L, 3L, 4L, 1L, 5L, 6L, 7L, 8L, 9L, 10L, 11L, 12L, 2L, 3L, 
4L, 1L, 5L, 6L, 7L, 8L, 9L, 10L, 11L, 12L, 2L, 3L, 4L, 1L, 5L, 
6L, 7L, 8L, 9L, 10L, 11L, 12L, 2L, 3L, 4L, 1L, 5L, 6L, 7L, 8L, 
9L, 10L, 11L, 12L, 2L, 3L, 4L, 1L, 5L, 6L, 7L, 8L, 9L, 10L, 11L, 
12L, 2L, 3L, 4L, 1L, 5L, 6L, 7L, 8L, 9L, 10L, 11L, 12L, 2L, 3L, 
4L, 1L, 5L, 6L, 7L, 8L, 9L, 10L, 11L, 12L, 2L, 3L, 4L), .Label = c("p1", 
"p10", "p11", "p12", "p2", "p3", "p4", "p5", "p6", "p7", "p8", 
"p9"), class = "factor"), veget = structure(c(2L, 1L, 2L, 1L, 
1L, 2L, 1L, 2L, 1L, 1L, 2L, 2L, 2L, 1L, 2L, 1L, 1L, 2L, 1L, 2L, 
1L, 1L, 2L, 2L, 2L, 1L, 2L, 1L, 1L, 2L, 1L, 2L, 1L, 1L, 2L, 2L, 
2L, 1L, 2L, 1L, 1L, 2L, 1L, 2L, 1L, 1L, 2L, 2L, 2L, 1L, 2L, 1L, 
1L, 2L, 1L, 2L, 1L, 1L, 2L, 2L, 2L, 1L, 2L, 1L, 1L, 2L, 1L, 2L, 
1L, 1L, 2L, 2L, 2L, 1L, 2L, 1L, 1L, 2L, 1L, 2L, 1L, 1L, 2L, 2L, 
2L, 1L, 2L, 1L, 1L, 2L, 1L, 2L, 1L, 1L, 2L, 2L, 2L, 1L, 2L, 1L, 
1L, 2L, 1L, 2L, 1L, 1L, 2L, 2L, 2L, 1L, 2L, 1L, 1L, 2L, 1L, 2L, 
1L, 1L, 2L, 2L, 2L, 1L, 2L, 1L, 1L, 2L, 1L, 2L, 1L, 1L, 2L, 2L, 
2L, 1L, 2L, 1L, 1L, 2L, 1L, 2L, 1L, 1L, 2L, 2L, 2L, 1L, 2L, 1L, 
1L, 2L, 1L, 2L, 1L, 1L, 2L, 2L, 2L, 1L, 2L, 1L, 1L, 2L, 1L, 2L, 
1L, 1L, 2L, 2L), .Label = c("Vintact", "Vremoved"), class = "factor"), 
    fruit = structure(c(1L, 1L, 2L, 2L, 2L, 1L, 1L, 2L, 1L, 2L, 
    1L, 2L, 1L, 1L, 2L, 2L, 2L, 1L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 
    1L, 2L, 2L, 2L, 1L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 1L, 2L, 2L, 
    2L, 1L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 1L, 2L, 2L, 2L, 1L, 1L, 
    2L, 1L, 2L, 1L, 2L, 1L, 1L, 2L, 2L, 2L, 1L, 1L, 2L, 1L, 2L, 
    1L, 2L, 1L, 1L, 2L, 2L, 2L, 1L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 
    1L, 2L, 2L, 2L, 1L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 1L, 2L, 2L, 
    2L, 1L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 1L, 2L, 2L, 2L, 1L, 1L, 
    2L, 1L, 2L, 1L, 2L, 1L, 1L, 2L, 2L, 2L, 1L, 1L, 2L, 1L, 2L, 
    1L, 2L, 1L, 1L, 2L, 2L, 2L, 1L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 
    1L, 2L, 2L, 2L, 1L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 1L, 2L, 2L, 
    2L, 1L, 1L, 2L, 1L, 2L, 1L, 2L), .Label = c("Fintact", "Fremoved"
    ), class = "factor"), time = structure(c(1L, 1L, 1L, 1L, 
    1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
    2L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 
    3L, 3L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 5L, 
    5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 6L, 6L, 6L, 6L, 
    6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 7L, 7L, 7L, 7L, 7L, 7L, 7L, 
    7L, 7L, 7L, 7L, 7L, 8L, 8L, 8L, 8L, 8L, 8L, 8L, 8L, 8L, 8L, 
    8L, 8L, 9L, 9L, 9L, 9L, 9L, 9L, 9L, 9L, 9L, 9L, 9L, 9L, 10L, 
    10L, 10L, 10L, 10L, 10L, 10L, 10L, 10L, 10L, 10L, 10L, 11L, 
    11L, 11L, 11L, 11L, 11L, 11L, 11L, 11L, 11L, 11L, 11L, 12L, 
    12L, 12L, 12L, 12L, 12L, 12L, 12L, 12L, 12L, 12L, 12L, 13L, 
    13L, 13L, 13L, 13L, 13L, 13L, 13L, 13L, 13L, 13L, 13L, 14L, 
    14L, 14L, 14L, 14L, 14L, 14L, 14L, 14L, 14L, 14L, 14L), .Label = c("1", 
    "2", "3", "4", "5", "6", "7", "8", "9", "10", "11", "12", 
    "13", "14"), class = "factor"), block = structure(c(1L, 1L, 
    1L, 1L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 1L, 1L, 1L, 1L, 2L, 
    2L, 2L, 2L, 3L, 3L, 3L, 3L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 
    3L, 3L, 3L, 3L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 
    3L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 1L, 1L, 
    1L, 1L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 1L, 1L, 1L, 1L, 2L, 
    2L, 2L, 2L, 3L, 3L, 3L, 3L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 
    3L, 3L, 3L, 3L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 
    3L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 1L, 1L, 
    1L, 1L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 1L, 1L, 1L, 1L, 2L, 
    2L, 2L, 2L, 3L, 3L, 3L, 3L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 
    3L, 3L, 3L, 3L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 
    3L), .Label = c("b1", "b2", "b3"), class = "factor"), rawcoon = c(0, 
    0.063, 0.125, 0.063, 0.25, 0.125, 0, 0, 0, 0.063, 0, 0, 0.125, 
    0, 0.375, 0.375, 0.188, 0.25, 0.063, 0.063, 0, 0.188, 0.125, 
    0.188, 0.063, 0, 0, 0, 0, 0, 0, 0, 0, 0.063, 0.063, 0, 0.063, 
    0, 0.188, 0.125, 0, 0.125, 0.063, 0, 0.125, 0, 0.188, 0, 
    0, 0, 0.063, 0.063, 0, 0, 0.063, 0, 0, 0.125, 0, 0, 0, 0, 
    0, 0.125, 0.125, 0.125, 0, 0, 0, 0, 0.063, 0, 0, 0.125, 0.063, 
    0.063, 0.188, 0, 0.063, 0.125, 0, 0, 0, 0.125, 0, 0, 0, 0.063, 
    0.063, 0, 0.188, 0.125, 0, 0.063, 0.063, 0, 0, 0.125, 0, 
    0.25, 0.188, 0, 0.063, 0.563, 0.188, 0.063, 0.188, 0, 0, 
    0.063, 0, 0, 0.063, 0, 0.063, 0.063, 0, 0.125, 0.063, 0, 
    0, 0.125, 0.063, 0.125, 0.125, 0.063, 0.063, 0.125, 0.125, 
    0.063, 0, 0.125, 0.125, 0.438, 0.125, 0.25, 0.188, 0.188, 
    0.375, 0.125, 0.25, 0.313, 0.063, 0.063, 0, 0.063, 0, 0, 
    0.063, 0.063, 0, 0, 0.125, 0.063, 0.063, 0.063, 0.063, 0.063, 
    0, 0.063, 0.125, 0, 0.313, 0, 0, 0.063, 0, 0.125), coon = c(8.605276544, 
    16.64929378, 22.12023855, 16.64929378, 30.73522506, 22.12023855, 
    8.605276544, 8.605276544, 8.605276544, 16.64929378, 8.605276544, 
    8.605276544, 22.12023855, 8.605276544, 38.09196076, 38.09196076, 
    26.67179738, 30.73522506, 16.64929378, 16.64929378, 8.605276544, 
    26.67179738, 22.12023855, 26.67179738, 16.64929378, 8.605276544, 
    8.605276544, 8.605276544, 8.605276544, 8.605276544, 8.605276544, 
    8.605276544, 8.605276544, 16.64929378, 16.64929378, 8.605276544, 
    16.64929378, 8.605276544, 26.67179738, 22.12023855, 8.605276544, 
    22.12023855, 16.64929378, 8.605276544, 22.12023855, 8.605276544, 
    26.67179738, 8.605276544, 8.605276544, 8.605276544, 16.64929378, 
    16.64929378, 8.605276544, 8.605276544, 16.64929378, 8.605276544, 
    8.605276544, 22.12023855, 8.605276544, 8.605276544, 8.605276544, 
    8.605276544, 8.605276544, 22.12023855, 22.12023855, 22.12023855, 
    8.605276544, 8.605276544, 8.605276544, 8.605276544, 16.64929378, 
    8.605276544, 8.605276544, 22.12023855, 16.64929378, 16.64929378, 
    26.67179738, 8.605276544, 16.64929378, 22.12023855, 8.605276544, 
    8.605276544, 8.605276544, 22.12023855, 8.605276544, 8.605276544, 
    8.605276544, 16.64929378, 16.64929378, 8.605276544, 26.67179738, 
    22.12023855, 8.605276544, 16.64929378, 16.64929378, 8.605276544, 
    8.605276544, 22.12023855, 8.605276544, 30.73522506, 26.67179738, 
    8.605276544, 16.64929378, 48.42882421, 26.67179738, 16.64929378, 
    26.67179738, 8.605276544, 8.605276544, 16.64929378, 8.605276544, 
    8.605276544, 16.64929378, 8.605276544, 16.64929378, 16.64929378, 
    8.605276544, 22.12023855, 16.64929378, 8.605276544, 8.605276544, 
    22.12023855, 16.64929378, 22.12023855, 22.12023855, 16.64929378, 
    16.64929378, 22.12023855, 22.12023855, 16.64929378, 8.605276544, 
    22.12023855, 22.12023855, 41.57117579, 22.12023855, 30.73522506, 
    26.67179738, 26.67179738, 38.09196076, 22.12023855, 30.73522506, 
    34.50487875, 16.64929378, 16.64929378, 8.605276544, 16.64929378, 
    8.605276544, 8.605276544, 16.64929378, 16.64929378, 8.605276544, 
    8.605276544, 22.12023855, 16.64929378, 16.64929378, 16.64929378, 
    16.64929378, 16.64929378, 8.605276544, 16.64929378, 22.12023855, 
    8.605276544, 34.50487875, 8.605276544, 8.605276544, 16.64929378, 
    8.605276544, 22.12023855), mice = c(26.67179738, 45, 26.67179738, 
    41.57117579, 38.09196076, 26.67179738, 8.605276544, 8.605276544, 
    16.64929378, 8.605276544, 30.73522506, 16.64929378, 41.57117579, 
    48.42882421, 16.64929378, 38.09196076, 30.73522506, 16.64929378, 
    16.64929378, 8.605276544, 16.64929378, 26.67179738, 16.64929378, 
    22.12023855, 51.90803924, 38.09196076, 34.50487875, 48.42882421, 
    8.605276544, 16.64929378, 26.67179738, 8.605276544, 8.605276544, 
    26.67179738, 22.12023855, 34.50487875, 41.57117579, 30.73522506, 
    38.09196076, 16.64929378, 34.50487875, 30.73522506, 8.605276544, 
    8.605276544, 51.90803924, 22.12023855, 30.73522506, 8.605276544, 
    34.50487875, 41.57117579, 26.67179738, 38.09196076, 22.12023855, 
    16.64929378, 16.64929378, 22.12023855, 30.73522506, 26.67179738, 
    22.12023855, 34.50487875, 26.67179738, 34.50487875, 22.12023855, 
    55.49512125, 16.64929378, 16.64929378, 34.50487875, 16.64929378, 
    41.57117579, 26.67179738, 16.64929378, 22.12023855, 59.26477494, 
    45, 41.57117579, 51.90803924, 41.57117579, 16.64929378, 34.50487875, 
    41.57117579, 34.50487875, 30.73522506, 22.12023855, 26.67179738, 
    30.73522506, 51.90803924, 51.90803924, 51.90803924, 41.57117579, 
    22.12023855, 38.09196076, 22.12023855, 26.67179738, 34.50487875, 
    26.67179738, 22.12023855, 41.57117579, 63.32820262, 26.67179738, 
    45, 34.50487875, 16.64929378, 30.73522506, 34.50487875, 30.73522506, 
    41.57117579, 30.73522506, 38.09196076, 38.09196076, 38.09196076, 
    22.12023855, 26.67179738, 34.50487875, 26.67179738, 38.09196076, 
    34.50487875, 34.50487875, 38.09196076, 16.64929378, 22.12023855, 
    38.09196076, 30.73522506, 51.90803924, 38.09196076, 16.64929378, 
    30.73522506, 16.64929378, 30.73522506, 59.26477494, 38.09196076, 
    38.09196076, 45, 34.50487875, 38.09196076, 26.67179738, 34.50487875, 
    22.12023855, 26.67179738, 45, 22.12023855, 34.50487875, 26.67179738, 
    8.605276544, 30.73522506, 22.12023855, 30.73522506, 16.64929378, 
    22.12023855, 26.67179738, 8.605276544, 30.73522506, 16.64929378, 
    26.67179738, 16.64929378, 8.605276544, 26.67179738, 51.90803924, 
    63.32820262, 30.73522506, 59.26477494, 48.42882421, 45, 59.26477494, 
    38.09196076, 48.42882421, 51.90803924, 51.90803924, 51.90803924
    ), rawmice = c(0.188, 0.5, 0.188, 0.438, 0.375, 0.188, 0, 
    0, 0.063, 0, 0.25, 0.063, 0.438, 0.563, 0.063, 0.375, 0.25, 
    0.063, 0.063, 0, 0.063, 0.188, 0.063, 0.125, 0.625, 0.375, 
    0.313, 0.563, 0, 0.063, 0.188, 0, 0, 0.188, 0.125, 0.313, 
    0.438, 0.25, 0.375, 0.063, 0.313, 0.25, 0, 0, 0.625, 0.125, 
    0.25, 0, 0.313, 0.438, 0.188, 0.375, 0.125, 0.063, 0.063, 
    0.125, 0.25, 0.188, 0.125, 0.313, 0.188, 0.313, 0.125, 0.688, 
    0.063, 0.063, 0.313, 0.063, 0.438, 0.188, 0.063, 0.125, 0.75, 
    0.5, 0.438, 0.625, 0.438, 0.063, 0.313, 0.438, 0.313, 0.25, 
    0.125, 0.188, 0.25, 0.625, 0.625, 0.625, 0.438, 0.125, 0.375, 
    0.125, 0.188, 0.313, 0.188, 0.125, 0.438, 0.813, 0.188, 0.5, 
    0.313, 0.063, 0.25, 0.313, 0.25, 0.438, 0.25, 0.375, 0.375, 
    0.375, 0.125, 0.188, 0.313, 0.188, 0.375, 0.313, 0.313, 0.375, 
    0.063, 0.125, 0.375, 0.25, 0.625, 0.375, 0.063, 0.25, 0.063, 
    0.25, 0.75, 0.375, 0.375, 0.5, 0.313, 0.375, 0.188, 0.313, 
    0.125, 0.188, 0.5, 0.125, 0.313, 0.188, 0, 0.25, 0.125, 0.25, 
    0.063, 0.125, 0.188, 0, 0.25, 0.063, 0.188, 0.063, 0, 0.188, 
    0.625, 0.813, 0.25, 0.75, 0.563, 0.5, 0.75, 0.375, 0.563, 
    0.625, 0.625, 0.625)), .Names = c("plot", "veget", "fruit", 
"time", "block", "rawcoon", "coon", "mice", "rawmice"), row.names = c(NA, 
-168L), class = "data.frame")

        
Richard M. Heiberger wrote:

            

      
      
      
    

        

      
      
    

        

  
    

      
      
    

  


  


      

    

    

      
      

        


  

        

      Richard M. Heiberger
    

    

      
      #
    

  


  

      
Your email program truncated the model.  It would not run,
hence was not reproducible.
The last few characters of the first line are

    +Error(block/plot,

which is syntactically impossible because there is no
closing parenthesis before the comma.
Try executing your email and see the difficulty.

Please resend the model statement from your source file,
not from the R output window,
with no more than 60 characters per line.

I am asking for the following from the trailer of all R-help email
"and provide commented, minimal, self-contained, reproducible code."

I have the dataset now

 > dim(track)
[1] 168   9
 > sapply(track, class)
      plot     veget     fruit      time     block   rawcoon      coon 
     mice
  "factor"  "factor"  "factor"  "factor"  "factor" "numeric" "numeric" 
"numeric"
   rawmice
"numeric"
 > sapply(track, function(x) length(levels(x)))
    plot   veget   fruit    time   block rawcoon    coon    mice rawmice
      12       2       2      14       3       0       0       0       0
 >


Rich

  
  


  


      

    

    

      
      

        


  

        

      hpdutra
    

    

      
      #
    

  


  

      
I can actually run the code from my post.
I used the nabble for my list server
http://www.nabble.com/ANCOVA-with-defined-error-terms-td25055311.html#a25100032
I don't know which server you use, but that one is not truncated, I can copy
the code just fine and run it.
Anyway, here it is again

summary(mymodel<-(aov(mice~veget*fruit*time+(time*block)+
(time*block*veget)+(time*block*fruit)+coon+Error(block/plot, data =
track))))

I have a question. After I use your suggestion of using 
sapply(track, function(x) length(levels(x))
I get a completely wrong structure for the DF. I get just 1 DF for time,
instead of 13. 
And I cannot run your model 
ricsmodel<-aov(mice ~ coon + block*veget*fruit*time - block:veget:fruit:time
+ Error(block/plot), data = track))

Any insights are helpful. Thank you for your patience.
Humberto

        
Richard M. Heiberger wrote:

            

      
      
      
    

        

      
      
    

        

  
    

      
      
    

  


  


      

    

    

      
      

        


  

        

      Richard M. Heiberger
    

    

      
      #
    

  


  

      
Your model formula cannot be correct.

The phrase
Error(block/plot, data = track)

is wrong.

It has to be something like this
Error(block/plot), data = track

The Error function requires a well-defined formula.
The "," character cannot be inside the Error function.



You misunderstood my use of sapply.
I am illustrating the number of levels for each of the factors.
That is information about the data.frame.
It is not part of a model specification.



I need reproducible code to do anything further.
I need exactly three complete statements of the form

    track <- structure(....)    ## from your earlier email
    mymodel <- aov(   formula + Error(errorformula), data=track)
    summary(mymodel)

Do not attach() anything.  I am wondering now if that use of the
attach() function in early email was part of your difficulty.

When you have those three statements ready, open a fresh copy
of R with the MSDOS command statement

    c:\progra~1\R\R-2.9.1\bin\Rgui --vanilla

Then paste the three statements in.  If that works, send the three
statements.

  
  


  


      

    

    

      
      

        


  

        

      hpdutra
    

    

      
      #
    

  


  

      
What does the command 
c:\progra~1\R\R-2.9.1\bin\Rgui --vanilla
do?
At first I thought that I could run it from R.  But then it did not work.
I'm using a mac, and I don't know how access MSDOS in a mac. I actually
installed R on the pc and tried to run the command from MSDOS and it also
didn't work. 
could you please tell me what were you trying to do? I'm guessing that you
were trying to get an R console that is absolutely clean, I think that what
you meant by fresh copy of R.
Anyway, I ran the same model on the pc,  this is a fresh copy, there is
nothing attached or stored. 

mymodel<-(aov(mice~veget*fruit*time+(time*block)+(time*block*veget)+(time*block*fruit)+coon+Error(block/plot),
data = track))
Warning message:
In aov(mice ~ veget * fruit * time + (time * block) + (time * block *  :
  Error() model is singular

summary(mymodel)

Error: block
      Df Sum Sq Mean Sq
block  2 4581.5  2290.7

Error: block:plot
            Df  Sum Sq Mean Sq  F value  Pr(>F)  
veget        1 1602.91 1602.91 379.1653 0.03267 *
fruit        1   80.07   80.07  18.9411 0.14378  
coon         1  193.53  193.53  45.7791 0.09341 .
veget:fruit  1   40.26   40.26   9.5245 0.19948  
veget:block  2  355.53  177.76  42.0495 0.10840  
fruit:block  2   15.39    7.69   1.8200 0.46423  
Residuals    1    4.23    4.23                   
---
Signif. codes:  0 ?***? 0.001 ?**? 0.01 ?*? 0.05 ?.? 0.1 ? ? 1 

Error: Within
                 Df Sum Sq Mean Sq F value    Pr(>F)    
time             13 8451.5   650.1  8.7462 2.384e-06 ***
coon              1 1070.7  1070.7 14.4047 0.0008364 ***
veget:time       13 1214.4    93.4  1.2568 0.3006216    
fruit:time       13  474.7    36.5  0.4912 0.9090773    
time:block       26 2943.4   113.2  1.5230 0.1482787    
veget:fruit:time 13  726.1    55.9  0.7514 0.6992585    
veget:time:block 26 2194.9    84.4  1.1357 0.3762584    
fruit:time:block 26 2672.5   102.8  1.3828 0.2103929    
Residuals        25 1858.3    74.3                      
---
Signif. codes:  0 ?***? 0.001 ?**? 0.01 ?*? 0.05 ?.? 0.1 ? ? 1 

or using yours
ricsmodel<-aov(mice ~ coon + block*veget*fruit*time -
block:veget:fruit:time+ Error(block/plot), data = track)
Warning message:
In aov(mice ~ coon + block * veget * fruit * time - block:veget:fruit:time + 
:
  Error() model is singular

            

      
      
      
    

        
Error: block
      Df Sum Sq Mean Sq
coon   1 1104.4  1104.4
block  1 3477.1  3477.1

Error: block:plot
                  Df  Sum Sq Mean Sq
coon               1  240.13  240.13
veget              1 1629.40 1629.40
fruit              1    6.98    6.98
block:veget        2  247.59  123.79
block:fruit        2   36.97   18.49
veget:fruit        1  126.62  126.62
block:veget:fruit  1    4.23    4.23

Error: Within
                 Df Sum Sq Mean Sq F value    Pr(>F)    
coon              1  904.4   904.4 12.1673  0.001819 ** 
time             13 8617.8   662.9  8.9183 1.986e-06 ***
block:time       26 2945.9   113.3  1.5243  0.147803    
veget:time       13 1209.2    93.0  1.2513  0.303764    
fruit:time       13  477.5    36.7  0.4941  0.907294    
block:veget:time 26 2201.5    84.7  1.1392  0.373409    
block:fruit:time 26 2706.7   104.1  1.4006  0.201399    
veget:fruit:time 13  685.2    52.7  0.7091  0.737650    
Residuals        25 1858.3    74.3                      
---
Signif. codes:  0 ?***? 0.001 ?**? 0.01 ?*? 0.05 ?.? 0.1 ? ? 1

            

      
      
      
    

        
Thank you
Humberto

        
Richard M. Heiberger wrote:

            

      
      
      
    

        

      
      
    

        

  
    

      
      
    

  


  


      

    

    

      
      

        


  

        

      Richard M. Heiberger
    

    

      
      #
    

  


  

      
Ok, now we can talk.

1. covariates: coon

Your model specification put coon sequentially last, effectively testing
the hypothesis that the slope associated with coon, after
adjusting for all the factors, is zero.

My model place coon sequentially first, effectively testing the
hypotheses that the factors are significant after adjusting to a common
level of coon.

Which of those is appropriate for your study is something for
you, as subject-matter investigator, to determine.  It is not a
computational question.

That is the only difference between your specification and mine.
As you can see, the degrees of freedom and Within stratum
residuals sums of squares are the same for both specifications.

The block:plot Residual in your specification is the same as the
block:veget:fruit SumSq in mine.


2.  Error() model is singular

The error model block/plot expands to block + block:plot
This model specifies 2 + 33 degrees of freedom.

There are only 9 df in the block:plot stratum.

Probably block does not belong in the treatment model at all,
it should appear in only the Error model.


3. The next obvious model to attempt is

nextmodel <- aov(mice ~ coon + veget*fruit*time + Error(block/plot),
                  data = track)
summary(nextmodel)

This also shows a singular Error().  We look at the data and see that
plot is identical to the three-way veget:fruit:block interaction.


4. Therefore rename the plots to be 1,2,3,4 within each block.
The new factor name begins with an uppercase P.

track$Plot <- with(track, interaction(veget, fruit))
## then use the model
modelPlot <- aov(mice ~ coon + veget*fruit*time + Error(block/Plot),
                  data = track)
summary(modelPlot)


5. modelPlot is something that you can study.  Now you can interpret it,
by looking for example at plots and tables of means, and decide what to
do next.


6.  --vanilla
Please see ?Startup
for a discussion of the meaning of the startup arguments.
I gave the Windows call using that argument.  The starting call on Mac
needs to follow Mac rules.


7. Error() notation.  Numerators and Denominators for F tests are
specified indirectly by the relation of the linear subspaces in the
treatment model and the Error model.

Rich

  
  


  


      

    

    

      
      

        


  

        

      hpdutra
    

    

      
      #
    

  


  

      
When we ran a regular ANOVA we showed that both mice and raccoons respond
positively to cover. So we decided to use raccoons as a co-variate because
those are mouse predators and their presence per se could explain part of
the variation on mice activity. That is the reason why I'm running this
ANCOVA. Thus your model that places coon sequentially first makes more sense
for this question. 

 "Error() model is singular" Some authors claim that one should test for the
interactions of block versus treatment,  to identify any spatial correlation
with the data. That is the reason why I tried to incorporate the
interactions of factors with block in my model. But like you said I don't
have enough degrees of freedom for it. Is there an R function that
calculates the degrees of freedom for complex designs?

I did not understand what you meant by

"Error() notation.  Numerators and Denominators for F tests are
specified indirectly by the relation of the linear subspaces in the
treatment model and the Error model." 

If this is last question is outside the scope of R-help, please just give me
a reference that covers it and I'll do my homework.

Thank you very much Richard, this was wonderful help.

  
    

      
      
    

  


  


      

    

    
2 days later

    

      
      

        


  

        

      hpdutra
    

    

      
      #
    

  


  

      
Could someone or Richard explain to me what he meant by 

"This also shows a singular Error().  We look at the data and see that
plot is identical to the three-way veget:fruit:block interaction."
It seems to me that I just needed to recoded the plots, in order to get rid
of the Error message. If that is true, then why I can't go back to the
original model proposed by Richard
track.aov <- aov(mice ~ coon
+ block * veget * fruit * time - block:veget:fruit:time
+ Error(block/Plot), data = track) 

and calculate F values myself? 

Thank you very much for the feedback.