[R] variance estimates in lme biased?

[Gary Allison]

Peter Dalgaard wrote:
[snip]
> 
> Or, try looking at a smaller example where things can be worked out
> explicitly: One-way ANOVA with random btw.group variation. Say 5
> groups and 3 obs per group. If I got this right (please do check!),
> the estimate of the between-group variance is 1/3 times the difference
> between two chi^2/f distributed variables with 4 and 10 DF
> respectively. This will become negative about half the time, and lme
> (and similar code) will set it to zero in that case. Now
> 
> 
>>sd.sim <- sqrt(pmax(1/3*(rchisq(1000,4)/4 - rchisq(1000,10)/10),0))
>>summary(sd.sim)
> 
>    Min. 1st Qu.  Median    Mean 3rd Qu.    Max.
>  0.0000  0.0000  0.0000  0.2017  0.3955  1.2160
> 
> does not seem to be too far from what Gary has been experiencing. 
> 
> Obviously, the fact that the estimator is censored at zero will make
> it biased, but an extended estimator (allowing negative values) is
> unbiased. 
> 
> 
>>var.sim <- 1/3*(rchisq(1000,4)/4 - rchisq(1000,10)/10)
>>summary(var.sim)
> 
>      Min.   1st Qu.    Median      Mean   3rd Qu.      Max.
> -0.796100 -0.196900 -0.042720 -0.007667  0.138500  1.104000
> 
> 
> 
As Peter implies, when I did the simulation in SAS with PROC VARCOMP, 
using the model option 'Type1' which is an extended estimator, estimates 
for level 1 average to 0:
http://david.science.oregonstate.edu/~allisong/R/sas_Type1_5.pdf
and the range is still quite large.

Thanks to the help of Drs. Lumley, Bates and Dalgaard, I'm beginning to 
understand what's going on.  But now it seems that much harder to 
interpret my experimental results -- high variance in both level 1 and 
level 2 at relatively low sample sizes -- given all this.  If variance 
at one level can 'bleed' into another, sometimes to a large degree, is 
this analysis even useful at relatively small sample sizes?

thanks again,
Gary