Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

Multiple glm in for loop

Tags:

for-loop

r

glm

I have an R dataframe, strongly simplified as:

id <- rep(1:2, c(6,8))
correct <- sample(0:1,14,TRUE)
phase <- c(rep("discr",3),rep("rev",3), rep("discr",4),rep("rev",4))
dat <- data.frame(id,correct,phase)

with id as my subjects (in reality I have a lot more than 2), correct = responses coded as incorrect (0) or correct (1), and the phases Discrimination and Reversal (within-subjects factor).

I want to perform a logistic regression in the form of

glm(correct~phase, dat, family="binomial")

later possibly adding additional predictors. However, since I have a varying amount of data for each subject, I would like to perform glm() seperately for each subject and later compare the coefficients with ANOVA for group effects. I would like to do this in a for loop in the form of

for(i in seq_along(dat$id)){
   my_glm[i] <- glm(correct~list,dat[dat$id==i,],family="binomial")
}

but keep receiving the error message

>Error in `contrasts<-`(`*tmp*`, value = contr.funs[1 + isOF[nn]]) :
 contrasts can be applied only to factors with 2 or more levels.

I have checked my data and there is no factor which contains only one level. All subjects gave at least one incorrect and one correct response, and all took part in Discrimination and Reversal. The function works outside the loop when I specify a particular subject.

like image 420
Lucy Vanes Avatar asked Feb 17 '23 18:02

Lucy Vanes


1 Answers

Here's an R Base solution:

> lapply(split(dat, dat$id), function(x) coef(summary(glm(correct~phase,family="binomial",data=x))))
$`1`
                 Estimate Std. Error       z value  Pr(>|z|)
(Intercept) -6.931472e-01   1.224745 -5.659524e-01 0.5714261
phaserev    -3.845925e-16   1.732050 -2.220446e-16 1.0000000

$`2`
                Estimate Std. Error      z value Pr(>|z|)
(Intercept) 3.356998e-16   1.000000 3.356998e-16 1.000000
phaserev    1.098612e+00   1.527524 7.192109e-01 0.472011
like image 112
Jilber Urbina Avatar answered Feb 27 '23 06:02

Jilber Urbina