ecmerkle / blavaan

An R package for Bayesian structural equation modeling
https://ecmerkle.github.io/blavaan
86 stars 23 forks source link

how to set latent mean equal to zero in one group? #92

Open jpritikin opened 1 day ago

jpritikin commented 1 day ago

lavaan 0.6-19 (cran), blavaan 0.5-6.1311 (via github)

The line phenom ~ c(0,pmM)*1 seems to not work correctly? The intercept is set to pmM in both groups.

library(rstan)
options(mc.cores = parallel::detectCores()/2)
library(blavaan)

future::plan("multicore")

NoYesLevels <- c('No','Not sure','Yes')
EuphoricLevels <- c('Euphoric', 'Mildly euphoric', 'No change', 'Mild discomfort', 'Severe discomfort')
ImpactLevels <- c('Helped','Slightly helped','No impact','Slightly hindered','Hindered')
AbsItemNames <- c('boredom','frustration','anxiety', 'abort','bodyAcute','bodyAfter','sleep')

ImportStudyItems <- function(df) {
  if (any(is.na(match(AbsItemNames, colnames(df))))) stop("Items are mising")

  df[['boredom']] <- ordered(df[['boredom']], levels=NoYesLevels)
  df[['frustration']] <- ordered(df[['frustration']], levels=NoYesLevels)
  df[['anxiety']] <- ordered(df[['anxiety']], levels=NoYesLevels)
  df[['abort']] <- ordered(df[['abort']], levels=NoYesLevels)
  df[['bodyAcute']] <- ordered(df[['bodyAcute']], levels=EuphoricLevels)
  df[['bodyAfter']] <- ordered(df[['bodyAfter']], levels=EuphoricLevels)
  df[['sleep']] <- ordered(df[['sleep']], levels=ImpactLevels)
  df
}

psiData <- ImportStudyItems(read.csv("psi.csv"))
pmData <- ImportStudyItems(read.csv("pm.csv"))
combinedData <- rbind(cbind(psiData, group="psi"),
                      cbind(pmData, group="pm"))
combinedData$group <- factor(combinedData$group, levels=c('psi','pm'))

spec <- '
phenom =~ boredomL*boredom + frustrationL*frustration + anxietyL*anxiety + abortL*abort + bodyAcuteL*bodyAcute + bodyAfterL*bodyAfter + sleepL*sleep

phenom ~ c(0,pmM)*1
boredom ~ boredomM*1
frustration ~ frustrationM*1
anxiety ~ anxietyM*1
abort ~ abortM*1
bodyAcute ~ bodyAcuteM*1
bodyAfter ~ bodyAfterM*1
sleep ~ sleepM*1

phenom ~~ 1*phenom
boredom ~~ 1*boredom
frustration ~~ 1*frustration
anxiety ~~ 1*anxiety
abort ~~ 1*abort
bodyAcute ~~ 1*bodyAcute
bodyAfter ~~ 1*bodyAfter
sleep ~~ 1*sleep

boredom | -.431*t1 + .431*t2
frustration | -.431*t1 + .431*t2
anxiety | -.431*t1 + .431*t2
abort | -.431*t1 + .431*t2
bodyAcute | -.842*t1 + -.253*t2 + .253*t3 + .842*t4
bodyAfter | -.842*t1 + -.253*t2 + .253*t3 + .842*t4
sleep | -.842*t1 + -.253*t2 + .253*t3 + .842*t4
'

fit <- blavaan(spec, data=combinedData,
               n.chains = 4,
               ordered = AbsItemNames,
               dp = dpriors(nu="normal(0,1)", alpha="normal(0,1)", lambda="normal(0,1)"),
               group = "group",
               allow.empty.cell = TRUE, test = "none")

summary(fit)

Data:

pm.csv psi.csv

jpritikin commented 1 day ago

If I use phenom ~ c(0,NA)*1 then does what is expected. The parameter is fixed in group 1 and free in group 2.

ecmerkle commented 1 hour ago

Thanks for the reports, and I will take a look at the others soon. For this one, I believe that the lavaan parser does not like parameter labels mixed with numerical constraints inside c(). You can see this by running

fit2 <- lavaan(spec, data=combinedData, ordered=AbsItemNames, group="group", allow.empty.cell=TRUE, cat.wls.w=FALSE)
lavInspect(fit2, 'free')

It produces a bunch of warnings, but the nonzero entries in the output of the second command signify free parameters. The $psi$alpha and $pm$alpha matrices show whether or not your phenom means are free.

If you want parameter labels and one mean fixed to 0, I think the line below would work. Though then a different warning comes from lavaan, that I believe can be ignored:

phenom ~ c(0,NA)*1 + c(psM,pmM)*1