Tuesday, January 23, 2007

F lmer

Hi Marie,
do you really need to use the quasibinomial distribution? Is there evidence of overdispersion in the maximal model you fit?

If not, then you can use the binomial distribution.
That allows you to use likelihood-ratio tests for the fixed and random effects.
using
anova(model,modelsimplified,test="Chisq")

If there is clear overdispersion, then I guess the best thing to do would be resampling according to the model, if that works for the quasibinomial "distribution". You need to use mcmcsamp() for that.
Otherwise, it will become hand-work to get some of those hotly debated F-ratios out....
Cheers, Tom

lmer P-value

hi everybody
Do you know the simplest and correct way to produce P-values of lmer and of the ANOVA (F test) on the lmer model in R?
My lmer model is of this form:
model<-lmer(pref1~brood.size*logtarsus*cohort+(1birth.nest/foster.brood),family=quasibinomial,method='Laplace')
summary(model)
anova(model,test="F")
I read an extensive eplanation by Douglas Bates on R help archive named "[R] lmer, p-values and all that". However, nothing is mentionned regarding Laplace. So i do not know what to do.
Could someone give me advice?
many thanks
Marie

Tuesday, January 16, 2007

Longitudinal Data

Finally, the participants in the course Longitudinal and Incomplete Data have been invited. Transcripts for R of the SAS code presented during the course will be added to the weblog soon.
Cheers, Tom

Loglinear models

Hi Nicola,
the best manual is http://www.statslab.cam.ac.uk/~pat/Splusdiscrete2.pdf, based on the book by Agresti.
You can basically use glm with a poisson link for log-linear models, I guess it is also discussed in the Venables and Ripley book.
Good luck, Tom

Loglin models and contingency tables

Hi Tom,

I'm trying to use loglinear models to investigate contingency tables with frequencies, but I can't figure out how to get the standardized residuals, odds ratios etc. Also the help function doesn't mention them, and the R-manuals I have don't even mention loglin models at all.
Do I have to calculate them by hand or is there perhaps a good manual on loglin models?
Thanks in advance!
Nicola

ps I have a 3-way table with one binary dependent factor and 2 explanatory factors with 3 levels each.

Friday, January 05, 2007

lmer, random effect with four levels

Dear Amir,
the "place" random effect in your dataset has only four levels, and they are nested within plants. The error which is returned states that the variance covariance matrix is singular, meaning that variation in intercept, slope and their covariance cannot be fitted. If you fit the random effect as (1|place), you will read that the random effect variance is effectively zero.
It does not make sense to include place both as fixed and as random, as in you second to fourth models. Since the number of levels of the random effects is very low, you could even just model them as fixed, using glm. The idea of random effects is that they are a representative sample from a large population, which is hard to defend with just two levels sampled per group!
More importantly, the distribution of your data does not look poisson at al. There are many zero observations. I would convert the counts to presence/absence data and analyse those as bernoulli probabilities. That works to some extent:

b<-ifelse(Predators.Allspecies.>0,1,0)
model<-lmer(b~julian.date+offset(log.time)+(1|place),family=binomial)

I've added you models to your post, so other visitors to this blog can also read and respond to your question,
all the best, Tom

Thursday, January 04, 2007

Mixed effects model

Beste Tom,
Tijdens de laatste GLM cursus hebben we samen een model gemaakt voor de
analyse van een dataset. Ik probeer het model nu te gebruiken, maar krijg
steeds fout meldingen in R. Helaas heb ik de R file van de analyse die we
samen gedaan hebben op een verkeerde wijze opgeslagen, waardoor ik niet
meer terug kan halen wat ik fout doe. Kan je misschien kijken wat er
mis gaat?

Ik stuur de dataset en de Tin R file met het model die we gebruikt hebben
als attachment naar je e-mail account. De dataset gaat over voorkomen van insecten in 4 plantages, twee van ieder plant soort (Guava en Eucalyptus), gedurende 2
jaar. De data is verzameld door bomen in de plantages te bemonsteren en
daarom nemen we de zoektijd per boom als weegfactor voor de datapunt. Ik
probeer erachter te komen of de plantsoort effect heeft op het voorkomen
van de verschillende insecten groepen.

Alvast bedankt,
Amir.


Here's Amir's code:

setwd("D:/Amir-18-9-06/field-data/All data/analysis")
library(lme4)
data<-read.table("pdata.txt",header=T)
attach(data)
model<-lmer(Predators.Allspecies.~plant*julian.date+offset(log.time)
+(1+julian.date|place),family=poisson,method="Laplace")
model<-lmer(Predators.Allspecies.~place*julian.date+offset(log.time)
+(1+julian.date|place),family=poisson,method="Laplace")
model<-lmer(Predators.Allspecies.~place*julian.date+offset(log.time)
+(1+julian.date|place),family=quasipoisson,method="Laplace")
model<-lmer(Predators.Allspecies.~place*julian.date+offset(log.time)
+(1+julian.date|place/tree.code),family=poisson,method="Laplace")