From: Benjamin Goehry Date: Tue, 10 Jan 2017 10:31:01 +0000 (+0100) Subject: typo + changement de clustering initial (mclust -> hclust) X-Git-Url: https://git.auder.net/%7B%7B%20asset%28%27mixstore/css/user/assets/rpsls.css?a=commitdiff_plain;h=4725af564bfdbeb40d0caf7f124ab94ac05a97eb;p=valse.git typo + changement de clustering initial (mclust -> hclust) --- diff --git a/R/initSmallEM.R b/R/initSmallEM.R index b70eea9..c836523 100644 --- a/R/initSmallEM.R +++ b/R/initSmallEM.R @@ -13,7 +13,7 @@ initSmallEM = function(k,X,Y,tau) m = ncol(Y) p = ncol(X) - Zinit1 = array(0, dim=c(n,20)) #doute sur la taille + Zinit1 = array(0, dim=c(n,20)) betaInit1 = array(0, dim=c(p,m,k,20)) sigmaInit1 = array(0, dim = c(m,m,k,20)) phiInit1 = array(0, dim = c(p,m,k,20)) @@ -24,12 +24,12 @@ initSmallEM = function(k,X,Y,tau) LLFinit1 = list() require(MASS) #Moore-Penrose generalized inverse of matrix - require(mclust) # K-means with selection of K for(repet in 1:20) { - clusters = Mclust(X,k) #default distance : euclidean #Mclust(matrix(c(X,Y)),k) - Zinit1[,repet] = clusters$classification - + distance_clus = dist(X) + tree_hier = hclust(distance_clus) + Zinit1[,repet] = cutree(tree_hier, k) + for(r in 1:k) { Z = Zinit1[,repet] @@ -39,7 +39,7 @@ initSmallEM = function(k,X,Y,tau) betaInit1[,,r,repet] = ginv( crossprod(X[Z_indice,]) ) %*% crossprod(X[Z_indice,], Y[Z_indice,]) sigmaInit1[,,r,repet] = diag(m) - phiInit1[,,r,repet] = betaInit1[,,r,repet]/sigmaInit1[,,r,repet] + phiInit1[,,r,repet] = betaInit1[,,r,repet] #/ sigmaInit1[,,r,repet] rhoInit1[,,r,repet] = solve(sigmaInit1[,,r,repet]) piInit1[repet,r] = sum(Z_vec)/n } @@ -48,7 +48,7 @@ initSmallEM = function(k,X,Y,tau) { for(r in 1:k) { - dotProduct = (Y[i,]%*%rhoInit1[,,r,repet]-X[i,]%*%phiInit1[,,r,repet]) %*% (Y[i,]%*%rhoInit1[,,r,repet]-X[i,]%*%phiInit1[,,r,repet]) + dotProduct = tcrossprod(Y[i,]%*%rhoInit1[,,r,repet]-X[i,]%*%phiInit1[,,r,repet]) Gam[i,r] = piInit1[repet,r]*det(rhoInit1[,,r,repet])*exp(-0.5*dotProduct) } sumGamI = sum(Gam[i,]) @@ -58,8 +58,7 @@ initSmallEM = function(k,X,Y,tau) miniInit = 10 maxiInit = 11 - new_EMG = .Call("EMGLLF_core",phiInit1[,,,repet],rhoInit1[,,,repet],piInit1[repet,], - gamInit1[,,repet],miniInit,maxiInit,1,0,X,Y,tau) + new_EMG = .Call("EMGLLF_core",phiInit1[,,,repet],rhoInit1[,,,repet],piInit1[repet,],gamInit1[,,repet],miniInit,maxiInit,1,0,X,Y,tau) LLFEessai = new_EMG$LLF LLFinit1[repet] = LLFEessai[length(LLFEessai)] }