[R-sig-ME] Different estimates for mixed effects Logistic regression and pwrssUpdate Error message with binomial glmer

[David Duffy]

Hi.

> https://stats.stackexchange.com/questions/499269/different-estimates-for-mixed-effects-logistic-regression-and-pwrssupdate-error

> The data are shared below.
> https://drive.google.com/file/d/1ZTiDUhTcoyOWUCa2vjXR95VlcpyzxGiU/view?usp=sharing

I'm having a little trouble understanding exactly what you want. Are these the actual data you want to analyse? Since you have such a simple random effects art of the model, there are many alternatives, including nonparametric models where you don't have to specify a distribution for the random effects (the default you are using here is Gaussian).  These programs will *all* give slightly different answers. Check out the different published results different programs give for the dataset from Crowder, "Beta-binomial ANOVA for proportions." Applied statistics 1978: 34-37. 

I tried out the glmmML package on your data, as I find it well behaved - see its documentation. One advantage (for me) is that it gives significance tests for the random effects.

glmmML(formula = y ~ x1 + x2, data = x, cluster = group, control = list(maxit = 1000),  method = "ghq") 

                          coef se(coef)      z Pr(>|z|)
(Intercept) -1.1521   0.8386 -1.374    0.169
x1                 0.6044   0.4626  1.306    0.191
x2                1.1099   0.8141  1.363    0.173
Scale parameter in mixing distribution:  2.12 gaussian 
Std. Error:                              2.599 
        LR p-value for H_0: sigma = 0:  0.2082 

Call:  glmmML(formula = y ~ x1 + x2, data = x[-55, ], cluster = group,  control = list(maxit = 1000), method = "ghq") 
               coef se(coef)      z Pr(>|z|)
(Intercept) -1.1491   0.8369 -1.373    0.170
x1                 0.5981   0.4630  1.292    0.196
x2                1.0992   0.8145  1.350    0.177
Scale parameter in mixing distribution:  2.101 gaussian 
Std. Error:                              2.615 
        LR p-value for H_0: sigma = 0:  0.2116 

Call:  glmmML(formula = y ~ x1 + x2, data = x, cluster = group, prior = "logistic",  control = list(maxit = 1000), method = "ghq") 
               coef se(coef)      z Pr(>|z|)
(Intercept) -1.4260   1.1504 -1.240    0.215
x1           0.8034   0.5658  1.420    0.156
x2           1.4835   1.0384  1.429    0.153
Scale parameter in mixing distribution:  1.76 logistic 
Std. Error:                              1.548 
        LR p-value for H_0: sigma = 0:  0.1227 

Call:  glmmML(formula = y ~ x1 + x2, data = x, cluster = group, prior = "cauchy",  control = list(maxit = 1000), method = "ghq") 

               coef se(coef)      z Pr(>|z|)
(Intercept) -1.5690   1.0233 -1.533  0.12500
x1           0.9551   0.3442  2.775  0.00553
x2           1.7864   0.6710  2.662  0.00776
Scale parameter in mixing distribution:  1.44 cauchy 
Std. Error:                              0.4347 
        LR p-value for H_0: sigma = 0:  0.5 
Residual deviance: 60.08 on 71 degrees of freedom 	AIC: 68.08

These look comparable to results from non-R programs.