### Simple structural equation model: Latent variable multiple regression 
# Parametric and with bootstrapping

### First, set the working directory where R will find the data

setwd("H:/data/User/Classes/multiv_13/YTs/YTSplit/Week 11 wmv - 2016/Wk11.R")
getwd()
## [1] "H:/data/User/Classes/multiv_13/YTs/YTSplit/Week 11 wmv - 2016/Wk11.R"
### Read the data
act<- read.csv("class12_practice_2014.csv")
head(act)
##   ï..const Cannabis   Heroin RecreationalPrescription  Cocaine
## 1 1 39.27526 37.96159 38.80631 35.61966
## 2 1 42.23284 46.17213 51.21147 57.59864
## 3 1 44.83347 46.99985 59.76539 49.09941
## 4 1 34.47021 41.91081 40.29767 43.06850
## 5 1 43.92398 41.56159 42.61687 38.01203
## 6 1 51.50516 44.16697 53.44381 44.51995
## InternalizedAnger ExternalizedAnger StateAnger TraitAnger ComorbidityIndex
## 1 53.69821 51.99119 53.46528 53.72893 38.54430
## 2 43.56667 51.99295 54.43380 44.36288 54.13281
## 3 63.13547 50.43772 48.84715 44.00640 42.28565
## 4 67.15893 48.21949 62.03162 51.11746 51.59931
## 5 60.95380 54.93101 53.09180 51.58114 40.10570
## 6 54.31116 59.33450 45.01612 56.69322 50.30999
## SubjectiveHealth FunctioningIndex FitnessIndex
## 1 46.61205 53.26231 48.87486
## 2 41.07006 52.10228 47.31790
## 3 58.08440 54.45158 50.20508
## 4 51.76131 48.43218 61.44567
## 5 55.22229 47.41978 58.48749
## 6 51.60633 45.14889 54.91742
library(lavaan)
## This is lavaan 0.6-7
## lavaan is BETA software! Please report any bugs.

act.mod4<- ' 
# measurement model
drugbelief =~ Cannabis + Heroin + RecreationalPrescription + Cocaine
anger =~ InternalizedAnger + ExternalizedAnger + StateAnger + TraitAnger
health =~ ComorbidityIndex + SubjectiveHealth + FunctioningIndex + FitnessIndex
#Correlated uniqueness
ExternalizedAnger ~~ TraitAnger
Heroin ~~ Cocaine
#correlated factor
drugbelief ~~ anger
# regressions
health~drugbelief+anger'

