State Eviction Policy Models 2

Author

Brian Surratt

Published

March 21, 2025

Reading in the data

Late renter households.

Code
latehh <- read.csv("./data/latehh.csv")

nrow(latehh)
[1] 16660
Code
latehh <- latehh %>%
  mutate(inclvl = case_when(income == 1 ~ "lessthan25",
                            income %in% c(2,3) ~ "25to49",
                            TRUE ~ "50andover"
                            )
         )

latehh$inclvl <- factor(latehh$inclvl, levels=c("25to49",  "lessthan25",  "50andover"))

tabyl(latehh$inclvl)
 latehh$inclvl    n   percent
        25to49 6158 0.3696279
    lessthan25 6894 0.4138055
     50andover 3608 0.2165666

Low late renter hh

Code
low_latehh <- read.csv("./data/low_latehh.csv")

nrow(low_latehh)
[1] 5136

Med late renter hh

Code
med_latehh <- read.csv("./data/medium_latehh.csv")

nrow(med_latehh)
[1] 5713

High late renter hh

Code
hi_latehh <- read.csv("./data/high_latehh.csv")

nrow(hi_latehh)
[1] 5811

Recode income as a categorical variable in all four dfs

Code
tabyl(latehh$income)
 latehh$income    n     percent
             1 6894 0.413805522
             2 3551 0.213145258
             3 2607 0.156482593
             4 1999 0.119987995
             5  829 0.049759904
             6  521 0.031272509
             7  126 0.007563025
             8  133 0.007983193
Code
nrow(latehh)
[1] 16660
Code
tabyl(low_latehh$income)
 low_latehh$income    n     percent
                 1 2353 0.458138629
                 2 1117 0.217484424
                 3  766 0.149143302
                 4  556 0.108255452
                 5  209 0.040693146
                 6  100 0.019470405
                 7   16 0.003115265
                 8   19 0.003699377
Code
nrow(low_latehh)
[1] 5136
Code
low_latehh <- low_latehh %>%
  mutate(inclvl = case_when(income == 1 ~ "lessthan25",
                            income %in% c(2,3) ~ "25to49",
                            TRUE ~ "50andover"
                            )
         )

low_latehh$inclvl <- factor(low_latehh$inclvl, levels=c("25to49",  "lessthan25",  "50andover"))

tabyl(low_latehh$inclvl)
 low_latehh$inclvl    n   percent
            25to49 1883 0.3666277
        lessthan25 2353 0.4581386
         50andover  900 0.1752336
Code
med_latehh <- med_latehh %>%
  mutate(inclvl = case_when(income == 1 ~ "lessthan25",
                            income %in% c(2,3) ~ "25to49",
                            TRUE ~ "50andover"
                            )
         )

med_latehh$inclvl <- factor(med_latehh$inclvl, levels=c("25to49",  "lessthan25",  "50andover"))

tabyl(med_latehh$inclvl)
 med_latehh$inclvl    n   percent
            25to49 2114 0.3700333
        lessthan25 2252 0.3941887
         50andover 1347 0.2357781
Code
hi_latehh <- hi_latehh %>%
  mutate(inclvl = case_when(income == 1 ~ "lessthan25",
                            income %in% c(2,3) ~ "25to49",
                            TRUE ~ "50andover"
                            )
         )

hi_latehh$inclvl <- factor(hi_latehh$inclvl, levels=c("25to49",  "lessthan25",  "50andover"))

tabyl(hi_latehh$inclvl)
 hi_latehh$inclvl    n   percent
           25to49 2161 0.3718809
       lessthan25 2289 0.3939081
        50andover 1361 0.2342110

Recode education

Code
tabyl(latehh$eeduc)
 latehh$eeduc    n    percent
            1  378 0.02268908
            2  733 0.04399760
            3 3548 0.21296519
            4 5230 0.31392557
            5 2126 0.12761104
            6 2770 0.16626651
            7 1875 0.11254502
Code
latehh <- latehh %>%
  mutate(educ = case_when(eeduc %in% c(1:3) ~ "HSorless",
                            eeduc == 4 ~ "somecollege",
                            TRUE ~ "BAorhigher"
                            )
         )

latehh$educ <- factor(latehh$educ, levels=c("somecollege",  "HSorless",  "BAorhigher"))

tabyl(latehh$educ)
 latehh$educ    n   percent
 somecollege 5230 0.3139256
    HSorless 4659 0.2796519
  BAorhigher 6771 0.4064226
Code
tabyl(low_latehh$eeduc)
 low_latehh$eeduc    n    percent
                1  114 0.02219626
                2  234 0.04556075
                3 1156 0.22507788
                4 1754 0.34151090
                5  738 0.14369159
                6  699 0.13609813
                7  441 0.08586449
Code
low_latehh <- low_latehh %>%
  mutate(educ = case_when(eeduc %in% c(1:3) ~ "HSorless",
                            eeduc == 4 ~ "somecollege",
                            TRUE ~ "BAorhigher"
                            )
         )

low_latehh$educ <- factor(low_latehh$educ, levels=c("somecollege",  "HSorless",  "BAorhigher"))

tabyl(low_latehh$educ)
 low_latehh$educ    n   percent
     somecollege 1754 0.3415109
        HSorless 1504 0.2928349
      BAorhigher 1878 0.3656542
Code
med_latehh <- med_latehh %>%
  mutate(educ = case_when(eeduc %in% c(1:3) ~ "HSorless",
                            eeduc == 4 ~ "somecollege",
                            TRUE ~ "BAorhigher"
                            )
         )

med_latehh$educ <- factor(med_latehh$educ, levels=c("somecollege",  "HSorless",  "BAorhigher"))

tabyl(med_latehh$educ)
 med_latehh$educ    n  percent
     somecollege 1701 0.297742
        HSorless 1550 0.271311
      BAorhigher 2462 0.430947
Code
hi_latehh <- hi_latehh %>%
  mutate(educ = case_when(eeduc %in% c(1:3) ~ "HSorless",
                            eeduc == 4 ~ "somecollege",
                            TRUE ~ "BAorhigher"
                            )
         )

hi_latehh$educ <- factor(hi_latehh$educ, levels=c("somecollege",  "HSorless",  "BAorhigher"))

tabyl(hi_latehh$educ)
 hi_latehh$educ    n   percent
    somecollege 1775 0.3054552
       HSorless 1605 0.2762003
     BAorhigher 2431 0.4183445

Relevel race_eth

Code
tabyl(latehh$race_eth)
 latehh$race_eth    n    percent
        hispanic 2960 0.17767107
        nh_asian 1085 0.06512605
        nh_black 3854 0.23133253
        nh_white 7599 0.45612245
           other 1162 0.06974790
Code
class(latehh$race_eth)
[1] "character"
Code
latehh$race_eth <- factor(latehh$race_eth, levels = c("nh_white", "nh_black", "hispanic", "other", "nh_asian"))

tabyl(latehh$race_eth)
 latehh$race_eth    n    percent
        nh_white 7599 0.45612245
        nh_black 3854 0.23133253
        hispanic 2960 0.17767107
           other 1162 0.06974790
        nh_asian 1085 0.06512605
Code
low_latehh$race_eth <- factor(low_latehh$race_eth, levels = c("nh_white", "nh_black", "hispanic", "other", "nh_asian"))

tabyl(low_latehh$race_eth)
 low_latehh$race_eth    n    percent
            nh_white 2611 0.50837227
            nh_black 1292 0.25155763
            hispanic  601 0.11701713
               other  452 0.08800623
            nh_asian  180 0.03504673
Code
med_latehh$race_eth <- factor(med_latehh$race_eth, levels = c("nh_white", "nh_black", "hispanic", "other", "nh_asian"))

tabyl(med_latehh$race_eth)
 med_latehh$race_eth    n   percent
            nh_white 2243 0.3926133
            nh_black 1316 0.2303518
            hispanic 1367 0.2392788
               other  297 0.0519867
            nh_asian  490 0.0857693
Code
hi_latehh$race_eth <- factor(hi_latehh$race_eth, levels = c("nh_white", "nh_black", "hispanic", "other", "nh_asian"))

tabyl(hi_latehh$race_eth)
 hi_latehh$race_eth    n    percent
           nh_white 2745 0.47237997
           nh_black 1246 0.21442093
           hispanic  992 0.17071072
              other  413 0.07107210
           nh_asian  415 0.07141628

3 state databases, categorical variables

Model 1: States with low scores

m1a, no controls

Code
m1a <- glmer(highrisk ~ score +
              (1 | abbv), 
            data = low_latehh, 
            family = binomial, 
            control = glmerControl(optimizer = "bobyqa"), 
            nAGQ = 10
            )

summary(m1a)
Generalized linear mixed model fit by maximum likelihood (Adaptive
  Gauss-Hermite Quadrature, nAGQ = 10) [glmerMod]
 Family: binomial  ( logit )
Formula: highrisk ~ score + (1 | abbv)
   Data: low_latehh
Control: glmerControl(optimizer = "bobyqa")

     AIC      BIC   logLik deviance df.resid 
  7077.5   7097.1  -3535.7   7071.5     5133 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.0737 -0.9368 -0.8191  1.0364  1.2444 

Random effects:
 Groups Name        Variance Std.Dev.
 abbv   (Intercept) 0.03164  0.1779  
Number of obs: 5136, groups:  abbv, 21

Fixed effects:
            Estimate Std. Error z value Pr(>|z|)
(Intercept) -0.07756    0.07384  -1.050    0.294
score       -0.14613    0.11851  -1.233    0.218

Correlation of Fixed Effects:
      (Intr)
score -0.753

m1b, with controls

Code
m1b <- glmer(highrisk ~ score +
               age +
               female +
               race_eth +
               inclvl +
               educ +
               unmarried +
               children +
               (1 | abbv),
             data = low_latehh,
             family = binomial,
             control = glmerControl(optimizer = "bobyqa"),
             nAGQ = 10
             )
Warning in checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, :
Model failed to converge with max|grad| = 0.002421 (tol = 0.002, component 1)
Warning in checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, : Model is nearly unidentifiable: very large eigenvalue
 - Rescale variables?
Code
summary(m1b)
Generalized linear mixed model fit by maximum likelihood (Adaptive
  Gauss-Hermite Quadrature, nAGQ = 10) [glmerMod]
 Family: binomial  ( logit )
Formula: highrisk ~ score + age + female + race_eth + inclvl + educ +  
    unmarried + children + (1 | abbv)
   Data: low_latehh
Control: glmerControl(optimizer = "bobyqa")

     AIC      BIC   logLik deviance df.resid 
  6887.8   6985.9  -3428.9   6857.8     5121 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.4897 -0.9350 -0.5894  0.9775  2.5846 

Random effects:
 Groups Name        Variance Std.Dev.
 abbv   (Intercept) 0.01879  0.1371  
Number of obs: 5136, groups:  abbv, 21

Fixed effects:
                  Estimate Std. Error z value   Pr(>|z|)    
(Intercept)      -0.193838   0.157686  -1.229   0.218974    
score            -0.129841   0.102754  -1.264   0.206369    
age              -0.001608   0.002390  -0.673   0.501029    
female           -0.009745   0.064883  -0.150   0.880611    
race_ethnh_black  0.205607   0.074120   2.774   0.005538 ** 
race_ethhispanic -0.235408   0.097013  -2.427   0.015242 *  
race_ethother     0.052634   0.105932   0.497   0.619283    
race_ethnh_asian -0.947954   0.196230  -4.831 0.00000136 ***
inclvllessthan25  0.292332   0.064451   4.536 0.00000574 ***
inclvl50andover  -0.409139   0.088824  -4.606 0.00000410 ***
educHSorless     -0.060477   0.072567  -0.833   0.404616    
educBAorhigher   -0.272336   0.070099  -3.885   0.000102 ***
unmarried         0.191951   0.067669   2.837   0.004559 ** 
children          0.194277   0.061592   3.154   0.001609 ** 
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Correlation matrix not shown by default, as p = 14 > 12.
Use print(x, correlation=TRUE)  or
    vcov(x)        if you need it
optimizer (bobyqa) convergence code: 0 (OK)
Model failed to converge with max|grad| = 0.002421 (tol = 0.002, component 1)
Model is nearly unidentifiable: very large eigenvalue
 - Rescale variables?

Model 2: States with med scores

m2a, no controls

Code
m2a <- glmer(highrisk ~ score +
              (1 | abbv), 
            data = med_latehh, 
            family = binomial, 
            control = glmerControl(optimizer = "bobyqa"), 
            nAGQ = 10
            )

summary(m2a)
Generalized linear mixed model fit by maximum likelihood (Adaptive
  Gauss-Hermite Quadrature, nAGQ = 10) [glmerMod]
 Family: binomial  ( logit )
