library(car)
## Loading required package: carData
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(foreign)
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
brfss2017=read.xport("~/Desktop/DataFolder/brfss.XPT ")
#Marital Status; Are you? (marital status): 1 Married, 2 Divorced, 3 Widowed, 4 Seprarated, 5 Never married, 6 Partner, 9 Refuse
brfss2017$marital<-Recode(brfss2017$MARITAL, recodes="1=1;5=0;else=NA",as.factor=T)
#Sexual Orientation; Do you consider yourself to be? (sexual orientation): 1 Straight, 2 Gay, 3 Bisexual, 4 Other, 7 Don't know, 9 Refuse
brfss2017$sxorient<-Recode(brfss2017$SXORIENT,recodes="1=0;2=1;else=NA",as.factor=T)
#Income 
brfss2017$income<-Recode(brfss2017$INCOME2,recodes="1:3='<$20k';4:5='$20k<$35k';6='35k<50k';7='50k<75';8='>75k';else=NA",as.factor=T)
#Age
brfss2017$age<-Recode(brfss2017$X_AGE_G, recodes="1='18-24';2='25-34';3='35-44';4='45-54';5='55-64';6='65+';else=NA",as.factor=T)
#Education
brfss2017$educ<-Recode(brfss2017$EDUCA, recodes="1:2='0Prim'; 3='1somehs'; 4='2hsgrad'; 5='3somecol'; 6='4colgrad';9=NA", as.factor=T)
brfss2017$educ<-relevel(brfss2017$educ, ref='0Prim')
#Cholesterol; Have you EVER been told by a doctor that your cholesterol is high? 1 Yes, 2 No, 7 Dont know, 9 Refuse
brfss2017$choles<-Recode(brfss2017$TOLDHI2,recodes="1='Y';2='N';else=NA",as.factor=T)
#Blood Pressure; Are you currently taking medicine for blood Pressure? 1 Yes, 2 No, 7 Dont know, 9 Refuse
brfss2017$blpressure<-Recode(brfss2017$BPMEDS,recodes="1='Y';2='N';else=NA",as.factor=T)
#Health Insurance; Do you have health care coverage? 1 Yes, 2 No, 7 Dont know, 9 Refuse
brfss2017$hlthinsur<-Recode(brfss2017$HLTHPLN1,recodes="1='Y';2='N';else=NA",as.factor=T)
#DR Checkup; Had a doctor checkup with the last...?  1 year, 2 years, 5 years, 5+ years
brfss2017$checkup<-Recode(brfss2017$CHECKUP1,recodes="1='1 year';2='2 years';3='5 years';4='5+ years';else=NA",as.factor=T)
#Live in a Safe Neighborhood; How safe from crime do you consider your neighborhood to be? 1 Extremely Safe, 2 Safe, 3 Unsafe, 4 Extremely Unsafe, 7 Don't know, 9 Refuse
brfss2017$safeneigh<-Recode(brfss2017$HOWSAFE1,recodes="1:2=0;3:4=1;else=NA",as.factor=T)
#Mental Health; How many days during the past 30 days was your mental health not good? 1-30, 88 None, 77 Don't Know, 99 Refuse
brfss2017$menth<-Recode(brfss2017$MENTHLTH,recodes="88=0;77=NA;99=NA")
#Interaction between marital status and sexual orientation
brfss2017$int<-interaction(brfss2017$marital,brfss2017$sxorient)
brfss2017$int<-as.factor(brfss2017$int)
#Interaction between marital status and sexual orientation recoded to gay people only (gay single and gay married)
brfss2017$int1<-Recode(brfss2017$int,recodes="'M.G'='M.G';'S.G'='S.G';else=NA",as.factor=T)
#Interaction between marital status and sexual orientation recoded to married people only (married straight and married gay)
brfss2017$int2<-Recode(brfss2017$int,recodes="'M.G'='M.G';'M.S'='M.S';else=NA",as.factor=T)
#Homeowner; Do you own or rent your home? 1 Own, 2 Rent, 3 Other, 7 Don't know, 9 Refuse
brfss2017$ownhome<-Recode(brfss2017$RENTHOM1,recodes="1='Y';2:3='N';else=NA",as.factor=T)
#Smoking; Do you now smoke cigarettes every day, some days, or not at all? 1 Everyday, 2 Some days, 3 Not at all, 7 Don't know, 9 Refuse
brfss2017$smoke<-Recode(brfss2017$SMOKDAY2,recodes="1:2='Y';3='N';else=NA",as.factor=T)
#Was there a time in the past 12 months when you needed to see a doctor but could not because of cost? 1 Yes, 2 No, 7 Don't know, 9 Refuse
brfss2017$medcost<-Recode(brfss2017$MEDCOST,recodes="1='Y';2='N';else=NA",as.factor=T)
#During the past 30 days, on the days when you drank, about how many drinks did you drink on the average?
brfss2017$drinks<-Recode(brfss2017$AVEDRNK2,recodes="1:10=0;10:30=1;77=NA;99=NA;else=NA",as.factor=T)
#BMI; Four-categories of Body Mass Index (BMI): 1 Underweight, 2 Normal Weight, 3 Overweight, 4 Obese
brfss2017$bmi<-brfss2017$X_BMI5CAT
#Not good mental health days: 1 (0 days), 2 (1-13 days), 3 (14-30 days), 9 Don't know
brfss2017$badmental<-Recode(brfss2017$X_MENT14D,recodes="1:2=0;3=1;9=NA;else=NA",as.factor=T)
#Race/Ethnicity
brfss2017$black<-Recode(brfss2017$racegr3, recodes="2=1; 9=NA; else=0")
brfss2017$white<-Recode(brfss2017$racegr3, recodes="1=1; 9=NA; else=0")
brfss2017$other<-Recode(brfss2017$racegr3, recodes="3:4=1; 9=NA; else=0")
brfss2017$hispanic<-Recode(brfss2017$racegr3, recodes="5=1; 9=NA; else=0")
#Race/Ethnicity; 1 White, 2 Black, 3 Other, 4 Multiracial, 5 Hispanic, 9 Don't know
brfss2017$race<-Recode(brfss2017$X_RACEGR3, recodes="2='black'; 1='awhite'; 3:4='other';5='hispanic'; else=NA",as.factor=T)
brfss2017$race<-relevel(brfss2017$race, ref='awhite')
#Sex; 1 Male, 2 Female
brfss2017$sex<-Recode(brfss2017$SEX,recodes="1=1;2=0;9=NA;else=NA",as.factor=T)
#Create a survey design
brfss2017%>%
  filter(complete.cases(.))
