root/applications/bdmtoolbox/tutorial/userguide/mixef_basic.m @ 1024

clear all

com.class='ARX';
com.yrv = RV({'d'},[2],[0]);
com.rgr = RV({},[],[]);

Data = [randn(2,50) randn(2,50)+[5*ones(1,50); 10*ones(1,50)]];

[Mix0,P0]=mixef_init(Data,com,4);
% show predictor

Pred = bm_epredictor(Mix0);
figure(1);
hold off
plot(Data(1,:), Data(2,:),'.');
hold on
epdf_2dplot(Pred);
hold on
for i=1:length(Pred.pdfs)
    mea=epdf_mean(Pred.pdfs{i});
    plot(mea(1),mea(2),'+r');
end

% Do Quasi Bayes

tic
%use same or new data
MixQB = Mix0; 
n = length(MixQB.Coms);
for t=1:size(Data,2)
    log_w_nn = zeros(1,n);
    % get predictions
    for c=1:n
        yt = Data(:,t);
        [dummy, log_w_nn(c)]=bm_bayes(MixQB.Coms{c}, yt);
    end
    % normalize weights
    w = exp(log_w_nn-max(log_w_nn));
    w = w/sum(w);
    
    for c=1:n
        yt = Data(:,t);
        [MixQB.Coms{c}]=bm_bayesweighted(MixQB.Coms{c}, yt, [], w(c));
    end
    MixQB.weights = bm_bayes(MixQB.weights, w);
end
toc

tic
% should be the same as:
MixQBc = bm_bayes_batch(MixQB,Data);
toc

%% display results

PredQB = bm_epredictor(MixQB);
figure(2);
hold off
plot(Data(1,:), Data(2,:),'.');
hold on
epdf_2dplot(PredQB);

PredQBc = bm_epredictor(MixQBc);
figure(3);
hold off
plot(Data(1,:), Data(2,:),'.');
hold on
epdf_2dplot(PredQBc);

%% relations between components
Batta_dist=zeros(n);
for i=1:n
    for j=i:n
        Batta_dist(i,j) = enorm_bhattacharyya(PredQB.pdfs{i}, PredQB.pdfs{j});
    end
end
figure(4)
imagesc(Batta_dist);