fit.act4 <- sem(act.mod4, data=act[2:13])
summary(fit.act4, fit.measures=TRUE, standardized=TRUE)
fitmeasures(fit.act4)
## lavaan 0.6-7 ended normally after 156 iterations
##
## Estimator ML
## Optimization method NLMINB
## Number of free parameters 29
##
## Number of observations 200
##
## Model Test User Model:
##
## Test statistic 55.570
## Degrees of freedom 49
## P-value (Chi-square) 0.241
##
## Model Test Baseline Model:
##
## Test statistic 960.853
## Degrees of freedom 66
## P-value 0.000
##
## User Model versus Baseline Model:
##
## Comparative Fit Index (CFI) 0.993
## Tucker-Lewis Index (TLI) 0.990
##
## Loglikelihood and Information Criteria:
##
## Loglikelihood user model (H0) -7610.508
## Loglikelihood unrestricted model (H1) -7582.723
##
## Akaike (AIC) 15279.015
## Bayesian (BIC) 15374.666
## Sample-size adjusted Bayesian (BIC) 15282.791
##
## Root Mean Square Error of Approximation:
##
## RMSEA 0.026
## 90 Percent confidence interval - lower 0.000
## 90 Percent confidence interval - upper 0.055
## P-value RMSEA <= 0.05 0.909
##
## Standardized Root Mean Square Residual:
##
## SRMR 0.042
##
## Parameter Estimates:
##
## Standard errors Standard
## Information Expected
## Information saturated (h1) model Structured
##
## Latent Variables:
## Estimate Std.Err z-value P(>|z|) Std.lv Std.all
## drugbelief =~
## Cannabis 1.000 9.263 0.943
## Heroin 0.542 0.049 11.062 0.000 5.019 0.743
## RcrtnlPrscrptn 0.573 0.053 10.823 0.000 5.307 0.727
## Cocaine 0.310 0.051 6.123 0.000 2.867 0.441
## anger =~
## InternalzdAngr 1.000 10.275 0.979
## ExternalzdAngr 0.532 0.042 12.716 0.000 5.467 0.750
## StateAnger 0.580 0.042 13.924 0.000 5.955 0.800
## TraitAnger 0.225 0.040 5.600 0.000 2.315 0.386
## health =~
## ComorbidtyIndx 1.000 5.300 0.853
## SubjectiveHlth 0.616 0.117 5.271 0.000 3.267 0.503
## FunctionngIndx 0.330 0.084 3.929 0.000 1.749 0.339
## FitnessIndex 0.526 0.103 5.128 0.000 2.787 0.481
##
## Regressions:
## Estimate Std.Err z-value P(>|z|) Std.lv Std.all
## health ~
## drugbelief 0.216 0.047 4.623 0.000 0.377 0.377
## anger 0.152 0.040 3.786 0.000 0.294 0.294
##
## Covariances:
## Estimate Std.Err z-value P(>|z|) Std.lv Std.all
## .ExternalizedAnger ~~
## .TraitAnger 12.063 2.195 5.496 0.000 12.063 0.452
## .Heroin ~~
## .Cocaine 11.535 2.315 4.982 0.000 11.535 0.437
## drugbelief ~~
## anger -9.454 7.212 -1.311 0.190 -0.099 -0.099
##
## Variances:
## Estimate Std.Err z-value P(>|z|) Std.lv Std.all
## .Cannabis 10.619 5.481 1.937 0.053 10.619 0.110
## .Heroin 20.407 2.597 7.857 0.000 20.407 0.448
## .RcrtnlPrscrptn 25.056 3.087 8.116 0.000 25.056 0.471
## .Cocaine 34.140 3.523 9.690 0.000 34.140 0.806
## .InternalzdAngr 4.604 5.010 0.919 0.358 4.604 0.042
## .ExternalzdAngr 23.205 2.723 8.521 0.000 23.205 0.437
## .StateAnger 19.997 2.613 7.652 0.000 19.997 0.361
## .TraitAnger 30.635 3.093 9.906 0.000 30.635 0.851
## .ComorbidtyIndx 10.513 4.219 2.492 0.013 10.513 0.272
## .SubjectiveHlth 31.587 3.623 8.718 0.000 31.587 0.747
## .FunctionngIndx 23.584 2.458 9.595 0.000 23.584 0.885
## .FitnessIndex 25.838 2.905 8.896 0.000 25.838 0.769
## drugbelief 85.800 10.989 7.808 0.000 1.000 1.000
## anger 105.567 12.085 8.735 0.000 1.000 1.000
## .health 22.283 4.900 4.547 0.000 0.793 0.793
#Bootstrapping

