Skip to content

explaining lme variance component results

3 messages · Mike Dunbar, John Maindonald, David Afshartous

Original
Mike Dunbar
# 
Thanks to all, a couple more comments following up on Kevin's comments below, and also ones sent to me directly.

I have plotted the data in many different ways, having spent several years (yes!) trying to work out a suitable analysis for these data. The aim of this particular analysis is to try to keep things as simple as possible, I'm aware in particular that there are differences between the behaviour of the factors across the months (so one option is a month by month analysis - which I have done but it was vetoed by co-workers for this paper so long as there's a simple interpretation as well), and also that time is generally the most important factor overall (this is already documented by others - the data are of drifting macroinvertebrates in rivers in case anyones interested).

The structure of the nesting is designed to mirror our expected view of the correlations in the data based on spatial/temporal proximity, a bit as Kevin describes below: so four times were measured across a day and the experiment repeated across four months, and for each of the 16 occasions, we have five transects, within those four poles each, and not described previously, 1-3 measures at different heights on the poles. 

Regarding the zero values: yes the normality is an assumption, I hope to do better once this initial analysis is over. What I hoped to show is despite this, and despite the assumptions of the variance components analysis, there is evidence of an effect of TRANSECT and / or POLE, once MONTH and TIME are accounted for. 

What is very pertinent (thanks John) is the fact that in the data as described, there is no replication within the lowest stratum, POLE. There was one seeming replicate, but that must be an error. This may well be the source of the problem that the POLE variance component was large but not significant.

I had thought that despite the lack of replication within POLE that it would still be possible to estimate a variance component for POLE separately from the residual. The very wooly reasoning being that the POLE component represents consistency in drift density between POLEs across TRANSECT, TIME and MONTH, and residual represents lack of consistency. 

If my reasoning above is flawed, I really don't want to ditch the POLE component, as its fairly central to the analysis, and I could bring in HEIGHT to give replication within POLE (previous data is for one height only). I'd prefer to do this as a fixed effect and I've posted below some example data/code: can anyone comment if this is valid?

Regarding the issue of magnitude of variance component/random effect vs significance, I wonder if there is more too it than that, certainly in this case we know that TIME is more important than MONTH, despite being nested, but more critically, I can show some data where the magnitude of the component doesn't seem to relate to its significance. I'll post this in a separate mail to avoid confusion, once again any comments are welcome. This gives me a real headache explaining my results to my co-workers, let alone reviewers. I ought to add that there could easily still be mistakes where, as one regarding a non-replicate 
has already been identified.

All the best again - hope this is interesting to others struggling with similar issues??

Mike




varcor.2h.insects.hf <- lme(log(insectdens+1) ~ HEIGHT, random=~1|MONTH/TIME/TRANSECT/POLE, data=temp2)
# introduce HEIGHT as a fixed effect, there are two heights per pole for some poles: hence unbalanced
VarCorr(varcor.2h.insects.hf)
# variances: MONTH - 0.639, TIME: 1.248, TRANSECT: 0.013, POLE: 0.160, Residual: 1.016

varcor.2h.insects.nospat.hf <- lme(log(insectdens+1) ~ HEIGHT, random=~1|MONTH/TIME, data=temp2)

anova(varcor.2h.insects.hf,varcor.2h.insects.nospat.hf)
# two spatial factors together marginally signficant: p=0.06, but test likely conservative 
# simulation approach for null distribution (Faraway) probably too difficult at this depth of nesting
intervals(varcor.2h.insects.hf)
# again some evidence for significance of TRANSECT, but POLE lower bound close to 0.

# delete transect term and just compare models with and without pole term
varcor.2h.insects.pole.hf <- lme(log(insectdens+1) ~ HEIGHT, random=~1|MONTH/TIME/POLE, data=temp2)
# test pole factor on its own. This is possible as pole is coded as a combination of transect and pole within transect
anova(varcor.2h.insects.pole.hf,varcor.2h.insects.nospat.hf)
# p=0.019. This would be great if analysis is valid



# read in data: this time with one or two heights per pole

