%%{init: {"flowchart": {"htmlLabels": false}} }%%
graph TD
subgraph "Measurement Procedure"
subgraph Modeling
direction LR
Start --> id1
id1([Specify model]) --> id2(["`Specify scale identification
method for latent variables`"])
id2 --> id3([Estimate model])
id3 --> id4{Adequate fit indices}
id4 -- No --> id1
end
Modeling -- Model fit acceptable --> one
subgraph one["Evaluation: Measurement Model with Auxiliary Components"]
direction LR
id5(["Score estimation
and secondary analyses with scores"]) --> id6([Item evaluation])
id6 --> id7([Scale construction])
id7 --> id8([Equating])
id8 --> id9([Measure Invariance/Differential item functioning])
id9 --> End
end
end
style Start fill:#f9f,stroke:#333,stroke-width:5px
style End fill:#bbf,stroke:#333,stroke-width:5px
Lecture 06
Generalized Measurement Models: An Introduction
Bayesian view: latent variables
Latent variables in Bayesian are built by following specification:
-
What are their distributions? (normal distribution or others)
For example, \(\theta_i\) values for one person and \(\theta\) values for samples. Factor score \(\theta\) is a mixture distribution of distributions of each individual’s factor score \(\theta_i\)
But, in MLE/WSLMV, we do not estimate mean and sd of each individual’s factor score for model to be converged
Code
library(tidyverse)
set.seed(12)
N = 15
ndraws = 200
FS = matrix(NA, nrow = N, ndraws)
FS_p = rnorm(N)
FS_p = FS_p[order(FS_p)]
for (i in 1:N) {
FS[i,] = rnorm(ndraws, mean = FS_p[i], sd = 1)
}
FS_plot <- as.data.frame(t(FS))
colnames(FS_plot) <- paste0("Person", 1:N)
FS_plot <- FS_plot |> pivot_longer(everything(), names_to = "Person", values_to = "Factor Score")
FS_plot$Person <- factor(FS_plot$Person, levels = paste0("Person", 1:N))
ggplot() +
geom_density(aes(x = `Factor Score`, fill = Person, col = Person ), alpha = .5, data = FS_plot) +
geom_density(aes(x = FS_p))
Simulation Study 1
- Let’s perform a small simulation study to see how to perform factor analysis in naive Stan.
- The model specification is a two-factor with each measured by 3 items. In total, there are 6 items with continuous responses. Sample size is 1000.
Model Specification
MCMC Result > Model diagnostic > PPMC
All items have PPP of item mean close to 0.5, suggesting great model-data fitting.
Item_Mean_rep_mat <- fit_cfa_twofactor$draws("Item_Mean_rep", format = 'matrix')
Item_Mean_obs <- colMeans(Y)
PPP <- rep(NA, J)
# colMeans(Item_Mean_rep_mat)
for (item in 1:J) {
PPP[item] <- mean(as.numeric(Item_Mean_rep_mat[,item]) > Item_Mean_obs[item])
}
PPP[1] 0.491250 0.505750 0.497500 0.499875 0.503250 0.486500data.frame(
Item = factor(1:J, levels = 1:J),
PPP = PPP
) |>
ggplot() +
geom_col(aes(x = Item, y = PPP)) +
geom_hline(aes(yintercept = .5), col = 'red', size = 1.3) +
theme_classic()
Posterior predictive distribution of item means
obs_mean <- data.frame(
Item = paste0("Item_Mean_rep[", 1:6,"]"),
y = Item_Mean_obs
)
Item_Mean_rep_mat |> as.data.frame() |>
pivot_longer(everything(), names_to = "Item", values_to = "Posterior") |>
ggplot() +
geom_density(aes(x = Posterior)) +
geom_vline(aes(xintercept = y), data = obs_mean, col ='red') +
facet_wrap(~ Item, scales = "free")
MCMC Results > Visualization > Factor Loadings
MCMC Results > Visualization > Factor Correlation
- We can also visual inspect trace plots of MCMC draws of
L - What type of error leads to the deviation between posterior values of
Land true value.
MCMC Results > Estimation > Item Intercepts
Item intercepts are set to .1. Let’s see what Bayesian model recovers the true values
# A tibble: 6 × 10
variable mean median sd mad q5 q95 rhat ess_bulk ess_tail
<chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
1 mu[1] 0.106 0.106 0.0247 0.0231 0.0673 0.147 1.03 83.2 122.
2 mu[2] 0.106 0.106 0.0177 0.0167 0.0777 0.135 1.03 85.8 144.
3 mu[3] 0.101 0.101 0.0109 0.0102 0.0836 0.119 1.03 88.6 130.
4 mu[4] 0.0842 0.0844 0.0208 0.0217 0.0488 0.118 1.02 128. 328.
5 mu[5] 0.0966 0.0968 0.0151 0.0157 0.0707 0.121 1.02 129. 348.
6 mu[6] 0.0989 0.0990 0.00952 0.00983 0.0828 0.114 1.02 138. 401.We notices item 4 is little off with true values.
MCMC Results > Estimation > Residual variances
# A tibble: 6 × 10
variable mean median sd mad q5 q95 rhat ess_bulk ess_tail
<chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
1 sigma[1] 0.107 0.108 0.00688 0.00674 0.0961 0.118 1.01 747. 1009.
2 sigma[2] 0.0995 0.0995 0.00412 0.00412 0.0928 0.106 1.00 1178. 2166.
3 sigma[3] 0.0961 0.0960 0.00254 0.00251 0.0919 0.100 1.00 5741. 6069.
4 sigma[4] 0.0991 0.0993 0.00676 0.00644 0.0877 0.110 1.00 738. 1131.
5 sigma[5] 0.0962 0.0962 0.00408 0.00416 0.0894 0.103 1.00 1221. 2231.
6 sigma[6] 0.0998 0.0998 0.00256 0.00260 0.0957 0.104 1.00 6468. 5918.
Simulation Study 2
- To illustrate how to model a more complex CFA with cross-loadings, let’s generate a new simulation data with test length 7 and two latent variables
- The mean structure of model is like this, each latent variable was measured by 4 items. The other parameters are as same as Model 1.
Mean Structure of Model 2
Strategies for Model 2
The general idea is that we assign a uninformative prior distribution for “main” factor loadings, and assign informative “close-to-zero” prior distribution for “zero” factor loadings:
Reference
