library(car)
## Loading required package: carData
library(mice)
## 
## Attaching package: 'mice'
## The following object is masked from 'package:stats':
## 
##     filter
## The following objects are masked from 'package:base':
## 
##     cbind, rbind
library(ggplot2)
library(dplyr)
## 
## Attaching package: 'dplyr'
## The following object is masked from 'package:car':
## 
##     recode
## The following objects are masked from 'package:stats':
## 
##     filter, lag
## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union
library(haven)
NSDUH_2019 <- read_sav("NSDUH_2019.SAV")
View(NSDUH_2019)
nams<-names(NSDUH_2019)
head(nams, n=10)
##  [1] "QUESTID2" "FILEDATE" "CIGEVER"  "CIGOFRSM" "CIGWILYR" "CIGTRY"  
##  [7] "CIGYFU"   "CIGMFU"   "CIGREC"   "CIG30USE"
newnames<-tolower(gsub(pattern = "_",replacement =  "",x =  nams))
names(NSDUH_2019)<-newnames
NSDUH_2019$attempt_suicide<-Recode(NSDUH_2019$adwrsatp, recodes="1=1; 2=0;else=NA")

## Main variable hopeless last thirty days
NSDUH_2019$hopelesslst30days<-Recode(NSDUH_2019$dsthop30,
                               recodes="1:2=3;3:4=2;5=1; else=NA")
NSDUH_2019$fltnervouslst30days<-Recode(NSDUH_2019$dstnrv30,
                               recodes="1:2=3;3:4=2;5=1; else=NA")
NSDUH_2019$fltrestlesslst30days<-Recode(NSDUH_2019$dstrst30,
                               recodes="1:2=3;3:4=2;5=1; else=NA")
NSDUH_2019$fltsadlst30days<-Recode(NSDUH_2019$dstchr30,
                               recodes="1:2=3;3:4=2;5=1; else=NA")
NSDUH_2019$effortlst30days<-Recode(NSDUH_2019$dsteff30,
                               recodes="1:2=3;3:4=2;5=1; else=NA")
NSDUH_2019$fltdwnlst30days<-Recode(NSDUH_2019$dstngd30,
                               recodes="1:2=3;3:4=2;5=1; else=NA")
## marital status
NSDUH_2019$marst<-Recode(NSDUH_2019$irmarit, recodes="1='married'; 2='divorced'; 3='widowed'; 4='separated'; else=NA", as.factor=T)
NSDUH_2019$marst<-relevel(NSDUH_2019$marst, ref='married')

## education recodes
NSDUH_2019$educ<-Recode(NSDUH_2019$ireduhighst2, recodes="1:7='LssThnHgh'; 8='highschool'; 9='someCollege'; 10='associates'; 11='colgrad';else=NA", as.factor=T)
NSDUH_2019$educ<-relevel(NSDUH_2019$educ, ref='colgrad')

## sexuality recodes
NSDUH_2019$sexuality<-Recode(NSDUH_2019$sexident, recodes="1='Heterosexual'; 2='Les/Gay'; 3='Bisexual';else=NA", as.factor=T)
NSDUH_2019$sexuality<-relevel(NSDUH_2019$sexuality, ref='Heterosexual')

## gender recodes
NSDUH_2019$male<-as.factor(ifelse(NSDUH_2019$irsex==1, "Male", "Female"))

## Race recoded items
NSDUH_2019$black<-Recode(NSDUH_2019$newrace2, recodes="2=1; 9=NA; else=0")
NSDUH_2019$white<-Recode(NSDUH_2019$newrace2, recodes="1=1; 9=NA; else=0")
NSDUH_2019$other<-Recode(NSDUH_2019$newrace2, recodes="3:4=1; 9=NA; else=0")
NSDUH_2019$mult_race<-Recode(NSDUH_2019$newrace2, recodes="6=1; 9=NA; else=0")
NSDUH_2019$asian<-Recode(NSDUH_2019$newrace2, recodes="5=1; 9=NA; else=0")
NSDUH_2019$hispanic<-Recode(NSDUH_2019$newrace2, recodes="7=1; 9=NA; else=0")
NSDUH_2019$race_eth<-Recode(NSDUH_2019$newrace2,
                          recodes="1='white'; 2='black'; 3='other'; 4='asian'; 5='mult_race'; 6='hispanic'; else=NA",
                          as.factor = T)