temp2 <-
structure(list(MONTH = structure(c(1L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 
3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 
3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 
3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 
3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 
3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 
3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 
3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 
3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 
4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 
4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 
4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 
4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 
4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 
4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 
4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 
4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L
), .Label = c("4", "5", "6", "7"), class = "factor"), TRANSECT = structure(c(1L, 
1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 3L, 3L, 
3L, 3L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 5L, 5L, 5L, 5L, 5L, 5L, 
1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 3L, 3L, 
3L, 3L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 5L, 5L, 5L, 5L, 5L, 5L, 
1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 3L, 
3L, 3L, 3L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 5L, 5L, 5L, 5L, 5L, 
5L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 
3L, 3L, 3L, 3L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 5L, 5L, 5L, 5L, 
5L, 5L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 
3L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 5L, 5L, 5L, 
5L, 5L, 5L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 
3L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 5L, 5L, 5L, 
5L, 5L, 5L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 
3L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 5L, 5L, 5L, 
5L, 5L, 5L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 3L, 3L, 
3L, 3L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 4L, 4L, 4L, 5L, 5L, 5L, 5L, 
5L, 5L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 
3L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 5L, 5L, 5L, 
5L, 5L, 5L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 3L, 3L, 
3L, 3L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 5L, 5L, 
5L, 5L, 5L, 5L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 3L, 
3L, 3L, 3L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 5L, 
5L, 5L, 5L, 5L, 5L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 
3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 
5L, 5L, 5L, 5L, 5L, 5L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 
3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 
5L, 5L, 5L, 5L, 5L, 5L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 
3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 
5L, 5L, 5L, 5L, 5L, 5L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 
3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 
5L, 5L, 5L, 5L, 5L, 5L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 
3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 
5L, 5L, 5L, 5L, 5L, 5L), .Label = c("1", "2", "3", "4", "5"), class = "factor"), 
    POLE = structure(c(1L, 1L, 2L, 2L, 3L, 4L, 4L, 5L, 5L, 6L, 
    6L, 7L, 7L, 8L, 8L, 9L, 9L, 10L, 10L, 11L, 11L, 12L, 12L, 
    13L, 13L, 14L, 14L, 15L, 16L, 17L, 17L, 18L, 18L, 1L, 2L, 
    2L, 3L, 4L, 4L, 5L, 5L, 5L, 6L, 7L, 7L, 8L, 8L, 9L, 9L, 10L, 
    10L, 11L, 11L, 12L, 12L, 13L, 13L, 14L, 14L, 15L, 16L, 17L, 
    17L, 18L, 18L, 1L, 1L, 2L, 2L, 3L, 4L, 4L, 5L, 5L, 6L, 6L, 
    7L, 7L, 8L, 8L, 9L, 9L, 10L, 10L, 11L, 11L, 12L, 12L, 13L, 
    13L, 14L, 14L, 15L, 16L, 17L, 17L, 18L, 18L, 1L, 1L, 2L, 
    2L, 3L, 4L, 4L, 5L, 5L, 6L, 6L, 7L, 7L, 8L, 8L, 9L, 9L, 10L, 
    10L, 11L, 11L, 12L, 12L, 13L, 13L, 14L, 14L, 15L, 16L, 17L, 
    17L, 18L, 18L, 1L, 1L, 2L, 2L, 3L, 4L, 4L, 5L, 5L, 6L, 6L, 
    7L, 7L, 8L, 8L, 9L, 9L, 10L, 10L, 11L, 11L, 12L, 12L, 13L, 
    13L, 14L, 14L, 15L, 16L, 17L, 17L, 18L, 18L, 1L, 2L, 2L, 
    3L, 4L, 4L, 5L, 5L, 6L, 6L, 7L, 7L, 8L, 8L, 9L, 9L, 10L, 
    10L, 11L, 11L, 12L, 12L, 13L, 13L, 14L, 14L, 15L, 16L, 17L, 
    17L, 18L, 18L, 1L, 2L, 2L, 3L, 4L, 4L, 5L, 5L, 6L, 6L, 7L, 
    7L, 8L, 8L, 9L, 9L, 10L, 10L, 11L, 11L, 12L, 12L, 13L, 13L, 
    14L, 14L, 15L, 16L, 17L, 17L, 18L, 18L, 1L, 1L, 2L, 2L, 3L, 
    4L, 4L, 5L, 5L, 6L, 6L, 7L, 7L, 8L, 8L, 9L, 9L, 10L, 10L, 
    11L, 11L, 13L, 13L, 14L, 14L, 15L, 16L, 17L, 17L, 18L, 18L, 
    1L, 1L, 2L, 2L, 3L, 4L, 4L, 5L, 5L, 6L, 6L, 7L, 7L, 8L, 8L, 
    9L, 9L, 10L, 10L, 11L, 11L, 12L, 12L, 13L, 13L, 14L, 14L, 
    15L, 16L, 17L, 17L, 18L, 18L, 1L, 1L, 2L, 2L, 3L, 4L, 4L, 
    5L, 5L, 6L, 6L, 7L, 7L, 8L, 8L, 9L, 9L, 10L, 10L, 11L, 11L, 
    12L, 12L, 13L, 13L, 14L, 14L, 15L, 16L, 17L, 17L, 18L, 18L, 
    1L, 1L, 2L, 2L, 3L, 4L, 4L, 5L, 5L, 6L, 6L, 7L, 7L, 8L, 8L, 
    9L, 9L, 10L, 10L, 11L, 11L, 12L, 12L, 13L, 13L, 14L, 14L, 
    15L, 16L, 17L, 17L, 18L, 18L, 1L, 1L, 2L, 2L, 3L, 4L, 4L, 
    5L, 5L, 6L, 6L, 7L, 7L, 8L, 8L, 9L, 9L, 10L, 10L, 11L, 11L, 
    12L, 12L, 13L, 13L, 14L, 14L, 15L, 16L, 17L, 17L, 18L, 18L, 
    1L, 2L, 2L, 3L, 4L, 4L, 5L, 5L, 6L, 6L, 7L, 7L, 8L, 8L, 9L, 
    9L, 10L, 10L, 11L, 11L, 12L, 12L, 13L, 13L, 14L, 14L, 15L, 
    16L, 17L, 17L, 18L, 18L, 1L, 2L, 2L, 3L, 4L, 4L, 5L, 5L, 
    6L, 6L, 7L, 7L, 8L, 8L, 9L, 9L, 10L, 10L, 11L, 11L, 12L, 
    12L, 13L, 13L, 14L, 14L, 15L, 16L, 17L, 17L, 18L, 18L, 1L, 
    2L, 2L, 3L, 4L, 4L, 5L, 5L, 6L, 6L, 7L, 7L, 8L, 8L, 9L, 9L, 
    10L, 10L, 11L, 11L, 12L, 12L, 13L, 13L, 14L, 14L, 15L, 16L, 
    17L, 17L, 18L, 18L, 1L, 2L, 2L, 3L, 4L, 4L, 5L, 5L, 6L, 6L, 
    7L, 7L, 8L, 8L, 9L, 9L, 10L, 10L, 11L, 11L, 12L, 12L, 13L, 
    13L, 14L, 14L, 15L, 16L, 17L, 17L, 18L, 18L), .Label = c("11", 
    "12", "13", "14", "23", "24", "31", "32", "33", "34", "41", 
    "42", "43", "44", "51", "52", "53", "54"), class = "factor"), 
    TIME = structure(c(1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
    1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
    1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
    2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
    2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 3L, 
    3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 
    3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 4L, 4L, 
    4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 
    4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 
    4L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
    1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
    1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
    2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
    2L, 2L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 
    3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 
    3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 
    4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 
    4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 1L, 1L, 1L, 1L, 1L, 1L, 
    1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
    1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 
    2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
    2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
    3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 
    3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 
    3L, 3L, 3L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 
    4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 
    4L, 4L, 4L, 4L, 4L, 4L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
    1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 
    1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
    2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 
    2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 3L, 
    3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 
    3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 4L, 4L, 4L, 
    4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 
    4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L), .Label = c("1", 
    "2", "3", "4"), class = "factor"), HEIGHT = structure(c(1L, 
    2L, 1L, 2L, 1L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 
    2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 1L, 1L, 2L, 
    1L, 2L, 2L, 1L, 2L, 1L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 
    2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 1L, 
    1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 1L, 2L, 1L, 2L, 1L, 2L, 
    1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 
    2L, 1L, 1L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 1L, 2L, 1L, 
    2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 
    1L, 2L, 1L, 2L, 1L, 1L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 
    1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 
    2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 1L, 1L, 2L, 1L, 2L, 1L, 1L, 
    2L, 1L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 
    2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 1L, 1L, 2L, 1L, 2L, 
    1L, 1L, 2L, 1L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 
    2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 1L, 1L, 2L, 
    1L, 2L, 1L, 2L, 1L, 2L, 1L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 
    1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 1L, 1L, 
    2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 
    2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 
    1L, 1L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 1L, 2L, 1L, 2L, 
    1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 
    2L, 1L, 2L, 1L, 1L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 1L, 
    2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 
    1L, 2L, 1L, 2L, 1L, 2L, 1L, 1L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 
    2L, 1L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 
    2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 1L, 1L, 2L, 1L, 2L, 
    1L, 1L, 2L, 1L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 
    2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 1L, 1L, 2L, 
    1L, 2L, 1L, 1L, 2L, 1L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 
    2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 1L, 
    1L, 2L, 1L, 2L, 1L, 1L, 2L, 1L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 
    2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 
    1L, 1L, 1L, 2L, 1L, 2L, 1L, 1L, 2L, 1L, 1L, 2L, 1L, 2L, 1L, 
    2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 1L, 2L, 
    1L, 2L, 1L, 1L, 1L, 2L, 1L, 2L), .Label = c("1", "2", "3"
    ), class = "factor"), insectdens = c(0, 0, 63.64, 11.99, 
    14.57, 22.5, 15.53, 0, 20.49, 107.6, 0, 87.16, 19.82, 22.24, 
    18.28, 51.92, 33.87, 42.1, 59.01, 0, 47.32, 15.78, 12.46, 
    43.02, 12.23, 9.98, 27.75, 7.47, 38.92, 11.78, 11.18, 0, 
    0, 120.6, 44.61, 24.02, 45.9, 26.78, 14.56, 80.2, 62.34, 
    37.4, 32.44, 17.58, 47.52, 8.94, 26.01, 54.7, 9.19, 141.89, 
    29.36, 10.39, 48.88, 14.6, 20.46, 158.34, 20.5, 9.52, 18.82, 
    14.36, 47.94, 12.26, 45.76, 31.44, 53.82, 104.37, 112, 74.4, 
    59.88, 73.38, 94.36, 73.78, 120.26, 305, 48.12, 129.45, 264.87, 
    53.88, 129.36, 87.9, 107.03, 57.33, 145.53, 90.48, 95.2, 
    110, 116.55, 110.44, 492, 50.7, 140.4, 68.16, 111.28, 104.8, 
    59.76, 75, 91.92, 68.4, 12.92, 19.94, 22.6, 17.38, 53.6, 
    102.6, 10.45, 151.92, 30.3, 0, 0, 0, 39.18, 34.96, 16.38, 
    21.38, 18.32, 60.4, 35.48, 16.9, 0, 24.96, 56.28, 263.76, 
    38.24, 37.12, 9.26, 30.76, 26.24, 25.88, 46.48, 7.2, 21.09, 
    48.87, 0, 28.1, 10.09, 44.28, 67.26, 0, 0, 29.72, 50.49, 
    63.92, 0, 0, 0, 18.28, 10.82, 7.5, 27.06, 21.48, 9.09, 21.94, 
    13.56, 10.4, 13.25, 46.6, 31.74, 8.57, 11.98, 12.08, 30.55, 
    12.46, 31.16, 27.27, 16.35, 78.15, 100.8, 13.54, 80.44, 69.35, 
    104.55, 83.6, 37.32, 0, 107.7, 91.55, 21.52, 50.76, 22.28, 
    17, 55.6, 52.85, 40.72, 15.76, 15.12, 41.08, 25.44, 10.79, 
    87.36, 19.58, 19.94, 78.32, 13.04, 39.54, 40.55, 74.08, 14.37, 
    34.68, 31.68, 69.4, 62.28, 13.13, 117.96, 41.02, 18.27, 72.66, 
    34.74, 30.2, 69.86, 17.4, 100.89, 16.72, 95.7, 43.92, 0, 
    27.6, 129.6, 73.64, 147.4, 107.82, 92.16, 46.9, 76.1, 52.78, 
    52.32, 60.57, 46.7, 48.65, 49.41, 0, 54.8, 30.18, 59.2, 0, 
    12.52, 0, 0, 15.89, 90.39, 35.42, 26.64, 8.54, 17.46, 52.98, 
    7.88, 48.81, 12.68, 49.85, 32.67, 64.6, 41.2, 20.2, 8.47, 
    80.29, 38.52, 17.28, 35.94, 41.55, 9.4, 237.25, 0, 38.88, 
    24.56, 25.69, 0, 15.42, 0, 0, 0, 0, 467.64, 25.82, 36, 11.64, 
    112.05, 31.54, 42.08, 0, 26.86, 79.74, 0, 27.18, 17.48, 0, 
    34.95, 14.45, 43.88, 33.76, 23.24, 32.2, 16.29, 72.84, 189.99, 
    436.05, 365.6, 259.98, 329.29, 228, 158.4, 140.91, 448.95, 
    433.84, 47.11, 228.9, 193.13, 130.3, 335.73, 609.9, 202.54, 
    371.88, 332, 360.36, 219.56, 338.91, 329.94, 139.15, 262.34, 
    285.9, 357.76, 253.68, 353.35, 839.16, 368, 717.42, 840.18, 
    2081.2, 900.15, 1052.03, 705.12, 1276.65, 512.25, 838.88, 
    614.46, 734.58, 479.52, 286.38, 3020.4, 750.6, 885.96, 796.8, 
    932.49, 824.67, 1476.09, 716.76, 576.46, 528.58, 568.8, 568.8, 
    712.53, 1168.86, 1864.56, 997.26, 792.05, 1807.52, 899.25, 
    939.03, 1487.7, 1121.12, 166.5, 84.96, 78.7, 31.98, 169.2, 
    99.35, 124.2, 176.85, 116.88, 104.6, 45.43, 0, 82.44, 193.05, 
    53.5, 204.49, 135.72, 201.9, 129.76, 49.71, 50.5, 93.06, 
    239.98, 75.72, 221.54, 207.79, 218.24, 73.26, 96.4, 227.63, 
    155.4, 141.7, 280.63, 98.25, 58.4, 16.6, 30.84, 141.72, 0, 
    277.16, 313.82, 534.19, 104.74, 508.04, 67.62, 68.44, 119.7, 
    215.37, 26.92, 0, 63.24, 48.68, 11.62, 81.36, 142.5, 65.07, 
    28.06, 133.5, 126.54, 70.28, 79.62, 107.73, 36.16, 30.14, 
    31.76, 407.76, 422.24, 274.24, 317.7, 241.5, 190.3, 644.49, 
    162.17, 1104.24, 324.78, 268.24, 214.2, 449.25, 363.22, 475.57, 
    197.12, 311.63, 154.28, 461.3, 352.52, 247.69, 382.65, 395.25, 
    270.63, 399.84, 338.4, 529.48, 440.82, 394.56, 270.48, 322, 
    441.22, 353.5, 452.4, 414.96, 699.72, 89.04, 173.7, 347.6, 
    10150.24, 563.67, 353.94, 456.88, 117.92, 513, 245.48, 440.37, 
    372.36, 398.86, 334.35, 428, 410.13, 398.06, 674.87, 438.75, 
    226.16, 367.9, 416.8, 501.48, 522.6, 616.11, 421.2, 309.96, 
    423.09, 232.08, 198.06, 48.66, 109.59, 49.59, 58.05, 152.08, 
    0, 617.83, 64.66, 372.75, 32.07, 66.81, 112.24, 68.28, 83.64, 
    157.48, 145.2, 46.24, 143, 99.18, 117.5, 158.05, 61.1, 91.68, 
    67.5, 112.62, 98.21, 117.54, 58.92, 77.3, 0)), .Names = c("MONTH", 
"TRANSECT", "POLE", "TIME", "HEIGHT", "insectdens"), class = "data.frame", row.names = c(1L, 
2L, 4L, 5L, 7L, 9L, 10L, 13L, 14L, 16L, 17L, 19L, 20L, 22L, 23L, 
25L, 26L, 28L, 29L, 31L, 32L, 34L, 35L, 37L, 38L, 40L, 41L, 43L, 
44L, 46L, 47L, 49L, 50L, 53L, 55L, 56L, 58L, 60L, 61L, 64L, 65L, 
67L, 68L, 70L, 71L, 73L, 74L, 76L, 77L, 79L, 80L, 82L, 83L, 85L, 
86L, 88L, 89L, 91L, 92L, 94L, 95L, 97L, 98L, 100L, 101L, 103L, 
104L, 106L, 107L, 109L, 111L, 112L, 115L, 116L, 118L, 119L, 121L, 
122L, 124L, 125L, 127L, 128L, 130L, 131L, 133L, 134L, 136L, 137L, 
139L, 140L, 142L, 143L, 145L, 146L, 148L, 149L, 151L, 152L, 154L, 
155L, 157L, 158L, 160L, 162L, 163L, 166L, 167L, 169L, 170L, 172L, 
173L, 175L, 176L, 178L, 179L, 181L, 182L, 184L, 185L, 187L, 188L, 
190L, 191L, 193L, 194L, 196L, 197L, 199L, 200L, 202L, 203L, 205L, 
206L, 208L, 209L, 211L, 213L, 214L, 217L, 218L, 220L, 221L, 223L, 
224L, 226L, 227L, 229L, 230L, 232L, 233L, 235L, 236L, 238L, 239L, 
241L, 242L, 244L, 245L, 247L, 248L, 250L, 251L, 253L, 254L, 256L, 
259L, 260L, 262L, 264L, 265L, 268L, 269L, 271L, 272L, 274L, 275L, 
277L, 278L, 280L, 281L, 283L, 284L, 286L, 287L, 289L, 290L, 292L, 
293L, 295L, 296L, 298L, 299L, 301L, 302L, 304L, 305L, 307L, 310L, 
311L, 313L, 315L, 316L, 319L, 320L, 322L, 323L, 325L, 326L, 328L, 
329L, 331L, 332L, 334L, 335L, 337L, 338L, 340L, 341L, 343L, 344L, 
346L, 347L, 349L, 350L, 352L, 353L, 355L, 356L, 358L, 359L, 361L, 
362L, 364L, 366L, 367L, 370L, 371L, 373L, 374L, 376L, 377L, 379L, 
380L, 382L, 383L, 385L, 386L, 388L, 389L, 394L, 395L, 397L, 398L, 
400L, 401L, 403L, 404L, 406L, 407L, 409L, 410L, 412L, 413L, 415L, 
417L, 418L, 421L, 422L, 424L, 425L, 427L, 428L, 430L, 431L, 433L, 
434L, 436L, 437L, 439L, 440L, 442L, 443L, 445L, 446L, 448L, 449L, 
451L, 452L, 454L, 455L, 457L, 458L, 460L, 461L, 463L, 464L, 466L, 
468L, 469L, 472L, 473L, 475L, 476L, 478L, 479L, 481L, 482L, 484L, 
485L, 487L, 488L, 490L, 491L, 493L, 494L, 496L, 497L, 499L, 500L, 
502L, 503L, 505L, 506L, 508L, 509L, 511L, 512L, 514L, 515L, 517L, 
519L, 520L, 523L, 524L, 526L, 527L, 529L, 530L, 532L, 533L, 535L, 
536L, 538L, 539L, 541L, 542L, 544L, 545L, 547L, 548L, 550L, 551L, 
553L, 554L, 556L, 557L, 559L, 560L, 562L, 563L, 565L, 566L, 568L, 
570L, 571L, 574L, 575L, 577L, 578L, 580L, 581L, 583L, 584L, 586L, 
587L, 589L, 590L, 592L, 593L, 595L, 596L, 598L, 599L, 601L, 602L, 
604L, 605L, 607L, 608L, 610L, 611L, 613L, 616L, 617L, 619L, 621L, 
622L, 625L, 626L, 628L, 629L, 631L, 632L, 634L, 635L, 637L, 638L, 
640L, 641L, 643L, 644L, 646L, 647L, 649L, 650L, 652L, 653L, 655L, 
656L, 658L, 659L, 661L, 662L, 664L, 667L, 668L, 670L, 672L, 673L, 
676L, 677L, 679L, 680L, 682L, 683L, 685L, 686L, 688L, 689L, 691L, 
692L, 694L, 695L, 697L, 698L, 700L, 701L, 703L, 704L, 706L, 707L, 
709L, 710L, 712L, 713L, 715L, 718L, 719L, 721L, 723L, 724L, 727L, 
728L, 730L, 731L, 733L, 734L, 736L, 737L, 739L, 740L, 742L, 743L, 
745L, 746L, 748L, 749L, 751L, 752L, 754L, 755L, 757L, 758L, 760L, 
761L, 763L, 764L, 766L, 769L, 770L, 772L, 774L, 775L, 778L, 779L, 
781L, 782L, 784L, 785L, 787L, 788L, 790L, 791L, 793L, 794L, 796L, 
797L, 799L, 800L, 802L, 803L, 805L, 806L, 808L, 809L, 811L, 812L, 
814L, 815L))
1 day later
John Maindonald
# 
While I was out of contact with the list for a couple of days, I put the
following together.  It is for the data that Mike gave initially, so  
that
it assumes no replicates within POLE. The mean square for
MONTH:TIME:TRANSECT:POLE is the Residual.

