Harmanjit,
Are you using an LLM to generate R code and then ask other people to validate it for you? (The code looks a lot like R code that LLMs generate)
The answer to your question is: No, this would not work. The code takes a latent variable model and then calculates the interaction term using scale scores. This approach produces biased estimates of the interaction effect. It also probably does not make sense to treat occupation as a continuous variable in the model.
This paper explains how you can estimate interaction models using semTools:
Schoemann, A. M., & Jorgensen, T. D. (2021). Testing and Interpreting Latent Variable Interactions Using the semTools Package. Psych, 3(3), Article 3. https://doi.org/10.3390/psych3030024
Best regards,
Mikko
From: lav...@googlegroups.com <lav...@googlegroups.com> on behalf of Harmanjit Singh <harman...@gmail.com>
Date: Saturday, 14. March 2026 at 11.34
To: lavaan <lav...@googlegroups.com>
Subject: Is following method of testing SEM and Moderation Suitable in Lavaan?
Dear all,
If following approach for testing SEM and Moderation Suitable?
The main doubt is, should moderator variable be added in the structural model part, considering it is being tested in Parts 3 and 4 of the code?
If yes why, if not why?
# ==============================================================================
# PART 1: SETUP & LIBRARIES
# ==============================================================================
library(lavaan)
library(semTools)
library(readxl)
# Load Data
Fulldata <- read_excel("your_data_path.xlsx")
# --- DATA CLEANING ---
Fulldata <- na.omit(Fulldata)
# Convert Demographics to Numeric (Generic Labels)
Fulldata$C1 <- as.numeric(as.factor(Fulldata$Gender)) - 1
Fulldata$C2 <- as.numeric(as.factor(Fulldata$Occupation))
Fulldata$C3 <- as.numeric(as.factor(Fulldata$Education))
Fulldata$C4 <- as.numeric(as.factor(Fulldata$Income))
Fulldata$C5 <- as.numeric(Fulldata$Age)
# ==============================================================================
# PART 2: THE MAIN EFFECTS MODEL
# ==============================================================================
Model_Main <- '
# Measurement Model
IV_1 =~ AB1 + AB2 + AB3
IV_2 =~ CD1 + CD2 + CD3
IV_3 =~ EF1 + EF2 + EF3
IV_4 =~ GH1 + GH2 + GH3
IV_5 =~ IJ1 + IJ2 + IJ3
IV_6 =~ KL1 + KL2 + KL3
Moderator =~ MN1 + MN2 + MN3
Dependent_Var =~ OP1 + OP2 + OP3
# Structural Model
Dependent_Var ~ IV_1 + IV_2 + IV_3 + IV_4 +
IV_5 + IV_6 + Moderator +
C1 + C2 + C3 + C4 + C5
'
# Using MLM for robust estimation
fit_main <- sem(Model_Main, estimator="MLM", data=Fulldata, meanstructure=TRUE)
print("--- MAIN EFFECTS MODEL RESULTS ---")
summary(fit_main, fit.measures=TRUE, standardized=TRUE, rsquare=TRUE)
print(fitmeasures(fit_main, c("chisq.scaled","df.scaled","pvalue.scaled","tli.robust","cfi.robust","rmsea.robust","srmr")))
# ==============================================================================
# PART 3: AUTOMATED MODERATION ANALYSIS
# ==============================================================================
# 1. VARIABLE DICTIONARY
item_dict <- list(
IV_1 = paste0("AB", 1:3),
IV_2 = paste0("CD", 1:3),
IV_3 = paste0("EF", 1:3),
IV_4 = paste0("GH", 1:3),
IV_5 = paste0("IJ", 1:3),
IV_6 = paste0("KL", 1:3),
Moderator = paste0("MN", 1:3),
Dependent_Var = paste0("OP", 1:3)
)
control_vars <- "C1 + C2 + C3 + C4 + C5"
# 2. MODERATION FUNCTION
run_moderation_robust <- function(df, iv_name, mod_name, dv_name) {
iv_items <- item_dict[[iv_name]]
mod_items <- item_dict[[mod_name]]
dv_items <- item_dict[[dv_name]]
# Calculate Mean Scores and Center for Interaction
iv_score <- rowMeans(df[, iv_items], na.rm = TRUE)
mod_score <- rowMeans(df[, mod_items], na.rm = TRUE)
iv_centered <- scale(iv_score, center=TRUE, scale=FALSE)
mod_centered <- scale(mod_score, center=TRUE, scale=FALSE)
# Add Interaction Term to df
df$Interaction_Observed <- as.numeric(iv_centered * mod_centered)
model_syntax <- paste0("
# Measurement Models
", iv_name, " =~ ", paste(iv_items, collapse=" + "), "
", mod_name, " =~ ", paste(mod_items, collapse=" + "), "
", dv_name, " =~ ", paste(dv_items, collapse=" + "), "
# Structural Regression
", dv_name, " ~ ", iv_name, " + ", mod_name, " + Interaction_Observed + ", control_vars, "
# Covariances
Interaction_Observed ~~ ", iv_name, "
Interaction_Observed ~~ ", mod_name, "
")
fit <- sem(model_syntax, data = df, std.lv = TRUE, estimator="MLM")
return(fit)
}
# ==============================================================================
# PART 4: EXECUTION LOOP
# ==============================================================================
IV_List <- c("IV_1", "IV_2", "IV_3", "IV_4", "IV_5", "IV_6")
Mod_List <- c("Moderator")
DV_List <- c("Dependent_Var")
results_table <- data.frame()
for (mod in Mod_List) {
for (dv in DV_List) {
for (iv in IV_List) {
tryCatch({
fit <- run_moderation_robust(Fulldata, iv, mod, dv)
if (inspect(fit, "converged")) {
params <- parameterEstimates(fit, standardized = TRUE)
interaction_row <- subset(params, op == "~" & rhs == "Interaction_Observed")
if(nrow(interaction_row) > 0) {
new_row <- data.frame(
Independent_Var = iv,
Moderator = mod,
Dependent_Var = dv,
Interaction_Est = round(interaction_row$est, 3),
P_Value = round(interaction_row$pvalue, 3),
Std_Beta = round(interaction_row$std.all, 3),
Significant = ifelse(interaction_row$pvalue < 0.05, "YES", "NO")
)
results_table <- rbind(results_table, new_row)
}
}
}, error = function(e) {
cat(paste("Error in model for IV:", iv, "\n"))
})
}
}
}
# ==============================================================================
# PART 5: VIEW RESULTS
# ==============================================================================
print("--- FINAL MODERATION RESULTS ---")
print(results_table)
--
You received this message because you are subscribed to the Google Groups "lavaan" group.
To unsubscribe from this group and stop receiving emails from it, send an email to
lavaan+un...@googlegroups.com.
To view this discussion visit
https://groups.google.com/d/msgid/lavaan/d79a7324-8b1f-4edb-8964-c385afcccb45n%40googlegroups.com.