Loading libraries and the RECS file

library(haven)
library(readr)
library(dplyr)
## 
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
## 
##     filter, lag
## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union
library(ggplot2)
library(broom)
library(car)
## 
## Attaching package: 'car'
## The following object is masked from 'package:dplyr':
## 
##     recode
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)

#importing RECS 2015 data
recs2015v2<-read.csv("https://raw.githubusercontent.com/demograf/_Thesis/b9101bf3e4a94c9437557c6a7361e380384ce365/recs2015_public_v2.csv")

Creating a subset without missing values:

myrecs2<-subset(recs2015v2,ZMICRO=1) #has a microwave
myrecs2<-subset(myrecs2,ZOVEN=1) #has an oven
myrecs2<-subset(myrecs2,ZEDUCATION=1) #eductaion question answered
myrecs2<-subset(myrecs2,ZSTOVE=1) #has a stove
myrecs2<-subset(myrecs2,ZOVEN=1) #has an oven
myrecs2<-subset(myrecs2,ZMONEYPY=1) #reported income

Recoding educational attainment and household income:

myrecs2$myedu<-recode(myrecs2$EDUCATION,recodes="
                        1='No High School';
                        2='High School';
                        3='Some College';
                        4='Bachelors';
                        5='Masters or PhD'
                        ",as.factor.result=T)
myrecs2$myincome<-recode(myrecs2$MONEYPY,recodes="
                        1='Less than $20,000';
                        2='$20,000 - $39,999';
                        3='$40,000 - $59,999';
                        4='$60,000 - $79,999';
                        5='$80,000 - $99,999';
                        6='$100,000 - $119,999';
                        7='$120,000 - $139,999';
                        8='$140,000 or more'
                        ",as.factor.result=T)

Recoding missing values as zeros for calcuation of cooking methods (-2 indicates “NA”, thus it is effectively a not used case, or 0)

myrecs2$COOKTUSE<-recode(myrecs2$COOKTUSE, recodes = "-2:0=0")
myrecs2$OVENUSE<-recode(myrecs2$OVENUSE, recodes = "-2:0=0")
myrecs2$SEPCOOKTUSE<-recode(myrecs2$SEPCOOKTUSE, recodes = "-2:0=0")
myrecs2$SEPOVENUSE<-recode(myrecs2$SEPOVENUSE, recodes = "-2:0=0")
myrecs2$AMTMICRO<-recode(myrecs2$AMTMICRO, recodes = "-2:0=0")

Calculating the value to be used for binary indicator

myrecs2$mycooks<-myrecs2$COOKTUSE
myrecs2$mycooks<-myrecs2$mycooks+myrecs2$OVENUSE
myrecs2$mycooks<-myrecs2$mycooks+myrecs2$SEPCOOKTUSE
myrecs2$mycooks<-myrecs2$mycooks+myrecs2$SEPOVENUSE
myrecs2$mycooks<-myrecs2$mycooks-myrecs2$AMTMICRO

Recoding as binary variable:

myrecs2$mycooks<-recode(myrecs2$mycooks, recodes="0:396=1;-99:-1=0")
## 1 = prefers cooking with stove / stovetop
## 0 = uses microwave at least as often

Survey Design

options(survey.lonely.psu = "adjust")
des<-svydesign(ids=~1, weights=~myrecs2$NWEIGHT, data = myrecs2[is.na(myrecs2$myedu)==F,])
cat<-svytable(~myrecs2$mycooks+myrecs2$myedu, design = des)
##prop.table(svytable(~myrecs$mycooks+myrecs$myedu, design = des), margin = 2)
sv.table<-svyby(formula = ~mycooks, by = ~myedu, design = des, FUN = svymean, na.rm=T)
sv.table
##                         myedu   mycooks         se
## Bachelors           Bachelors 0.5722880 0.01639110
## High School       High School 0.5415268 0.01548395
## Masters or PhD Masters or PhD 0.5338177 0.01867883
## No High School No High School 0.6145932 0.02795179
## Some College     Some College 0.5361810 0.01318545
fit.logit<-svyglm(mycooks~myedu+myincome, design = des, family=binomial)
## Warning in eval(family$initialize): non-integer #successes in a binomial
## glm!
myexponent<-function(x) (exp(x))

stargazer(fit.logit, type = "text", apply.coef = myexponent, style="demography",covariate.labels=c("No High School","High School","Some College","Bachelors","Masters or PhD","Less than $20,000","$20,000 - $39,999","$40,000 - $59,999","$60,000 - $79,999","$80,000 - $99,999","$100,000 - $119,999","$120,000 - $139,999","$140,000 or more"))
## 
## ---------------------------------
##                        mycooks   
## ---------------------------------
## No High School         0.833***  
##                        (0.098)   
## High School            0.868***  
##                        (0.101)   
## Some College           1.118***  
##                        (0.144)   
## Bachelors              0.832***  
##                        (0.088)   
## Masters or PhD         0.887***  
##                        (0.173)   
## Less than 20,000       1.123***  
##                        (0.145)   
## 20,000 - 39,999        1.265***  
##                        (0.129)   
## 40,000 - 59,999        1.171***  
##                        (0.135)   
## 60,000 - 79,999        1.169***  
##                        (0.137)   
## 80,000 - 99,999        1.170***  
##                        (0.148)   
## 100,000 - 119,999      1.203***  
##                        (0.137)   
## 120,000 - 139,999      1.183***  
##                        (0.121)   
## N                       5,686    
## Log Likelihood        -4,091.140 
## AIC                   8,206.281  
## ---------------------------------
## *p < .05; **p < .01; ***p < .001
dat<-expand.grid(myedu=levels(myrecs2$myedu),myincome=levels(myrecs2$myincome))
fit<-predict(fit.logit,newdata = dat, type = "response")
dat$fitted.prob.lrm<-round(fit,3)
head(dat, n=40)
##             myedu            myincome fitted.prob.lrm
## 1       Bachelors $100,000 - $119,999           0.542
## 2     High School $100,000 - $119,999           0.496
## 3  Masters or PhD $100,000 - $119,999           0.507
## 4  No High School $100,000 - $119,999           0.569
## 5    Some College $100,000 - $119,999           0.496
## 6       Bachelors $120,000 - $139,999           0.512
## 7     High School $120,000 - $139,999           0.466
## 8  Masters or PhD $120,000 - $139,999           0.477
## 9  No High School $120,000 - $139,999           0.540
## 10   Some College $120,000 - $139,999           0.466
## 11      Bachelors    $140,000 or more           0.571
## 12    High School    $140,000 or more           0.525
## 13 Masters or PhD    $140,000 or more           0.536
## 14 No High School    $140,000 or more           0.598
## 15   Some College    $140,000 or more           0.525
## 16      Bachelors   $20,000 - $39,999           0.599
## 17    High School   $20,000 - $39,999           0.555
## 18 Masters or PhD   $20,000 - $39,999           0.565
## 19 No High School   $20,000 - $39,999           0.626
## 20   Some College   $20,000 - $39,999           0.554
## 21      Bachelors   $40,000 - $59,999           0.581
## 22    High School   $40,000 - $59,999           0.536
## 23 Masters or PhD   $40,000 - $59,999           0.546
## 24 No High School   $40,000 - $59,999           0.608
## 25   Some College   $40,000 - $59,999           0.535
## 26      Bachelors   $60,000 - $79,999           0.580
## 27    High School   $60,000 - $79,999           0.535
## 28 Masters or PhD   $60,000 - $79,999           0.546
## 29 No High School   $60,000 - $79,999           0.607
## 30   Some College   $60,000 - $79,999           0.535
## 31      Bachelors   $80,000 - $99,999           0.581
## 32    High School   $80,000 - $99,999           0.536
## 33 Masters or PhD   $80,000 - $99,999           0.546
## 34 No High School   $80,000 - $99,999           0.607
## 35   Some College   $80,000 - $99,999           0.535
## 36      Bachelors   Less than $20,000           0.587
## 37    High School   Less than $20,000           0.543
## 38 Masters or PhD   Less than $20,000           0.553
## 39 No High School   Less than $20,000           0.614
## 40   Some College   Less than $20,000           0.542
dat[which(dat$myedu=="Masters or PhD"&dat$myincome=="$120,000 - $139,999"),]
##            myedu            myincome fitted.prob.lrm
## 8 Masters or PhD $120,000 - $139,999           0.477
dat[which(dat$myedu=="No High School"&dat$myincome=="$20,000 - $39,999"),]
##             myedu          myincome fitted.prob.lrm
## 19 No High School $20,000 - $39,999           0.626

Question:

Are households with a maximum attainment level of Masters or PhD and higher income ($120k-139k) less likely to use a stove than a microwave for cooking compares to households with no High School and low income ($20k-29k)?

dat[which(dat$myedu=="Masters or PhD"&dat$myincome=="$120,000 - $139,999"),]
##            myedu            myincome fitted.prob.lrm
## 8 Masters or PhD $120,000 - $139,999           0.477
dat[which(dat$myedu=="No High School"&dat$myincome=="$20,000 - $39,999"),]
##             myedu          myincome fitted.prob.lrm
## 19 No High School $20,000 - $39,999           0.626

Answer:

Yes, only 48% of the time, the former use stoves and 63% of the latter use stoves.

(This is for future use:)

##fit.logit1<-svyglm(mycooks~myincome, design=des,family=binomial) #income only
##fit.logit2<-svyglm(mycooks~myeduc, design=des,family=binomial) #education only
##fit.logit3<-svyglm(mycooks~myedu+myincome, design=des,family=binomial) #both