// OUTPUTS //
        /////////////
 
-       SEXP phi, rho, pi, LLF, S, dimPhiS, dimRho;
+       SEXP phi, rho, pi, LLF, S, affec, dimPhiS, dimRho;
        PROTECT(dimPhiS = allocVector(INTSXP, 3));
        int* pDimPhiS = INTEGER(dimPhiS);
        pDimPhiS[0] = p; pDimPhiS[1] = m; pDimPhiS[2] = k;
        PROTECT(pi = allocVector(REALSXP, k));
        PROTECT(LLF = allocVector(REALSXP, maxi-mini+1));
        PROTECT(S = allocArray(REALSXP, dimPhiS));
+       PROTECT(affec = allocVector(INTSXP, n));
        double *pPhi=REAL(phi), *pRho=REAL(rho), *pPi=REAL(pi), *pLLF=REAL(LLF), *pS=REAL(S);
+       int *pAffec=INTEGER(affec);
 
        ////////////////////
        // Call to EMGLLF //
        ////////////////////
 
        EMGLLF_core(phiInit,rhoInit,piInit,gamInit,mini,maxi,gamma,lambda,X,Y,tau,
-               pPhi,pRho,pPi,pLLF,pS,
+               pPhi,pRho,pPi,pLLF,pS,pAffec,
                n,p,m,k);
 
        // Build list from OUT params and return it
        SEXP listParams, listNames;
-       PROTECT(listParams = allocVector(VECSXP, 5));
-       char* lnames[5] = {"phi", "rho", "pi", "LLF", "S"}; //lists labels
-       PROTECT(listNames = allocVector(STRSXP,5));
-       for (int i=0; i<5; i++)
+       int nouts = 6;
+       PROTECT(listParams = allocVector(VECSXP, nouts));
+       char* lnames[6] = {"phi", "rho", "pi", "LLF", "S", "affec"}; //lists labels
+       PROTECT(listNames = allocVector(STRSXP,nouts));
+       for (int i=0; i<nouts; i++)
                SET_STRING_ELT(listNames,i,mkChar(lnames[i]));
        setAttrib(listParams, R_NamesSymbol, listNames);
        SET_VECTOR_ELT(listParams, 0, phi);
        SET_VECTOR_ELT(listParams, 2, pi);
        SET_VECTOR_ELT(listParams, 3, LLF);
        SET_VECTOR_ELT(listParams, 4, S);
+       SET_VECTOR_ELT(listParams, 5, affec);
 
-       UNPROTECT(9);
+       UNPROTECT(10);
        return listParams;
 }
 
        Real* pi, // parametre des proportions renormalisé, calculé par l'EM
        Real* LLF, // log vraisemblance associée à cet échantillon, pour les valeurs estimées des paramètres
        Real* S,
+       int* affec,
        // additional size parameters
        int n, // nombre d'echantillons
        int p, // nombre de covariables
                ite++;
        }
 
-       //TODO: affec = apply(gam, 1,which.max) à traduire, et retourner affec aussi.
+       //affec = apply(gam, 1, which.max)
+       for (int i=0; i<n; i++)
+       {
+               Real rowMax = 0.;
+               affec[i] = 0;
+               for (int j=0; j<k; j++)
+               {
+                       if (gam[mi(i,j,n,k)] > rowMax)
+                       {
+                               affec[i] = j+1; //R indices start at 1
+                               rowMax = gam[mi(i,j,n,k)];
+                       }
+               }
+       }
 
        //free memory
        free(b);
 
        write.table(as.double(res$pi), paste(testFolder,"pi",sep=""), row.names=F, col.names=F)
        write.table(as.double(res$LLF), paste(testFolder,"LLF",sep=""), row.names=F, col.names=F)
        write.table(as.double(res$S), paste(testFolder,"S",sep=""), row.names=F, col.names=F)
+       write.table(as.integer(res$affec), paste(testFolder,"affec",sep=""), row.names=F, col.names=F)
 }
 
        Real* pi = (Real*)malloc(k*sizeof(Real));
        Real* LLF = (Real*)malloc(maxi*sizeof(Real));
        Real* S = (Real*)malloc(p*m*k*sizeof(Real));
+       int* affec = (int*)malloc(n*sizeof(int));
        /////////////
 
        ////////////////////
        // Call to EMGLLF //
        EMGLLF_core(phiInit,rhoInit,piInit,gamInit,mini,maxi,gamma,lambda,X,Y,tau,
-               phi,rho,pi,LLF,S,
+               phi,rho,pi,LLF,S,affec,
                n,p,m,k);
        ////////////////////
 
        free(S);
        free(ref_S);
 
+       int* ref_affec = readArray_int("affec");
+       compareArray_int("affec", affec, ref_affec, n);
+       free(affec);
+       free(ref_affec);
+
        return 0;
 }