diff --git a/DESCRIPTION b/DESCRIPTION index cb12b1c3..ccfcc6c0 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,6 +1,6 @@ Package: lavaan Title: Latent Variable Analysis -Version: 0.6-19.2169 +Version: 0.6-19.2170 Authors@R: c(person(given = "Yves", family = "Rosseel", role = c("aut", "cre"), email = "Yves.Rosseel@UGent.be", diff --git a/R/ctr_pml_utils.R b/R/ctr_pml_utils.R index 6642460c..fd812cbb 100644 --- a/R/ctr_pml_utils.R +++ b/R/ctr_pml_utils.R @@ -98,6 +98,8 @@ uni_lik <- function(Y1, th.y1, eXo = NULL, PI.y1 = NULL) { uni_lik <- pnorm(th.y1.upper) - pnorm(th.y1.lower) uni_lik[is.na(uni_lik)] <- 0 + + uni_lik } ################################################################# diff --git a/R/lav_model_compute.R b/R/lav_model_compute.R index bada097f..1906763e 100644 --- a/R/lav_model_compute.R +++ b/R/lav_model_compute.R @@ -735,7 +735,6 @@ computeNU <- function(lavmodel = NULL, GLIST = NULL, NU } - # E(Y): expectation (mean) of observed variables # returns vector 1 x nvar computeEY <- function(lavmodel = NULL, GLIST = NULL, lavsamplestats = NULL, @@ -778,6 +777,48 @@ computeEY <- function(lavmodel = NULL, GLIST = NULL, lavsamplestats = NULL, EY } +# E(Y|x_i): conditional expectation (mean) of observed variables +# returns matrix N x nvar +computeEYx <- function(lavmodel = NULL, GLIST = NULL, lavsamplestats = NULL, + eXo = NULL, delta = TRUE) { + # state or final? + if (is.null(GLIST)) GLIST <- lavmodel@GLIST + + nblocks <- lavmodel@nblocks + nmat <- lavmodel@nmat + representation <- lavmodel@representation + + # return a list + EYx <- vector("list", length = nblocks) + + # compute E(Y) for each group + for (g in 1:nblocks) { + # which mm belong to group g? + mm.in.group <- 1:nmat[g] + cumsum(c(0, nmat))[g] + MLIST <- GLIST[mm.in.group] + + if (representation == "LISREL") { + EYx.g <- computeEYx.LISREL( + MLIST = MLIST, + eXo = eXo[[g]], + sample.mean = lavsamplestats@mean[[g]], + ov.y.dummy.ov.idx = lavmodel@ov.y.dummy.ov.idx[[g]], + ov.x.dummy.ov.idx = lavmodel@ov.x.dummy.ov.idx[[g]], + ov.y.dummy.lv.idx = lavmodel@ov.y.dummy.lv.idx[[g]], + ov.x.dummy.lv.idx = lavmodel@ov.x.dummy.lv.idx[[g]], + delta = delta + ) + } else { + lav_msg_stop(gettext( + "only representation LISREL has been implemented for now")) + } + + EYx[[g]] <- EYx.g + } + + EYx +} + # E(Y | ETA, x_i): conditional expectation (means) of observed variables # for a given value of x_i AND eta_i diff --git a/R/lav_objective.R b/R/lav_objective.R index 93045d08..7b6c0758 100644 --- a/R/lav_objective.R +++ b/R/lav_objective.R @@ -552,7 +552,7 @@ estimator.PML <- function(Sigma.hat = NULL, # model-based var/cov/cor pstar.idx <- PSTAR[i, j] # cat("pstar.idx =", pstar.idx, "i = ", i, " j = ", j, "\n") if (ov.types[i] == "numeric" && - ov.types[j] == "numeric") { + ov.types[j] == "numeric") { # ordinary pearson correlation LIK[, pstar.idx] <- lav_bvreg_lik( @@ -565,7 +565,7 @@ estimator.PML <- function(Sigma.hat = NULL, # model-based var/cov/cor rho = Cor.hat[i, j] ) } else if (ov.types[i] == "numeric" && - ov.types[j] == "ordered") { + ov.types[j] == "ordered") { # polyserial correlation ### FIXME: th.y2 should go into ps_lik!!! LIK[, pstar.idx] <- @@ -579,7 +579,7 @@ estimator.PML <- function(Sigma.hat = NULL, # model-based var/cov/cor rho = Cor.hat[i, j] ) } else if (ov.types[j] == "numeric" && - ov.types[i] == "ordered") { + ov.types[i] == "ordered") { # polyserial correlation ### FIXME: th.y1 should go into ps_lik!!! LIK[, pstar.idx] <- @@ -593,7 +593,7 @@ estimator.PML <- function(Sigma.hat = NULL, # model-based var/cov/cor rho = Cor.hat[i, j] ) } else if (ov.types[i] == "ordered" && - ov.types[j] == "ordered") { + ov.types[j] == "ordered") { LIK[, pstar.idx] <- pc_lik_PL_with_cov( Y1 = X[, i], diff --git a/R/lav_sam_step1.R b/R/lav_sam_step1.R index 60e6b092..f5df5148 100644 --- a/R/lav_sam_step1.R +++ b/R/lav_sam_step1.R @@ -482,7 +482,11 @@ lav_sam_step1_local <- function(STEP1 = NULL, FIT = NULL, FS[[b]] <- do.call("cbind", tmp) colnames(FS[[b]]) <- LABEL - # dummy lv's? + # dummy lv's? (both 'x' and 'y'!) + dummy.ov.idx <- c(FIT@Model@ov.y.dummy.ov.idx[[b]], + FIT@Model@ov.x.dummy.ov.idx[[b]]) + dummy.lv.idx <- c(FIT@Model@ov.y.dummy.lv.idx[[b]], + FIT@Model@ov.x.dummy.lv.idx[[b]]) if (length(dummy.lv.idx) > 0L) { FS.obs <- FIT@Data@X[[b]][, dummy.ov.idx, drop = FALSE] colnames(FS.obs) <- FIT@Data@ov.names[[b]][dummy.ov.idx]