Formula: highrisk ~ score + (1 | abbv)
   Data: med_latehh
Control: glmerControl(optimizer = "bobyqa")

     AIC      BIC   logLik deviance df.resid 
    7772     7792    -3883     7766     5710 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-0.9649 -0.8356 -0.8121  1.1366  1.3597 

Random effects:
 Groups Name        Variance Std.Dev.
 abbv   (Intercept) 0.04501  0.2122  
Number of obs: 5713, groups:  abbv, 13

Fixed effects:
             Estimate Std. Error z value Pr(>|z|)
(Intercept) -0.357460   0.278535  -1.283    0.199
score        0.009001   0.149640   0.060    0.952

Correlation of Fixed Effects:
      (Intr)
score -0.969

m2b, with controls

Code
m2b <- glmer(highrisk ~ score + 
              age + 
              female + 
              race_eth + 
              inclvl +
              educ +
              unmarried + 
              children +
              (1 | abbv), 
            data = med_latehh, 
            family = binomial, 
            control = glmerControl(optimizer = "bobyqa"), 
            nAGQ = 10
            )
Warning in checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, : Model is nearly unidentifiable: very large eigenvalue
 - Rescale variables?
Code
summary(m2b)
Generalized linear mixed model fit by maximum likelihood (Adaptive
  Gauss-Hermite Quadrature, nAGQ = 10) [glmerMod]
 Family: binomial  ( logit )
Formula: highrisk ~ score + age + female + race_eth + inclvl + educ +  
    unmarried + children + (1 | abbv)
   Data: med_latehh
Control: glmerControl(optimizer = "bobyqa")

     AIC      BIC   logLik deviance df.resid 
  7386.1   7485.8  -3678.0   7356.1     5698 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.6094 -0.8754 -0.5298  1.0051  3.6690 

Random effects:
 Groups Name        Variance Std.Dev.
 abbv   (Intercept) 0.04196  0.2048  
Number of obs: 5713, groups:  abbv, 13

Fixed effects:
                  Estimate Std. Error z value             Pr(>|z|)    
(Intercept)       0.194311   0.311100   0.625             0.532238    
score            -0.075875   0.147789  -0.513             0.607669    
age              -0.008848   0.002264  -3.908        0.00009288624 ***
female           -0.115774   0.061385  -1.886             0.059291 .  
race_ethnh_black  0.387834   0.075505   5.137        0.00000027984 ***
race_ethhispanic -0.178683   0.075743  -2.359             0.018321 *  
race_ethother     0.129151   0.127161   1.016             0.309796    
race_ethnh_asian -1.319651   0.143014  -9.227 < 0.0000000000000002 ***
inclvllessthan25  0.296591   0.063540   4.668        0.00000304522 ***
inclvl50andover  -0.456429   0.079148  -5.767        0.00000000808 ***
educHSorless     -0.080312   0.073276  -1.096             0.273065    
educBAorhigher   -0.239495   0.067967  -3.524             0.000426 ***
unmarried         0.184539   0.063541   2.904             0.003681 ** 
children          0.036908   0.060958   0.605             0.544866    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Correlation matrix not shown by default, as p = 14 > 12.
Use print(x, correlation=TRUE)  or
    vcov(x)        if you need it
optimizer (bobyqa) convergence code: 0 (OK)
Model is nearly unidentifiable: very large eigenvalue
 - Rescale variables?

Model 3: States with hi scores

m3a, no controls

Code
m3a <- glmer(highrisk ~ score +
              (1 | abbv), 
            data = hi_latehh, 
            family = binomial, 
            control = glmerControl(optimizer = "bobyqa"), 
            nAGQ = 10
            )

summary(m3a)
Generalized linear mixed model fit by maximum likelihood (Adaptive
  Gauss-Hermite Quadrature, nAGQ = 10) [glmerMod]
 Family: binomial  ( logit )
Formula: highrisk ~ score + (1 | abbv)
   Data: hi_latehh
Control: glmerControl(optimizer = "bobyqa")

     AIC      BIC   logLik deviance df.resid 
  7744.9   7764.9  -3869.5   7738.9     5808 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-0.9659 -0.8194 -0.7216  1.2040  1.4330 

Random effects:
 Groups Name        Variance Std.Dev.
 abbv   (Intercept) 0.03858  0.1964  
Number of obs: 5811, groups:  abbv, 17

Fixed effects:
            Estimate Std. Error z value Pr(>|z|)
(Intercept) -0.09486    0.38312  -0.248    0.804
score       -0.10196    0.10800  -0.944    0.345

Correlation of Fixed Effects:
      (Intr)
score -0.989

m3b, with controls

Code
m3b <- glmer(highrisk ~ score + 
              age + 
              female + 
              race_eth + 
              inclvl +
              educ +
              unmarried + 
              children +
              (1 | abbv), 
            data = hi_latehh, 
            family = binomial, 
            control = glmerControl(optimizer = "bobyqa"), 
            nAGQ = 10
            )
Warning in checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, :
Model failed to converge with max|grad| = 0.00225944 (tol = 0.002, component 1)
Warning in checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, : Model is nearly unidentifiable: very large eigenvalue
 - Rescale variables?
Code
summary(m3b)
Generalized linear mixed model fit by maximum likelihood (Adaptive
  Gauss-Hermite Quadrature, nAGQ = 10) [glmerMod]
 Family: binomial  ( logit )
Formula: highrisk ~ score + age + female + race_eth + inclvl + educ +  
    unmarried + children + (1 | abbv)
   Data: hi_latehh
Control: glmerControl(optimizer = "bobyqa")

     AIC      BIC   logLik deviance df.resid 
  7399.7   7499.7  -3684.9   7369.7     5796 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.4774 -0.8205 -0.5543  1.0435  3.5279 

Random effects:
 Groups Name        Variance Std.Dev.
 abbv   (Intercept) 0.03209  0.1791  
Number of obs: 5811, groups:  abbv, 17

Fixed effects:
                  Estimate Std. Error z value         Pr(>|z|)    
(Intercept)      -0.361746   0.386642  -0.936           0.3495    
score            -0.073030   0.102515  -0.712           0.4762    
age              -0.002073   0.002197  -0.944           0.3453    
female           -0.045320   0.060949  -0.744           0.4571    
race_ethnh_black  0.127168   0.075827   1.677           0.0935 .  
race_ethhispanic -0.036370   0.081234  -0.448           0.6544    
race_ethother     0.053942   0.111688   0.483           0.6291    
race_ethnh_asian -1.052080   0.151914  -6.926 0.00000000000434 ***
inclvllessthan25  0.416263   0.063262   6.580 0.00000000004705 ***
inclvl50andover  -0.336663   0.080453  -4.185 0.00002856654372 ***
educHSorless      0.033031   0.071936   0.459           0.6461    
educBAorhigher   -0.326489   0.068196  -4.787 0.00000168895613 ***
unmarried         0.283250   0.067319   4.208 0.00002581340082 ***
children          0.295783   0.060641   4.878 0.00000107383412 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Correlation matrix not shown by default, as p = 14 > 12.
Use print(x, correlation=TRUE)  or
    vcov(x)        if you need it
optimizer (bobyqa) convergence code: 0 (OK)
Model failed to converge with max|grad| = 0.00225944 (tol = 0.002, component 1)
Model is nearly unidentifiable: very large eigenvalue
 - Rescale variables?

Graphs of score’s effect on risk, without controls

Plot of low score states

Code
plot_low <- plot_predictions(m1a, condition = c("score")) +
  labs(title = "Eviction Risk",
       subtitle = "Renters Who Are Not Current On Rent",
       caption = "Source: HH Pulse and COVID State Housing Policy Score"
       ) +
  xlab(label = "State Score") +
  ylab(label = "Predicted Probability of High Risk of Eviction") # ylim(0, 0)
  
  #png('plot_low.png')

Plot of medium score states

Code
plot_med <- plot_predictions(m2a, condition = c("score")) +
  labs(title = "Eviction Risk",
       subtitle = "Renters Who Are Not Current On Rent",
       caption = "Source: HH Pulse and COVID State Housing Policy Score"
       ) +
  xlab(label = "State Score") +
  ylab(label = "Predicted Probability of High Risk of Eviction") # ylim(0, 0)
  
  #png('plot_med.png')

Plot of low score states

Code
plot_high <- plot_predictions(m3a, condition = c("score")) +
  labs(title = "Eviction Risk",
       subtitle = "Renters Who Are Not Current On Rent",
       caption = "Source: HH Pulse and COVID State Housing Policy Score"
       ) +
  xlab(label = "State Score") +
  ylab(label = "Predicted Probability of High Risk of Eviction") # ylim(0, 0)
  
  #png('plot_high.png')

##All plots

Code
plot_low

Code
plot_med

Code
plot_high

Model 10: latehh df, focal variable children, with categorical variables

m10, highrisk by score, children only

Code
m10 <- glmer(highrisk ~ score +
               children +
              (1 | abbv), 
            data = latehh, 
            family = binomial, 
            control = glmerControl(optimizer = "bobyqa"), 
            nAGQ = 10
            )

summary(m10)
Generalized linear mixed model fit by maximum likelihood (Adaptive
  Gauss-Hermite Quadrature, nAGQ = 10) [glmerMod]
 Family: binomial  ( logit )
Formula: highrisk ~ score + children + (1 | abbv)
   Data: latehh
Control: glmerControl(optimizer = "bobyqa")

     AIC      BIC   logLik deviance df.resid 
 22550.7  22581.6 -11271.4  22542.7    16656 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.1226 -0.8640 -0.7413  1.1045  1.4823 

Random effects:
 Groups Name        Variance Std.Dev.
 abbv   (Intercept) 0.03803  0.195   
Number of obs: 16660, groups:  abbv, 51

Fixed effects:
            Estimate Std. Error z value      Pr(>|z|)    
(Intercept) -0.21608    0.05689  -3.798      0.000146 ***
score       -0.09497    0.02346  -4.048 0.00005168659 ***
children     0.18839    0.03174   5.936 0.00000000292 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Correlation of Fixed Effects:
         (Intr) score 
score    -0.774       
children -0.289  0.024
Code
plot_10 <- plot_predictions(m10, condition = c("score", "children")) +
  labs(title = "Predicted Probability of High Eviction Risk, Children in Household",
       subtitle = "Renters Who Are Not Current On Rent",
       caption = "Source: HH Pulse and COVID State Housing Policy Score"
       ) +
  scale_color_discrete(name="Model 10",
                       breaks=c("1", "0"),
                       labels=c("Children in HH", "No Children in HH")
                       ) +
  scale_fill_discrete(name="Model 10",
                       breaks=c("1", "0"),
                       labels=c("Children in HH", "No Children in HH")
                       ) +
  xlab(label = "State Score") +
  ylab(label = "Predicted Probability of High Eviction Risk") #+ ylim(0, 0)
  # png('10.png')
Code
plot_10

m10a, latehh, highrisk ~ score, all covariates

Code
m10a <- glmer(highrisk ~ score +
                children +
                age +
                female +
                race_eth +
                inclvl +
                educ +
                married +
                (1 | abbv), 
            data = latehh, 
            family = binomial, 
            control = glmerControl(optimizer = "bobyqa"), 
            nAGQ = 10
            )
Warning in checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, : Model is nearly unidentifiable: very large eigenvalue
 - Rescale variables?
Code
summary(m10a)
Generalized linear mixed model fit by maximum likelihood (Adaptive
  Gauss-Hermite Quadrature, nAGQ = 10) [glmerMod]
 Family: binomial  ( logit )
Formula: highrisk ~ score + children + age + female + race_eth + inclvl +  
    educ + married + (1 | abbv)
   Data: latehh
Control: glmerControl(optimizer = "bobyqa")

     AIC      BIC   logLik deviance df.resid 
 21654.9  21770.7 -10812.5  21624.9    16645 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.5797 -0.8763 -0.5673  1.0066  3.4430 

Random effects:
 Groups Name        Variance Std.Dev.
 abbv   (Intercept) 0.03011  0.1735  
Number of obs: 16660, groups:  abbv, 51

Fixed effects:
                  Estimate Std. Error z value             Pr(>|z|)    
