I get an error message when fitting a full random two-level null model with the package R2MLwiN. My dataframe is a subset of a Multiple Indicator Cluster Survey developed by UNICEF for Mozambique. My response variable is agem.18 a binary ("Yes", "No") indicating whether the woman got married before the age of 18 y.o. or not.
> table(moz.20$agem.18, useNA = "ifany")
No Yes
5934 5405
My two levels are, in descending order, province and w.id. This is the code I run
# Random effect
F1 <- logit(agem.18) ~ 1 + (1 | province) + (1 | w.id)
rand.eff <- runMLwiN(Formula = F1, D = "Binomial", data = moz.20)
This is the error message I get
Invalid link function specified: NA
Error in runMLwiN(Formula = F1, D = "Binomial", data = df) :
I initially thought that the logit function was masked by the car package, but this is not the case. I also thought that the proble was in the denom call in the link function, but I read that R2MLwiN should automatically create the denom as a set of 1s. I do not get any error if I use the package lme4 with the same data, variable and levels:
(fit <- glmer(agem.18 ~ 1 + (1 | province), family = binomial("logit"), data = moz.20)
Generalized linear mixed model fit by maximum likelihood (Laplace
Approximation) [glmerMod]
Family: binomial ( logit )
Formula: agem.18 ~ 1 + (1 | province)
Data: moz.20
AIC BIC logLik deviance df.resid
14907.498 14922.170 -7451.749 14903.498 11337
Random effects:
Groups Name Std.Dev.
province (Intercept) 0.5366
Number of obs: 11339, groups: province, 11
Fixed Effects:
(Intercept)
-0.04992
I do not encounter the same problem if I use a very similar formula included in the demo UserGuide09.R for the package R2MLwiN.
(mymodel1 <- runMLwiN(logit(use) ~ 1 + lc, D = "Binomial", data = bang))
My only guess at the moment is that for some reasons, R2MLwiN fails to recognise my response variable agem.18 as binary.
Any suggestion?
Thanks
Manolo
Problem solved. It appears there is a bug in R2MLwiN, as noted by Chris Charlton in the R2MLwiN forum (see post here). I replicate the code written by Chris that shows the bug.
> library(R2MLwiN)
> data(bang, package = "R2MLwiN")
> (mymodel1 <- runMLwiN(logit(use) ~ 1 + lc, D = "Binomial", data = bang))
-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-
MLwiN (version: 3.02) multilevel model (Binomial)
Estimation algorithm: IGLS MQL1 Elapsed time : 0.6s
Number of obs: 2867 (from total 2867) The model converged after 4 iterations.
Log likelihood: NA
Deviance statistic: NA
---------------------------------------------------------------------------------------------------
The model formula:
logit(use) ~ 1 + lc
Level 1: l1id
---------------------------------------------------------------------------------------------------
The fixed part estimates:
Coef. Std. Err. z Pr(>|z|) [95% Conf. Interval]
Intercept -1.12288 0.08348 -13.45 3.05e-41 *** -1.28650 -0.95926
lcOne_child 0.93275 0.12156 7.67 1.676e-14 *** 0.69450 1.17100
lcTwo_children 1.09250 0.12509 8.73 2.466e-18 *** 0.84733 1.33768
lcThree_plus 0.87224 0.10302 8.47 2.523e-17 *** 0.67033 1.07416
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
---------------------------------------------------------------------------------------------------
The random part estimates at the l1id level:
Coef. Std. Err.
var_bcons_1 1.00000 0.00000
-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-
> bang$use.test <- bang$use
> (mymodel1 <- runMLwiN(logit(use.test) ~ 1 + lc, D = "Binomial", data = bang))
Invalid link function specified: NA
Error in runMLwiN(logit(use.test) ~ 1 + lc, D = "Binomial", data = bang) :