boot.act4<-(bootstrapLavaan(fit.act4))
summary(boot.act4)
#fitmeasures(boot.act4)
##  drugbelief=~Heroin drugbelief=~RecreationalPrescription drugbelief=~Cocaine
## Min. :0.3393 Min. :0.3898 Min. :0.1811
## 1st Qu.:0.5097 1st Qu.:0.5412 1st Qu.:0.2770
## Median :0.5420 Median :0.5740 Median :0.3079
## Mean :0.5424 Mean :0.5765 Mean :0.3112
## 3rd Qu.:0.5751 3rd Qu.:0.6148 3rd Qu.:0.3434
## Max. :0.7062 Max. :0.7584 Max. :0.4831
## anger=~ExternalizedAnger anger=~StateAnger anger=~TraitAnger
## Min. :0.4202 Min. :0.4299 Min. :0.08874
## 1st Qu.:0.5069 1st Qu.:0.5532 1st Qu.:0.19989
## Median :0.5323 Median :0.5794 Median :0.22639
## Mean :0.5324 Mean :0.5788 Mean :0.22648
## 3rd Qu.:0.5587 3rd Qu.:0.6076 3rd Qu.:0.25374
## Max. :0.6362 Max. :0.7156 Max. :0.36066
## health=~SubjectiveHealth health=~FunctioningIndex health=~FitnessIndex
## Min. :0.1142 Min. :0.0858 Min. :0.1251
## 1st Qu.:0.5088 1st Qu.:0.2715 1st Qu.:0.4336
## Median :0.6117 Median :0.3246 Median :0.5180
## Mean :0.6192 Mean :0.3251 Mean :0.5322
## 3rd Qu.:0.7249 3rd Qu.:0.3746 3rd Qu.:0.6170
## Max. :1.1953 Max. :0.6240 Max. :1.5448
## ExternalizedAnger~~TraitAnger Heroin~~Cocaine drugbelief~~anger
## Min. : 3.777 Min. : 4.564 Min. :-28.026
## 1st Qu.:10.049 1st Qu.: 9.741 1st Qu.:-13.716
## Median :11.823 Median :11.473 Median : -9.597
## Mean :11.814 Mean :11.470 Mean : -9.398
## 3rd Qu.:13.484 3rd Qu.:13.084 3rd Qu.: -5.013
## Max. :18.796 Max. :19.655 Max. : 10.330
## health~drugbelief health~anger Cannabis~~Cannabis Heroin~~Heroin
## Min. :0.05277 Min. :0.01711 Min. :-13.498 Min. :12.54
## 1st Qu.:0.18629 1st Qu.:0.11971 1st Qu.: 6.611 1st Qu.:18.21
## Median :0.21447 Median :0.14773 Median : 10.592 Median :20.12
## Mean :0.21404 Mean :0.14833 Mean : 10.359 Mean :20.10
## 3rd Qu.:0.24172 3rd Qu.:0.17534 3rd Qu.: 14.465 3rd Qu.:21.91
## Max. :0.40613 Max. :0.31175 Max. : 28.040 Max. :29.76
## RecreationalPrescription~~RecreationalPrescription Cocaine~~Cocaine
## Min. :16.73 Min. :21.84
## 1st Qu.:22.69 1st Qu.:30.59
## Median :24.68 Median :33.72
## Mean :24.81 Mean :33.95
## 3rd Qu.:26.85 3rd Qu.:36.96
## Max. :33.33 Max. :52.22
## InternalizedAnger~~InternalizedAnger ExternalizedAnger~~ExternalizedAnger
## Min. :-13.860 Min. :15.53
## 1st Qu.: 1.424 1st Qu.:21.03
## Median : 4.653 Median :22.75
## Mean : 4.433 Mean :22.84
## 3rd Qu.: 7.618 3rd Qu.:24.55
## Max. : 17.844 Max. :30.99
## StateAnger~~StateAnger TraitAnger~~TraitAnger
## Min. :13.55 Min. :20.86
## 1st Qu.:18.27 1st Qu.:27.75
## Median :19.78 Median :29.94
## Mean :19.83 Mean :30.16
## 3rd Qu.:21.34 3rd Qu.:32.36
## Max. :27.87 Max. :43.03
## ComorbidityIndex~~ComorbidityIndex SubjectiveHealth~~SubjectiveHealth
## Min. :-45.480 Min. :20.64
## 1st Qu.: 5.579 1st Qu.:28.65
## Median : 9.940 Median :31.06
## Mean : 9.229 Mean :31.10
## 3rd Qu.: 13.645 3rd Qu.:33.50
## Max. : 24.680 Max. :43.31
## FunctioningIndex~~FunctioningIndex FitnessIndex~~FitnessIndex
## Min. :15.27 Min. :15.12
## 1st Qu.:21.80 1st Qu.:22.67
## Median :23.52 Median :25.15
## Mean :23.47 Mean :25.37
## 3rd Qu.:25.09 3rd Qu.:27.81
## Max. :31.99 Max. :39.68
## drugbelief~~drugbelief anger~~anger health~~health
## Min. : 56.38 Min. : 68.10 Min. : 7.514
## 1st Qu.: 77.96 1st Qu.: 96.77 1st Qu.:18.922
## Median : 85.07 Median :104.08 Median :22.427
## Mean : 85.25 Mean :104.73 Mean :23.139
## 3rd Qu.: 92.15 3rd Qu.:112.34 3rd Qu.:26.539
## Max. :127.34 Max. :139.30 Max. :77.496
boot.act4.df<-as.data.frame(boot.act4)
head(boot.act4.df)
##   drugbelief=~Heroin drugbelief=~RecreationalPrescription drugbelief=~Cocaine
## 1 0.5306324 0.5820211 0.2739624
## 2 0.5324146 0.5541148 0.3545429
## 3 0.4483420 0.4213323 0.2480484
## 4 0.5102962 0.5707966 0.2270199
## 5 0.5719150 0.5633284 0.3177219
## 6 0.6493614 0.6834911 0.2890858
## anger=~ExternalizedAnger anger=~StateAnger anger=~TraitAnger
## 1 0.5509046 0.6096060 0.2935393
## 2 0.5457708 0.5890187 0.2722200
## 3 0.4536425 0.5909232 0.1595334
## 4 0.5327138 0.5081917 0.1895334
## 5 0.5698685 0.5535323 0.2649729
## 6 0.5290908 0.5531622 0.2754555
## health=~SubjectiveHealth health=~FunctioningIndex health=~FitnessIndex
## 1 0.7775360 0.3993354 0.7546237
## 2 0.5408158 0.3622083 0.7165513
## 3 0.3746239 0.1545994 0.4090513
## 4 0.8465778 0.2860563 0.5772059
## 5 0.5670046 0.3473295 0.4429461
## 6 0.4359502 0.3956468 0.4797628
## ExternalizedAnger~~TraitAnger Heroin~~Cocaine drugbelief~~anger
## 1 7.926916 13.360828 1.966133
## 2 13.201764 12.412088 -11.826336
## 3 16.511724 17.296815 -14.266460
## 4 14.575275 10.710551 -4.174472
## 5 14.122003 11.269143 -16.755728
## 6 10.265001 9.291838 -4.125708
## health~drugbelief health~anger Cannabis~~Cannabis Heroin~~Heroin
## 1 0.1665823 0.1426707 12.202182 18.87461
## 2 0.2606514 0.1697081 10.058974 21.92845
## 3 0.2193985 0.1155621 -13.498248 27.00665
## 4 0.2573788 0.1556157 8.200333 19.15658
## 5 0.1765568 0.1382384 7.435513 20.58540
## 6 0.2162261 0.2267753 14.720383 13.51219
## RecreationalPrescription~~RecreationalPrescription Cocaine~~Cocaine
## 1 21.28774 39.18077
## 2 31.47606 36.71030
## 3 25.09081 44.72043
## 4 20.69933 29.38347
## 5 23.05223 31.43890
## 6 28.00981 34.97539
## InternalizedAnger~~InternalizedAnger ExternalizedAnger~~ExternalizedAnger
## 1 6.62312680 19.71793
## 2 6.41807469 22.97717
## 3 3.03143902 28.37988
## 4 3.31161393 23.24621
## 5 -0.04889051 24.24528
## 6 6.14797658 22.31119
## StateAnger~~StateAnger TraitAnger~~TraitAnger
## 1 16.41461 26.96738
## 2 20.37653 30.41898
## 3 22.50772 32.71626
## 4 23.44050 33.59271
## 5 20.62584 29.47410
## 6 19.37756 33.86184
## ComorbidityIndex~~ComorbidityIndex SubjectiveHealth~~SubjectiveHealth
## 1 14.908360 25.23323
## 2 13.383856 24.18982
## 3 -6.139919 37.44216
## 4 13.969106 23.34736
## 5 7.958794 36.41906
## 6 12.816144 29.30958
## FunctioningIndex~~FunctioningIndex FitnessIndex~~FitnessIndex
## 1 24.08676 30.07860
## 2 25.05585 21.99100
## 3 23.86151 31.32080
## 4 24.86028 21.30004
## 5 22.92639 28.49262
## 6 21.78320 26.82356
## drugbelief~~drugbelief anger~~anger health~~health
## 1 79.09153 90.19545 16.77512
## 2 68.77235 91.01809 20.75482
## 3 120.61872 117.42422 32.54614
## 4 89.02119 123.26487 22.73888
## 5 96.01150 108.02994 27.17884
## 6 80.41765 90.67876 25.38392
colnames(boot.act4.df) <- c("lx.drug_her","lx.drug_rec","lx.drug_coc","lx.ang_ext",
"lx.ang_state", "lx.ang_trait", "lx.htlh_subj", "lx.hlth_func","lx.hlth_fit",
"corres.extang_traitang", "corres.her_coc", "ph.drug_ang", "ph.health_drug",
"ph.health_ang", "res_canna", "res.her", "res.recpresc", "res.coc",
"res.intang", "res.extang", "res.stateang", "res.traitang", "res.comorb",
"res.subjhlth", "res.func", "res.fit", "var.drug", "var.ang", "var.hlth")

