X-Git-Url: https://git.auder.net/?a=blobdiff_plain;f=pkg%2FR%2FcomputeMu.R;h=f961fdebe40bb81fca5409ae1c4cc8aeb9099385;hb=4b2f17bb108bab0f263619cfe00eabfb1e9b8860;hp=f7f82ad20a3565522f3e6f69fb1d3da90f7de911;hpb=6dd5c2acccd10635449230faa824b7e8906911bf;p=morpheus.git diff --git a/pkg/R/computeMu.R b/pkg/R/computeMu.R index f7f82ad..f961fde 100644 --- a/pkg/R/computeMu.R +++ b/pkg/R/computeMu.R @@ -23,6 +23,7 @@ #' @examples #' io = generateSampleIO(10000, 1/2, matrix(c(1,0,0,1),ncol=2), c(0,0), "probit") #' μ = computeMu(io$X, io$Y, list(K=2)) #or just X and Y for estimated K +#' #' @export computeMu = function(X, Y, optargs=list()) { @@ -45,6 +46,8 @@ computeMu = function(X, Y, optargs=list()) large_ratio <- ( abs(Σ[-d] / Σ[-1]) > 3 ) K <- if (any(large_ratio)) max(2, which.min(large_ratio)) else d } + else if (K > d) + stop("K: integer >= 2, <= d") # Step 1: generate a family of d matrices to joint-diagonalize to increase robustness d = ncol(X) @@ -66,9 +69,8 @@ computeMu = function(X, Y, optargs=list()) jd_method = ifelse(!is.null(optargs$jd_method), optargs$jd_method, "uwedge") V = if (jd_nvects > 1) { - #NOTE: increasing itermax does not help to converge, thus we suppress warnings + # NOTE: increasing itermax does not help to converge, thus we suppress warnings suppressWarnings({jd = jointDiag::ajd(M2_t, method=jd_method)}) -# if (jd_method=="uwedge") jd$B else solve(jd$A) if (jd_method=="uwedge") jd$B else MASS::ginv(jd$A) } else @@ -79,7 +81,6 @@ computeMu = function(X, Y, optargs=list()) for (i in seq_len(K)) M2_t[,,i] = .T_I_I_w(M[[3]],V[,i]) suppressWarnings({jd = jointDiag::ajd(M2_t, method=jd_method)}) -# U = if (jd_method=="uwedge") solve(jd$B) else jd$A U = if (jd_method=="uwedge") MASS::ginv(jd$B) else jd$A μ = normalize(U[,1:K])