Skip to content
Snippets Groups Projects
Commit c2f62fe5 authored by adeline.paiement's avatar adeline.paiement
Browse files

simplification pour estimer tous les X toujours

parent 6fcf549f
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
parmorceauxMax-git/Dynamics_Model.m
+ 37
− 128
View file @ c2f62fe5
......@@ -12,13 +12,13 @@ classdef Dynamics_Model < matlab.mixin.Copyable %handle
X
T
Y
previous_N_fixed
N_fixed
nb_times_fixed
%previous_N_fixed
%N_fixed
%nb_times_fixed
current_alpha
alpha_variable
confidence
%alpha_variable
%confidence
window_size
display_screen
......@@ -53,12 +53,12 @@ classdef Dynamics_Model < matlab.mixin.Copyable %handle
obj.T = [];
obj.X = [];
obj.Y = [];
obj.N_fixed = 0;
%obj.N_fixed = 0;
obj.current_alpha = alpha_ini;
obj.alpha_variable = [];
%obj.alpha_variable = [];
obj.window_size = [];
obj.confidence = [];
obj.nb_times_fixed = [];
%obj.confidence = [];
%obj.nb_times_fixed = [];
obj.initial_estimate_X = initial_X;
obj.alpha_initial = alpha_ini;
......@@ -74,33 +74,33 @@ classdef Dynamics_Model < matlab.mixin.Copyable %handle
alpha = obj.estimate_alpha();
sprintf('estimated alpha_variable: %f', alpha);
obj.current_alpha = alpha;
obj.window_size = [obj.window_size, frame - obj.N_fixed];
obj.window_size = [obj.window_size, frame];
% estimate X inside the window (and update N_fixed)
obj.previous_N_fixed = obj.N_fixed;
obj.estimate_X(alpha);
sprintf('"normalised" log likelihood of X (or confidence): %f', obj.confidence(frame));
%obj.previous_N_fixed = obj.N_fixed;
obj.estimate_X();
%sprintf('"normalised" log likelihood of X (or confidence): %f', obj.confidence(frame));
if obj.display_screen ~= -1
obj.plot_estimated_X();
end
% save the value of alpha for the frames that have converged
if obj.previous_N_fixed < obj.N_fixed
for i=obj.previous_N_fixed+1:obj.N_fixed
obj.alpha_variable(i) = alpha;
end
end
%if obj.previous_N_fixed < obj.N_fixed
% for i=obj.previous_N_fixed+1:obj.N_fixed
% obj.alpha_variable(i) = alpha;
% end
%end
end
function alpha = estimate_alpha(obj)
N = size(obj.X,1);
if N < 5
if N < 3
alpha = obj.alpha_initial;
return
end
first_index = max(obj.N_fixed, 1);
nb_points = N - first_index + 1;
first_index = 1;
nb_points = N;
if nb_points < 3 % we need at least 2 intervals to compute the mean
alpha = obj.alpha_initial;
......@@ -132,34 +132,17 @@ classdef Dynamics_Model < matlab.mixin.Copyable %handle
% end
end
end
function estimate_X(obj, alpha)
if obj.N_fixed > 0
last_fixed_X = obj.X(obj.N_fixed);
else
last_fixed_X = obj.initial_estimate_X;
end
function estimate_X(obj)
%function to be minimised over new_X, with extra parameters Y, alpha and T
f = @(variable_X)(-obj.logP_x_knowing_y_partialUpdate(last_fixed_X,variable_X,alpha,obj.sigma_estimation));
f = @(variable_X)(-obj.logP_x_knowing_y_paquet(variable_X,obj.current_alpha,obj.sigma_estimation));
N = size(obj.X,1);
if N > obj.N_fixed
X_ini = obj.X(obj.N_fixed+1:N);
else
X_ini = [];
end
if N == 0
new_X_ini = obj.initial_estimate_X;
elseif N == 1
new_X_ini = obj.X(N) + alpha;
X_ini = obj.initial_estimate_X;
else
new_X_ini = obj.X(N) + alpha * (obj.T(N+1) - obj.T(N));
X_ini = [obj.X; obj.X(N) + obj.current_alpha * (obj.T(N+1) - obj.T(N))];
end
% voir pb donc effacer X_ini=obj.normalise_X(X_ini)
X_ini = [X_ini; new_X_ini];
X_ini=obj.normalise_X(X_ini)
lb = repmat(0,1,length(X_ini));
ub = repmat(1,1,length(X_ini));
......@@ -168,66 +151,9 @@ classdef Dynamics_Model < matlab.mixin.Copyable %handle
Aeq = [];
beq = [];
[new_X_variable,loglikelihood] = fmincon(f,X_ini,A,b,Aeq,beq,lb,ub);
obj.confidence = [obj.confidence, -loglikelihood / size(new_X_variable,1)];
new_X_variable = obj.normalise_X(new_X_variable);
new_X = [obj.X(1:obj.N_fixed); new_X_variable];
% comparer X_estimated et X_0 et noter ceux qui ne bougent plus
% et ceux qui bougent encore
candidates_new_N_fixed = [];
for i = obj.N_fixed + 1 : N
while i > size(obj.nb_times_fixed,1)
obj.nb_times_fixed = [obj.nb_times_fixed; 0];
end
if abs(obj.X(i) - new_X(i)) < 1e-3%1e-2
obj.nb_times_fixed(i) = obj.nb_times_fixed(i)+1;
else
obj.nb_times_fixed(i) = 0;
end
if obj.nb_times_fixed(i) >= 2
candidates_new_N_fixed = [candidates_new_N_fixed i];
end
end
new_N_fixed = obj.N_fixed;
i = obj.N_fixed + 1;
for j=1:size(candidates_new_N_fixed,2)
if i ~= candidates_new_N_fixed(j)
break;
end
new_N_fixed = i;
i = i+1;
end
if new_N_fixed > N-3
new_N_fixed = N-3;
end
%
if new_N_fixed < 0
new_N_fixed = 0;
end
%
if new_N_fixed < N-13 %N-8
new_N_fixed = N-13; %N-8;
end
%obj.confidence = [obj.confidence, -loglikelihood / size(new_X_variable,1)];
if size(new_X, 1) < 4 %12
new_N_fixed = size(new_X, 1) - 2;
if new_N_fixed < 0
new_N_fixed = 0;
end
elseif size(new_X, 1) == 4 %12
new_N_fixed = 0;
end
obj.N_fixed = new_N_fixed;
obj.X = new_X;
obj.X = new_X_variable;
end
function spread_X(obj)
N = size(obj.X,1);
......@@ -296,26 +222,21 @@ classdef Dynamics_Model < matlab.mixin.Copyable %handle
obj.plot_estimated_X()
obj.compute_loglikelihood_dynamics_forced(1:N)
end
function loglikelihood = logP_x_knowing_y_partialUpdate(obj, last_fixed_X, variable_X, alpha, sigma)
function loglikelihood = logP_x_knowing_y_paquet(obj, variable_X, alpha, sigma)
N_variable = size(variable_X,1);
% variable_X = obj.normalise_X(variable_X);
% if N == 1, we only maximise P_cond_y(Y(1),X(1))
% if N == 1, we maximise P_cond_y(Y(1),X(1)) * P_x_Markov_simple(X(1),X(0)) avec X(0) virtuel
% if N == 2, we maximise P_cond_y(Y(2),X(2)) * P_x_Markov_simple(X(2),X(1)) * P_cond_y(Y(1),X(1)) with P_x_Markov_simple(X(2),X(1)) = 0 if delta_X is negative, and uniform proba otherwise
if obj.N_fixed == 0
loglikelihood = obj.logP_cond_y(obj.Y(1,:),variable_X(1));
else
loglikelihood = obj.logP_cond_y(obj.Y(obj.N_fixed+1,:),variable_X(1)) + obj.logP_x_Markov(variable_X(1),last_fixed_X,obj.T(obj.N_fixed+1),obj.T(obj.N_fixed),alpha,sigma);
end
loglikelihood = obj.logP_cond_y(obj.Y(1,:),variable_X(1));
loglikelihood = loglikelihood + obj.logP_x_Markov(variable_X(1), variable_X(1)-alpha, 2, 1, alpha, sigma);
if N_variable > 1
Y_2_to_N = obj.Y(obj.N_fixed+2:obj.N_fixed+N_variable,:);
Y_2_to_N = obj.Y(2:N_variable,:);
X_2_to_N = variable_X(2:N_variable);
X_1_to_Nminus1 = variable_X(1:N_variable-1);
T_2_to_N = obj.T(obj.N_fixed+2:obj.N_fixed+N_variable);
T_1_to_Nminus1 = obj.T(obj.N_fixed+1:obj.N_fixed+N_variable-1);
T_2_to_N = obj.T(2:N_variable);
T_1_to_Nminus1 = obj.T(1:N_variable-1);
ll_proba_cond_y = obj.logP_cond_y(Y_2_to_N, X_2_to_N);
ll_proba_x_Markov = obj.logP_x_Markov(X_2_to_N,X_1_to_Nminus1,T_2_to_N,T_1_to_Nminus1,alpha,sigma);
......@@ -376,28 +297,16 @@ classdef Dynamics_Model < matlab.mixin.Copyable %handle
else
loglikelihood = log(normpdf(difference,0,sigma));
index2 = loglikelihood == -Inf;%a été supprimé
index2 = loglikelihood < -100;%a été supprimé
loglikelihood(index2) = -100;%a été supprimé
end
end
function list_frames = get_converged_frames(obj)
if obj.previous_N_fixed < obj.N_fixed
list_frames = [obj.previous_N_fixed+1:obj.N_fixed];
else
list_frames = [];
end
end
function llh_dynamics = compute_loglikelihood_dynamics_forced(obj, frames)
function llh_dynamics = compute_loglikelihood_dynamics(obj, frames)
%llh_dynamics(frames <= 1) = nan;
indexes = frames > 1 & frames <= length(obj.T);
% if ~isempty(obj.alpha_variable)
% alpha = obj.alpha_variable(end);
% else
% alpha = obj.current_alpha;
% end
alpha = obj.current_alpha;
terme1 = obj.logP_cond_y(obj.Y(frames(indexes),:), obj.X(frames(indexes)));
......
parmorceauxMax-git/Interface.m
+ 7
− 7
View file @ c2f62fe5
......@@ -105,8 +105,8 @@ classdef Interface < handle
%% calcul les scores pour les frames du dernier morceau
% mise à jour des x pour bien les étaler
obj.meilleure_sequence_models(end).spread_X();
new_scores = obj.meilleure_sequence_models(end).compute_loglikelihood_dynamics_forced([1:obj.N-obj.first_frame_of_models(indice)+1]);
%obj.meilleure_sequence_models(end).spread_X();
new_scores = obj.meilleure_sequence_models(end).compute_loglikelihood_dynamics([1:obj.N-obj.first_frame_of_models(indice)+1]);
obj.loglikelihood = [obj.scores_finis{indice} new_scores];
%% à faire : afficher la séquence des modèles et la séquence de scores comme dans plot
......@@ -167,7 +167,7 @@ classdef Interface < handle
new_first_frame_of_models(morceau) = obj.first_frame_of_models(morceau);
new_models{morceau}.process_new_frame(Yt, Tt);
proba_obs = new_models{morceau}.compute_loglikelihood_dynamics_forced(1);
proba_obs = new_models{morceau}.compute_loglikelihood_dynamics(1);
new_proba_morceaux(morceau) = log(obj.proba_morceaux(morceau)) + proba_obs;
new_morceaux_finis{morceau} = [];
......@@ -206,7 +206,7 @@ classdef Interface < handle
%print branche et morceau
BRANCHE=branche
MORCEAU=morceau
proba_obs_branche(branche) = modele_branche{branche}.compute_loglikelihood_dynamics_forced(frame_number-first_frame_of_model_branche(branche)+1);
proba_obs_branche(branche) = modele_branche{branche}.compute_loglikelihood_dynamics(frame_number-first_frame_of_model_branche(branche)+1);
%% on garde en mémoire les morceaux qui sont finis, pour l'affichage final des X
%%%% %%%%%%%%%%%%%%%% a faire : mettre à jour les x pour bien les étaler
......@@ -215,17 +215,17 @@ classdef Interface < handle
else
nouveau_morceau_fini_branche{branche} = copy(obj.dynamics_models{branche}); % ce modèle n'a pas encore vu la frame courante, on peut donc le sauvegader dans son état à la frame précédente
%%%% %%%%%% mettre à jour la copie ici
nouveau_morceau_fini_branche{branche}.spread_X();
%nouveau_morceau_fini_branche{branche}.spread_X();
end
%% on récupère les scores des frames vues par le modèle
if branche == morceau
%% on récupère les scores pour les frames précédentes, qui ont pu bouger avec le argmax
proba_morceau_prec(branche) = sum(modele_branche{branche}.compute_loglikelihood_dynamics_forced([1:frame_number-first_frame_of_model_branche(branche)]));
proba_morceau_prec(branche) = sum(modele_branche{branche}.compute_loglikelihood_dynamics([1:frame_number-first_frame_of_model_branche(branche)]));
proba_morceau_prec(branche) = proba_morceau_prec(branche) + log(trans) * (frame_number-first_frame_of_model_branche(branche)-1);
else
%% on récupère les scores pour les frames précédentes, qui ont pu bouger avec le argmax
new_scores_branche{branche} = nouveau_morceau_fini_branche{branche}.compute_loglikelihood_dynamics_forced([1:frame_number-obj.first_frame_of_models(branche)]);
new_scores_branche{branche} = nouveau_morceau_fini_branche{branche}.compute_loglikelihood_dynamics([1:frame_number-obj.first_frame_of_models(branche)]);
proba_morceau_prec(branche) = sum(new_scores_branche{branche});
trans_morceau_prec = obj.transition(branche, branche);
proba_morceau_prec(branche) = proba_morceau_prec(branche) + log(trans_morceau_prec) * (frame_number-obj.first_frame_of_models(branche)-1);
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment