From 7b2720733e9aebe177c211119a9ec160c7e7117b Mon Sep 17 00:00:00 2001
From: Benjamin Auder <benjamin.auder@somewhere>
Date: Wed, 11 Jan 2017 17:49:31 +0100
Subject: [PATCH] add test folder

---
 src/test/.test.EMGLLF.c.swp                   | Bin 0 -> 12288 bytes
 src/test/.utils.c.swp                         | Bin 0 -> 12288 bytes
 src/test/Makefile                             |  37 +++++
 src/test/OLD_TEST_MATLAB/TODO                 |   4 +
 .../OLD_TEST_MATLAB/basicInitParameters.m     |  19 +++
 src/test/OLD_TEST_MATLAB/checkOutput.m        |  11 ++
 src/test/OLD_TEST_MATLAB/covariance.m         |   9 ++
 src/test/OLD_TEST_MATLAB/generateIO.m         |  37 +++++
 src/test/OLD_TEST_MATLAB/generateIOdefault.m  |  23 +++
 .../generateRunSaveTest_EMGLLF.m              |  47 ++++++
 ...eRunSaveTest_constructionModelesLassoMLE.m |  62 ++++++++
 .../generateRunSaveTest_selectiontotale.m     |  49 ++++++
 .../testConstructionModelesLassoMLE.m         |  46 ++++++
 src/test/OLD_TEST_MATLAB/testEMGLLF.m         |  44 ++++++
 .../OLD_TEST_MATLAB/testSelectiontotale.m     |  44 ++++++
 src/test/test.ConstructionModelesLassoMLE.c   | 143 ++++++++++++++++++
 src/test/test.EMGLLF.c                        |  81 ++++++++++
 src/test/test.EMGrank.c                       |  92 +++++++++++
 src/test/test.constructionModelesLassoRank.c  | 107 +++++++++++++
 src/test/test.selectiontotale.c               | 134 ++++++++++++++++
 src/test/utils.c                              |  81 ++++++++++
 21 files changed, 1070 insertions(+)
 create mode 100644 src/test/.test.EMGLLF.c.swp
 create mode 100644 src/test/.utils.c.swp
 create mode 100644 src/test/Makefile
 create mode 100644 src/test/OLD_TEST_MATLAB/TODO
 create mode 100644 src/test/OLD_TEST_MATLAB/basicInitParameters.m
 create mode 100644 src/test/OLD_TEST_MATLAB/checkOutput.m
 create mode 100644 src/test/OLD_TEST_MATLAB/covariance.m
 create mode 100644 src/test/OLD_TEST_MATLAB/generateIO.m
 create mode 100644 src/test/OLD_TEST_MATLAB/generateIOdefault.m
 create mode 100644 src/test/OLD_TEST_MATLAB/generateRunSaveTest_EMGLLF.m
 create mode 100644 src/test/OLD_TEST_MATLAB/generateRunSaveTest_constructionModelesLassoMLE.m
 create mode 100644 src/test/OLD_TEST_MATLAB/generateRunSaveTest_selectiontotale.m
 create mode 100644 src/test/OLD_TEST_MATLAB/testConstructionModelesLassoMLE.m
 create mode 100644 src/test/OLD_TEST_MATLAB/testEMGLLF.m
 create mode 100644 src/test/OLD_TEST_MATLAB/testSelectiontotale.m
 create mode 100644 src/test/test.ConstructionModelesLassoMLE.c
 create mode 100644 src/test/test.EMGLLF.c
 create mode 100644 src/test/test.EMGrank.c
 create mode 100644 src/test/test.constructionModelesLassoRank.c
 create mode 100644 src/test/test.selectiontotale.c
 create mode 100644 src/test/utils.c

