|
wolffd@0
|
1 function [mu, Sigma, weights, mask] = cwr_predict(cwr, X, mask_data)
|
|
wolffd@0
|
2 % CWR_PREDICT cluster weighted regression: predict Y given X
|
|
wolffd@0
|
3 % function [mu, Sigma] = cwr_predict(cwr, X)
|
|
wolffd@0
|
4 %
|
|
wolffd@0
|
5 % mu(:,t) = E[Y|X(:,t)] = sum_c P(c | X(:,t)) E[Y|c, X(:,t)]
|
|
wolffd@0
|
6 % Sigma(:,:,t) = Cov[Y|X(:,t)]
|
|
wolffd@0
|
7 %
|
|
wolffd@0
|
8 % [mu, Sigma, weights, mask] = cwr_predict(cwr, X, mask_data)
|
|
wolffd@0
|
9 % mask(i) = sum_t sum_c p(mask_data(:,i) | X(:,t), c) P(c|X(:,t))
|
|
wolffd@0
|
10 % This evaluates the predictive density on a set of points
|
|
wolffd@0
|
11 % (This is only sensible if T=1, ie. X is a single vector)
|
|
wolffd@0
|
12
|
|
wolffd@0
|
13 [nx T] = size(X);
|
|
wolffd@0
|
14 [ny nx nc] = size(cwr.weightsY);
|
|
wolffd@0
|
15 mu = zeros(ny, T);
|
|
wolffd@0
|
16 Sigma = zeros(ny, ny, T);
|
|
wolffd@0
|
17
|
|
wolffd@0
|
18 if nargout == 4
|
|
wolffd@0
|
19 comp_mask = 1;
|
|
wolffd@0
|
20 N = size(mask_data,2);
|
|
wolffd@0
|
21 mask = zeros(N,1);
|
|
wolffd@0
|
22 else
|
|
wolffd@0
|
23 comp_mask = 0;
|
|
wolffd@0
|
24 end
|
|
wolffd@0
|
25
|
|
wolffd@0
|
26 if nc==1
|
|
wolffd@0
|
27 if isempty(cwr.weightsY)
|
|
wolffd@0
|
28 mu = repmat(cwr.muY, 1, T);
|
|
wolffd@0
|
29 Sigma = repmat(cwr.SigmaY, [1 1 T]);
|
|
wolffd@0
|
30 else
|
|
wolffd@0
|
31 mu = repmat(cwr.muY, 1, T) + cwr.weightsY * X;
|
|
wolffd@0
|
32 Sigma = repmat(cwr.SigmaY, [1 1 T]);
|
|
wolffd@0
|
33 %for t=1:T
|
|
wolffd@0
|
34 % mu(:,t) = cwr.muY + cwr.weightsY*X(:,t);
|
|
wolffd@0
|
35 % Sigma(:,:,t) = cwr.SigmaY;
|
|
wolffd@0
|
36 %end
|
|
wolffd@0
|
37 end
|
|
wolffd@0
|
38 if comp_mask, mask = gaussian_prob(mask_data, mu, Sigma); end
|
|
wolffd@0
|
39 weights = [];
|
|
wolffd@0
|
40 return;
|
|
wolffd@0
|
41 end
|
|
wolffd@0
|
42
|
|
wolffd@0
|
43
|
|
wolffd@0
|
44 % likX(c,t) = p(x(:,t) | c)
|
|
wolffd@0
|
45 likX = mixgauss_prob(X, cwr.muX, cwr.SigmaX);
|
|
wolffd@0
|
46 weights = normalize(repmat(cwr.priorC, 1, T) .* likX, 1);
|
|
wolffd@0
|
47 for t=1:T
|
|
wolffd@0
|
48 mut = zeros(ny, nc);
|
|
wolffd@0
|
49 for c=1:nc
|
|
wolffd@0
|
50 mut(:,c) = cwr.muY(:,c) + cwr.weightsY(:,:,c)*X(:,t);
|
|
wolffd@0
|
51 if comp_mask
|
|
wolffd@0
|
52 mask = mask + gaussian_prob(mask_data, mut(:,c), cwr.SigmaY(:,:,c)) * weights(c);
|
|
wolffd@0
|
53 end
|
|
wolffd@0
|
54 end
|
|
wolffd@0
|
55 %w = normalise(cwr.priorC(:) .* likX(:,t));
|
|
wolffd@0
|
56 [mu(:,t), Sigma(:,:,t)] = collapse_mog(mut, cwr.SigmaY, weights(:,t));
|
|
wolffd@0
|
57 end
|
