library(psych)
## Warning: package 'psych' was built under R version 4.5.3
library(tidyverse)
## Warning: package 'ggplot2' was built under R version 4.5.3
## Warning: package 'tibble' was built under R version 4.5.3
## Warning: package 'tidyr' was built under R version 4.5.3
## Warning: package 'readr' was built under R version 4.5.3
## Warning: package 'purrr' was built under R version 4.5.3
## Warning: package 'dplyr' was built under R version 4.5.3
## Warning: package 'stringr' was built under R version 4.5.3
## Warning: package 'forcats' was built under R version 4.5.3
## Warning: package 'lubridate' was built under R version 4.5.3
## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr 1.2.1 ✔ readr 2.2.0
## ✔ forcats 1.0.1 ✔ stringr 1.6.0
## ✔ ggplot2 4.0.3 ✔ tibble 3.3.1
## ✔ lubridate 1.9.5 ✔ tidyr 1.3.2
## ✔ purrr 1.2.2
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ ggplot2::%+%() masks psych::%+%()
## ✖ ggplot2::alpha() masks psych::alpha()
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag() masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
library (ggplot2)
library(lavaan)
## Warning: package 'lavaan' was built under R version 4.5.3
## This is lavaan 0.7-2
## lavaan is FREE software! Please report any bugs.
##
## Attaching package: 'lavaan'
##
## The following object is masked from 'package:psych':
##
## cor2cov
library(lavaanPlot)
## Warning: package 'lavaanPlot' was built under R version 4.5.2
library(stats)
pov <- read.csv("260724_Final_Cleaned_21_3720_Blanks.csv", header = TRUE)
pov.model <- '
income =~ conrinc + coninc + satfin
hlth =~ health + hlthstrt
education =~ educ + degree
econ =~ weekswrk + partfull + wrkslf + prestg105plus
pol =~ vote20 + lkelyvot + polint
soc =~ socommun + socfrend + socrel
income ~~ education
income ~ hlth
income ~ econ
income ~ pol
income ~ soc'
#Code taken from https://slesche.github.io/r-for-empra/instructor/factor-analysis-02.html
#fiml code taken from https://www.geeksforgeeks.org/machine-learning/full-information-maximum-likelihood-for-missing-data-in-r/
fit_pov <- cfa(model = pov.model, data = pov, missing = "fiml")
## Warning: lavaan->lav_data_full():
## some observed variances are (at least) a factor 1000 times larger than
## others; use varTable(fit) to investigate
## Warning: lavaan->lav_data_full():
## due to missing values, some pairwise combinations have less than 10%
## coverage; use lavInspect(fit, "coverage") to investigate.
summary(fit_pov, fit.measures = TRUE, standardized = TRUE)
## lavaan 0.7-2 ended normally after 359 iterations
##
## Estimator ML
## Optimization method NLMINB
## Number of model parameters 66
##
## Number of observations 3720
## Number of missing patterns 205
##
## Model Test User Model:
##
## Test statistic 2805.812
## Degrees of freedom 104
## P-value (Chi-square) 0.000
##
## Model Test Baseline Model:
##
## Test statistic 12549.577
## Degrees of freedom 136
## P-value 0.000
##
## User Model versus Baseline Model:
##
## Comparative Fit Index (CFI) 0.782
## Tucker-Lewis Index (TLI) 0.715
##
## Robust Comparative Fit Index (CFI) 0.704
## Robust Tucker-Lewis Index (TLI) 0.612
##
## Loglikelihood and Information Criteria:
##
## Loglikelihood user model (H0) -76730.891
## Loglikelihood unrestricted model (H1) -75327.985
##
## Akaike (AIC) 153593.781
## Bayesian (BIC) 154004.399
## Sample-size adjusted Bayesian (SABIC) 153794.683
##
## Root Mean Square Error of Approximation:
##
## RMSEA 0.084
## 90 Percent confidence interval - lower 0.081
## 90 Percent confidence interval - upper 0.086
## P-value H_0: RMSEA <= 0.050 0.000
## P-value H_0: RMSEA >= 0.080 0.986
##
## Robust RMSEA 0.113
## 90 Percent confidence interval - lower 0.109
## 90 Percent confidence interval - upper 0.117
## P-value H_0: Robust RMSEA <= 0.050 0.000
## P-value H_0: Robust RMSEA >= 0.080 1.000
##
## Standardized Root Mean Square Residual:
##
## SRMR 0.090
##
## Goodness of Fit Index:
##
## Goodness of Fit Index (GFI) 0.921
## 90 Percent confidence interval - lower 0.917
## 90 Percent confidence interval - upper 0.926
##
## Robust GFI 0.865
## 90 Percent confidence interval - lower 0.857
## 90 Percent confidence interval - upper 0.874
##
## Parameter Estimates:
##
## Standard errors Standard
## Information Observed
## Observed information based on Hessian
##
## Latent Variables:
## Estimate Std.Err z-value P(>|z|) Std.lv Std.all
## income =~
## conrinc 1.000 0.284 0.236
## coninc 2.231 0.241 9.252 0.000 0.633 0.573
## satfin 1.144 0.220 5.194 0.000 0.325 0.442
## hlth =~
## health 1.000 0.543 0.724
## hlthstrt 0.526 0.126 4.171 0.000 0.286 0.399
## education =~
## educ 1.000 0.910 0.886
## degree 1.318 0.023 57.880 0.000 1.199 0.949
## econ =~
## weekswrk 1.000 3.372 0.146
## partfull 0.023 0.005 4.908 0.000 0.076 0.184
## wrkslf 0.000 0.002 0.102 0.919 0.001 0.002
## prestg105plus 0.277 0.045 6.097 0.000 0.933 0.688
## pol =~
## vote20 1.000 0.323 0.767
## lkelyvot 3.607 0.131 27.434 0.000 1.164 0.868
## polint 2.001 0.104 19.205 0.000 0.646 0.506
## soc =~
## socommun 1.000 0.837 0.427
## socfrend 1.391 0.162 8.604 0.000 1.164 0.723
## socrel 0.829 0.069 11.949 0.000 0.694 0.404
##
## Regressions:
## Estimate Std.Err z-value P(>|z|) Std.lv Std.all
## income ~
## hlth 0.136 0.042 3.262 0.001 0.261 0.261
## econ 0.037 0.008 4.454 0.000 0.440 0.440
## pol 0.109 0.033 3.328 0.001 0.124 0.124
## soc -0.003 0.016 -0.186 0.853 -0.009 -0.009
##
## Covariances:
## Estimate Std.Err z-value P(>|z|) Std.lv Std.all
## .income ~~
## education -0.009 0.011 -0.861 0.389 -0.046 -0.046
## hlth ~~
## education 0.155 0.012 12.944 0.000 0.314 0.314
## econ 0.561 0.118 4.748 0.000 0.306 0.306
## pol 0.035 0.005 7.458 0.000 0.197 0.197
## soc 0.139 0.019 7.453 0.000 0.305 0.305
## education ~~
## econ 2.453 0.407 6.031 0.000 0.799 0.799
## pol 0.101 0.007 14.878 0.000 0.344 0.344
## soc 0.053 0.021 2.481 0.013 0.069 0.069
## econ ~~
## pol 0.360 0.063 5.691 0.000 0.331 0.331
## soc 0.256 0.122 2.102 0.036 0.091 0.091
## pol ~~
## soc 0.042 0.008 4.960 0.000 0.154 0.154
##
## Intercepts:
## Estimate Std.Err z-value P(>|z|) Std.lv Std.all
## .conrinc 2.813 0.020 142.828 0.000 2.813 2.342
## .coninc 2.336 0.018 128.791 0.000 2.336 2.112
## .satfin 1.902 0.012 157.363 0.000 1.902 2.587
## .health 2.828 0.012 229.480 0.000 2.828 3.772
## .hlthstrt 3.144 0.027 117.560 0.000 3.144 4.386
## .educ 4.097 0.017 243.165 0.000 4.097 3.992
## .degree 1.895 0.021 91.443 0.000 1.895 1.500
## .weekswrk 30.991 0.386 80.273 0.000 30.991 1.338
## .partfull 1.776 0.008 213.044 0.000 1.776 4.294
## .wrkslf 1.113 0.005 210.833 0.000 1.113 3.513
## .prestg105plus 2.911 0.023 127.986 0.000 2.911 2.145
## .vote20 1.761 0.007 247.102 0.000 1.761 4.188
## .lkelyvot 4.094 0.023 178.826 0.000 4.094 3.054
## .polint 2.998 0.029 104.196 0.000 2.998 2.348
## .socommun 2.987 0.040 75.345 0.000 2.987 1.526
## .socfrend 3.689 0.032 113.641 0.000 3.689 2.293
## .socrel 4.320 0.035 124.234 0.000 4.320 2.515
##
## Variances:
## Estimate Std.Err z-value P(>|z|) Std.lv Std.all
## .conrinc 1.362 0.037 36.333 0.000 1.362 0.944
## .coninc 0.822 0.060 13.749 0.000 0.822 0.672
## .satfin 0.435 0.017 25.096 0.000 0.435 0.805
## .health 0.267 0.070 3.804 0.000 0.267 0.475
## .hlthstrt 0.432 0.030 14.226 0.000 0.432 0.841
## .educ 0.226 0.013 17.508 0.000 0.226 0.214
## .degree 0.159 0.021 7.651 0.000 0.159 0.100
## .weekswrk 524.886 12.531 41.886 0.000 524.886 0.979
## .partfull 0.165 0.005 34.286 0.000 0.165 0.966
## .wrkslf 0.100 0.002 42.444 0.000 0.100 1.000
## .prestg105plus 0.971 0.097 9.977 0.000 0.971 0.527
## .vote20 0.073 0.004 18.766 0.000 0.073 0.411
## .lkelyvot 0.442 0.045 9.723 0.000 0.442 0.246
## .polint 1.214 0.045 27.163 0.000 1.214 0.744
## .socommun 3.134 0.121 25.856 0.000 3.134 0.817
## .socfrend 1.235 0.153 8.052 0.000 1.235 0.477
## .socrel 2.470 0.092 26.931 0.000 2.470 0.837
## .income 0.049 0.017 2.890 0.004 0.606 0.606
## hlth 0.295 0.071 4.160 0.000 1.000 1.000
## education 0.828 0.027 31.016 0.000 1.000 1.000
## econ 11.373 3.260 3.489 0.000 1.000 1.000
## pol 0.104 0.005 19.761 0.000 1.000 1.000
## soc 0.700 0.103 6.791 0.000 1.000 1.000
lavaanPlot(model = fit_pov, node_options = list(shape = "box", fontname = "Helvetica"), edge_options = list(color = "red"), coefs = TRUE, stars = "latent")