Correction d'une erreur dans la génération des données manquantes

MATLAB/calib_meth/EMNeNMF/emnenmf.m

@@ -3,8 +3,8 @@ function [ T , RMSE ] = emnenmf( W , X , G , F , Omega_G, Omega_F, Phi_G, Phi_F
 X0 =X;
 Omega_G = (Omega_G == 1); % Logical mask is faster than indexing in matlab.
 Omega_F = (Omega_F == 1); % Logical mask is faster than indexing in matlab.
-% nOmega_G = ~Omega_G; % Logical mask is faster than indexing in matlab.
-% nOmega_F = ~Omega_F; % Logical mask is faster than indexing in matlab.
+nOmega_G = ~Omega_G; % Logical mask is faster than indexing in matlab.
+nOmega_F = ~Omega_F; % Logical mask is faster than indexing in matlab.
 
 [~, num_sensor] = size(F);
 num_sensor = num_sensor-1;
@@ -38,7 +38,7 @@ while toc<Tmax
         XF = X*F' - Phi_G*FF;
 
         G(Omega_G) = 0; % Convert G to \Delta G
-        [ G , iterG , ~ ] = MaJ_G_EM_NeNMF( G , FF , XF , InnerMaxIter , StoppingCritG , Omega_G); % Update \Delta G
+        [ G , iterG ] = MaJ_G_EM_NeNMF( FF , XF , G , InnerMinIter , InnerMaxIter , StoppingCritG , nOmega_G); % Update \Delta G
         G(Omega_G) = Phi_G(Omega_G); % Convert \Delta G to G
         niter = niter + iterG;
         if(iterG<=InnerMinIter)
@@ -49,7 +49,7 @@ while toc<Tmax
         GX = G'*X-GG*Phi_F;
 
         F(Omega_F) = 0; % Convert F to \Delta F
-        [ F , iterF ] = MaJ_F_EM_NeNMF( GG , GX , F , InnerMaxIter , StoppingCritF , Omega_F); % Update \Delta F
+        [ F , iterF ] = MaJ_F_EM_NeNMF( GG , GX , F , InnerMinIter , InnerMaxIter , StoppingCritF , nOmega_F); % Update \Delta F
         F(Omega_F) = Phi_F(Omega_F); % Convert \Delta F to F
         niter = niter + iterF;
         if(iterF<=InnerMinIter)
@@ -65,9 +65,9 @@ while toc<Tmax
             RMSE(:,i) = vecnorm(F(:,1:end-1) - F_theo(:,1:end-1),2,2)/sqrt(num_sensor);
         end
     end
-    
+  
 end
-
+niter  
 end
 
 

MATLAB/calib_meth/REMNeNMF/remnenmf.m

@@ -4,8 +4,8 @@ r = 0;
 X0 = X;
 Omega_G = (Omega_G == 1); % Logical mask is faster than indexing in matlab.
 Omega_F = (Omega_F == 1); % Logical mask is faster than indexing in matlab.
-% nOmega_G = ~Omega_G; % Logical mask is faster than indexing in matlab.
-% nOmega_F = ~Omega_F; % Logical mask is faster than indexing in matlab.
+nOmega_G = ~Omega_G; % Logical mask is faster than indexing in matlab.
+nOmega_F = ~Omega_F; % Logical mask is faster than indexing in matlab.
 
 [~, num_sensor] = size(F);
 num_sensor = num_sensor-1;
@@ -54,9 +54,9 @@ while toc<Tmax
         XF = X_R * F_R' - Phi_G * FF;
 
         G(Omega_G) = 0; % Convert G to \Delta G
-        [ G , iterG , ~ ] = MaJ_G_EM_NeNMF( G , FF , XF , InnerMaxIter , StoppingCritG , Omega_G); % Update \Delta G
+        [ G , iterG ] = MaJ_G_EM_NeNMF( FF , XF , G , InnerMinIter , InnerMaxIter , StoppingCritG , nOmega_G); % Update \Delta G
         G(Omega_G) = Phi_G(Omega_G); % Convert \Delta G to G
-
+        niter = niter + iterG;
         if(iterG<=InnerMinIter)
             StoppingCritG = 1.e-1*StoppingCritG;
         end
@@ -66,9 +66,9 @@ while toc<Tmax
         GX = G_L' * X_L - GG * Phi_F;
 
         F(Omega_F) = 0; % Convert F to \Delta F
-        [ F , iterF ] = MaJ_F_EM_NeNMF( GG , GX , F , InnerMaxIter , StoppingCritF , Omega_F); % Update \Delta F
+        [ F , iterF ] = MaJ_F_EM_NeNMF( GG , GX , F , InnerMinIter , InnerMaxIter , StoppingCritF , nOmega_F); % Update \Delta F
         F(Omega_F) = Phi_F(Omega_F); % Convert \Delta F to F
-
+        niter = niter + iterF;
         if(iterF<=InnerMinIter)
             StoppingCritF = 1.e-1*StoppingCritF;
         end
@@ -82,7 +82,7 @@ while toc<Tmax
             RMSE(:,i) = vecnorm(F(:,1:end-1) - F_theo(:,1:end-1),2,2)/sqrt(num_sensor);
         end
     end
-    niter = niter + 1;
+    
 end
 niter
 end
@@ -97,7 +97,7 @@ function [ L,R ] = RSI_compression(X,r)
 %             Date: 13/04/2018
 %
 
-compressionLevel=10;
+compressionLevel=20;
 [m,n]=size(X);
 
 l = min(n, max(compressionLevel, r + 10));

MATLAB/calib_meth/common_nesterov/MaJ_F_EM_NeNMF.m

@@ -1,4 +1,4 @@
-function [ F1 , i ] = MaJ_F_EM_NeNMF( GG , GX , F1 , InnerMaxIter , StoppingCritF , Omega_F)
+function [ F1 , i ] = MaJ_F_EM_NeNMF( GG , GX , F1 , InnerMinIter , InnerMaxIter , StoppingCritF , nOmega_F)
 
 Y = F1;
 
@@ -11,8 +11,8 @@ L = norm(GG);
 
 for i = 1 : InnerMaxIter
     
-    F2 = max(Y-(1/L)*Grad_y,0);
-    F2(Omega_F) = 0; % Projection on not(Omega_F)
+    F2 = max(Y-(1/L)*Grad_y,0).*nOmega_F;
+%     F2 = F2.*nOmega_F; % Projection on not(Omega_F)
     alpha(2) = (1+sqrt(4*alpha(1)^2+1))/2;
     
     Y = F2 + ((alpha(1)-1)/alpha(2))*(F2-F1);
@@ -21,10 +21,11 @@ for i = 1 : InnerMaxIter
     F1 = F2;
     alpha(1) = alpha(2);
     
-    if(Grad_P(Grad_y , Y)<=StoppingCritF)
-        break
+    if mod(i , InnerMinIter) == 0
+        if(Grad_P(Grad_y , Y)<=StoppingCritF)
+            break
+        end
     end
     
 end
-
 end

MATLAB/calib_meth/common_nesterov/MaJ_G_EM_NeNMF.m

@@ -1,4 +1,4 @@
-function [ G1 , i , GradG ] = MaJ_G_EM_NeNMF( G1 , FF , XF ,  InnerMaxIter , StoppingCritG , Omega_G)
+function [ G1 , i ] = MaJ_G_EM_NeNMF( FF , XF , G1 , InnerMinIter , InnerMaxIter , StoppingCritG , nOmega_G)
 
 Y = G1;
 
@@ -11,8 +11,8 @@ Grad_y = Y*FF-XF;
 
 for i = 1 : InnerMaxIter
     
-    G2 = max(Y-(1/L)*Grad_y,0);
-    G2(Omega_G) = 0;
+    G2 = max(Y-(1/L)*Grad_y,0).*nOmega_G;
+%     G2 = G2.*nOmega_G; % Projection on not(Omega_G)
     alpha(2) = (1+sqrt(4*alpha(1)^2+1))/2;
     
     Y = G2 + ((alpha(1)-1)/alpha(2))*(G2-G1);
@@ -21,10 +21,12 @@ for i = 1 : InnerMaxIter
     G1 = G2;
     alpha(1) = alpha(2);
     
-    if(Grad_P(Grad_y , Y)<=StoppingCritG)
-        break
+    if mod(i , InnerMinIter) == 0
+        if(Grad_P(Grad_y , Y)<=StoppingCritG)
+            break
+        end
     end
     
 end
-GradG = G1*FF-XF;
-end
+end
+

MATLAB/common/data_gen.m

@@ -50,7 +50,7 @@ for i = 1 : N_RV
     W(idx_RefRV(i),idx_CptRV(i)) = 1;
 end
 
-N_data = round((1-MV)*(N_Cpt)*(100-N_Ref)); % Nb. of measurements in data matrix X
+N_data = round((1-MV)*(N_Cpt)*(s_n-N_Ref)); % Nb. of measurements in data matrix X
 xCpt = 1 : s_n;
 xCpt(idx_Ref) = []; % Reference free locations
 [xx,yy] = meshgrid(xCpt,1:N_Cpt); % Possibly sensed locations

MATLAB/main_emnenmf_vs_remnenmf.m

@@ -0,0 +1,113 @@
+clear
+close all
+clc
+
+%% Adding every files in the path
+
+addpath(genpath(pwd))
+
+%% Simulation parameters
+
+s_width = 20;  % Scene width
+s_length = 20; % Scene length
+N_Ref = 4; % Nb. of reference measurements
+N_Cpt = 25; % Nb. of mobile sensors
+Mu_beta = .9; % Mean sensors gain
+Mu_alpha = 5; % Mean sensors offset
+Bound_beta = [.01;1.5]; % Gain boundaries
+Bound_alpha = [3.5;6.5]; % Offset boundaries
+
+MV = .5; % Missing Value prop.
+RV = 0.3; % RendezVous prop.
+var_n = 0; % Noise variance
+M_loop = 50; % Number of M loop per E step
+runs = 1; % Total number of runs
+
+%% Nesterov parameters
+
+InnerMaxIter = 300;
+InnerMinIter = floor(InnerMaxIter/10);
+Tmax = 30;
+
+%% 
+delta_measure = 1;
+iter_max = round(Tmax / delta_measure);
+
+
+%% Allocation for the RMSE values
+% remnenmf
+RMSE_offset_remnenmf = nan(runs, iter_max);
+RMSE_gain_remnenmf = nan(runs, iter_max);
+% emnenmf
+RMSE_offset_emnenmf = nan(runs, iter_max);
+RMSE_gain_emnenmf = nan(runs, iter_max);
+
+for run = 1:runs
+    % data generation
+    [X, X_theo, W, F_theo, Omega_G, Omega_F, Phi_G, Phi_F, Ginit, Finit] = data_gen(s_width, s_length, run, N_Ref, N_Cpt, Mu_beta, Mu_alpha, Bound_beta, Bound_alpha, MV, RV, var_n);
+    % remnenmf
+    [T_remnenmf, RMSE] = remnenmf( W , X , Ginit , Finit, Omega_G, Omega_F, Phi_G, Phi_F , InnerMinIter , InnerMaxIter , Tmax , M_loop, F_theo, delta_measure);
+    RMSE_offset_remnenmf(run,:) = RMSE(1,:);
+    RMSE_gain_remnenmf(run,:) = RMSE(2,:);
+    % emnenmf
+    [T_emnenmf, RMSE] = emnenmf( W , X , Ginit , Finit, Omega_G, Omega_F, Phi_G, Phi_F , InnerMinIter , InnerMaxIter , Tmax , M_loop, F_theo, delta_measure);
+    RMSE_offset_emnenmf(run,:) = RMSE(1,:);
+    RMSE_gain_emnenmf(run,:) = RMSE(2,:);    
+    
+end
+
+% remnenmf
+min_offset_remnenmf = min(RMSE_offset_remnenmf,[],1,'omitnan');
+med_offset_remnenmf = median(RMSE_offset_remnenmf,1,'omitnan');
+max_offset_remnenmf = max(RMSE_offset_remnenmf,[],1,'omitnan');
+
+min_gain_remnenmf = min(RMSE_gain_remnenmf,[],1,'omitnan');
+med_gain_remnenmf = median(RMSE_gain_remnenmf,1,'omitnan');
+max_gain_remnenmf = max(RMSE_gain_remnenmf,[],1,'omitnan');
+
+subplot(121)
+semilogy(T_remnenmf,min_offset_remnenmf,'b')
+hold on
+semilogy(T_remnenmf,med_offset_remnenmf,'b')
+o_e = semilogy(T_remnenmf,max_offset_remnenmf,'b');
+hold off
+
+subplot(122)
+semilogy(T_remnenmf,min_gain_remnenmf,'b')
+hold on
+semilogy(T_remnenmf,med_gain_remnenmf,'b')
+g_e = semilogy(T_remnenmf,max_gain_remnenmf,'b');
+hold off
+
+% emnenmf
+min_offset_emnenmf = min(RMSE_offset_emnenmf,[],1,'omitnan');
+med_offset_emnenmf = median(RMSE_offset_emnenmf,1,'omitnan');
+max_offset_emnenmf = max(RMSE_offset_emnenmf,[],1,'omitnan');
+
+min_gain_emnenmf = min(RMSE_gain_emnenmf,[],1,'omitnan');
+med_gain_emnenmf = median(RMSE_gain_emnenmf,1,'omitnan');
+max_gain_emnenmf = max(RMSE_gain_emnenmf,[],1,'omitnan');
+
+subplot(121)
+hold on
+semilogy(T_emnenmf,min_offset_emnenmf,'r')
+semilogy(T_emnenmf,med_offset_emnenmf,'r')
+o_i = semilogy(T_emnenmf,max_offset_emnenmf,'r');
+hold off
+
+subplot(122)
+hold on
+semilogy(T_emnenmf,min_gain_emnenmf,'r')
+semilogy(T_emnenmf,med_gain_emnenmf,'r')
+g_i = semilogy(T_emnenmf,max_gain_emnenmf,'r');
+hold off
+
+% adding title and labels
+subplot(121)
+title('RMSE offset')
+legend([o_e o_i],'REMNeNMF','EMNeNMF')
+
+subplot(122)
+title('RMSE gain')
+legend([g_e g_i],'REMNeNMF','EMNeNMF')
+

MATLAB/main_incal_vs_emnenmf.m

@@ -8,8 +8,8 @@ addpath(genpath(pwd))
 
 %% Simulation parameters
 
-s_width = 100;  % Scene width
-s_length = 100; % Scene length
+s_width = 50;  % Scene width
+s_length = 50; % Scene length
 N_Ref = 4; % Nb. of reference measurements
 N_Cpt = 25; % Nb. of mobile sensors
 Mu_beta = .9; % Mean sensors gain
@@ -24,9 +24,9 @@ M_loop = 5; % Number of M loop per E step
 runs = 1; % Total number of runs
 
 %% Nesterov parameters
-InnerMinIter = 5;
+InnerMinIter = 10;
 InnerMaxIter = 100;
-Tmax = 300;
+Tmax = 30;
 
 %% 
 delta_measure = 1;
@@ -110,3 +110,5 @@ subplot(122)
 title('RMSE gain')
 legend([g_e g_i],'EMNeNMF','IN\_Cal')
 
+med_offset_incal(end)
+med_offset_emnenmf(end)

MATLAB/main_incal_vs_remnenmf.m

@@ -8,8 +8,8 @@ addpath(genpath(pwd))
 
 %% Simulation parameters
 
-s_width = 100;  % Scene width
-s_length = 100; % Scene length
+s_width = 50;  % Scene width
+s_length = 50; % Scene length
 N_Ref = 4; % Nb. of reference measurements
 N_Cpt = 25; % Nb. of mobile sensors
 Mu_beta = .9; % Mean sensors gain
@@ -20,13 +20,13 @@ Bound_alpha = [3.5;6.5]; % Offset boundaries
 MV = .5; % Missing Value prop.
 RV = 0.3; % RendezVous prop.
 var_n = 0; % Noise variance
-M_loop = 20; % Number of M loop per E step
+M_loop = 50; % Number of M loop per E step
 runs = 1; % Total number of runs
 
 %% Nesterov parameters
 InnerMinIter = 5;
 InnerMaxIter = 100;
-Tmax = 300;
+Tmax = 30;
 
 %% 
 delta_measure = 1;