The variance component for TRANSECT is effectively zero.  For
purposes of exposition, it can be set to zero.

The terms formula MONTH/TIME then suffices. lme() gives the numbers of
groups as:

Number of Observations: 286
Number of Groups:
           MONTH TIME %in% MONTH
               4              16

Observe that the design is then very nearly a complete balanced design:

 > with(temp, table(MONTH,TIME))
      TIME
MONTH  1  2  3  4
     4 18 17 18 18
     5 18 18 18 17
     6 18 18 18 18
     7 18 18 18 18

Observe that
   286/16=17.875
less than 18 because two of the TIME:MONTH combinations have
only 17 values.  The harmonic mean, which is 17.87, seems however
preferable.

If the design were balanced, with replicates at
MONTH/TIME/POLE level equal to n1=4/n2=4/n3=17.87
then the mean squares would estimate, respectively

4 x 17.87 s1^2 +  17.87 s2^2 + s3^2  = 4 x 17.87 (s1^2 +  s2^2 /4 +  
s3^2 / (4*17.87)
17.87 s2^2 + s3^2  = 17.87 ( s2^2 + s3^2 / 17.87)
s3^2

In the first two cases, a second version of the formula is given
that makes (to me, at least) better intuitive sense.

We have
          Variance     Anova        DF     No. of repeats
          Component    mean square         per 'parent'
MONTH    0.75287      74.10751      3     4
TIME     1.07395      20.52433     12     4
Residual 1.20354      1.20356     270     17.87 (harmonic mean)

Observe that the residual mean square estimates the residual variance.
The TIME variance component is calculated pretty much as
(20.52433 - 1.20356)/17.87 = 1.0813, which does not quite agree
with the ML estimate (nor should it; equating mean squares to expected
mean squares is not the same as REML!)).

