Commit | Line | Data |
---|---|---|
2279a641 BA |
1 | #' constructionModelesLassoRank |
2 | #' | |
43d76c49 | 3 | #' Construct a collection of models with the Lasso-Rank procedure. |
4 | #' | |
5 | #' @param S output of selectVariables.R | |
6 | #' @param k number of components | |
7 | #' @param mini integer, minimum number of iterations in the EM algorithm, by default = 10 | |
8 | #' @param maxi integer, maximum number of iterations in the EM algorithm, by default = 100 | |
9 | #' @param X matrix of covariates (of size n*p) | |
10 | #' @param Y matrix of responses (of size n*m) | |
11 | #' @param eps real, threshold to say the EM algorithm converges, by default = 1e-4 | |
12 | #' @param rank.min integer, minimum rank in the low rank procedure, by default = 1 | |
13 | #' @param rank.max integer, maximum rank in the low rank procedure, by default = 5 | |
14 | #' @param ncores Number of cores, by default = 3 | |
15 | #' @param fast TRUE to use compiled C code, FALSE for R code only | |
16 | #' @param verbose TRUE to show some execution traces | |
17 | #' | |
18 | #' @return a list with several models, defined by phi, rho, pi, llh | |
2279a641 | 19 | #' |
43d76c49 | 20 | #' @export |
21 | constructionModelesLassoRank = function(S, k, mini, maxi, X, Y, eps, rank.min, | |
22 | rank.max, ncores, fast=TRUE, verbose=FALSE) | |
ef67d338 | 23 | { |
3f145e9a BG |
24 | n = dim(X)[1] |
25 | p = dim(X)[2] | |
43d76c49 | 26 | m = dim(Y)[2] |
27 | L = length(S) | |
28 | ||
29 | # Possible interesting ranks | |
30 | deltaRank = rank.max - rank.min + 1 | |
3f145e9a | 31 | Size = deltaRank^k |
43d76c49 | 32 | RankLambda = matrix(0, nrow=Size*L, ncol=k+1) |
0eb161e3 | 33 | for (r in 1:k) |
43d76c49 | 34 | { |
35 | # On veut le tableau de toutes les combinaisons de rangs possibles, et des lambdas | |
36 | # Dans la première colonne : on répète (rank.max-rank.min)^(k-1) chaque chiffre : | |
37 | # ça remplit la colonne | |
38 | # Dans la deuxieme : on répète (rank.max-rank.min)^(k-2) chaque chiffre, | |
39 | # et on fait ça (rank.max-rank.min)^2 fois | |
40 | # ... | |
41 | # Dans la dernière, on répète chaque chiffre une fois, | |
42 | # et on fait ça (rank.min-rank.max)^(k-1) fois. | |
43 | RankLambda[,r] = rep(rank.min + rep(0:(deltaRank-1), deltaRank^(r-1), each=deltaRank^(k-r)), each = L) | |
ef67d338 | 44 | } |
43d76c49 | 45 | RankLambda[,k+1] = rep(1:L, times = Size) |
46 | ||
0eb161e3 | 47 | if (ncores > 1) |
43d76c49 | 48 | { |
b9b0b72a | 49 | cl = parallel::makeCluster(ncores, outfile='') |
0eb161e3 | 50 | parallel::clusterExport( cl, envir=environment(), |
43d76c49 | 51 | varlist=c("A1","Size","Pi","Rho","mini","maxi","X","Y","eps", |
52 | "Rank","m","phi","ncores","verbose") ) | |
53 | } | |
54 | ||
55 | computeAtLambda <- function(index) | |
56 | { | |
57 | lambdaIndex = RankLambda[index,k+1] | |
58 | rankIndex = RankLambda[index,1:k] | |
59 | if (ncores > 1) | |
60 | require("valse") #workers start with an empty environment | |
61 | ||
62 | # 'relevant' will be the set of relevant columns | |
63 | selected = S[[lambdaIndex]]$selected | |
64 | relevant = c() | |
65 | for (j in 1:p){ | |
66 | if (length(selected[[j]])>0){ | |
67 | relevant = c(relevant,j) | |
3f145e9a BG |
68 | } |
69 | } | |
43d76c49 | 70 | if (max(rankIndex)<length(relevant)){ |
71 | phi = array(0, dim=c(p,m,k)) | |
72 | if (length(relevant) > 0) | |
73 | { | |
74 | res = EMGrank(S[[lambdaIndex]]$Pi, S[[lambdaIndex]]$Rho, mini, maxi, | |
75 | X[,relevant], Y, eps, rankIndex, fast) | |
76 | llh = c( res$LLF, sum(rankIndex * (length(relevant)- rankIndex + m)) ) | |
77 | phi[relevant,,] = res$phi | |
78 | } | |
79 | list("llh"=llh, "phi"=phi, "pi" = S[[lambdaIndex]]$Pi, "rho" = S[[lambdaIndex]]$Rho) | |
80 | ||
81 | } | |
82 | } | |
83 | ||
84 | #For each lambda in the grid we compute the estimators | |
0eb161e3 | 85 | out = |
43d76c49 | 86 | if (ncores > 1) |
87 | parLapply(cl, seq_len(length(S)*Size), computeAtLambda) | |
88 | else | |
89 | lapply(seq_len(length(S)*Size), computeAtLambda) | |
90 | ||
91 | if (ncores > 1) | |
0eb161e3 | 92 | parallel::stopCluster(cl) |
43d76c49 | 93 | |
94 | out | |
9ade3f1b | 95 | } |