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(magrittr)
library(car)
## Loading required package: carData
##
## Attaching package: 'car'
## The following object is masked from 'package:dplyr':
##
## recode
library(foreign)
#download datasets - AddHealth
wave1<-read.dta("/projects/add_health/data/11613344/ICPSR_27021/DS0001/da27021p1.dta")
wave3<-read.dta("/projects/add_health/data/11613344/ICPSR_27021/DS0003/da27021p3.dta")
wave4<-read.dta("/projects/add_health/data/11613344/ICPSR_27021/DS0012/da27021p12.dta")
wave4a<-read.dta("/projects/add_health/data/11613344/ICPSR_27021/DS0013/da27021p13.dta")
wave5<-read.xport("/projects/add_health/data/wave_5/wave5.xpt")
wave5age<-read.xport("/projects/add_health/data/wave_5/Constructed Files/Wave V Constructed Age/agew5.xpt")
#wave3section19<-read.dta("/projects/add_health/data/11613344/ICPSR_27021/DS0006/da27021p6.dta")
#download biomarkers Wave 4
w4glu <- read.xport("/projects/add_health/data/wave_4/biomarkers/glu_a1c.xpt")
w4crp <- read.xport("/projects/add_health/data/wave_4/biomarkers/crp_ebv.xpt")
w4lipids <- read.xport("/projects/add_health/data/wave_4/biomarkers/lipids.xpt")
#weightW4 <- read.dta("/projects/add_health/data/11613344/ICPSR_27021/DS0028/da27021p28.dta")
#weightW4<-read.xport("/projects/add_health/data/wave_4/Wave IV Grand Sample Weights/weights4/weights4.xpt")
#download biomarkers wave 5
w5brenal<-read.xport("/projects/add_health/data/wave_5/brenal5.xpt")
w5blipids<-read.xport("/projects/add_health/data/wave_5/blipids5.xpt")
w5bcardio<-read.xport("/projects/add_health/data/wave_5/bcardio5.xpt")
w5bglua1c<-read.xport("/projects/add_health/data/wave_5/bglua1c5.xpt")
w5banthro<-read.xport("/projects/add_health/data/wave_5/banthro5.xpt")
w5bcrp<-read.xport("/projects/add_health/data/wave_5/bcrp5.xpt")
w5bdemo<-read.xport("/projects/add_health/data/wave_5/bdemo5.xpt")
w5bweight<-read.xport("/projects/add_health/data/wave_5/bweight5.xpt")
#w5bweight=read.dta("/projects/add_health/data/wave_5/wave5weight.dta")
friendship <- read.dta("/projects/add_health/data/11613344/ICPSR_27033/DS0001/schoolnetworkvariables.dta")
w1school<-read.dta("/projects/add_health/data/11613344/ICPSR_27021/DS0019/da27021p19.dta")
w1school$s11<-w1school$S11
w1school$s17<-w1school$S17
w1school$s25<-w1school$S25
w1school$s26<-w1school$S26
w1school$s44<-w1school$S44
w1school$s46a<-w1school$S46A
w1school$s46d<-w1school$S46D
w1school2<-w1school%>%
dplyr::select(aid,s11,s17,s25,s26,s44,s46a,s46d)
library(data.table)
##
## Attaching package: 'data.table'
## The following objects are masked from 'package:dplyr':
##
## between, first, last
names(wave5)[1] <- "aid"
names(wave5age)[1] <- "aid"
names(w4glu)[1] <- "aid"
names(w4crp)[1] <- "aid"
names(w4lipids)[1] <- "aid"
#names(weightW4)[1] <- "aid"
names(w5brenal)[1] <- "aid"
names(w5blipids)[1] <- "aid"
names(w5bcardio)[1] <- "aid"
names(w5bglua1c)[1] <- "aid"
names(w5banthro)[1] <- "aid"
names(w5bcrp)[1] <- "aid"
names(w5bdemo)[1] <- "aid"
names(w5bweight)[1] <- "aid"
#change uppercase to lower case
#wave5$aid<-wave5$AID
#wave5age$aid<-wave5age$AID
#w4glu$aid<-w4glu$AID
#w4crp$aid<-w4crp$AID
#w4lipids$aid<-w4lipids$AID
#remove variables
#wave5 = subset(wave5, select = -c(AID) )
#w4glu = subset(w4glu, select = -c(AID) )
#w4crp = subset(w4crp, select = -c(AID) )
#w4lipids = subset(w4lipids, select = -c(AID) )
#merge datasets
a<-merge(wave1,wave3,by="aid")
bb<-merge(a,wave4,by="aid")
b<-merge(bb,wave4a,by="aid")
c<-merge(b,wave5,by="aid")
d<-merge(c,w4glu,by="aid")
e<-merge(d,w4crp,by="aid")
f<-merge(e,w4lipids,by="aid")
addhealtha<-merge(f,w5bweight,by="aid")
addhealthb<-merge(addhealtha,wave5age,by="aid")
addhealthc<-merge(addhealthb,w5brenal,by="aid")
addhealthd<-merge(addhealthc,w5blipids,by="aid")
addhealthe<-merge(addhealthd,w5bcardio,by="aid")
addhealthf<-merge(addhealthe,w5bglua1c,by="aid")
addhealthg<-merge(addhealthf,w5banthro,by="aid")
addhealthh<-merge(addhealthg,w5bcrp,by="aid")
addhealthi<-merge(addhealthh,w5bdemo,by="aid")
addhealthj<-merge(addhealthi,friendship,by="aid")
addhealth<-merge(addhealthj,w1school2,by="aid")
addhealth$livesMother<-Recode(addhealth$s11,recodes="0=0;1=1;else=NA")
addhealth$livesFather<-Recode(addhealth$s17,recodes="0=0;1=1;else=NA")
addhealth$livesParent<-ifelse(addhealth$livesMother==1&addhealth$livesFather==1,"lives_2_parents",
ifelse(addhealth$livesMother==1&addhealth$livesFather==0,"lives_mom",
ifelse(addhealth$livesMother==0&addhealth$livesFather==1,"lives_dad",
ifelse(addhealth$livesMother==0&addhealth$livesFather==0,"lives_other",NA))))
addhealth$livesParent2<-ifelse(addhealth$livesMother==1&addhealth$livesFather==1,"lives_2_parents",
ifelse(addhealth$livesMother==1&addhealth$livesFather==0,"lives_1_parent",
ifelse(addhealth$livesMother==0&addhealth$livesFather==1,"lives_1_parent",
ifelse(addhealth$livesMother==0&addhealth$livesFather==0,"lives_other",NA))))
addhealth$adopt<-ifelse(addhealth$s25==0,"no",
ifelse(addhealth$s25==1&addhealth$s26==0,"adoptNONbioparent",
ifelse(addhealth$s25==1&addhealth$s26==1,"adoptbioparent",NA)))
addhealth$schoolClubJoin<-addhealth$s44
addhealth$teacherTrouble<-Recode(addhealth$s46a,recodes="0:1=0;2:4=1;else=NA")
addhealth$studentTrouble<-Recode(addhealth$s46d,recodes="0:1=0;2:4=1;else=NA")
table(addhealth$schoolClubJoin)
##
## 0 1
## 4980 907
library(zoo)
##
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numeric
#age wave 4
#birth month wave 4
#addhealth$birthmonthw4<-Recode(addhealth$H4OD1M,recodes="1=1;2=2;3=3;4=4;5=5;6=6;7=7;8=8;9=9;10=10;11=11;12=12")
#birth year wave 4
#addhealth$birthyearw4<-addhealth$H4OD1Y
#addhealth$birthyearw4 <- Recode(addhealth$H4OD1Y,recodes="1974=1974;1975=1975;1976=1976;1977=1977;1978=1978;1979=1979;1980=1980;1981=1981;1982=1982;1983=1983")
#combine year and month for birth date wave 4
#addhealth$birthdatew4 <- as.yearmon(paste(addhealth$birthyearw4, addhealth$birthmonthw4), "%Y %m")
#interview date w4
#addhealth$interviewmonthw4 <- Recode(addhealth$IMONTH4,recodes="1=1;2=2;3=3;4=4;5=5;6=6;7=7;8=8;9=9;10=10;11=11;12=12")
#addhealth$interviewyearw4 <- Recode(addhealth$IYEAR4,recodes="2007=2007;2008=2008;2009=2009")
#addhealth$interviewdatew4<-as.yearmon(paste(addhealth$interviewyearw4,addhealth$interviewmonthw4),"%Y %m")
#age - derived from date of birth subtracted from wave interview date
#addhealth$ageW4<-(as.numeric(round(addhealth$interviewdatew4 - addhealth$birthdatew4)))
#table(addhealth$ageW4)
#Date of birth variables; Wave 1: H1GI1M (month), H1GI1Y (year); Wave 3: H3OD1M (month), H3OD1Y (year); Wave 4: H4OD1M (month), H4OD1Y (year)
addhealth$birthmonthw4 <- Recode(addhealth$H4OD1M,recodes="1=1;2=2;3=3;4=4;5=5;6=6;7=7;8=8;9=9;10=10;11=11;12=12")
addhealth$birthyearw1 <- factor(ifelse(addhealth$H1GI1Y=="74",1974,
ifelse(addhealth$H1GI1Y=="75",1975,
ifelse(addhealth$H1GI1Y=="76",1976,
ifelse(addhealth$H1GI1Y=="77",1977,
ifelse(addhealth$H1GI1Y=="78",1978,
ifelse(addhealth$H1GI1Y=="79",1979,
ifelse(addhealth$H1GI1Y=="80",1980,
ifelse(addhealth$H1GI1Y=="81",1981,
ifelse(addhealth$H1GI1Y=="82",1982,
ifelse(addhealth$H1GI1Y=="83",1983,1979)))))))))))
addhealth$birthmonthw1 <- Recode(addhealth$H1GI1M,recodes="1=1;2=2;3=3;4=4;5=5;6=6;7=7;8=8;9=9;10=10;11=11;12=12;else=NA")
addhealth$birthyearw4<-addhealth$H4OD1Y
addhealth$birthyearw4 <- Recode(addhealth$H4OD1Y,recodes="1974=1974;1975=1975;1976=1976;1977=1977;1978=1978;1979=1979;1980=1980;1981=1981;1982=1982;1983=1983")
addhealth$interviewmonthw4 <- Recode(addhealth$IMONTH4,recodes="1=1;2=2;3=3;4=4;5=5;6=6;7=7;8=8;9=9;10=10;11=11;12=12")
table(addhealth$interviewmonthw4)
##
## 1 2 3 4 5 6 7 8 9 10 11 12
## 277 1056 1141 1147 870 450 332 241 147 98 81 47
addhealth$interviewyearw4 <- Recode(addhealth$IYEAR4,recodes="2007=2007;2008=2008;2009=2009")
addhealth$iyearfix <- Recode(addhealth$iyear,recodes="'94'='1994';'95'='1995'")
addhealth$birthdatew4 <- as.yearmon(paste(addhealth$birthyearw4, addhealth$birthmonthw4), "%Y %m")
addhealth$interviewdatew4<-as.yearmon(paste(addhealth$interviewyearw4,addhealth$interviewmonthw4),"%Y %m")
addhealth$birthdatew3<-as.yearmon(paste(addhealth$H3OD1Y,addhealth$H3OD1M),"%Y %m")
addhealth$interviewdatew3<-as.yearmon(paste(addhealth$IYEAR3,addhealth$IMONTH3 ),"%Y %m")
addhealth$birthdatew1 <- as.yearmon(paste(addhealth$birthyearw1, addhealth$birthmonthw1), "%Y %m")
addhealth$interviewdatew1 <- as.yearmon(paste(addhealth$iyearfix, addhealth$imonth), "%Y %m")
addhealth$ageW4<-(as.numeric(round(addhealth$interviewdatew4 - addhealth$birthdatew4)))
addhealth$ageW3<-(as.numeric(round(addhealth$interviewdatew3 - addhealth$birthdatew3)))
addhealth$ageW1<-(as.numeric(round(addhealth$interviewdatew1 - addhealth$birthdatew1)))
addhealth$ageW5<-addhealth$AGE
#sex assigned at birth
addhealth$sex<-Recode(addhealth$H5OD2A,recodes="1='Male';2='Female';else=NA")
addhealth$sexW5<-factor(addhealth$sex,levels=c("Male","Female"))
#Sexual Orientation
# Wave 3
addhealth$sxorW3<-Recode(addhealth$H3SE13, recodes="1:2='Straight'; 3='Bisexual';4:5='Gay';else=NA")
addhealth$sexualOrientationW3 <- factor(addhealth$sxorW3, levels = c("Straight", "Gay", "Bisexual"))
# Wave 4
addhealth$sxorW4<-Recode(addhealth$H4SE31, recodes="1:2='Straight'; 3='Bisexual';4:5='Gay';else=NA")
addhealth$sexualOrientationW4 <- factor(addhealth$sxorW4, levels = c("Straight", "Gay", "Bisexual"))
# Wave 5
addhealth$sxorW5<-Recode(addhealth$H5SE15,recodes="1:2='Straight';3='Bisexual';4:5='Gay';else=NA")
addhealth$sexualOrientationW5 <- factor(addhealth$sxorW5, levels = c("Straight", "Gay", "Bisexual"))
#education wave 5
addhealth$educW5<-Recode(addhealth$H5OD11,recodes="1:2=0;3:9=1;10:16=2;else=NA")
addhealth$educW5<-as.factor(addhealth$educW5)
#employment status wave 1
#addhealth$H1EE3[addhealth$H1EE3 == 6] <- NA
#addhealth$H1EE3[addhealth$H1EE3 == 8] <- NA
#addhealth$employmentW1<-addhealth$H1EE3
#employment status wave 3
#addhealth$H3LM1[addhealth$H3LM1 == 6] <- NA
#addhealth$H3LM1[addhealth$H3LM1 == 8] <- NA
#addhealth$H3LM2[addhealth$H3LM2 == 6] <- NA
#addhealth$H3LM2[addhealth$H3LM2 == 7] <- NA
#addhealth$H3LM2[addhealth$H3LM2 == 8] <- NA
#addhealth$H3LM2[addhealth$H3LM2 == 9] <- NA
#addhealth$H3LM7[addhealth$H3LM7 == 7] <- NA
#addhealth$employmentW3<-ifelse(addhealth$H3LM1<2&addhealth$H3LM7>0,1,
# ifelse(addhealth$H3LM1<2&addhealth$H3LM7<1,0,
# ifelse(addhealth$H3LM2<2&addhealth$H3LM7>0,1,
# ifelse(addhealth$H3LM2<2&addhealth$H3LM7<1,0,NA))))
#employment status wave 4
#addhealth$H4LM1[addhealth$H4LM1 == 6] <- NA
#addhealth$H4LM2[addhealth$H4LM2 == 6] <- NA
#addhealth$H4LM2[addhealth$H4LM2 == 7] <- NA
#addhealth$employmentW4<-ifelse(addhealth$H4LM1<2&addhealth$H4LM14<11,0,
# ifelse(addhealth$H4LM2<2&addhealth$H4LM14<11,0,1))
#employment status Wave 5
#addhealth$employmentW5<-ifelse(addhealth$H5LM5==1,1,0)
#addhealth$employmentW1<-as.factor(addhealth$employmentW1)
#addhealth$employmentW3<-as.factor(addhealth$employmentW3)
#addhealth$employmentW4<-as.factor(addhealth$employmentW4)
#addhealth$employmentW5<-as.factor(addhealth$employmentW5)
#marital status Wave 1
addhealth$H1GI15[addhealth$H1GI15 == 6] <- NA
addhealth$H1GI15[addhealth$H1GI15 == 7] <- NA
addhealth$H1GI15[addhealth$H1GI15 == 8] <- NA
addhealth$maritalW1<-ifelse(addhealth$ageW1<15,0,
ifelse(addhealth$ageW1>14&addhealth$H1GI15<1,0,
ifelse(addhealth$ageW1>14&addhealth$H1GI15>0,1,NA)))
#marital status Wave 3 NOTE 0=never married, 1=married, 2=widow, 3=divorce
addhealth$H3MR1[addhealth$H3MR1 == 6] <- NA
addhealth$H3MR1[addhealth$H3MR1 == 8] <- NA
addhealth$H3MR1[addhealth$H3MR1 == 9] <- NA
addhealth$H3MR3_A[addhealth$H3MR3_A == 7] <- NA
addhealth$H3MR3_A[addhealth$H3MR3_A == 8] <- NA
addhealth$maritalW3<-ifelse(addhealth$H3MR1==0,0,
ifelse(addhealth$H3MR1>0&addhealth$H3MR3_A==1,1,2))
#marital status wave 4
addhealth$H4TR1[addhealth$H4TR1 == 6] <- NA
addhealth$H4TR1[addhealth$H4TR1 == 8] <- NA
addhealth$maritalW4<-ifelse(addhealth$H4TR1==0,0,
ifelse(addhealth$H4TR1>0&addhealth$H4TR13==1,1,
ifelse(addhealth$H4TR1>0&addhealth$H4TR13!=1,2,NA)))
#marital status wave 5
addhealth$maritalW5<-Recode(addhealth$H5HR1,recodes="1=1;2:4=2;5=0;else=NA")
#addhealth$maritalW1<-as.factor(addhealth$maritalW1)
#addhealth$maritalW3<-as.factor(addhealth$maritalW3)
#addhealth$maritalW4<-as.factor(addhealth$maritalW4)
#addhealth$maritalW5<-as.factor(addhealth$maritalW5)
#table(addhealth$maritalW5)
#combine race to one variable. This is for white, black, hispanic, asian, native american, pacific islander, and other.
addhealth$raceW5all <- factor(ifelse(addhealth$H5OD4A==1,"a_white",
ifelse(addhealth$H5OD4B==1, "b_black",
ifelse(addhealth$H5OD4C==1, "c_hispanic",
ifelse(addhealth$H5OD4D==1,"d_asian",
ifelse(addhealth$H5OD4F==1,"e_nativeamer",
ifelse(addhealth$H5OD4E==1,"f_pacificisl",
ifelse(addhealth$H5OD4G==1,"g_other",NA))))))))
#combine race to one variable. This is for white, black, hispanic, asian, native american, pacific islander, and other.
addhealth$raceW5min <- factor(ifelse(addhealth$H5OD4A==1,"a_white",
ifelse(addhealth$H5OD4B==1, "b_black",
ifelse(addhealth$H5OD4C==1, "c_hispanic",NA))))
library(StatMeasures)
#biomarkers wave 5
#1 systolic bp
addhealth$H5SBP<-Recode(addhealth$H5SBP,recodes="9999=NA;9997=NA;9996=NA")
addhealth$systolicbpW5<-ifelse(addhealth$H5SBP>=140,1,0)
#2 diastolic bp
addhealth$H5DBP<-Recode(addhealth$H5DBP,recodes="9999=NA;9997=NA;9996=NA")
addhealth$diastolicbpW5<-ifelse(addhealth$H5DBP>=90,1,0)
#3 pulse rate
addhealth$H5PR<-Recode(addhealth$H5PR,recodes="9999=NA;9997=NA;9996=NA")
addhealth$pulserateW5<-ifelse(addhealth$H5PR>=85,1,0)
#4 pulse pressure
addhealth$H5PP<-Recode(addhealth$H5PP,recodes="9999=NA;9997=NA;9996=NA")
addhealth$pulsepressureW5<-ifelse(addhealth$H5PP>=60,1,0)
#5 TOTAL CHOLESTEROL W5
addhealth$tcW5<-decile(addhealth$H5TC)
addhealth$totalcholW5<-ifelse(addhealth$tcW5>=8,1,0)
table(addhealth$totalcholW5)
##
## 0 1
## 3741 1519
#6 HDL, HIGH-DENSITY LIPOPROTEIN CHOLESTEROL
addhealth$H5HDL<-Recode(addhealth$H5HDL,recodes="9999=NA")
addhealth$hdlW5<-decile(addhealth$H5HDL)
addhealth$hidensity_lipW5<-ifelse(addhealth$hdlW5<=2,1,0)
#7 TRIGLYCERIDES
addhealth$TG<-decile(addhealth$H5TG)
addhealth$triglycerideW5<-ifelse(addhealth$TG>=8,1,0)
#8 HEMOGLOBIN
addhealth$hemoglobinW5<-ifelse(addhealth$H5HBA1C>=6.4,1,0)
table(addhealth$hemoglobinW5)
##
## 0 1
## 4901 246
#9 glucose
addhealth$GLUCOSE<-Recode(addhealth$H5GLUCOS,recodes="9999=NA")
addhealth$glucoseW5<-ifelse(addhealth$GLUCOSE>=200,1,0)
table(addhealth$glucoseW5)
##
## 0 1
## 5123 132
#10 BMI
addhealth$H5BMI<-Recode(addhealth$H5BMI,recodes="9999=NA;9996=NA;9997=NA;9994=NA")
addhealth$bmiW5<-ifelse(addhealth$H5BMI>=30,1,0)
table(addhealth$bmiW5)
##
## 0 1
## 2982 2780
#11 waist circumference
addhealth$waistW5<-Recode(addhealth$H5WSTCLS,recodes="1=0;2=1")
#12 high-sensitivity C-reactive protein test.
addhealth$hsCRPW5<-ifelse(addhealth$H5CRP>=3,1,0)
table(addhealth$hsCRPW5)
##
## 0 1
## 3124 1938
# COMBINE BIOMARKERS
addhealth$scoreW5<-addhealth$systolicbpW5+addhealth$diastolicbpW5+addhealth$pulserateW5+addhealth$pulsepressureW5+addhealth$totalcholW5+addhealth$hidensity_lipW5+addhealth$triglycerideW5+addhealth$hemoglobinW5+addhealth$glucoseW5+addhealth$bmiW5+addhealth$waistW5+addhealth$hsCRPW5
table(addhealth$scoreW5)
##
## 0 1 2 3 4 5 6 7 8 9 10
## 784 687 761 788 608 418 337 169 76 32 20
#biomarkers Wave 4
#1 systolic bp
addhealth$H4SBP<-Recode(addhealth$H4SBP,recodes="999=NA;997=NA;996=NA")
addhealth$systolicbp<-ifelse(addhealth$H4SBP>=140,1,0)
#2 diastolic bp
addhealth$H4DBP<-Recode(addhealth$H4DBP,recodes="999=NA;997=NA;996=NA")
addhealth$diastolicbp<-ifelse(addhealth$H4DBP>=90,1,0)
#3 pulse rate
addhealth$H4PR<-Recode(addhealth$H4PR,recodes="999=NA;997=NA;996=NA")
addhealth$pulserate<-ifelse(addhealth$H4PR>=85,1,0)
#4 pulse pressure
addhealth$H4PP<-Recode(addhealth$H4PP,recodes="999=NA;997=NA;996=NA")
addhealth$pulsepressure<-ifelse(addhealth$H4PP>=60,1,0)
#5 TOTAL CHOLESTEROL
addhealth$TC<-Recode(addhealth$TC,recodes="99=NA")
addhealth$totalchol<-ifelse(addhealth$TC>=8,1,0)
#6 HDL, HIGH-DENSITY LIPOPROTEIN CHOLESTEROL
addhealth$HDL<-Recode(addhealth$HDL,recodes="99=NA")
addhealth$hidensity_lip<-ifelse(addhealth$HDL<=2,1,0)
#7 TRIGLYCERIDES
addhealth$TG<-Recode(addhealth$TG,recodes="99=NA")
addhealth$triglyceride<-ifelse(addhealth$TG>=8,1,0)
#8 HEMOGLOBIN
addhealth$HBA1C<-Recode(addhealth$HBA1C,recodes="99=NA")
addhealth$hemoglobin<-ifelse(addhealth$HBA1C>=6.4,1,0)
#9 glucose
addhealth$GLUCOSE<-Recode(addhealth$GLUCOSE,recodes="999=NA")
addhealth$glucose<-ifelse(addhealth$GLUCOSE>=200,1,0)
#10 BMI
addhealth$H4BMI<-Recode(addhealth$H4BMI,recodes="999=NA;888=NA;889=NA;996=NA;997=NA")
addhealth$BMI<-ifelse(addhealth$H4BMI>=30,1,0)
# 11 high-sensitivity C-reactive protein test.
addhealth$CRPa<-Recode(addhealth$CRP,recodes="999=NA;998=NA")
addhealth$hsCRP<-ifelse(addhealth$CRPa>=3,1,0)
#EPSTEIN BARR
#quantile(addhealth$EBVa, c(.80), na.rm = T)
#addhealth$EBVa<-Recode(addhealth$EBV,recodes="9999=NA;9998=NA")
#addhealth$EpsteinBar<-ifelse(addhealth$EBVa>=243,1,0)
#12 waist circumference
addhealth$H4WAIST[addhealth$H4WAIST == 996] <- NA
addhealth$H4WAIST[addhealth$H4WAIST == 997] <- NA
addhealth$H4WAIST[addhealth$H4WAIST == 998] <- NA
addhealth$waist<-ifelse(addhealth$H5OD2A==2&addhealth$H4WAIST>88,1,
ifelse(addhealth$H5OD2A==2&addhealth$H4WAIST<89,0,
ifelse(addhealth$H5OD2A==1&addhealth$H4WAIST>102,1,
ifelse(addhealth$H5OD2A==1&addhealth$H4WAIST<103,0,NA))))
#combine biomarkers
addhealth$scoreW4<-addhealth$hsCRP+addhealth$BMI+addhealth$triglyceride+addhealth$hemoglobin+addhealth$hidensity_lip+addhealth$totalchol+addhealth$pulserate+addhealth$systolicbp+addhealth$pulsepressure+addhealth$diastolicbp+addhealth$glucose+addhealth$waist
summary(addhealth$HHincomeW3)
## Length Class Mode
## 0 NULL NULL
#household income wave 3
#addhealth$HHincomeW3<-Recode(addhealth$H3EC6)
addhealth$HHincomeW3<-ifelse(addhealth$H3EC6>0&addhealth$H3EC6<30000,'a_<30K',
ifelse(addhealth$H3EC6>29999&addhealth$H3EC6<50000,'b_30k<50k',
ifelse(addhealth$H3EC6>49999&addhealth$H3EC6<100000,'c_50k>100k',
ifelse(addhealth$H3EC6>99999,'d_>100k',NA))))
addhealth$HHincomeW3<-as.factor(addhealth$HHincomeW3)
#household income wave 5
addhealth$HHincomeW5<- Recode(addhealth$H5EC2,recodes="1:6='a_<30K';7:8='b_30k<50k';9:10='c_50k>100k';11:13='d_>100k';else=NA")
#household income wave 4
addhealth$HHincomeW4<- Recode(addhealth$H4EC1,recodes="1:6='a_<30K';7:8='b_30k<50k';9:10='c_50k>100k';11:12='d_>100k';else=NA")
#percieved closeness mom WAVE 1
addhealth$close.momW1<-Recode(addhealth$H1WP9,recodes="1=1;2=2;3=3;4=4;5=5;else=NA")
#percieved closeness mom WAVE 3
addhealth$close.momW3<-Recode(addhealth$H3WP20,recodes="5=1;4=2;3=3;2=4;1=5;else=NA")
#percieved closeness mom WAVE 4
addhealth$close.momW4<-Recode(addhealth$H4WP24,recodes="1=1;2=2;3=3;4=4;5=5;else=NA")
#percieved closeness mom WAVE 5
addhealth$close.momW5<-Recode(addhealth$H5WP14,recodes="1=1;2=2;3=3;4=4;5=5;else=NA")
#parent SUPPORT/SUPERVISION involvement
addhealth$curfew<-Recode(addhealth$H1WP1,recodes="0=0;1=1;else=NA")
addhealth$hang<-Recode(addhealth$H1WP2,recodes="0=0;1=1;else=NA")
addhealth$wear<-Recode(addhealth$H1WP3,recodes="0=0;1=1;else=NA")
addhealth$tv<-Recode(addhealth$H1WP4,recodes="0=0;1=1;else=NA")
addhealth$prog<-Recode(addhealth$H1WP5,recodes="0=0;1=1;else=NA")
addhealth$bed<-Recode(addhealth$H1WP6,recodes="0=0;1=1;else=NA")
addhealth$eat<-Recode(addhealth$H1WP7,recodes="0=0;1=1;else=NA")
addhealth$parentSupervise<-addhealth$curfew+addhealth$hang+addhealth$wear+addhealth$tv+addhealth$prog+addhealth$bed+addhealth$eat
#mothers education
addhealth$mom.educ<-Recode(addhealth$H1RM1,recodes="1:3=0;4:7=1;8:9=2;else=NA")
addhealth$mom.educ2<-Recode(addhealth$H1RM1,recodes="1:3=0;4:9=1;else=NA")
#mothers nativity
addhealth$mom.nativity<-Recode(addhealth$H1RM2,recodes="0=0;1=1;else=NA")
#mother employed
addhealth$mom.employed<-Recode(addhealth$H1RM4,recodes="16=0;1:15=1;else=NA")
#mother accecpt public assistance
addhealth$mom.pubass<-Recode(addhealth$H1RM9,recodes="0=0;1=1;else=NA")
addhealth$motherInequity<-addhealth$mom.educ2+addhealth$mom.nativity+addhealth$mom.employed+addhealth$mom.pubass
addhealth$relw1<-Recode(addhealth$H1RE1,recodes="0=0;1:28=1;else=NA")
addhealth$religionW1<-ifelse(addhealth$relw1==0,0,
ifelse(addhealth$relw1<2&addhealth$H1RE3==4,0,
ifelse(addhealth$relw1<2&addhealth$H1RE3==3,0,
ifelse(addhealth$relw1<2&addhealth$H1RE3==2,0,
ifelse(addhealth$relw1<2&addhealth$H1RE3==1,1,NA)))))
addhealth$religionW3<-Recode(addhealth$H3RE24,recodes="0:4=0;5:6=1;else=NA")
addhealth$religionW4<-Recode(addhealth$H4RE7,recodes="0:3=0;4:5=1;else=NA")
addhealth$religionW5<-Recode(addhealth$H5RE2,recodes="1:4=0;5:6=1")
#school prestige friends
addhealth$pre<-addhealth$prxprest
addhealth$schoolSeg<-addhealth$SEG1RCE5
#education level of 0=less high school, 1=highschool,
#education wave 3
addhealth$educW3<-ifelse(addhealth$H3ED2==1,1,
ifelse(addhealth$H3ED3==1,1,
ifelse(addhealth$H3ED5==1,2,
ifelse(addhealth$H3ED6==1&addhealth$H3ED5!=1,2,
ifelse(addhealth$H3ED7==1&addhealth$H3ED5!=1,2,
ifelse(addhealth$H3ED8==1&addhealth$H3ED5!=1,2,
ifelse(addhealth$H3ED9==1,0,NA)))))))
#education wave 4
addhealth$educW4<-Recode(addhealth$H4ED2,recodes="1:2=0;3:6=1;7:13=2;else=NA")
#education wave 5
addhealth$educW5<-Recode(addhealth$H5OD11,recodes="1:2=0;3:9=1;10:16=2;else=NA")
#education wave 1
addhealth$educW1<-addhealth$H1GI18
addhealth$educW1<-as.factor(addhealth$educW1)
addhealth$educW3<-as.factor(addhealth$educW3)
addhealth$educW4<-as.factor(addhealth$educW4)
addhealth$educW5<-as.factor(addhealth$educW5)
#employment
#employment binary scaled
addhealth$currently.employed<-Recode(addhealth$H5LM5,recodes="1=1;2:3=0;else=NA")
#employer provides insurance, benefits, sick leave
addhealth$H5LM13A<-Recode(addhealth$H5LM13A,recodes="0=0;1=1;else=NA")
addhealth$H5LM13B<-Recode(addhealth$H5LM13B,recodes="0=0;1=1;else=NA")
addhealth$H5LM13C<-Recode(addhealth$H5LM13C,recodes="0=0;1=1;else=NA")
#scale benefits
addhealth$benefitscore<-addhealth$H5LM13A+addhealth$H5LM13B+addhealth$H5LM13C
addhealth$beneScore<-Recode(addhealth$benefitscore,recodes="0=0;1:3=1")
table(addhealth$benefitscore)
##
## 0 1 2 3
## 604 380 351 3406
table(addhealth$beneScore)
##
## 0 1
## 604 4137
addhealth$employment<-ifelse(addhealth$currently.employed==0,"unemployed",
ifelse(addhealth$currently.employed==1&addhealth$beneScore==0,"empNoBene",
ifelse(addhealth$currently.employed==1&addhealth$beneScore==1,"empBene",NA)))
addhealth$employment<-factor(addhealth$employment,levels=c("unemployed","empNoBene","empBene"))
table(addhealth$employment)
##
## unemployed empNoBene empBene
## 1084 604 4137
#install.packages("magrittr")
#library(data.table)
addhealthSub<-addhealth%>%
dplyr::select(aid,W5BIOWGT,raceW5min,sexW5,scoreW4,scoreW5,ageW1,ageW3,ageW4,ageW5,sexualOrientationW3,sexualOrientationW4,sexualOrientationW5,educW1,educW3,educW4,educW5,employment,maritalW1,maritalW3,maritalW4,maritalW5,religionW1,religionW3,religionW4,religionW5,HHincomeW3,HHincomeW4,HHincomeW5,livesParent2,close.momW1,close.momW3,close.momW4,close.momW5,parentSupervise,mom.educ,mom.educ2,mom.nativity,mom.employed,mom.pubass,pre,teacherTrouble,studentTrouble,schoolClubJoin,schoolSeg,scoreW4,scoreW5)
out<-melt(setDT(addhealthSub), id = "aid",
measure.vars = list(age=c("ageW1","ageW3","ageW4","ageW5"),
sexualorientation=c("sexualOrientationW3","sexualOrientationW3", "sexualOrientationW4", "sexualOrientationW5"),
education=c("educW1","educW3", "educW4", "educW5"),
hhincome=c("HHincomeW3","HHincomeW4", "HHincomeW4", "HHincomeW5"),
marriage=c("maritalW1","maritalW3","maritalW4","maritalW5"),
religion=c("religionW1","religionW3","religionW4","religionW5"),
closeMom=c("close.momW1","close.momW3","close.momW4","close.momW5"),
score=c("scoreW4","scoreW4","scoreW4","scoreW5")))%>%
setorder(aid)
addhealthLong<-addhealthSub%>%
dplyr::select(aid, W5BIOWGT,raceW5min,sexW5,livesParent2,parentSupervise,mom.educ,mom.educ2,mom.nativity,mom.employed,mom.pubass,pre,teacherTrouble,studentTrouble,schoolClubJoin,schoolSeg,employment)%>%
left_join(., out, "aid")
addhealthLong$wave<-Recode(addhealthLong$variable,recodes="1=1;2=3;3=4;4=5")
addhealthS1long <- addhealthLong[,c(27,1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,19,20,21,22,23,24,25,26)]
addhealthS1long$sexualorientation <- factor(addhealthS1long$sexualorientation, levels = c("Straight", "Gay", "Bisexual"))
addhealthS1long$raceXsexo<-interaction(addhealthS1long$raceW5min,addhealthS1long$sexualorientation,drop=TRUE)
addhealthS1long$education<-as.factor(addhealthS1long$education)
addhealthS1long$hhincome<-as.factor(addhealthS1long$hhincome)
#addhealthS1long$employment<-as.factor(addhealthS1long$employment)
addhealthS1long$marriage<-as.factor(addhealthS1long$marriage)
addhealthS1long$mom.educ<-as.factor(addhealthS1long$mom.educ)
addhealthS1long$mom.educ2<-as.factor(addhealthS1long$mom.educ2)
addhealthS1long$mom.employed<-as.factor(addhealthS1long$mom.employed)
addhealthS1long$mom.pubass<-as.factor(addhealthS1long$mom.pubass)
addhealthS1long$mom.nativity<-as.factor(addhealthS1long$mom.nativity)
addhealthS1long$religion<-as.factor(addhealthS1long$religion)
addhealthS1long$livesParent2<-as.factor(addhealthS1long$livesParent2)
addhealthS1long$schoolClubJoin<-as.factor(addhealthS1long$schoolClubJoin)
addhealthS1long$teacherTrouble<-as.factor(addhealthS1long$teacherTrouble)
addhealthS1long$studentTrouble<-as.factor(addhealthS1long$studentTrouble)
addhealthS1long$parentSupervise<-as.factor(addhealthS1long$parentSupervise)
addhealthS1long$closeMom<-as.factor(addhealthS1long$closeMom)
addhealthexp<-addhealthS1long%>%
dplyr::select(wave,aid,W5BIOWGT,age,sexW5,education,hhincome,employment,marriage,livesParent2,parentSupervise,mom.educ2,mom.educ,mom.employed,mom.nativity,mom.pubass,pre,teacherTrouble,studentTrouble,schoolClubJoin,schoolSeg,sexualorientation,religion,closeMom,score,raceXsexo)
#library(tidyr)
addhealthBOMB<-addhealthexp%>%
filter(complete.cases(.))
addhealthBOMB$livesParent2<-factor(addhealthBOMB$livesParent2,levels=c("lives_2_parents","lives_1_parent","lives_other"))
addhealthBOMB$closeMom<-as.numeric(addhealthBOMB$closeMom)
addhealthBOMB$closeMom2<-Recode(addhealthBOMB$closeMom,recodes="1:2='a_notClose';3:5='b_close'")
addhealthBOMB$parentSupervise<-as.numeric(addhealthBOMB$parentSupervise)
#table(addhealthBOMB$parentSupervise)
addhealthBOMB$parentSupervise2<-Recode(addhealthBOMB$parentSupervise,recodes="1:3='a_easy';4:6='b_moderate';7:8='c_strict'")
addhealthBOMB$teacherTrouble2<-as.numeric(addhealthBOMB$teacherTrouble)
addhealthBOMB$studentTrouble2<-as.numeric(addhealthBOMB$studentTrouble)
addhealthBOMB$schoolconflict<-addhealthBOMB$studentTrouble2+addhealthBOMB$teacherTrouble2
addhealthBOMB$schoolconflict<-as.factor(addhealthBOMB$schoolconflict)
addhealthBOMB$parentSupervise5<-as.factor(Recode(addhealthBOMB$parentSupervise,recodes="1:6='a';7='b';8='c'"))
addhealthBOMB$closeMom2<-Recode(addhealthBOMB$closeMom,recodes="1:2='a_notClose';3:5='b_close'")
#install.packages("geepack")
library(geepack)
#wave,aid,W5BIOWGT,age,sexW5,education,hhincome,employment,marriage,,,,,,,,pre,teacherTrouble,studentTrouble,schoolClubJoin,schoolSeg,
fit<-geeglm(score~age+I(age^2)+sexW5+raceXsexo+education+employment+hhincome+marriage,id=aid , wave = wave, corstr ="exchangeable", data=addhealthBOMB, weights=W5BIOWGT/mean(W5BIOWGT,na.rm=T))
summary(fit)
##
## Call:
## geeglm(formula = score ~ age + I(age^2) + sexW5 + raceXsexo +
## education + employment + hhincome + marriage, data = addhealthBOMB,
## weights = W5BIOWGT/mean(W5BIOWGT, na.rm = T), id = aid, waves = wave,
## corstr = "exchangeable")
##
## Coefficients:
## Estimate Std.err Wald Pr(>|W|)
## (Intercept) 3.9306210 0.6364701 38.139 6.59e-10 ***
## age -0.0876183 0.0379190 5.339 0.02085 *
## I(age^2) 0.0020084 0.0007521 7.130 0.00758 **
## sexW5Female 0.1637235 0.1214367 1.818 0.17759
## raceXsexob_black.Straight 0.3667933 0.1538609 5.683 0.01713 *
## raceXsexoc_hispanic.Straight 0.2203025 0.2448025 0.810 0.36816
## raceXsexoa_white.Gay -0.3916751 0.1585820 6.100 0.01352 *
## raceXsexob_black.Gay 0.1604503 0.2207732 0.528 0.46737
## raceXsexoc_hispanic.Gay 0.8854806 0.7063889 1.571 0.21001
## raceXsexoa_white.Bisexual 0.3046398 0.1869142 2.656 0.10314
## raceXsexob_black.Bisexual -0.2025419 0.3979780 0.259 0.61080
## raceXsexoc_hispanic.Bisexual 0.2313986 0.3067633 0.569 0.45066
## education1 -0.1446085 0.2038388 0.503 0.47806
## education2 -0.4330877 0.2213278 3.829 0.05037 .
## employmentempNoBene -0.3708867 0.2502034 2.197 0.13825
## employmentempBene -0.4375551 0.1802724 5.891 0.01522 *
## hhincomeb_30k<50k -0.3013741 0.1581239 3.633 0.05666 .
## hhincomec_50k>100k -0.0192379 0.1552076 0.015 0.90136
## hhincomed_>100k -0.0960963 0.1867552 0.265 0.60686
## marriage1 0.0168895 0.1052325 0.026 0.87249
## marriage2 0.2330924 0.2060255 1.280 0.25790
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation structure = exchangeable
## Estimated Scale Parameters:
##
## Estimate Std.err
## (Intercept) 3.598 0.262
## Link = identity
##
## Estimated Correlation Parameters:
## Estimate Std.err
## alpha 0.7861 0.1519
## Number of clusters: 2063 Maximum cluster size: 450
#education+employment+hhincome+marriage
#fit2<-geeglm(score~age+I(age^2)+sexW5+raceXsexo+education+employment+hhincome+marriage+livesParent2+mom.educ+mom.pubass+closeMom2+parentSupervise2+religion,id=aid , wave = wave, corstr ="exchangeable", data=addhealthBOMB, weights=W5BIOWGT/mean(W5BIOWGT,na.rm=T))
#summary(fit2)
fit3<-geeglm(score~age+I(age^2)+sexW5+raceXsexo+education+employment+hhincome+marriage+livesParent2+mom.educ+mom.pubass+closeMom2+parentSupervise5+religion+schoolconflict+schoolClubJoin+scale(schoolSeg)+scale(pre),id=aid , wave = wave, corstr ="exchangeable", data=addhealthBOMB, weights=W5BIOWGT/mean(W5BIOWGT,na.rm=T))
summary(fit3)
##
## Call:
## geeglm(formula = score ~ age + I(age^2) + sexW5 + raceXsexo +
## education + employment + hhincome + marriage + livesParent2 +
## mom.educ + mom.pubass + closeMom2 + parentSupervise5 + religion +
## schoolconflict + schoolClubJoin + scale(schoolSeg) + scale(pre),
## data = addhealthBOMB, weights = W5BIOWGT/mean(W5BIOWGT, na.rm = T),
## id = aid, waves = wave, corstr = "exchangeable")
##
## Coefficients:
## Estimate Std.err Wald Pr(>|W|)
## (Intercept) 4.34458 0.69666 38.89 4.5e-10 ***
## age -0.09504 0.04086 5.41 0.0200 *
## I(age^2) 0.00210 0.00079 7.10 0.0077 **
## sexW5Female 0.13201 0.11942 1.22 0.2690
## raceXsexob_black.Straight 0.20509 0.16204 1.60 0.2056
## raceXsexoc_hispanic.Straight -0.04934 0.26740 0.03 0.8536
## raceXsexoa_white.Gay -0.39123 0.16138 5.88 0.0153 *
## raceXsexob_black.Gay -0.05195 0.22722 0.05 0.8191
## raceXsexoc_hispanic.Gay 0.48481 0.63630 0.58 0.4461
## raceXsexoa_white.Bisexual 0.31378 0.17731 3.13 0.0768 .
## raceXsexob_black.Bisexual -0.42928 0.40777 1.11 0.2925
## raceXsexoc_hispanic.Bisexual -0.13422 0.31742 0.18 0.6724
## education1 -0.17261 0.20664 0.70 0.4035
## education2 -0.45428 0.22753 3.99 0.0459 *
## employmentempNoBene -0.34157 0.24475 1.95 0.1628
## employmentempBene -0.40643 0.17670 5.29 0.0214 *
## hhincomeb_30k<50k -0.28508 0.14960 3.63 0.0567 .
## hhincomec_50k>100k -0.02468 0.15512 0.03 0.8736
## hhincomed_>100k -0.09408 0.17884 0.28 0.5988
## marriage1 0.02606 0.10822 0.06 0.8097
## marriage2 0.24479 0.20325 1.45 0.2284
## livesParent2lives_1_parent 0.23396 0.16038 2.13 0.1446
## livesParent2lives_other -0.10902 0.32912 0.11 0.7405
## mom.educ1 -0.36376 0.16828 4.67 0.0306 *
## mom.educ2 -0.49191 0.18630 6.97 0.0083 **
## mom.pubass1 0.36861 0.23363 2.49 0.1146
## closeMom2b_close -0.05646 0.18158 0.10 0.7558
## parentSupervise5b 0.21650 0.12240 3.13 0.0769 .
## parentSupervise5c 0.33844 0.15702 4.65 0.0311 *
## religion1 -0.15121 0.12327 1.50 0.2199
## schoolconflict3 0.04708 0.13267 0.13 0.7227
## schoolconflict4 -0.03325 0.15327 0.05 0.8283
## schoolClubJoin1 0.28519 0.18497 2.38 0.1231
## scale(schoolSeg) 0.09753 0.06014 2.63 0.1049
## scale(pre) -0.08613 0.04387 3.85 0.0496 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation structure = exchangeable
## Estimated Scale Parameters:
##
## Estimate Std.err
## (Intercept) 3.59 0.255
## Link = identity
##
## Estimated Correlation Parameters:
## Estimate Std.err
## alpha 0.824 0.133
## Number of clusters: 2063 Maximum cluster size: 450
#fit4<-geeglm(score~age+I(age^2)+sexW5+raceXsexo+livesParent2+mom.educ+mom.pubass+closeMom2+parentSupervise2+religion+schoolconflict+schoolClubJoin+schoolSeg+pre+education+employment+hhincome+marriage,id=aid , wave = wave, corstr ="exchangeable", data=addhealthBOMB, weights=W5BIOWGT/mean(W5BIOWGT,na.rm=T))
#summary(fit4)
library(huxtable)
##
## Attaching package: 'huxtable'
## The following object is masked from 'package:StatMeasures':
##
## contents
## The following object is masked from 'package:dplyr':
##
## add_rownames
allmodels<-huxreg(fit,fit3,statistics = character(0))
allmodels
| (1) | (2) |
| (Intercept) | 3.931 *** | 4.345 *** |
| (0.636) | (0.697) |
| age | -0.088 * | -0.095 * |
| (0.038) | (0.041) |
| I(age^2) | 0.002 ** | 0.002 ** |
| (0.001) | (0.001) |
| sexW5Female | 0.164 | 0.132 |
| (0.121) | (0.119) |
| raceXsexob_black.Straight | 0.367 * | 0.205 |
| (0.154) | (0.162) |
| raceXsexoc_hispanic.Straight | 0.220 | -0.049 |
| (0.245) | (0.267) |
| raceXsexoa_white.Gay | -0.392 * | -0.391 * |
| (0.159) | (0.161) |
| raceXsexob_black.Gay | 0.160 | -0.052 |
| (0.221) | (0.227) |
| raceXsexoc_hispanic.Gay | 0.885 | 0.485 |
| (0.706) | (0.636) |
| raceXsexoa_white.Bisexual | 0.305 | 0.314 |
| (0.187) | (0.177) |
| raceXsexob_black.Bisexual | -0.203 | -0.429 |
| (0.398) | (0.408) |
| raceXsexoc_hispanic.Bisexual | 0.231 | -0.134 |
| (0.307) | (0.317) |
| education1 | -0.145 | -0.173 |
| (0.204) | (0.207) |
| education2 | -0.433 | -0.454 * |
| (0.221) | (0.228) |
| employmentempNoBene | -0.371 | -0.342 |
| (0.250) | (0.245) |
| employmentempBene | -0.438 * | -0.406 * |
| (0.180) | (0.177) |
| hhincomeb_30k<50k | -0.301 | -0.285 |
| (0.158) | (0.150) |
| hhincomec_50k>100k | -0.019 | -0.025 |
| (0.155) | (0.155) |
| hhincomed_>100k | -0.096 | -0.094 |
| (0.187) | (0.179) |
| marriage1 | 0.017 | 0.026 |
| (0.105) | (0.108) |
| marriage2 | 0.233 | 0.245 |
| (0.206) | (0.203) |
| livesParent2lives_1_parent | | 0.234 |
| | (0.160) |
| livesParent2lives_other | | -0.109 |
| | (0.329) |
| mom.educ1 | | -0.364 * |
| | (0.168) |
| mom.educ2 | | -0.492 ** |
| | (0.186) |
| mom.pubass1 | | 0.369 |
| | (0.234) |
| closeMom2b_close | | -0.056 |
| | (0.182) |
| parentSupervise5b | | 0.216 |
| | (0.122) |
| parentSupervise5c | | 0.338 * |
| | (0.157) |
| religion1 | | -0.151 |
| | (0.123) |
| schoolconflict3 | | 0.047 |
| | (0.133) |
| schoolconflict4 | | -0.033 |
| | (0.153) |
| schoolClubJoin1 | | 0.285 |
| | (0.185) |
| scale(schoolSeg) | | 0.098 |
| | (0.060) |
| scale(pre) | | -0.086 * |
| | (0.044) |
| *** p < 0.001; ** p < 0.01; * p < 0.05. |
nobs(fit3)
## [1] 24657
addhealthwh<-addhealthS1long%>%
dplyr::select(wave,aid,W5BIOWGT,age,sexW5,education,hhincome,employment,marriage,livesParent2,parentSupervise,mom.educ2,mom.educ,mom.employed,mom.nativity,mom.pubass,pre,teacherTrouble,studentTrouble,schoolClubJoin,schoolSeg,sexualorientation,religion,closeMom,score,raceXsexo,raceW5min)
addhealthwhi<-addhealthwh%>%
subset(raceW5min=='a_white')
addhealthWhite<-addhealthwhi%>%
filter(complete.cases(.))
addhealthWhite$livesParent2<-factor(addhealthWhite$livesParent2,levels=c("lives_2_parents","lives_1_parent","lives_other"))
addhealthWhite$closeMom<-as.numeric(addhealthWhite$closeMom)
addhealthWhite$closeMom2<-Recode(addhealthWhite$closeMom,recodes="1:2='a_notClose';3:5='b_close'")
addhealthWhite$parentSupervise<-as.numeric(addhealthWhite$parentSupervise)
#table(addhealthWhite$parentSupervise)
addhealthWhite$parentSupervise2<-Recode(addhealthWhite$parentSupervise,recodes="1:3='a_easy';4:6='b_moderate';7:8='c_strict'")
addhealthWhite$parentSupervise5<-as.factor(Recode(addhealthWhite$parentSupervise,recodes="1:6='a';7='b';8='c'"))
addhealthWhite$teacherTrouble2<-as.numeric(addhealthWhite$teacherTrouble)
addhealthWhite$studentTrouble2<-as.numeric(addhealthWhite$studentTrouble)
addhealthWhite$schoolconflict<-addhealthWhite$studentTrouble2+addhealthWhite$teacherTrouble2
addhealthWhite$schoolconflict<-as.factor(addhealthWhite$schoolconflict)
addhealthWhite$raceXsexo<-as.factor(addhealthWhite$raceXsexo)
addhealthWhite$raceXsexo <- droplevels(addhealthWhite$raceXsexo)
table(addhealthWhite$raceXsexo)
##
## a_white.Straight a_white.Gay a_white.Bisexual
## 15312 490 609
fitWh<-geeglm(score~age+I(age^2)+sexW5+raceXsexo+education+employment+hhincome+marriage,id=aid , wave = wave, corstr ="exchangeable", data=addhealthWhite, weights=W5BIOWGT/mean(W5BIOWGT,na.rm=T))
summary(fitWh)
##
## Call:
## geeglm(formula = score ~ age + I(age^2) + sexW5 + raceXsexo +
## education + employment + hhincome + marriage, data = addhealthWhite,
## weights = W5BIOWGT/mean(W5BIOWGT, na.rm = T), id = aid, waves = wave,
## corstr = "exchangeable")
##
## Coefficients:
## Estimate Std.err Wald Pr(>|W|)
## (Intercept) 3.254979 0.599160 29.51 5.6e-08 ***
## age -0.031402 0.030696 1.05 0.306
## I(age^2) 0.000900 0.000648 1.93 0.165
## sexW5Female 0.107083 0.135351 0.63 0.429
## raceXsexoa_white.Gay -0.340044 0.175145 3.77 0.052 .
## raceXsexoa_white.Bisexual 0.225125 0.188360 1.43 0.232
## education1 -0.181065 0.270147 0.45 0.503
## education2 -0.537815 0.286648 3.52 0.061 .
## employmentempNoBene -0.359582 0.283543 1.61 0.205
## employmentempBene -0.414350 0.200245 4.28 0.039 *
## hhincomeb_30k<50k -0.150612 0.180683 0.69 0.405
## hhincomec_50k>100k 0.159975 0.145703 1.21 0.272
## hhincomed_>100k 0.070060 0.180286 0.15 0.698
## marriage1 -0.021531 0.127882 0.03 0.866
## marriage2 0.127726 0.259544 0.24 0.623
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation structure = exchangeable
## Estimated Scale Parameters:
##
## Estimate Std.err
## (Intercept) 3.62 0.32
## Link = identity
##
## Estimated Correlation Parameters:
## Estimate Std.err
## alpha 0.558 0.084
## Number of clusters: 1532 Maximum cluster size: 392
#education+employment+hhincome+marriage
#fit2<-geeglm(score~age+I(age^2)+sexW5+raceXsexo+livesParent2+mom.educ+mom.pubass+closeMom2+parentSupervise2+religion,id=aid , wave = wave, corstr ="exchangeable", data=addhealthWhite, weights=W5BIOWGT/mean(W5BIOWGT,na.rm=T))
#summary(fit2)
fitWh3<-geeglm(score~age+I(age^2)+sexW5+raceXsexo+education+employment+hhincome+marriage+livesParent2+mom.educ+mom.pubass+closeMom2+parentSupervise5+religion+schoolconflict+schoolClubJoin+scale(schoolSeg)+scale(pre),id=aid , wave = wave, corstr ="exchangeable", data=addhealthWhite, weights=W5BIOWGT/mean(W5BIOWGT,na.rm=T))
summary(fitWh3)
##
## Call:
## geeglm(formula = score ~ age + I(age^2) + sexW5 + raceXsexo +
## education + employment + hhincome + marriage + livesParent2 +
## mom.educ + mom.pubass + closeMom2 + parentSupervise5 + religion +
## schoolconflict + schoolClubJoin + scale(schoolSeg) + scale(pre),
## data = addhealthWhite, weights = W5BIOWGT/mean(W5BIOWGT,
## na.rm = T), id = aid, waves = wave, corstr = "exchangeable")
##
## Coefficients:
## Estimate Std.err Wald Pr(>|W|)
## (Intercept) 3.520402 0.625665 31.66 1.8e-08 ***
## age -0.031474 0.030882 1.04 0.3081
## I(age^2) 0.000884 0.000654 1.83 0.1764
## sexW5Female 0.098991 0.133591 0.55 0.4587
## raceXsexoa_white.Gay -0.354489 0.179700 3.89 0.0485 *
## raceXsexoa_white.Bisexual 0.227230 0.185544 1.50 0.2207
## education1 -0.168051 0.271921 0.38 0.5366
## education2 -0.494809 0.288364 2.94 0.0862 .
## employmentempNoBene -0.295717 0.278414 1.13 0.2882
## employmentempBene -0.389425 0.197561 3.89 0.0487 *
## hhincomeb_30k<50k -0.154855 0.182291 0.72 0.3956
## hhincomec_50k>100k 0.160387 0.148895 1.16 0.2814
## hhincomed_>100k 0.071771 0.175426 0.17 0.6824
## marriage1 -0.021781 0.128078 0.03 0.8650
## marriage2 0.121351 0.250760 0.23 0.6284
## livesParent2lives_1_parent 0.137031 0.196634 0.49 0.4859
## livesParent2lives_other 0.122458 0.387336 0.10 0.7519
## mom.educ1 -0.435967 0.195605 4.97 0.0258 *
## mom.educ2 -0.574066 0.210143 7.46 0.0063 **
## mom.pubass1 0.551120 0.322643 2.92 0.0876 .
## closeMom2b_close -0.076188 0.196802 0.15 0.6987
## parentSupervise5b 0.202431 0.136497 2.20 0.1381
## parentSupervise5c 0.427318 0.178261 5.75 0.0165 *
## religion1 -0.000505 0.076848 0.00 0.9948
## schoolconflict3 -0.026227 0.151591 0.03 0.8626
## schoolconflict4 0.025872 0.173540 0.02 0.8815
## schoolClubJoin1 0.285860 0.208287 1.88 0.1699
## scale(schoolSeg) 0.118013 0.061045 3.74 0.0532 .
## scale(pre) -0.107946 0.040339 7.16 0.0075 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation structure = exchangeable
## Estimated Scale Parameters:
##
## Estimate Std.err
## (Intercept) 3.49 0.292
## Link = identity
##
## Estimated Correlation Parameters:
## Estimate Std.err
## alpha 0.549 0.112
## Number of clusters: 1532 Maximum cluster size: 392
#fit4<-geeglm(score~age+I(age^2)+sexW5+raceXsexo+livesParent2+mom.educ+mom.pubass+closeMom2+parentSupervise2+religion+schoolconflict+schoolClubJoin+schoolSeg+pre+education+employment+hhincome+marriage,id=aid , wave = wave, corstr ="exchangeable", data=addhealthWhite, weights=W5BIOWGT/mean(W5BIOWGT,na.rm=T))
#summary(fit4)
#library(huxtable)
allmodelsWh<-huxreg(fitWh,fitWh3,statistics = character(0))
allmodelsWh
| (1) | (2) |
| (Intercept) | 3.255 *** | 3.520 *** |
| (0.599) | (0.626) |
| age | -0.031 | -0.031 |
| (0.031) | (0.031) |
| I(age^2) | 0.001 | 0.001 |
| (0.001) | (0.001) |
| sexW5Female | 0.107 | 0.099 |
| (0.135) | (0.134) |
| raceXsexoa_white.Gay | -0.340 | -0.354 * |
| (0.175) | (0.180) |
| raceXsexoa_white.Bisexual | 0.225 | 0.227 |
| (0.188) | (0.186) |
| education1 | -0.181 | -0.168 |
| (0.270) | (0.272) |
| education2 | -0.538 | -0.495 |
| (0.287) | (0.288) |
| employmentempNoBene | -0.360 | -0.296 |
| (0.284) | (0.278) |
| employmentempBene | -0.414 * | -0.389 * |
| (0.200) | (0.198) |
| hhincomeb_30k<50k | -0.151 | -0.155 |
| (0.181) | (0.182) |
| hhincomec_50k>100k | 0.160 | 0.160 |
| (0.146) | (0.149) |
| hhincomed_>100k | 0.070 | 0.072 |
| (0.180) | (0.175) |
| marriage1 | -0.022 | -0.022 |
| (0.128) | (0.128) |
| marriage2 | 0.128 | 0.121 |
| (0.260) | (0.251) |
| livesParent2lives_1_parent | | 0.137 |
| | (0.197) |
| livesParent2lives_other | | 0.122 |
| | (0.387) |
| mom.educ1 | | -0.436 * |
| | (0.196) |
| mom.educ2 | | -0.574 ** |
| | (0.210) |
| mom.pubass1 | | 0.551 |
| | (0.323) |
| closeMom2b_close | | -0.076 |
| | (0.197) |
| parentSupervise5b | | 0.202 |
| | (0.136) |
| parentSupervise5c | | 0.427 * |
| | (0.178) |
| religion1 | | -0.001 |
| | (0.077) |
| schoolconflict3 | | -0.026 |
| | (0.152) |
| schoolconflict4 | | 0.026 |
| | (0.174) |
| schoolClubJoin1 | | 0.286 |
| | (0.208) |
| scale(schoolSeg) | | 0.118 |
| | (0.061) |
| scale(pre) | | -0.108 ** |
| | (0.040) |
| *** p < 0.001; ** p < 0.01; * p < 0.05. |
addhealthM<-addhealthS1long%>%
dplyr::select(wave,aid,W5BIOWGT,age,sexW5,education,hhincome,employment,marriage,livesParent2,parentSupervise,mom.educ2,mom.educ,mom.employed,mom.nativity,mom.pubass,pre,teacherTrouble,studentTrouble,schoolClubJoin,schoolSeg,sexualorientation,religion,closeMom,score,raceXsexo,raceW5min)
addhealthMa<-addhealthM%>%
subset(sexW5=='Male')
addhealthMale<-addhealthMa%>%
filter(complete.cases(.))
addhealthMale$livesParent2<-factor(addhealthMale$livesParent2,levels=c("lives_2_parents","lives_1_parent","lives_other"))
addhealthMale$closeMom<-as.numeric(addhealthMale$closeMom)
addhealthMale$closeMom2<-Recode(addhealthMale$closeMom,recodes="1:2='a_notClose';3:5='b_close'")
addhealthMale$parentSupervise<-as.numeric(addhealthMale$parentSupervise)
#table(addhealthMale$parentSupervise)
addhealthMale$parentSupervise2<-Recode(addhealthMale$parentSupervise,recodes="1:3='a_easy';4:6='b_moderate';7:8='c_strict'")
addhealthMale$parentSupervise5<-as.factor(Recode(addhealthMale$parentSupervise,recodes="1:6='a';7='b';8='c'"))
addhealthMale$teacherTrouble2<-as.numeric(addhealthMale$teacherTrouble)
addhealthMale$studentTrouble2<-as.numeric(addhealthMale$studentTrouble)
addhealthMale$schoolconflict<-addhealthMale$studentTrouble2+addhealthMale$teacherTrouble2
addhealthMale$schoolconflict<-as.factor(addhealthMale$schoolconflict)
addhealthMale$raceXsexo<-factor(addhealthMale$raceXsexo)
table(addhealthMale$raceXsexo)
##
## a_white.Straight b_black.Straight c_hispanic.Straight a_white.Gay
## 5846 2887 643 188
## b_black.Gay c_hispanic.Gay a_white.Bisexual
## 92 13 48
fitM<-geeglm(score~age+I(age^2)+raceXsexo+education+employment+hhincome+marriage,id=aid , wave = wave, corstr ="exchangeable", data=addhealthMale, weights=W5BIOWGT/mean(W5BIOWGT,na.rm=T))
summary(fitM)
##
## Call:
## geeglm(formula = score ~ age + I(age^2) + raceXsexo + education +
## employment + hhincome + marriage, data = addhealthMale, weights = W5BIOWGT/mean(W5BIOWGT,
## na.rm = T), id = aid, waves = wave, corstr = "exchangeable")
##
## Coefficients:
## Estimate Std.err Wald Pr(>|W|)
## (Intercept) 4.36895 0.97254 20.18 7e-06 ***
## age -0.15694 0.06020 6.80 0.0091 **
## I(age^2) 0.00344 0.00117 8.64 0.0033 **
## raceXsexob_black.Straight -0.01338 0.30724 0.00 0.9653
## raceXsexoc_hispanic.Straight 0.14892 0.41489 0.13 0.7196
## raceXsexoa_white.Gay -0.50385 0.35320 2.03 0.1537
## raceXsexob_black.Gay -0.09229 0.71761 0.02 0.8977
## raceXsexoc_hispanic.Gay 1.00825 0.93751 1.16 0.2822
## raceXsexoa_white.Bisexual 0.31402 0.48853 0.41 0.5204
## education1 0.19862 0.13961 2.02 0.1548
## education2 -0.14957 0.15733 0.90 0.3418
## employmentempNoBene -0.30049 0.48242 0.39 0.5334
## employmentempBene -0.13805 0.39970 0.12 0.7298
## hhincomeb_30k<50k -0.42720 0.23158 3.40 0.0651 .
## hhincomec_50k>100k -0.10488 0.24240 0.19 0.6653
## hhincomed_>100k -0.37208 0.25951 2.06 0.1516
## marriage1 -0.09105 0.18995 0.23 0.6317
## marriage2 0.11402 0.26340 0.19 0.6651
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation structure = exchangeable
## Estimated Scale Parameters:
##
## Estimate Std.err
## (Intercept) 3.2 0.303
## Link = identity
##
## Estimated Correlation Parameters:
## Estimate Std.err
## alpha 0.702 0.0993
## Number of clusters: 737 Maximum cluster size: 450
#education+employment+hhincome+marriage
#fitM2<-geeglm(score~age+I(age^2)+raceXsexo+livesParent2+mom.educ+mom.pubass+closeMom2+parentSupervise2+religion,id=aid , wave = wave, corstr ="exchangeable", data=addhealthMale, weights=W5BIOWGT/mean(W5BIOWGT,na.rm=T))
#summary(fitM2)
fitM3<-geeglm(score~age+I(age^2)+raceXsexo+education+employment+hhincome+marriage+livesParent2+mom.educ+mom.pubass+closeMom2+parentSupervise5+religion+schoolconflict+schoolClubJoin+scale(schoolSeg)+scale(pre),id=aid , wave = wave, corstr ="exchangeable", data=addhealthMale, weights=W5BIOWGT/mean(W5BIOWGT,na.rm=T))
summary(fitM3)
##
## Call:
## geeglm(formula = score ~ age + I(age^2) + raceXsexo + education +
## employment + hhincome + marriage + livesParent2 + mom.educ +
## mom.pubass + closeMom2 + parentSupervise5 + religion + schoolconflict +
## schoolClubJoin + scale(schoolSeg) + scale(pre), data = addhealthMale,
## weights = W5BIOWGT/mean(W5BIOWGT, na.rm = T), id = aid, waves = wave,
## corstr = "exchangeable")
##
## Coefficients:
## Estimate Std.err Wald Pr(>|W|)
## (Intercept) 5.08124 1.09924 21.37 3.8e-06 ***
## age -0.17154 0.06205 7.64 0.0057 **
## I(age^2) 0.00362 0.00117 9.59 0.0020 **
## raceXsexob_black.Straight -0.09864 0.31089 0.10 0.7510
## raceXsexoc_hispanic.Straight 0.09820 0.42940 0.05 0.8191
## raceXsexoa_white.Gay -0.53348 0.33001 2.61 0.1060
## raceXsexob_black.Gay -0.12039 0.65959 0.03 0.8552
## raceXsexoc_hispanic.Gay 0.81713 0.91298 0.80 0.3708
## raceXsexoa_white.Bisexual 0.22202 0.40438 0.30 0.5830
## education1 0.13642 0.16558 0.68 0.4100
## education2 -0.24925 0.19818 1.58 0.2085
## employmentempNoBene -0.35080 0.45140 0.60 0.4371
## employmentempBene -0.19020 0.35998 0.28 0.5972
## hhincomeb_30k<50k -0.32153 0.20730 2.41 0.1209
## hhincomec_50k>100k -0.09296 0.22516 0.17 0.6797
## hhincomed_>100k -0.37086 0.23256 2.54 0.1108
## marriage1 -0.03812 0.20484 0.03 0.8524
## marriage2 0.18105 0.27948 0.42 0.5171
## livesParent2lives_1_parent -0.07147 0.28426 0.06 0.8015
## livesParent2lives_other -0.29656 0.69107 0.18 0.6678
## mom.educ1 -0.60474 0.33442 3.27 0.0706 .
## mom.educ2 -0.59771 0.35453 2.84 0.0918 .
## mom.pubass1 -0.19662 0.33361 0.35 0.5556
## closeMom2b_close 0.05834 0.15141 0.15 0.7000
## parentSupervise5b 0.40053 0.21412 3.50 0.0614 .
## parentSupervise5c 0.57323 0.25224 5.16 0.0231 *
## religion1 -0.28687 0.19624 2.14 0.1438
## schoolconflict3 -0.04313 0.21454 0.04 0.8407
## schoolconflict4 -0.19254 0.25326 0.58 0.4471
## schoolClubJoin1 0.31650 0.31837 0.99 0.3202
## scale(schoolSeg) 0.18931 0.09286 4.16 0.0415 *
## scale(pre) 0.04190 0.09322 0.20 0.6531
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation structure = exchangeable
## Estimated Scale Parameters:
##
## Estimate Std.err
## (Intercept) 3.14 0.34
## Link = identity
##
## Estimated Correlation Parameters:
## Estimate Std.err
## alpha 0.743 0.172
## Number of clusters: 737 Maximum cluster size: 450
#fitM4<-geeglm(score~age+I(age^2)+raceXsexo+livesParent2+mom.educ+mom.pubass+closeMom2+parentSupervise2+religion+schoolconflict+schoolClubJoin+schoolSeg+pre+education+employment+hhincome+marriage,id=aid , wave = wave, corstr ="exchangeable", data=addhealthMale, weights=W5BIOWGT/mean(W5BIOWGT,na.rm=T))
#summary(fitM4)
#library(huxtable)
allmodelsM<-huxreg(fitM,fitM3,statistics = character(0))
allmodelsM
| (1) | (2) |
| (Intercept) | 4.369 *** | 5.081 *** |
| (0.973) | (1.099) |
| age | -0.157 ** | -0.172 ** |
| (0.060) | (0.062) |
| I(age^2) | 0.003 ** | 0.004 ** |
| (0.001) | (0.001) |
| raceXsexob_black.Straight | -0.013 | -0.099 |
| (0.307) | (0.311) |
| raceXsexoc_hispanic.Straight | 0.149 | 0.098 |
| (0.415) | (0.429) |
| raceXsexoa_white.Gay | -0.504 | -0.533 |
| (0.353) | (0.330) |
| raceXsexob_black.Gay | -0.092 | -0.120 |
| (0.718) | (0.660) |
| raceXsexoc_hispanic.Gay | 1.008 | 0.817 |
| (0.938) | (0.913) |
| raceXsexoa_white.Bisexual | 0.314 | 0.222 |
| (0.489) | (0.404) |
| education1 | 0.199 | 0.136 |
| (0.140) | (0.166) |
| education2 | -0.150 | -0.249 |
| (0.157) | (0.198) |
| employmentempNoBene | -0.300 | -0.351 |
| (0.482) | (0.451) |
| employmentempBene | -0.138 | -0.190 |
| (0.400) | (0.360) |
| hhincomeb_30k<50k | -0.427 | -0.322 |
| (0.232) | (0.207) |
| hhincomec_50k>100k | -0.105 | -0.093 |
| (0.242) | (0.225) |
| hhincomed_>100k | -0.372 | -0.371 |
| (0.260) | (0.233) |
| marriage1 | -0.091 | -0.038 |
| (0.190) | (0.205) |
| marriage2 | 0.114 | 0.181 |
| (0.263) | (0.279) |
| livesParent2lives_1_parent | | -0.071 |
| | (0.284) |
| livesParent2lives_other | | -0.297 |
| | (0.691) |
| mom.educ1 | | -0.605 |
| | (0.334) |
| mom.educ2 | | -0.598 |
| | (0.355) |
| mom.pubass1 | | -0.197 |
| | (0.334) |
| closeMom2b_close | | 0.058 |
| | (0.151) |
| parentSupervise5b | | 0.401 |
| | (0.214) |
| parentSupervise5c | | 0.573 * |
| | (0.252) |
| religion1 | | -0.287 |
| | (0.196) |
| schoolconflict3 | | -0.043 |
| | (0.215) |
| schoolconflict4 | | -0.193 |
| | (0.253) |
| schoolClubJoin1 | | 0.316 |
| | (0.318) |
| scale(schoolSeg) | | 0.189 * |
| | (0.093) |
| scale(pre) | | 0.042 |
| | (0.093) |
| *** p < 0.001; ** p < 0.01; * p < 0.05. |
addhealthFem<-addhealthS1long%>%
dplyr::select(wave,aid,W5BIOWGT,age,sexW5,education,hhincome,employment,marriage,livesParent2,parentSupervise,mom.educ2,mom.educ,mom.employed,mom.nativity,mom.pubass,pre,teacherTrouble,studentTrouble,schoolClubJoin,schoolSeg,sexualorientation,religion,closeMom,score,raceXsexo,raceW5min)
addhealthFema<-addhealthFem%>%
subset(sexW5=='Female')
addhealthFemale<-addhealthFema%>%
filter(complete.cases(.))
addhealthFemale$livesParent2<-factor(addhealthFemale$livesParent2,levels=c("lives_2_parents","lives_1_parent","lives_other"))
addhealthFemale$closeMom<-as.numeric(addhealthFemale$closeMom)
addhealthFemale$closeMom2<-Recode(addhealthFemale$closeMom,recodes="1:2='a_notClose';3:5='b_close'")
addhealthFemale$parentSupervise<-as.numeric(addhealthFemale$parentSupervise)
#table(addhealthFemale$parentSupervise)
addhealthFemale$parentSupervise2<-Recode(addhealthFemale$parentSupervise,recodes="1:3='a_easy';4:6='b_moderate';7:8='c_strict'")
addhealthFemale$parentSupervise5<-as.factor(Recode(addhealthFemale$parentSupervise,recodes="1:6='a';7='b';8='c'"))
addhealthFemale$teacherTrouble2<-as.numeric(addhealthFemale$teacherTrouble)
addhealthFemale$studentTrouble2<-as.numeric(addhealthFemale$studentTrouble)
addhealthFemale$schoolconflict<-addhealthFemale$studentTrouble2+addhealthFemale$teacherTrouble2
addhealthFemale$schoolconflict<-as.factor(addhealthFemale$schoolconflict)
addhealthFemale$raceXsexo<-factor(addhealthFemale$raceXsexo)
table(addhealthFemale$raceXsexo)
##
## a_white.Straight b_black.Straight c_hispanic.Straight a_white.Gay
## 9466 3282 1081 302
## b_black.Gay c_hispanic.Gay a_white.Bisexual b_black.Bisexual
## 120 30 561 70
## c_hispanic.Bisexual
## 28
fitF<-geeglm(score~age+I(age^2)+raceXsexo+education+employment+hhincome+marriage,id=aid , wave = wave, corstr ="exchangeable", data=addhealthFemale, weights=W5BIOWGT/mean(W5BIOWGT,na.rm=T))
summary(fitF)
##
## Call:
## geeglm(formula = score ~ age + I(age^2) + raceXsexo + education +
## employment + hhincome + marriage, data = addhealthFemale,
## weights = W5BIOWGT/mean(W5BIOWGT, na.rm = T), id = aid, waves = wave,
## corstr = "exchangeable")
##
## Coefficients:
## Estimate Std.err Wald Pr(>|W|)
## (Intercept) 3.764520 0.655961 32.94 9.5e-09 ***
## age -0.026754 0.032346 0.68 0.40817
## I(age^2) 0.000721 0.000679 1.13 0.28875
## raceXsexob_black.Straight 0.614303 0.161815 14.41 0.00015 ***
## raceXsexoc_hispanic.Straight 0.222312 0.261436 0.72 0.39513
## raceXsexoa_white.Gay -0.387073 0.284263 1.85 0.17330
## raceXsexob_black.Gay 0.628323 0.302702 4.31 0.03792 *
## raceXsexoc_hispanic.Gay 2.176097 0.450513 23.33 1.4e-06 ***
## raceXsexoa_white.Bisexual 0.270645 0.207083 1.71 0.19123
## raceXsexob_black.Bisexual 0.117985 0.386653 0.09 0.76026
## raceXsexoc_hispanic.Bisexual 0.161091 0.272065 0.35 0.55378
## education1 -0.640717 0.371165 2.98 0.08431 .
## education2 -0.891014 0.409477 4.73 0.02956 *
## employmentempNoBene -0.357644 0.291190 1.51 0.21937
## employmentempBene -0.618127 0.189547 10.63 0.00111 **
## hhincomeb_30k<50k -0.077400 0.173188 0.20 0.65494
## hhincomec_50k>100k 0.201465 0.141598 2.02 0.15480
## hhincomed_>100k 0.258890 0.193251 1.79 0.18036
## marriage1 0.102578 0.117153 0.77 0.38125
## marriage2 0.313272 0.302870 1.07 0.30097
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation structure = exchangeable
## Estimated Scale Parameters:
##
## Estimate Std.err
## (Intercept) 3.82 0.364
## Link = identity
##
## Estimated Correlation Parameters:
## Estimate Std.err
## alpha 0.658 0.0908
## Number of clusters: 1326 Maximum cluster size: 324
#education+employment+hhincome+marriage
#fitF2<-geeglm(score~age+I(age^2)+raceXsexo+livesParent2+mom.educ+mom.pubass+closeMom2+parentSupervise2+religion,id=aid , wave = wave, corstr ="exchangeable", data=addhealthFemale, weights=W5BIOWGT/mean(W5BIOWGT,na.rm=T))
#summary(fitF2)
fitF3<-geeglm(score~age+I(age^2)+raceXsexo+education+employment+hhincome+marriage+livesParent2+mom.educ+mom.pubass+closeMom2+parentSupervise5+religion+schoolconflict+schoolClubJoin+scale(schoolSeg)+scale(pre),id=aid , wave = wave, corstr ="exchangeable", data=addhealthFemale, weights=W5BIOWGT/mean(W5BIOWGT,na.rm=T))
summary(fitF3)
##
## Call:
## geeglm(formula = score ~ age + I(age^2) + raceXsexo + education +
## employment + hhincome + marriage + livesParent2 + mom.educ +
## mom.pubass + closeMom2 + parentSupervise5 + religion + schoolconflict +
## schoolClubJoin + scale(schoolSeg) + scale(pre), data = addhealthFemale,
## weights = W5BIOWGT/mean(W5BIOWGT, na.rm = T), id = aid, waves = wave,
## corstr = "exchangeable")
##
## Coefficients:
## Estimate Std.err Wald Pr(>|W|)
## (Intercept) 3.787340 0.690230 30.11 4.1e-08 ***
## age -0.015342 0.032607 0.22 0.63799
## I(age^2) 0.000469 0.000675 0.48 0.48676
## raceXsexob_black.Straight 0.420823 0.180691 5.42 0.01986 *
## raceXsexoc_hispanic.Straight -0.055680 0.315804 0.03 0.86005
## raceXsexoa_white.Gay -0.471602 0.319660 2.18 0.14013
## raceXsexob_black.Gay 0.731547 0.271859 7.24 0.00713 **
## raceXsexoc_hispanic.Gay 1.838609 0.845220 4.73 0.02961 *
## raceXsexoa_white.Bisexual 0.338356 0.217134 2.43 0.11917
## raceXsexob_black.Bisexual -0.531550 0.449047 1.40 0.23652
## raceXsexoc_hispanic.Bisexual -0.143264 0.402364 0.13 0.72180
## education1 -0.680048 0.383219 3.15 0.07597 .
## education2 -0.833118 0.426110 3.82 0.05056 .
## employmentempNoBene -0.201227 0.289495 0.48 0.48699
## employmentempBene -0.509880 0.197982 6.63 0.01001 *
## hhincomeb_30k<50k -0.075537 0.178720 0.18 0.67255
## hhincomec_50k>100k 0.240976 0.146423 2.71 0.09982 .
## hhincomed_>100k 0.334657 0.196250 2.91 0.08815 .
## marriage1 0.132927 0.132725 1.00 0.31657
## marriage2 0.350602 0.321328 1.19 0.27523
## livesParent2lives_1_parent 0.422354 0.170475 6.14 0.01323 *
## livesParent2lives_other -0.239029 0.377104 0.40 0.52618
## mom.educ1 -0.179599 0.194494 0.85 0.35579
## mom.educ2 -0.470781 0.219356 4.61 0.03186 *
## mom.pubass1 0.747363 0.265228 7.94 0.00484 **
## closeMom2b_close -0.284597 0.298189 0.91 0.33987
## parentSupervise5b 0.113404 0.146511 0.60 0.43892
## parentSupervise5c 0.139097 0.190235 0.53 0.46466
## religion1 0.048661 0.086365 0.32 0.57314
## schoolconflict3 0.198417 0.169994 1.36 0.24313
## schoolconflict4 0.100896 0.181824 0.31 0.57896
## schoolClubJoin1 0.273934 0.209711 1.71 0.19147
## scale(schoolSeg) -0.014625 0.073975 0.04 0.84327
## scale(pre) -0.151608 0.044301 11.71 0.00062 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation structure = exchangeable
## Estimated Scale Parameters:
##
## Estimate Std.err
## (Intercept) 3.94 0.391
## Link = identity
##
## Estimated Correlation Parameters:
## Estimate Std.err
## alpha 1.09 0.25
## Number of clusters: 1326 Maximum cluster size: 324
#fitF4<-geeglm(score~age+I(age^2)+raceXsexo+livesParent2+mom.educ+mom.pubass+closeMom2+parentSupervise2+religion+schoolconflict+schoolClubJoin+schoolSeg+pre+education+employment+hhincome+marriage,id=aid , wave = wave, corstr ="exchangeable", data=addhealthFemale, weights=W5BIOWGT/mean(W5BIOWGT,na.rm=T))
#summary(fitF4)
#library(huxtable)
allmodelsF<-huxreg(fitF,fitF3,statistics = character(0))
allmodelsF
| (1) | (2) |
| (Intercept) | 3.765 *** | 3.787 *** |
| (0.656) | (0.690) |
| age | -0.027 | -0.015 |
| (0.032) | (0.033) |
| I(age^2) | 0.001 | 0.000 |
| (0.001) | (0.001) |
| raceXsexob_black.Straight | 0.614 *** | 0.421 * |
| (0.162) | (0.181) |
| raceXsexoc_hispanic.Straight | 0.222 | -0.056 |
| (0.261) | (0.316) |
| raceXsexoa_white.Gay | -0.387 | -0.472 |
| (0.284) | (0.320) |
| raceXsexob_black.Gay | 0.628 * | 0.732 ** |
| (0.303) | (0.272) |
| raceXsexoc_hispanic.Gay | 2.176 *** | 1.839 * |
| (0.451) | (0.845) |
| raceXsexoa_white.Bisexual | 0.271 | 0.338 |
| (0.207) | (0.217) |
| raceXsexob_black.Bisexual | 0.118 | -0.532 |
| (0.387) | (0.449) |
| raceXsexoc_hispanic.Bisexual | 0.161 | -0.143 |
| (0.272) | (0.402) |
| education1 | -0.641 | -0.680 |
| (0.371) | (0.383) |
| education2 | -0.891 * | -0.833 |
| (0.409) | (0.426) |
| employmentempNoBene | -0.358 | -0.201 |
| (0.291) | (0.289) |
| employmentempBene | -0.618 ** | -0.510 * |
| (0.190) | (0.198) |
| hhincomeb_30k<50k | -0.077 | -0.076 |
| (0.173) | (0.179) |
| hhincomec_50k>100k | 0.201 | 0.241 |
| (0.142) | (0.146) |
| hhincomed_>100k | 0.259 | 0.335 |
| (0.193) | (0.196) |
| marriage1 | 0.103 | 0.133 |
| (0.117) | (0.133) |
| marriage2 | 0.313 | 0.351 |
| (0.303) | (0.321) |
| livesParent2lives_1_parent | | 0.422 * |
| | (0.170) |
| livesParent2lives_other | | -0.239 |
| | (0.377) |
| mom.educ1 | | -0.180 |
| | (0.194) |
| mom.educ2 | | -0.471 * |
| | (0.219) |
| mom.pubass1 | | 0.747 ** |
| | (0.265) |
| closeMom2b_close | | -0.285 |
| | (0.298) |
| parentSupervise5b | | 0.113 |
| | (0.147) |
| parentSupervise5c | | 0.139 |
| | (0.190) |
| religion1 | | 0.049 |
| | (0.086) |
| schoolconflict3 | | 0.198 |
| | (0.170) |
| schoolconflict4 | | 0.101 |
| | (0.182) |
| schoolClubJoin1 | | 0.274 |
| | (0.210) |
| scale(schoolSeg) | | -0.015 |
| | (0.074) |
| scale(pre) | | -0.152 *** |
| | (0.044) |
| *** p < 0.001; ** p < 0.01; * p < 0.05. |
addhealthBl<-addhealthS1long%>%
dplyr::select(wave,aid,W5BIOWGT,age,sexW5,education,hhincome,employment,marriage,livesParent2,parentSupervise,mom.educ2,mom.educ,mom.employed,mom.nativity,mom.pubass,pre,teacherTrouble,studentTrouble,schoolClubJoin,schoolSeg,sexualorientation,religion,closeMom,score,raceXsexo,raceW5min)
addhealthBla<-addhealthBl%>%
subset(raceW5min=='b_black')
addhealthBlack<-addhealthBla%>%
filter(complete.cases(.))
addhealthBlack$livesParent2<-factor(addhealthBlack$livesParent2,levels=c("lives_2_parents","lives_1_parent","lives_other"))
addhealthBlack$closeMom<-as.numeric(addhealthBlack$closeMom)
addhealthBlack$closeMom2<-Recode(addhealthBlack$closeMom,recodes="1:2='a_notClose';3:5='b_close'")
addhealthBlack$parentSupervise<-as.numeric(addhealthBlack$parentSupervise)
#table(addhealthBlack$parentSupervise)
addhealthBlack$parentSupervise2<-Recode(addhealthBlack$parentSupervise,recodes="1:3='a_easy';4:6='b_moderate';7:8='c_strict'")
addhealthBlack$parentSupervise5<-as.factor(Recode(addhealthBlack$parentSupervise,recodes="1:6='a';7='b';8='c'"))
addhealthBlack$teacherTrouble2<-as.numeric(addhealthBlack$teacherTrouble)
addhealthBlack$studentTrouble2<-as.numeric(addhealthBlack$studentTrouble)
addhealthBlack$schoolconflict<-addhealthBlack$studentTrouble2+addhealthBlack$teacherTrouble2
addhealthBlack$schoolconflict<-as.factor(addhealthBlack$schoolconflict)
addhealthBlack$raceXsexo<-factor(addhealthBlack$raceXsexo)
table(addhealthBlack$raceXsexo)
##
## b_black.Straight b_black.Gay b_black.Bisexual
## 6169 212 70
fitB<-geeglm(score~age+I(age^2)+sexW5+raceXsexo+education+employment+hhincome+marriage,id=aid , wave = wave, corstr ="exchangeable", data=addhealthBlack, weights=W5BIOWGT/mean(W5BIOWGT,na.rm=T))
summary(fitB)
##
## Call:
## geeglm(formula = score ~ age + I(age^2) + sexW5 + raceXsexo +
## education + employment + hhincome + marriage, data = addhealthBlack,
## weights = W5BIOWGT/mean(W5BIOWGT, na.rm = T), id = aid, waves = wave,
## corstr = "exchangeable")
##
## Coefficients:
## Estimate Std.err Wald Pr(>|W|)
## (Intercept) 6.48836 1.50728 18.53 1.7e-05 ***
## age -0.25977 0.08810 8.69 0.0032 **
## I(age^2) 0.00530 0.00164 10.50 0.0012 **
## sexW5Female 0.65808 0.31346 4.41 0.0358 *
## raceXsexob_black.Gay -0.33636 0.19240 3.06 0.0804 .
## raceXsexob_black.Bisexual -0.55858 0.47048 1.41 0.2351
## education1 -0.17873 0.10221 3.06 0.0804 .
## education2 -0.14973 0.16904 0.78 0.3757
## employmentempNoBene -0.72440 0.58904 1.51 0.2188
## employmentempBene -0.79163 0.46242 2.93 0.0869 .
## hhincomeb_30k<50k -0.44998 0.26088 2.98 0.0846 .
## hhincomec_50k>100k -0.44943 0.28273 2.53 0.1119
## hhincomed_>100k -0.44801 0.37622 1.42 0.2337
## marriage1 0.20712 0.23745 0.76 0.3831
## marriage2 0.40459 0.38784 1.09 0.2969
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation structure = exchangeable
## Estimated Scale Parameters:
##
## Estimate Std.err
## (Intercept) 3.56 0.637
## Link = identity
##
## Estimated Correlation Parameters:
## Estimate Std.err
## alpha 0.98 0.262
## Number of clusters: 367 Maximum cluster size: 450
#education+employment+hhincome+marriage
#fitB2<-geeglm(score~age+I(age^2)+sexW5+raceXsexo+livesParent2+mom.educ+mom.pubass+closeMom2+parentSupervise2+religion,id=aid , wave = wave, corstr ="exchangeable", data=addhealthBlack, weights=W5BIOWGT/mean(W5BIOWGT,na.rm=T))
#summary(fitB2)
fitB3<-geeglm(score~age+I(age^2)+sexW5+raceXsexo+education+employment+hhincome+marriage+livesParent2+mom.educ+mom.pubass+closeMom2+parentSupervise5+religion+schoolconflict+schoolClubJoin+scale(schoolSeg)+scale(pre),id=aid , wave = wave, corstr ="exchangeable", data=addhealthBlack, weights=W5BIOWGT/mean(W5BIOWGT,na.rm=T))
summary(fitB3)
##
## Call:
## geeglm(formula = score ~ age + I(age^2) + sexW5 + raceXsexo +
## education + employment + hhincome + marriage + livesParent2 +
## mom.educ + mom.pubass + closeMom2 + parentSupervise5 + religion +
## schoolconflict + schoolClubJoin + scale(schoolSeg) + scale(pre),
## data = addhealthBlack, weights = W5BIOWGT/mean(W5BIOWGT,
## na.rm = T), id = aid, waves = wave, corstr = "exchangeable")
##
## Coefficients:
## Estimate Std.err Wald Pr(>|W|)
## (Intercept) 7.12998 1.56602 20.73 5.3e-06 ***
## age -0.26384 0.08170 10.43 0.00124 **
## I(age^2) 0.00525 0.00147 12.71 0.00036 ***
## sexW5Female 0.55015 0.30201 3.32 0.06851 .
## raceXsexob_black.Gay -0.46160 0.21136 4.77 0.02896 *
## raceXsexob_black.Bisexual -0.57173 0.47352 1.46 0.22728
## education1 -0.39118 0.18115 4.66 0.03081 *
## education2 -0.36418 0.29761 1.50 0.22108
## employmentempNoBene -0.83843 0.62586 1.79 0.18036
## employmentempBene -0.61087 0.49310 1.53 0.21540
## hhincomeb_30k<50k -0.48345 0.23335 4.29 0.03829 *
## hhincomec_50k>100k -0.47813 0.26680 3.21 0.07312 .
## hhincomed_>100k -0.45092 0.34976 1.66 0.19732
## marriage1 0.22175 0.24134 0.84 0.35819
## marriage2 0.43789 0.35472 1.52 0.21703
## livesParent2lives_1_parent 0.36361 0.28610 1.62 0.20376
## livesParent2lives_other -0.83942 0.52699 2.54 0.11119
## mom.educ1 -0.46395 0.42682 1.18 0.27704
## mom.educ2 -0.27551 0.46798 0.35 0.55605
## mom.pubass1 0.27562 0.37902 0.53 0.46711
## closeMom2b_close -0.33629 0.22423 2.25 0.13367
## parentSupervise5b 0.14587 0.28806 0.26 0.61258
## parentSupervise5c 0.11811 0.34656 0.12 0.73324
## religion1 -0.27829 0.21839 1.62 0.20258
## schoolconflict3 0.39032 0.32743 1.42 0.23323
## schoolconflict4 0.13093 0.34278 0.15 0.70249
## schoolClubJoin1 0.95095 0.46283 4.22 0.03992 *
## scale(schoolSeg) -0.07700 0.14300 0.29 0.59025
## scale(pre) -0.03373 0.14740 0.05 0.81900
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation structure = exchangeable
## Estimated Scale Parameters:
##
## Estimate Std.err
## (Intercept) 3.65 0.712
## Link = identity
##
## Estimated Correlation Parameters:
## Estimate Std.err
## alpha 0.967 0.41
## Number of clusters: 367 Maximum cluster size: 450
#fitB4<-geeglm(score~age+I(age^2)+sexW5+raceXsexo+livesParent2+mom.educ+mom.pubass+closeMom2+parentSupervise2+religion+schoolconflict+schoolClubJoin+schoolSeg+pre+education+employment+hhincome+marriage,id=aid , wave = wave, corstr ="exchangeable", data=addhealthBlack, weights=W5BIOWGT/mean(W5BIOWGT,na.rm=T))
#summary(fitB4)
#library(huxtable)
allmodelsB<-huxreg(fitB,fitB3,statistics = character(0))
allmodelsB
| (1) | (2) |
| (Intercept) | 6.488 *** | 7.130 *** |
| (1.507) | (1.566) |
| age | -0.260 ** | -0.264 ** |
| (0.088) | (0.082) |
| I(age^2) | 0.005 ** | 0.005 *** |
| (0.002) | (0.001) |
| sexW5Female | 0.658 * | 0.550 |
| (0.313) | (0.302) |
| raceXsexob_black.Gay | -0.336 | -0.462 * |
| (0.192) | (0.211) |
| raceXsexob_black.Bisexual | -0.559 | -0.572 |
| (0.470) | (0.474) |
| education1 | -0.179 | -0.391 * |
| (0.102) | (0.181) |
| education2 | -0.150 | -0.364 |
| (0.169) | (0.298) |
| employmentempNoBene | -0.724 | -0.838 |
| (0.589) | (0.626) |
| employmentempBene | -0.792 | -0.611 |
| (0.462) | (0.493) |
| hhincomeb_30k<50k | -0.450 | -0.483 * |
| (0.261) | (0.233) |
| hhincomec_50k>100k | -0.449 | -0.478 |
| (0.283) | (0.267) |
| hhincomed_>100k | -0.448 | -0.451 |
| (0.376) | (0.350) |
| marriage1 | 0.207 | 0.222 |
| (0.237) | (0.241) |
| marriage2 | 0.405 | 0.438 |
| (0.388) | (0.355) |
| livesParent2lives_1_parent | | 0.364 |
| | (0.286) |
| livesParent2lives_other | | -0.839 |
| | (0.527) |
| mom.educ1 | | -0.464 |
| | (0.427) |
| mom.educ2 | | -0.276 |
| | (0.468) |
| mom.pubass1 | | 0.276 |
| | (0.379) |
| closeMom2b_close | | -0.336 |
| | (0.224) |
| parentSupervise5b | | 0.146 |
| | (0.288) |
| parentSupervise5c | | 0.118 |
| | (0.347) |
| religion1 | | -0.278 |
| | (0.218) |
| schoolconflict3 | | 0.390 |
| | (0.327) |
| schoolconflict4 | | 0.131 |
| | (0.343) |
| schoolClubJoin1 | | 0.951 * |
| | (0.463) |
| scale(schoolSeg) | | -0.077 |
| | (0.143) |
| scale(pre) | | -0.034 |
| | (0.147) |
| *** p < 0.001; ** p < 0.01; * p < 0.05. |
addhealthHi<-addhealthS1long%>%
dplyr::select(wave,aid,W5BIOWGT,age,sexW5,education,hhincome,employment,marriage,livesParent2,parentSupervise,mom.educ2,mom.educ,mom.employed,mom.nativity,mom.pubass,pre,teacherTrouble,studentTrouble,schoolClubJoin,schoolSeg,sexualorientation,religion,closeMom,score,raceXsexo,raceW5min)
addhealthHis<-addhealthHi%>%
subset(raceW5min=='c_hispanic')
addhealthHispanic<-addhealthHis%>%
filter(complete.cases(.))
addhealthHispanic$livesParent2<-factor(addhealthHispanic$livesParent2,levels=c("lives_2_parents","lives_1_parent","lives_other"))
addhealthHispanic$closeMom<-as.numeric(addhealthHispanic$closeMom)
addhealthHispanic$closeMom2<-Recode(addhealthHispanic$closeMom,recodes="1:2='a_notClose';3:5='b_close'")
addhealthHispanic$parentSupervise<-as.numeric(addhealthHispanic$parentSupervise)
#table(addhealthHispanic$parentSupervise)
addhealthHispanic$parentSupervise2<-Recode(addhealthHispanic$parentSupervise,recodes="1:3='a_easy';4:6='b_moderate';7:8='c_strict'")
addhealthHispanic$parentSupervise5<-as.factor(Recode(addhealthHispanic$parentSupervise,recodes="1:6='a';7='b';8='c'"))
addhealthHispanic$teacherTrouble2<-as.numeric(addhealthHispanic$teacherTrouble)
addhealthHispanic$studentTrouble2<-as.numeric(addhealthHispanic$studentTrouble)
addhealthHispanic$schoolconflict<-addhealthHispanic$studentTrouble2+addhealthHispanic$teacherTrouble2
addhealthHispanic$schoolconflict<-as.factor(addhealthHispanic$schoolconflict)
addhealthHispanic$raceXsexo<-as.factor(addhealthHispanic$raceXsexo)
addhealthHispanic$raceXsexo <- droplevels(addhealthHispanic$raceXsexo)
fitH<-geeglm(score~age+I(age^2)+sexW5+raceXsexo+education+employment+hhincome+marriage,id=aid , wave = wave, corstr ="exchangeable", data=addhealthHispanic, weights=W5BIOWGT/mean(W5BIOWGT,na.rm=T))
summary(fitH)
##
## Call:
## geeglm(formula = score ~ age + I(age^2) + sexW5 + raceXsexo +
## education + employment + hhincome + marriage, data = addhealthHispanic,
## weights = W5BIOWGT/mean(W5BIOWGT, na.rm = T), id = aid, waves = wave,
## corstr = "exchangeable")
##
## Coefficients:
## Estimate Std.err Wald Pr(>|W|)
## (Intercept) 3.153229 1.257950 6.28 0.012 *
## age 0.023047 0.045035 0.26 0.609
## I(age^2) -0.000471 0.000915 0.26 0.607
## sexW5Female 0.048591 0.514096 0.01 0.925
## raceXsexoc_hispanic.Gay 0.620930 0.556958 1.24 0.265
## raceXsexoc_hispanic.Bisexual -0.083040 0.190454 0.19 0.663
## education1 0.268039 0.179261 2.24 0.135
## education2 -0.025925 0.290457 0.01 0.929
## employmentempNoBene -0.840258 1.117429 0.57 0.452
## employmentempBene -0.754035 0.844692 0.80 0.372
## hhincomeb_30k<50k -0.188787 0.222819 0.72 0.397
## hhincomec_50k>100k -0.101185 0.137636 0.54 0.462
## hhincomed_>100k -0.318719 0.175774 3.29 0.070 .
## marriage1 0.165478 0.383819 0.19 0.666
## marriage2 0.204743 0.262861 0.61 0.436
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation structure = exchangeable
## Estimated Scale Parameters:
##
## Estimate Std.err
## (Intercept) 3.17 0.437
## Link = identity
##
## Estimated Correlation Parameters:
## Estimate Std.err
## alpha 0.826 0.234
## Number of clusters: 164 Maximum cluster size: 144
#education+employment+hhincome+marriage
#fitH2<-geeglm(score~age+I(age^2)+sexW5+raceXsexo+livesParent2+mom.educ+mom.pubass+closeMom2+parentSupervise2+religion,id=aid , wave = wave, corstr ="exchangeable", data=addhealthHispanic, weights=W5BIOWGT/mean(W5BIOWGT,na.rm=T))
#summary(fitH2)
fitH3<-geeglm(score~age+I(age^2)+sexW5+raceXsexo+education+employment+hhincome+marriage+livesParent2+mom.educ+mom.pubass+closeMom2+parentSupervise5+religion+schoolconflict+schoolClubJoin+scale(schoolSeg)+scale(pre),id=aid , wave = wave, corstr ="exchangeable", data=addhealthHispanic, weights=W5BIOWGT/mean(W5BIOWGT,na.rm=T))
summary(fitH3)
##
## Call:
## geeglm(formula = score ~ age + I(age^2) + sexW5 + raceXsexo +
## education + employment + hhincome + marriage + livesParent2 +
## mom.educ + mom.pubass + closeMom2 + parentSupervise5 + religion +
## schoolconflict + schoolClubJoin + scale(schoolSeg) + scale(pre),
## data = addhealthHispanic, weights = W5BIOWGT/mean(W5BIOWGT,
## na.rm = T), id = aid, waves = wave, corstr = "exchangeable")
##
## Coefficients:
## Estimate Std.err Wald Pr(>|W|)
## (Intercept) 3.429169 1.073394 10.21 0.00140 **
## age 0.020969 0.046117 0.21 0.64934
## I(age^2) -0.000449 0.000948 0.22 0.63548
## sexW5Female -0.569474 0.435996 1.71 0.19150
## raceXsexoc_hispanic.Gay 1.028967 0.872719 1.39 0.23838
## raceXsexoc_hispanic.Bisexual 0.167226 0.444432 0.14 0.70672
## education1 0.344015 0.191559 3.23 0.07252 .
## education2 -0.046019 0.292435 0.02 0.87496
## employmentempNoBene -0.722340 0.905050 0.64 0.42480
## employmentempBene -1.251005 0.694028 3.25 0.07146 .
## hhincomeb_30k<50k -0.117370 0.209638 0.31 0.57557
## hhincomec_50k>100k -0.141256 0.147347 0.92 0.33773
## hhincomed_>100k -0.374078 0.169799 4.85 0.02759 *
## marriage1 0.194149 0.372788 0.27 0.60250
## marriage2 0.214566 0.256627 0.70 0.40310
## livesParent2lives_1_parent -0.341037 0.480267 0.50 0.47764
## livesParent2lives_other -2.301388 0.657622 12.25 0.00047 ***
## mom.educ1 0.470467 0.541161 0.76 0.38465
## mom.educ2 -0.720301 0.586555 1.51 0.21944
## mom.pubass1 -0.830110 0.593385 1.96 0.16183
## closeMom2b_close 0.096473 0.217064 0.20 0.65672
## parentSupervise5b 0.950843 0.507405 3.51 0.06094 .
## parentSupervise5c 0.091437 0.622866 0.02 0.88329
## religion1 -0.055191 0.130975 0.18 0.67347
## schoolconflict3 1.140959 0.487440 5.48 0.01925 *
## schoolconflict4 -0.417558 0.608157 0.47 0.49234
## schoolClubJoin1 -0.197599 0.491380 0.16 0.68759
## scale(schoolSeg) -0.095323 0.259892 0.13 0.71378
## scale(pre) 0.223660 0.139838 2.56 0.10973
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation structure = exchangeable
## Estimated Scale Parameters:
##
## Estimate Std.err
## (Intercept) 2.6 0.569
## Link = identity
##
## Estimated Correlation Parameters:
## Estimate Std.err
## alpha 0.533 0.396
## Number of clusters: 164 Maximum cluster size: 144
#fitH4<-geeglm(score~age+I(age^2)+sexW5+raceXsexo+livesParent2+mom.educ+mom.pubass+closeMom2+parentSupervise2+religion+schoolconflict+schoolClubJoi#n+schoolSeg+pre+education+employment+hhincome+marriage,id=aid , wave = wave, corstr ="exchangeable", data=addhealthHispanic, #weights=W5BIOWGT/mean(W5BIOWGT,na.rm=T))
#summary(fitH4)
#library(huxtable)
allmodelsH<-huxreg(fitH,fitH3,statistics = character(0))
allmodelsH
| (1) | (2) |
| (Intercept) | 3.153 * | 3.429 ** |
| (1.258) | (1.073) |
| age | 0.023 | 0.021 |
| (0.045) | (0.046) |
| I(age^2) | -0.000 | -0.000 |
| (0.001) | (0.001) |
| sexW5Female | 0.049 | -0.569 |
| (0.514) | (0.436) |
| raceXsexoc_hispanic.Gay | 0.621 | 1.029 |
| (0.557) | (0.873) |
| raceXsexoc_hispanic.Bisexual | -0.083 | 0.167 |
| (0.190) | (0.444) |
| education1 | 0.268 | 0.344 |
| (0.179) | (0.192) |
| education2 | -0.026 | -0.046 |
| (0.290) | (0.292) |
| employmentempNoBene | -0.840 | -0.722 |
| (1.117) | (0.905) |
| employmentempBene | -0.754 | -1.251 |
| (0.845) | (0.694) |
| hhincomeb_30k<50k | -0.189 | -0.117 |
| (0.223) | (0.210) |
| hhincomec_50k>100k | -0.101 | -0.141 |
| (0.138) | (0.147) |
| hhincomed_>100k | -0.319 | -0.374 * |
| (0.176) | (0.170) |
| marriage1 | 0.165 | 0.194 |
| (0.384) | (0.373) |
| marriage2 | 0.205 | 0.215 |
| (0.263) | (0.257) |
| livesParent2lives_1_parent | | -0.341 |
| | (0.480) |
| livesParent2lives_other | | -2.301 *** |
| | (0.658) |
| mom.educ1 | | 0.470 |
| | (0.541) |
| mom.educ2 | | -0.720 |
| | (0.587) |
| mom.pubass1 | | -0.830 |
| | (0.593) |
| closeMom2b_close | | 0.096 |
| | (0.217) |
| parentSupervise5b | | 0.951 |
| | (0.507) |
| parentSupervise5c | | 0.091 |
| | (0.623) |
| religion1 | | -0.055 |
| | (0.131) |
| schoolconflict3 | | 1.141 * |
| | (0.487) |
| schoolconflict4 | | -0.418 |
| | (0.608) |
| schoolClubJoin1 | | -0.198 |
| | (0.491) |
| scale(schoolSeg) | | -0.095 |
| | (0.260) |
| scale(pre) | | 0.224 |
| | (0.140) |
| *** p < 0.001; ** p < 0.01; * p < 0.05. |
addhealthL<-addhealthS1long%>%
dplyr::select(wave,aid,W5BIOWGT,age,sexW5,education,hhincome,employment,marriage,livesParent2,parentSupervise,mom.educ2,mom.educ,mom.employed,mom.nativity,mom.pubass,pre,teacherTrouble,studentTrouble,schoolClubJoin,schoolSeg,sexualorientation,religion,closeMom,score,raceXsexo,raceW5min)
#library(tidyr)
addhealthLG<-addhealthL%>%
subset(sexualorientation=='Gay'|sexualorientation=="Bisexual")
addhealthLG$sexualorientation<-factor(addhealthLG$sexualorientation)
addhealthLG$raceXsexo<-interaction(addhealthLG$raceW5min,addhealthLG$sexualorientation)
addhealthLG$raceXsexo<-factor(addhealthLG$raceXsexo)
addhealthLGBTQ<-addhealthLG%>%
filter(complete.cases(.))
addhealthLGBTQ$livesParent2<-factor(addhealthLGBTQ$livesParent2,levels=c("lives_2_parents","lives_1_parent","lives_other"))
addhealthLGBTQ$parentSupervise5<-as.factor(Recode(addhealthLGBTQ$parentSupervise,recodes="1:6='a';7='b';8='c'"))
addhealthLGBTQ$closeMom<-as.numeric(addhealthLGBTQ$closeMom)
addhealthLGBTQ$closeMom2<-Recode(addhealthLGBTQ$closeMom,recodes="1:2='a_notClose';3:5='b_close'")
addhealthLGBTQ$parentSupervise<-as.numeric(addhealthLGBTQ$parentSupervise)
#table(addhealthLGBTQ$parentSupervise)
addhealthLGBTQ$parentSupervise2<-Recode(addhealthLGBTQ$parentSupervise,recodes="1:3='a_easy';4:6='b_moderate';7:8='c_strict'")
addhealthLGBTQ$teacherTrouble2<-as.numeric(addhealthLGBTQ$teacherTrouble)
addhealthLGBTQ$studentTrouble2<-as.numeric(addhealthLGBTQ$studentTrouble)
addhealthLGBTQ$schoolconflict<-addhealthLGBTQ$studentTrouble2+addhealthLGBTQ$teacherTrouble2
addhealthLGBTQ$schoolconflict<-as.factor(addhealthLGBTQ$schoolconflict)
#install.packages("geepack")
#library(geepack)
#wave,aid,W5BIOWGT,age,sexW5,education,hhincome,employment,marriage,,,,,,,,pre,teacherTrouble,studentTrouble,schoolClubJoin,schoolSeg,
fitL<-geeglm(score~age+I(age^2)+sexW5+raceXsexo+education+employment+hhincome+marriage,id=aid , wave = wave, corstr ="exchangeable", data=addhealthLGBTQ, weights=W5BIOWGT/mean(W5BIOWGT,na.rm=T))
summary(fitL)
##
## Call:
## geeglm(formula = score ~ age + I(age^2) + sexW5 + raceXsexo +
## education + employment + hhincome + marriage, data = addhealthLGBTQ,
## weights = W5BIOWGT/mean(W5BIOWGT, na.rm = T), id = aid, waves = wave,
## corstr = "exchangeable")
##
## Coefficients:
## Estimate Std.err Wald Pr(>|W|)
## (Intercept) -0.29051 2.34022 0.02 0.9012
## age 0.08179 0.14578 0.31 0.5748
## I(age^2) -0.00143 0.00256 0.31 0.5773
## sexW5Female -0.06095 0.66340 0.01 0.9268
## raceXsexob_black.Gay 1.57149 0.90320 3.03 0.0819 .
## raceXsexoc_hispanic.Gay 1.16287 1.71219 0.46 0.4970
## raceXsexoa_white.Bisexual 0.54619 0.29989 3.32 0.0686 .
## raceXsexob_black.Bisexual 3.24041 1.22784 6.96 0.0083 **
## raceXsexoc_hispanic.Bisexual 1.57958 1.36590 1.34 0.2475
## education1 -0.38501 1.03418 0.14 0.7097
## education2 -0.97528 1.12450 0.75 0.3858
## employmentempNoBene 0.95903 1.03201 0.86 0.3527
## employmentempBene 1.83500 0.74471 6.07 0.0137 *
## hhincomeb_30k<50k -0.25509 0.34698 0.54 0.4622
## hhincomec_50k>100k 0.65592 0.25141 6.81 0.0091 **
## hhincomed_>100k 0.08902 0.37506 0.06 0.8124
## marriage1 0.47095 0.54630 0.74 0.3886
## marriage2 -0.40852 0.41658 0.96 0.3268
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation structure = exchangeable
## Estimated Scale Parameters:
##
## Estimate Std.err
## (Intercept) 3.4 0.757
## Link = identity
##
## Estimated Correlation Parameters:
## Estimate Std.err
## alpha 0.491 0.192
## Number of clusters: 107 Maximum cluster size: 108
#education+employment+hhincome+marriage
#fitL2<-geeglm(score~age+I(age^2)+sexW5+raceXsexo+livesParent2+mom.educ+mom.pubass+closeMom2+parentSupervise2+religion,id=aid , wave = wave, corstr ="exchangeable", data=addhealthLGBTQ, weights=W5BIOWGT/mean(W5BIOWGT,na.rm=T))
#summary(fitL2)
fitL3<-geeglm(score~age+I(age^2)+sexW5+raceXsexo+education+employment+hhincome+marriage+livesParent2+mom.educ+mom.pubass+closeMom2+parentSupervise5+religion+schoolconflict+schoolClubJoin+scale(schoolSeg)+scale(pre),id=aid , wave = wave, corstr ="exchangeable", data=addhealthLGBTQ, weights=W5BIOWGT/mean(W5BIOWGT,na.rm=T))
summary(fitL3)
##
## Call:
## geeglm(formula = score ~ age + I(age^2) + sexW5 + raceXsexo +
## education + employment + hhincome + marriage + livesParent2 +
## mom.educ + mom.pubass + closeMom2 + parentSupervise5 + religion +
## schoolconflict + schoolClubJoin + scale(schoolSeg) + scale(pre),
## data = addhealthLGBTQ, weights = W5BIOWGT/mean(W5BIOWGT,
## na.rm = T), id = aid, waves = wave, corstr = "exchangeable")
##
## Coefficients:
## Estimate Std.err Wald Pr(>|W|)
## (Intercept) -0.53155 2.08816 0.06 0.79907
## age 0.06533 0.13660 0.23 0.63245
## I(age^2) -0.00107 0.00249 0.19 0.66709
## sexW5Female 0.20711 0.56299 0.14 0.71297
## raceXsexob_black.Gay 0.88428 0.93691 0.89 0.34526
## raceXsexoc_hispanic.Gay 0.45256 1.30760 0.12 0.72927
## raceXsexoa_white.Bisexual 0.11702 0.39531 0.09 0.76721
## raceXsexob_black.Bisexual 3.31633 1.12497 8.69 0.00320 **
## raceXsexoc_hispanic.Bisexual 1.12533 1.17185 0.92 0.33690
## education1 -0.51324 1.23090 0.17 0.67671
## education2 -0.97128 1.31229 0.55 0.45921
## employmentempNoBene -0.36317 0.95431 0.14 0.70353
## employmentempBene 1.51128 0.82353 3.37 0.06649 .
## hhincomeb_30k<50k -0.06739 0.45293 0.02 0.88171
## hhincomec_50k>100k 0.74160 0.20840 12.66 0.00037 ***
## hhincomed_>100k -0.22845 0.32106 0.51 0.47673
## marriage1 0.46829 0.39073 1.44 0.23072
## marriage2 -0.27348 0.39764 0.47 0.49160
## livesParent2lives_1_parent 0.08858 0.60717 0.02 0.88400
## livesParent2lives_other 2.06178 0.76572 7.25 0.00709 **
## mom.educ1 1.15207 0.64725 3.17 0.07508 .
## mom.educ2 0.90072 0.71703 1.58 0.20905
## mom.pubass1 0.46264 1.00910 0.21 0.64662
## closeMom2b_close -0.35612 0.37611 0.90 0.34372
## parentSupervise5b 1.73244 0.67039 6.68 0.00976 **
## religion1 1.07006 0.20416 27.47 1.6e-07 ***
## schoolconflict3 -0.38854 0.61513 0.40 0.52762
## schoolconflict4 0.47232 0.69053 0.47 0.49397
## schoolClubJoin1 -0.60010 0.83882 0.51 0.47436
## scale(schoolSeg) -0.05870 0.32895 0.03 0.85838
## scale(pre) -0.82332 0.25195 10.68 0.00108 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation structure = exchangeable
## Estimated Scale Parameters:
##
## Estimate Std.err
## (Intercept) 1.82 0.505
## Link = identity
##
## Estimated Correlation Parameters:
## Estimate Std.err
## alpha 0.256 0.128
## Number of clusters: 107 Maximum cluster size: 108
#fitL4<-geeglm(score~age+I(age^2)+sexW5+raceXsexo+livesParent2+mom.educ+mom.pubass+closeMom2+parentSupervise2+religion+schoolconflict+schoolClubJoin+schoolSeg+pre+education+employment+hhincome+marriage,id=aid , wave = wave, corstr ="exchangeable", data=addhealthLGBTQ, weights=W5BIOWGT/mean(W5BIOWGT,na.rm=T))
#summary(fitL4)
#library(huxtable)
allmodelsL<-huxreg(fitL,fitL3,statistics = character(0))
allmodelsL
| (1) | (2) |
| (Intercept) | -0.291 | -0.532 |
| (2.340) | (2.088) |
| age | 0.082 | 0.065 |
| (0.146) | (0.137) |
| I(age^2) | -0.001 | -0.001 |
| (0.003) | (0.002) |
| sexW5Female | -0.061 | 0.207 |
| (0.663) | (0.563) |
| raceXsexob_black.Gay | 1.571 | 0.884 |
| (0.903) | (0.937) |
| raceXsexoc_hispanic.Gay | 1.163 | 0.453 |
| (1.712) | (1.308) |
| raceXsexoa_white.Bisexual | 0.546 | 0.117 |
| (0.300) | (0.395) |
| raceXsexob_black.Bisexual | 3.240 ** | 3.316 ** |
| (1.228) | (1.125) |
| raceXsexoc_hispanic.Bisexual | 1.580 | 1.125 |
| (1.366) | (1.172) |
| education1 | -0.385 | -0.513 |
| (1.034) | (1.231) |
| education2 | -0.975 | -0.971 |
| (1.125) | (1.312) |
| employmentempNoBene | 0.959 | -0.363 |
| (1.032) | (0.954) |
| employmentempBene | 1.835 * | 1.511 |
| (0.745) | (0.824) |
| hhincomeb_30k<50k | -0.255 | -0.067 |
| (0.347) | (0.453) |
| hhincomec_50k>100k | 0.656 ** | 0.742 *** |
| (0.251) | (0.208) |
| hhincomed_>100k | 0.089 | -0.228 |
| (0.375) | (0.321) |
| marriage1 | 0.471 | 0.468 |
| (0.546) | (0.391) |
| marriage2 | -0.409 | -0.273 |
| (0.417) | (0.398) |
| livesParent2lives_1_parent | | 0.089 |
| | (0.607) |
| livesParent2lives_other | | 2.062 ** |
| | (0.766) |
| mom.educ1 | | 1.152 |
| | (0.647) |
| mom.educ2 | | 0.901 |
| | (0.717) |
| mom.pubass1 | | 0.463 |
| | (1.009) |
| closeMom2b_close | | -0.356 |
| | (0.376) |
| parentSupervise5b | | 1.732 ** |
| | (0.670) |
| religion1 | | 1.070 *** |
| | (0.204) |
| schoolconflict3 | | -0.389 |
| | (0.615) |
| schoolconflict4 | | 0.472 |
| | (0.691) |
| schoolClubJoin1 | | -0.600 |
| | (0.839) |
| scale(schoolSeg) | | -0.059 |
| | (0.329) |
| scale(pre) | | -0.823 ** |
| | (0.252) |
| *** p < 0.001; ** p < 0.01; * p < 0.05. |
summary(fitL3)
##
## Call:
## geeglm(formula = score ~ age + I(age^2) + sexW5 + raceXsexo +
## education + employment + hhincome + marriage + livesParent2 +
## mom.educ + mom.pubass + closeMom2 + parentSupervise5 + religion +
## schoolconflict + schoolClubJoin + scale(schoolSeg) + scale(pre),
## data = addhealthLGBTQ, weights = W5BIOWGT/mean(W5BIOWGT,
## na.rm = T), id = aid, waves = wave, corstr = "exchangeable")
##
## Coefficients:
## Estimate Std.err Wald Pr(>|W|)
## (Intercept) -0.53155 2.08816 0.06 0.79907
## age 0.06533 0.13660 0.23 0.63245
## I(age^2) -0.00107 0.00249 0.19 0.66709
## sexW5Female 0.20711 0.56299 0.14 0.71297
## raceXsexob_black.Gay 0.88428 0.93691 0.89 0.34526
## raceXsexoc_hispanic.Gay 0.45256 1.30760 0.12 0.72927
## raceXsexoa_white.Bisexual 0.11702 0.39531 0.09 0.76721
## raceXsexob_black.Bisexual 3.31633 1.12497 8.69 0.00320 **
## raceXsexoc_hispanic.Bisexual 1.12533 1.17185 0.92 0.33690
## education1 -0.51324 1.23090 0.17 0.67671
## education2 -0.97128 1.31229 0.55 0.45921
## employmentempNoBene -0.36317 0.95431 0.14 0.70353
## employmentempBene 1.51128 0.82353 3.37 0.06649 .
## hhincomeb_30k<50k -0.06739 0.45293 0.02 0.88171
## hhincomec_50k>100k 0.74160 0.20840 12.66 0.00037 ***
## hhincomed_>100k -0.22845 0.32106 0.51 0.47673
## marriage1 0.46829 0.39073 1.44 0.23072
## marriage2 -0.27348 0.39764 0.47 0.49160
## livesParent2lives_1_parent 0.08858 0.60717 0.02 0.88400
## livesParent2lives_other 2.06178 0.76572 7.25 0.00709 **
## mom.educ1 1.15207 0.64725 3.17 0.07508 .
## mom.educ2 0.90072 0.71703 1.58 0.20905
## mom.pubass1 0.46264 1.00910 0.21 0.64662
## closeMom2b_close -0.35612 0.37611 0.90 0.34372
## parentSupervise5b 1.73244 0.67039 6.68 0.00976 **
## religion1 1.07006 0.20416 27.47 1.6e-07 ***
## schoolconflict3 -0.38854 0.61513 0.40 0.52762
## schoolconflict4 0.47232 0.69053 0.47 0.49397
## schoolClubJoin1 -0.60010 0.83882 0.51 0.47436
## scale(schoolSeg) -0.05870 0.32895 0.03 0.85838
## scale(pre) -0.82332 0.25195 10.68 0.00108 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation structure = exchangeable
## Estimated Scale Parameters:
##
## Estimate Std.err
## (Intercept) 1.82 0.505
## Link = identity
##
## Estimated Correlation Parameters:
## Estimate Std.err
## alpha 0.256 0.128
## Number of clusters: 107 Maximum cluster size: 108