Dear All

I am currently analysing ecological data looking at the effect of
agricultural intensity on a rnage of biodiversity measures. I have data from
26 sites grouped into 9 clusters. I am therefore running a linear mixed
model using the nlme package. I am aware that there is issues in mixed
models of determing the number of degrees of freedom and hence the p values
but as these are reported by nlme I have a query I was hoping I might be
able to resolve.

The model I have run is below:

beerichmix1<-lme(MeanBeeSppRich~cultivate, data=bees, random=~1|cluster,
 method="ML")

Linear mixed-effects model fit by maximum likelihood
Data: bees
      AIC      BIC    logLik
 165.0555 170.0879 -78.52774

Random effects:
Formula: ~1 | cluster
       (Intercept) Residual
StdDev:    1.899494 4.636802

Fixed effects: MeanBeeSppRich ~ cultivate
              Value Std.Error DF   t-value p-value
(Intercept) 38.19658  4.267833 16  8.949878  0.0000
cultivate   -0.22248  0.059514 16 -3.738322  0.0018
Correlation:
         (Intr)
cultivate -0.962

Standardized Within-Group Residuals:
      Min         Q1        Med         Q3        Max
-1.8576001 -0.5935824  0.0859121  0.7933137  1.9665311

Number of Observations: 26
Number of Groups: 9

I can see that the df has been calculated as 26 sites - 9 clusters - 1 = 16.
However, this is clearly a misrepresentation of the data structure in the
model and hence so is the p value. Using some sort of correction (i.e.
Kenward-Roger correction), even if not precise would offer a slightly better
solution in my circumstance than not correcting at all there fore I was
wondering if anybody has addressed this in R or if people are just moving
more and more away from reporting df and p values.

Any help would be much appreciated.

Many thanks
AR




-- 
Dr Anna R. Renwick
Research Ecologist
British Trust for Ornithology,
The Nunnery,
Thetford,
Norfolk,
IP24 2PU,
UK
Tel: +44 (0)1842 750050; Fax: +44 (0)1842 750030

	[[alternative HTML version deleted]]