The statistical error is affected both by the statistical error in  
the TIME
anova mean square, and by the statistical error in the Residual mean
square.  The variance formula that is given below for MONTH can be
readily adapted for this case also.

The MONTH variance component is calculated pretty much as: (74.10751 -
20.52433)/(4*17.87) = 0.7497 Again the agreement with the REML
estimate is, quite properly, not perfect.

The estimate for s1^2 (MONTH} has a statistical error that is a
compound of the errors in the ANOVA mean squares for both TIME and
MONTH.  The variances of the two anova mean squares (both sample
values of chi-squared statistics, under the usual assumptions) add,
while the quantities themselves are subtracted.  The SE (sqrt of
variance), for the estimate s1^2 of sig1^2, is

sqrt( (n2 * n3 * sig1^2 + n3 * sig2^2 + sig3^2)^2 / nu1
      + (n3 * sig2^2 + sig3^2)^2 / nu2 ) / (n2*n3)

[The n's are the numbers of repeats.  The nu's are degrees of freedom.]

Variances are not, for quantities that are differences multiples of
chi-squared statistics, a good basis for inference. (Here I am tempted
to make rude comments about over-reliance on variances in much
sample  survey work!). The variance calculations may however be useful
in giving an idea of the relative contributions of the different  
sources of
statistical noise.

I'd expect, though I have not gone through the arithmetic in detail,
that the estimate for SE[s1^2] will increase, if we allow for the
possibility that the TRANSECT component of variance, even though
estimated as zero, may actually be greater than zero.

John Maindonald             email: john.maindonald at anu.edu.au
phone : +61 2 (6125)3473    fax  : +61 2(6125)5549
Centre for Mathematics & Its Applications, Room 1194,
John Dedman Mathematical Sciences Building (Building 27)
Australian National University, Canberra ACT 0200.
On 5 Sep 2007, at 4:54 AM, Mike Dunbar wrote:

            

      
      
      
    

  
  


  


      

    

    
7 days later

    

      
      

        


  

        

      David Afshartous
    

    

      
      #
    

  


  

      
All,
I don't have the initial e-mails from this thread and when I searched the
mail archives for "explaining lme variance component results" it doesn't
come up.  Do the SIGs need to be searched differently  then r-help?
Thanks,
David

        
On 9/6/07 5:07 AM, "John Maindonald" <john.maindonald at anu.edu.au> wrote: