KalmanFilter trueKf1 = new KalmanFilter(model1, modelCovariance1, measurementCovariance);
trueKf1.setCurrentInput(VectorFactory.getDefault().copyValues(truePsis.get(0).getElement(0)));
KalmanFilter trueKf2 = new KalmanFilter(model2, modelCovariance2, measurementCovariance);
trueKf2.setCurrentInput(VectorFactory.getDefault().copyValues(truePsis.get(1).getElement(0)));
Vector initialClassProbs = VectorFactory.getDefault()
.copyArray(new double[] { 0.5d, 0.5d });
Matrix classTransProbs = MatrixFactory.getDefault().copyArray(
new double[][] { { 0.5d, 0.5d },
{ 0.5d, 0.5d } });
DlmHiddenMarkovModel trueHmm1 = new DlmHiddenMarkovModel(
Lists.newArrayList(trueKf1, trueKf2),
initialClassProbs, classTransProbs);
final double sigmaPriorMean = Math.pow(0.4, 2);
final double sigmaPriorShape = 2d;
final double sigmaPriorScale = sigmaPriorMean*(sigmaPriorShape + 1d);
final InverseGammaDistribution sigmaPrior = new InverseGammaDistribution(sigmaPriorShape,
sigmaPriorScale);
final Vector phiMean1 = VectorFactory.getDefault().copyArray(new double[] {
0d, 0.8d
});
final Matrix phiCov1 = MatrixFactory.getDefault().copyArray(new double[][] {
{2d + 4d * sigmaPriorMean, 0d},
{ 0d, 4d * sigmaPriorMean}
});
final MultivariateGaussian phiPrior1 = new MultivariateGaussian(phiMean1, phiCov1);
final Vector phiMean2 = VectorFactory.getDefault().copyArray(new double[] {
0d, 0.1d
});
final Matrix phiCov2 = MatrixFactory.getDefault().copyArray(new double[][] {
{ 1d + 4d * sigmaPriorMean, 0d},
{ 0d, 4d * sigmaPriorMean}
});
final MultivariateGaussian phiPrior2 = new MultivariateGaussian(phiMean2, phiCov2);
List<MultivariateGaussian> priorPhis = Lists.newArrayList(phiPrior1, phiPrior2);
final HmmPlFilter<DlmHiddenMarkovModel, GaussianArHpTransitionState, Vector> wfFilter =
new GaussianArHpHmmPLFilter(trueHmm1, sigmaPrior, priorPhis, random, true);
final int K = 3;
final int T = 200;
final int N = 1000;
/*
* Note: replications are over the same set of simulated observations.
*/
List<SimHmmObservedValue<Vector, Vector>> simulation = trueHmm1.sample(random, T);
wfFilter.setNumParticles(N);
wfFilter.setResampleOnly(false);
log.info("rep\tt\tfilter.type\tmeasurement.type\tresample.type\tmeasurement");
GaussianArHmmClassEvaluator wfClassEvaluator = new GaussianArHmmClassEvaluator("wf-pl",
null);
GaussianArHmmRmseEvaluator wfRmseEvaluator = new GaussianArHmmRmseEvaluator("wf-pl",
null);
GaussianArHmmPsiLearningEvaluator wfPsiEvaluator = new GaussianArHmmPsiLearningEvaluator("wf-pl",
truePsis, null);
RingAccumulator<MutableDouble> wfLatency =
new RingAccumulator<MutableDouble>();
Stopwatch wfWatch = new Stopwatch();
for (int k = 0; k < K; k++) {
log.info("Processing replication " + k);
CountedDataDistribution<GaussianArHpTransitionState> wfDistribution =
(CountedDataDistribution<GaussianArHpTransitionState>) wfFilter.getUpdater().createInitialParticles(N);
final long numPreRuns = -1l;//wfDistribution.getMaxValueKey().getTime();
/*
* Recurse through the particle filter
*/
for (int i = 0; i < T; i++) {
final double x = simulation.get(i).getClassId();
final Vector y = simulation.get(i).getObservedValue();
if (i > numPreRuns) {
if (i > 0) {
wfWatch.reset();