(Intercept)       0.091885   0.091066   1.009              0.31298    
score            -0.070948   0.021847  -3.247              0.00116 ** 
children          0.177098   0.035179   5.034     0.00000047971852 ***
age              -0.004195   0.001312  -3.197              0.00139 ** 
female           -0.056044   0.035905  -1.561              0.11854    
race_ethnh_black  0.239295   0.042895   5.579     0.00000002424333 ***
race_ethhispanic -0.142958   0.047487  -3.010              0.00261 ** 
race_ethother     0.075269   0.065623   1.147              0.25139    
race_ethnh_asian -1.154767   0.091586 -12.609 < 0.0000000000000002 ***
inclvllessthan25  0.330804   0.036716   9.010 < 0.0000000000000002 ***
inclvl50andover  -0.403344   0.047477  -8.496 < 0.0000000000000002 ***
educHSorless     -0.036378   0.041843  -0.869              0.38464    
educBAorhigher   -0.279033   0.039552  -7.055     0.00000000000173 ***
married          -0.224558   0.038067  -5.899     0.00000000365723 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Correlation matrix not shown by default, as p = 14 > 12.
Use print(x, correlation=TRUE)  or
    vcov(x)        if you need it
optimizer (bobyqa) convergence code: 0 (OK)
Model is nearly unidentifiable: very large eigenvalue
 - Rescale variables?
Code
plot_10a <- plot_predictions(m10a, condition = c("score", "children")) +
  labs(title = "Predicted Probability of High Eviction Risk, Children in Household",
       subtitle = "Renters Who Are Not Current On Rent",
       caption = "Source: HH Pulse and COVID State Housing Policy Score"
       ) +
  scale_color_discrete(name="Model 10a",
                       breaks=c("1", "0"),
                       labels=c("Children in HH", "No Children in HH")
                       ) +
  scale_fill_discrete(name="Model 10a",
                       breaks=c("1", "0"),
                       labels=c("Children in HH", "No Children in HH")
                       ) +
  xlab(label = "State Score") +
  ylab(label = "Predicted Probability of High Eviction Risk") #+ ylim(0, 0)
  # png('10a.png')
Code
plot_10a

m10b, latehh, highrisk ~ score, all covariates, score:children

Code
m10b <- glmer(highrisk ~ score +
                score:children +
                age +
                female +
                race_eth +
                inclvl +
                educ +
                married +
                (1 | abbv), 
            data = latehh, 
            family = binomial, 
            control = glmerControl(optimizer = "bobyqa"), 
            nAGQ = 10
            )
Warning in checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, :
Model failed to converge with max|grad| = 0.00407356 (tol = 0.002, component 1)
Warning in checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, : Model is nearly unidentifiable: very large eigenvalue
 - Rescale variables?
Code
summary(m10b)
Generalized linear mixed model fit by maximum likelihood (Adaptive
  Gauss-Hermite Quadrature, nAGQ = 10) [glmerMod]
 Family: binomial  ( logit )
Formula: highrisk ~ score + score:children + age + female + race_eth +  
    inclvl + educ + married + (1 | abbv)
   Data: latehh
Control: glmerControl(optimizer = "bobyqa")

     AIC      BIC   logLik deviance df.resid 
 21657.6  21773.4 -10813.8  21627.6    16645 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.5162 -0.8773 -0.5668  1.0075  3.4713 

Random effects:
 Groups Name        Variance Std.Dev.
 abbv   (Intercept) 0.02997  0.1731  
Number of obs: 16660, groups:  abbv, 51

Fixed effects:
                  Estimate Std. Error z value             Pr(>|z|)    
(Intercept)       0.190439   0.087746   2.170             0.029982 *  
score            -0.105346   0.022838  -4.613     0.00000397527289 ***
age              -0.004620   0.001297  -3.563             0.000367 ***
female           -0.049590   0.035739  -1.388             0.165274    
race_ethnh_black  0.244927   0.042815   5.721     0.00000001061256 ***
race_ethhispanic -0.139563   0.047442  -2.942             0.003263 ** 
race_ethother     0.081736   0.065581   1.246             0.212640    
race_ethnh_asian -1.154301   0.091582 -12.604 < 0.0000000000000002 ***
inclvllessthan25  0.329989   0.036707   8.990 < 0.0000000000000002 ***
inclvl50andover  -0.402719   0.047473  -8.483 < 0.0000000000000002 ***
educHSorless     -0.035897   0.041834  -0.858             0.390849    
educBAorhigher   -0.280387   0.039541  -7.091     0.00000000000133 ***
married          -0.216856   0.037871  -5.726     0.00000001027417 ***
score:children    0.071395   0.014997   4.761     0.00000192982423 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Correlation matrix not shown by default, as p = 14 > 12.
Use print(x, correlation=TRUE)  or
    vcov(x)        if you need it
optimizer (bobyqa) convergence code: 0 (OK)
Model failed to converge with max|grad| = 0.00407356 (tol = 0.002, component 1)
Model is nearly unidentifiable: very large eigenvalue
 - Rescale variables?
Code
plot_10b <- plot_predictions(m10b, condition = c("score", "children")) +
  labs(title = "Predicted Probability of High Eviction Risk, Children in Household",
       subtitle = "Renters Who Are Not Current On Rent",
       caption = "Source: HH Pulse and COVID State Housing Policy Score"
       ) +
  scale_color_discrete(name="Model 10b",
                       breaks=c("1", "0"),
                       labels=c("Children in HH", "No Children in HH")
                       ) +
  scale_fill_discrete(name="Model 10b",
                       breaks=c("1", "0"),
                       labels=c("Children in HH", "No Children in HH")
                       ) +
  xlab(label = "State Score") +
  ylab(label = "Predicted Probability of High Eviction Risk") #+ ylim(0, 0)
  # png('10b.png')

plot_10b

m10c, latehh, highrisk ~ score, all covariates, score:children

Code
m10c <- glmer(highrisk ~ 
                score*(children +
                age +
                female +
                race_eth +
                inclvl +
                educ +
                married
                ) +
                (1 | abbv), 
            data = latehh, 
            family = binomial, 
            control = glmerControl(optimizer = "bobyqa"), 
            nAGQ = 10
            )
Warning in checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, :
Model failed to converge with max|grad| = 0.0558651 (tol = 0.002, component 1)
Warning in checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, : Model is nearly unidentifiable: very large eigenvalue
 - Rescale variables?
Code
summary(m10c)
Generalized linear mixed model fit by maximum likelihood (Adaptive
  Gauss-Hermite Quadrature, nAGQ = 10) [glmerMod]
 Family: binomial  ( logit )
Formula: highrisk ~ score * (children + age + female + race_eth + inclvl +  
    educ + married) + (1 | abbv)
   Data: latehh
Control: glmerControl(optimizer = "bobyqa")

     AIC      BIC   logLik deviance df.resid 
 21664.2  21872.7 -10805.1  21610.2    16633 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.5785 -0.8786 -0.5642  1.0056  3.5183 

Random effects:
 Groups Name        Variance Std.Dev.
 abbv   (Intercept) 0.02914  0.1707  
Number of obs: 16660, groups:  abbv, 51

Fixed effects:
                         Estimate Std. Error z value        Pr(>|z|)    
(Intercept)             0.2018341  0.1419646   1.422        0.155107    
score                  -0.1266596  0.0611364  -2.072        0.038288 *  
children                0.1138942  0.0607914   1.874        0.060996 .  
age                    -0.0053742  0.0023150  -2.321        0.020263 *  
female                 -0.0436234  0.0630844  -0.692        0.489246    
race_ethnh_black        0.3191940  0.0728970   4.379 0.0000119390936 ***
race_ethhispanic       -0.2785820  0.0881250  -3.161        0.001571 ** 
race_ethother           0.0471319  0.1107578   0.426        0.670443    
race_ethnh_asian       -1.1728778  0.1725570  -6.797 0.0000000000107 ***
inclvllessthan25        0.2630142  0.0635813   4.137 0.0000352403417 ***
inclvl50andover        -0.4679409  0.0845541  -5.534 0.0000000312617 ***
educHSorless           -0.0823847  0.0722975  -1.140        0.254485    
educBAorhigher         -0.2360528  0.0686646  -3.438        0.000587 ***
married                -0.2092736  0.0660089  -3.170        0.001522 ** 
score:children          0.0329777  0.0265672   1.241        0.214497    
score:age               0.0006016  0.0009885   0.609        0.542802    
score:female           -0.0080514  0.0272241  -0.296        0.767425    
score:race_ethnh_black -0.0442036  0.0319099  -1.385        0.165972    
score:race_ethhispanic  0.0670706  0.0378075   1.774        0.076062 .  
score:race_ethother     0.0155116  0.0487545   0.318        0.750366    
score:race_ethnh_asian  0.0079672  0.0733980   0.109        0.913561    
score:inclvllessthan25  0.0367850  0.0278164   1.322        0.186027    
score:inclvl50andover   0.0346523  0.0361250   0.959        0.337441    
score:educHSorless      0.0245139  0.0315886   0.776        0.437729    
score:educBAorhigher   -0.0230497  0.0300300  -0.768        0.442750    
score:married          -0.0073946  0.0292440  -0.253        0.800378    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Correlation matrix not shown by default, as p = 26 > 12.
Use print(x, correlation=TRUE)  or
    vcov(x)        if you need it
optimizer (bobyqa) convergence code: 0 (OK)
Model failed to converge with max|grad| = 0.0558651 (tol = 0.002, component 1)
Model is nearly unidentifiable: very large eigenvalue
 - Rescale variables?
Code
plot_10c <- plot_predictions(m10c, condition = c("score", "children")) +
  labs(title = "Predicted Probability of High Eviction Risk, Children in Household",
       subtitle = "Renters Who Are Not Current On Rent",
       caption = "Source: HH Pulse and COVID State Housing Policy Score"
       ) +
  scale_color_discrete(name="Model 10c",
                       breaks=c("1", "0"),
                       labels=c("Children in HH", "No Children in HH")
                       ) +
  scale_fill_discrete(name="Model 10c",
                       breaks=c("1", "0"),
                       labels=c("Children in HH", "No Children in HH")
                       ) +
  xlab(label = "State Score") +
  ylab(label = "Predicted Probability of High Eviction Risk") #+ ylim(0, 0)
  # png('10b.png')

plot_10c

Model 4: unmarried households with children

Make dummy variable for unmarried households with children

Code
latehh <- latehh %>%
  mutate(unwch = case_when(married == 0 & children == 1 ~ 1,
                              TRUE ~ 0
                              )
        )

tabyl(latehh$unwch)
 latehh$unwch     n   percent
            0 11658 0.6997599
            1  5002 0.3002401

Eviction risk by state score for Unmarried with Children

Code
scoreprobs <- latehh %>%
  group_by(unwch, score) %>%
  summarize(p=mean(highrisk), n=n())
`summarise()` has grouped output by 'unwch'. You can override using the
`.groups` argument.
Code
scoreprobs$num <- scoreprobs$p*(1-scoreprobs$p)

scoreprobs$sep <- sqrt(scoreprobs$num/scoreprobs$n)

scoreprobs$me <- 2*scoreprobs$sep

scoreprobs <- scoreprobs %>% 
  mutate(umwc = case_when(unwch == 1 ~ 'Unmarried with Children',
                           unwch == 0 ~ 'Married or No Children'
                           )
         )

ggplot(scoreprobs, aes(x = score, y = p, group=umwc)) +
  geom_point() +
  geom_line(aes(color = umwc)) + 
  labs(x = "Score", y = "Predicted Probability of High Eviction Risk", title = "Predicted Probability of High Eviction Risk for Unmarried with Children") +
  scale_color_discrete(name="Probability of High Eviction Risk",
                       breaks=c("Unmarried with Children", "Married or No Children"),
                       labels=c("Unmarried with Children", "Married or No Children")
                       )

m4, unmarried with children the only covariate

Code
m4 <- glmer(highrisk ~ score +
               unwch +
              (1 | abbv), 
            data = latehh, 
            family = binomial, 
            control = glmerControl(optimizer = "bobyqa"), 
            nAGQ = 10
            )

summary(m4)
Generalized linear mixed model fit by maximum likelihood (Adaptive
  Gauss-Hermite Quadrature, nAGQ = 10) [glmerMod]
 Family: binomial  ( logit )
Formula: highrisk ~ score + unwch + (1 | abbv)
   Data: latehh
Control: glmerControl(optimizer = "bobyqa")

     AIC      BIC   logLik deviance df.resid 
 22414.1  22445.0 -11203.1  22406.1    16656 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.2442 -0.8474 -0.7100  1.1024  1.5162 

Random effects:
 Groups Name        Variance Std.Dev.
 abbv   (Intercept) 0.03743  0.1935  
Number of obs: 16660, groups:  abbv, 51

Fixed effects:
            Estimate Std. Error z value             Pr(>|z|)    
(Intercept) -0.27159    0.05548  -4.895          0.000000982 ***
score       -0.09287    0.02335  -3.977          0.000069782 ***
unwch        0.45132    0.03444  13.106 < 0.0000000000000002 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Correlation of Fixed Effects:
      (Intr) score 
