From 5deedb4ce4b94294f07f9f604ba37a5cbf3425af Mon Sep 17 00:00:00 2001 From: Benjamin Auder Date: Sat, 14 Jan 2017 03:54:13 +0100 Subject: [PATCH] OLD_MATLAB_CODE not required anymore --- .../InputParameters/basicInitParameters.m | 19 -- OLD_MATLAB/InputParameters/compileMex.m | 7 - OLD_MATLAB/InputParameters/generateIO.m | 37 ---- .../InputParameters/generateIOdefault.m | 23 --- OLD_MATLAB/InputParameters/grillelambda.m | 17 -- OLD_MATLAB/InputParameters/initSmallEM.m | 33 ---- OLD_MATLAB/InputParameters/selectiontotale.m | 54 ------ OLD_MATLAB/ProcLassoMLE/EMGLLF.m | 174 ------------------ OLD_MATLAB/ProcLassoMLE/compileMex.m | 9 - .../constructionModelesLassoMLE.m | 58 ------ OLD_MATLAB/ProcLassoRank/EMGrank.m | 69 ------- OLD_MATLAB/ProcLassoRank/compileMex.m | 9 - .../constructionModelesLassoRank.m | 40 ---- .../SelectModel/selectionModelesLassoMLE.m | 12 -- .../SelectModel/selectionModelesLassoRank.m | 11 -- OLD_MATLAB/SelectModel/selectiondindice.m | 19 -- OLD_MATLAB/SelectModel/selectionmodele.m | 20 -- .../SelectModel/suppressionmodelesegaux.m | 18 -- .../SelectModel/suppressionmodelesegaux2.m | 17 -- 19 files changed, 646 deletions(-) delete mode 100644 OLD_MATLAB/InputParameters/basicInitParameters.m delete mode 100644 OLD_MATLAB/InputParameters/compileMex.m delete mode 100644 OLD_MATLAB/InputParameters/generateIO.m delete mode 100644 OLD_MATLAB/InputParameters/generateIOdefault.m delete mode 100644 OLD_MATLAB/InputParameters/grillelambda.m delete mode 100644 OLD_MATLAB/InputParameters/initSmallEM.m delete mode 100644 OLD_MATLAB/InputParameters/selectiontotale.m delete mode 100644 OLD_MATLAB/ProcLassoMLE/EMGLLF.m delete mode 100644 OLD_MATLAB/ProcLassoMLE/compileMex.m delete mode 100644 OLD_MATLAB/ProcLassoMLE/constructionModelesLassoMLE.m delete mode 100644 OLD_MATLAB/ProcLassoRank/EMGrank.m delete mode 100644 OLD_MATLAB/ProcLassoRank/compileMex.m delete mode 100644 OLD_MATLAB/ProcLassoRank/constructionModelesLassoRank.m delete mode 100644 OLD_MATLAB/SelectModel/selectionModelesLassoMLE.m delete mode 100644 OLD_MATLAB/SelectModel/selectionModelesLassoRank.m delete mode 100644 OLD_MATLAB/SelectModel/selectiondindice.m delete mode 100644 OLD_MATLAB/SelectModel/selectionmodele.m delete mode 100644 OLD_MATLAB/SelectModel/suppressionmodelesegaux.m delete mode 100644 OLD_MATLAB/SelectModel/suppressionmodelesegaux2.m diff --git a/OLD_MATLAB/InputParameters/basicInitParameters.m b/OLD_MATLAB/InputParameters/basicInitParameters.m deleted file mode 100644 index 50410f5..0000000 --- a/OLD_MATLAB/InputParameters/basicInitParameters.m +++ /dev/null @@ -1,19 +0,0 @@ -function[phiInit,rhoInit,piInit,gamInit] = basicInitParameters(n,p,m,k) - - phiInit = zeros(p,m,k); - - piInit = (1.0/k) * ones(1,k); - - rhoInit = zeros(m,m,k); - for r=1:k - rhoInit(:,:,r) = eye(m,m); - end - - gamInit = 0.1 * ones(n,k); - R = random('unid',k,n,1); - for i=1:n - gamInit(i,R(i)) = 0.9; - end - gamInit = gamInit / (sum(gamInit(1,:))); - -end diff --git a/OLD_MATLAB/InputParameters/compileMex.m b/OLD_MATLAB/InputParameters/compileMex.m deleted file mode 100644 index 9e5f64b..0000000 --- a/OLD_MATLAB/InputParameters/compileMex.m +++ /dev/null @@ -1,7 +0,0 @@ -%compile C code (for MATLAB or Octave) -if exist('octave_config_info') - setenv('CFLAGS','-O2 -std=gnu99 -fopenmp') - mkoctfile --mex -DOctave -I../Util -I../ProcLassoMLE selectiontotale.c selectiontotale_interface.c ../Util/ioutils.c ../ProcLassoMLE/EMGLLF.c -o selectiontotale -lm -lgsl -lgslcblas -lgomp -else - mex CFLAGS="\$CFLAGS -O2 -std=gnu99 -fopenmp" -I../Util -I../ProcLassoMLE selectiontotale.c selectiontotale_interface.c ../Util/ioutils.c ../ProcLassoMLE/EMGLLF.c -output selectiontotale -lm -lgsl -lgslcblas -lgomp -end diff --git a/OLD_MATLAB/InputParameters/generateIO.m b/OLD_MATLAB/InputParameters/generateIO.m deleted file mode 100644 index c677fd2..0000000 --- a/OLD_MATLAB/InputParameters/generateIO.m +++ /dev/null @@ -1,37 +0,0 @@ -%X is generated following a gaussian mixture \sum pi_r N(meanX_k, covX_r) -%Y is generated then, with Y_i ~ \sum pi_r N(Beta_r.X_i, covY_r) -function[X,Y,Z] = generateIO(meanX, covX, covY, pi, beta, n) - - [p, ~, k] = size(covX); - [m, ~, ~] = size(covY); - if exist('octave_config_info') - %Octave statistics package doesn't have gmdistribution() - X = zeros(n, p); - Z = zeros(n); - cs = cumsum(pi); - for i=1:n - %TODO: vectorize ? http://stackoverflow.com/questions/2977497/weighted-random-numbers-in-matlab - tmpRand01 = rand(); - [~,Z(i)] = min(cs - tmpRand01 >= 0); - X(i,:) = mvnrnd(meanX(Z(i),:), covX(:,:,Z(i)), 1); - end - else - gmDistX = gmdistribution(meanX, covX, pi); - [X, Z] = random(gmDistX, n); - end - - Y = zeros(n, m); - BX = zeros(n,m,k); - for i=1:n - for r=1:k - %compute beta_r . X_i - BXir = zeros(1, m); - for mm=1:m - BXir(mm) = dot(X(i,:), beta(:,mm,r)); - end - %add pi(r) * N(beta_r . X_i, covY) to Y_i - Y(i,:) = Y(i,:) + pi(r) * mvnrnd(BXir, covY(:,:,r), 1); - end - end - -end diff --git a/OLD_MATLAB/InputParameters/generateIOdefault.m b/OLD_MATLAB/InputParameters/generateIOdefault.m deleted file mode 100644 index f4c3c1f..0000000 --- a/OLD_MATLAB/InputParameters/generateIOdefault.m +++ /dev/null @@ -1,23 +0,0 @@ -%call generateIO with default parameters (random means, covariances = identity, equirepartition) -function[X,Y,Z] = generateIOdefault(n, p, m, k) - - rangeX = 100; - meanX = rangeX * (1 - 2*rand(k, p)); - covX = zeros(p,p,k); - covY = zeros(m,m,k); - for r=1:k - covX(:,:,r) = eye(p); - covY(:,:,r) = eye(m); - end - pi = (1/k) * ones(1,k); - - %initialize beta to a random number of non-zero random value - beta = zeros(p,m,k); - for j=1:p - nonZeroCount = ceil(m*rand(1)); - beta(j,1:nonZeroCount,:) = rand(nonZeroCount, k); - end - - [X,Y,Z] = generateIO(meanX, covX, covY, pi, beta, n); - -end diff --git a/OLD_MATLAB/InputParameters/grillelambda.m b/OLD_MATLAB/InputParameters/grillelambda.m deleted file mode 100644 index 8d90665..0000000 --- a/OLD_MATLAB/InputParameters/grillelambda.m +++ /dev/null @@ -1,17 +0,0 @@ -function[gridLambda] = grillelambda(phiInit,rhoInit,piInit,gamInit,X,Y,gamma,mini,maxi,tau) - - n = size(X, 1); - [p,m,k] = size(phiInit); - [phi,rho,pi,~,S] = EMGLLF(phiInit,rhoInit,piInit,gamInit,mini,maxi,1,0,X,Y,tau); - - gridLambda = zeros(p,m,k); - for j=1:p - for mm=1:m - gridLambda(j,mm,:) = abs(reshape(S(j,mm,:),[k,1])) ./ (n*pi.^gamma); - end - end - - gridLambda = unique(gridLambda); - gridLambda(gridLambda()>1) = []; - -end diff --git a/OLD_MATLAB/InputParameters/initSmallEM.m b/OLD_MATLAB/InputParameters/initSmallEM.m deleted file mode 100644 index 6a2ccc9..0000000 --- a/OLD_MATLAB/InputParameters/initSmallEM.m +++ /dev/null @@ -1,33 +0,0 @@ -function[phiInit,rhoInit,piInit,gamInit] = initSmallEM(k,x,y,tau) -[n,m]=size(y); -gamInit1=zeros(n,k,20); -for repet=1:20 - Zinit1(:,repet)=clusterdata(y,k); - for r=1:k - betaInit1(:,:,r,repet)=pinv(transpose(x(Zinit1(:,repet)==r,:))*x(Zinit1(:,repet)==r,:))*transpose(x(Zinit1(:,repet)==r,:))*y(Zinit1(:,repet)==r,:); - sigmaInit1(:,:,r,repet)=eye(m); - phiInit1(:,:,r,repet)=betaInit1(:,:,r,repet)/sigmaInit1(:,:,r,repet); - rhoInit1(:,:,r,repet)=inv(sigmaInit1(:,:,r,repet)); - piInit1(repet,r)=sum(Zinit1(:,repet)==r)/n; - end - for i=1:n - for r=1:k - dotProduct = (y(i,:)*rhoInit1(:,:,r,repet)-x(i,:)*phiInit1(:,:,r,repet)) * transpose(y(i,:)*rhoInit1(:,:,r,repet)-x(i,:)*phiInit1(:,:,r,repet)); - Gam(i,r) = piInit1(repet,r)*det(rhoInit1(:,:,r,repet))*exp(-0.5*dotProduct); - end - sumGamI = sum(Gam(i,:)); - gamInit1(i,:,repet) = Gam(i,:) / sumGamI; - end - miniInit=int64(10); - maxiInit=int64(11); - - [~,~,~,LLFEssai,~] = EMGLLF(phiInit1(:,:,:,repet),rhoInit1(:,:,:,repet),piInit1(repet,:),gamInit1(:,:,repet),miniInit,maxiInit,1,0,x,y,tau); - LLFinit1(repet)=LLFEssai(end); -end -[~,b]=max(LLFinit1); - -phiInit=phiInit1(:,:,:,b); -rhoInit=rhoInit1(:,:,:,b); -piInit=piInit1(b,:); -gamInit=gamInit1(:,:,b); -end \ No newline at end of file diff --git a/OLD_MATLAB/InputParameters/selectiontotale.m b/OLD_MATLAB/InputParameters/selectiontotale.m deleted file mode 100644 index 6c24d05..0000000 --- a/OLD_MATLAB/InputParameters/selectiontotale.m +++ /dev/null @@ -1,54 +0,0 @@ -function[A1,A2,Rho,Pi] = selectiontotale(phiInit,rhoInit,piInit,gamInit,mini,maxi,gamma,glambda,X,Y,seuil,tau) - - [p,m,k] = size(phiInit); - L = length(glambda); - A1 = zeros(p,m+1,L,'int64'); - A2 = zeros(p,m+1,L,'int64'); - Rho = zeros(m,m,k,L); - Pi = zeros(k,L); - - %Pour chaque lambda de la grille, on calcule les coefficients - for lambdaIndex=1:L - [phi,rho,pi,~,~] = EMGLLF(phiInit,rhoInit,piInit,gamInit,mini,maxi,gamma,glambda(lambdaIndex),X,Y,tau); - - %Si un des coefficients est supérieur au seuil, on garde cette variable - selectedVariables = zeros(p,m); - discardedVariables = zeros(p,m); - atLeastOneSelectedVariable = false; - for j=1:p - cpt=1; - cpt2=1; - for mm=1:m - if max(abs(phi(j,mm,:))) > seuil - selectedVariables(j,cpt) = mm; - cpt = cpt+1; - atLeastOneSelectedVariable = true; - else - discardedVariables(j,cpt2) = mm; - cpt2 = cpt2+1; - end - end - end - - %Si aucun des coefficients n'a été gardé on renvoit la matrice nulle - %Et si on enlevait ces colonnes de zéro ??? Indices des colonnes vides - if atLeastOneSelectedVariable - vec = []; - for j=1:p - if selectedVariables(j,1) ~= 0 - vec = [vec;j]; - end - end - - %Sinon on renvoit les numéros des coefficients utiles - A1(:,1,lambdaIndex) = [vec;zeros(p-length(vec),1)]; - A1(1:length(vec),2:m+1,lambdaIndex) = selectedVariables(vec,:); - A2(:,1,lambdaIndex) = 1:p; - A2(:,2:m+1,lambdaIndex) = discardedVariables; - Rho(:,:,:,lambdaIndex) = rho; - Pi(:,lambdaIndex) = pi; - end - - end - -end diff --git a/OLD_MATLAB/ProcLassoMLE/EMGLLF.m b/OLD_MATLAB/ProcLassoMLE/EMGLLF.m deleted file mode 100644 index 1be6ba0..0000000 --- a/OLD_MATLAB/ProcLassoMLE/EMGLLF.m +++ /dev/null @@ -1,174 +0,0 @@ -function[phi,rho,pi,LLF,S] = EMGLLF(phiInit,rhoInit,piInit,gamInit,mini,maxi,gamma,lambda,X,Y,tau) - - %Get matrices dimensions - PI = 4.0 * atan(1.0); - n = size(X, 1); - [p,m,k] = size(phiInit); - - %Initialize outputs - phi = phiInit; - rho = rhoInit; - pi = piInit; - LLF = zeros(maxi,1); - S = zeros(p,m,k); - - %Other local variables - %NOTE: variables order is always n,p,m,k - gam = gamInit; - Gram2 = zeros(p,p,k); - ps2 = zeros(p,m,k); - b = zeros(k,1); - pen = zeros(maxi,k); - X2 = zeros(n,p,k); - Y2 = zeros(n,m,k); - dist = 0; - dist2 = 0; - ite = 1; - pi2 = zeros(k,1); - ps = zeros(m,k); - nY2 = zeros(m,k); - ps1 = zeros(n,m,k); - nY21 = zeros(n,m,k); - Gam = zeros(n,k); - EPS = 1e-15; - - while ite<=mini || (ite<=maxi && (dist>=tau || dist2>=sqrt(tau))) - - Phi = phi; - Rho = rho; - Pi = pi; - - %Calculs associés à Y et X - for r=1:k - for mm=1:m - Y2(:,mm,r) = sqrt(gam(:,r)) .* Y(:,mm); - end - for i=1:n - X2(i,:,r) = X(i,:) .* sqrt(gam(i,r)); - end - for mm=1:m - ps2(:,mm,r) = transpose(X2(:,:,r)) * Y2(:,mm,r); - end - for j=1:p - for s=1:p - Gram2(j,s,r) = dot(X2(:,j,r), X2(:,s,r)); - end - end - end - - %%%%%%%%%% - %Etape M % - %%%%%%%%%% - - %Pour pi - for r=1:k - b(r) = sum(sum(abs(phi(:,:,r)))); - end - gam2 = sum(gam,1); - a = sum(gam*transpose(log(pi))); - - %tant que les proportions sont negatives - kk = 0; - pi2AllPositive = false; - while ~pi2AllPositive - pi2 = pi + 0.1^kk * ((1/n)*gam2 - pi); - pi2AllPositive = true; - for r=1:k - if pi2(r) < 0 - pi2AllPositive = false; - break; - end - end - kk = kk+1; - end - - %t(m) la plus grande valeur dans la grille O.1^k tel que ce soit - %décroissante ou constante - while (-1/n*a+lambda*((pi.^gamma)*b))<(-1/n*gam2*transpose(log(pi2))+lambda.*(pi2.^gamma)*b) && kk<1000 - pi2 = pi+0.1^kk*(1/n*gam2-pi); - kk = kk+1; - end - t = 0.1^(kk); - pi = (pi+t*(pi2-pi)) / sum(pi+t*(pi2-pi)); - - %Pour phi et rho - for r=1:k - for mm=1:m - for i=1:n - ps1(i,mm,r) = Y2(i,mm,r) * dot(X2(i,:,r), phi(:,mm,r)); - nY21(i,mm,r) = (Y2(i,mm,r))^2; - end - ps(mm,r) = sum(ps1(:,mm,r)); - nY2(mm,r) = sum(nY21(:,mm,r)); - rho(mm,mm,r) = ((ps(mm,r)+sqrt(ps(mm,r)^2+4*nY2(mm,r)*(gam2(r))))/(2*nY2(mm,r))); - end - end - for r=1:k - for j=1:p - for mm=1:m - S(j,mm,r) = -rho(mm,mm,r)*ps2(j,mm,r) + dot(phi(1:j-1,mm,r),Gram2(j,1:j-1,r)')... - + dot(phi(j+1:p,mm,r),Gram2(j,j+1:p,r)'); - if abs(S(j,mm,r)) <= n*lambda*(pi(r)^gamma) - phi(j,mm,r)=0; - else - if S(j,mm,r)> n*lambda*(pi(r)^gamma) - phi(j,mm,r)=(n*lambda*(pi(r)^gamma)-S(j,mm,r))/Gram2(j,j,r); - else - phi(j,mm,r)=-(n*lambda*(pi(r)^gamma)+S(j,mm,r))/Gram2(j,j,r); - end - end - end - end - end - - %%%%%%%%%% - %Etape E % - %%%%%%%%%% - - sumLogLLF2 = 0.0; - for i=1:n - %precompute dot products to numerically adjust their values - dotProducts = zeros(k,1); - for r=1:k - dotProducts(r)= (Y(i,:)*rho(:,:,r)-X(i,:)*phi(:,:,r)) * transpose(Y(i,:)*rho(:,:,r)-X(i,:)*phi(:,:,r)); - end - shift = 0.5*min(dotProducts); - - %compute Gam(:,:) using shift determined above - sumLLF1 = 0.0; - for r=1:k - Gam(i,r) = pi(r)*det(rho(:,:,r))*exp(-0.5*dotProducts(r) + shift); - sumLLF1 = sumLLF1 + Gam(i,r)/(2*PI)^(m/2); - end - sumLogLLF2 = sumLogLLF2 + log(sumLLF1); - sumGamI = sum(Gam(i,:)); - if sumGamI > EPS - gam(i,:) = Gam(i,:) / sumGamI; - else - gam(i,:) = zeros(k,1); - end - end - - sumPen = 0.0; - for r=1:k - sumPen = sumPen + pi(r).^gamma .* b(r); - end - LLF(ite) = -(1/n)*sumLogLLF2 + lambda*sumPen; - - if ite == 1 - dist = LLF(ite); - else - dist = (LLF(ite)-LLF(ite-1))/(1+abs(LLF(ite))); - end - - Dist1=max(max(max((abs(phi-Phi))./(1+abs(phi))))); - Dist2=max(max(max((abs(rho-Rho))./(1+abs(rho))))); - Dist3=max(max((abs(pi-Pi))./(1+abs(Pi)))); - dist2=max([Dist1,Dist2,Dist3]); - - ite=ite+1; - end - - pi = transpose(pi); - -end diff --git a/OLD_MATLAB/ProcLassoMLE/compileMex.m b/OLD_MATLAB/ProcLassoMLE/compileMex.m deleted file mode 100644 index e595409..0000000 --- a/OLD_MATLAB/ProcLassoMLE/compileMex.m +++ /dev/null @@ -1,9 +0,0 @@ -%compile C code (for MATLAB or Octave) -if exist('octave_config_info') - setenv('CFLAGS','-O2 -std=gnu99 -fopenmp') - mkoctfile --mex -DOctave -I../Util EMGLLF.c EMGLLF_interface.c ../Util/ioutils.c -o EMGLLF -lm -lgsl -lgslcblas -lgomp - mkoctfile --mex -DOctave -I../Util constructionModelesLassoMLE.c constructionModelesLassoMLE_interface.c EMGLLF.c ../Util/ioutils.c -o constructionModelesLassoMLE -lm -lgsl -lgslcblas -lgomp -else - mex CFLAGS="\$CFLAGS -O2 -std=gnu99 -fopenmp" -I../Util EMGLLF.c EMGLLF_interface.c ../Util/ioutils.c -output EMGLLF -lm -lgsl -lgslcblas -lgomp - mex CFLAGS="\$CFLAGS -O2 -std=gnu99 -fopenmp" -I../Util constructionModelesLassoMLE.c constructionModelesLassoMLE_interface.c EMGLLF.c ../Util/ioutils.c -output constructionModelesLassoMLE -lm -lgsl -lgslcblas -lgomp -end diff --git a/OLD_MATLAB/ProcLassoMLE/constructionModelesLassoMLE.m b/OLD_MATLAB/ProcLassoMLE/constructionModelesLassoMLE.m deleted file mode 100644 index 9f977ac..0000000 --- a/OLD_MATLAB/ProcLassoMLE/constructionModelesLassoMLE.m +++ /dev/null @@ -1,58 +0,0 @@ -function[phi,rho,pi,lvraisemblance] = constructionModelesLassoMLE(... - phiInit,rhoInit,piInit,gamInit,mini,maxi,gamma,glambda,X,Y,seuil,tau,A1,A2) - - PI = 4.0 * atan(1.0); - - %get matrix sizes - n = size(X, 1); - [p,m,k] = size(phiInit); - L = length(glambda); - - %output parameters - phi = zeros(p,m,k,L); - rho = zeros(m,m,k,L); - pi = zeros(k,L); - lvraisemblance = zeros(L,2); - - for lambdaIndex=1:L - % Procedure Lasso-MLE - a = A1(:,1,lambdaIndex); - a(a==0) = []; - if length(a) == 0 - continue; - end - [phiLambda,rhoLambda,piLambda,~,~] = EMGLLF(... - phiInit(a,:,:),rhoInit,piInit,gamInit,mini,maxi,gamma,0,X(:,a),Y,tau); - - for j=1:length(a) - phi(a(j),:,:,lambdaIndex) = phiLambda(j,:,:); - end - rho(:,:,:,lambdaIndex) = rhoLambda; - pi(:,lambdaIndex) = piLambda; - - dimension = 0; - for j=1:p - b = A2(j,2:end,lambdaIndex); - b(b==0) = []; - if length(b) > 0 - phi(A2(j,1,lambdaIndex),b,:,lambdaIndex) = 0.0; - end - c = A1(j,2:end,lambdaIndex); - c(c==0) = []; - dimension = dimension + length(c); - end - - %on veut calculer l'EMV avec toutes nos estimations - densite = zeros(n,L); - for i=1:n - for r=1:k - delta = Y(i,:)*rho(:,:,r,lambdaIndex) - (X(i,a)*(phi(a,:,r,lambdaIndex))); - densite(i,lambdaIndex) = densite(i,lambdaIndex) +... - pi(r,lambdaIndex)*det(rho(:,:,r,lambdaIndex))/(sqrt(2*PI))^m*exp(-dot(delta,delta)/2.0); - end - end - lvraisemblance(lambdaIndex,1) = sum(log(densite(:,lambdaIndex))); - lvraisemblance(lambdaIndex,2) = (dimension+m+1)*k-1; - end - -end diff --git a/OLD_MATLAB/ProcLassoRank/EMGrank.m b/OLD_MATLAB/ProcLassoRank/EMGrank.m deleted file mode 100644 index 6ffc313..0000000 --- a/OLD_MATLAB/ProcLassoRank/EMGrank.m +++ /dev/null @@ -1,69 +0,0 @@ -function[phi,LLF] = EMGrank(Pi,Rho,mini,maxi,X,Y,tau,rank) - - % get matrix sizes - [~,m,k] = size(Rho); - [n,p] = size(X); - - % allocate output matrices - phi = zeros(p,m,k); - Z = ones(n,1,'int64'); - LLF = 0.0; - - % local variables - Phi = zeros(p,m,k); - deltaPhi = 0.0; - deltaPhi = []; - sumDeltaPhi = 0.0; - deltaPhiBufferSize = 20; - - %main loop (at least mini iterations) - ite = int64(1); - while ite<=mini || (ite<=maxi && sumDeltaPhi>tau) - - %M step: Mise à jour de Beta (et donc phi) - for r=1:k - if (sum(Z==r) == 0) - continue; - end - %U,S,V = SVD of (t(Xr)Xr)^{-1} * t(Xr) * Yr - [U,S,V] = svd(pinv(transpose(X(Z==r,:))*X(Z==r,:))*transpose(X(Z==r,:))*Y(Z==r,:)); - %Set m-rank(r) singular values to zero, and recompose - %best rank(r) approximation of the initial product - S(rank(r)+1:end,:) = 0; - phi(:,:,r) = U * S * transpose(V) * Rho(:,:,r); - end - - %Etape E et calcul de LLF - sumLogLLF2 = 0.0; - for i=1:n - sumLLF1 = 0.0; - maxLogGamIR = -Inf; - for r=1:k - dotProduct = (Y(i,:)*Rho(:,:,r)-X(i,:)*phi(:,:,r)) * transpose(Y(i,:)*Rho(:,:,r)-X(i,:)*phi(:,:,r)); - logGamIR = log(Pi(r)) + log(det(Rho(:,:,r))) - 0.5*dotProduct; - %Z(i) = index of max (gam(i,:)) - if logGamIR > maxLogGamIR - Z(i) = r; - maxLogGamIR = logGamIR; - end - sumLLF1 = sumLLF1 + exp(logGamIR) / (2*pi)^(m/2); - end - sumLogLLF2 = sumLogLLF2 + log(sumLLF1); - end - - LLF = -1/n * sumLogLLF2; - - % update distance parameter to check algorithm convergence (delta(phi, Phi)) - deltaPhi = [ deltaPhi, max(max(max((abs(phi-Phi))./(1+abs(phi))))) ]; - if length(deltaPhi) > deltaPhiBufferSize - deltaPhi = deltaPhi(2:length(deltaPhi)); - end - sumDeltaPhi = sum(abs(deltaPhi)); - - % update other local variables - Phi = phi; - ite = ite+1; - - end - -end diff --git a/OLD_MATLAB/ProcLassoRank/compileMex.m b/OLD_MATLAB/ProcLassoRank/compileMex.m deleted file mode 100644 index 044f442..0000000 --- a/OLD_MATLAB/ProcLassoRank/compileMex.m +++ /dev/null @@ -1,9 +0,0 @@ -%compile C code (for MATLAB or Octave) -if exist('octave_config_info') - setenv('CFLAGS','-O2 -std=gnu99 -fopenmp') - mkoctfile --mex -DOctave -I../Util EMGrank.c EMGrank_interface.c ../Util/ioutils.c -o EMGrank -lm -lgsl -lgslcblas -lgomp - mkoctfile --mex -DOctave -I../Util constructionModelesLassoRank.c constructionModelesLassoRank_interface.c EMGrank.c ../Util/ioutils.c -o constructionModelesLassoRank -lm -lgsl -lgslcblas -lgomp -else - mex CFLAGS="\$CFLAGS -O2 -std=gnu99 -fopenmp" -I../Util EMGrank.c EMGrank_interface.c ../Util/ioutils.c -output EMGrank -lm -lgsl -lgslcblas -lgomp - mex CFLAGS="\$CFLAGS -O2 -std=gnu99 -fopenmp" -I../Util constructionModelesLassoRank.c constructionModelesLassoRank_interface.c EMGrank.c ../Util/ioutils.c -output constructionModelesLassoRank -lm -lgsl -lgslcblas -lgomp -end diff --git a/OLD_MATLAB/ProcLassoRank/constructionModelesLassoRank.m b/OLD_MATLAB/ProcLassoRank/constructionModelesLassoRank.m deleted file mode 100644 index ae5e34e..0000000 --- a/OLD_MATLAB/ProcLassoRank/constructionModelesLassoRank.m +++ /dev/null @@ -1,40 +0,0 @@ -function[phi,lvraisemblance] = constructionModelesLassoRank(Pi,Rho,mini,maxi,X,Y,tau,A1,rangmin,rangmax) - - PI = 4.0 * atan(1.0); - - %get matrix sizes - [n,p] = size(X); - [~,m,k,~] = size(Rho); - L = size(A1, 2); %A1 est p x m+1 x L ou p x L ?! - - %On cherche les rangs possiblement intéressants - deltaRank = rangmax - rangmin + 1; - Size = deltaRank^k; - Rank = zeros(Size,k,'int64'); - for r=1:k - %On veut le tableau de toutes les combinaisons de rangs possibles - %Dans la première colonne : on répète (rangmax-rangmin)^(k-1) chaque chiffre : ca remplit la colonne - %Dans la deuxieme : on répète (rangmax-rangmin)^(k-2) chaque chiffre, et on fait ca (rangmax-rangmin)^2 fois - %... - %Dans la dernière, on répète chaque chiffre une fois, et on fait ca (rangmin-rangmax)^(k-1) fois. - Rank(:,r) = rangmin + reshape(repmat(0:(deltaRank-1), deltaRank^(k-r), deltaRank^(r-1)), Size, 1); - end - - %output parameters - phi = zeros(p,m,k,L*Size); - lvraisemblance = zeros(L*Size,2); - for lambdaIndex=1:L - %On ne garde que les colonnes actives - %active sera l'ensemble des variables informatives - active = A1(:,lambdaIndex); - active(active==0) = []; - if length(active) > 0 - for j=1:Size - [phiLambda,LLF] = EMGrank(Pi(:,lambdaIndex),Rho(:,:,:,lambdaIndex),mini,maxi,X(:,active),Y,tau,Rank(j,:)); - lvraisemblance((lambdaIndex-1)*Size+j,:) = [LLF, sum(Rank(j,:) .* (length(active)-Rank(j,:)+m))]; - phi(active,:,:,(lambdaIndex-1)*Size+j) = phiLambda; - end - end - end - -end diff --git a/OLD_MATLAB/SelectModel/selectionModelesLassoMLE.m b/OLD_MATLAB/SelectModel/selectionModelesLassoMLE.m deleted file mode 100644 index ab39371..0000000 --- a/OLD_MATLAB/SelectModel/selectionModelesLassoMLE.m +++ /dev/null @@ -1,12 +0,0 @@ -%Selection de modèle dans la procédure Lasso-MLE -% Selection de modele : voici les parametres selectionnes pour chaque -% dimension, et l'indice associe -[indiceLassoMLE,D1LassoMLE] = selectionmodele(vraisLassoMLE); - -%on veut tracer la courbe -%figure(2) -%plot(D1LassoMLE/n,vraisLassoMLE(indiceLassoMLE),'.') -%on veut appliquer Capushe : il nous faut un tableau de données -tableauLassoMLE = enregistrerdonnees(D1LassoMLE,vraisLassoMLE(indiceLassoMLE,:),n); -save tableauLassoMLE.mat tableauLassoMLE -%On conclut avec Capushe ! diff --git a/OLD_MATLAB/SelectModel/selectionModelesLassoRank.m b/OLD_MATLAB/SelectModel/selectionModelesLassoRank.m deleted file mode 100644 index aacc460..0000000 --- a/OLD_MATLAB/SelectModel/selectionModelesLassoRank.m +++ /dev/null @@ -1,11 +0,0 @@ -%Selection de modèle dans la procedure Lasso Rank -vraisLassoRank=vraisLassoRank'; -% Selection de modele : voici les parametres selectionnes pour chaque -% dimension, et l'indice associe -[indiceLassoRank,D1LassoRank]=selectionmodele(vraisLassoRank); -tableauLassoRank=enregistrerdonnees(D1LassoRank,vraisLassoRank(indiceLassoRank,:),n); -%On veut tracer la courbe des pentes -%figure(2) -%plot(D1LassoRank/n,vraisLassoRank(indiceLassoRank),'.') -save tableauLassoRank.mat tableauLassoRank -%On conclut avec Capushe ! diff --git a/OLD_MATLAB/SelectModel/selectiondindice.m b/OLD_MATLAB/SelectModel/selectiondindice.m deleted file mode 100644 index 4ef8b96..0000000 --- a/OLD_MATLAB/SelectModel/selectiondindice.m +++ /dev/null @@ -1,19 +0,0 @@ -%utile dans selectiontotale.m, remarque quels coefficients sont nuls et lesquels ne le sont pas -function[A,B]=selectiondindice(phi,seuil) - - [pp,m,~]=size(phi); - A=zeros(pp,m); - B=zeros(pp,m); - for j=1:pp - cpt=0;cpt2=0; - for mm=1:m - if (max(phi(j,mm,:)) > seuil) - cpt=cpt+1; - A(j,cpt)=mm; - else cpt2=cpt2+1; - B(j,cpt2)=mm; - end - end - end - -end diff --git a/OLD_MATLAB/SelectModel/selectionmodele.m b/OLD_MATLAB/SelectModel/selectionmodele.m deleted file mode 100644 index 33828ab..0000000 --- a/OLD_MATLAB/SelectModel/selectionmodele.m +++ /dev/null @@ -1,20 +0,0 @@ -function [indice,D1]=selectionmodele(vraisemblance) - - D=vraisemblance(:,2); - [D1]=unique(D); - indice=ones(1,length(D1)); - %On ne sélectionne que celui qui maximise : l'EMV - if length(D1)>2 - for i=1:length(D1) - a=[]; - for j=1:length(D) - if D(j)==D1(i) - a=[a,vraisemblance(j,1)]; - end - end - b=max(a); - indice(i)=find(vraisemblance(:,1)==b,1); - end - end - -end diff --git a/OLD_MATLAB/SelectModel/suppressionmodelesegaux.m b/OLD_MATLAB/SelectModel/suppressionmodelesegaux.m deleted file mode 100644 index 3d2a3a9..0000000 --- a/OLD_MATLAB/SelectModel/suppressionmodelesegaux.m +++ /dev/null @@ -1,18 +0,0 @@ -function[B1,B2,glambda,ind,rho,pi]=suppressionmodelesegaux(B1,B2,glambda,rho,pi) - - ind=[]; - for l=1:length(glambda) - for ll=1:l-1 - if B1(:,:,l)==B1(:,:,ll) - ind=[ind l]; - end - end - end - ind=unique(ind); - B1(:,:,ind)=[]; - glambda(ind)=[]; - B2(:,:,ind)=[]; - rho(:,:,:,ind)=[]; - pi(:,ind)=[]; - -end diff --git a/OLD_MATLAB/SelectModel/suppressionmodelesegaux2.m b/OLD_MATLAB/SelectModel/suppressionmodelesegaux2.m deleted file mode 100644 index 3158d36..0000000 --- a/OLD_MATLAB/SelectModel/suppressionmodelesegaux2.m +++ /dev/null @@ -1,17 +0,0 @@ -function[B1,ind,rho,pi]=suppressionmodelesegaux2(B1,rho,pi) - - ind=[]; - nombreLambda=size(B1,2); - for l=1:nombreLambda - for ll=1:l-1 - if B1(:,l)==B1(:,ll) - ind=[ind l]; - end - end - end - ind=unique(ind); - B1(:,ind)=[]; - rho(:,:,:,ind)=[]; - pi(:,ind)=[]; - -end -- 2.44.0