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=0;2: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) -0.14876    0.01506  -9.877   <2e-16 ***
## marital1    -0.75125    0.01809 -41.528   <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) -0.65585    0.01057 -62.020   <2e-16 ***
## sxorient1    0.70461    0.07080   9.952   <2e-16 ***
## ---
## 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) -0.18154    0.02400  -7.565 3.89e-14 ***
## int1.0      -0.75810    0.02828 -26.808  < 2e-16 ***
## int0.1       0.48298    0.10399   4.644 3.41e-06 ***
## int1.1      -0.39163    0.13634  -2.872  0.00407 ** 
## ---
## 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: 4
#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)      0.54830    1.25089   0.438 0.661181    
## int1.0          -0.16169    0.24001  -0.674 0.500556    
## int0.1           0.37280    0.58463   0.638 0.523737    
## int1.1           1.13669    0.98863   1.150 0.250345    
## sex1            -0.45078    0.15885  -2.838 0.004576 ** 
## age25-34        -0.57685    0.94754  -0.609 0.542716    
## age35-44        -1.12223    0.96125  -1.167 0.243122    
## age45-54        -0.79643    0.94962  -0.839 0.401719    
## age55-64        -1.11845    0.95800  -1.167 0.243121    
## age65+          -1.74400    0.95456  -1.827 0.067807 .  
## raceblack       -0.11213    0.31341  -0.358 0.720531    
## racehispanic     0.14686    0.47654   0.308 0.757972    
## raceother       -0.06902    0.42516  -0.162 0.871045    
## income>75k      -0.72867    0.35497  -2.053 0.040190 *  
## income$20k<$35k -0.62400    0.35315  -1.767 0.077347 .  
## income35k<50k   -0.50571    0.36094  -1.401 0.161301    
## income50k<75    -0.79036    0.36071  -2.191 0.028528 *  
## educ1somehs     -0.04017    0.79056  -0.051 0.959474    
## educ2hsgrad      0.19865    0.74582   0.266 0.789989    
## educ3somecol     0.23813    0.76295   0.312 0.754971    
## educ4colgrad     0.41696    0.76405   0.546 0.585305    
## smokeY           0.10633    0.18426   0.577 0.563952    
## drinks1         -0.96876    0.54821  -1.767 0.077322 .  
## cholesY          0.60338    0.15683   3.847 0.000122 ***
## blpressureY     -0.15075    0.22526  -0.669 0.503415    
## bmi              0.04895    0.10158   0.482 0.629935    
## medcostY         0.59763    0.24935   2.397 0.016608 *  
## hlthinsurY       0.62896    0.36508   1.723 0.085036 .  
## checkup2 years  -0.10696    0.27540  -0.388 0.697756    
## checkup5 years   0.32157    0.37913   0.848 0.396415    
## checkup5+ years  0.81412    0.45872   1.775 0.076050 .  
## ownhomeY        -0.49428    0.24551  -2.013 0.044184 *  
## safeneigh1       0.93952    0.42573   2.207 0.027409 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1.187565)
## 
## Number of Fisher Scoring iterations: 4