score -0.787       
unwch -0.211  0.015
Code
plot_m4 <- plot_predictions(m4, condition = c("score", "unwch")) +
  labs(title = "Predicted Probability of High Eviction Risk, Unmarried with Children",
       subtitle = "Renters Who Are Not Current On Rent",
       caption = "Source: HH Pulse and COVID State Housing Policy Score"
       ) +
  scale_color_discrete(name="Model 4",
                       breaks=c("1", "0"),
                       labels=c("Unmarried with Chilren", "Married or No Children")
                       ) +
  scale_fill_discrete(name="Model 4",
                       breaks=c("1", "0"),
                       labels=c("Unmarried with Chilren", "Married or No Children")
                      ) +
  xlab(label = "State Score") +
  ylab(label = "Predicted Probability of High Eviction Risk") #+ ylim(0, 0)
  # png('m4.png')
Code
plot_m4

m4a, unmarried with children with covariates

Code
m4a <- glmer(highrisk ~ score +
              unwch +
              age + 
              female + 
              race_eth + 
              inclvl +
              educ +
              (1 | abbv), 
            data = latehh, 
            family = binomial, 
            control = glmerControl(optimizer = "bobyqa"), 
            nAGQ = 10
            )
Warning in checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, :
Model failed to converge with max|grad| = 0.00681388 (tol = 0.002, component 1)
Warning in checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, : Model is nearly unidentifiable: very large eigenvalue
 - Rescale variables?
Code
summary(m4a)
Generalized linear mixed model fit by maximum likelihood (Adaptive
  Gauss-Hermite Quadrature, nAGQ = 10) [glmerMod]
 Family: binomial  ( logit )
Formula: highrisk ~ score + unwch + age + female + race_eth + inclvl +  
    educ + (1 | abbv)
   Data: latehh
Control: glmerControl(optimizer = "bobyqa")

     AIC      BIC   logLik deviance df.resid 
 21649.8  21757.9 -10810.9  21621.8    16646 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.6156 -0.8730 -0.5689  1.0102  3.2630 

Random effects:
 Groups Name        Variance Std.Dev.
 abbv   (Intercept) 0.02986  0.1728  
Number of obs: 16660, groups:  abbv, 51

Fixed effects:
                  Estimate Std. Error z value             Pr(>|z|)    
(Intercept)       0.033850   0.088372   0.383             0.701689    
score            -0.069643   0.021779  -3.198             0.001385 ** 
unwch             0.275451   0.036998   7.445   0.0000000000000969 ***
age              -0.004305   0.001290  -3.337             0.000846 ***
female           -0.063767   0.035859  -1.778             0.075357 .  
race_ethnh_black  0.237965   0.042879   5.550   0.0000000286290402 ***
race_ethhispanic -0.150553   0.047366  -3.178             0.001481 ** 
race_ethother     0.076339   0.065627   1.163             0.244744    
race_ethnh_asian -1.171410   0.091268 -12.835 < 0.0000000000000002 ***
inclvllessthan25  0.342820   0.036428   9.411 < 0.0000000000000002 ***
inclvl50andover  -0.413240   0.047196  -8.756 < 0.0000000000000002 ***
educHSorless     -0.042148   0.041795  -1.008             0.313246    
educBAorhigher   -0.275301   0.039560  -6.959   0.0000000000034270 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Correlation matrix not shown by default, as p = 13 > 12.
Use print(x, correlation=TRUE)  or
    vcov(x)        if you need it
optimizer (bobyqa) convergence code: 0 (OK)
Model failed to converge with max|grad| = 0.00681388 (tol = 0.002, component 1)
Model is nearly unidentifiable: very large eigenvalue
 - Rescale variables?
Code
plot_m4a <- plot_predictions(m4a, condition = c("score", "unwch")) +
  labs(title = "Predicted Probability of High Eviction Risk, Unmarried with Children",
       subtitle = "Renters Who Are Not Current On Rent",
       caption = "Source: HH Pulse and COVID State Housing Policy Score"
       ) +
  scale_color_discrete(name="Model 4a",
                       breaks=c("1", "0"),
                       labels=c("Unmarried with Chilren", "Married or No Children")
                       ) +
  scale_fill_discrete(name="Model 4a",
                       breaks=c("1", "0"),
                       labels=c("Unmarried with Chilren", "Married or No Children")
                      ) +
  xlab(label = "State Score") +
  ylab(label = "Predicted Probability of High Eviction Risk") #+ ylim(0, 0)
  # png('m4a.png')


plot_m4a

m4b, unmarried with children, interacted with score

Code
m4b <- glmer(highrisk ~ score +
               age +
               female +
               race_eth +
               inclvl +
               educ +
               score:unwch +
               (1 | abbv),
             data = latehh,
             family = binomial,
             control = glmerControl(optimizer = "bobyqa"),
             nAGQ = 10
             )
Warning in checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, :
Model failed to converge with max|grad| = 0.00328838 (tol = 0.002, component 1)
Warning in checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, : Model is nearly unidentifiable: very large eigenvalue
 - Rescale variables?
Code
summary(m4b)
Generalized linear mixed model fit by maximum likelihood (Adaptive
  Gauss-Hermite Quadrature, nAGQ = 10) [glmerMod]
 Family: binomial  ( logit )
Formula: highrisk ~ score + age + female + race_eth + inclvl + educ +  
    score:unwch + (1 | abbv)
   Data: latehh
Control: glmerControl(optimizer = "bobyqa")

     AIC      BIC   logLik deviance df.resid 
 21660.1  21768.2 -10816.0  21632.1    16646 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.5741 -0.8779 -0.5664  1.0078  3.3497 

Random effects:
 Groups Name        Variance Std.Dev.
 abbv   (Intercept) 0.02949  0.1717  
Number of obs: 16660, groups:  abbv, 51

Fixed effects:
                  Estimate Std. Error z value             Pr(>|z|)    
(Intercept)       0.133640   0.086726   1.541             0.123332    
score            -0.104537   0.022262  -4.696     0.00000265488432 ***
age              -0.004743   0.001284  -3.694             0.000221 ***
female           -0.051251   0.035663  -1.437             0.150695    
race_ethnh_black  0.248686   0.042761   5.816     0.00000000603755 ***
race_ethhispanic -0.147017   0.047336  -3.106             0.001898 ** 
race_ethother     0.083923   0.065566   1.280             0.200552    
race_ethnh_asian -1.174493   0.091287 -12.866 < 0.0000000000000002 ***
inclvllessthan25  0.345787   0.036402   9.499 < 0.0000000000000002 ***
inclvl50andover  -0.417776   0.047174  -8.856 < 0.0000000000000002 ***
educHSorless     -0.042523   0.041774  -1.018             0.308718    
educBAorhigher   -0.278480   0.039533  -7.044     0.00000000000187 ***
score:unwch       0.107172   0.015940   6.723     0.00000000001775 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Correlation matrix not shown by default, as p = 13 > 12.
Use print(x, correlation=TRUE)  or
    vcov(x)        if you need it
optimizer (bobyqa) convergence code: 0 (OK)
Model failed to converge with max|grad| = 0.00328838 (tol = 0.002, component 1)
Model is nearly unidentifiable: very large eigenvalue
 - Rescale variables?
Code
plot_m4b <- plot_predictions(m4b, condition = c("score","unwch")) +
  labs(title = "Eviction Risk by Score x unmarried with children",
       subtitle = "Renters Who Are Not Current On Rent",
       caption = "Source: HH Pulse and COVID State Housing Policy Score"
       ) +
   scale_color_discrete(name="Model 4b",
                       breaks=c("1", "0"),
                       labels=c("Unmarried with Chilren", "Married or No Children")
                       ) +
  scale_fill_discrete(name="Model 4b",
                       breaks=c("1", "0"),
                       labels=c("Unmarried with Chilren", "Married or No Children")
                      ) +
  xlab(label = "State Score") +
  ylab(label = "Predicted Probability of High Risk of Eviction") #+ # ylim(0, 0)
  #png('m4b.png')

plot_m4b

Model 5: Black/unmarried/children

Code
tabyl(latehh, unwch, race_eth)
 unwch nh_white nh_black hispanic other nh_asian
     0     5620     2324     1978   757      979
     1     1979     1530      982   405      106

We have 1979 unwed nh_white hh with children and 1530 unwed nh_black hh with children, so let’s make dummy variables for those.

Make dummy variable for black unmarried households with children

Code
latehh <- latehh %>%
  mutate(bunwch = case_when(race_eth == "nh_black" & married == 0 & children == 1 ~ 1,
                              TRUE ~ 0
                              )
        )

tabyl(latehh$bunwch)
 latehh$bunwch     n    percent
             0 15130 0.90816327
             1  1530 0.09183673

m5, black unmarried with children the only covariate

Code
m5 <- glmer(highrisk ~ score +
               bunwch +
              (1 | abbv), 
            data = latehh, 
            family = binomial, 
            control = glmerControl(optimizer = "bobyqa"), 
            nAGQ = 10
            )

summary(m5)
Generalized linear mixed model fit by maximum likelihood (Adaptive
  Gauss-Hermite Quadrature, nAGQ = 10) [glmerMod]
 Family: binomial  ( logit )
Formula: highrisk ~ score + bunwch + (1 | abbv)
   Data: latehh
Control: glmerControl(optimizer = "bobyqa")

     AIC      BIC   logLik deviance df.resid 
 22467.3  22498.2 -11229.7  22459.3    16656 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.3665 -0.8374 -0.7375  1.1040  1.4980 

Random effects:
 Groups Name        Variance Std.Dev.
 abbv   (Intercept) 0.03408  0.1846  
Number of obs: 16660, groups:  abbv, 51

Fixed effects:
            Estimate Std. Error z value             Pr(>|z|)    
(Intercept) -0.18056    0.05275  -3.423              0.00062 ***
score       -0.09500    0.02258  -4.207            0.0000259 ***
bunwch       0.60489    0.05578  10.843 < 0.0000000000000002 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Correlation of Fixed Effects:
       (Intr) score 
score  -0.799       
bunwch -0.106  0.011
Code
plot_m5 <- plot_predictions(m5, condition = c("score", "bunwch")) +
  labs(title = "Predicted Probability of High Eviction Risk, Black Unmarried with Children",
       subtitle = "Renters Who Are Not Current On Rent",
       caption = "Source: HH Pulse and COVID State Housing Policy Score"
       ) +
  scale_color_discrete(name="Model 5",
                       breaks=c("1", "0"),
                       labels=c("Black Unmarried with Chilren", "Not Black Unmarried with Chilren")
                       ) +
  scale_fill_discrete(name="Model 5",
                       breaks=c("1", "0"),
                       labels=c("Black Unmarried with Chilren", "Not Black Unmarried with Chilren")
                      ) +
  xlab(label = "State Score") +
  ylab(label = "Predicted Probability of High Eviction Risk") #+ ylim(0, 0)
  # png('m5.png')
Code
plot_m5

m5a, black unmarried with children with covariates

Code
m5a <- glmer(highrisk ~ score +
              bunwch +
              age + 
              female + 
              race_eth + 
              inclvl +
              educ +
              (1 | abbv), 
            data = latehh, 
            family = binomial, 
            control = glmerControl(optimizer = "bobyqa"), 
            nAGQ = 10
            )
Warning in checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, : Model is nearly unidentifiable: very large eigenvalue
 - Rescale variables?
Code
summary(m5a)
Generalized linear mixed model fit by maximum likelihood (Adaptive
  Gauss-Hermite Quadrature, nAGQ = 10) [glmerMod]
 Family: binomial  ( logit )
Formula: highrisk ~ score + bunwch + age + female + race_eth + inclvl +  
    educ + (1 | abbv)
   Data: latehh
Control: glmerControl(optimizer = "bobyqa")

     AIC      BIC   logLik deviance df.resid 
 21686.3  21794.4 -10829.1  21658.3    16646 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.6262 -0.8778 -0.5776  1.0162  3.3424 

Random effects:
 Groups Name        Variance Std.Dev.
 abbv   (Intercept) 0.02854  0.1689  
Number of obs: 16660, groups:  abbv, 51

Fixed effects:
                  Estimate Std. Error z value             Pr(>|z|)    