NSDUH_2019$race_eth<-relevel(NSDUH_2019$race_eth, ref='white')

NSDUH_2019$dep_year2<-Recode(NSDUH_2019$amdeyr, recodes="1=1; 2=0;else=NA")

NSDUH_2019$daysalc<-Recode(NSDUH_2019$alcdays, recodes = "85=NA; 91=NA;  93=NA; 94=NA; 97=NA; 98=NA ")
imp<-mice(data = NSDUH_2019[,c("daysalc",  "race_eth", "attempt_suicide", "dep_year2", "male")], m = 5)
## 
##  iter imp variable
##   1   1  daysalc  race_eth  attempt_suicide  dep_year2
##   1   2  daysalc  race_eth  attempt_suicide  dep_year2
##   1   3  daysalc  race_eth  attempt_suicide  dep_year2
##   1   4  daysalc  race_eth  attempt_suicide  dep_year2
##   1   5  daysalc  race_eth  attempt_suicide  dep_year2
##   2   1  daysalc  race_eth  attempt_suicide  dep_year2
##   2   2  daysalc  race_eth  attempt_suicide  dep_year2
##   2   3  daysalc  race_eth  attempt_suicide  dep_year2
##   2   4  daysalc  race_eth  attempt_suicide  dep_year2
##   2   5  daysalc  race_eth  attempt_suicide  dep_year2
##   3   1  daysalc  race_eth  attempt_suicide  dep_year2
##   3   2  daysalc  race_eth  attempt_suicide  dep_year2
##   3   3  daysalc  race_eth  attempt_suicide  dep_year2
##   3   4  daysalc  race_eth  attempt_suicide  dep_year2
##   3   5  daysalc  race_eth  attempt_suicide  dep_year2
##   4   1  daysalc  race_eth  attempt_suicide  dep_year2
##   4   2  daysalc  race_eth  attempt_suicide  dep_year2
##   4   3  daysalc  race_eth  attempt_suicide  dep_year2
##   4   4  daysalc  race_eth  attempt_suicide  dep_year2
##   4   5  daysalc  race_eth  attempt_suicide  dep_year2
##   5   1  daysalc  race_eth  attempt_suicide  dep_year2
##   5   2  daysalc  race_eth  attempt_suicide  dep_year2
##   5   3  daysalc  race_eth  attempt_suicide  dep_year2
##   5   4  daysalc  race_eth  attempt_suicide  dep_year2
##   5   5  daysalc  race_eth  attempt_suicide  dep_year2
print(imp)
## Class: mids
## Number of multiple imputations:  5 
## Imputation methods:
##         daysalc        race_eth attempt_suicide       dep_year2            male 
##           "pmm"       "polyreg"           "pmm"           "pmm"              "" 
## PredictorMatrix:
##                 daysalc race_eth attempt_suicide dep_year2 male
## daysalc               0        1               1         1    1
## race_eth              1        0               1         1    1
## attempt_suicide       1        1               0         1    1
## dep_year2             1        1               1         0    1
## male                  1        1               1         1    0
library(car)
library(stargazer)
## 
## Please cite as:
##  Hlavac, Marek (2018). stargazer: Well-Formatted Regression and Summary Statistics Tables.
##  R package version 5.2.2. https://CRAN.R-project.org/package=stargazer
library(survey)
## Loading required package: grid
## Loading required package: Matrix
## Loading required package: survival
## 
## Attaching package: 'survey'
## The following object is masked from 'package:graphics':
## 
##     dotchart
library(questionr)
library(dplyr)
library(questionr)
library(tableone)
## Warning: package 'tableone' was built under R version 4.0.4
library(boot)
## 
## Attaching package: 'boot'
## The following object is masked from 'package:survival':
## 
##     aml
## The following object is masked from 'package:car':
## 
##     logit
library(table1)
## Warning: package 'table1' was built under R version 4.0.4
## 
## Attaching package: 'table1'
## The following objects are masked from 'package:base':
## 
##     units, units<-
dat.imp<-complete(imp, action = 1)
head(dat.imp, n=10)
##    daysalc  race_eth attempt_suicide dep_year2   male
## 1       20     white               0         1   Male
## 2        5 mult_race               0         0 Female
## 3        1     white               0         0   Male
## 4       15     white               0         0 Female
## 5        6     white               0         0   Male
## 6        5     black               0         1   Male
## 7        6     white               0         0 Female
## 8        3     white               0         0   Male
## 9       10 mult_race               1         0 Female
## 10       1     black               1         0 Female
glm1<-glm(attempt_suicide~factor(race_eth),
          data=dat.imp,
          family = binomial)
