From: Benjamin Auder Date: Fri, 21 Apr 2017 15:14:10 +0000 (+0200) Subject: Merge branch 'master' of auder.net:valse X-Git-Url: https://git.auder.net/?p=valse.git;a=commitdiff_plain;h=c3bf2821bce67c75504e303fae23dd41c00f06c8;hp=465c0e07fbf8363a2625f26681fa7ba31bad82b7 Merge branch 'master' of auder.net:valse --- diff --git a/pkg/R/constructionModelesLassoMLE.R b/pkg/R/constructionModelesLassoMLE.R index c55035e..78e11ad 100644 --- a/pkg/R/constructionModelesLassoMLE.R +++ b/pkg/R/constructionModelesLassoMLE.R @@ -40,10 +40,10 @@ constructionModelesLassoMLE <- function(phiInit, rhoInit, piInit, gamInit, mini, if (verbose) print(paste("Computations for lambda=", lambda)) - n <- dim(X)[1] - p <- dim(phiInit)[1] - m <- dim(phiInit)[2] - k <- dim(phiInit)[3] + n <- nrow(X) + p <- ncol(X) + m <- ncol(Y) + k <- length(piInit) sel.lambda <- S[[lambda]]$selected # col.sel = which(colSums(sel.lambda)!=0) #if boolean matrix col.sel <- which(sapply(sel.lambda, length) > 0) #if list of selected vars @@ -51,8 +51,8 @@ constructionModelesLassoMLE <- function(phiInit, rhoInit, piInit, gamInit, mini, return(NULL) # lambda == 0 because we compute the EMV: no penalization here - res <- EMGLLF(array(phiInit[col.sel, , ],dim=c(length(col.sel),m,k)), rhoInit, - piInit, gamInit, mini, maxi, gamma, 0, as.matrix(X[, col.sel]), Y, eps, fast) + res <- EMGLLF(array(phiInit,dim=c(p,m,k))[col.sel, , ], rhoInit, piInit, gamInit, + mini, maxi, gamma, 0, as.matrix(X[, col.sel]), Y, eps, fast) # Eval dimension from the result + selected phiLambda2 <- res$phi @@ -65,7 +65,7 @@ constructionModelesLassoMLE <- function(phiInit, rhoInit, piInit, gamInit, mini, ## Computation of the loglikelihood # Precompute det(rhoLambda[,,r]) for r in 1...k - detRho <- sapply(1:k, function(r) det(rhoLambda[, , r])) + detRho <- sapply(1:k, function(r) gdet(rhoLambda[, , r])) sumLogLLH <- 0 for (i in 1:n) { @@ -81,17 +81,6 @@ constructionModelesLassoMLE <- function(phiInit, rhoInit, piInit, gamInit, mini, sumLogLLH <- sumLogLLH + log(norm_fact) - log((2 * base::pi)^(m/2)) } llhLambda <- c(sumLogLLH/n, (dimension + m + 1) * k - 1) - # densite <- vector("double", n) - # for (r in 1:k) - # { - # if (length(col.sel) == 1) - # { - # delta <- (Y %*% rhoLambda[, , r] - (X[, col.sel] %*% t(phiLambda[col.sel, , r]))) - # } else delta <- (Y %*% rhoLambda[, , r] - (X[, col.sel] %*% phiLambda[col.sel, , r])) - # densite <- densite + piLambda[r] * det(rhoLambda[, , r])/(sqrt(2 * base::pi))^m * - # exp(-rowSums(delta^2)/2) - # } - # llhLambda <- c(mean(log(densite)), (dimension + m + 1) * k - 1) list(phi = phiLambda, rho = rhoLambda, pi = piLambda, llh = llhLambda) } diff --git a/pkg/R/main.R b/pkg/R/main.R index d710b7e..0a3487b 100644 --- a/pkg/R/main.R +++ b/pkg/R/main.R @@ -17,8 +17,7 @@ #' @param ncores_outer Number of cores for the outer loop on k #' @param ncores_inner Number of cores for the inner loop on lambda #' @param thresh real, threshold to say a variable is relevant, by default = 1e-8 -#' @param compute_grid_lambda, TRUE to compute the grid, FALSE if known (in arguments) -#' @param grid_lambda, a vector with regularization parameters if known, by default 0 +#' @param grid_lambda, a vector with regularization parameters if known, by default numeric(0) #' @param size_coll_mod (Maximum) size of a collection of models #' @param fast TRUE to use compiled C code, FALSE for R code only #' @param verbose TRUE to show some execution traces @@ -30,8 +29,8 @@ #' @export valse <- function(X, Y, procedure = "LassoMLE", selecMod = "DDSE", gamma = 1, mini = 10, maxi = 50, eps = 1e-04, kmin = 2, kmax = 3, rank.min = 1, rank.max = 5, ncores_outer = 1, - ncores_inner = 1, thresh = 1e-08, compute_grid_lambda = TRUE, grid_lambda = 0, size_coll_mod = 10, fast = TRUE, verbose = FALSE, - plot = TRUE) + ncores_inner = 1, thresh = 1e-08, grid_lambda = numeric(0), size_coll_mod = 10, + fast = TRUE, verbose = FALSE, plot = TRUE) { n <- nrow(X) p <- ncol(X) @@ -60,7 +59,7 @@ valse <- function(X, Y, procedure = "LassoMLE", selecMod = "DDSE", gamma = 1, mi # component, doing this 20 times, and keeping the values maximizing the # likelihood after 10 iterations of the EM algorithm. P <- initSmallEM(k, X, Y, fast) - if (compute_grid_lambda == TRUE) + if (length(grid_lambda) == 0) { grid_lambda <- computeGridLambda(P$phiInit, P$rhoInit, P$piInit, P$gamInit, X, Y, gamma, mini, maxi, eps, fast) @@ -125,9 +124,7 @@ valse <- function(X, Y, procedure = "LassoMLE", selecMod = "DDSE", gamma = 1, mi })) tableauRecap <- tableauRecap[which(tableauRecap[, 4] != Inf), ] if (verbose == TRUE) - { print(tableauRecap) - } modSel <- capushe::capushe(tableauRecap, n) indModSel <- if (selecMod == "DDSE") as.numeric(modSel@DDSE@model) else if (selecMod == "Djump")