##thoughts of suicide wave 1
mergeddata4$sui_1<-Recode(mergeddata4$H1SU1, recodes="1=1; 0=0; else=NA")
## age, if friend attempted suicide, if family attempted, if live alone, education, thoughts of suicide, sex, race, at wave 3
mergeddata4$age_3<-(mergeddata4$IYEAR3 - mergeddata4$H3OD1Y)
mergeddata4$fri_sui3<-Recode(mergeddata4$H3TO133, recodes="1=1; 0=0; else=NA")
mergeddata4$fam_sui3 <-Recode(mergeddata4$H3TO135, recodes="1=1; 0=0; else=NA")
mergeddata4$live_03<-Recode(mergeddata4$H3HR5, recodes="1=1; 2=0; else=NA")
mergeddata4$H3ED1 <- as.numeric(mergeddata4$H3ED1)
mergeddata4$educ_re3<-Recode(mergeddata4$H3ED1, recodes= "6:11='lesshigh';
12='high'; 13:15='somecol'; 16:22='col'; else=NA", as.factor=T)
mergeddata4$sui_3<-Recode(mergeddata4$H3TO130, recodes="1=1; 0=0; else=NA")
mergeddata4$male_3<-Recode(mergeddata4$BIO_SEX3, recodes="1=1; 2=0; else=NA")
mergeddata4$H3IR4 <- as.numeric(mergeddata4$H3IR4)
mergeddata4$race_3<-Recode(mergeddata4$H3IR4, recodes= "1='nwhite';
2='nhblack'; 3:5='other'; else=NA", as.factor=T)
## Adverse childhood experiences grade 6
## how often left alone as a kid
mergeddata4$H3MA1 <- as.numeric(mergeddata4$H3MA1)
mergeddata4$gr6la<-Recode(mergeddata4$H3MA1, recodes= "1:2='oneortwotimes';
3='thrtofivetimes'; 4:5='sixtimesormore'; 6='never'; else=NA", as.factor=T)
## needs not met grade 6
mergeddata4$H3MA2 <- as.numeric(mergeddata4$H3MA2)
mergeddata4$gr6nnm<-Recode(mergeddata4$H3MA2, recodes= "1:2='oneortwotimes';
3='thrtofivetimes'; 4:5='sixtimesormore'; 6='never'; else=NA", as.factor=T)
## times hit grade 6
mergeddata4$H3MA3 <- as.numeric(mergeddata4$H3MA3)
mergeddata4$gr6th<-Recode(mergeddata4$H3MA3, recodes= "1:2='oneortwotimes';
3='thrtofivetimes'; 4:5='sixtimesormore'; 6='never'; else=NA", as.factor=T)
## touched sexually
mergeddata4$H3MA4 <- as.numeric(mergeddata4$H3MA4)
mergeddata4$gr6ts<-Recode(mergeddata4$H3MA4, recodes= "1:2='oneortwotimes';
3='thrtofivetimes'; 4:5='sixtimesormore'; 6='never'; else=NA", as.factor=T)
## thoughts of suicide at wave four
mergeddata4$sui_4<-Recode(mergeddata4$H4SE1, recodes="1=1; 0=0; else=NA")
## use scale ( ) function for continuous, makes a z score
## in order to make an interaction term you do the following (variable name * variable name 2) within the regression model, do not attempt to make a unique variable
## no strata value, does not effect standard errors that much
##Selecting my variables
myvars<-c( "AID", "age_3", "male_3", "fri_sui3", "educ_re3", "race_3", "fam_sui3", "sui_1", "sui_3", "sui_4",
"gr6th", "gr6ts", "GSWGT134", "CLUSTER2")
mergeddata4<-mergeddata4[,myvars]
##subsetting data
sam <- mergeddata4 %>%
filter(complete.cases(.))
##suicide transition variable
sam<-sam %>%
mutate(suitran_1 =ifelse(sui_1==0 & sui_3 ==0, 0,1),
sui_tran2 =ifelse(suitran_1==1, NA,
ifelse(sui_3==0 & sui_4==0,0,1)))
##suicide thought transition based on self reported education at wave 3
##You must form a person-period data set
subpp<-survSplit(Surv(age_3, sui_3)~.,
data = sam, cut=seq(19,28,1),
episode ="age")
## Survey design
options(survey.lonely.psu = "adjust")
des<-svydesign(ids =~CLUSTER2, weights = ~GSWGT134,
data = subpp, nest = T)
##Fit the discrete time hazard model to your outcome
##Consider both the general model and other time specifications
##Include all main effects in the model
##Test for an interaction between at least two of the predictors
## Basic discrete time model
m0<-svyglm(sui_3~age_3+male_3+as.factor(educ_re3)+fri_sui3+fam_sui3+as.factor(gr6th)+as.factor(gr6ts),
design=des, family=binomial (link="cloglog"))
## Warning in eval(family$initialize): non-integer #successes in a binomial glm!
summary(m0)
##
## Call:
## svyglm(formula = sui_3 ~ age_3 + male_3 + as.factor(educ_re3) +
## fri_sui3 + fam_sui3 + as.factor(gr6th) + as.factor(gr6ts),
## design = des, family = binomial(link = "cloglog"))
##
## Survey design:
## svydesign(ids = ~CLUSTER2, weights = ~GSWGT134, data = subpp,
## nest = T)
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -10.1372 0.9257 -10.951 < 2e-16 ***
## age_3 0.2361 0.0411 5.745 7.31e-08 ***
## male_3 -0.1934 0.1503 -1.287 0.20062
## as.factor(educ_re3)high 0.7775 0.2782 2.794 0.00607 **
## as.factor(educ_re3)lesshigh 0.6725 0.2862 2.350 0.02043 *
## as.factor(educ_re3)somecol 0.7678 0.2536 3.028 0.00302 **
## fri_sui3 1.0682 0.2109 5.064 1.53e-06 ***
## fam_sui3 0.8810 0.2533 3.477 0.00071 ***
## as.factor(gr6th)oneortwotimes 0.3931 0.2419 1.625 0.10685
## as.factor(gr6th)sixtimesormore 0.6353 0.2198 2.891 0.00458 **
## as.factor(gr6th)thrtofivetimes 0.6381 0.2822 2.261 0.02560 *
## as.factor(gr6ts)oneortwotimes 0.4047 0.3120 1.297 0.19724
## as.factor(gr6ts)sixtimesormore 1.0682 0.3724 2.869 0.00488 **
## as.factor(gr6ts)thrtofivetimes -0.2938 0.6346 -0.463 0.64421
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 0.9663934)
##
## Number of Fisher Scoring iterations: 7
## With interaction terms
m1<-svyglm(suitran_1~age_3+male_3+as.factor(educ_re3)+male_3*as.factor(educ_re3)+male_3*fri_sui3+male_3*fam_sui3+
as.factor(gr6th)+as.factor(gr6ts),
design=des, family=binomial (link="cloglog"))
## Warning in eval(family$initialize): non-integer #successes in a binomial glm!
summary(m1)
##
## Call:
## svyglm(formula = suitran_1 ~ age_3 + male_3 + as.factor(educ_re3) +
## male_3 * as.factor(educ_re3) + male_3 * fri_sui3 + male_3 *
## fam_sui3 + as.factor(gr6th) + as.factor(gr6ts), design = des,
## family = binomial(link = "cloglog"))
##
## Survey design:
## svydesign(ids = ~CLUSTER2, weights = ~GSWGT134, data = subpp,
## nest = T)
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -2.208089 0.337114 -6.550 1.77e-09 ***
## age_3 -0.005725 0.013827 -0.414 0.679631
## male_3 0.318522 0.277852 1.146 0.254063
## as.factor(educ_re3)high 0.780871 0.162947 4.792 5.06e-06 ***
## as.factor(educ_re3)lesshigh 1.008308 0.208610 4.833 4.26e-06 ***
## as.factor(educ_re3)somecol 0.640042 0.192738 3.321 0.001209 **
## fri_sui3 0.165989 0.235567 0.705 0.482486
## fam_sui3 0.336336 0.223911 1.502 0.135861
## as.factor(gr6th)oneortwotimes 0.167244 0.148407 1.127 0.262162
## as.factor(gr6th)sixtimesormore 0.705600 0.154436 4.569 1.26e-05 ***
## as.factor(gr6th)thrtofivetimes 0.638014 0.169702 3.760 0.000271 ***
## as.factor(gr6ts)oneortwotimes 0.718940 0.344392 2.088 0.039085 *
## as.factor(gr6ts)sixtimesormore 0.660536 0.398266 1.659 0.099984 .
## as.factor(gr6ts)thrtofivetimes 0.314391 0.389518 0.807 0.421288
## male_3:as.factor(educ_re3)high -0.854779 0.302438 -2.826 0.005570 **
## male_3:as.factor(educ_re3)lesshigh -1.214610 0.403727 -3.008 0.003237 **
## male_3:as.factor(educ_re3)somecol -0.840957 0.327137 -2.571 0.011450 *
## male_3:fri_sui3 0.413105 0.328052 1.259 0.210529
## male_3:fam_sui3 0.641296 0.432985 1.481 0.141361
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 0.9950363)
##
## Number of Fisher Scoring iterations: 5
## linear model
m2<-svyglm(suitran_1~age_3+male_3+as.factor(educ_re3)+fri_sui3+fam_sui3+as.factor(gr6th)+as.factor(gr6ts),
design=des, family=binomial (link="cloglog"))
## Warning in eval(family$initialize): non-integer #successes in a binomial glm!
summary(m2)
##
## Call:
## svyglm(formula = suitran_1 ~ age_3 + male_3 + as.factor(educ_re3) +
## fri_sui3 + fam_sui3 + as.factor(gr6th) + as.factor(gr6ts),
## design = des, family = binomial(link = "cloglog"))
##
## Survey design:
## svydesign(ids = ~CLUSTER2, weights = ~GSWGT134, data = subpp,
## nest = T)
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -1.893780 0.308835 -6.132 1.19e-08 ***
## age_3 -0.005306 0.013520 -0.392 0.695425
## male_3 -0.372764 0.101839 -3.660 0.000378 ***
## as.factor(educ_re3)high 0.380923 0.162689 2.341 0.020888 *
## as.factor(educ_re3)lesshigh 0.437184 0.175501 2.491 0.014126 *
## as.factor(educ_re3)somecol 0.255153 0.157512 1.620 0.107924
## fri_sui3 0.362631 0.178476 2.032 0.044418 *
## fam_sui3 0.602398 0.183274 3.287 0.001335 **
## as.factor(gr6th)oneortwotimes 0.176657 0.148643 1.188 0.237037
## as.factor(gr6th)sixtimesormore 0.695561 0.156256 4.451 1.95e-05 ***
## as.factor(gr6th)thrtofivetimes 0.654549 0.161624 4.050 9.21e-05 ***
## as.factor(gr6ts)oneortwotimes 0.719430 0.357050 2.015 0.046185 *
## as.factor(gr6ts)sixtimesormore 0.718086 0.364966 1.968 0.051467 .
## as.factor(gr6ts)thrtofivetimes 0.317569 0.392951 0.808 0.420622
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 0.9900493)
##
## Number of Fisher Scoring iterations: 5
exp(coef(m0))
## (Intercept) age_3
## 3.958161e-05 1.266348e+00
## male_3 as.factor(educ_re3)high
## 8.241676e-01 2.176063e+00
## as.factor(educ_re3)lesshigh as.factor(educ_re3)somecol
## 1.959144e+00 2.154960e+00
## fri_sui3 fam_sui3
## 2.910013e+00 2.413247e+00
## as.factor(gr6th)oneortwotimes as.factor(gr6th)sixtimesormore
## 1.481487e+00 1.887612e+00
## as.factor(gr6th)thrtofivetimes as.factor(gr6ts)oneortwotimes
## 1.892806e+00 1.498791e+00
## as.factor(gr6ts)sixtimesormore as.factor(gr6ts)thrtofivetimes
## 2.910088e+00 7.454159e-01
exp(coef(m1))
## (Intercept) age_3
## 0.1099105 0.9942913
## male_3 as.factor(educ_re3)high
## 1.3750938 2.1833727
## as.factor(educ_re3)lesshigh as.factor(educ_re3)somecol
## 2.7409589 1.8965612
## fri_sui3 fam_sui3
## 1.1805602 1.3998099
## as.factor(gr6th)oneortwotimes as.factor(gr6th)sixtimesormore
## 1.1820426 2.0250610
## as.factor(gr6th)thrtofivetimes as.factor(gr6ts)oneortwotimes
## 1.8927183 2.0522567
## as.factor(gr6ts)sixtimesormore as.factor(gr6ts)thrtofivetimes
## 1.9358301 1.3694252
## male_3:as.factor(educ_re3)high male_3:as.factor(educ_re3)lesshigh
## 0.4253773 0.2968258
## male_3:as.factor(educ_re3)somecol male_3:fri_sui3
## 0.4312974 1.5115036
## male_3:fam_sui3
## 1.8989404
exp(coef(m2))
## (Intercept) age_3
## 0.1505018 0.9947078
## male_3 as.factor(educ_re3)high
## 0.6888274 1.4636355
## as.factor(educ_re3)lesshigh as.factor(educ_re3)somecol
## 1.5483408 1.2906595
## fri_sui3 fam_sui3
## 1.4371054 1.8264934
## as.factor(gr6th)oneortwotimes as.factor(gr6th)sixtimesormore
## 1.1932218 2.0048332
## as.factor(gr6th)thrtofivetimes as.factor(gr6ts)oneortwotimes
## 1.9242736 2.0532621
## as.factor(gr6ts)sixtimesormore as.factor(gr6ts)thrtofivetimes
## 2.0505053 1.3737843
## For the present research the shift was focused away from protective risk factors to those of adverse childhood experiences during grade 6. In this case if a child was hit, touched sexually, neglected, or left alone at the age of six, with each variable being that of a categorical one.
## Although a general discrete time model was fitted to the data, the only alternative time specification used was that of a linear model. Unfortunately, the nature of the Add Health survey only allows for the linear specification. As we only make use of three periods, and time is not continuous it would not be recommended to use any of the other models. Additionally, the survey waves take place in durations of at least 6 years prior to each wave being measured, so this also factors into why it is only recommended to use the linear model.
## The models were rerun after going over the summary statistics and finding in terms of adverse childhood experience, children being left alone and their needs not being met at grade 6 were not significant risk factors for predicting thoughts of suicide. Those two variables were then removed from all three models.
##M0 indicates the general model with all risk factors associated with increased risk of having thoughts of suicide. Having less education than college, friends and family recent suicide attempt, and being hit or touched sexually in sixth grade was associated with increase risk of a suicide attempt compared to those who were never. was associate with increased risk of suicidal thoughts
##M1 indicates the model with interaction terms being sex multiplied by family/friend suicide attempt, and the effect on education. With largely similar statistical significance however the predictive power of family/friend suicide disappears after interacting it with sex. Contrasting this with sex and education indicates that males with less than college education are less likely to think of suicide compared to females of similar educational standings.
##M2 The linear model presents some puzzling results. Only in this model the risk of thoughts of suicide is significantly less for males compared to females. The association between some college compared to college, and friend suicide has also lost predictive power. But perhaps the most egregious example is being touched sexually once or twice at grade six compared to those who had never was a significant risk of having thoughts of suicide, compared to any other category in that group.
##Generate hazard plots for interesting cases highlighting the significant predictors in your analysis
## Model AIC
AIC0<-AIC(m0)
## Warning in eval(family$initialize): non-integer #successes in a binomial glm!
AIC1<-AIC(m1)
## Warning in eval(family$initialize): non-integer #successes in a binomial glm!
AIC2<-AIC(m2)
## Warning in eval(family$initialize): non-integer #successes in a binomial glm!
AICS<-data.frame(AIC = c(AIC0["AIC"],AIC1["AIC"],AIC2["AIC"]),
mod = factor(c("general", "interaction", "linear")) )
AICS$mod<-forcats::fct_relevel(AICS$mod,
c("general", "interaction" , "linear"))
AICS$deltaAIC<-AICS$AIC - AICS$AIC[AICS$mod=="general"]
knitr::kable(AICS[, c("mod", "AIC", "deltaAIC")],
caption = "Relative AIC for alternative time specifications")
Relative AIC for alternative time specifications
|
mod
|
AIC
|
deltaAIC
|
|
general
|
2463.929
|
0.00
|
|
interaction
|
14321.956
|
11858.03
|
|
linear
|
14370.942
|
11907.01
|
##In terms of model fit the general model appears to be the best fit but rather interestingly the interaction and linear based models have relatively close aics.