summary(glm1)
## 
## Call:
## glm(formula = attempt_suicide ~ factor(race_eth), family = binomial, 
##     data = dat.imp)
## 
## Deviance Residuals: 
##    Min      1Q  Median      3Q     Max  
## -1.051  -0.780  -0.780   1.376   1.781  
## 
## Coefficients:
##                           Estimate Std. Error z value Pr(>|z|)    
## (Intercept)               -1.03398    0.01142 -90.575  < 2e-16 ***
## factor(race_eth)asian      0.72966    0.10859   6.719 1.83e-11 ***
## factor(race_eth)black      0.57885    0.02440  23.725  < 2e-16 ***
## factor(race_eth)hispanic   0.56470    0.04102  13.766  < 2e-16 ***
## factor(race_eth)mult_race -0.32311    0.04404  -7.337 2.19e-13 ***
## factor(race_eth)other      0.30619    0.06988   4.382 1.18e-05 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 67317  on 56135  degrees of freedom
## Residual deviance: 66493  on 56130  degrees of freedom
## AIC: 66505
## 
## Number of Fisher Scoring iterations: 4
glm2<-glm(attempt_suicide~factor(race_eth)+factor(male),
          data=dat.imp,
          family = binomial)
summary(glm2)
## 
## Call:
## glm(formula = attempt_suicide ~ factor(race_eth) + factor(male), 
##     family = binomial, data = dat.imp)
## 
## Deviance Residuals: 
##     Min       1Q   Median       3Q      Max  
## -1.2108  -0.9015  -0.6370   1.2200   1.9818  
## 
## Coefficients:
##                           Estimate Std. Error z value Pr(>|z|)    
## (Intercept)               -0.69046    0.01377 -50.138  < 2e-16 ***
## factor(race_eth)asian      0.76870    0.11072   6.942 3.85e-12 ***
## factor(race_eth)black      0.59087    0.02485  23.782  < 2e-16 ***
## factor(race_eth)hispanic   0.60270    0.04181  14.415  < 2e-16 ***
## factor(race_eth)mult_race -0.32072    0.04460  -7.191 6.41e-13 ***
## factor(race_eth)other      0.32155    0.07109   4.523 6.08e-06 ***
## factor(male)Male          -0.80136    0.01957 -40.938  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 67317  on 56135  degrees of freedom
## Residual deviance: 64754  on 56129  degrees of freedom
## AIC: 64768
## 
## Number of Fisher Scoring iterations: 4
glm3<-glm(attempt_suicide~factor(race_eth)+scale(daysalc)+factor(male),
          data=dat.imp,
          family = binomial)
