Skip to content

Predicted probabilites with CIs for multilevel logistic regression with prior weights

13 messages · Mollie Brooks, d@iuedecke m@iii@g oii uke@de, Rose, Charles E. (CDC/DDNID/NCBDDD/OD) +1 more

Original
Sam Crawley
# 
Hello all,

I am doing a multilevel binomial logistic regression using lme4, and the survey data I am using requires weights to be used. I would like to calculate various predicted probabilities with confidence intervals based on the estimated model. The predict function obviously doesn't give me standard errors, and the recommended method to get these is to use the bootMer function.

However, in my case, the weights provided are not integers, and the bootMer function exits with an error if the weights are not integers (I raised a GitHub issue about this, and was pointed to this list: https://github.com/lme4/lme4/issues/524 ).

So my question is, what is the best way to calculate the predicted probabilities (with confidence intervals) in my case?

I would appreciate any help you can give me, and I'm happy to provide more details if required.

Thanks,
Sam Crawley.
d@iuedecke m@iii@g oii uke@de
# 
Hi Sam,

you could the "ggeffects" package
(https://strengejacke.github.io/ggeffects/), and there is also an example
for a logistic mixed effects model
(https://strengejacke.github.io/ggeffects/articles/practical_logisticmixedmo
del.html), which might help you.

For binomial models, using weights often results in the following warning:
#> non-integer #successes in a binomial glm!

However, CIs for the predicted probabilities can be calculated nevertheless
(at least in my quick example). Note that afaik, mixed models in R do
correctly not account for sampling weights. However, Thomas Lumley, author
of the survey-package, works on a survey-function for mixed models
(https://github.com/tslumley/svylme), probably the GitHub version is quite
stable (haven't tested yet).

An alternative would be the "scale_weights()" function from the
sjstats-package
(https://strengejacke.github.io/sjstats/articles/mixedmodels-statistics.html
#rescale-model-weights-for-complex-samples ), which rescales sampling
weights so they can be used as "weights" for the mixed models function you
have in R (lme4, lme, ...).

Based on that function, I have a small example that demonstrates how to
compute predicted probabilities for mixed models with (sampling) weights
(ignore the warnings, this is just for demonstration purposes):

library(lme4)
library(sjstats) # for scale_weights() and sample data
library(ggeffects) # for ggpredict()

data(nhanes_sample)
set.seed(123)
nhanes_sample$bin <- rbinom(nrow(nhanes_sample), 1, prob = .3)
nhanes_sample <- scale_weights(nhanes_sample, SDMVSTRA, WTINT2YR)

m <- glmer(
  bin ~ factor(RIAGENDR) * age + factor(RIDRETH1) + (1 | SDMVPSU),
  family = binomial(),
  data = nhanes_sample,
  weights = svywght_a
)

ggpredict(m, c("age", "RIAGENDR")) %>% plot()


Best
Daniel

-----Urspr?ngliche Nachricht-----
Von: R-sig-mixed-models <r-sig-mixed-models-bounces at r-project.org> Im
Auftrag von Sam Crawley
Gesendet: Montag, 10. Juni 2019 10:36
An: r-sig-mixed-models at r-project.org
Betreff: [R-sig-ME] Predicted probabilites with CIs for multilevel logistic
regression with prior weights

Hello all,

I am doing a multilevel binomial logistic regression using lme4, and the
survey data I am using requires weights to be used. I would like to
calculate various predicted probabilities with confidence intervals based on
the estimated model. The predict function obviously doesn't give me standard
errors, and the recommended method to get these is to use the bootMer
function.

However, in my case, the weights provided are not integers, and the bootMer
function exits with an error if the weights are not integers (I raised a
GitHub issue about this, and was pointed to this list:
https://github.com/lme4/lme4/issues/524 ).

So my question is, what is the best way to calculate the predicted
probabilities (with confidence intervals) in my case?

I would appreciate any help you can give me, and I'm happy to provide more
details if required.

Thanks,
Sam Crawley.

_______________________________________________
R-sig-mixed-models at r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models

--

_____________________________________________________________________

Universit?tsklinikum Hamburg-Eppendorf; K?rperschaft des ?ffentlichen Rechts; Gerichtsstand: Hamburg | www.uke.de
Vorstandsmitglieder: Prof. Dr. Burkhard G?ke (Vorsitzender), Prof. Dr. Dr. Uwe Koch-Gromus, Joachim Pr?l?, Marya Verdel
_____________________________________________________________________

SAVE PAPER - THINK BEFORE PRINTING
d@iuedecke m@iii@g oii uke@de
# 
Should be: mixed models in R do *currently* not account for sampling weights

-----Urspr?ngliche Nachricht-----
Von: R-sig-mixed-models <r-sig-mixed-models-bounces at r-project.org> Im
Auftrag von d.luedecke at uke.de
Gesendet: Montag, 10. Juni 2019 17:31
An: 'Sam Crawley' <sam_crawley at warpmail.net>;
r-sig-mixed-models at r-project.org
Betreff: Re: [R-sig-ME] Predicted probabilites with CIs for multilevel
logistic regression with prior weights

Hi Sam,

you could the "ggeffects" package
(https://strengejacke.github.io/ggeffects/), and there is also an example
for a logistic mixed effects model
(https://strengejacke.github.io/ggeffects/articles/practical_logisticmixedmo
del.html), which might help you.

For binomial models, using weights often results in the following warning:
#> non-integer #successes in a binomial glm!

However, CIs for the predicted probabilities can be calculated nevertheless
(at least in my quick example). Note that afaik, mixed models in R do
correctly not account for sampling weights. However, Thomas Lumley, author
of the survey-package, works on a survey-function for mixed models
(https://github.com/tslumley/svylme), probably the GitHub version is quite
stable (haven't tested yet).

An alternative would be the "scale_weights()" function from the
sjstats-package
(https://strengejacke.github.io/sjstats/articles/mixedmodels-statistics.html
#rescale-model-weights-for-complex-samples ), which rescales sampling
weights so they can be used as "weights" for the mixed models function you
have in R (lme4, lme, ...).

Based on that function, I have a small example that demonstrates how to
compute predicted probabilities for mixed models with (sampling) weights
(ignore the warnings, this is just for demonstration purposes):

library(lme4)
library(sjstats) # for scale_weights() and sample data
library(ggeffects) # for ggpredict()

data(nhanes_sample)
set.seed(123)
nhanes_sample$bin <- rbinom(nrow(nhanes_sample), 1, prob = .3)
nhanes_sample <- scale_weights(nhanes_sample, SDMVSTRA, WTINT2YR)

m <- glmer(
  bin ~ factor(RIAGENDR) * age + factor(RIDRETH1) + (1 | SDMVPSU),
  family = binomial(),
  data = nhanes_sample,
  weights = svywght_a
)

ggpredict(m, c("age", "RIAGENDR")) %>% plot()


Best
Daniel

-----Urspr?ngliche Nachricht-----
Von: R-sig-mixed-models <r-sig-mixed-models-bounces at r-project.org> Im
Auftrag von Sam Crawley
Gesendet: Montag, 10. Juni 2019 10:36
An: r-sig-mixed-models at r-project.org
Betreff: [R-sig-ME] Predicted probabilites with CIs for multilevel logistic
regression with prior weights

Hello all,

I am doing a multilevel binomial logistic regression using lme4, and the
survey data I am using requires weights to be used. I would like to
calculate various predicted probabilities with confidence intervals based on
the estimated model. The predict function obviously doesn't give me standard
errors, and the recommended method to get these is to use the bootMer
function.

However, in my case, the weights provided are not integers, and the bootMer
function exits with an error if the weights are not integers (I raised a
GitHub issue about this, and was pointed to this list:
https://github.com/lme4/lme4/issues/524 ).

So my question is, what is the best way to calculate the predicted
probabilities (with confidence intervals) in my case?

I would appreciate any help you can give me, and I'm happy to provide more
details if required.

Thanks,
Sam Crawley.

_______________________________________________
R-sig-mixed-models at r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models

--

_____________________________________________________________________

Universit?tsklinikum Hamburg-Eppendorf; K?rperschaft des ?ffentlichen
Rechts; Gerichtsstand: Hamburg | www.uke.de
Vorstandsmitglieder: Prof. Dr. Burkhard G?ke (Vorsitzender), Prof. Dr. Dr.
Uwe Koch-Gromus, Joachim Pr?l?, Marya Verdel
_____________________________________________________________________

SAVE PAPER - THINK BEFORE PRINTING

_______________________________________________
R-sig-mixed-models at r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models

--

_____________________________________________________________________

Universit?tsklinikum Hamburg-Eppendorf; K?rperschaft des ?ffentlichen Rechts; Gerichtsstand: Hamburg | www.uke.de
Vorstandsmitglieder: Prof. Dr. Burkhard G?ke (Vorsitzender), Prof. Dr. Dr. Uwe Koch-Gromus, Joachim Pr?l?, Marya Verdel
_____________________________________________________________________

SAVE PAPER - THINK BEFORE PRINTING
Mollie Brooks
# 
I?m still trying to get a handle of the different definitions of "weights" but I believe we implemented sampling weights in glmmTMB. We do this by weighting the log-likelihood contribution of each observation. I think this is different from prior weights if you mean Bayesian priors. There has been some discussion of the different implementations of "weights" in different R functions (link below) and we still need to update the documentation for glmmTMB 
https://github.com/glmmTMB/glmmTMB/issues/285

Here?s a binomial example:

library(glmmTMB)
set.seed(123)
n=100
dat=data.frame(trials=rpois(n, lambda=50), rownum=1:n)
dat$success=rbinom(n, dat$trials, prob=.3)
dat$rep=sample(1:5, size=n, replace=TRUE) #each observation is repeated 1 to 5 times
rows=rep(dat$rownum, each=1, times=dat$rep)
dat_disaggregated=dat[rows, ]

summary(glmmTMB(cbind(success, trials-success)~1, weights=rep, dat, family=binomial))
summary(glmmTMB(cbind(success, trials-success)~1, dat_disaggregated, family=binomial))

and it works with non-integer weights

summary(glmmTMB(cbind(success, trials-success)~1, weights=rep/5, dat, family=binomial))

cheers,
Mollie

  
  


  


      

    

    

      
      

        


  

        

      d@iuedecke m@iii@g oii uke@de
    

    

      
      #
    

  


  

      
I think that Sam is talking about ?sampling? or ?survey? weights (as compared to analytical or frequency weights, used by ?normal? regression models).

 

The issue you?re referring to is referenced by another issue (https://github.com/glmmTMB/glmmTMB/issues/440), which in turn shows an example from Cross Validated:

https://stats.stackexchange.com/questions/57107/use-of-weights-in-svyglm-vs-glm

 

If I use that example, and add a third model fitted with glmmTMB, I get following result when comparing the weights from the fitted objects:

 

library(glmmTMB)

glm2 <- glmmTMB(re78 ~ treat, weights = w , data = lalonde)

cbind(glm1$weights, glm11$weights, glm2$frame$`(weights)`)

#>        [,1]     [,2]     [,3]

#> 1 1.4682453 2.108394 2.108394

#> 2 0.9593877 1.377677 1.377677

#> 3 0.7489954 1.075554 1.075554

#> 4 0.7319955 1.051143 1.051143

#> 5 0.7283328 1.045883 1.045883

#> 6 0.7244569 1.040317 1.040317

 

As you can see, ?glm? and ?glmmTMB? produce the same weights, while the survey-package has different weights? I?m not sure that the weights implemented in glmmTMB are actually ?sampling? weights (for surveys, as implemented in the survey package), or how to reproduce such weights using glmmTMB.

 

Von: Mollie Brooks <mollieebrooks at gmail.com> 
Gesendet: Montag, 10. Juni 2019 19:04
An: Sam Crawley <sam_crawley at warpmail.net>; Help Mixed Models <r-sig-mixed-models at r-project.org>
Cc: d.luedecke at uke.de
Betreff: Re: [R-sig-ME] Predicted probabilites with CIs for multilevel logistic regression with prior weights

        
On 10Jun 2019, at 17:33, <d.luedecke at uke.de <mailto:d.luedecke at uke.de> > <d.luedecke at uke.de <mailto:d.luedecke at uke.de> > wrote:

        
mixed models in R do correctly not account for sampling weights


Should be: mixed models in R do *currently* not account for sampling weights

 

I?m still trying to get a handle of the different definitions of "weights" but I believe we implemented sampling weights in glmmTMB. We do this by weighting the log-likelihood contribution of each observation. I think this is different from prior weights if you mean Bayesian priors. There has been some discussion of the different implementations of "weights" in different R functions (link below) and we still need to update the documentation for glmmTMB 

https://github.com/glmmTMB/glmmTMB/issues/285

 

Here?s a binomial example:

 

library(glmmTMB)

set.seed(123)

n=100

dat=data.frame(trials=rpois(n, lambda=50), rownum=1:n)

dat$success=rbinom(n, dat$trials, prob=.3)

dat$rep=sample(1:5, size=n, replace=TRUE) #each observation is repeated 1 to 5 times

rows=rep(dat$rownum, each=1, times=dat$rep)

dat_disaggregated=dat[rows, ]





summary(glmmTMB(cbind(success, trials-success)~1, weights=rep, dat, family=binomial))

summary(glmmTMB(cbind(success, trials-success)~1, dat_disaggregated, family=binomial))





and it works with non-integer weights





summary(glmmTMB(cbind(success, trials-success)~1, weights=rep/5, dat, family=binomial))





cheers,

Mollie






-----Urspr?ngliche Nachricht-----
Von: R-sig-mixed-models <r-sig-mixed-models-bounces at r-project.org <mailto:r-sig-mixed-models-bounces at r-project.org> > Im
Auftrag von d.luedecke at uke.de <mailto:d.luedecke at uke.de> 
Gesendet: Montag, 10. Juni 2019 17:31
An: 'Sam Crawley' <sam_crawley at warpmail.net <mailto:sam_crawley at warpmail.net> >;
r-sig-mixed-models at r-project.org <mailto:r-sig-mixed-models at r-project.org> 
Betreff: Re: [R-sig-ME] Predicted probabilites with CIs for multilevel
logistic regression with prior weights

Hi Sam,

you could the "ggeffects" package
(https://strengejacke.github.io/ggeffects/), and there is also an example
for a logistic mixed effects model
(https://strengejacke.github.io/ggeffects/articles/practical_logisticmixedmo
del.html), which might help you.

For binomial models, using weights often results in the following warning:
#> non-integer #successes in a binomial glm!

However, CIs for the predicted probabilities can be calculated nevertheless
(at least in my quick example). Note that afaik, mixed models in R do
correctly not account for sampling weights. However, Thomas Lumley, author
of the survey-package, works on a survey-function for mixed models
(https://github.com/tslumley/svylme), probably the GitHub version is quite
stable (haven't tested yet).

An alternative would be the "scale_weights()" function from the
sjstats-package
(https://strengejacke.github.io/sjstats/articles/mixedmodels-statistics.html
#rescale-model-weights-for-complex-samples ), which rescales sampling
weights so they can be used as "weights" for the mixed models function you
have in R (lme4, lme, ...).

Based on that function, I have a small example that demonstrates how to
compute predicted probabilities for mixed models with (sampling) weights
(ignore the warnings, this is just for demonstration purposes):

library(lme4)
library(sjstats) # for scale_weights() and sample data
library(ggeffects) # for ggpredict()

data(nhanes_sample)
set.seed(123)
nhanes_sample$bin <- rbinom(nrow(nhanes_sample), 1, prob = .3)
nhanes_sample <- scale_weights(nhanes_sample, SDMVSTRA, WTINT2YR)

m <- glmer(
 bin ~ factor(RIAGENDR) * age + factor(RIDRETH1) + (1 | SDMVPSU),
 family = binomial(),
 data = nhanes_sample,
 weights = svywght_a
)

ggpredict(m, c("age", "RIAGENDR")) %>% plot()


Best
Daniel

-----Urspr?ngliche Nachricht-----
Von: R-sig-mixed-models <r-sig-mixed-models-bounces at r-project.org <mailto:r-sig-mixed-models-bounces at r-project.org> > Im
Auftrag von Sam Crawley
Gesendet: Montag, 10. Juni 2019 10:36
An: r-sig-mixed-models at r-project.org <mailto:r-sig-mixed-models at r-project.org> 
Betreff: [R-sig-ME] Predicted probabilites with CIs for multilevel logistic
regression with prior weights

Hello all,

I am doing a multilevel binomial logistic regression using lme4, and the
survey data I am using requires weights to be used. I would like to
calculate various predicted probabilities with confidence intervals based on
the estimated model. The predict function obviously doesn't give me standard
errors, and the recommended method to get these is to use the bootMer
function.

However, in my case, the weights provided are not integers, and the bootMer
function exits with an error if the weights are not integers (I raised a
GitHub issue about this, and was pointed to this list:
https://github.com/lme4/lme4/issues/524 ).

So my question is, what is the best way to calculate the predicted
probabilities (with confidence intervals) in my case?

I would appreciate any help you can give me, and I'm happy to provide more
details if required.

Thanks,
Sam Crawley.

_______________________________________________
R-sig-mixed-models at r-project.org <mailto:R-sig-mixed-models at r-project.org>  mailing list
https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models

--

_____________________________________________________________________

Universit?tsklinikum Hamburg-Eppendorf; K?rperschaft des ?ffentlichen
Rechts; Gerichtsstand: Hamburg | www.uke.de <http://www.uke.de> 
Vorstandsmitglieder: Prof. Dr. Burkhard G?ke (Vorsitzender), Prof. Dr. Dr.
Uwe Koch-Gromus, Joachim Pr?l?, Marya Verdel
_____________________________________________________________________

SAVE PAPER - THINK BEFORE PRINTING

_______________________________________________
R-sig-mixed-models at r-project.org <mailto:R-sig-mixed-models at r-project.org>  mailing list
https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models

--

_____________________________________________________________________

Universit?tsklinikum Hamburg-Eppendorf; K?rperschaft des ?ffentlichen Rechts; Gerichtsstand: Hamburg | www.uke.de <http://www.uke.de> 
Vorstandsmitglieder: Prof. Dr. Burkhard G?ke (Vorsitzender), Prof. Dr. Dr. Uwe Koch-Gromus, Joachim Pr?l?, Marya Verdel
_____________________________________________________________________

SAVE PAPER - THINK BEFORE PRINTING

_______________________________________________
R-sig-mixed-models at r-project.org <mailto:R-sig-mixed-models at r-project.org>  mailing list
https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models

 

--

_____________________________________________________________________

Universit?tsklinikum Hamburg-Eppendorf; K?rperschaft des ?ffentlichen Rechts; Gerichtsstand: Hamburg | www.uke.de
Vorstandsmitglieder: Prof. Dr. Burkhard G?ke (Vorsitzender), Prof. Dr. Dr. Uwe Koch-Gromus, Joachim Pr?l?, Marya Verdel
_____________________________________________________________________

SAVE PAPER - THINK BEFORE PRINTING

  
  


  


      

    

    

      
      

        


  

        

      Mollie Brooks
    

    

      
      #
    

  


  

          

      
      
      
    

        
Yes, I (mebrooks) am the one who referenced it and the user (mmeierer) said it fit their needs for "sample weights".

            

      
      
      
    

        
Ok. I don?t know the survey package and don?t have time to look into it now. Feel free to follow up on issue 285 if you have more insight. 

cheers,
Mollie

            

      
      
      
    

        

      
      
    

            

      
      
      
    

        

  
  


  


      

    

    

      
      

        


  

        

      d@iuedecke m@iii@g oii uke@de
    

    

      
      #
    

  


  

      
Short addition, I also fitted two more glmmTMB-models, using sjstats::scale_weights(), which seems to come close to survey-weights (not sure if this is always the case that results match that good)

library(sjstats)
lalonde$id <- 1
lalonde <- scale_weights(lalonde, cluster.id = id, pweight = w)

glm3 <- glmmTMB(re78 ~ treat, weights=svywght_a , data=lalonde)
glm4 <- glmmTMB(re78 ~ treat, weights=svywght_b , data=lalonde)

# glm3 with svywght_a are similar to survey-weights
 head(cbind(glm1$weights, glm11$weights, glm2$frame$`(weights)`, glm3$frame$`(weights)`, glm4$frame$`(weights)`))

       [,1]     [,2]     [,3]      [,4]      [,5]
1 1.4682453 2.108394 2.108394 1.4682453 0.7302372
2 0.9593877 1.377677 1.377677 0.9593877 0.4771550
3 0.7489954 1.075554 1.075554 0.7489954 0.3725156
4 0.7319955 1.051143 1.051143 0.7319955 0.3640607
5 0.7283328 1.045883 1.045883 0.7283328 0.3622390
6 0.7244569 1.040317 1.040317 0.7244569 0.3603113


-----Urspr?ngliche Nachricht-----
Von: R-sig-mixed-models <r-sig-mixed-models-bounces at r-project.org> Im Auftrag von d.luedecke at uke.de
Gesendet: Montag, 10. Juni 2019 19:41
An: 'Mollie Brooks' <mollieebrooks at gmail.com>; 'Sam Crawley' <sam_crawley at warpmail.net>; 'Help Mixed Models' <r-sig-mixed-models at r-project.org>
Betreff: Re: [R-sig-ME] Predicted probabilites with CIs for multilevel logistic regression with prior weights

I think that Sam is talking about ?sampling? or ?survey? weights (as compared to analytical or frequency weights, used by ?normal? regression models).

 The issue you?re referring to is referenced by another issue (https://github.com/glmmTMB/glmmTMB/issues/440), which in turn shows an example from Cross Validated:
https://stats.stackexchange.com/questions/57107/use-of-weights-in-svyglm-vs-glm

 If I use that example, and add a third model fitted with glmmTMB, I get following result when comparing the weights from the fitted objects:

 library(glmmTMB)

glm2 <- glmmTMB(re78 ~ treat, weights = w , data = lalonde)
cbind(glm1$weights, glm11$weights, glm2$frame$`(weights)`)

#>        [,1]     [,2]     [,3]
#> 1 1.4682453 2.108394 2.108394
#> 2 0.9593877 1.377677 1.377677
#> 3 0.7489954 1.075554 1.075554
#> 4 0.7319955 1.051143 1.051143
#> 5 0.7283328 1.045883 1.045883
#> 6 0.7244569 1.040317 1.040317

As you can see, ?glm? and ?glmmTMB? produce the same weights, while the survey-package has different weights? I?m not sure that the weights implemented in glmmTMB are actually ?sampling? weights (for surveys, as implemented in the survey package), or how to reproduce such weights using glmmTMB.
 

Von: Mollie Brooks <mollieebrooks at gmail.com> 
Gesendet: Montag, 10. Juni 2019 19:04
An: Sam Crawley <sam_crawley at warpmail.net>; Help Mixed Models <r-sig-mixed-models at r-project.org>
Cc: d.luedecke at uke.de
Betreff: Re: [R-sig-ME] Predicted probabilites with CIs for multilevel logistic regression with prior weights

        
On 10Jun 2019, at 17:33, <d.luedecke at uke.de <mailto:d.luedecke at uke.de> > <d.luedecke at uke.de <mailto:d.luedecke at uke.de> > wrote:

        
mixed models in R do correctly not account for sampling weights


Should be: mixed models in R do *currently* not account for sampling weights

 

I?m still trying to get a handle of the different definitions of "weights" but I believe we implemented sampling weights in glmmTMB. We do this by weighting the log-likelihood contribution of each observation. I think this is different from prior weights if you mean Bayesian priors. There has been some discussion of the different implementations of "weights" in different R functions (link below) and we still need to update the documentation for glmmTMB 

https://github.com/glmmTMB/glmmTMB/issues/285

 

Here?s a binomial example:

 

library(glmmTMB)

set.seed(123)

n=100

dat=data.frame(trials=rpois(n, lambda=50), rownum=1:n)

dat$success=rbinom(n, dat$trials, prob=.3)

dat$rep=sample(1:5, size=n, replace=TRUE) #each observation is repeated 1 to 5 times

rows=rep(dat$rownum, each=1, times=dat$rep)

dat_disaggregated=dat[rows, ]





summary(glmmTMB(cbind(success, trials-success)~1, weights=rep, dat, family=binomial))

summary(glmmTMB(cbind(success, trials-success)~1, dat_disaggregated, family=binomial))





and it works with non-integer weights





summary(glmmTMB(cbind(success, trials-success)~1, weights=rep/5, dat, family=binomial))





cheers,

Mollie






-----Urspr?ngliche Nachricht-----
Von: R-sig-mixed-models <r-sig-mixed-models-bounces at r-project.org <mailto:r-sig-mixed-models-bounces at r-project.org> > Im
Auftrag von d.luedecke at uke.de <mailto:d.luedecke at uke.de> 
Gesendet: Montag, 10. Juni 2019 17:31
An: 'Sam Crawley' <sam_crawley at warpmail.net <mailto:sam_crawley at warpmail.net> >;
r-sig-mixed-models at r-project.org <mailto:r-sig-mixed-models at r-project.org> 
Betreff: Re: [R-sig-ME] Predicted probabilites with CIs for multilevel
logistic regression with prior weights

Hi Sam,

you could the "ggeffects" package
(https://strengejacke.github.io/ggeffects/), and there is also an example
for a logistic mixed effects model
(https://strengejacke.github.io/ggeffects/articles/practical_logisticmixedmo
del.html), which might help you.

For binomial models, using weights often results in the following warning:
#> non-integer #successes in a binomial glm!

However, CIs for the predicted probabilities can be calculated nevertheless
(at least in my quick example). Note that afaik, mixed models in R do
correctly not account for sampling weights. However, Thomas Lumley, author
of the survey-package, works on a survey-function for mixed models
(https://github.com/tslumley/svylme), probably the GitHub version is quite
stable (haven't tested yet).

An alternative would be the "scale_weights()" function from the
sjstats-package
(https://strengejacke.github.io/sjstats/articles/mixedmodels-statistics.html
#rescale-model-weights-for-complex-samples ), which rescales sampling
weights so they can be used as "weights" for the mixed models function you
have in R (lme4, lme, ...).

Based on that function, I have a small example that demonstrates how to
compute predicted probabilities for mixed models with (sampling) weights
(ignore the warnings, this is just for demonstration purposes):

library(lme4)
library(sjstats) # for scale_weights() and sample data
library(ggeffects) # for ggpredict()

data(nhanes_sample)
set.seed(123)
nhanes_sample$bin <- rbinom(nrow(nhanes_sample), 1, prob = .3)
nhanes_sample <- scale_weights(nhanes_sample, SDMVSTRA, WTINT2YR)

m <- glmer(
 bin ~ factor(RIAGENDR) * age + factor(RIDRETH1) + (1 | SDMVPSU),
 family = binomial(),
 data = nhanes_sample,
 weights = svywght_a
)

ggpredict(m, c("age", "RIAGENDR")) %>% plot()


Best
Daniel

-----Urspr?ngliche Nachricht-----
Von: R-sig-mixed-models <r-sig-mixed-models-bounces at r-project.org <mailto:r-sig-mixed-models-bounces at r-project.org> > Im
Auftrag von Sam Crawley
Gesendet: Montag, 10. Juni 2019 10:36
An: r-sig-mixed-models at r-project.org <mailto:r-sig-mixed-models at r-project.org> 
Betreff: [R-sig-ME] Predicted probabilites with CIs for multilevel logistic
regression with prior weights

Hello all,

I am doing a multilevel binomial logistic regression using lme4, and the
survey data I am using requires weights to be used. I would like to
calculate various predicted probabilities with confidence intervals based on
the estimated model. The predict function obviously doesn't give me standard
errors, and the recommended method to get these is to use the bootMer
function.

However, in my case, the weights provided are not integers, and the bootMer
function exits with an error if the weights are not integers (I raised a
GitHub issue about this, and was pointed to this list:
https://github.com/lme4/lme4/issues/524 ).

So my question is, what is the best way to calculate the predicted
probabilities (with confidence intervals) in my case?

I would appreciate any help you can give me, and I'm happy to provide more
details if required.

Thanks,
Sam Crawley.

_______________________________________________
R-sig-mixed-models at r-project.org <mailto:R-sig-mixed-models at r-project.org>  mailing list
https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models

--

_____________________________________________________________________

Universit?tsklinikum Hamburg-Eppendorf; K?rperschaft des ?ffentlichen
Rechts; Gerichtsstand: Hamburg | www.uke.de <http://www.uke.de> 
Vorstandsmitglieder: Prof. Dr. Burkhard G?ke (Vorsitzender), Prof. Dr. Dr.
Uwe Koch-Gromus, Joachim Pr?l?, Marya Verdel
_____________________________________________________________________

SAVE PAPER - THINK BEFORE PRINTING

_______________________________________________
R-sig-mixed-models at r-project.org <mailto:R-sig-mixed-models at r-project.org>  mailing list
https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models

--

_____________________________________________________________________

Universit?tsklinikum Hamburg-Eppendorf; K?rperschaft des ?ffentlichen Rechts; Gerichtsstand: Hamburg | www.uke.de <http://www.uke.de> 
Vorstandsmitglieder: Prof. Dr. Burkhard G?ke (Vorsitzender), Prof. Dr. Dr. Uwe Koch-Gromus, Joachim Pr?l?, Marya Verdel
_____________________________________________________________________

SAVE PAPER - THINK BEFORE PRINTING

_______________________________________________
R-sig-mixed-models at r-project.org <mailto:R-sig-mixed-models at r-project.org>  mailing list
https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models

 

--

_____________________________________________________________________

Universit?tsklinikum Hamburg-Eppendorf; K?rperschaft des ?ffentlichen Rechts; Gerichtsstand: Hamburg | www.uke.de
Vorstandsmitglieder: Prof. Dr. Burkhard G?ke (Vorsitzender), Prof. Dr. Dr. Uwe Koch-Gromus, Joachim Pr?l?, Marya Verdel
_____________________________________________________________________

SAVE PAPER - THINK BEFORE PRINTING


_______________________________________________
R-sig-mixed-models at r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models

--

_____________________________________________________________________

Universit?tsklinikum Hamburg-Eppendorf; K?rperschaft des ?ffentlichen Rechts; Gerichtsstand: Hamburg | www.uke.de
Vorstandsmitglieder: Prof. Dr. Burkhard G?ke (Vorsitzender), Prof. Dr. Dr. Uwe Koch-Gromus, Joachim Pr?l?, Marya Verdel
_____________________________________________________________________

SAVE PAPER - THINK BEFORE PRINTING

  
  


  


      

    

    

      
      

        


  

        

      d@iuedecke m@iii@g oii uke@de
    

    

      
      #
    

  


  

          

      
      
      
    

        
At least from the technical side, I don?t have more insights, I guess. I already noticed the discussion in #285 some time ago, so I?m lurking, but not actively following ?

 

Terminology used in different papers or from method reports of different surveys also doesn?t seem always consistent to me. I think, ?post-stratification weights? were requested by Sam, which are weights based on group (or stratum) characteristics (like the distribution of age or gender proportions). Ben Bolker also mentioned sample weights in #285 (?in addition to these two cases, there's also the case of sampling weights, which is difficult/a mess for complex regression models but worth discussing at least ...?).

 

The difference between the weights-argument in typical regression model functions and the survey-package is ?The survey package not only allows for adjusting the composition of a sample to the characteristics of the general population. Most base packages would allow you to do that by specifying a weights argument. The survey package goes further by correcting the design effect introduced by the application of post-stratification weights.? (https://tophcito.blogspot.com/2014/04/social-science-goes-r-weighted-survey.html).

 

This of course only applies if you actually have survey-data.

 

Von: Mollie Brooks <mollieebrooks at gmail.com> 
Gesendet: Montag, 10. Juni 2019 19:46
An: d.luedecke at uke.de
Cc: Sam Crawley <sam_crawley at warpmail.net>; Help Mixed Models <r-sig-mixed-models at r-project.org>
Betreff: Re: [R-sig-ME] Predicted probabilites with CIs for multilevel logistic regression with prior weights

        
On 10Jun 2019, at 19:40, <d.luedecke at uke.de <mailto:d.luedecke at uke.de> > <d.luedecke at uke.de <mailto:d.luedecke at uke.de> > wrote:

        
I think that Sam is talking about ?sampling? or ?survey? weights (as compared to analytical or frequency weights, used by ?normal? regression models).

 

The issue you?re referring to is referenced by another issue ( <https://github.com/glmmTMB/glmmTMB/issues/440> https://github.com/glmmTMB/glmmTMB/issues/440), 

 

Yes, I (mebrooks) am the one who referenced it and the user (mmeierer) said it fit their needs for "sample weights".





which in turn shows an example from Cross Validated:

 <https://stats.stackexchange.com/questions/57107/use-of-weights-in-svyglm-vs-glm> https://stats.stackexchange.com/questions/57107/use-of-weights-in-svyglm-vs-glm

 

If I use that example, and add a third model fitted with glmmTMB, I get following result when comparing the weights from the fitted objects:

 

library(glmmTMB)

glm2 <- glmmTMB(re78 ~ treat, weights = w , data = lalonde)

cbind(glm1$weights, glm11$weights, glm2$frame$`(weights)`)

#>        [,1]     [,2]     [,3]

#> 1 1.4682453 2.108394 2.108394

#> 2 0.9593877 1.377677 1.377677

#> 3 0.7489954 1.075554 1.075554

#> 4 0.7319955 1.051143 1.051143

#> 5 0.7283328 1.045883 1.045883

#> 6 0.7244569 1.040317 1.040317

 

As you can see, ?glm? and ?glmmTMB? produce the same weights, while the survey-package has different weights? I?m not sure that the weights implemented in glmmTMB are actually ?sampling? weights (for surveys, as implemented in the survey package), 

 

Ok. I don?t know the survey package and don?t have time to look into it now. Feel free to follow up on issue 285 if you have more insight. 

 

cheers,

Mollie





or how to reproduce such weights using glmmTMB.

 

Von: Mollie Brooks <mollieebrooks at gmail.com <mailto:mollieebrooks at gmail.com> > 
Gesendet: Montag, 10. Juni 2019 19:04
An: Sam Crawley <sam_crawley at warpmail.net <mailto:sam_crawley at warpmail.net> >; Help Mixed Models <r-sig-mixed-models at r-project.org <mailto:r-sig-mixed-models at r-project.org> >
Cc: d.luedecke at uke.de <mailto:d.luedecke at uke.de> 
Betreff: Re: [R-sig-ME] Predicted probabilites with CIs for multilevel logistic regression with prior weights

        
On 10Jun 2019, at 17:33, < <mailto:d.luedecke at uke.de> d.luedecke at uke.de> < <mailto:d.luedecke at uke.de> d.luedecke at uke.de> wrote:

        
mixed models in R do correctly not account for sampling weights


Should be: mixed models in R do *currently* not account for sampling weights

 

I?m still trying to get a handle of the different definitions of "weights" but I believe we implemented sampling weights in glmmTMB. We do this by weighting the log-likelihood contribution of each observation. I think this is different from prior weights if you mean Bayesian priors. There has been some discussion of the different implementations of "weights" in different R functions (link below) and we still need to update the documentation for glmmTMB 

 <https://github.com/glmmTMB/glmmTMB/issues/285> https://github.com/glmmTMB/glmmTMB/issues/285

 

Here?s a binomial example:

 

library(glmmTMB)

set.seed(123)

n=100

dat=data.frame(trials=rpois(n, lambda=50), rownum=1:n)

dat$success=rbinom(n, dat$trials, prob=.3)

dat$rep=sample(1:5, size=n, replace=TRUE) #each observation is repeated 1 to 5 times

rows=rep(dat$rownum, each=1, times=dat$rep)

dat_disaggregated=dat[rows, ]






summary(glmmTMB(cbind(success, trials-success)~1, weights=rep, dat, family=binomial))

summary(glmmTMB(cbind(success, trials-success)~1, dat_disaggregated, family=binomial))






and it works with non-integer weights






summary(glmmTMB(cbind(success, trials-success)~1, weights=rep/5, dat, family=binomial))






cheers,

Mollie







-----Urspr?ngliche Nachricht-----
Von: R-sig-mixed-models < <mailto:r-sig-mixed-models-bounces at r-project.org> r-sig-mixed-models-bounces at r-project.org> Im
Auftrag von  <mailto:d.luedecke at uke.de> d.luedecke at uke.de
Gesendet: Montag, 10. Juni 2019 17:31
An: 'Sam Crawley' < <mailto:sam_crawley at warpmail.net> sam_crawley at warpmail.net>;
 <mailto:r-sig-mixed-models at r-project.org> r-sig-mixed-models at r-project.org
Betreff: Re: [R-sig-ME] Predicted probabilites with CIs for multilevel
logistic regression with prior weights

Hi Sam,

you could the "ggeffects" package
( <https://strengejacke.github.io/ggeffects/> https://strengejacke.github.io/ggeffects/), and there is also an example
for a logistic mixed effects model
( <https://strengejacke.github.io/ggeffects/articles/practical_logisticmixedmo> https://strengejacke.github.io/ggeffects/articles/practical_logisticmixedmo
del.html), which might help you.

For binomial models, using weights often results in the following warning:
#> non-integer #successes in a binomial glm!

However, CIs for the predicted probabilities can be calculated nevertheless
(at least in my quick example). Note that afaik, mixed models in R do
correctly not account for sampling weights. However, Thomas Lumley, author
of the survey-package, works on a survey-function for mixed models
( <https://github.com/tslumley/svylme> https://github.com/tslumley/svylme), probably the GitHub version is quite
stable (haven't tested yet).

An alternative would be the "scale_weights()" function from the
sjstats-package
( <https://strengejacke.github.io/sjstats/articles/mixedmodels-statistics.html> https://strengejacke.github.io/sjstats/articles/mixedmodels-statistics.html
#rescale-model-weights-for-complex-samples ), which rescales sampling
weights so they can be used as "weights" for the mixed models function you
have in R (lme4, lme, ...).

Based on that function, I have a small example that demonstrates how to
compute predicted probabilities for mixed models with (sampling) weights
(ignore the warnings, this is just for demonstration purposes):

library(lme4)
library(sjstats) # for scale_weights() and sample data
library(ggeffects) # for ggpredict()

data(nhanes_sample)
set.seed(123)
nhanes_sample$bin <- rbinom(nrow(nhanes_sample), 1, prob = .3)
nhanes_sample <- scale_weights(nhanes_sample, SDMVSTRA, WTINT2YR)

m <- glmer(
 bin ~ factor(RIAGENDR) * age + factor(RIDRETH1) + (1 | SDMVPSU),
 family = binomial(),
 data = nhanes_sample,
 weights = svywght_a
)

ggpredict(m, c("age", "RIAGENDR")) %>% plot()


Best
Daniel

-----Urspr?ngliche Nachricht-----
Von: R-sig-mixed-models < <mailto:r-sig-mixed-models-bounces at r-project.org> r-sig-mixed-models-bounces at r-project.org> Im
Auftrag von Sam Crawley
Gesendet: Montag, 10. Juni 2019 10:36
An:  <mailto:r-sig-mixed-models at r-project.org> r-sig-mixed-models at r-project.org
Betreff: [R-sig-ME] Predicted probabilites with CIs for multilevel logistic
regression with prior weights

Hello all,

I am doing a multilevel binomial logistic regression using lme4, and the
survey data I am using requires weights to be used. I would like to
calculate various predicted probabilities with confidence intervals based on
the estimated model. The predict function obviously doesn't give me standard
errors, and the recommended method to get these is to use the bootMer
function.

However, in my case, the weights provided are not integers, and the bootMer
function exits with an error if the weights are not integers (I raised a
GitHub issue about this, and was pointed to this list:
 <https://github.com/lme4/lme4/issues/524> https://github.com/lme4/lme4/issues/524 ).

So my question is, what is the best way to calculate the predicted
probabilities (with confidence intervals) in my case?

I would appreciate any help you can give me, and I'm happy to provide more
details if required.

Thanks,
Sam Crawley.

_______________________________________________
 <mailto:R-sig-mixed-models at r-project.org> R-sig-mixed-models at r-project.org mailing list
 <https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models

--

_____________________________________________________________________

Universit?tsklinikum Hamburg-Eppendorf; K?rperschaft des ?ffentlichen
Rechts; Gerichtsstand: Hamburg |  <http://www.uke.de/> www.uke.de
Vorstandsmitglieder: Prof. Dr. Burkhard G?ke (Vorsitzender), Prof. Dr. Dr.
Uwe Koch-Gromus, Joachim Pr?l?, Marya Verdel
_____________________________________________________________________

SAVE PAPER - THINK BEFORE PRINTING

_______________________________________________
 <mailto:R-sig-mixed-models at r-project.org> R-sig-mixed-models at r-project.org mailing list
 <https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models

--

_____________________________________________________________________

Universit?tsklinikum Hamburg-Eppendorf; K?rperschaft des ?ffentlichen Rechts; Gerichtsstand: Hamburg |  <http://www.uke.de/> www.uke.de
Vorstandsmitglieder: Prof. Dr. Burkhard G?ke (Vorsitzender), Prof. Dr. Dr. Uwe Koch-Gromus, Joachim Pr?l?, Marya Verdel
_____________________________________________________________________

SAVE PAPER - THINK BEFORE PRINTING

_______________________________________________
 <mailto:R-sig-mixed-models at r-project.org> R-sig-mixed-models at r-project.org mailing list
 <https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models

 

 


  _____  


Universit?tsklinikum Hamburg-Eppendorf; K?rperschaft des ?ffentlichen Rechts; Gerichtsstand: Hamburg |  <http://www.uke.de/> www.uke.de
Vorstandsmitglieder: Prof. Dr. Burkhard G?ke (Vorsitzender), Prof. Dr. Dr. Uwe Koch-Gromus, Joachim Pr?l?, Marya Verdel


  _____  


SAVE PAPER - THINK BEFORE PRINTING

 

--

_____________________________________________________________________

Universit?tsklinikum Hamburg-Eppendorf; K?rperschaft des ?ffentlichen Rechts; Gerichtsstand: Hamburg | www.uke.de
Vorstandsmitglieder: Prof. Dr. Burkhard G?ke (Vorsitzender), Prof. Dr. Dr. Uwe Koch-Gromus, Joachim Pr?l?, Marya Verdel
_____________________________________________________________________

SAVE PAPER - THINK BEFORE PRINTING

  
  


  


      

    

    

      
      

        


  

        

      Rose, Charles E. (CDC/DDNID/NCBDDD/OD)
    

    

      
      #
    

  


  

      
I don?t know the package glmmTMB but I'm wondering if the use of survey weights in the weight statement results in the appropriate SEs as one would expect from the Lumley survey package; any comments about this issue are appreciated, thanks, chuck


-----Original Message-----
From: R-sig-mixed-models <r-sig-mixed-models-bounces at r-project.org> On Behalf Of d.luedecke at uke.de
Sent: Monday, June 10, 2019 1:47 PM
To: 'Mollie Brooks' <mollieebrooks at gmail.com>; 'Sam Crawley' <sam_crawley at warpmail.net>; 'Help Mixed Models' <r-sig-mixed-models at r-project.org>
Subject: Re: [R-sig-ME] Predicted probabilites with CIs for multilevel logistic regression with prior weights

Short addition, I also fitted two more glmmTMB-models, using sjstats::scale_weights(), which seems to come close to survey-weights (not sure if this is always the case that results match that good)

library(sjstats)
lalonde$id <- 1
lalonde <- scale_weights(lalonde, cluster.id = id, pweight = w)

glm3 <- glmmTMB(re78 ~ treat, weights=svywght_a , data=lalonde)
glm4 <- glmmTMB(re78 ~ treat, weights=svywght_b , data=lalonde)

# glm3 with svywght_a are similar to survey-weights  head(cbind(glm1$weights, glm11$weights, glm2$frame$`(weights)`, glm3$frame$`(weights)`, glm4$frame$`(weights)`))

       [,1]     [,2]     [,3]      [,4]      [,5]
1 1.4682453 2.108394 2.108394 1.4682453 0.7302372
2 0.9593877 1.377677 1.377677 0.9593877 0.4771550
3 0.7489954 1.075554 1.075554 0.7489954 0.3725156
4 0.7319955 1.051143 1.051143 0.7319955 0.3640607
5 0.7283328 1.045883 1.045883 0.7283328 0.3622390
6 0.7244569 1.040317 1.040317 0.7244569 0.3603113


-----Urspr?ngliche Nachricht-----
Von: R-sig-mixed-models <r-sig-mixed-models-bounces at r-project.org> Im Auftrag von d.luedecke at uke.de
Gesendet: Montag, 10. Juni 2019 19:41
An: 'Mollie Brooks' <mollieebrooks at gmail.com>; 'Sam Crawley' <sam_crawley at warpmail.net>; 'Help Mixed Models' <r-sig-mixed-models at r-project.org>
Betreff: Re: [R-sig-ME] Predicted probabilites with CIs for multilevel logistic regression with prior weights

I think that Sam is talking about ?sampling? or ?survey? weights (as compared to analytical or frequency weights, used by ?normal? regression models).

 The issue you?re referring to is referenced by another issue (https://github.com/glmmTMB/glmmTMB/issues/440), which in turn shows an example from Cross Validated:
https://stats.stackexchange.com/questions/57107/use-of-weights-in-svyglm-vs-glm

 If I use that example, and add a third model fitted with glmmTMB, I get following result when comparing the weights from the fitted objects:

 library(glmmTMB)

glm2 <- glmmTMB(re78 ~ treat, weights = w , data = lalonde) cbind(glm1$weights, glm11$weights, glm2$frame$`(weights)`)

#>        [,1]     [,2]     [,3]
#> 1 1.4682453 2.108394 2.108394
#> 2 0.9593877 1.377677 1.377677
#> 3 0.7489954 1.075554 1.075554
#> 4 0.7319955 1.051143 1.051143
#> 5 0.7283328 1.045883 1.045883
#> 6 0.7244569 1.040317 1.040317

As you can see, ?glm? and ?glmmTMB? produce the same weights, while the survey-package has different weights? I?m not sure that the weights implemented in glmmTMB are actually ?sampling? weights (for surveys, as implemented in the survey package), or how to reproduce such weights using glmmTMB.
 

Von: Mollie Brooks <mollieebrooks at gmail.com>
Gesendet: Montag, 10. Juni 2019 19:04
An: Sam Crawley <sam_crawley at warpmail.net>; Help Mixed Models <r-sig-mixed-models at r-project.org>
Cc: d.luedecke at uke.de
Betreff: Re: [R-sig-ME] Predicted probabilites with CIs for multilevel logistic regression with prior weights

        
On 10Jun 2019, at 17:33, <d.luedecke at uke.de <mailto:d.luedecke at uke.de> > <d.luedecke at uke.de <mailto:d.luedecke at uke.de> > wrote:

        
mixed models in R do correctly not account for sampling weights


Should be: mixed models in R do *currently* not account for sampling weights

 

I?m still trying to get a handle of the different definitions of "weights" but I believe we implemented sampling weights in glmmTMB. We do this by weighting the log-likelihood contribution of each observation. I think this is different from prior weights if you mean Bayesian priors. There has been some discussion of the different implementations of "weights" in different R functions (link below) and we still need to update the documentation for glmmTMB 

https://github.com/glmmTMB/glmmTMB/issues/285

 

Here?s a binomial example:

 

library(glmmTMB)

set.seed(123)

n=100

dat=data.frame(trials=rpois(n, lambda=50), rownum=1:n)

dat$success=rbinom(n, dat$trials, prob=.3)

dat$rep=sample(1:5, size=n, replace=TRUE) #each observation is repeated 1 to 5 times

rows=rep(dat$rownum, each=1, times=dat$rep)

dat_disaggregated=dat[rows, ]





summary(glmmTMB(cbind(success, trials-success)~1, weights=rep, dat, family=binomial))

summary(glmmTMB(cbind(success, trials-success)~1, dat_disaggregated, family=binomial))





and it works with non-integer weights





summary(glmmTMB(cbind(success, trials-success)~1, weights=rep/5, dat, family=binomial))





cheers,

Mollie






-----Urspr?ngliche Nachricht-----
Von: R-sig-mixed-models <r-sig-mixed-models-bounces at r-project.org <mailto:r-sig-mixed-models-bounces at r-project.org> > Im Auftrag von d.luedecke at uke.de <mailto:d.luedecke at uke.de>
Gesendet: Montag, 10. Juni 2019 17:31
An: 'Sam Crawley' <sam_crawley at warpmail.net <mailto:sam_crawley at warpmail.net> >; r-sig-mixed-models at r-project.org <mailto:r-sig-mixed-models at r-project.org>
Betreff: Re: [R-sig-ME] Predicted probabilites with CIs for multilevel logistic regression with prior weights

Hi Sam,

you could the "ggeffects" package
(https://strengejacke.github.io/ggeffects/), and there is also an example for a logistic mixed effects model (https://strengejacke.github.io/ggeffects/articles/practical_logisticmixedmo
del.html), which might help you.

For binomial models, using weights often results in the following warning:
#> non-integer #successes in a binomial glm!

However, CIs for the predicted probabilities can be calculated nevertheless (at least in my quick example). Note that afaik, mixed models in R do correctly not account for sampling weights. However, Thomas Lumley, author of the survey-package, works on a survey-function for mixed models (https://github.com/tslumley/svylme), probably the GitHub version is quite stable (haven't tested yet).

An alternative would be the "scale_weights()" function from the sjstats-package (https://strengejacke.github.io/sjstats/articles/mixedmodels-statistics.html
#rescale-model-weights-for-complex-samples ), which rescales sampling weights so they can be used as "weights" for the mixed models function you have in R (lme4, lme, ...).

Based on that function, I have a small example that demonstrates how to compute predicted probabilities for mixed models with (sampling) weights (ignore the warnings, this is just for demonstration purposes):

library(lme4)
library(sjstats) # for scale_weights() and sample data
library(ggeffects) # for ggpredict()

data(nhanes_sample)
set.seed(123)
nhanes_sample$bin <- rbinom(nrow(nhanes_sample), 1, prob = .3) nhanes_sample <- scale_weights(nhanes_sample, SDMVSTRA, WTINT2YR)

m <- glmer(
 bin ~ factor(RIAGENDR) * age + factor(RIDRETH1) + (1 | SDMVPSU),  family = binomial(),  data = nhanes_sample,  weights = svywght_a
)

ggpredict(m, c("age", "RIAGENDR")) %>% plot()


Best
Daniel

-----Urspr?ngliche Nachricht-----
Von: R-sig-mixed-models <r-sig-mixed-models-bounces at r-project.org <mailto:r-sig-mixed-models-bounces at r-project.org> > Im Auftrag von Sam Crawley
Gesendet: Montag, 10. Juni 2019 10:36
An: r-sig-mixed-models at r-project.org <mailto:r-sig-mixed-models at r-project.org>
Betreff: [R-sig-ME] Predicted probabilites with CIs for multilevel logistic regression with prior weights

Hello all,

I am doing a multilevel binomial logistic regression using lme4, and the survey data I am using requires weights to be used. I would like to calculate various predicted probabilities with confidence intervals based on the estimated model. The predict function obviously doesn't give me standard errors, and the recommended method to get these is to use the bootMer function.

However, in my case, the weights provided are not integers, and the bootMer function exits with an error if the weights are not integers (I raised a GitHub issue about this, and was pointed to this list:
https://github.com/lme4/lme4/issues/524 ).

So my question is, what is the best way to calculate the predicted probabilities (with confidence intervals) in my case?

I would appreciate any help you can give me, and I'm happy to provide more details if required.

Thanks,
Sam Crawley.

_______________________________________________
R-sig-mixed-models at r-project.org <mailto:R-sig-mixed-models at r-project.org>  mailing list https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models

--

_____________________________________________________________________

Universit?tsklinikum Hamburg-Eppendorf; K?rperschaft des ?ffentlichen Rechts; Gerichtsstand: Hamburg | www.uke.de <http://www.uke.de>
Vorstandsmitglieder: Prof. Dr. Burkhard G?ke (Vorsitzender), Prof. Dr. Dr.
Uwe Koch-Gromus, Joachim Pr?l?, Marya Verdel _____________________________________________________________________

SAVE PAPER - THINK BEFORE PRINTING

_______________________________________________
R-sig-mixed-models at r-project.org <mailto:R-sig-mixed-models at r-project.org>  mailing list https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models

--

_____________________________________________________________________

Universit?tsklinikum Hamburg-Eppendorf; K?rperschaft des ?ffentlichen Rechts; Gerichtsstand: Hamburg | www.uke.de <http://www.uke.de>
Vorstandsmitglieder: Prof. Dr. Burkhard G?ke (Vorsitzender), Prof. Dr. Dr. Uwe Koch-Gromus, Joachim Pr?l?, Marya Verdel _____________________________________________________________________

SAVE PAPER - THINK BEFORE PRINTING

_______________________________________________
R-sig-mixed-models at r-project.org <mailto:R-sig-mixed-models at r-project.org>  mailing list https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models

 

--

_____________________________________________________________________

Universit?tsklinikum Hamburg-Eppendorf; K?rperschaft des ?ffentlichen Rechts; Gerichtsstand: Hamburg | www.uke.de
Vorstandsmitglieder: Prof. Dr. Burkhard G?ke (Vorsitzender), Prof. Dr. Dr. Uwe Koch-Gromus, Joachim Pr?l?, Marya Verdel _____________________________________________________________________

SAVE PAPER - THINK BEFORE PRINTING


_______________________________________________
R-sig-mixed-models at r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models

--

_____________________________________________________________________

Universit?tsklinikum Hamburg-Eppendorf; K?rperschaft des ?ffentlichen Rechts; Gerichtsstand: Hamburg | www.uke.de
Vorstandsmitglieder: Prof. Dr. Burkhard G?ke (Vorsitzender), Prof. Dr. Dr. Uwe Koch-Gromus, Joachim Pr?l?, Marya Verdel _____________________________________________________________________

SAVE PAPER - THINK BEFORE PRINTING
_______________________________________________
R-sig-mixed-models at r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models

  
  


  


      

    

    

      
      

        


  

        

      Sam Crawley
    

    

      
      #
    

  


  

      
Yes, I am looking to apply post-stratification weights (the dataset I am using is from Eurobarometer, which provides the weights).

Thanks very much all for your help so far. It looks like I have a bit of reading to do, so I will continue to work on this and let you know if I have further questions.

Cheers,
Sam Crawley.

        
On Tue, 11 Jun 2019, at 06:15, d.luedecke at uke.de wrote:

            

      
      
      
    

        

      
      
    

        

            

      
      
      
    

        

  
  


  


      

    

    
7 days later

    

      
      

        


  

        

      Sam Crawley
    

    

      
      #
    

  


  

      
Hi again Daniel (and list),

Thanks again for the below. I have been using the ggpredict() function, and it works well. However, should I be using the type = "re" parameter? Or is this only required when attempting to predict values for each group? (I have read the ggeffects documentation on this, but it's still not entirely clear to me).

When adding type="re", the confidence intervals become very wide, which is obviously not ideal.

Thanks,
Sam Crawley.

        
On Tue, 11 Jun 2019, at 03:30, d.luedecke at uke.de wrote:

            

      
      
      
    

        

      
      
    

  
  


  


      

    

    

      
      

        


  

        

      d@iuedecke m@iii@g oii uke@de
    

    

      
      #
    

  


  

          

      
      
      
    

        
This depends on whether you would like to have confidence or prediction intervals. I'm not sure if there's a clear definition, or a consensus on how to obtain prediction intervals. There's a section about this in Ben Bolker's GLMM-FAQ:
http://bbolker.github.io/mixedmodels-misc/glmmFAQ.html#predictions-andor-confidence-or-prediction-intervals-on-predictions

He states that you want to add the residual variance to compute prediction intervals. However, there are also some caveats: https://stackoverflow.com/questions/14358811/extract-prediction-band-from-lme-fit

Ben also has a comment in his examples: "## must be adapted for more complex models". I'm not sure if this refers only to make sure to compute sigma properly, or if you also need to take the random effect variances into account. ggeffects currently computes the standard error for type = "re" based on the *random effects variances* (see details here: https://easystats.github.io/insight/reference/get_variance.html). This seems to be even a bit more conservative than just taking the residual variance. But again, any comments on this are welcome. I've been in email exchange with Russell Lenth, maintainer of the emmeans-package, and we were also talking about if it's actually straightforward to obtain prediction intervals, or not. In the latest emmeans-release, there is a vignette on this topic (https://cran.r-project.org/web/packages/emmeans/vignettes/predictions.html).

In short: you don't need to specify type = "re" if you want to get predictions from your mixed models. Predictions for both type = "fe" and type = "re" are always on a population-level. However, you need type = "re" if you want to make predictions at each level of the random effects (group factor), see https://strengejacke.github.io/ggeffects/articles/randomeffects.html#marginal-effects-for-each-level-of-random-effects (though here you don't get any intervals).

Best
Daniel

-----Urspr?ngliche Nachricht-----
Von: Sam Crawley <sam_crawley at warpmail.net> 
Gesendet: Dienstag, 18. Juni 2019 03:44
An: d.luedecke <d.luedecke at uke.de>; r-sig-mixed-models at r-project.org
Betreff: Re: [R-sig-ME] Predicted probabilites with CIs for multilevel logistic regression with prior weights

Hi again Daniel (and list),

Thanks again for the below. I have been using the ggpredict() function, and it works well. However, should I be using the type = "re" parameter? Or is this only required when attempting to predict values for each group? (I have read the ggeffects documentation on this, but it's still not entirely clear to me).

When adding type="re", the confidence intervals become very wide, which is obviously not ideal.

Thanks,
Sam Crawley.

        
On Tue, 11 Jun 2019, at 03:30, d.luedecke at uke.de wrote:

            

      
      
      
    

        

      
      
    

        
--

_____________________________________________________________________

Universit?tsklinikum Hamburg-Eppendorf; K?rperschaft des ?ffentlichen Rechts; Gerichtsstand: Hamburg | www.uke.de
Vorstandsmitglieder: Prof. Dr. Burkhard G?ke (Vorsitzender), Prof. Dr. Dr. Uwe Koch-Gromus, Joachim Pr?l?, Marya Verdel
_____________________________________________________________________

SAVE PAPER - THINK BEFORE PRINTING

  
  


  


      

    

    

      
      

        


  

        

      Sam Crawley
    

    

      
      #
    

  


  

      
Hi Daniel,

Thanks very much, that clears it up.

Cheers,
Sam.

        
On Tue, 18 Jun 2019, at 18:30, d.luedecke at uke.de wrote: