|
wolffd@0
|
1 % Example of fixed lag smoothing
|
|
wolffd@0
|
2
|
|
wolffd@0
|
3 rand('state', 1);
|
|
wolffd@0
|
4 S = 2;
|
|
wolffd@0
|
5 O = 2;
|
|
wolffd@0
|
6 T = 7;
|
|
wolffd@0
|
7 data = sample_discrete([0.5 0.5], 1, T);
|
|
wolffd@0
|
8 transmat = mk_stochastic(rand(S,S));
|
|
wolffd@0
|
9 obsmat = mk_stochastic(rand(S,O));
|
|
wolffd@0
|
10 obslik = multinomial_prob(data, obsmat);
|
|
wolffd@0
|
11 prior = [0.5 0.5]';
|
|
wolffd@0
|
12
|
|
wolffd@0
|
13
|
|
wolffd@0
|
14 [alpha0, beta0, gamma0, ll0, xi0] = fwdback(prior, transmat, obslik);
|
|
wolffd@0
|
15
|
|
wolffd@0
|
16 w = 3;
|
|
wolffd@0
|
17 alpha1 = zeros(S, T);
|
|
wolffd@0
|
18 gamma1 = zeros(S, T);
|
|
wolffd@0
|
19 xi1 = zeros(S, S, T-1);
|
|
wolffd@0
|
20 t = 1;
|
|
wolffd@0
|
21 b = obsmat(:, data(t));
|
|
wolffd@0
|
22 olik_win = b; % window of conditional observation likelihoods
|
|
wolffd@0
|
23 alpha_win = normalise(prior .* b);
|
|
wolffd@0
|
24 alpha1(:,t) = alpha_win;
|
|
wolffd@0
|
25 for t=2:T
|
|
wolffd@0
|
26 [alpha_win, olik_win, gamma_win, xi_win] = ...
|
|
wolffd@0
|
27 fixed_lag_smoother(w, alpha_win, olik_win, obsmat(:, data(t)), transmat);
|
|
wolffd@0
|
28 alpha1(:,max(1,t-w+1):t) = alpha_win;
|
|
wolffd@0
|
29 gamma1(:,max(1,t-w+1):t) = gamma_win;
|
|
wolffd@0
|
30 xi1(:,:,max(1,t-w+1):t-1) = xi_win;
|
|
wolffd@0
|
31 end
|
|
wolffd@0
|
32
|
|
wolffd@0
|
33 e = 1e-1;
|
|
wolffd@0
|
34 %assert(approxeq(alpha0, alpha1, e));
|
|
wolffd@0
|
35 assert(approxeq(gamma0(:, T-w+1:end), gamma1(:, T-w+1:end), e));
|
|
wolffd@0
|
36 %assert(approxeq(xi0(:,:,T-w+1:end), xi1(:,:,T-w+1:end), e));
|
|
wolffd@0
|
37
|
|
wolffd@0
|
38
|