summary(glm3)
## 
## Call:
## glm(formula = attempt_suicide ~ factor(race_eth) + scale(daysalc) + 
##     factor(male), family = binomial, data = dat.imp)
## 
## Deviance Residuals: 
##     Min       1Q   Median       3Q      Max  
## -1.7928  -0.8303  -0.6310   1.1734   2.1490  
## 
## Coefficients:
##                           Estimate Std. Error z value Pr(>|z|)    
## (Intercept)               -0.64595    0.01405 -45.978  < 2e-16 ***
## factor(race_eth)asian      0.84621    0.11285   7.499 6.45e-14 ***
## factor(race_eth)black      0.63650    0.02536  25.100  < 2e-16 ***
## factor(race_eth)hispanic   0.72192    0.04257  16.960  < 2e-16 ***
## factor(race_eth)mult_race -0.21687    0.04522  -4.796 1.62e-06 ***
## factor(race_eth)other      0.43848    0.07212   6.080 1.20e-09 ***
## scale(daysalc)             0.37578    0.00955  39.350  < 2e-16 ***
## factor(male)Male          -0.99427    0.02067 -48.106  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 67317  on 56135  degrees of freedom
## Residual deviance: 63216  on 56128  degrees of freedom
## AIC: 63232
## 
## Number of Fisher Scoring iterations: 4
glm4<-glm(attempt_suicide~factor(race_eth)+factor(male)+factor(dep_year2),
          data=dat.imp,
          family = binomial)
summary(glm4)
## 
## Call:
## glm(formula = attempt_suicide ~ factor(race_eth) + factor(male) + 
##     factor(dep_year2), family = binomial, data = dat.imp)
## 
## Deviance Residuals: 
##     Min       1Q   Median       3Q      Max  
## -1.2144  -0.9058  -0.6395   1.2163   2.0188  
## 
## Coefficients:
##                           Estimate Std. Error z value Pr(>|z|)    
## (Intercept)               -0.67903    0.01431 -47.441  < 2e-16 ***
## factor(race_eth)asian      0.76551    0.11071   6.915 4.69e-12 ***
## factor(race_eth)black      0.58793    0.02487  23.640  < 2e-16 ***
## factor(race_eth)hispanic   0.60361    0.04182  14.435  < 2e-16 ***
## factor(race_eth)mult_race -0.32450    0.04462  -7.273 3.51e-13 ***
## factor(race_eth)other      0.32354    0.07107   4.552 5.31e-06 ***
## factor(male)Male          -0.80418    0.01960 -41.027  < 2e-16 ***
## factor(dep_year2)1        -0.09044    0.03124  -2.895  0.00379 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 67317  on 56135  degrees of freedom
## Residual deviance: 64745  on 56128  degrees of freedom
## AIC: 64761
## 
## Number of Fisher Scoring iterations: 4
glm5<-glm(attempt_suicide~factor(race_eth)+factor(male)+scale(daysalc)+factor(dep_year2),
          data=dat.imp,
          family = binomial)
summary(glm5)
## 
## Call:
## glm(formula = attempt_suicide ~ factor(race_eth) + factor(male) + 
##     scale(daysalc) + factor(dep_year2), family = binomial, data = dat.imp)
## 
## Deviance Residuals: 
##     Min       1Q   Median       3Q      Max  
## -1.7991  -0.8169  -0.6344   1.1662   2.1966  
## 
## Coefficients:
##                            Estimate Std. Error z value Pr(>|z|)    
## (Intercept)               -0.630758   0.014597 -43.213  < 2e-16 ***
## factor(race_eth)asian      0.841964   0.112850   7.461 8.59e-14 ***
## factor(race_eth)black      0.632484   0.025390  24.910  < 2e-16 ***
## factor(race_eth)hispanic   0.723100   0.042575  16.984  < 2e-16 ***
## factor(race_eth)mult_race -0.221947   0.045240  -4.906 9.29e-07 ***
## factor(race_eth)other      0.441454   0.072105   6.122 9.22e-10 ***
## factor(male)Male          -0.997952   0.020691 -48.231  < 2e-16 ***
## scale(daysalc)             0.376796   0.009557  39.427  < 2e-16 ***
## factor(dep_year2)1        -0.119573   0.031613  -3.782 0.000155 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 67317  on 56135  degrees of freedom
## Residual deviance: 63201  on 56127  degrees of freedom
## AIC: 63219
## 
## Number of Fisher Scoring iterations: 4
library(car)
library(stargazer)
library(survey)
library(questionr)
library(dplyr)
library(tidyverse)
## -- Attaching packages --------------------------------------- tidyverse 1.3.0 --
## v tibble  3.0.3     v purrr   0.3.4
## v tidyr   1.1.1     v stringr 1.4.0
## v readr   1.3.1     v forcats 0.5.0
## -- Conflicts ------------------------------------------ tidyverse_conflicts() --
## x tidyr::expand() masks Matrix::expand()
## x dplyr::filter() masks mice::filter(), stats::filter()
## x dplyr::lag()    masks stats::lag()
## x tidyr::pack()   masks Matrix::pack()
## x dplyr::recode() masks car::recode()
## x purrr::some()   masks car::some()
## x tidyr::unpack() masks Matrix::unpack()
library(emmeans)
rg<-ref_grid(glm3)
marg_logit<-emmeans(object = rg,
              specs = c( "race_eth",  "daysalc"),
             type="response" ,
              data=sub)