##   [1] X_STATE    FMONTH     IDATE      IMONTH     IDAY       IYEAR     
##   [7] DISPCODE   SEQNO      X_PSU      CTELENM1   PVTRESD1   COLGHOUS  
##  [13] STATERE1   CELLFON4   LADULT     NUMADULT   NUMMEN     NUMWOMEN  
##  [19] SAFETIME   CTELNUM1   CELLFON5   CADULT     PVTRESD3   CCLGHOUS  
##  [25] CSTATE1    LANDLINE   HHADULT    GENHLTH    PHYSHLTH   MENTHLTH  
##  [31] POORHLTH   HLTHPLN1   PERSDOC2   MEDCOST    CHECKUP1   BPHIGH4   
##  [37] BPMEDS     CHOLCHK1   TOLDHI2    CHOLMED1   CVDINFR4   CVDCRHD4  
##  [43] CVDSTRK3   ASTHMA3    ASTHNOW    CHCSCNCR   CHCOCNCR   CHCCOPD1  
##  [49] HAVARTH3   ADDEPEV2   CHCKIDNY   DIABETE3   DIABAGE2   LMTJOIN3  
##  [55] ARTHDIS2   ARTHSOCL   JOINPAI1   SEX        MARITAL    EDUCA     
##  [61] RENTHOM1   NUMHHOL2   NUMPHON2   CPDEMO1A   VETERAN3   EMPLOY1   
##  [67] CHILDREN   INCOME2    INTERNET   WEIGHT2    HEIGHT3    PREGNANT  
##  [73] DEAF       BLIND      DECIDE     DIFFWALK   DIFFDRES   DIFFALON  
##  [79] SMOKE100   SMOKDAY2   STOPSMK2   LASTSMK2   USENOW3    ECIGARET  
##  [85] ECIGNOW    ALCDAY5    AVEDRNK2   DRNK3GE5   MAXDRNKS   FRUIT2    
##  [91] FRUITJU2   FVGREEN1   FRENCHF1   POTATOE1   VEGETAB2   EXERANY2  
##  [97] EXRACT11   EXEROFT1   EXERHMM1   EXRACT21   EXEROFT2   EXERHMM2  
## [103] STRENGTH   SEATBELT   FLUSHOT6   FLSHTMY2   PNEUVAC3   SHINGLE2  
## [109] HIVTST6    HIVTSTD3   HIVRISK5   PDIABTST   PREDIAB1   INSULIN   
## [115] BLDSUGAR   FEETCHK2   DOCTDIAB   CHKHEMO3   FEETCHK    EYEEXAM   
## [121] DIABEYE    DIABEDU    COPDCOGH   COPDFLEM   COPDBRTH   COPDBTST  
## [127] COPDSMOK   HAREHAB1   STREHAB1   CVDASPRN   ASPUNSAF   RLIVPAIN  
## [133] RDUCHART   RDUCSTRK   BPEATHBT   BPSALT     BPALCHOL   BPEXER    
## [139] BPEATADV   BPSLTADV   BPALCADV   BPEXRADV   BPMEDADV   BPHI2MR   
## [145] ARTTODAY   ARTHWGT    ARTHEXER   ARTHEDU    ASTHMAGE   ASATTACK  
## [151] ASERVIST   ASDRVIST   ASRCHKUP   ASACTLIM   ASYMPTOM   ASNOSLEP  
## [157] ASTHMED3   ASINHALR   PAINACT2   QLMENTL2   QLSTRES2   QLHLTH2   
## [163] SLEPTIM1   ADSLEEP    SLEPDAY1   SLEPSNO2   SLEPBRTH   MEDICARE  
## [169] HLTHCVR1   DELAYMED   DLYOTHER   NOCOV121   LSTCOVRG   DRVISITS  
## [175] MEDSCOS1   CARERCVD   MEDBILL1   ASBIALCH   ASBIDRNK   ASBIBING  
## [181] ASBIADVC   ASBIRDUC   CNCRDIFF   CNCRAGE    CNCRTYP1   CSRVTRT2  
## [187] CSRVDOC1   CSRVSUM    CSRVRTRN   CSRVINST   CSRVINSR   CSRVDEIN  
## [193] CSRVCLIN   CSRVPAIN   CSRVCTL1   SSBSUGR2   SSBFRUT3   WTCHSALT  
## [199] DRADVISE   MARIJANA   USEMRJN1   RSNMRJNA   PFPPRVN2   TYPCNTR7  
## [205] NOBCUSE6   IMFVPLAC   HPVADVC2   HPVADSHT   TETANUS    LCSFIRST  
## [211] LCSLAST    LCSNUMCG   LCSCTSCN   CAREGIV1   CRGVREL2   CRGVLNG1  
## [217] CRGVHRS1   CRGVPRB2   CRGVPERS   CRGVHOUS   CRGVMST2   CRGVEXPT  
## [223] CIMEMLOS   CDHOUSE    CDASSIST   CDHELP     CDSOCIAL   CDDISCUS  
## [229] EMTSUPRT   LSATISFY   SDHBILLS   SDHMOVE    HOWSAFE1   SDHFOOD   
## [235] SDHMEALS   SDHMONEY   SDHSTRES   SXORIENT   TRNSGNDR   FIREARM4  
## [241] GUNLOAD    LOADULK2   RCSGENDR   RCSRLTN2   CASTHDX2   CASTHNO2  
## [247] QSTVER     QSTLANG    MSCODE     X_STSTR    X_STRWT    X_RAWRAKE 
## [253] X_WT2RAKE  X_IMPRACE  X_CHISPNC  X_CRACE1   X_CPRACE   X_CLLCPWT 
## [259] X_DUALUSE  X_DUALCOR  X_LLCPWT2  X_LLCPWT   X_RFHLTH   X_PHYS14D 
## [265] X_MENT14D  X_HCVU651  X_RFHYPE5  X_CHOLCH1  X_RFCHOL1  X_MICHD   
## [271] X_LTASTH1  X_CASTHM1  X_ASTHMS1  X_DRDXAR1  X_LMTACT1  X_LMTWRK1 
## [277] X_LMTSCL1  X_PRACE1   X_MRACE1   X_HISPANC  X_RACE     X_RACEG21 
## [283] X_RACEGR3  X_RACE_G1  X_AGEG5YR  X_AGE65YR  X_AGE80    X_AGE_G   
## [289] HTIN4      HTM4       WTKG3      X_BMI5     X_BMI5CAT  X_RFBMI5  
## [295] X_CHLDCNT  X_EDUCAG   X_INCOMG   X_SMOKER3  X_RFSMOK3  X_ECIGSTS 
## [301] X_CURECIG  DRNKANY5   DROCDY3_   X_RFBING5  X_DRNKWEK  X_RFDRHV5 
## [307] FTJUDA2_   FRUTDA2_   GRENDA1_   FRNCHDA_   POTADA1_   VEGEDA2_  
## [313] X_MISFRT1  X_MISVEG1  X_FRTRES1  X_VEGRES1  X_FRUTSU1  X_VEGESU1 
## [319] X_FRTLT1A  X_VEGLT1A  X_FRT16A   X_VEG23A   X_FRUITE1  X_VEGETE1 
## [325] X_TOTINDA  METVL11_   METVL21_   MAXVO2_    FC60_      ACTIN11_  
## [331] ACTIN21_   PADUR1_    PADUR2_    PAFREQ1_   PAFREQ2_   X_MINAC11 
## [337] X_MINAC21  STRFREQ_   PAMISS1_   PAMIN11_   PAMIN21_   PA1MIN_   
## [343] PAVIG11_   PAVIG21_   PA1VIGM_   X_PACAT1   X_PAINDX1  X_PA150R2 
## [349] X_PA300R2  X_PA30021  X_PASTRNG  X_PAREC1   X_PASTAE1  X_RFSEAT2 
## [355] X_RFSEAT3  X_FLSHOT6  X_PNEUMO2  X_AIDTST3  marital    sxorient  
## [361] income     age        educ       choles     blpressure hlthinsur 
## [367] checkup    safeneigh  menth      int        int1       int2      
## [373] ownhome    smoke      medcost    drinks     bmi        badmental 
## [379] black      white      other      hispanic   race       sex       
## <0 rows> (or 0-length row.names)
options(survey.lonely.psu = "adjust")
des<-svydesign(ids = ~1,strata=~X_STRWT,weights = ~X_LLCPWT,data = brfss2017)
#Mental Health by Marital Status
fit.logit1<-svyglm(badmental~marital,design=des,family=binomial)
## Warning in eval(family$initialize): non-integer #successes in a binomial
## glm!
summary(fit.logit1)
## 
## Call:
## svyglm(formula = badmental ~ marital, design = des, family = binomial)
## 
## Survey design:
## svydesign(ids = ~1, strata = ~X_STRWT, weights = ~X_LLCPWT, data = brfss2017)
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -1.66807    0.02009  -83.02   <2e-16 ***
## marital1    -0.66855    0.02577  -25.94   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1.088388)
## 
## Number of Fisher Scoring iterations: 4
#Mental Health by Sexual Orientation
fit.logit2<-svyglm(badmental~sxorient,design=des,family=binomial)
## Warning in eval(family$initialize): non-integer #successes in a binomial
## glm!
summary(fit.logit2)
## 
## Call:
## svyglm(formula = badmental ~ sxorient, design = des, family = binomial)
## 
## Survey design:
## svydesign(ids = ~1, strata = ~X_STRWT, weights = ~X_LLCPWT, data = brfss2017)
## 
## Coefficients:
##             Estimate Std. Error  t value Pr(>|t|)    
## (Intercept) -2.04262    0.01562 -130.809   <2e-16 ***
## sxorient1    0.72776    0.09511    7.651    2e-14 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1.295568)
## 
## Number of Fisher Scoring iterations: 4
table(brfss2017$sxorient)
## 
##      0      1 
## 190310   3167
#Mental Health by Interaction between Mental Health and Sexual Orientation
fit.logit3<-svyglm(badmental~int,design=des,family=binomial)
## Warning in eval(family$initialize): non-integer #successes in a binomial
## glm!
summary(fit.logit3)
## 
## Call:
## svyglm(formula = badmental ~ int, design = des, family = binomial)
## 
## Survey design:
## svydesign(ids = ~1, strata = ~X_STRWT, weights = ~X_LLCPWT, data = brfss2017)
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -1.75010    0.03266 -53.592  < 2e-16 ***
## int1.0      -0.64897    0.04087 -15.879  < 2e-16 ***
## int0.1       0.64021    0.13675   4.682 2.85e-06 ***
## int1.1      -0.32492    0.17163  -1.893   0.0583 .  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1.425193)
## 
## Number of Fisher Scoring iterations: 5
#Mental Health by a slew of variables
fit.logit4<-svyglm(badmental~int+sex+age+race+income+educ+smoke+drinks+choles+blpressure+bmi+medcost+hlthinsur+checkup+ownhome+safeneigh,design=des,family=binomial)
## Warning in eval(family$initialize): non-integer #successes in a binomial
## glm!
summary(fit.logit4)
## 
## Call:
## svyglm(formula = badmental ~ int + sex + age + race + income + 
##     educ + smoke + drinks + choles + blpressure + bmi + medcost + 
##     hlthinsur + checkup + ownhome + safeneigh, design = des, 
##     family = binomial)
## 
## Survey design:
## svydesign(ids = ~1, strata = ~X_STRWT, weights = ~X_LLCPWT, data = brfss2017)
## 
## Coefficients:
##                 Estimate Std. Error t value Pr(>|t|)    
## (Intercept)     -2.99171    1.26825  -2.359 0.018398 *  
## int1.0          -0.48930    0.31087  -1.574 0.115606    
## int0.1          -1.24947    0.89501  -1.396 0.162817    
## int1.1           1.11477    0.67300   1.656 0.097754 .  
## sex1            -0.08402    0.22308  -0.377 0.706472    
## age25-34         1.81098    0.82636   2.192 0.028499 *  
## age35-44         2.14937    0.82924   2.592 0.009594 ** 
## age45-54         1.43176    0.80551   1.777 0.075604 .  
## age55-64         1.53957    0.82029   1.877 0.060645 .  
## age65+           1.52439    0.85831   1.776 0.075840 .  
## raceblack       -0.25435    0.46092  -0.552 0.581115    
## racehispanic    -0.36336    0.58387  -0.622 0.533779    
## raceother        0.34168    0.45264   0.755 0.450394    
## income>75k      -0.33530    0.38270  -0.876 0.381026    
## income$20k<$35k -0.52783    0.40804  -1.294 0.195921    
## income35k<50k   -0.29414    0.41532  -0.708 0.478873    
## income50k<75    -0.77445    0.40058  -1.933 0.053299 .  
## educ1somehs     -0.61253    0.77422  -0.791 0.428920    
## educ2hsgrad     -0.94169    0.71563  -1.316 0.188325    
## educ3somecol    -1.45708    0.73365  -1.986 0.047126 *  
## educ4colgrad    -1.16338    0.74718  -1.557 0.119579    
## smokeY           0.89896    0.25989   3.459 0.000551 ***
## drinks1         -0.17957    0.80586  -0.223 0.823686    
## cholesY          0.52870    0.23297   2.269 0.023322 *  
## blpressureY     -0.08530    0.32237  -0.265 0.791340    
## bmi              0.06629    0.14762   0.449 0.653457    
## medcostY         1.14382    0.28042   4.079 4.65e-05 ***
## hlthinsurY       0.82625    0.43084   1.918 0.055242 .  
## checkup2 years   0.17113    0.40743   0.420 0.674504    
## checkup5 years   0.19548    0.61833   0.316 0.751921    
## checkup5+ years  0.16349    0.44959   0.364 0.716161    
## ownhomeY        -0.56904    0.29674  -1.918 0.055258 .  
## safeneigh1       0.45019    0.49237   0.914 0.360627    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1.250957)
## 
## Number of Fisher Scoring iterations: 5