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(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(forcats)

brfss20<- readRDS(url("https://github.com/coreysparks/DEM7283/blob/master/data/brfss20sm.rds?raw=true"))
names(brfss20) <- tolower(gsub(pattern = "_",replacement =  "",x =  names(brfss20)))

Sample

Data for this analysis come from the Behavioral Risk Factor Surveillance System (BRFSS) survey, published in 2021.

Dependent Variable

The outcome of interest is whether respondents avoided going to the doctor when they needed to in the past twelve month due to cost. This outcome was measured as a binary variable, yes people avoided going to the doctor due to cost or no they did not avoid going to the doctor due to cost. I coded the values as 0 for no and 1 for yes. I removed all “unknown” and blank responses from the analysis.

Research Question

Does having health insurance or a primary care physician influence whether a person had to skip going to the doctor due to cost?

Independent Variables

I included a variable measuring whether respondents had health insurance, defined as “any kind of health coverage, including health insurance, prepaid plans such as HMOs, or government plans such as Medicare, or Indian Health Service.” I coded the values as 0 for no and 1 for yes. I also included a variable measuring whether the respondent has a primary care physician, defined as one or more people respondents “think of as your personal doctor or health care provider.” I coded responses of no as 0 and responses of “yes, only one” and “more than one” as 1. For both variables, I removed all “unknown” and blank responses from the analysis.

#recode
brfss20$ins<-Recode(brfss20$hlthpln1,
                    recodes ="7:9=NA; 1=1;2=0")

brfss20$doctor<-Recode(brfss20$persdoc2,
                    recodes ="7:9=NA; 1:2=1;3=0")

brfss20$skip<-Recode(brfss20$medcost,
                    recodes ="7:9=NA; 1=1;2=0")

#filter out NAs
brfss20<-brfss20%>%
  filter(is.na(ins)==F, is.na(doctor)==F, is.na(skip)==F)

Analysis

I examined whether having insurance and/or having a primary care physician influenced whether a person skipped going to the doctor when they needed to because of cost. I used R (version 4.1.2) to conduct the analysis.

I conducted chi-square tests of significance between the two independent variables and each independent variable and the dependent variable. All chi-square tests were significant (P < 0.001). Next I examined the proportion of people who reported skipping going to the doctor due to cost who have some form of health insurance and who have a primary care doctor. I conducted chi-square tests and found proportions first for the unweighted data then with the weighted data.

Unweighted Survey Design and Results

#skipped dr bc of cost x have primary care doctor
100*prop.table(table(brfss20$doctor, brfss20$skip), margin=2)
##    
##            0        1
##   0 15.90592 35.90432
##   1 84.09408 64.09568
#skipped dr bc of cost x insurance
100*prop.table(table(brfss20$ins, brfss20$skip), margin=2)
##    
##             0         1
##   0  6.115479 31.167979
##   1 93.884521 68.832021
#install.packages("tableone")
library(tableone)

pstable<-CreateTableOne(vars=c("doctor"),strata="skip",data=brfss20)
pstable
##                     Stratified by skip
##                      0             1            p      test
##   n                  176274        16764                   
##   doctor (mean (SD))   0.84 (0.37)  0.64 (0.48) <0.001

84% of people who reported not skipping going to the doctor when they needed to in the last twelve months due to the cost have a primary care doctor while 64% of people who reported skipping doctor when they needed to in the last twelve months due to the cost have a primary care doctor (P < 0.001).

istable<-CreateTableOne(vars=c("ins"),strata="skip",data=brfss20)
istable
##                  Stratified by skip
##                   0             1            p      test
##   n               176274        16764                   
##   ins (mean (SD))   0.94 (0.24)  0.69 (0.46) <0.001

94% of people who reported not skipping going to the doctor when they needed to in the last twelve months due to the cost reported having insurance while 64% of people who reported skipping doctor when they needed to in the last twelve months due to the cost reported having insurance (P < 0.001).

The percentage of people who reported skipping going to the doctor due to cost is higher if they do not have a primary care physician or if they do not have insurance. However, the percentage of people who have insurance and still reported skipping going to the doctor due to cost is high.

#basic chi square test of independence (signif = variables are NOT independent)
#skipped dr bc of cost and primary care dr; x2=4170.4, p<2.2e-16
chisq.test(table(brfss20$doctor, brfss20$skip))
## 
##  Pearson's Chi-squared test with Yates' continuity correction
## 
## data:  table(brfss20$doctor, brfss20$skip)
## X-squared = 4212.2, df = 1, p-value < 2.2e-16
#skipped dr bc of cost and insurance, x2=12352, p<2.2e-16
chisq.test(table(brfss20$ins, brfss20$skip))
## 
##  Pearson's Chi-squared test with Yates' continuity correction
## 
## data:  table(brfss20$ins, brfss20$skip)
## X-squared = 12632, df = 1, p-value < 2.2e-16

Weighted Survey Design and Results

library(srvyr)
## 
## Attaching package: 'srvyr'
## The following object is masked from 'package:stats':
## 
##     filter
options(survey.lonely.psu = "adjust") #good idea to run whether you have lonely psu or not
surveydesign<-svydesign(ids=~1, strata=~ststr, weights=~mmsawt, data = brfss20 )
surveydesign
## Stratified Independent Sampling design (with replacement)
## svydesign(ids = ~1, strata = ~ststr, weights = ~mmsawt, data = brfss20)
#have primary care dr and skipping dr
doctorskip<-wtd.table(brfss20$doctor,
               brfss20$skip,
               weights = brfss20$mmsawt)
doctorskip
##           0         1
## 0  27451107   7050394
## 1 104852172   9034661
100*prop.table(
  doctorskip,
  margin=2)
##          0        1
## 0 20.74862 43.83195
## 1 79.25138 56.16805
#insurance and skipping dr
insuranceskip<-wtd.table(brfss20$ins,
               brfss20$skip,
               weights = brfss20$mmsawt)

100*prop.table(
  insuranceskip,
  margin=2)
##           0         1
## 0  9.176328 38.116124
## 1 90.823672 61.883876

Weighted Proportions with Test Statistics

#primary doctor and skipped dr
primaryskiptable<-svyby(formula = ~doctor,
           by=~skip,
           design = surveydesign,
           FUN=svymean)
primaryskiptable
##   skip    doctor          se
## 0    0 0.7925138 0.002349790
## 1    1 0.5616805 0.008403185
svychisq(~doctor+skip, design = surveydesign)
## 
##  Pearson's X^2: Rao & Scott adjustment
## 
## data:  svychisq(~doctor + skip, design = surveydesign)
## F = 928.57, ndf = 1, ddf = 191856, p-value < 2.2e-16
#insurance and skipped dr
insuranceskiptable<-svyby(formula = ~ins,
           by=~skip,
           design = surveydesign,
           FUN=svymean)
insuranceskiptable
##   skip       ins          se
## 0    0 0.9082367 0.001852145
## 1    1 0.6188388 0.008362045
svychisq(~skip+ins, design = surveydesign)
## 
##  Pearson's X^2: Rao & Scott adjustment
## 
## data:  svychisq(~skip + ins, design = surveydesign)
## F = 2235.5, ndf = 1, ddf = 191856, p-value < 2.2e-16
library(gtsummary)

brfss20%>%
  as_survey_design( strata =ststr,
                    weights = mmsawt)%>%
  select(doctor, skip)%>%
  tbl_svysummary(by = skip, 
              label = list(doctor = "Have a Primary Care Doctor"))%>%
  add_p()%>%
  add_n()
Characteristic N 0, N = 132,303,2801 1, N = 16,085,0551 p-value2
Have a Primary Care Doctor 148,388,335 104,852,172 (79%) 9,034,661 (56%) <0.001

1 n (%)

2 chi-squared test with Rao & Scott's second-order correction 

brfss20%>%
  as_survey_design( strata =ststr,
                    weights = mmsawt)%>%
  select(ins, skip)%>%
  tbl_svysummary(by = skip, 
              label = list(ins = "Have Insurance"))%>%
  add_p()%>%
  add_n()
Characteristic N 0, N = 132,303,2801 1, N = 16,085,0551 p-value2
Have Insurance 148,388,335 120,162,697 (91%) 9,954,056 (62%) <0.001

1 n (%)

2 chi-squared test with Rao & Scott's second-order correction

Differences Between Unweighted and Weighted Results

There are differences between the unweighted and weighted percentages of people who reported skipping going to the doctor for both the primary care and insurance variables. For the primary care variable, the weighted percentages of people with no primary care doctor who reported skipping going to the doctor and not skipping going to the doctor increased from the unweighted percentages. Similarly, the weighted percentages of people who reported not having insurance increased for both people who skipped going to the doctor and who did not skip going to the doctor compared to the unweighted analysis.