(Intercept)       0.131660   0.086493   1.522              0.12796    
score            -0.070340   0.021452  -3.279              0.00104 ** 
bunwch            0.298368   0.068734   4.341    0.000014191671619 ***
age              -0.005259   0.001280  -4.110    0.000039542086568 ***
female           -0.027641   0.035312  -0.783              0.43377    
race_ethnh_black  0.152418   0.050230   3.034              0.00241 ** 
race_ethhispanic -0.130110   0.047156  -2.759              0.00580 ** 
race_ethother     0.094955   0.065428   1.451              0.14670    
race_ethnh_asian -1.193447   0.091160 -13.092 < 0.0000000000000002 ***
inclvllessthan25  0.347694   0.036379   9.557 < 0.0000000000000002 ***
inclvl50andover  -0.427625   0.047080  -9.083 < 0.0000000000000002 ***
educHSorless     -0.035412   0.041734  -0.849              0.39615    
educBAorhigher   -0.284332   0.039485  -7.201    0.000000000000598 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Correlation matrix not shown by default, as p = 13 > 12.
Use print(x, correlation=TRUE)  or
    vcov(x)        if you need it
optimizer (bobyqa) convergence code: 0 (OK)
Model is nearly unidentifiable: very large eigenvalue
 - Rescale variables?
Code
plot_m5a <- plot_predictions(m5a, condition = c("score", "bunwch")) +
  labs(title = "Predicted Probability of High Eviction Risk, Black Unmarried with Children",
       subtitle = "Renters Who Are Not Current On Rent",
       caption = "Source: HH Pulse and COVID State Housing Policy Score"
       ) +
  scale_color_discrete(name="Model 5a",
                       breaks=c("1", "0"),
                       labels=c("Black Unmarried with Chilren", "Not Black Unmarried with Chilren")
                       ) +
  scale_fill_discrete(name="Model 5a",
                       breaks=c("1", "0"),
                       labels=c("Black Unmarried with Chilren", "Not Black Unmarried with Chilren")
                      ) +
  xlab(label = "State Score") +
  ylab(label = "Predicted Probability of High Eviction Risk") #+ ylim(0, 0)
  # png('m5a.png')
Code
plot_m5a

m5b, black unmarried with children, interacted with score

Code
m5b <- glmer(highrisk ~ score +
               age +
               female +
               race_eth +
               inclvl +
               educ +
               score:bunwch +
               (1 | abbv),
             data = latehh,
             family = binomial,
             control = glmerControl(optimizer = "bobyqa"),
             nAGQ = 10
             )
Warning in checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, : Model is nearly unidentifiable: very large eigenvalue
 - Rescale variables?
Code
summary(m5b)
Generalized linear mixed model fit by maximum likelihood (Adaptive
  Gauss-Hermite Quadrature, nAGQ = 10) [glmerMod]
 Family: binomial  ( logit )
Formula: highrisk ~ score + age + female + race_eth + inclvl + educ +  
    score:bunwch + (1 | abbv)
   Data: latehh
Control: glmerControl(optimizer = "bobyqa")

     AIC      BIC   logLik deviance df.resid 
 21699.7  21807.7 -10835.8  21671.7    16646 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.5348 -0.8814 -0.5758  1.0145  3.4012 

Random effects:
 Groups Name        Variance Std.Dev.
 abbv   (Intercept) 0.02949  0.1717  
Number of obs: 16660, groups:  abbv, 51

Fixed effects:
                  Estimate Std. Error z value             Pr(>|z|)    
(Intercept)       0.159868   0.086589   1.846              0.06485 .  
score            -0.077614   0.021839  -3.554              0.00038 ***
age              -0.005692   0.001275  -4.464    0.000008051232182 ***
female           -0.019199   0.035246  -0.545              0.58595    
race_ethnh_black  0.222149   0.046818   4.745    0.000002085900085 ***
race_ethhispanic -0.129841   0.047195  -2.751              0.00594 ** 
race_ethother     0.092857   0.065474   1.418              0.15613    
race_ethnh_asian -1.188984   0.091188 -13.039 < 0.0000000000000002 ***
inclvllessthan25  0.349770   0.036360   9.620 < 0.0000000000000002 ***
inclvl50andover  -0.430810   0.047069  -9.153 < 0.0000000000000002 ***
educHSorless     -0.035716   0.041712  -0.856              0.39185    
educBAorhigher   -0.286458   0.039466  -7.258    0.000000000000392 ***
score:bunwch      0.066740   0.028354   2.354              0.01858 *  
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Correlation matrix not shown by default, as p = 13 > 12.
Use print(x, correlation=TRUE)  or
    vcov(x)        if you need it
optimizer (bobyqa) convergence code: 0 (OK)
Model is nearly unidentifiable: very large eigenvalue
 - Rescale variables?
Code
plot_m5b <- plot_predictions(m5b, condition = c("score","bunwch")) +
  labs(title = "Eviction Risk by Score : black unmarried with children",
       subtitle = "Renters Who Are Not Current On Rent",
       caption = "Source: HH Pulse and COVID State Housing Policy Score"
       ) +
  scale_color_discrete(name="Model 5b",
                       breaks=c("1", "0"),
                       labels=c("Black Unmarried with Chilren", "Not Black Unmarried with Chilren")
                       ) +
  scale_fill_discrete(name="Model 5b",
                       breaks=c("1", "0"),
                       labels=c("Black Unmarried with Chilren", "Not Black Unmarried with Chilren")
                      ) +
  xlab(label = "State Score") +
  ylab(label = "Predicted Probability of High Risk of Eviction") #+ # ylim(0, 0)
  #png('m5b.png')
Code
plot_m5b

Model 6: White/unmarried/children

Make dummy variable for white unmarried households with children

Code
latehh <- latehh %>%
  mutate(wunwch = case_when(race_eth == "nh_white" & married == 0 & children == 1 ~ 1,
                              TRUE ~ 0
                              )
        )

tabyl(latehh$wunwch)
 latehh$wunwch     n   percent
             0 14681 0.8812125
             1  1979 0.1187875

m6, white unmarried with children the only covariate

Code
m6 <- glmer(highrisk ~ score +
               wunwch +
              (1 | abbv), 
            data = latehh, 
            family = binomial, 
            control = glmerControl(optimizer = "bobyqa"), 
            nAGQ = 10
            )

summary(m6)
Generalized linear mixed model fit by maximum likelihood (Adaptive
  Gauss-Hermite Quadrature, nAGQ = 10) [glmerMod]
 Family: binomial  ( logit )
Formula: highrisk ~ score + wunwch + (1 | abbv)
   Data: latehh
Control: glmerControl(optimizer = "bobyqa")

     AIC      BIC   logLik deviance df.resid 
 22558.2  22589.0 -11275.1  22550.2    16656 

Scaled residuals: 
   Min     1Q Median     3Q    Max 
-1.211 -0.846 -0.738  1.096  1.445 

Random effects:
 Groups Name        Variance Std.Dev.
 abbv   (Intercept) 0.04043  0.2011  
Number of obs: 16660, groups:  abbv, 51

Fixed effects:
            Estimate Std. Error z value    Pr(>|z|)    
(Intercept) -0.15915    0.05619  -2.832     0.00462 ** 
score       -0.09617    0.02398  -4.011 0.000060577 ***
wunwch       0.25892    0.04903   5.280 0.000000129 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Correlation of Fixed Effects:
       (Intr) score 
score  -0.796       
wunwch -0.137  0.018
Code
plot_m6 <- plot_predictions(m6, condition = c("score", "wunwch")) +
  labs(title = "Predicted Probability of High Eviction Risk, White Unmarried with Children",
       subtitle = "Renters Who Are Not Current On Rent",
       caption = "Source: HH Pulse and COVID State Housing Policy Score"
       ) +
  scale_color_discrete(name="Model 6",
                       breaks=c("1", "0"),
                       labels=c("White Unmarried with Chilren", "Not White Unmarried with Chilren")
                       ) +
  scale_fill_discrete(name="Model 6",
                       breaks=c("1", "0"),
                       labels=c("White Unmarried with Chilren", "Not White Unmarried with Chilren")
                      ) +
  xlab(label = "State Score") +
  ylab(label = "Predicted Probability of High Eviction Risk") #+ ylim(0, 0)
  # png('m6.png')
Code
plot_m6

m6a, white unmarried with children with covariates

Code
m6a <- glmer(highrisk ~ score +
              wunwch +
              age + 
              female + 
              race_eth + 
              inclvl +
              educ +
              (1 | abbv), 
            data = latehh, 
            family = binomial, 
            control = glmerControl(optimizer = "bobyqa"), 
            nAGQ = 10
            )
Warning in checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, : Model is nearly unidentifiable: very large eigenvalue
 - Rescale variables?
Code
summary(m6a)
Generalized linear mixed model fit by maximum likelihood (Adaptive
  Gauss-Hermite Quadrature, nAGQ = 10) [glmerMod]
 Family: binomial  ( logit )
Formula: highrisk ~ score + wunwch + age + female + race_eth + inclvl +  
    educ + (1 | abbv)
   Data: latehh
Control: glmerControl(optimizer = "bobyqa")

     AIC      BIC   logLik deviance df.resid 
 21689.4  21797.5 -10830.7  21661.4    16646 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.5095 -0.8795 -0.5730  1.0117  3.3481 

Random effects:
 Groups Name        Variance Std.Dev.
 abbv   (Intercept) 0.02976  0.1725  
Number of obs: 16660, groups:  abbv, 51

Fixed effects:
                  Estimate Std. Error z value             Pr(>|z|)    
(Intercept)       0.083785   0.088714   0.944              0.34495    
score            -0.070762   0.021746  -3.254              0.00114 ** 
wunwch            0.217350   0.054632   3.978    0.000069360265750 ***
age              -0.005451   0.001277  -4.268    0.000019698930097 ***
female           -0.031255   0.035462  -0.881              0.37811    
race_ethnh_black  0.327082   0.045193   7.237    0.000000000000457 ***
race_ethhispanic -0.073970   0.049384  -1.498              0.13417    
race_ethother     0.151872   0.067188   2.260              0.02380 *  
race_ethnh_asian -1.138936   0.092027 -12.376 < 0.0000000000000002 ***
inclvllessthan25  0.350885   0.036364   9.649 < 0.0000000000000002 ***
inclvl50andover  -0.424999   0.047118  -9.020 < 0.0000000000000002 ***
educHSorless     -0.037086   0.041724  -0.889              0.37409    
educBAorhigher   -0.283860   0.039490  -7.188    0.000000000000657 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Correlation matrix not shown by default, as p = 13 > 12.
Use print(x, correlation=TRUE)  or
    vcov(x)        if you need it
optimizer (bobyqa) convergence code: 0 (OK)
Model is nearly unidentifiable: very large eigenvalue
 - Rescale variables?
Code
plot_m6a <- plot_predictions(m6a, condition = c("score", "wunwch")) +
  labs(title = "Predicted Probability of High Eviction Risk, White Unmarried with Children",
       subtitle = "Renters Who Are Not Current On Rent",
       caption = "Source: HH Pulse and COVID State Housing Policy Score"
       ) +
  scale_color_discrete(name="Model 6a",
                       breaks=c("1", "0"),
                       labels=c("White Unmarried with Chilren", "Not White Unmarried with Chilren")
                       ) +
  scale_fill_discrete(name="Model 6a",
                       breaks=c("1", "0"),
                       labels=c("White Unmarried with Chilren", "Not White Unmarried with Chilren")
                      ) +
  xlab(label = "State Score") +
  ylab(label = "Predicted Probability of High Eviction Risk") #+ ylim(0, 0)
  # png('m6a.png')
Code
plot_m6a

m6b, white unmarried with children, interacted with score

Code
m6b <- glmer(highrisk ~ score +
               age +
               female +
               race_eth +
               inclvl +
               educ +
               score:wunwch +
               (1 | abbv),
             data = latehh,
             family = binomial,
             control = glmerControl(optimizer = "bobyqa"),
             nAGQ = 10
             )
Warning in checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, : Model is nearly unidentifiable: very large eigenvalue
 - Rescale variables?
Code
summary(m6b)
Generalized linear mixed model fit by maximum likelihood (Adaptive
  Gauss-Hermite Quadrature, nAGQ = 10) [glmerMod]
 Family: binomial  ( logit )
Formula: highrisk ~ score + age + female + race_eth + inclvl + educ +  
    score:wunwch + (1 | abbv)
   Data: latehh
Control: glmerControl(optimizer = "bobyqa")

     AIC      BIC   logLik deviance df.resid 
 21688.3  21796.4 -10830.2  21660.3    16646 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.5238 -0.8794 -0.5715  1.0115  3.3713 

Random effects:
 Groups Name        Variance Std.Dev.
 abbv   (Intercept) 0.02868  0.1693  
Number of obs: 16660, groups:  abbv, 51

Fixed effects:
                  Estimate Std. Error z value             Pr(>|z|)    
