}
# Get summary "tableauRecap" from models
- tableauRecap = do.call( rbind, lapply( models_list, function(models) {
+ tableauRecap = do.call( rbind, lapply( seq_along(models_list), function(i) {
+ models <- models_list[[i]]
#Pour un groupe de modeles (même k, différents lambda):
- llh = matrix(ncol = 2)
- for (l in seq_along(models))
- llh = rbind(llh, models[[l]]$llh) #TODO: LLF? harmonize between EMGLLF and EMGrank?
- LLH = llh[-1,1]
- D = llh[-1,2]
- k = length(models[[1]]$pi)
- cbind(LLH, D, rep(k, length(models)), 1:length(models))
+ LLH <- sapply( models, function(model) model$llh )
+ k == length(models[[1]]$pi)
+ # TODO: chuis pas sûr du tout des lignes suivantes...
+ # J'ai l'impression qu'il manque des infos
+ sumPen = sapply( models, function(model)
+ sum( model$pi^gamma * sapply(1:k, function(r) sum(abs(model$phi[,,r]))) ) )
+ data.frame(model=paste(i,".",seq_along(models),sep=""),
+ pen=sumPen/1000, complexity=sumPen, contrast=LLH)
} ) )
- tableauRecap = tableauRecap[rowSums(tableauRecap[, 2:4])!=0,]
- tableauRecap = tableauRecap[(tableauRecap[,1])!=Inf,]
- data = cbind(1:dim(tableauRecap)[1], tableauRecap[,2], tableauRecap[,2], tableauRecap[,1])
-browser()
+
modSel = capushe::capushe(data, n)
indModSel <-
if (selecMod == 'DDSE')
modSel@BIC_capushe$model
else if (selecMod == 'AIC')
modSel@AIC_capushe$model
+
models_list[[tableauRecap[indModSel,3]]][[tableauRecap[indModSel,4]]]
}
#include <gsl/gsl_linalg.h>
// TODO: don't recompute indexes ai(...) and mi(...) when possible
-void EMGLLH_core(
+void EMGLLF_core(
// IN parameters
const Real* phiInit, // parametre initial de moyenne renormalisé
const Real* rhoInit, // parametre initial de variance renormalisé
/////////////
// Precompute det(rho[,,r]) for r in 1...k
+ int signum;
for (int r=0; r<k; r++)
{
for (int u=0; u<m; u++)
detRho[r] = gsl_linalg_LU_det(matrix, signum);
}
- int signum;
Real sumLogLLH = 0.;
for (int i=0; i<n; i++)
{