Nothing
set_phantoms <- function(partable, ov.names, lv.names, ov.names.x, lv.names.x, ov.cp, lv.cp, lv.x.wish, ngroups) {
## Add phantom lvs for covariance parameters
## first: parameter matrices + indexing
partable <- lavMatrixRepresentation(partable, add.attributes = TRUE)
## for defined parameters
defpar <- which(partable$op == ":=")
if(length(defpar) > 0){
partable$mat[defpar] <- "def"
partable$row[defpar] <- 1:length(defpar)
partable$col[defpar] <- 1
partable$group[defpar] <- 1
}
## add prior column if it doesn't exist
if(is.na(match("prior", names(partable)))) partable$prior <- rep("", length(partable$id))
## exclude lv.x if we are using dmnorm/dwish:
if(lv.x.wish){
covpars <- which(partable$op == "~~" &
partable$lhs != partable$rhs &
partable$group == 1 &
!(partable$lhs %in% ov.names.x &
partable$free == 0) &
!(partable$lhs %in% lv.names.x))
} else {
if(lv.cp == "srs"){
covpars <- which(partable$op == "~~" &
partable$lhs != partable$rhs &
partable$group == 1 &
!(partable$lhs %in% ov.names.x &
partable$free == 0))
}
if(lv.cp == "fa"){
covpars <- which(partable$op == "~~" &
partable$lhs != partable$rhs &
partable$group == 1 &
!(partable$lhs %in% ov.names.x &
partable$free == 0) &
!(partable$lhs %in% lv.names &
partable$free == 0))
}
}
## remove covariances fixed to zero under fa
## if(lv.cp == "fa"){
## fcovs <- which(partable$free[covpars] == 0)
## zcovs <- which(partable$ustart[covpars][fcovs] == 0)
## covpars <- covpars[-fcovs[zcovs]]
## }
blkrow <- rep(NA, ncol(partable)) #length(partable$id))
facovs <- NULL
## Only do this if covpars exist
if(length(covpars) > 0){
## add to model matrices
## added entries in lambda vs in beta
nmvcovs <- sum(partable$lhs[covpars] %in% ov.names)
nlvcovs <- length(covpars) - nmvcovs
patts <- attributes(partable)
for(k in 1:ngroups){
if(!("lambda" %in% patts$mmNames[[k]]) & nmvcovs > 0){
lcolstart <- 0
attributes(partable)$mmNames[[k]] <- c(patts$mmNames[[k]],
"lambda")
attributes(partable)$mmRows[[k]] <- c(patts$mmRows[[k]], lambda=length(ov.names))
attributes(partable)$mmCols[[k]] <- c(patts$mmCols[[k]], lambda=nmvcovs)
} else {
lcolstart <- patts$mmCols[[k]]["lambda"]
attributes(partable)$mmCols[[k]]["lambda"] <- patts$mmCols[[k]]["lambda"] + nmvcovs
}
if(!("beta" %in% patts$mmNames[[k]]) & nlvcovs > 0){
bcolstart <- 0
attributes(partable)$mmNames[[k]] <- c(patts$mmNames[[k]], "beta")
attributes(partable)$mmRows[[k]] <- c(patts$mmRows[[k]], beta=
nlvcovs)
attributes(partable)$mmCols[[k]] <- c(patts$mmCols[[k]], beta=nlvcovs)
} else {
bcolstart <- patts$mmCols[[k]]["beta"]
attributes(partable)$mmRows[[k]]["beta"] <- patts$mmRows[[k]]["beta"] + nlvcovs
attributes(partable)$mmCols[[k]]["beta"] <- patts$mmCols[[k]]["beta"] + nlvcovs
}
if(!("psi" %in% patts$mmNames[[k]])){
psicolstart <- 0
attributes(partable)$mmNames[[k]] <- c(patts$mmNames[[k]], "psi")
attributes(partable)$mmRows[[k]] <- c(patts$mmRows[[k]], psi=length(covpars))
attributes(partable)$mmCols[[k]] <- c(patts$mmCols[[k]], psi=length(covpars))
} else {
psicolstart <- patts$mmCols[[k]]["psi"]
attributes(partable)$mmRows[[k]]["psi"] <- patts$mmRows[[k]]["psi"] + length(covpars)
attributes(partable)$mmCols[[k]]["psi"] <- patts$mmCols[[k]]["psi"] + length(covpars)
}
}
cprm <- NULL
## which covariances are under srs?
ridx <- 1:length(covpars)
for(k in 1:ngroups){
tlcs <- lcolstart
tbcs <- bcolstart
tpcs <- psicolstart
for(i in 1:length(covpars)){
## Find the row for group k (needed for multiple groups)
covparg <- which(partable$op == "~~" &
partable$lhs == partable$lhs[covpars[i]] &
partable$rhs == partable$rhs[covpars[i]] &
partable$group == k)
## Is this constrained equal to a previous parameter?
eq.const <- FALSE
grp.idx <- k
eq.idx <- which(partable$op == "==" & partable$rhs == partable$plabel[covparg])
if(length(eq.idx) > 0){
eq.const <- TRUE
## TODO? assumes it is equal to another covariance; do any models
## restrict covariances to be equal to other types of parameters?
full.idx <- which(partable$plabel == partable$lhs[eq.idx])
old.idx <- which(partable$lhs[covpars] == partable$lhs[full.idx[1]] &
partable$rhs[covpars] == partable$rhs[full.idx[1]])
old.ridx <- ridx[old.idx]
grp.idx <- partable$group[full.idx[1]]
}
tmprows <- nrow(partable) + 1:3
phname <- paste(".phant", i, sep="")
partable <- rbind(partable, blkrow, blkrow, blkrow)
## TODO? should 'block' ever differ from 'group'?
partable$group[tmprows] <- partable$block[tmprows] <- k
partable$rhs[tmprows[1]] <- partable$lhs[covpars[i]]
partable$rhs[tmprows[2]] <- partable$rhs[covpars[i]]
## Decide on =~ (ov) vs ~ (lv)
if(partable$lhs[covpars[i]] %in% ov.names){
partable$lhs[tmprows[1]] <- phname
partable$op[tmprows[1]] <- "=~"
partable$rhs[tmprows[1]] <- partable$lhs[covpars[i]]
partable$mat[tmprows[1]] <- "lambda"
partable$row[tmprows[1]] <- match(partable$lhs[covpars[i]],
patts$mmDimNames[[k]]$lambda[[1]])
tlcs <- tlcs + 1
partable$col[tmprows[1]] <- tlcs
v1var <- which(partable$lhs == partable$lhs[covpars[i]] &
partable$rhs == partable$lhs[covpars[i]] &
partable$group == k &
partable$op == "~~")
tmpv1 <- paste(partable$mat[v1var], "[", partable$row[v1var], ",", partable$col[v1var], ",", k,
"]", sep="")
if(eq.const){
oldr <- match(partable$lhs[full.idx], patts$mmDimNames[[k]]$lambda[[1]])
oldv1 <- paste(partable$mat[v1var], "[", oldr , ",", oldr, ",", grp.idx, "]", sep="")
}
ctype <- "ov"
} else {
partable$lhs[tmprows[1]] <- partable$lhs[covpars[i]]
partable$op[tmprows[1]] <- "~"
partable$rhs[tmprows[1]] <- phname
partable$mat[tmprows[1]] <- "beta"
tbcs <- tbcs + 1
partable$row[tmprows[1]] <- tbcs
partable$col[tmprows[1]] <- match(partable$lhs[covpars[i]],
patts$mmDimNames[[k]]$psi[[1]])
tmpv1 <- paste("psi[", partable$col[tmprows[1]], ",", partable$col[tmprows[1]],
",", k, "]", sep="")
if(eq.const){
oldr <- match(partable$lhs[full.idx], patts$mmDimNames[[k]]$psi[[1]])
oldv1 <- paste("psi[", oldr, ",", oldr, ",", grp.idx, "]", sep="")
}
ctype <- "lv"
}
if(partable$rhs[covpars[i]] %in% ov.names){
partable$lhs[tmprows[2]] <- phname
partable$op[tmprows[2]] <- "=~"
partable$rhs[tmprows[2]] <- partable$rhs[covpars[i]]
partable$mat[tmprows[2]] <- "lambda"
partable$row[tmprows[2]] <- match(partable$rhs[covpars[i]],
patts$mmDimNames[[k]]$lambda[[1]])
tlcs <- tlcs + 1
partable$col[tmprows[2]] <- tlcs
v2var <- which(partable$lhs == partable$rhs[covpars[i]] &
partable$rhs == partable$rhs[covpars[i]] &
partable$group == k &
partable$op == "~~")
tmpv2 <- paste(partable$mat[v2var], "[", partable$row[v2var], ",", partable$col[v2var], ",", k, "]", sep="")
if(eq.const){
oldr <- match(partable$rhs[full.idx], patts$mmDimNames[[k]]$lambda[[1]])
oldv2 <- paste(partable$mat[v2var], "[", oldr, ",", oldr, ",", grp.idx, "]", sep="")
}
} else {
partable$lhs[tmprows[2]] <- partable$rhs[covpars[i]]
partable$op[tmprows[2]] <- "~"
partable$rhs[tmprows[2]] <- phname
partable$mat[tmprows[2]] <- "beta"
tbcs <- tbcs + 1
partable$row[tmprows[2]] <- tbcs
partable$col[tmprows[2]] <- match(partable$rhs[covpars[i]],
patts$mmDimNames[[k]]$psi[[1]])
tmpv2 <- paste("psi[", partable$col[tmprows[2]], ",", partable$col[tmprows[2]],
",", k, "]", sep="")
if(eq.const){
oldr <- match(partable$lhs[full.idx], patts$mmDimNames[[k]]$psi[[1]])
oldv2 <- paste("psi[", oldr, ",", oldr, ",", grp.idx, "]", sep="")
}
}
## TODO tmprows[3] should get some prior associated with covariance param for fa priors?
partable$prior[tmprows[1:3]] <- ""
## Decide what priors to use
tpcs <- tpcs + 1
if((ctype == "ov" & ov.cp == "srs") | (ctype == "lv" & lv.cp == "srs")){
rhomat <- "rho"
if(partable$mat[covpars[i]] == "psi") rhomat <- "lvrho"
rhoind <- paste(partable$row[covpars[i]], ",", partable$col[covpars[i]], sep="")
## srs priors
partable$free[tmprows[1:3]] <- 0L
partable$exo[tmprows[1:3]] <- 0L
partable$ustart[tmprows[1]] <- paste("sqrt(abs(", rhomat, "[", rhoind, ",", k, "])*", tmpv1, ")", sep="")
partable$ustart[tmprows[2]] <- paste("(-1 + 2*step(", rhomat, "[", rhoind, ",", k,
"]))*sqrt(abs(", rhomat, "[", rhoind, ",", k, "])*", tmpv2, ")", sep="")
partable$ustart[tmprows[3]] <- 1
partable$mat[tmprows[3]] <- "psi"
partable$row[tmprows[3]] <- partable$col[tmprows[3]] <- tpcs
partable$id[covparg] <- paste(rhomat, "[", rhoind, ",", k, "]", sep="")
partable$plabel[tmprows] <- paste(".p", tmprows,
".", sep="")
partable$lhs[tmprows[3]] <- partable$rhs[tmprows[3]] <- phname
partable$op[tmprows[3]] <- "~~"
} else {
## factor analysis priors
partable$mat[tmprows[3]] <- "psi"
partable$row[tmprows[3]] <- partable$col[tmprows[3]] <- tpcs
if(partable$free[covparg] == 0){
if(partable$ustart[covparg] != 0) stop("blavaan ERROR: Cannot fix covariances to nonzero values under fa priors.\n")
partable$free[tmprows[1:3]] <- 0
partable$ustart[tmprows[3]] <- 1
partable$ustart[tmprows[1:2]] <- 0
partable$exo[tmprows[1:3]] <- 0
} else {
partable$plabel[tmprows[3]] <- partable$plabel[covparg]
partable$free[tmprows[3]] <- partable$free[covparg]
partable$exo[tmprows[1:3]] <- 0
partable$plabel[tmprows[1:2]] <- paste(".p", tmprows[1:2], ".", sep="")
partable$free[tmprows[1:2]] <- tmprows[1:2]
}
cprm <- c(cprm, covparg)
partable$lhs[tmprows[3]] <- partable$rhs[tmprows[3]] <- phname
partable$op[tmprows[3]] <- "~~"
}
## equality constraints
if(eq.const){
## set loading/regression parameters equal to others
partable$free[tmprows[1:2]] <- tmprows[1:2]
old.labels <- which(partable$op %in% c("=~", "~") &
partable$group == grp.idx &
(grepl(paste(".phant", old.ridx, sep=""),
partable$lhs) |
grepl(paste(".phant", old.ridx, sep=""),
partable$rhs)))
partable <- rbind(partable, blkrow, blkrow)
nr <- nrow(partable)
partable$lhs[(nr-1):nr] <- partable$plabel[old.labels[1:2]]
partable$op[(nr-1):nr] <- "=="
partable$rhs[(nr-1):nr] <- partable$plabel[tmprows[1:2]]
partable$mat[(nr-1):nr] <- ""
partable$user[(nr-1):nr] <- 2
partable$free[(nr-1):nr] <- 0
partable$group[(nr-1):nr] <- 0
partable$exo[(nr-1):nr] <- 0
if((ctype == "ov" & ov.cp == "fa") | (ctype == "lv" & lv.cp == "fa")){
old.label <- which(partable$op == "~~" &
partable$group == grp.idx &
grepl(paste(".phant", old.ridx, sep=""), partable$lhs))
partable <- rbind(partable, blkrow)
nr <- nrow(partable)
partable$lhs[nr] <- partable$plabel[old.label]
partable$op[nr] <- "=="
partable$rhs[nr] <- partable$plabel[tmprows[3]]
partable$mat[(nr-1):nr] <- ""
partable$user[nr] <- 2
partable$free[nr] <- 0
partable$group[nr] <- 0
partable$exo[nr] <- 0
}
}
}
}
## Now remove rows of fa covariance parameters
if(!is.null(cprm)){
facovs <- partable[cprm,]
partable <- partable[-cprm,]
}
}
## FIXME?
## Remove covariances associated with fixed x
covpars <- which(partable$op == "~~" &
partable$lhs != partable$rhs &
partable$group == 1 &
partable$lhs %in% ov.names.x &
partable$free == 0)
if(length(covpars) > 0) partable <- partable[-covpars,]
## to avoid errors in lavaan 0.6-1
if(any(is.na(partable$exo))) partable$exo[is.na(partable$exo)] <- 0
## Add parameter numbers now that we have phantoms
parnums <- rep(NA, nrow(partable))
parrows <- which(!(partable$op == "=="))
parnums[parrows] <- 1:length(parrows)
partable$parnums <- parnums
list(partable = partable, facovs = facovs)
}
set_mv0 <- function(partable, ov.names, ngroups) {
## If any mvs have fixed 0 variance (single indicator lv),
## move the fixed 0 variance to the lv
mv0 <- which(partable$op == "~~" &
partable$lhs %in% ov.names &
partable$rhs == partable$lhs &
partable$group == 1 &
partable$free == 0 &
partable$ustart == 0)
if(length(mv0) > 0){
ovn <- partable$lhs[mv0]
for(i in 1:length(ovn)){
for(j in 1:ngroups){
mvloc <- which(partable$op == "~~" &
partable$lhs == ovn[i] &
partable$rhs == partable$lhs &
partable$group == j &
partable$free == 0 &
partable$ustart == 0)
lvloc <- which(partable$op == "=~" &
partable$rhs == ovn[i] &
partable$group == j)
lvreg <- which(partable$op == "~" &
partable$lhs %in% c(ovn[i], partable$lhs[lvloc]) &
partable$group == j)
lvcov <- which(partable$op == "~~" &
(partable$lhs %in% partable$lhs[lvloc] |
partable$rhs %in% partable$lhs[lvloc]) &
partable$lhs != partable$rhs)
## If this is an indicator of multiple lvs or
## has covariances attached, we cannot handle it
## in conditional (on lv) form:
if(length(lvloc) + length(lvreg) + length(lvcov) > 1){
if(length(mvloc) > 1){
stop("blavaan ERROR: Problem with ov variances fixed to 0.")
}
partable$ustart[mvloc] <- .001
message(paste("blavaan NOTE: The variance of variable", ovn[i],
"in group", j, "has been fixed to .001 instead of 0 (necessary for conditional model specification).\n"))
} else {
lvname <- partable$lhs[lvloc]
## TODO? check for covariances with the lv?
lvvar <- which(partable$lhs == lvname &
partable$rhs == lvname &
partable$op == "~~" &
partable$group == j)
tmpfree <- partable$free[lvvar]
tmpustart <- partable$ustart[lvvar]
tmpplabel <- partable$plabel[lvvar]
tmpstart <- partable$start[lvvar]
partable$free[lvvar] <- partable$free[mvloc]
partable$ustart[lvvar] <- partable$ustart[mvloc]
partable$plabel[lvvar] <- partable$plabel[mvloc]
partable$start[lvvar] <- partable$start[mvloc]
partable$free[mvloc] <- tmpfree
partable$ustart[mvloc] <- tmpustart
partable$plabel[mvloc] <- tmpplabel
partable$start[mvloc] <- tmpstart
}
}
}
}
partable
}
set_phanvars <- function(partable, ov.names, lv.names, ov.cp, lv.cp, ngroups){
## once we have defined all phantoms, go back to set equality
## constraints on phantom variances
vnames <- c(ov.names, lv.names)
for(i in 1:length(vnames)){
for(k in 1:ngroups){
if(vnames[i] %in% ov.names){
phanlvs <- which(partable$rhs == vnames[i] &
grepl(".phant", partable$lhs) &
partable$op == "=~" &
partable$group == k)
} else {
phanlvs <- which(partable$lhs == vnames[i] &
grepl(".phant", partable$rhs) &
partable$op == "~" &
partable$group == k)
}
if(length(phanlvs) > 0){
vvar <- which(partable$lhs == vnames[i] &
partable$lhs == partable$rhs &
partable$op == "~~" &
partable$group == k)
if(ov.cp == "srs"){
eqconst <- paste(partable$mat[vvar], "[", partable$row[vvar],
",", partable$col[vvar], ",", k, "]", sep="")
for(j in 1:length(phanlvs)){
eqconst <- paste(eqconst, " - (", partable$ustart[phanlvs[j]],
")^2", sep="")
}
partable <- rbind(partable, partable[vvar,])
partable$parnums[nrow(partable)] <- max(partable$parnums, na.rm=TRUE) + 1
partable$mat[vvar] <- paste(partable$mat[vvar], "star", sep="")
} else {
## fa priors TODO change around so it is faster
eqconst <- paste(partable$mat[vvar], "star[", partable$row[vvar],
",", partable$col[vvar], ",", k, "]", sep="")
for(j in 1:length(phanlvs)){
if(vnames[i] %in% ov.names){
phname <- partable$lhs[phanlvs[j]]
} else {
phname <- partable$rhs[phanlvs[j]]
}
phanvar <- which(partable$lhs == phname &
partable$lhs == partable$rhs &
partable$op == "~~" &
partable$group == k)
eqconst <- paste(eqconst, " + (", partable$mat[phanlvs[j]], "[",
partable$row[phanlvs[j]], ",",
partable$col[phanlvs[j]], ",", k, "]^2*",
partable$mat[phanvar], "[",
partable$row[phanvar], ",", partable$col[phanvar],
",", k, "])", sep="")
}
partable <- rbind(partable, partable[vvar,])
partable$mat[nrow(partable)] <- paste(partable$mat[vvar], "star",
sep="")
partable$parnums[nrow(partable)] <- max(partable$parnums, na.rm=TRUE) + 1
}
partable$free[vvar] <- 0
partable$label[vvar] <- ""
partable$ustart[vvar] <- eqconst
}
}
}
partable
}
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.