(Intercept)       0.122236   0.086801   1.408             0.159059    
score            -0.083782   0.021699  -3.861             0.000113 ***
age              -0.005527   0.001275  -4.335    0.000014577389466 ***
female           -0.028912   0.035381  -0.817             0.413834    
race_ethnh_black  0.314830   0.044147   7.131    0.000000000000994 ***
race_ethhispanic -0.083297   0.048643  -1.712             0.086819 .  
race_ethother     0.139549   0.066498   2.099             0.035857 *  
race_ethnh_asian -1.146341   0.091759 -12.493 < 0.0000000000000002 ***
inclvllessthan25  0.352830   0.036364   9.703 < 0.0000000000000002 ***
inclvl50andover  -0.425864   0.047106  -9.041 < 0.0000000000000002 ***
educHSorless     -0.037259   0.041722  -0.893             0.371843    
educBAorhigher   -0.284253   0.039485  -7.199    0.000000000000607 ***
score:wunwch      0.095078   0.023105   4.115    0.000038712321240 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Correlation matrix not shown by default, as p = 13 > 12.
Use print(x, correlation=TRUE)  or
    vcov(x)        if you need it
optimizer (bobyqa) convergence code: 0 (OK)
Model is nearly unidentifiable: very large eigenvalue
 - Rescale variables?
Code
plot_m6b <- plot_predictions(m6b, condition = c("score","wunwch")) +
  labs(title = "Eviction Risk by Score : white unmarried with children",
       subtitle = "Renters Who Are Not Current On Rent",
       caption = "Source: HH Pulse and COVID State Housing Policy Score"
       ) +
  scale_color_discrete(name="Model 6b",
                       breaks=c("1", "0"),
                       labels=c("White Unmarried with Chilren", "Not White Unmarried with Chilren")
                       ) +
  scale_fill_discrete(name="Model 6b",
                       breaks=c("1", "0"),
                       labels=c("White Unmarried with Chilren", "Not White Unmarried with Chilren")
                      ) +
  xlab(label = "State Score") +
  ylab(label = "Predicted Probability of High Risk of Eviction") #+ # ylim(0, 0)
  #png('m6b.png')
Code
plot_m6b

Model 7: White and black, unmarried with children, in same model

m7, unmarried with children

Code
m7 <- glmer(highrisk ~ score +
              wunwch +
              bunwch +
              (1 | abbv), 
            data = latehh, 
            family = binomial, 
            control = glmerControl(optimizer = "bobyqa"), 
            nAGQ = 10
            )

summary(m7)
Generalized linear mixed model fit by maximum likelihood (Adaptive
  Gauss-Hermite Quadrature, nAGQ = 10) [glmerMod]
 Family: binomial  ( logit )
Formula: highrisk ~ score + wunwch + bunwch + (1 | abbv)
   Data: latehh
Control: glmerControl(optimizer = "bobyqa")

     AIC      BIC   logLik deviance df.resid 
 22426.4  22465.0 -11208.2  22416.4    16655 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.3715 -0.8415 -0.7254  1.1032  1.5169 

Random effects:
 Groups Name        Variance Std.Dev.
 abbv   (Intercept) 0.03495  0.1869  
Number of obs: 16660, groups:  abbv, 51

Fixed effects:
            Estimate Std. Error z value             Pr(>|z|)    
(Intercept) -0.23523    0.05389  -4.365      0.0000127269040 ***
score       -0.09187    0.02279  -4.030      0.0000556597003 ***
wunwch       0.32395    0.04933   6.568      0.0000000000511 ***
bunwch       0.64691    0.05618  11.516 < 0.0000000000000002 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Correlation of Fixed Effects:
       (Intr) score  wunwch
score  -0.792              
wunwch -0.156  0.020       
bunwch -0.121  0.013  0.115

m7a, white and black unmarried with children with covariates

Code
m7a <- glmer(highrisk ~ score +
               wunwch +
              bunwch +
              age + 
              female + 
              race_eth + 
              inclvl +
              educ +
              (1 | abbv), 
            data = latehh, 
            family = binomial, 
            control = glmerControl(optimizer = "bobyqa"), 
            nAGQ = 10
            )
Warning in checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, :
Model failed to converge with max|grad| = 0.00657554 (tol = 0.002, component 1)
Warning in checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, : Model is nearly unidentifiable: very large eigenvalue
 - Rescale variables?
Code
summary(m7a)
Generalized linear mixed model fit by maximum likelihood (Adaptive
  Gauss-Hermite Quadrature, nAGQ = 10) [glmerMod]
 Family: binomial  ( logit )
Formula: highrisk ~ score + wunwch + bunwch + age + female + race_eth +  
    inclvl + educ + (1 | abbv)
   Data: latehh
Control: glmerControl(optimizer = "bobyqa")

     AIC      BIC   logLik deviance df.resid 
 21671.2  21787.0 -10820.6  21641.2    16645 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.6296 -0.8758 -0.5736  1.0135  3.2736 

Random effects:
 Groups Name        Variance Std.Dev.
 abbv   (Intercept) 0.02953  0.1718  
Number of obs: 16660, groups:  abbv, 51

Fixed effects:
                  Estimate Std. Error z value             Pr(>|z|)    
(Intercept)       0.053300   0.088904   0.600             0.548824    
score            -0.069638   0.021696  -3.210             0.001329 ** 
wunwch            0.225751   0.054645   4.131     0.00003608037258 ***
bunwch            0.308080   0.068757   4.481     0.00000743905262 ***
age              -0.004683   0.001288  -3.636             0.000277 ***
female           -0.046210   0.035624  -1.297             0.194582    
race_ethnh_black  0.210364   0.052202   4.030     0.00005581749894 ***
race_ethhispanic -0.070482   0.049364  -1.428             0.153351    
race_ethother     0.157158   0.067162   2.340             0.019284 *  
race_ethnh_asian -1.142930   0.091998 -12.423 < 0.0000000000000002 ***
inclvllessthan25  0.346748   0.036402   9.525 < 0.0000000000000002 ***
inclvl50andover  -0.419830   0.047142  -8.906 < 0.0000000000000002 ***
educHSorless     -0.037413   0.041761  -0.896             0.370320    
educBAorhigher   -0.279826   0.039525  -7.080     0.00000000000144 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Correlation matrix not shown by default, as p = 14 > 12.
Use print(x, correlation=TRUE)  or
    vcov(x)        if you need it
optimizer (bobyqa) convergence code: 0 (OK)
Model failed to converge with max|grad| = 0.00657554 (tol = 0.002, component 1)
Model is nearly unidentifiable: very large eigenvalue
 - Rescale variables?

Model 8: Combination of married and children

Make dummy variable for dependent variable

Code
latehh <- latehh %>%
  mutate(marrchild = case_when(married == 0 & children == 1 ~ "UnWCh",
                               married == 1 & children == 1 ~ "MWCh",
                               married == 0 & children == 0 ~ "UnNoCh",
                               married == 1 & children == 0 ~ "MNoCh"
                              )
        )

tabyl(latehh$marrchild)
 latehh$marrchild    n   percent
            MNoCh 2111 0.1267107
             MWCh 2863 0.1718487
           UnNoCh 6684 0.4012005
            UnWCh 5002 0.3002401
Code
latehh$marrchild <- factor(latehh$marrchild, levels=c("UnNoCh",  "UnWCh",  "MNoCh", "MWCh"))

tabyl(latehh$marrchild)
 latehh$marrchild    n   percent
           UnNoCh 6684 0.4012005
            UnWCh 5002 0.3002401
            MNoCh 2111 0.1267107
             MWCh 2863 0.1718487

m8, marrchild dependent variable with no covariates

Code
m8 <- glmer(highrisk ~ score +
               marrchild +
              (1 | abbv), 
            data = latehh, 
            family = binomial, 
            control = glmerControl(optimizer = "bobyqa"), 
            nAGQ = 10
            )

summary(m8)
Generalized linear mixed model fit by maximum likelihood (Adaptive
  Gauss-Hermite Quadrature, nAGQ = 10) [glmerMod]
 Family: binomial  ( logit )
Formula: highrisk ~ score + marrchild + (1 | abbv)
   Data: latehh
Control: glmerControl(optimizer = "bobyqa")

     AIC      BIC   logLik deviance df.resid 
 22355.4  22401.7 -11171.7  22343.4    16654 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.2456 -0.8466 -0.7045  1.0669  1.7090 

Random effects:
 Groups Name        Variance Std.Dev.
 abbv   (Intercept) 0.03826  0.1956  
Number of obs: 16660, groups:  abbv, 51

Fixed effects:
               Estimate Std. Error z value             Pr(>|z|)    
(Intercept)    -0.14142    0.05831  -2.425               0.0153 *  
score          -0.09507    0.02355  -4.037       0.000054113154 ***
marrchildUnWCh  0.32351    0.03801   8.510 < 0.0000000000000002 ***
marrchildMNoCh -0.34172    0.05254  -6.504       0.000000000078 ***
marrchildMWCh  -0.27882    0.04660  -5.984       0.000000002184 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Correlation of Fixed Effects:
            (Intr) score  mrrUWC mrrMNC
score       -0.759                     
mrrchldUnWC -0.302  0.020              
mrrchldMNCh -0.203  0.004  0.312       
mrrchldMWCh -0.242  0.019  0.354  0.257
Code
plot_m8 <- plot_predictions(m8, condition = c("score", "marrchild")) +
  labs(title = "Predicted Probability of High Eviction Risk, Model 8",
       subtitle = "Renters Who Are Not Current On Rent",
       caption = "Source: HH Pulse and COVID State Housing Policy Score"
       ) +
  scale_color_discrete(name="Model 8",
                       breaks=c("UnWCh", "UnNoCh", "MWCh", "MNoCh"),
                       labels=c("Unmarried, With Chilren", "Unmarried, No Children", "Married, With Children", "Married, No Children")
                       ) +
  scale_fill_discrete(name="Model 8",
                       breaks=c("UnWCh", "UnNoCh", "MWCh", "MNoCh"),
                       labels=c("Unmarried, With Chilren", "Unmarried, No Children", "Married, With Children", "Married, No Children")
                       ) +
  xlab(label = "State Score") +
  ylab(label = "Predicted Probability of High Eviction Risk") #+ ylim(0, 0)
  # png('m8.png')
Code
plot_m8

m8a, marrchild dependent variable with covariates

Code
m8a <- glmer(highrisk ~ score +
              marrchild +
              age + 
              female + 
              race_eth + 
              inclvl +
              educ +
              (1 | abbv), 
            data = latehh, 
            family = binomial, 
            control = glmerControl(optimizer = "bobyqa"), 
            nAGQ = 10
            )
Warning in checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, :
Model failed to converge with max|grad| = 0.0137291 (tol = 0.002, component 1)
Warning in checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, : Model is nearly unidentifiable: very large eigenvalue
 - Rescale variables?
Code
summary(m8a)
Generalized linear mixed model fit by maximum likelihood (Adaptive
  Gauss-Hermite Quadrature, nAGQ = 10) [glmerMod]
 Family: binomial  ( logit )
Formula: highrisk ~ score + marrchild + age + female + race_eth + inclvl +  
    educ + (1 | abbv)
   Data: latehh
Control: glmerControl(optimizer = "bobyqa")

     AIC      BIC   logLik deviance df.resid 
 21648.3  21771.8 -10808.2  21616.3    16644 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.6090 -0.8727 -0.5701  1.0071  3.2785 

Random effects:
 Groups Name        Variance Std.Dev.
 abbv   (Intercept) 0.03023  0.1739  
Number of obs: 16660, groups:  abbv, 51

Fixed effects:
                  Estimate Std. Error z value             Pr(>|z|)    
(Intercept)       0.077299   0.091231   0.847             0.396837    
score            -0.070836   0.021877  -3.238             0.001204 ** 
marrchildUnWCh    0.237051   0.040683   5.827     0.00000000565054 ***
marrchildMNoCh   -0.107856   0.054813  -1.968             0.049099 *  
marrchildMWCh    -0.086955   0.049601  -1.753             0.079587 .  
age              -0.004392   0.001313  -3.344             0.000824 ***
female           -0.061470   0.035956  -1.710             0.087343 .  
race_ethnh_black  0.237230   0.042918   5.528     0.00000003247917 ***
race_ethhispanic -0.144175   0.047500  -3.035             0.002403 ** 
race_ethother     0.076241   0.065662   1.161             0.245599    
race_ethnh_asian -1.153975   0.091561 -12.603 < 0.0000000000000002 ***
inclvllessthan25  0.331822   0.036731   9.034 < 0.0000000000000002 ***
inclvl50andover  -0.401560   0.047481  -8.457 < 0.0000000000000002 ***
educHSorless     -0.037271   0.041863  -0.890             0.373306    
educBAorhigher   -0.276365   0.039570  -6.984     0.00000000000287 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Correlation matrix not shown by default, as p = 15 > 12.
Use print(x, correlation=TRUE)  or
    vcov(x)        if you need it
