# Fit interaction models
rm(list=ls())
# Load libraries ---------------------------
library(foreign)
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
# Set wd for generated output ---------------------------
setwd("/Volumes/caas/CMHIV/WISE/AK-R-output")
# Load data ---------------------------
load(file="table2.RData")
options(na.action='na.pass')
## Adjustment variables: drugp1, smhis1_6m, scrn73_b.out
levels(vuln_dt$drugp1)
## [1] "No" "Yes" "Don't Know" "Refuse to Answer"
## [5] "Not Applicable"
levels(vuln_dt$smhis1_6m)
## [1] "6 months or less" "More than 6 months"
table(scrn7r_b.out, exclude = NULL)
## scrn7r_b.out
## 0 1 <NA>
## 126 120 1
drugp1 <- na_if(vuln_dt$drugp1, "Not Applicable")
levels(drugp1)
## [1] "No" "Yes" "Don't Know" "Refuse to Answer"
## [5] "Not Applicable"
## Outcomes
outcomes <- cbind(Final_QuitStatus.out)
## Future Vulnerability
FPV.logreg <- apply(outcomes, 2, function (x) {
glm(x ~ FPV * vuln_dt$trt_grp+
drugp1+ #drug treatment program in the ACI
vuln_dt$smhis1_6m+ #
scrn7r_b.out,
family="binomial",
na.action = 'na.omit'
)
}
)
lapply(FPV.logreg, summary)
## $Final_QuitStatus.out
##
## Call:
## glm(formula = x ~ FPV * vuln_dt$trt_grp + drugp1 + vuln_dt$smhis1_6m +
## scrn7r_b.out, family = "binomial", na.action = "na.omit")
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.3314 -0.5460 -0.3566 -0.1621 2.7358
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -2.00545 1.26783 -1.582 0.1137
## FPV -0.15830 0.08236 -1.922 0.0546 .
## vuln_dt$trt_grpIntervention -1.40072 1.49097 -0.939 0.3475
## drugp1Yes 0.66156 0.42038 1.574 0.1155
## vuln_dt$smhis1_6mMore than 6 months 1.41559 0.46639 3.035 0.0024 **
## scrn7r_b.out 0.52926 0.41383 1.279 0.2009
## FPV:vuln_dt$trt_grpIntervention 0.21329 0.10141 2.103 0.0355 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 207.26 on 232 degrees of freedom
## Residual deviance: 165.97 on 226 degrees of freedom
## (14 observations deleted due to missingness)
## AIC: 179.97
##
## Number of Fisher Scoring iterations: 6
lapply(FPV.logreg, function (x) exp(coef(x)))
## $Final_QuitStatus.out
## (Intercept) FPV
## 0.1346002 0.8535960
## vuln_dt$trt_grpIntervention drugp1Yes
## 0.2464200 1.9378193
## vuln_dt$smhis1_6mMore than 6 months scrn7r_b.out
## 4.1188954 1.6976699
## FPV:vuln_dt$trt_grpIntervention
## 1.2377379
## Future Precaution
FUTPREC.logreg <- apply(outcomes, 2, function (x) {
glm(x ~ FUTPREC * vuln_dt$trt_grp+
drugp1+
vuln_dt$smhis1_6m+
scrn7r_b.out,
family="binomial",
na.action = 'na.omit'
)
}
)
lapply(FUTPREC.logreg, summary)
## $Final_QuitStatus.out
##
## Call:
## glm(formula = x ~ FUTPREC * vuln_dt$trt_grp + drugp1 + vuln_dt$smhis1_6m +
## scrn7r_b.out, family = "binomial", na.action = "na.omit")
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.2715 -0.5490 -0.3695 -0.1857 3.0631
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -2.6760 1.5783 -1.695 0.08999 .
## FUTPREC -0.5015 0.4546 -1.103 0.26988
## vuln_dt$trt_grpIntervention 0.5809 1.9438 0.299 0.76504
## drugp1Yes 0.4361 0.4102 1.063 0.28776
## vuln_dt$smhis1_6mMore than 6 months 1.4596 0.4813 3.033 0.00242 **
## scrn7r_b.out 0.6954 0.4148 1.676 0.09367 .
## FUTPREC:vuln_dt$trt_grpIntervention 0.3627 0.5523 0.657 0.51137
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 199.23 on 228 degrees of freedom
## Residual deviance: 162.65 on 222 degrees of freedom
## (18 observations deleted due to missingness)
## AIC: 176.65
##
## Number of Fisher Scoring iterations: 6
lapply(FUTPREC.logreg, function (x) exp(coef(x)))
## $Final_QuitStatus.out
## (Intercept) FUTPREC
## 0.06884061 0.60559913
## vuln_dt$trt_grpIntervention drugp1Yes
## 1.78771255 1.54662362
## vuln_dt$smhis1_6mMore than 6 months scrn7r_b.out
## 4.30410203 2.00457259
## FUTPREC:vuln_dt$trt_grpIntervention
## 1.43723255
## Relative Pessimism
RELPESS.logreg <- apply(outcomes, 2, function (x) {
glm(x ~ RELPESS * vuln_dt$trt_grp+
drugp1+
vuln_dt$smhis1_6m+
scrn7r_b.out,
family="binomial",
na.action = 'na.omit'
)
}
)
lapply(RELPESS.logreg, summary)
## $Final_QuitStatus.out
##
## Call:
## glm(formula = x ~ RELPESS * vuln_dt$trt_grp + drugp1 + vuln_dt$smhis1_6m +
## scrn7r_b.out, family = "binomial", na.action = "na.omit")
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.3775 -0.5411 -0.3863 -0.1951 2.8986
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -3.33963 1.24013 -2.693 0.00708 **
## RELPESS -0.31730 0.39653 -0.800 0.42361
## vuln_dt$trt_grpIntervention 1.65144 1.29418 1.276 0.20194
## drugp1Yes 0.42238 0.41247 1.024 0.30582
## vuln_dt$smhis1_6mMore than 6 months 1.35747 0.47022 2.887 0.00389 **
## scrn7r_b.out 0.65706 0.41465 1.585 0.11306
## RELPESS:vuln_dt$trt_grpIntervention 0.02747 0.45110 0.061 0.95144
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 198.20 on 225 degrees of freedom
## Residual deviance: 160.66 on 219 degrees of freedom
## (21 observations deleted due to missingness)
## AIC: 174.66
##
## Number of Fisher Scoring iterations: 6
lapply(RELPESS.logreg, function (x) exp(coef(x)))
## $Final_QuitStatus.out
## (Intercept) RELPESS
## 0.03545017 0.72811532
## vuln_dt$trt_grpIntervention drugp1Yes
## 5.21449523 1.52558877
## vuln_dt$smhis1_6mMore than 6 months scrn7r_b.out
## 3.88636373 1.92911227
## RELPESS:vuln_dt$trt_grpIntervention
## 1.02785410
## Current Vulnerability
ropbais5.logreg <- apply(outcomes, 2, function (x) {
glm(x ~ ropbais5 * vuln_dt$trt_grp+
drugp1+ #drug treatment program in the ACI
vuln_dt$smhis1_6m+ #
scrn7r_b.out,
family="binomial",
na.action = 'na.omit'
)
}
)
lapply(ropbais5.logreg, summary)
## $Final_QuitStatus.out
##
## Call:
## glm(formula = x ~ ropbais5 * vuln_dt$trt_grp + drugp1 + vuln_dt$smhis1_6m +
## scrn7r_b.out, family = "binomial", na.action = "na.omit")
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.1939 -0.5519 -0.3678 -0.1720 2.7675
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -3.8074 0.6912 -5.508 3.62e-08 ***
## ropbais51 -1.0504 0.7763 -1.353 0.1760
## vuln_dt$trt_grpIntervention 1.1919 0.5975 1.995 0.0461 *
## drugp1Yes 0.4945 0.4036 1.225 0.2204
## vuln_dt$smhis1_6mMore than 6 months 1.5083 0.4623 3.262 0.0011 **
## scrn7r_b.out 0.6516 0.4156 1.568 0.1169
## ropbais51:vuln_dt$trt_grpIntervention 1.0061 0.8989 1.119 0.2630
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 207.97 on 234 degrees of freedom
## Residual deviance: 169.18 on 228 degrees of freedom
## (12 observations deleted due to missingness)
## AIC: 183.18
##
## Number of Fisher Scoring iterations: 6
lapply(ropbais5.logreg, function (x) exp(coef(x)))
## $Final_QuitStatus.out
## (Intercept) ropbais51
## 0.02220505 0.34978303
## vuln_dt$trt_grpIntervention drugp1Yes
## 3.29336968 1.63969619
## vuln_dt$smhis1_6mMore than 6 months scrn7r_b.out
## 4.51890556 1.91851418
## ropbais51:vuln_dt$trt_grpIntervention
## 2.73502916
## Future Pessimism
OPBIAS4R.logreg <- apply(outcomes, 2, function (x) {
glm(x ~ OPBIAS4R_2way * vuln_dt$trt_grp +
drugp1+
vuln_dt$smhis1_6m+
scrn7r_b.out,
family="binomial",
na.action = 'na.omit'
)
}
)
lapply(OPBIAS4R.logreg, summary)
## $Final_QuitStatus.out
##
## Call:
## glm(formula = x ~ OPBIAS4R_2way * vuln_dt$trt_grp + drugp1 +
## vuln_dt$smhis1_6m + scrn7r_b.out, family = "binomial", na.action = "na.omit")
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.2140 -0.5404 -0.3485 -0.1885 2.8424
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -4.0218 0.6459 -6.226 4.78e-10
## OPBIAS4R_2way1 -1.3633 1.1132 -1.225 0.22070
## vuln_dt$trt_grpIntervention 1.5819 0.5188 3.049 0.00229
## drugp1Yes 0.4623 0.4088 1.131 0.25807
## vuln_dt$smhis1_6mMore than 6 months 1.4733 0.4639 3.176 0.00149
## scrn7r_b.out 0.5898 0.4074 1.448 0.14767
## OPBIAS4R_2way1:vuln_dt$trt_grpIntervention 1.0322 1.2175 0.848 0.39652
##
## (Intercept) ***
## OPBIAS4R_2way1
## vuln_dt$trt_grpIntervention **
## drugp1Yes
## vuln_dt$smhis1_6mMore than 6 months **
## scrn7r_b.out
## OPBIAS4R_2way1:vuln_dt$trt_grpIntervention
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 203.61 on 231 degrees of freedom
## Residual deviance: 164.51 on 225 degrees of freedom
## (15 observations deleted due to missingness)
## AIC: 178.51
##
## Number of Fisher Scoring iterations: 6
lapply(OPBIAS4R.logreg, function (x) exp(coef(x)))
## $Final_QuitStatus.out
## (Intercept)
## 0.01792132
## OPBIAS4R_2way1
## 0.25582705
## vuln_dt$trt_grpIntervention
## 4.86432645
## drugp1Yes
## 1.58774443
## vuln_dt$smhis1_6mMore than 6 months
## 4.36380987
## scrn7r_b.out
## 1.80362149
## OPBIAS4R_2way1:vuln_dt$trt_grpIntervention
## 2.80733284
# Save object ---------------------------
save.image(file="table3.RData")