-
Notifications
You must be signed in to change notification settings - Fork 0
Commit
This commit does not belong to any branch on this repository, and may belong to a fork outside of the repository.
- Loading branch information
1 parent
7d324ba
commit 9831c63
Showing
13 changed files
with
521 additions
and
132 deletions.
There are no files selected for viewing
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| Original file line number | Diff line number | Diff line change |
|---|---|---|
| @@ -0,0 +1,106 @@ | ||
| #' Compute a Jacobian Matrix of the Implied Covariance/Correlation Matrix based on a RAM model. | ||
| #' | ||
| #' It computes a symbolic Jacobian matrix of the model-implied covariance (or correlation) matrix in | ||
| #' SEM using the RAM specification. | ||
| #' | ||
| #' @param RAM A RAM object including a list of matrices of the model returned | ||
| #' from \code{\link[metaSEM]{lavaan2RAM}} | ||
| #' @param vars A vector of characters of the random variables. If the random | ||
| #' variables are not listed in `vars`, they are treated as constants. If `vars` | ||
| #' is missing, all names in `RAM` are treated as random variables. | ||
| #' @param corr Whether the model implied matrix is covariance (default) or | ||
| #' correlation structure. | ||
| #' @return A Jacobian matrix. | ||
| #' @author Mike W.-L. Cheung <mikewlcheung@@nus.edu.sg> | ||
| #' @export | ||
| #' @examples | ||
| #' | ||
| #' #### A mediation model | ||
| #' model1 <- "y ~ c*x + b*m | ||
| #' m ~ a*x | ||
| #' ## Means | ||
| #' y ~ b0*1 | ||
| #' m ~ m0*1 | ||
| #' x ~ x0*1" | ||
| #' | ||
| #' RAM1 <- metaSEM::lavaan2RAM(model1) | ||
| #' | ||
| #' ## Model-implied covariance matrix and mean structure | ||
| #' JacobianRAM(RAM1, corr=FALSE) | ||
| #' | ||
| #' ## Model-implied correlation matrix | ||
| #' JacobianRAM(RAM1, corr=TRUE) | ||
| #' | ||
| #' #### A CFA model | ||
| #' model2 <- "f =~ x1 + x2 + x3 + x4#' | ||
| #' ## Mean | ||
| #' f ~ fmean*1" | ||
| #' | ||
| #' RAM2 <- metaSEM::lavaan2RAM(model2) | ||
| #' | ||
| #' ## Model-implied covariance matrix | ||
| #' JacobianRAM(RAM2, corr=FALSE) | ||
| #' | ||
| #' ## Model-implied correlation matrix | ||
| #' JacobianRAM(RAM2, corr=TRUE) | ||
| #' | ||
| JacobianRAM <- function(RAM, vars, corr=FALSE) { | ||
| RAM.old <- RAM | ||
|
|
||
| ## Convert reserved words to random strings | ||
| RAM$A <- .convert(RAM$A, type="sympy2random") | ||
| RAM$S <- .convert(RAM$S, type="sympy2random") | ||
| RAM$M <- .convert(RAM$M, type="sympy2random") | ||
| dim(RAM$M) <- dim(RAM.old$M) | ||
| dimnames(RAM$M) <- dimnames(RAM.old$M) | ||
|
|
||
| ## Keep the random strings as sympy reserved words are not allowed in | ||
| ## Jacobian in Sympy | ||
| x <- impliedS(RAM, corr=corr, convert=FALSE) | ||
|
|
||
| vecM <- x$Mu | ||
| Sigma <- x$Sigma | ||
| labels.mu <- paste0("Mu_", colnames(vecM)) | ||
|
|
||
| p <- ncol(Sigma) | ||
| labels.sigma <- outer(seq_len(p), seq_len(p), function(x, y) paste0("Cov", x, "_", y)) | ||
|
|
||
| ## Vector of correlation matrix only | ||
| if (corr) { | ||
| vecS <- matrix(OpenMx::vechs(Sigma), ncol=1) | ||
| labels <- OpenMx::vechs(labels.sigma) | ||
| } else { | ||
| ## Vector of covariance matrix and then vector of means | ||
| vecS <- matrix(c(OpenMx::vech(Sigma), vecM), ncol=1) | ||
| labels <- c(OpenMx::vech(labels.sigma), labels.mu) | ||
| } | ||
|
|
||
| varlist <- sort(all.vars(parse(text=vecS))) | ||
| vecS <- caracas::as_sym(vecS) | ||
|
|
||
| ## Ensure the column order of the Jmatrix is sorted according to the parameters | ||
| if (missing(vars)) { | ||
| vars <- sort(varlist) | ||
| } else { | ||
| vars <- .convert(sort(vars), type="sympy2random") | ||
| if (any(!vars %in% varlist)) { | ||
| stop("Some of \"vars\" do not agree with those names in \"x\".\n") | ||
| } | ||
| } | ||
|
|
||
| ## Jmatrix <- caracas::jacobian(vecS, vars) | ||
| ## It cannot handle many vars, say 15 | ||
| ## Error in py_call_impl(callable, call_args$unnamed, call_args$named) : | ||
| ## ValueError: too many values to unpack (expected 2) | ||
|
|
||
| ## Adhoc solution to split vars to batches | ||
| vars.batch <- split(x=vars, | ||
| f=rep(1:length(vars), each=10, length.out = length(vars)) ) | ||
| J.batch <- lapply(vars.batch, function(x) caracas::jacobian(vecS, x)) | ||
| Jmatrix <- do.call(cbind, J.batch) | ||
| Jmatrix <- caracas::as_character_matrix(Jmatrix) | ||
| Jmatrix <- .convert(Jmatrix, type="random2sympy") | ||
| dimnames(Jmatrix) <- list(labels, .convert(vars, type="random2sympy")) | ||
| Jmatrix | ||
| } | ||
|
|
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| Original file line number | Diff line number | Diff line change |
|---|---|---|
| @@ -0,0 +1,35 @@ | ||
| .convert <- function(x, type=c("sympy2random", "random2sympy")) { | ||
| x.old <- x | ||
|
|
||
| ## Conversion the reserved words in sympy with some random strings | ||
| words <- matrix( c( | ||
| ## random sympy | ||
| "a4dbepz", "I", | ||
| "ajljael", "lambda", | ||
| "akr0dsf", "E", | ||
| "data_", "data.", | ||
| "jeljsxe", "oo" | ||
| ), byrow = TRUE, ncol = 2) | ||
| colnames(words) <- c("random", "sympy") | ||
|
|
||
| type <- match.arg(type) | ||
| switch (type, | ||
| sympy2random = for (i in 1:nrow(words)) { | ||
| x <- gsub(words[i, 2], words[i, 1], x, fixed=TRUE)}, | ||
| random2sympy = for (i in 1:nrow(words)) { | ||
| x <- gsub(words[i, 1], words[i, 2], x, fixed=TRUE)} | ||
| ) | ||
|
|
||
| ## Try to convert these matrices from string to numeric, if possible | ||
| if (suppressWarnings(all(!is.na(as.numeric(x))))) { | ||
| x <- as.numeric(x) | ||
| } | ||
|
|
||
| if (is.matrix(x.old)) { | ||
| x <- matrix(x, nrow=nrow(x.old), ncol=ncol(x.old), | ||
| dimnames=dimnames(x.old)) | ||
| } else if (is.matrix(x)) { | ||
| names(x) <- names(x.old) | ||
| } | ||
| x | ||
| } |
Oops, something went wrong.