optimizer (bobyqa) convergence code: 0 (OK)
Model failed to converge with max|grad| = 0.0137291 (tol = 0.002, component 1)
Model is nearly unidentifiable: very large eigenvalue
 - Rescale variables?
Code
plot_m8a <- plot_predictions(m8a, condition = c("score", "marrchild")) +
  labs(title = "Predicted Probability of High Eviction Risk, Model 8",
       subtitle = "Renters Who Are Not Current On Rent",
       caption = "Source: HH Pulse and COVID State Housing Policy Score"
       ) +
  scale_color_discrete(name="Model 8",
                       breaks=c("UnWCh", "UnNoCh", "MWCh", "MNoCh"),
                       labels=c("Unmarried, With Chilren", "Unmarried, No Children", "Married, With Children", "Married, No Children")
                       ) +
  scale_fill_discrete(name="Model 8",
                       breaks=c("UnWCh", "UnNoCh", "MWCh", "MNoCh"),
                       labels=c("Unmarried, With Chilren", "Unmarried, No Children", "Married, With Children", "Married, No Children")
                       ) +
  xlab(label = "State Score") +
  ylab(label = "Predicted Probability of High Eviction Risk") #+ ylim(0, 0)
  # png('m8a.png')
Code
plot_m8a

m8b, marrchild dependent variable with covariates

Code
m8b <- glmer(highrisk ~ score +
               age +
               female +
               race_eth +
               inclvl +
               educ +
               score:marrchild +
               (1 | abbv),
             data = latehh,
             family = binomial,
             control = glmerControl(optimizer = "bobyqa"),
             nAGQ = 10
             )
Warning in checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, :
Model failed to converge with max|grad| = 0.0181597 (tol = 0.002, component 1)
Warning in checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, : Model is nearly unidentifiable: very large eigenvalue
 - Rescale variables?
Code
summary(m8b)
Generalized linear mixed model fit by maximum likelihood (Adaptive
  Gauss-Hermite Quadrature, nAGQ = 10) [glmerMod]
 Family: binomial  ( logit )
Formula: highrisk ~ score + age + female + race_eth + inclvl + educ +  
    score:marrchild + (1 | abbv)
   Data: latehh
Control: glmerControl(optimizer = "bobyqa")

     AIC      BIC   logLik deviance df.resid 
 21661.2  21784.7 -10814.6  21629.2    16644 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.5739 -0.8796 -0.5662  1.0076  3.3643 

Random effects:
 Groups Name        Variance Std.Dev.
 abbv   (Intercept) 0.02979  0.1726  
Number of obs: 16660, groups:  abbv, 51

Fixed effects:
                      Estimate Std. Error z value             Pr(>|z|)    
(Intercept)           0.136591   0.087070   1.569             0.116705    
score                -0.093532   0.023240  -4.025     0.00005708323742 ***
age                  -0.004803   0.001297  -3.702             0.000214 ***
female               -0.049369   0.035736  -1.381             0.167135    
race_ethnh_black      0.248031   0.042785   5.797     0.00000000674547 ***
race_ethhispanic     -0.143213   0.047431  -3.019             0.002533 ** 
race_ethother         0.083171   0.065579   1.268             0.204707    
race_ethnh_asian     -1.162037   0.091526 -12.696 < 0.0000000000000002 ***
inclvllessthan25      0.339000   0.036598   9.263 < 0.0000000000000002 ***
inclvl50andover      -0.410764   0.047354  -8.674 < 0.0000000000000002 ***
educHSorless         -0.039168   0.041825  -0.936             0.349034    
educBAorhigher       -0.278673   0.039538  -7.048     0.00000000000181 ***
score:marrchildUnWCh  0.095549   0.017341   5.510     0.00000003590281 ***
score:marrchildMNoCh -0.032294   0.023699  -1.363             0.172983    
score:marrchildMWCh  -0.029227   0.022154  -1.319             0.187076    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Correlation matrix not shown by default, as p = 15 > 12.
Use print(x, correlation=TRUE)  or
    vcov(x)        if you need it
optimizer (bobyqa) convergence code: 0 (OK)
Model failed to converge with max|grad| = 0.0181597 (tol = 0.002, component 1)
Model is nearly unidentifiable: very large eigenvalue
 - Rescale variables?
Code
plot_m8b <- plot_predictions(m8b, condition = c("score", "marrchild")) +
  labs(title = "Predicted Probability of High Eviction Risk, Model 8b",
       subtitle = "Renters Who Are Not Current On Rent",
       caption = "Source: HH Pulse and COVID State Housing Policy Score"
       ) +
  scale_color_discrete(name="Model 8b",
                       breaks=c("UnWCh", "UnNoCh", "MWCh", "MNoCh"),
                       labels=c("Unmarried, With Chilren", "Unmarried, No Children", "Married, With Children", "Married, No Children")
                       ) +
  scale_fill_discrete(name="Model 8b",
                       breaks=c("UnWCh", "UnNoCh", "MWCh", "MNoCh"),
                       labels=c("Unmarried, With Chilren", "Unmarried, No Children", "Married, With Children", "Married, No Children")
                       ) +
  xlab(label = "State Score") +
  ylab(label = "Predicted Probability of High Eviction Risk") #+ ylim(0, 0)
  # png('m8b.png')
Code
plot_m8b

Model 9: Sample of only unmarried with children, then test by race

Code
tabyl(latehh, married, children)
 married    0    1
       0 6684 5002
       1 2111 2863

Filter for unmarried with children late households

Code
umwc_latehh <- latehh %>% 
  filter(married == 0 & children == 1)

nrow(umwc_latehh)
[1] 5002

Write out csv for unmarried with children late households

Code
# Export to csv files

fwrite(umwc_latehh, "./data/umwc_latehh.csv")
Code
tabyl(umwc_latehh$race_eth)
 umwc_latehh$race_eth    n    percent
             nh_white 1979 0.39564174
             nh_black 1530 0.30587765
             hispanic  982 0.19632147
                other  405 0.08096761
             nh_asian  106 0.02119152

m9, race_eth with no covariages, umwc sample

Code
m9 <- glmer(highrisk ~ score +
               race_eth +
              (1 | abbv), 
            data = umwc_latehh, 
            family = binomial, 
            control = glmerControl(optimizer = "bobyqa"), 
            nAGQ = 10
            )

summary(m9)
Generalized linear mixed model fit by maximum likelihood (Adaptive
  Gauss-Hermite Quadrature, nAGQ = 10) [glmerMod]
 Family: binomial  ( logit )
Formula: highrisk ~ score + race_eth + (1 | abbv)
   Data: umwc_latehh
Control: glmerControl(optimizer = "bobyqa")

     AIC      BIC   logLik deviance df.resid 
  6885.1   6930.7  -3435.5   6871.1     4995 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.3026 -0.9757  0.7677  0.9843  1.5148 

Random effects:
 Groups Name        Variance Std.Dev.
 abbv   (Intercept) 0.02958  0.172   
Number of obs: 5002, groups:  abbv, 51

Fixed effects:
                 Estimate Std. Error z value   Pr(>|z|)    
(Intercept)       0.01566    0.07350   0.213     0.8312    
score            -0.05881    0.02826  -2.081     0.0374 *  
race_ethnh_black  0.34159    0.07142   4.783 0.00000173 ***
race_ethhispanic -0.06422    0.08157  -0.787     0.4311    
race_ethother     0.14653    0.11099   1.320     0.1868    
race_ethnh_asian -0.51212    0.21056  -2.432     0.0150 *  
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Correlation of Fixed Effects:
            (Intr) score  rc_thnh_b rc_thh rc_tht
score       -0.704                               
rc_thnh_blc -0.403  0.004                        
rac_thhspnc -0.303 -0.049  0.369                 
race_eththr -0.268  0.015  0.264     0.233       
rac_thnh_sn -0.109 -0.030  0.142     0.151  0.094
Code
plot_m9 <- plot_predictions(m9, condition = c("score", "race_eth")) +
  labs(title = "Predicted Probability of High Eviction Risk, Model 9",
       subtitle = "Renters Who Are Not Current On Rent",
       caption = "Source: HH Pulse and COVID State Housing Policy Score"
       ) +
  scale_color_discrete(name="Model 9",
                       breaks=c("nh_black", "other", "nh_white", "hispanic", "nh_asian"),
                       labels=c("Non-Hispanic Black", "Other", "Non-Hispanic White", "Hispanic", "Non-Hispanic Asian")
                       ) +
  scale_fill_discrete(name="Model 9",
                      breaks=c("nh_black", "other", "nh_white", "hispanic", "nh_asian"),
                      labels=c("Non-Hispanic Black", "Other", "Non-Hispanic White", "Hispanic", "Non-Hispanic Asian")
                      ) +
  xlab(label = "State Score") +
  ylab(label = "Predicted Probability of High Eviction Risk") #+ ylim(0, 0)
  # png('m9.png')
Code
plot_m9

m9a, race_eth with covariages, umwc sample

Code
m9a <- glmer(highrisk ~ score +
              race_eth +
              age + 
              female + 
              race_eth + 
              inclvl +
              educ +
              (1 | abbv), 
            data = umwc_latehh, 
            family = binomial, 
            control = glmerControl(optimizer = "bobyqa"), 
            nAGQ = 10
            )
Warning in checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, :
Model failed to converge with max|grad| = 0.0230086 (tol = 0.002, component 1)
Warning in checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, : Model is nearly unidentifiable: very large eigenvalue
 - Rescale variables?
Code
summary(m9a)
Generalized linear mixed model fit by maximum likelihood (Adaptive
  Gauss-Hermite Quadrature, nAGQ = 10) [glmerMod]
 Family: binomial  ( logit )
Formula: highrisk ~ score + race_eth + age + female + race_eth + inclvl +  
    educ + (1 | abbv)
   Data: umwc_latehh
Control: glmerControl(optimizer = "bobyqa")

     AIC      BIC   logLik deviance df.resid 
  6816.5   6901.2  -3395.2   6790.5     4989 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.7416 -0.9762  0.6614  0.9553  1.7411 

Random effects:
 Groups Name        Variance Std.Dev.
 abbv   (Intercept) 0.02635  0.1623  
Number of obs: 5002, groups:  abbv, 51

Fixed effects:
                  Estimate Std. Error z value      Pr(>|z|)    
(Intercept)       0.488812   0.158742   3.079       0.00207 ** 
score            -0.052242   0.027794  -1.880       0.06017 .  
race_ethnh_black  0.322464   0.071996   4.479 0.00000750351 ***
race_ethhispanic -0.115321   0.083088  -1.388       0.16516    
race_ethother     0.111075   0.111854   0.993       0.32069    
race_ethnh_asian -0.469729   0.213317  -2.202       0.02766 *  
age              -0.008185   0.002709  -3.021       0.00252 ** 
female           -0.240267   0.077308  -3.108       0.00188 ** 
inclvllessthan25  0.382573   0.062934   6.079 0.00000000121 ***
inclvl50andover  -0.144421   0.095545  -1.512       0.13065    
educHSorless     -0.116092   0.070436  -1.648       0.09931 .  
educBAorhigher   -0.199756   0.072947  -2.738       0.00617 ** 
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Correlation of Fixed Effects:
            (Intr) score  rc_thnh_b rc_thh rc_tht rc_thnh_s age    female
score       -0.317                                                       
rc_thnh_blc -0.210  0.003                                                
rac_thhspnc -0.187 -0.050  0.368                                         
race_eththr -0.123  0.015  0.264     0.232                               
rac_thnh_sn -0.057 -0.030  0.137     0.143  0.094                        
age         -0.724 -0.006  0.052     0.061  0.004  0.004                 
female      -0.480  0.007 -0.036     0.041  0.005  0.047     0.102       
inclvllss25 -0.147  0.016  0.003    -0.005 -0.035 -0.021    -0.045 -0.064
inclvl50ndv -0.136 -0.020  0.035     0.035  0.003 -0.021    -0.036  0.090
educHSorlss -0.192 -0.023 -0.010    -0.102  0.002  0.006    -0.007  0.055
educBArhghr -0.100 -0.005 -0.012     0.001  0.000 -0.069    -0.138 -0.001
            incl25 incl50 edcHSr
