Commit | Line | Data |
---|---|---|
2279a641 BA |
1 | #' constructionModelesLassoMLE |
2 | #' | |
5965d116 | 3 | #' Construct a collection of models with the Lasso-MLE procedure. |
4 | #' | |
2279a641 BA |
5 | #' |
6 | #' @param ... | |
7 | #' | |
8 | #' @return ... | |
9 | #' | |
10 | #' export | |
11 | constructionModelesLassoMLE = function(phiInit, rhoInit, piInit, gamInit, mini, maxi, | |
bb11d873 | 12 | gamma, X, Y, thresh, tau, S, ncores=3, fast=TRUE, verbose=FALSE) |
46a2e676 | 13 | { |
08f4604c BA |
14 | if (ncores > 1) |
15 | { | |
b9b0b72a | 16 | cl = parallel::makeCluster(ncores, outfile='') |
08f4604c BA |
17 | parallel::clusterExport( cl, envir=environment(), |
18 | varlist=c("phiInit","rhoInit","gamInit","mini","maxi","gamma","X","Y","thresh", | |
19 | "tau","S","ncores","verbose") ) | |
20 | } | |
21 | ||
22 | # Individual model computation | |
23 | computeAtLambda <- function(lambda) | |
24 | { | |
25 | if (ncores > 1) | |
26 | require("valse") #nodes start with an empty environment | |
27 | ||
28 | if (verbose) | |
29 | print(paste("Computations for lambda=",lambda)) | |
30 | ||
31 | n = dim(X)[1] | |
32 | p = dim(phiInit)[1] | |
33 | m = dim(phiInit)[2] | |
34 | k = dim(phiInit)[3] | |
08f4604c BA |
35 | sel.lambda = S[[lambda]]$selected |
36 | # col.sel = which(colSums(sel.lambda)!=0) #if boolean matrix | |
37 | col.sel <- which( sapply(sel.lambda,length) > 0 ) #if list of selected vars | |
08f4604c BA |
38 | if (length(col.sel) == 0) |
39 | return (NULL) | |
40 | ||
41 | # lambda == 0 because we compute the EMV: no penalization here | |
42 | res = EMGLLF(phiInit[col.sel,,],rhoInit,piInit,gamInit,mini,maxi,gamma,0, | |
aa480ac1 | 43 | X[,col.sel], Y, tau, fast) |
08f4604c BA |
44 | |
45 | # Eval dimension from the result + selected | |
46 | phiLambda2 = res$phi | |
47 | rhoLambda = res$rho | |
48 | piLambda = res$pi | |
49 | phiLambda = array(0, dim = c(p,m,k)) | |
50 | for (j in seq_along(col.sel)) | |
fb6e49cb | 51 | phiLambda[col.sel[j],sel.lambda[[j]],] = phiLambda2[j,sel.lambda[[j]],] |
08f4604c BA |
52 | dimension = length(unlist(sel.lambda)) |
53 | ||
54 | # Computation of the loglikelihood | |
55 | densite = vector("double",n) | |
56 | for (r in 1:k) | |
57 | { | |
fb6e49cb | 58 | if (length(col.sel)==1){ |
59 | delta = (Y%*%rhoLambda[,,r] - (X[, col.sel]%*%t(phiLambda[col.sel,,r]))) | |
60 | } else delta = (Y%*%rhoLambda[,,r] - (X[, col.sel]%*%phiLambda[col.sel,,r])) | |
08f4604c | 61 | densite = densite + piLambda[r] * |
bb11d873 | 62 | det(rhoLambda[,,r])/(sqrt(2*base::pi))^m * exp(-diag(tcrossprod(delta))/2.0) |
08f4604c | 63 | } |
bb11d873 | 64 | llhLambda = c( sum(log(densite)), (dimension+m+1)*k-1 ) |
08f4604c BA |
65 | list("phi"= phiLambda, "rho"= rhoLambda, "pi"= piLambda, "llh" = llhLambda) |
66 | } | |
67 | ||
68 | # For each lambda, computation of the parameters | |
69 | out = | |
70 | if (ncores > 1) | |
71 | parLapply(cl, 1:length(S), computeAtLambda) | |
b9b0b72a BA |
72 | else |
73 | lapply(1:length(S), computeAtLambda) | |
08f4604c BA |
74 | |
75 | if (ncores > 1) | |
76 | parallel::stopCluster(cl) | |
77 | ||
78 | out | |
c3bc4705 | 79 | } |