diff --git a/src/test/.test.EMGLLF.c.swp b/src/test/.test.EMGLLF.c.swp
new file mode 100644
index 0000000000000000000000000000000000000000..176ac00789fa190350c8e35d63ce68e5f447f061
GIT binary patch
literal 12288
zcmeI2&x;&I6vrz<NXEEP^x(xnIZ>D$%k4>4A!JbqYQkz_)Y(mD5LtSAYBR-7|FF7e
zHHwJnCE`JE>RAK@Z{k7y3j{CX5A+gGo}$s4r{K4K%uM%gFH0^$;5B@@tLnY#SMOV|
zdl|OZI~SJE@e>OjhT|^A?p@pJo!@_hwYC@w#{;R&_#VgbawfG(lRey6FdPg5EweNa
zw!%>^gIq^JA@jnX7M?nH=Iq(i3(;P;Av-BR3f!y$<6`6C6LW0w_%ZvaA8p>x@40*R
zX757&qyQ;E3XlS%04YEUkOHItDR2ubP-Oeqn^^yCK8NT1b>FtD@6v@7AO%PPQh*d7
z1xNu>fD|AFNC8rS6d(m|K?TAAWBLGNKODs3@&EtH@BeFeGWHAj8GHr41fPNr!MorJ
z=!3_>QE(9ad6==^z^~v7@D{iT9s+aV`W$25gU`SfpujWWX|MzYI0Ehi*Y05KTkru`
z2akY7a0DCx{~TiMZ}16tAG`zJ204hp3V0Hn1RW57!{CqGkq`JBd<;GUFN0IyQScym
z0NiLW_7nIHd;`7)?|}$h1ZTi$&;<WNyT8Ei;41hL7(K58I!FOhfD|AFNC8sd|4`r!
z*3hyT>x6ff7;9{3DO*}@taay?e9yL`G|oaTPiY-qUehuhwVK_g;9bG3v*1}fZf}}7
zFTmJ#qx>;k!$EZtd`e;P3F9=7Z`Xu^1;KE3;k@H45(a8wxmj`es;b7!Oy6bdnl%wK
zg@BoI@yvBwZkehHh)LB*dYDc_`F2e#SP%?n7tJ=}opNHlCZ?HO4Vn#0`MGvYP<WXP
zhm#adqh=4`P!rTlHi-oRf5H_8FH*DTQp+TgJRKL=xX33P-Jgoko3<{6@p7VyX+Kky
zF4C&>5z+TqqmhVA#KP`K0|sFf!ycLX=nB?wU&gLQsg^B2&$yWRaI+U5UBoI;!c>UB
zNF0h$7_SdP(Gz`9gkymOf||o69*#yPYZr|N!7jd&C!z}-Tef89V@yq>mX*7mtJh?@
zVJ-7U4Oi8!qO1O9ox8Z(a;$o9?&1n1cC5r!+b;f4aEE{G%D3fjMLNH_^4#i5cc(f~
z)+EZIvac`yR>6QREhj_Y8_Tbde#zdevG=BS?zXcUz>Up(?8aAJ?2z?XP6}>hpV(DW
zw3=pMW<>>6tr(czz0G~4d6#k?-cs^)#aZL(3?)xr&NZ&yQu4(Hnwqci^%mpn--iPg
z%OuCg>A5Mz_?0e-lW$_`%<$skOEXhkH61URva0D=$&^%09b<yh?NStY=;@yhtNLsO
zZoTE_($?DYv$dthUT2t(T~tZSG(vqT5<D7)y4^ONt$F+;S4lL&=li^=@aj^dd|}uu
K4f;>*G4>x6ig=a)

literal 0
HcmV?d00001

diff --git a/src/test/.utils.c.swp b/src/test/.utils.c.swp
new file mode 100644
index 0000000000000000000000000000000000000000..61f8f0a0a69a776daa26eed84cb892c58b0f0ff6
GIT binary patch
literal 12288
zcmeHN&yN&E7%fOdTR@E-Jb18J1KYDZyEBVolx4()g*f4t2r6MASkqm-)1jw&QeC~v
z2!wd^;?2awc=iwQpaJFTQ4bzGXgnZSLt;3YnE1Zx?%A1z4T*3esZJiV^L2e+y?XV&
z>fO+t<_GO#{EgWr!}T0v--wknKYDkWzrvWvL#ZkLfU7&`OC6``0mE3M8-@)n`)Sr#
z5lJQ+nGPC8W~PywILT&%2NiJlodQmQhfrWUdv*7mZS_|scJLRTfBzwBxm2fsQ@|<U
z6mSYS1)Ks-0jGdd;BiyH^f$3<AicRHezN@BwC=O~>OPzTP64NYQ@|<U6mSYS1)Ks-
z0jGddz$xGqcnlQ~A!BEsVC=>gkjMZ3SAYNCdX}-zfUCe|;5_g)Fa>M{o&xSX!`MyW
zQy>M-0v`ekz(3m=y9@jQ{0{sEd;(kpjsR`o9pKJ3#{L4X1K$DP0vc!m^T1}{##Y9@
z1O~uKU>`6KYyoaR&Dfv7ufTQSE8q%{0B3<Sz>C0bkFlSD8^9O9=RgE(1O9%Bv0s26
zfDrfycoTRTcnSFTNpJ(c1_poz7J&s|H_!xbLBF4Xo51$~>3JP+ms7wg;1qBQJaPrD
zGEYmBYsEWyL2EIXXI`YG^n1c2qGLlf#8%QcoT|p~7B6m}P7f4Wqr>YwjC)dLc$qD$
zbJbald-#mAwlZ>AYO01gf+WqPAH|8J0BCrg@4qKSGF2O$dMQ7T*@Q_WzmrFi)W`E4
zC3q6=z7<KRyFpI8E+>w-zs3XEjRrSbs4Pl#k7q{5YMIB%q@0M<xlkdZsg^u8;J`({
zuhYJcg^@f-17YZWJ<l5wTPl8Pt59Ad+D#&AmlovuU{`3wqGf4<;kq|So}Q-kGKq>@
z23ETei~Bnw2{Ej2eAljO(2mw5pB(>wPOT4^b%FY14Rx_#K0>_W4Q;d(Ay;{?BQ^I$
z#zm<SWqBJe>Ha+Za(v4Yy>g`@j_*WtJOH9&0T87*K?(0Y#w0t(CMkKzaFwc=O!kF_
zM>G#0vIq>`xFpH2t=#Xeo{B$~Ja0AUdA`r2MkH_~ex9C2O-d+N4`(gC{@QuSZG&s;
zgR%$~9Wdfe2-On_eBf6_r(`1Gd_xnGNh=94isNO?T%r%iQ$ea$q&D1iWserp8O-Q(
zHKV4|aUe5JeiVlK#M{d=RKt#wC9U)KR}mN=E1=sK3?G>q#D*?Mg^C@YYK?T}QcNqg
z3XuU~k|*9*bSr+c;?Ru{)#m1zT}*?n*?TX3`o!W1J0gm8W_Vv@8LjK0N>_Oo=m|Bb
z7}`5lH?$Gz+S4Y=ljJb$eO&aca);YT4^440Ii`=&zEu9`kUF2x6C?(D7Hzoh9;H=<
z?mJ*26h<^`0orOzG%`IV^lh{j6uY};8buOfP}QBDtCyQX18e10D=9--bsul%O6>CF
z>}<nMG|7FH8tSGb8E}Oih66$g9c%0|D||(*W)E8*=75AHuX08lpR`BABb32e@I%N2
z9og^?3MP3^WwuRY#g8qVKDuzAX5Chu@1@D4!sIyU@Mp>J><`i&I!G1-+KX^0)vbea
zWMg$}xjLW2QVYhIa4qEm*p2a~62&@-)s^K&Lu;jJX(VUu|1fH9z=*_tRFpg%`N~B8
z1jlva$T?-5Q$B!3@3oZ(aBPTRuosrU<cm(`ORZB_1g_e-%-^DakzdE%Lv*8l=-pFf
zA9^xXF&ahT@~Josx45l^(RR4~LSw4|NTA~)qFJj)sTID4UgtBsidD)}nCAAp{Z&$%
z9;M?sHZ+!k40tao*kD~lTi$c#mWK*$Zgr&C%!)CRo3r$gtVBH6mBB@N2=B~DrS>B6
zvndb+qvzzzo_XGApv8iie!+kn?elBMp=SGO%$*N!gF5w_<;Jm5Rv+GK4ZDN=2a{1Z
A1^@s6

literal 0
HcmV?d00001

diff --git a/src/test/Makefile b/src/test/Makefile
new file mode 100644
index 0000000..c4931ad
--- /dev/null
+++ b/src/test/Makefile
@@ -0,0 +1,37 @@
+CC = gcc
+CFLAGS = -g -std=gnu99 -Wno-implicit-function-declaration
+LDFLAGS = -Lobj/ -lvalse_core
+LIB = valse_core.so
+LIB_SRC = $(wildcard ../sources/*.c)
+LIB_OBJ = $(LIB_SRC:.c=.o)
+
+all: $(LIB) test.EMGLLF test.EMGrank test.constructionModelesLassoMLE test.EMGrank test.constructionModelesLassoRank test.selectionTotale
+
+$(LIB): $(LIB_OBJ)
+	$(CC) -o $@ $^
+
+test.EMGLLF: test.EMGLLF.o
+	$(CC) -o $@ $^ $(LDFLAGS)
+
+test.constructionModelesLassoMLE: test.constructionModelesLassoMLE.o
+	$(CC) -o $@ $^ $(LDFLAGS)
+
+test.EMGrank: test.EMGrank.o
+	$(CC) -o $@ $^ $(LDFLAGS)
+
+test.constructionModelesLassoRank: test.constructionModelesLassoRank.o
+	$(CC) -o $@ $^ $(LDFLAGS)
+
+test.selectionTotale: test.selectionTotale.o
+	$(CC) -o $@ $^ $(LDFLAGS)
+
+%.o: %.c
+	$(CC) -o $@ -c $< $(CFLAGS)
+
+clean:
+	rm -f *.o ../sources/*.o
+
+cclean: clean
+	rm -f $(LIB) $(TEST)
+
+.PHONY: all clean cclean
diff --git a/src/test/OLD_TEST_MATLAB/TODO b/src/test/OLD_TEST_MATLAB/TODO
new file mode 100644
index 0000000..d4b05cf
--- /dev/null
+++ b/src/test/OLD_TEST_MATLAB/TODO
@@ -0,0 +1,4 @@
+1) transformer les .m en .R
+2) sauvegarder les résultats + entrée/sorties à partir du code R uniquement
+3) dans test.Truc.c, prendre en entrée les paramètres de data/ juste sauvegardés,
+   puis comparer les sorties aux sorties enregistrées
diff --git a/src/test/OLD_TEST_MATLAB/basicInitParameters.m b/src/test/OLD_TEST_MATLAB/basicInitParameters.m
new file mode 100644
index 0000000..50410f5
--- /dev/null
+++ b/src/test/OLD_TEST_MATLAB/basicInitParameters.m
@@ -0,0 +1,19 @@
+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/src/test/OLD_TEST_MATLAB/checkOutput.m b/src/test/OLD_TEST_MATLAB/checkOutput.m
new file mode 100644
index 0000000..c56ed0b
--- /dev/null
+++ b/src/test/OLD_TEST_MATLAB/checkOutput.m
@@ -0,0 +1,11 @@
+function[]=checkOutput(varName, matrix, refMatrix, tol)
+
+	fprintf('Checking %s\n',varName);
+	maxError = max(max(max(max(abs(matrix - refMatrix)))));
+	if maxError >= tol
+		fprintf('    Inaccuracy: max(abs(error)) = %g >= %g\n',maxError,tol);
+	else
+		fprintf('    OK\n');
+	end
+
+end
diff --git a/src/test/OLD_TEST_MATLAB/covariance.m b/src/test/OLD_TEST_MATLAB/covariance.m
new file mode 100644
index 0000000..839c33c
--- /dev/null
+++ b/src/test/OLD_TEST_MATLAB/covariance.m
@@ -0,0 +1,9 @@
+%covariance matrix for tests on synthetic data: A(i,j) = a^|i-j|
+function[A] = covariance(p,a)
+
+	A = a*ones(p,p);
+	for i=1:p
+		A(i,:) = A(i,:) .^ abs(i-(1:p));
+	end
+
+end
diff --git a/src/test/OLD_TEST_MATLAB/generateIO.m b/src/test/OLD_TEST_MATLAB/generateIO.m
new file mode 100644
index 0000000..c677fd2
--- /dev/null
+++ b/src/test/OLD_TEST_MATLAB/generateIO.m
@@ -0,0 +1,37 @@
+%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/src/test/OLD_TEST_MATLAB/generateIOdefault.m b/src/test/OLD_TEST_MATLAB/generateIOdefault.m
new file mode 100644
index 0000000..f4c3c1f
--- /dev/null
+++ b/src/test/OLD_TEST_MATLAB/generateIOdefault.m
@@ -0,0 +1,23 @@
+%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/src/test/OLD_TEST_MATLAB/generateRunSaveTest_EMGLLF.m b/src/test/OLD_TEST_MATLAB/generateRunSaveTest_EMGLLF.m
new file mode 100644
index 0000000..bf0badf
--- /dev/null
+++ b/src/test/OLD_TEST_MATLAB/generateRunSaveTest_EMGLLF.m
@@ -0,0 +1,47 @@
+function[] = generateRunSaveTest_EMGLLF(n, p, m, k, mini, maxi, gamma, lambda, varargin)
+
+	%set defaults for optional inputs
+	optargs = {200 15 10 3 5 10 1.0 0.5};
+	%replace defaults by user parameters
+	optargs(1:length(varargin)) = varargin;
+	[n, p, m, k, mini, maxi, gamma, lambda] = optargs{:};
+	tau = 1e-6;
+	mini = int64(mini);
+	maxi = int64(maxi);
+
+	%Generate phiInit,piInit,...
+	[phiInit,rhoInit,piInit,gamInit] = basicInitParameters(n, p, m, k);
+
+	%Generate X and Y
+	[X, Y, ~] = generateIOdefault(n, p, m, k);
+
+	%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+	testFolder = 'data/';
+	mkdir(testFolder);
+	delimiter = '\n';
+
+	%save inputs
+	dlmwrite(strcat(testFolder,'phiInit'), reshape(phiInit,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'rhoInit'), reshape(rhoInit,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'piInit'), piInit, delimiter);
+	dlmwrite(strcat(testFolder,'gamInit'), reshape(gamInit,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'mini'), mini, delimiter);
+	dlmwrite(strcat(testFolder,'maxi'), maxi, delimiter);
+	dlmwrite(strcat(testFolder,'gamma'), gamma, delimiter);
+	dlmwrite(strcat(testFolder,'lambda'), lambda, delimiter);
+	dlmwrite(strcat(testFolder,'X'), reshape(X,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'Y'), reshape(Y,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'tau'), tau, delimiter);
+	dlmwrite(strcat(testFolder,'dimensions'), [n,p,m,k], delimiter);
+
+	[phi,rho,pi,LLF,S] = EMGLLF(phiInit,rhoInit,piInit,gamInit,mini,maxi,gamma,lambda,X,Y,tau);
+
+	%save output
+	dlmwrite(strcat(testFolder,'phi'), reshape(phi,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'rho'), reshape(rho,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'pi'), pi, delimiter);
+	dlmwrite(strcat(testFolder,'LLF'), LLF, delimiter);
+	dlmwrite(strcat(testFolder,'S'), reshape(S,1,[]), delimiter);
+
+end
diff --git a/src/test/OLD_TEST_MATLAB/generateRunSaveTest_constructionModelesLassoMLE.m b/src/test/OLD_TEST_MATLAB/generateRunSaveTest_constructionModelesLassoMLE.m
new file mode 100644
index 0000000..6e48d45
--- /dev/null
+++ b/src/test/OLD_TEST_MATLAB/generateRunSaveTest_constructionModelesLassoMLE.m
@@ -0,0 +1,62 @@
+function[] = generateRunSaveTest_constructionModelesLassoMLE(n, p, m, k, mini, maxi, gamma, glambda, varargin)
+	
+	%set defaults for optional inputs
+	optargs = {200 15 10 3 5 10 1.0 [0.0,0.01,0.02,0.03,0.05,0.1,0.2,0.3,0.5,0.7,0.85,0.99]};
+	%replace defaults by user parameters
+	optargs(1:length(varargin)) = varargin;
+	[n, p, m, k, mini, maxi, gamma, glambda] = optargs{:};
+	tau = 1e-6;
+	seuil = 1e-15;
+	mini = int64(mini);
+	maxi = int64(maxi);
+	L = length(glambda);
+
+	%Generate phiInit,piInit,...
+	[phiInit,rhoInit,piInit,gamInit] = basicInitParameters(n, p, m, k);
+
+	%Generate X and Y
+	[X, Y, ~] = generateIOdefault(n, p, m, k);
+
+	A2 = zeros(p,m+1,L,'int64');
+    for i=1:L
+        A2(:,1,i) = 1:p;
+        A2(1:5,2,i) = 1:5;
+    end
+	A1 = zeros(p,m+1,L,'int64');
+    for i=1:L
+        A1(:,1,i) = 1:p;
+        A1(1:5,2,i) = 1:5;
+    end
+    
+	%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+	testFolder = 'data/';
+	mkdir(testFolder);
+	delimiter = '\n';
+
+	%save inputs
+	dlmwrite(strcat(testFolder,'phiInit'), reshape(phiInit,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'rhoInit'), reshape(rhoInit,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'piInit'), piInit, delimiter);
+	dlmwrite(strcat(testFolder,'gamInit'), reshape(gamInit,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'mini'), mini, delimiter);
+	dlmwrite(strcat(testFolder,'maxi'), maxi, delimiter);
+	dlmwrite(strcat(testFolder,'gamma'), gamma, delimiter);
+	dlmwrite(strcat(testFolder,'glambda'), glambda, delimiter);
+	dlmwrite(strcat(testFolder,'X'), reshape(X,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'Y'), reshape(Y,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'seuil'), seuil, delimiter);
+	dlmwrite(strcat(testFolder,'tau'), tau, delimiter);
+	dlmwrite(strcat(testFolder,'A1'), reshape(A1,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'A2'), reshape(A2,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'dimensions'), [n,p,m,k,L], delimiter);
+
+	[phi,rho,pi,lvraisemblance] = constructionModelesLassoMLE(phiInit,rhoInit,piInit,gamInit,mini,maxi,gamma,glambda,X,Y,seuil,tau,A1,A2);
+
+	%save output
+	dlmwrite(strcat(testFolder,'phi'), reshape(phi,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'rho'), reshape(rho,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'pi'), reshape(pi,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'lvraisemblance'), reshape(lvraisemblance,1,[]), delimiter);
+
+end
diff --git a/src/test/OLD_TEST_MATLAB/generateRunSaveTest_selectiontotale.m b/src/test/OLD_TEST_MATLAB/generateRunSaveTest_selectiontotale.m
new file mode 100644
index 0000000..6caee15
--- /dev/null
+++ b/src/test/OLD_TEST_MATLAB/generateRunSaveTest_selectiontotale.m
@@ -0,0 +1,49 @@
+function[] = generateRunSaveTest_selectiontotale(n, p, m, k, mini, maxi, gamma, glambda, varargin)
+
+	%set defaults for optional inputs
+	optargs = {200 15 10 3 5 10 1.0 [0.0,0.01,0.02,0.03,0.05,0.1,0.2,0.3,0.5,0.7,0.85,0.99]};
+	%replace defaults by user parameters
+	optargs(1:length(varargin)) = varargin;
+	[n, p, m, k, mini, maxi, gamma, glambda] = optargs{:};
+	tau = 1e-6;
+	seuil = 1e-15;
+	mini = int64(mini);
+	maxi = int64(maxi);
+	L = length(glambda);
+
+	%Generate phiInit,piInit,...
+	[phiInit,rhoInit,piInit,gamInit] = basicInitParameters(n, p, m, k);
+	
+	%Generate X and Y
+	[X, Y, ~] = generateIOdefault(n, p, m, k);
+
+	%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+	testFolder = 'data/';
+	mkdir(testFolder);
+	delimiter = '\n';
+
+	%save inputs
+	dlmwrite(strcat(testFolder,'phiInit'), reshape(phiInit,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'rhoInit'), reshape(rhoInit,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'piInit'), piInit, delimiter);
+	dlmwrite(strcat(testFolder,'gamInit'), reshape(gamInit,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'mini'), mini, delimiter);
+	dlmwrite(strcat(testFolder,'maxi'), maxi, delimiter);
+	dlmwrite(strcat(testFolder,'gamma'), gamma, delimiter);
+	dlmwrite(strcat(testFolder,'glambda'), glambda, delimiter);
+	dlmwrite(strcat(testFolder,'X'), reshape(X,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'Y'), reshape(Y,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'seuil'), seuil, delimiter);
+	dlmwrite(strcat(testFolder,'tau'), tau, delimiter);
+	dlmwrite(strcat(testFolder,'dimensions'), [n,p,m,k,L], delimiter);
+
+	[A1,A2,Rho,Pi] = selectiontotale(phiInit,rhoInit,piInit,gamInit,mini,maxi,gamma,glambda,X,Y,seuil,tau);
+
+	%save output
+	dlmwrite(strcat(testFolder,'A1'), reshape(A1,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'A2'), reshape(A2,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'Rho'), reshape(Rho,1,[]), delimiter);
+	dlmwrite(strcat(testFolder,'Pi'), reshape(Pi,1,[]), delimiter);
+
+end
diff --git a/src/test/OLD_TEST_MATLAB/testConstructionModelesLassoMLE.m b/src/test/OLD_TEST_MATLAB/testConstructionModelesLassoMLE.m
new file mode 100644
index 0000000..27c1208
--- /dev/null
+++ b/src/test/OLD_TEST_MATLAB/testConstructionModelesLassoMLE.m
@@ -0,0 +1,46 @@
+function[] = testConstructionModelesLassoMLE()
+
+	testFolder = 'data/';
+	delimiter = '\n';
+
+	%get dimensions
+	dimensions = dlmread(strcat(testFolder,'dimensions'), delimiter);
+	n = dimensions(1);
+	p = dimensions(2);
+	m = dimensions(3);
+	k = dimensions(4);
+	L = dimensions(5);
+
+	%get all input arrays
+	phiInit = reshape(dlmread(strcat(testFolder,'phiInit'), delimiter), p, m, k);
+	rhoInit = reshape(dlmread(strcat(testFolder,'rhoInit'), delimiter), m, m, k);
+	piInit = transpose(dlmread(strcat(testFolder,'piInit'), delimiter));
+	gamInit = reshape(dlmread(strcat(testFolder,'gamInit'), delimiter), n, k);
+	mini = int64(dlmread(strcat(testFolder,'mini'), delimiter));
+	maxi = int64(dlmread(strcat(testFolder,'maxi'), delimiter));
+	gamma = dlmread(strcat(testFolder,'gamma'), delimiter);
+	glambda = dlmread(strcat(testFolder,'glambda'), delimiter);
+	X = reshape(dlmread(strcat(testFolder,'X'), delimiter), n, p);
+	Y = reshape(dlmread(strcat(testFolder,'Y'), delimiter), n, m);
+	seuil = dlmread(strcat(testFolder,'seuil'), delimiter);
+	tau = dlmread(strcat(testFolder,'tau'), delimiter);
+	A1 = int64(reshape(dlmread(strcat(testFolder,'A1'), delimiter), p, m+1, L));
+	A2 = int64(reshape(dlmread(strcat(testFolder,'A2'), delimiter), p, m+1, L));
+
+	%run constructionModelesLassoMLE.m
+	[phi,rho,pi,lvraisemblance] = constructionModelesLassoMLE(phiInit,rhoInit,piInit,gamInit,mini,maxi,gamma,glambda,X,Y,seuil,tau,A1,A2);
+
+	%get all stored outputs
+	ref_phi = reshape(dlmread(strcat(testFolder,'phi'), delimiter), p, m, k, L);
+	ref_rho = reshape(dlmread(strcat(testFolder,'rho'), delimiter), m, m, k, L);
+	ref_pi = reshape(dlmread(strcat(testFolder,'pi'), delimiter), k, L);
+	ref_lvraisemblance = reshape(dlmread(strcat(testFolder,'lvraisemblance'), delimiter), L, 2);
+
+	%check that output correspond to stored output
+	tol = 1e-5;
+	checkOutput('phi',phi,ref_phi,tol);
+	checkOutput('rho',rho,ref_rho,tol);
+	checkOutput('pi',pi,ref_pi,tol);
+	checkOutput('lvraisemblance',lvraisemblance,ref_lvraisemblance,tol);
+
+end
diff --git a/src/test/OLD_TEST_MATLAB/testEMGLLF.m b/src/test/OLD_TEST_MATLAB/testEMGLLF.m
new file mode 100644
index 0000000..0dd9db8
--- /dev/null
+++ b/src/test/OLD_TEST_MATLAB/testEMGLLF.m
@@ -0,0 +1,44 @@
+function[] = testEMGLLF()
+
+	testFolder = 'data/';
+	delimiter = '\n';
+
+	%get dimensions
+	dimensions = dlmread(strcat(testFolder,'dimensions'), delimiter);
+	n = dimensions(1);
+	p = dimensions(2);
+	m = dimensions(3);
+	k = dimensions(4);
+
+	%get all input arrays
+	phiInit = reshape(dlmread(strcat(testFolder,'phiInit'), delimiter), p, m, k);
+	rhoInit = reshape(dlmread(strcat(testFolder,'rhoInit'), delimiter), m, m, k);
+	piInit = transpose(dlmread(strcat(testFolder,'piInit'), delimiter));
+	gamInit = reshape(dlmread(strcat(testFolder,'gamInit'), delimiter), n, k);
+	mini = int64(dlmread(strcat(testFolder,'mini'), delimiter));
+	maxi = int64(dlmread(strcat(testFolder,'maxi'), delimiter));
+	gamma = dlmread(strcat(testFolder,'gamma'), delimiter);
+	lambda = dlmread(strcat(testFolder,'lambda'), delimiter);
+	X = reshape(dlmread(strcat(testFolder,'X'), delimiter), n, p);
+	Y = reshape(dlmread(strcat(testFolder,'Y'), delimiter), n, m);
+	tau = dlmread(strcat(testFolder,'tau'), delimiter);
+
+	%run EMGLLF.m
+	[phi,rho,pi,LLF,S] = EMGLLF(phiInit,rhoInit,piInit,gamInit,mini,maxi,gamma,lambda,X,Y,tau);
+
+	%get all stored outputs
+	ref_phi = reshape(dlmread(strcat(testFolder,'phi'), delimiter), p, m, k);
+	ref_rho = reshape(dlmread(strcat(testFolder,'rho'), delimiter), m, m, k);
+	ref_pi = dlmread(strcat(testFolder,'pi'), delimiter);
+	ref_LLF = dlmread(strcat(testFolder,'LLF'), delimiter);
+	ref_S = reshape(dlmread(strcat(testFolder,'S'), delimiter), p, m, k);
+
+	%check that output correspond to stored output
+	tol = 1e-5;
+	checkOutput('phi',phi,ref_phi,tol);
+	checkOutput('rho',rho,ref_rho,tol);
+	checkOutput('pi',pi,ref_pi,tol);
+	checkOutput('LLF',LLF,ref_LLF,tol);
+	checkOutput('S',S,ref_S,tol);
+
+end
diff --git a/src/test/OLD_TEST_MATLAB/testSelectiontotale.m b/src/test/OLD_TEST_MATLAB/testSelectiontotale.m
new file mode 100644
index 0000000..996017f
--- /dev/null
+++ b/src/test/OLD_TEST_MATLAB/testSelectiontotale.m
@@ -0,0 +1,44 @@
+function[] = testSelectiontotale()
+
+	testFolder = 'data/';
+	delimiter = '\n';
+
+	%get dimensions
+	dimensions = dlmread(strcat(testFolder,'dimensions'), delimiter);
+	n = dimensions(1);
+	p = dimensions(2);
+	m = dimensions(3);
+	k = dimensions(4);
+	L = dimensions(5);
+
+	%get all input arrays
+	phiInit = reshape(dlmread(strcat(testFolder,'phiInit'), delimiter), p, m, k);
+	rhoInit = reshape(dlmread(strcat(testFolder,'rhoInit'), delimiter), m, m, k);
+	piInit = transpose(dlmread(strcat(testFolder,'piInit'), delimiter));
+	gamInit = reshape(dlmread(strcat(testFolder,'gamInit'), delimiter), n, k);
+	mini = int64(dlmread(strcat(testFolder,'mini'), delimiter));
+	maxi = int64(dlmread(strcat(testFolder,'maxi'), delimiter));
+	gamma = dlmread(strcat(testFolder,'gamma'), delimiter);
+	glambda = dlmread(strcat(testFolder,'glambda'), delimiter);
+	X = reshape(dlmread(strcat(testFolder,'X'), delimiter), n, p);
+	Y = reshape(dlmread(strcat(testFolder,'Y'), delimiter), n, m);
+	seuil = dlmread(strcat(testFolder,'seuil'), delimiter);
+	tau = dlmread(strcat(testFolder,'tau'), delimiter);
+
+	%run constructionModelesLassoMLE.m
+	[A1,A2,Rho,Pi] = selectiontotale(phiInit,rhoInit,piInit,gamInit,mini,maxi,gamma,glambda,X,Y,seuil,tau);
+
+	%get all stored outputs
+	ref_A1 = int64(reshape(dlmread(strcat(testFolder,'A1'), delimiter), p, m+1, L));
+	ref_A2 = int64(reshape(dlmread(strcat(testFolder,'A2'), delimiter), p, m+1, L));
+	ref_Rho = reshape(dlmread(strcat(testFolder,'Rho'), delimiter), m, m, k, L);
+	ref_Pi = reshape(dlmread(strcat(testFolder,'Pi'), delimiter), k, L);
+
+	%check that output correspond to stored output
+	tol = 1e-5;
+	checkOutput('A1',A1,ref_A1,tol);
+	checkOutput('A2',A2,ref_A2,tol);
+	checkOutput('Rho',Rho,ref_Rho,tol);
+	checkOutput('Pi',Pi,ref_Pi,tol);
+
+end
diff --git a/src/test/test.ConstructionModelesLassoMLE.c b/src/test/test.ConstructionModelesLassoMLE.c
new file mode 100644
index 0000000..a6c0cce
--- /dev/null
+++ b/src/test/test.ConstructionModelesLassoMLE.c
@@ -0,0 +1,143 @@
+#include "constructionModelesLassoMLE.h"
+#include "ioutils.h"
+
+int main(int argc, char** argv)
+{
+	// read dimensions
+	const Int nbDims = 5;
+	Int* dimensions = readArray_int("dimensions",&nbDims,1);
+	mwSize n = dimensions[0];
+	mwSize p = dimensions[1];
+	mwSize m = dimensions[2];
+	mwSize k = dimensions[3];
+	mwSize L = dimensions[4];
+	free(dimensions);
+	mwSize lengthOne = 1;
+
+	////////////
+	// INPUTS //
+	////////////
+
+	// phiInit
+	const mwSize dimPhiInit[] = {p, m, k};
+	Real* phiInit = readArray_real("phiInit",dimPhiInit,3);
+
+	// rhoInit
+	const mwSize dimRhoInit[] = {m, m, k};
+	Real* rhoInit = readArray_real("rhoInit",dimRhoInit,3);
+
+	// piInit
+	Real* piInit = readArray_real("piInit",&k,1);
+
+	// gamInit
+	const mwSize dimGamInit[] = {n, k};
+	Real* gamInit = readArray_real("gamInit",dimGamInit,2);
+
+	// min number of iterations
+	Int* pmini = readArray_int("mini",&lengthOne,1);
+	Int mini = *pmini;
+	free(pmini);
+	
+	// max number of iterations
+	Int* pmaxi = readArray_int("maxi",&lengthOne,1);
+	Int maxi = *pmaxi;
+	free(pmaxi);
+
+	// gamma
+	Real* pgamma = readArray_real("gamma",&lengthOne,1);
+	Real gamma = *pgamma;
+	free(pgamma);
+
+	// lambda
+	Real* glambda = readArray_real("glambda",&L,1);
+
+	// X
+	const mwSize dimX[] = {n, p};
+	Real* X = readArray_real("X",dimX,2);
+
+	// Y
+	const mwSize dimY[] = {n, m};
+	Real* Y = readArray_real("Y",dimY,2);
+
+	// seuil
+	Real* pseuil = readArray_real("seuil",&lengthOne,1);
+	Real seuil = *pseuil;
+	free(pseuil);
+
+	// tau
+	Real* ptau = readArray_real("tau",&lengthOne,1);
+	Real tau = *ptau;
+	free(ptau);
+	
+	// A1
+	const mwSize dimA[] = {p, m+1, L};
+	Int* A1 = readArray_int("A1",dimA,3);
+	
+	// A2
+	Int* A2 = readArray_int("A2",dimA,3);
+	
+	/////////////
+	// OUTPUTS //
+	/////////////
+
+	// phi
+	const mwSize dimPhi[] = {dimPhiInit[0], dimPhiInit[1], dimPhiInit[2], L};
+	Real* phi = (Real*)malloc(dimPhi[0]*dimPhi[1]*dimPhi[2]*dimPhi[3]*sizeof(Real));
+
+	// rho
+	const mwSize dimRho[] = {dimRhoInit[0], dimRhoInit[1], dimRhoInit[2], L};
+	Real* rho = (Real*)malloc(dimRho[0]*dimRho[1]*dimRho[2]*dimRho[3]*sizeof(Real));
+
+	// pi
+	const mwSize dimPi[] = {k, L};
+	Real* pi = (Real*)malloc(dimPi[0]*dimPi[1]*sizeof(Real));
+
+	// lvraisemblance
+	const mwSize dimLvraisemblance[] = {L, 2};
+	Real* lvraisemblance = (Real*)malloc(dimLvraisemblance[0]*dimLvraisemblance[1]*sizeof(Real));
+
+	/////////////////////////////////////////
+	// Call to constructionModelesLassoMLE //
+	/////////////////////////////////////////
+
+	constructionModelesLassoMLE(
+		phiInit,rhoInit,piInit,gamInit,mini,maxi,gamma,glambda,X,Y,seuil,tau,A1,A2,
+		phi,rho,pi,lvraisemblance,
+		n,p,m,k,L);
+	
+	free(phiInit);
+	free(rhoInit);
+	free(piInit);
+	free(gamInit);
+	free(X);
+	free(Y);
+	free(A1);
+	free(A2);
+	free(glambda);
+	
+	// Compare to reference outputs
+	Real* ref_phi = readArray_real("phi",dimPhi,4);
+	compareArray_real("phi", phi, ref_phi, dimPhi[0]*dimPhi[1]*dimPhi[2]*dimPhi[3]);
+	free(phi);
+	free(ref_phi);
+	
+	// rho
+	Real* ref_rho = readArray_real("rho",dimRho,4);
+	compareArray_real("rho", rho, ref_rho, dimRho[0]*dimRho[1]*dimRho[2]*dimRho[3]);
+	free(rho);
+	free(ref_rho);
+	
+	// pi
+	Real* ref_pi = readArray_real("pi",dimPi,2);
+	compareArray_real("pi", pi, ref_pi, dimPi[0]*dimPi[1]);
+	free(pi);
+	free(ref_pi);
+	
+	// lvraisemblance
+	Real* ref_lvraisemblance = readArray_real("lvraisemblance",dimLvraisemblance,2);
+	compareArray_real("lvraisemblance", lvraisemblance, ref_lvraisemblance, dimLvraisemblance[0]*dimLvraisemblance[1]);
+	free(lvraisemblance);
+	free(ref_lvraisemblance);
+
+	return 0;
+}
diff --git a/src/test/test.EMGLLF.c b/src/test/test.EMGLLF.c
new file mode 100644
index 0000000..60516bc
--- /dev/null
+++ b/src/test/test.EMGLLF.c
@@ -0,0 +1,81 @@
+#include "EMGLLF.h"
+#include "ioutils.h"
+
+int main(int argc, char** argv)
+{
+	////////////
+	// INPUTS //
+	////////////
+
+	Int* dimensions = readArray_int("dimensions");
+	mwSize n = dimensions[0];
+	mwSize p = dimensions[1];
+	mwSize m = dimensions[2];
+	mwSize k = dimensions[3];
+	free(dimensions);
+
+	Real* phiInit = readArray_real("phiInit");
+	Real* rhoInit = readArray_real("rhoInit");
+	Real* piInit = readArray_real("piInit");
+	Real* gamInit = readArray_real("gamInit");
+	Int mini = read_int("mini");
+	Int maxi = read_int("maxi");
+	Real gamma = read_real("gamma");
+	Real lambda = read_real("lambda");
+	Real* X = readArray_real("X");
+	Real* Y = readArray_real("Y");
+	Real tau = read_real("tau");
+
+	/////////////
+	// OUTPUTS //
+	/////////////
+
+	Real* phi = (Real*)malloc(p*m*k*sizeof(Real));
+	Real* rho = (Real*)malloc(m*m*k*sizeof(Real));
+	Real* pi = (Real*)malloc(k*sizeof(Real));
+	Real* LLF = (Real*)malloc(maxi*sizeof(Real));
+	Real* S = (Real*)malloc(p*m*k*sizeof(Real));
+
+	////////////////////
+	// Call to EMGLLF //
+	////////////////////
+
+	EMGLLF_core(phiInit,rhoInit,piInit,gamInit,mini,maxi,gamma,lambda,X,Y,tau,
+		phi,rho,pi,LLF,S,
+		n,p,m,k);
+
+	free(phiInit);
+	free(rhoInit);
+	free(piInit);
+	free(gamInit);
+	free(X);
+	free(Y);
+
+	// Compare to reference outputs
+	Real* ref_phi = readArray_real("phi");
+	compareArray_real("phi", phi, ref_phi, p*m*k);
+	free(phi);
+	free(ref_phi);
+
+	Real* ref_rho = readArray_real("rho");
+	compareArray_real("rho", rho, ref_rho, m*m*k);
+	free(rho);
+	free(ref_rho);
+
+	Real* ref_pi = readArray_real("pi");
+	compareArray_real("pi", pi, ref_pi, k);
+	free(pi);
+	free(ref_pi);
+
+	Real* ref_LLF = readArray_real("LLF", maxi);
+	compareArray_real("LLF", LLF, ref_LLF);
+	free(LLF);
+	free(ref_LLF);
+
+	Real* ref_S = readArray_real("S");
+	compareArray_real("S", S, ref_S, p*m*k);
+	free(S);
+	free(ref_S);
+
+	return 0;
+}
diff --git a/src/test/test.EMGrank.c b/src/test/test.EMGrank.c
new file mode 100644
index 0000000..1227738
--- /dev/null
+++ b/src/test/test.EMGrank.c
@@ -0,0 +1,92 @@
+#include "EMGrank.h"
+#include "ioutils.h"
+
+int main(int argc, char** argv)
+{
+	// read dimensions
+	const Int nbDims = 4;
+	Int* dimensions = readArray_int("dimensions",&nbDims,1);
+	mwSize n = dimensions[0];
+	mwSize p = dimensions[1];
+	mwSize m = dimensions[2];
+	mwSize k = dimensions[3];
+	free(dimensions);
+	mwSize lengthOne = 1;
+
+	////////////
+	// INPUTS //
+	////////////
+
+	// Rho
+	const mwSize dimRho[] = {m, m, k};
+	Real* Rho = readArray_real("Rho",dimRho,3);
+
+	// Pi
+	Real* Pi = readArray_real("Pi",&k,1);
+
+	// min number of iterations
+	Int* pmini = readArray_int("mini",&lengthOne,1);
+	Int mini = *pmini;
+	free(pmini);
+
+	// max number of iterations
+	Int* pmaxi = readArray_int("maxi",&lengthOne,1);
+	Int maxi = *pmaxi;
+	free(pmaxi);
+
+	// X
+	const mwSize dimX[] = {n, p};
+	Real* X = readArray_real("X",dimX, 2);
+
+	// Y
+	const mwSize dimY[] = {n, m};
+	Real* Y = readArray_real("Y",dimY, 2);
+
+	// tau
+	Real* ptau = readArray_real("tau",&lengthOne,1);
+	Real tau = *ptau;
+	free(ptau);
+
+	// tau
+	Int* rank = readArray_int("rank",&k,1);
+
+	/////////////
+	// OUTPUTS //
+	/////////////
+
+	// phi
+	const mwSize dimPhi[] = {p, m, k};
+	Real* phi = (Real*)malloc(dimPhi[0]*dimPhi[1]*dimPhi[2]*sizeof(Real));
+
+	// LLF
+	Real* LLF = (Real*)malloc(1*sizeof(Real));
+
+	//////////////////////////
+	// Main call to EMGrank //
+	//////////////////////////
+
+	EMGrank(Pi,Rho,mini,maxi,X,Y,tau,rank,
+		phi,LLF,
+		n,p,m,k);
+	
+	// free input pointers
+	free(Rho);
+	free(Pi);
+	free(X);
+	free(Y);
+	free(rank);
+	
+	// Compare to reference outputs
+	Real* ref_phi = readArray_real("phi",dimPhi, 3);
+	compareArray_real("phi", phi, ref_phi, dimPhi[0]*dimPhi[1]*dimPhi[2]);
+	free(phi);
+	free(ref_phi);
+	
+	// LLF
+	Real* ref_LLF = readArray_real("LLF",&lengthOne,1);
+	compareArray_real("LLF", LLF, ref_LLF, 1);
+	free(LLF);
+	free(ref_LLF);
+	
+	return 0;
+}
diff --git a/src/test/test.constructionModelesLassoRank.c b/src/test/test.constructionModelesLassoRank.c
new file mode 100644
index 0000000..2232e9e
--- /dev/null
+++ b/src/test/test.constructionModelesLassoRank.c
@@ -0,0 +1,107 @@
+#include "constructionModelesLassoRank.h"
+#include "ioutils.h"
+
+int main(int argc, char** argv)
+{
+	// read dimensions
+	const Int nbDims = 5;
+	Int* dimensions = readArray_int("dimensions",&nbDims,1);
+	mwSize n = dimensions[0];
+	mwSize p = dimensions[1];
+	mwSize m = dimensions[2];
+	mwSize k = dimensions[3];
+	mwSize L = dimensions[4];
+	free(dimensions);
+	mwSize lengthOne = 1;
+
+	////////////
+	// INPUTS //
+	////////////
+
+	// piInit
+	const mwSize dimPi[] = {k, L};
+	Real* Pi = readArray_real("Pi",dimPi,2);
+
+	// rhoInit
+	const mwSize dimRho[] = {m, m, k, L};
+	Real* Rho = readArray_real("Rho",dimRho,4);
+
+	// min number of iterations
+	Int* pmini = readArray_int("mini",&lengthOne,1);
+	Int mini = *pmini;
+	free(pmini);
+
+	// max number of iterations
+	Int* pmaxi = readArray_int("maxi",&lengthOne,1);
+	Int maxi = *pmaxi;
+	free(pmaxi);
+
+	// X
+	const mwSize dimX[] = {n, p};
+	Real* X = readArray_real("X",dimX,2);
+
+	// Y
+	const mwSize dimY[] = {n, m};
+	Real* Y = readArray_real("Y",dimY,2);
+
+	// tau
+	Real* ptau = readArray_real("tau",&lengthOne,1);
+	Real tau = *ptau;
+	free(ptau);
+
+	// A1
+	const mwSize dimA[] = {p, L};
+	Int* A1 = readArray_int("A1",dimA,2);
+
+	// rangmin
+	Int* prangmin = readArray_int("rangmin",&lengthOne,1);
+	Int rangmin = *prangmin;
+	free(prangmin);
+
+	// rangmax
+	Int* prangmax = readArray_int("rangmax",&lengthOne,1);
+	Int rangmax = *prangmax;
+	free(prangmax);
+	
+	/////////////
+	// OUTPUTS //
+	/////////////
+
+	// phi
+	mwSize Size = (mwSize)pow(rangmax-rangmin+1, k);
+	const mwSize dimPhi[] = {p, m, k, L*Size};
+	Real* phi = (Real*)malloc(dimPhi[0]*dimPhi[1]*dimPhi[2]*dimPhi[3]*sizeof(Real));
+
+	// lvraisemblance
+	const mwSize dimLvraisemblance[] = {L*Size, 2};
+	Real* lvraisemblance = (Real*)malloc(dimLvraisemblance[0]*dimLvraisemblance[1]*sizeof(Real));
+
+	//////////////////////////////////////////////
+	// Main call to constructionModelesLassoMLE //
+	//////////////////////////////////////////////
+
+	constructionModelesLassoRank(
+		Pi,Rho,mini,maxi,X,Y,tau,A1,rangmin,rangmax,
+		phi,lvraisemblance,
+		n,p,m,k,L);
+	
+	free(Rho);
+	free(Pi);
+	free(X);
+	free(Y);
+	free(A1);
+	
+	// Compare to reference outputs
+	Real* ref_phi = readArray_real("phi",dimPhi, 4);
+	compareArray_real("phi", phi, ref_phi, dimPhi[0]*dimPhi[1]*dimPhi[2]*dimPhi[3]);
+	free(phi);
+	free(ref_phi);
+	
+	// lvraisemblance
+	Real* ref_lvraisemblance = readArray_real("lvraisemblance",dimLvraisemblance,2);
+	compareArray_real("lvraisemblance", lvraisemblance, ref_lvraisemblance, dimLvraisemblance[0]*dimLvraisemblance[1]);
+	free(lvraisemblance);
+	free(ref_lvraisemblance);
+	
+	return 0;
+}
diff --git a/src/test/test.selectiontotale.c b/src/test/test.selectiontotale.c
new file mode 100644
index 0000000..45cb461
--- /dev/null
+++ b/src/test/test.selectiontotale.c
@@ -0,0 +1,134 @@
+#include "selectiontotale.h"
+#include "ioutils.h"
+
+mwSize main(mwSize argc, char** argv)
+{
+	// read dimensions
+	const Int nbDims = 5;
+	Int* dimensions = readArray_int("dimensions",&nbDims,1);
+	mwSize n = dimensions[0];
+	mwSize p = dimensions[1];
+	mwSize m = dimensions[2];
+	mwSize k = dimensions[3];
+	mwSize L = dimensions[4];
+	free(dimensions);
+	mwSize lengthOne = 1;
+	
+	////////////
+	// INPUTS //
+	////////////
+
+	// phiInit
+	const mwSize dimPhiInit[] = {p, m, k};
+	Real* phiInit = readArray_real("phiInit",dimPhiInit,3);
+
+	// rhoInit
+	const mwSize dimRhoInit[] = {m, m, k};
+	Real* rhoInit = readArray_real("rhoInit",dimRhoInit,3);
+
+	// piInit
+	Real* piInit = readArray_real("piInit",&k,1);
+
+	// gamInit
+	const mwSize dimGamInit[] = {n, k};
+	Real* gamInit = readArray_real("gamInit",dimGamInit, 2);
+
+	// min number of iterations
+	Int* pmini = readArray_int("mini",&lengthOne,1);
+	Int mini = *pmini;
+	free(pmini);
+	
+	// max number of iterations
+	Int* pmaxi = readArray_int("maxi",&lengthOne,1);
+	Int maxi = *pmaxi;
+	free(pmaxi);
+	
+	// gamma
+	Real* pgamma = readArray_real("gamma",&lengthOne,1);
+	Real gamma = *pgamma;
+	free(pgamma);
+	
+	// lambda
+	Real* glambda = readArray_real("glambda",&L,1);
+
+	// X
+	const mwSize dimX[] = {n, p};
+	Real* X = readArray_real("X",dimX,2);
+
+	// Y
+	const mwSize dimY[] = {n, m};
+	Real* Y = readArray_real("Y",dimY,2);
+
+	// seuil
+	Real* pseuil = readArray_real("seuil",&lengthOne,1);
+	Real seuil = *pseuil;
+	free(pseuil);
+	
+	// tau
+	Real* ptau = readArray_real("tau",&lengthOne,1);
+	Real tau = *ptau;
+	free(ptau);
+	
+	/////////////
+	// OUTPUTS //
+	/////////////
+
+	// A1
+	const mwSize dimA[] = {p, m+1, L};
+	Int* A1 = (Int*)malloc(dimA[0]*dimA[1]*dimA[2]*sizeof(Int));
+
+	// A2
+	Int* A2 = (Int*)malloc(dimA[0]*dimA[1]*dimA[2]*sizeof(Int));
+
+	// rho
+	const mwSize dimRho[] = {m, m, k, L};
+	Real* Rho = (Real*)malloc(dimRho[0]*dimRho[1]*dimRho[2]*dimRho[3]*sizeof(Real));
+
+	// pi
+	const mwSize dimPi[] = {k, L};
+	Real* Pi = (Real*)malloc(dimPi[0]*dimPi[1]*sizeof(Real));
+
+	//////////////////////////////////////////////
+	// Main call to constructionModelesLassoMLE //
+	//////////////////////////////////////////////
+
+	selectiontotale(
+		phiInit,rhoInit,piInit,gamInit,mini,maxi,gamma,glambda,X,Y,seuil,tau,
+		A1,A2,Rho,Pi,
+		n,p,m,k,L);
+	
+	// free input pointers
+	free(phiInit);
+	free(rhoInit);
+	free(piInit);
+	free(gamInit);
+	free(glambda);
+	free(X);
+	free(Y);
+	
+	// Compare to reference outputs
+	Int* ref_A1 = readArray_int("A1",dimA, 3);
+	compareArray_int("A1", A1, ref_A1, dimA[0]*dimA[1]*dimA[2]);
+	free(A1);
+	free(ref_A1);
+	
+	// A2
+	Int* ref_A2 = readArray_int("A2",dimA, 3);
+	compareArray_int("A2", A2, ref_A2, dimA[0]*dimA[1]*dimA[2]);
+	free(A2);
+	free(ref_A2);
+	
+	// Rho
+	Real* ref_Rho = readArray_real("Rho",dimRho, 4);
+	compareArray_real("Rho", Rho, ref_Rho, dimRho[0]*dimRho[1]*dimRho[2]*dimRho[3]);
+	free(Rho);
+	free(ref_Rho);
+	
+	// Pi
+	Real* ref_Pi = readArray_real("Pi",dimPi, 2);
+	compareArray_real("Pi", Pi, ref_Pi, dimPi[0]*dimPi[1]);
+	free(Pi);
+	free(ref_Pi);
+	
+	return 0;
+}
diff --git a/src/test/utils.c b/src/test/utils.c
new file mode 100644
index 0000000..7cdf740
--- /dev/null
+++ b/src/test/utils.c
@@ -0,0 +1,81 @@
+// Check if array == refArray
+void compareArray(const char* ID, const void* array, const void* refArray, int size,
+	int isInteger)
+{
+	Real EPS = 1e-5; //precision
+	printf("Checking %s\n",ID);
+	Real maxError = 0.0;
+	for (mwSize i=0; i<size; i++)
+	{
+		Real error = isInteger
+			? fabs(((Int*)array)[i] - ((Int*)refArray)[i])
+			: fabs(((Real*)array)[i] - ((Real*)refArray)[i]);
+		if (error >= maxError)
+			maxError = error;
+	}
+	if (maxError >= EPS)
+		printf("    Inaccuracy: max(abs(error)) = %g >= %g\n",maxError,EPS);
+	else
+		printf("    OK\n");
+}
+
+void compareArray_real(const char* ID, const void* array, const void* refArray, int size)
+{
+	return compareArray(ID, array, refArray, size, 0);
+}
+
+void compareArray_int(const char* ID, const void* array, const void* refArray, int size)
+{
+	return compareArray(ID, array, refArray, size, 1);
+}
+
+// Read array by columns (as in MATLAB) and return by-rows encoding
+void* readArray(const char* fileName, int isInteger)
+{
+	// need to prepend '../data/' (not really nice code...)
+	char* fullFileName = (char*)calloc(5+strlen(fileName)+1,sizeof(char));
+	strcat(fullFileName, "data/");
+	strcat(fullFileName, fileName);
+	FILE* file = fopen(fullFileName, "r");
+	free(fullFileName);
+
+	// first pass to know how many elements to allocate
+	// /////................... TODO
+
+	int d = 1;
+	for (int i=0; i<nbDims; i++)
+		totalDim *= dimensions[i];
+	size_t elementSize = isInteger
+		? sizeof(Int)
+		: sizeof(Real);
+	
+	// read all values, and convert them to by-rows matrices format
+	void* matlabArray = malloc(totalDim*elementSize);
+	char curChar = ' ';
+	char bufferNum[64];
+	for (mwSize u=0; u<totalDim; u++)
+	{
+		// position to next non-separator character
+		while (!feof(file) && (curChar==' ' || curChar=='\n' || curChar=='\t' || curChar==','))
+			curChar = fgetc(file);
+		// read number (as a string)
+		int bufferIndex = 0;
+		while (!feof(file) && curChar!=' ' && curChar!='\n' && curChar!='\t' && curChar!=',')
+		{
+			bufferNum[bufferIndex++] = curChar;
+			curChar = fgetc(file);
+		}
+		bufferNum[bufferIndex] = 0;
+		// transform string into Real, and store it at appropriate location
+		if (isInteger)
+			((Int*)matlabArray)[u] = atoi(bufferNum);
+		else
+			((Real*)matlabArray)[u] = atof(bufferNum);
+	}
+	fclose(file);
+	
+	void* brArray = matlabToBrArray(matlabArray, dimensions, nbDims, isInteger);
+	free(matlabArray);
+	return brArray;
+}
+
-- 
2.44.0