knitr::kable(marg_logit,  digits = 4)
race_eth daysalc prob SE df asymp.LCL asymp.UCL
white 7.5837 0.2418 0.0022 Inf 0.2374 0.2461
asian 7.5837 0.4263 0.0274 Inf 0.3736 0.4808
black 7.5837 0.3760 0.0053 Inf 0.3657 0.3864
hispanic 7.5837 0.3962 0.0098 Inf 0.3773 0.4155
mult_race 7.5837 0.2042 0.0071 Inf 0.1907 0.2185
other 7.5837 0.3308 0.0157 Inf 0.3007 0.3623
rg2<-ref_grid(glm2)
marg_logit<-emmeans(object = rg2,
              specs = c( "race_eth",  "male"),
             type="response" ,
              data=sub)
knitr::kable(marg_logit,  digits = 4)
race_eth male prob SE df asymp.LCL asymp.UCL
white Female 0.3339 0.0031 Inf 0.3280 0.3400
asian Female 0.5195 0.0276 Inf 0.4655 0.5732
black Female 0.4751 0.0058 Inf 0.4637 0.4866
hispanic Female 0.4781 0.0103 Inf 0.4580 0.4982
mult_race Female 0.2667 0.0085 Inf 0.2503 0.2838
other Female 0.4088 0.0171 Inf 0.3758 0.4426
white Male 0.1836 0.0025 Inf 0.1788 0.1886
asian Male 0.3267 0.0243 Inf 0.2809 0.3761
black Male 0.2889 0.0051 Inf 0.2790 0.2989
hispanic Male 0.2913 0.0086 Inf 0.2748 0.3084
mult_race Male 0.1403 0.0054 Inf 0.1301 0.1513
other Male 0.2368 0.0128 Inf 0.2126 0.2629
rg3<-ref_grid(glm3)
marg_logit<-emmeans(object = rg3,
              specs = c( "race_eth",  "male", "daysalc"),
             type="response" ,
              data=sub)
knitr::kable(marg_logit,  digits = 4)
race_eth male daysalc prob SE df asymp.LCL asymp.UCL
white Female 7.5837 0.3439 0.0032 Inf 0.3377 0.3501
asian Female 7.5837 0.5499 0.0279 Inf 0.4949 0.6037
black Female 7.5837 0.4976 0.0060 Inf 0.4859 0.5094
hispanic Female 7.5837 0.5190 0.0105 Inf 0.4985 0.5394
mult_race Female 7.5837 0.2968 0.0092 Inf 0.2790 0.3152
other Female 7.5837 0.4483 0.0177 Inf 0.4139 0.4832
white Male 7.5837 0.1624 0.0024 Inf 0.1578 0.1672
asian Male 7.5837 0.3113 0.0242 Inf 0.2660 0.3605
black Male 7.5837 0.2682 0.0050 Inf 0.2585 0.2781
hispanic Male 7.5837 0.2853 0.0086 Inf 0.2687 0.3025
mult_race Male 7.5837 0.1350 0.0053 Inf 0.1250 0.1458
other Male 7.5837 0.2312 0.0128 Inf 0.2070 0.2572