score                           
rc_thnh_blc                     
rac_thhspnc                     
race_eththr                     
rac_thnh_sn                     
age                             
female                          
inclvllss25                     
inclvl50ndv  0.338              
educHSorlss -0.124  0.019       
educBArhghr  0.013 -0.121  0.450
optimizer (bobyqa) convergence code: 0 (OK)
Model failed to converge with max|grad| = 0.0230086 (tol = 0.002, component 1)
Model is nearly unidentifiable: very large eigenvalue
 - Rescale variables?
Code
plot_m9a <- plot_predictions(m9a, condition = c("score", "race_eth")) +
  labs(title = "Predicted Probability of High Eviction Risk, Model 9a",
       subtitle = "Renters Who Are Not Current On Rent",
       caption = "Source: HH Pulse and COVID State Housing Policy Score"
       ) +
  scale_color_discrete(name="Model 9a",
                       breaks=c("nh_black", "other", "nh_white", "hispanic", "nh_asian"),
                       labels=c("Non-Hispanic Black", "Other", "Non-Hispanic White", "Hispanic", "Non-Hispanic Asian")
                       ) +
  scale_fill_discrete(name="Model 9a",
                      breaks=c("nh_black", "other", "nh_white", "hispanic", "nh_asian"),
                      labels=c("Non-Hispanic Black", "Other", "Non-Hispanic White", "Hispanic", "Non-Hispanic Asian")
                      ) +
  xlab(label = "State Score") +
  ylab(label = "Predicted Probability of High Eviction Risk") #+ ylim(0, 0)
  # png('m9a.png')
Code
plot_m9a

m9b, race_eth with covariages and interaction, umwc sample

Code
m9b <- glmer(highrisk ~ score +
               age +
               female +
               inclvl +
               educ +
               score:race_eth +
               (1 | abbv),
             data = umwc_latehh,
             family = binomial,
             control = glmerControl(optimizer = "bobyqa"),
             nAGQ = 10
             )
Warning in checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, : Model is nearly unidentifiable: very large eigenvalue
 - Rescale variables?
Code
summary(m9b)
Generalized linear mixed model fit by maximum likelihood (Adaptive
  Gauss-Hermite Quadrature, nAGQ = 10) [glmerMod]
 Family: binomial  ( logit )
Formula: highrisk ~ score + age + female + inclvl + educ + score:race_eth +  
    (1 | abbv)
   Data: umwc_latehh
Control: glmerControl(optimizer = "bobyqa")

     AIC      BIC   logLik deviance df.resid 
  6842.4   6927.2  -3408.2   6816.4     4989 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.5995 -0.9817  0.7031  0.9541  1.8331 

Random effects:
 Groups Name        Variance Std.Dev.
 abbv   (Intercept) 0.03169  0.178   
Number of obs: 5002, groups:  abbv, 51

Fixed effects:
                        Estimate Std. Error z value      Pr(>|z|)    
(Intercept)             0.577444   0.154705   3.733       0.00019 ***
score                  -0.074469   0.032962  -2.259       0.02387 *  
age                    -0.008470   0.002703  -3.133       0.00173 ** 
female                 -0.218110   0.077009  -2.832       0.00462 ** 
inclvllessthan25        0.380297   0.062861   6.050 0.00000000145 ***
inclvl50andover        -0.155998   0.095382  -1.636       0.10194    
educHSorless           -0.125276   0.070197  -1.785       0.07432 .  
educBAorhigher         -0.205633   0.072750  -2.827       0.00471 ** 
score:race_ethnh_black  0.077916   0.032350   2.409       0.01602 *  
score:race_ethhispanic -0.033008   0.035170  -0.939       0.34798    
score:race_ethother     0.014547   0.051048   0.285       0.77566    
score:race_ethnh_asian -0.104550   0.086746  -1.205       0.22811    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Correlation of Fixed Effects:
              (Intr) score  age    female incl25 incl50 edcHSr edcBAr
score         -0.288                                                 
age           -0.729 -0.025                                          
female        -0.492  0.002  0.103                                   
inclvllss25   -0.153  0.030 -0.046 -0.063                            
inclvl50ndv   -0.132 -0.030 -0.037  0.092  0.337                     
educHSorlss   -0.207 -0.004 -0.006  0.056 -0.124  0.019              
educBArhghr   -0.105  0.003 -0.139  0.000  0.013 -0.122  0.452       
scr:rc_thnh_b -0.016 -0.412  0.037 -0.015 -0.020  0.022 -0.006 -0.012
scr:rc_thhs   -0.036 -0.367  0.053  0.040 -0.018  0.038 -0.087 -0.001
scr:rc_thth    0.011 -0.251 -0.007  0.009 -0.050  0.007  0.004 -0.006
scr:rc_thnh_s  0.002 -0.151  0.005  0.029 -0.027 -0.023  0.007 -0.057
              scr:rc_thnh_b scr:rc_thh scr:rc_tht
score                                            
age                                              
female                                           
inclvllss25                                      
inclvl50ndv                                      
educHSorlss                                      
educBArhghr                                      
scr:rc_thnh_b                                    
scr:rc_thhs    0.379                             
scr:rc_thth    0.255         0.232               
scr:rc_thnh_s  0.150         0.143      0.106    
optimizer (bobyqa) convergence code: 0 (OK)
Model is nearly unidentifiable: very large eigenvalue
 - Rescale variables?
Code
plot_m9b <- plot_predictions(m9b, condition = c("score", "race_eth")) +
  labs(title = "Predicted Probability of High Eviction Risk",
       subtitle = "Unmarried, Children in Household, Who Are Not Current On Rent (Model 9b)",
       caption = "Source: HH Pulse and COVID State Housing Policy Score"
       ) +
  scale_color_discrete(name="Model 9b",
                       breaks=c("nh_black", "other", "nh_white", "hispanic", "nh_asian"),
                       labels=c("Non-Hispanic Black", "Other", "Non-Hispanic White", "Hispanic", "Non-Hispanic Asian")
                       ) +
  scale_fill_discrete(name="Model 9b",
                      breaks=c("nh_black", "other", "nh_white", "hispanic", "nh_asian"),
                      labels=c("Non-Hispanic Black", "Other", "Non-Hispanic White", "Hispanic", "Non-Hispanic Asian")
                      ) +
  xlab(label = "State Score") +
  ylab(label = "Predicted Probability of High Eviction Risk") #+ ylim(0, 0)
  # png('m9b.png')
Code
plot_m9b

m9c, race_eth with covariages and interaction, umwc sample

Code
m9c <- glmer(highrisk ~ 
               score*(race_eth +
                        age +
                        female +
                        inclvl +
                        educ
                        ) +
               (1 | abbv),
             data = umwc_latehh,
             family = binomial,
             control = glmerControl(optimizer = "bobyqa"),
             nAGQ = 10
             )
Warning in checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, :
Model failed to converge with max|grad| = 0.0696614 (tol = 0.002, component 1)
Warning in checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, : Model is nearly unidentifiable: very large eigenvalue
 - Rescale variables?;Model is nearly unidentifiable: large eigenvalue ratio
 - Rescale variables?
Code
summary(m9c)
Generalized linear mixed model fit by maximum likelihood (Adaptive
  Gauss-Hermite Quadrature, nAGQ = 10) [glmerMod]
 Family: binomial  ( logit )
Formula: highrisk ~ score * (race_eth + age + female + inclvl + educ) +  
    (1 | abbv)
   Data: umwc_latehh
Control: glmerControl(optimizer = "bobyqa")

     AIC      BIC   logLik deviance df.resid 
  6818.1   6968.0  -3386.1   6772.1     4979 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.8541 -0.9749  0.6462  0.9425  1.6173 

Random effects:
 Groups Name        Variance Std.Dev.
 abbv   (Intercept) 0.02307  0.1519  
Number of obs: 5002, groups:  abbv, 51

Fixed effects:
                        Estimate Std. Error z value  Pr(>|z|)    
(Intercept)             0.785708   0.258003   3.045   0.00232 ** 
score                  -0.212282   0.113235  -1.875   0.06083 .  
race_ethnh_black        0.512516   0.121084   4.233 0.0000231 ***
race_ethhispanic       -0.178227   0.152456  -1.169   0.24239    
race_ethother           0.243760   0.182057   1.339   0.18060    
race_ethnh_asian       -1.002961   0.434774  -2.307   0.02106 *  
age                    -0.011476   0.004697  -2.443   0.01455 *  
female                 -0.282187   0.133707  -2.110   0.03482 *  
inclvllessthan25        0.155952   0.107323   1.453   0.14619    
inclvl50andover        -0.272201   0.168888  -1.612   0.10702    
educHSorless           -0.195295   0.120691  -1.618   0.10563    
educBAorhigher         -0.292957   0.124259  -2.358   0.01839 *  
score:race_ethnh_black -0.110429   0.054240  -2.036   0.04176 *  
score:race_ethhispanic  0.024189   0.064691   0.374   0.70847    
score:race_ethother    -0.082584   0.083069  -0.994   0.32015    
score:race_ethnh_asian  0.236387   0.173393   1.363   0.17279    
score:age               0.001807   0.002029   0.891   0.37313    
score:female            0.020783   0.057803   0.360   0.71919    
score:inclvllessthan25  0.125572   0.047401   2.649   0.00807 ** 
score:inclvl50andover   0.070778   0.072139   0.981   0.32653    
score:educHSorless      0.045472   0.053217   0.854   0.39284    
score:educBAorhigher    0.054533   0.055029   0.991   0.32169    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Correlation matrix not shown by default, as p = 22 > 12.
Use print(x, correlation=TRUE)  or
    vcov(x)        if you need it
optimizer (bobyqa) convergence code: 0 (OK)
Model failed to converge with max|grad| = 0.0696614 (tol = 0.002, component 1)
Model is nearly unidentifiable: very large eigenvalue
 - Rescale variables?
Model is nearly unidentifiable: large eigenvalue ratio
 - Rescale variables?
Code
plot_m9c <- plot_predictions(m9c, condition = c("score", "race_eth")) +
  labs(title = "Predicted Probability of High Eviction Risk",
       subtitle = "Unmarried, Children in Household, Who Are Not Current On Rent (Model 9b)",
       caption = "Source: HH Pulse and COVID State Housing Policy Score"
       ) +
  scale_color_discrete(name="Model 9c",
                       breaks=c("nh_black", "other", "nh_white", "hispanic", "nh_asian"),
                       labels=c("Non-Hispanic Black", "Other", "Non-Hispanic White", "Hispanic", "Non-Hispanic Asian")
                       ) +
  scale_fill_discrete(name="Model 9c",
                      breaks=c("nh_black", "other", "nh_white", "hispanic", "nh_asian"),
                      labels=c("Non-Hispanic Black", "Other", "Non-Hispanic White", "Hispanic", "Non-Hispanic Asian")
                      ) +
  xlab(label = "State Score") +
  ylab(label = "Predicted Probability of High Eviction Risk") #+ ylim(0, 0)
  # png('m9c.png')

plot_m9c

Unused code

m5: unmarried x children interaction

latehh <- latehh %>% mutate(inclvl = case_when(income == 1 ~ “lessthan25”, income %in% c(2,3) ~ “25to49”, TRUE ~ “50andover” ) )

latehhinclvl <- factor(latehhinclvl, levels=c(“25to49”, “lessthan25”, “50andover”))

tabyl(latehh$inclvl)

latehh <- latehh %>% mutate(educ = case_when(eeduc %in% c(1:3) ~ “HSorless”, eeduc == 4 ~ “somecollege”, TRUE ~ “BAorhigher” ) )

latehheduc <- factor(latehheduc, levels=c(“somecollege”, “HSorless”, “BAorhigher”))

tabyl(latehh$educ)

m5 <- glmer(highrisk ~ score + age + female + race_eth + inclvl + educ + unmarried * children + (1 | abbv), data = latehh, family = binomial, control = glmerControl(optimizer = “bobyqa”), nAGQ = 10 )

summary(m5)

Plot of Eviction Risk by Unmarried x Children

How do you plot an interaction in this scenario?

plot_uxc <- plot_predictions(m5, condition = c(“score”,“children”)) + labs(title = “Eviction Risk by Unmarried x Children”, subtitle = “Renters Who Are Not Current On Rent”, caption = “Source: HH Pulse and COVID State Housing Policy Score” ) + scale_color_discrete(name=“Unmarried x Children”, breaks=c(1,0), labels=c(“Unmarried with Children”, “Married or No Children”)) + scale_fill_discrete(name=“Unmarried x Children”, breaks=c(1,0), labels=c(“Unmarried with Children”, “Married or No Children”)) + xlab(label = “State Score”) + ylab(label = “Predicted Probability of High Risk of Eviction”) #+ # ylim(0, 0) #png(‘plotuxc.png’)

plot_uxc