library(car)
library(stargazer)
## 
## Please cite as:
##  Hlavac, Marek (2018). stargazer: Well-Formatted Regression and Summary Statistics Tables.
##  R package version 5.2.1. 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(haven)
dat<-read_xpt("C:/Users/Monica/Documents/CSparks_StatsII/CNTY05.xpt")

a) My outcome variable will be as the self-reporting of number of poor mental health days in the past month. I will recode the variable as a binary variable in the following manner: a value of 88==0 (for NO) and a value of 1 - 30== 1 (for YES)

b) My research question is as follows: Does the threat of violence affect a person’s mental health as much as the act of violence?

c) My two predictor variables are 1)the threat of violence and 2)the act of violence

Rename variables

nams<-names(dat)
head(nams, n=10)
##  [1] "_CNTYNAM" "_STATE"   "FMONTH"   "IDATE"    "IMONTH"   "IDAY"    
##  [7] "IYEAR"    "INTVID"   "DISPCODE" "SEQNO"
myvariables<-tolower(gsub(pattern = "_",replacement =  "",x =  nams))
names(dat)<-myvariables

Recode for poor mental health (outcome variable)

dat$poormenhlth<-recode (dat$menthlth,recodes= "1:30=1; 88=0; else=NA")

Recode for intimate partner threat (predictor variable)

dat$threat<-recode (dat$ipvthrat, recodes = "7:9=NA; 1='1YESthreat'; 2='0NOthreat'")

Recode for intimate partner hurt (predictor variable)

dat$hurt<-recode(dat$ipvphhrt, recodes = "7:9=NA; 1='1YEShurt'; 2='0NOhurt'")

Perform descriptive analysis of the outcome variable for intimate partner threat.

table(dat$poormenhlth, dat$threat)
##    
##     0NOthreat 1YESthreat
##   0     17346       2187
##   1      7784       2559

Column percentages of poor mental health and threat

prop.table(table(dat$poormenhlth, dat$threat), margin=2)
##    
##     0NOthreat 1YESthreat
##   0 0.6902507  0.4608091
##   1 0.3097493  0.5391909

Perform descriptive analysis of outcome variable for intimate partner hurt.

table(dat$poormenhlth, dat$hurt)
##    
##     0NOhurt 1YEShurt
##   0   17148     2354
##   1    7622     2714

Column percentages of poor mental health and hurt

prop.table(table(dat$poormenhlth, dat$hurt), margin=2)
##    
##       0NOhurt  1YEShurt
##   0 0.6922891 0.4644830
##   1 0.3077109 0.5355170

Create survey design

options(survey.lonely.psu = "adjust")
mydesign<-svydesign(ids = ~1, strata=~ststr, weights=~cntywt, data = dat[is.na(dat$cntywt)==F,])

Weighted analysis for poor mental health and threat

newwts<-wtd.table(dat$poormenhlth, dat$threat, weights= dat$cntywt)
prop.table(wtd.table(dat$poormenhlth, dat$threat, weights= dat$cntywt), margin = 2)
##   0NOthreat 1YESthreat
## 0 0.6688194  0.4560180
## 1 0.3311806  0.5439820

Non-weighted analysis for poor mental health and threat

prop.table(table(dat$poormenhlth, dat$threat), margin=2)
##    
##     0NOthreat 1YESthreat
##   0 0.6902507  0.4608091
##   1 0.3097493  0.5391909

THe results of the weighted analysis for poor mental health and threat is slightly higher for respondents who reported no poor mental health and no threats than the non-weighted analysis.

Weighted analysis for poor mental health and hurt

newwts<-wtd.table(dat$poormenhlth, dat$hurt, weights= dat$cntywt)
prop.table(wtd.table(dat$poormenhlth, dat$hurt, weights= dat$cntywt), margin = 2)
##     0NOhurt  1YEShurt
## 0 0.6695058 0.4680454
## 1 0.3304942 0.5319546

Non-weighted analyis for poor mental health and hurt

prop.table(table(dat$poormenhlth, dat$hurt), margin=2)
##    
##       0NOhurt  1YEShurt
##   0 0.6922891 0.4644830
##   1 0.3077109 0.5355170

THe results of the weighted analysis for poor mental health and hurt is slightly lower for respondents who reported no poor mental health and no hurt than the non-weighted analysis but slightly higher for respondents who reported poor mental health and no hurt.

Repeating the same process as above, we find that the outcome for weighted proportions for the individuals who were hurt compared to the unweighted proportions are essentially the same as those who were only threatened.

Determining the standard error of the percentages for poor mental health and threat

n<-table(is.na(dat$poormenhlth)==F)
n
## 
##  FALSE   TRUE 
##   2705 168978
p<-prop.table(wtd.table(dat$poormenhlth, dat$threat, weights = dat$cntywt), margin=2)
se<-sqrt((p*(1-p))/n[2])
data.frame(proportion=p, se=se)
##   proportion.Var1 proportion.Var2 proportion.Freq se.Var1    se.Var2
## 1               0       0NOthreat       0.6688194       0  0NOthreat
## 2               1       0NOthreat       0.3311806       1  0NOthreat
## 3               0      1YESthreat       0.4560180       0 1YESthreat
## 4               1      1YESthreat       0.5439820       1 1YESthreat
##       se.Freq
## 1 0.001144911
## 2 0.001144911
## 3 0.001211625
## 4 0.001211625

Determine the standard error of the percentages for poor mental health and hurt

n<-table(is.na(dat$poormenhlth)==F)
n
## 
##  FALSE   TRUE 
##   2705 168978
p<-prop.table(wtd.table(dat$poormenhlth, dat$hurt, weights = dat$cntywt), margin=2)
data.frame(proportion=p, se=se)
##   proportion.Var1 proportion.Var2 proportion.Freq se.Var1    se.Var2
## 1               0         0NOhurt       0.6695058       0  0NOthreat
## 2               1         0NOhurt       0.3304942       1  0NOthreat
## 3               0        1YEShurt       0.4680454       0 1YESthreat
## 4               1        1YEShurt       0.5319546       1 1YESthreat
##       se.Freq
## 1 0.001144911
## 2 0.001144911
## 3 0.001211625
## 4 0.001211625

Creating proper survey design analysis

newwts<-svytable(~poormenhlth+threat, design= mydesign)
prop.table(svytable(~poormenhlth+threat, design=mydesign), margin = 2)
##            threat
## poormenhlth 0NOthreat 1YESthreat
##           0 0.6688194  0.4560180
##           1 0.3311806  0.5439820

Calculate statistics

mysvytable<-svyby(formula=~poormenhlth, by=~threat, design= mydesign, FUN = svymean, na.rm = T)
mysvytable
##                threat poormenhlth          se
## 0NOthreat   0NOthreat   0.3311806 0.006387965
## 1YESthreat 1YESthreat   0.5439820 0.013579769
mysvytable<-svyby(formula=~poormenhlth, by=~hurt, design= mydesign, FUN = svymean, na.rm = T)
mysvytable
##              hurt poormenhlth         se
## 0NOhurt   0NOhurt   0.3304942 0.00644844
## 1YEShurt 1YEShurt   0.5319546 0.01311238