quantile(boot.act4.df$lx.drug_her, c(0.025,0.975))
##      2.5%     97.5% 
## 0.4453478 0.6416322
quantile(boot.act4.df$lx.drug_rec, c(0.025,0.975)) 
##      2.5%     97.5% 
## 0.4498283 0.6957687
quantile(boot.act4.df$lx.drug_coc, c(0.025,0.975)) 
##      2.5%     97.5% 
## 0.2138968 0.4233742
quantile(boot.act4.df$lx.ang_ext, c(0.025,0.975)) 
##      2.5%     97.5% 
## 0.4576703 0.6032620
quantile(boot.act4.df$lx.ang_state, c(0.025,0.975)) 
##      2.5%     97.5% 
## 0.4983117 0.6518977
quantile(boot.act4.df$lx.ang_trait, c(0.025,0.975)) 
##      2.5%     97.5% 
## 0.1380021 0.3110086
quantile(boot.act4.df$lx.htlh_subj, c(0.025,0.975)) 
##      2.5%     97.5% 
## 0.3229178 0.9420919
quantile(boot.act4.df$lx.hlth_func, c(0.025,0.975)) 
##      2.5%     97.5% 
## 0.1851500 0.4768042
quantile(boot.act4.df$lx.hlth_fit, c(0.025,0.975)) 
##      2.5%     97.5% 
## 0.2707904 0.8575554
quantile(boot.act4.df$corres.extang_traitang, c(0.025,0.975)) 
##      2.5%     97.5% 
## 6.852617 16.814963
quantile(boot.act4.df$corres.her_coc, c(0.025,0.975)) 
##      2.5%     97.5% 
## 6.587128 16.715327
quantile(boot.act4.df$ph.drug_ang, c(0.025,0.975)) 
##       2.5%      97.5% 
## -21.612454 3.469618
quantile(boot.act4.df$ph.health_drug, c(0.025,0.975)) 
##      2.5%     97.5% 
## 0.1282511 0.3021391
quantile(boot.act4.df$ph.health_ang, c(0.025,0.975)) 
##       2.5%      97.5% 
## 0.06631403 0.23864955
quantile(boot.act4.df$res_canna, c(0.025,0.975)) 
##      2.5%     97.5% 
## -1.780251 21.684686
quantile(boot.act4.df$res.her, c(0.025,0.975)) 
##     2.5%    97.5% 
## 14.97063 25.73912
quantile(boot.act4.df$res.recpresc, c(0.025,0.975)) 
##     2.5%    97.5% 
## 19.45648 30.89448
quantile(boot.act4.df$res.coc, c(0.025,0.975)) 
##     2.5%    97.5% 
## 25.89767 44.00118
quantile(boot.act4.df$res.intang, c(0.025,0.975)) 
##     2.5%    97.5% 
## -5.78297 13.30813
quantile(boot.act4.df$res.extang, c(0.025,0.975)) 
##     2.5%    97.5% 
## 18.31223 28.27189
quantile(boot.act4.df$res.stateang, c(0.025,0.975)) 
##     2.5%    97.5% 
## 15.70952 24.33129
quantile(boot.act4.df$res.traitang, c(0.025,0.975)) 
##     2.5%    97.5% 
## 23.72643 37.27873
quantile(boot.act4.df$res.comorb, c(0.025,0.975)) 
##      2.5%     97.5% 
## -4.939536 19.718825
quantile(boot.act4.df$res.subjhlth, c(0.025,0.975)) 
##     2.5%    97.5% 
## 24.43650 38.19684
quantile(boot.act4.df$res.func, c(0.025,0.975)) 
##    2.5%   97.5% 
## 18.6520 28.6282
quantile(boot.act4.df$res.fit, c(0.025,0.975)) 
##     2.5%    97.5% 
## 18.71365 32.79103
quantile(boot.act4.df$var.drug, c(0.025,0.975)) 
##      2.5%     97.5% 
## 66.35804 105.25409
quantile(boot.act4.df$var.ang, c(0.025,0.975)) 
##      2.5%     97.5% 
## 84.07243 127.90763
quantile(boot.act4.df$var.hlth, c(0.025,0.975)) 
##     2.5%    97.5% 
## 12.41092 37.95433
#install.packages("lavaanPlot")
library(lavaanPlot)
## Warning: package 'lavaanPlot' was built under R version 4.0.3

lavaanPlot(model = fit.act4, sig=.05,stars=c("regress","latent","covs"),
node_options = list(shape = "box", fontname = "Arial"),
edge_options = list(color = "black"),
coefs = TRUE, covs = TRUE, stand=TRUE)