pf
cf_filters.h
Go to the documentation of this file.
1 #ifndef CF_FILTERS_H
2 #define CF_FILTERS_H
3 
4 #include <Eigen/Dense> //linear algebra stuff
5 #include <math.h> /* log */
6 
7 #include "rv_eval.h"
8 
9 
11 
17 template<size_t dimstate, size_t dimobs, typename float_t>
18 class cf_filter{
19 
20 public:
21 
23  using ssv = Eigen::Matrix<float_t,dimstate,1>;
25  using osv = Eigen::Matrix<float_t,dimstate,1>;
26 
30  virtual ~cf_filter();
31 
32 
37  virtual float_t getLogCondLike() const = 0;
38 };
39 
40 
41 template<size_t dimstate, size_t dimobs, typename float_t>
43 
44 
46 
52 template<size_t dimstate, size_t dimobs, size_t diminput, typename float_t>
53 class kalman : public cf_filter<dimstate, dimobs, float_t> {
54 
55 public:
56 
58  using ssv = Eigen::Matrix<float_t,dimstate,1>;
59 
61  using osv = Eigen::Matrix<float_t,dimobs,1>;
62 
64  using isv = Eigen::Matrix<float_t,diminput,1>;
65 
67  using ssMat = Eigen::Matrix<float_t,dimstate,dimstate>;
68 
70  using osMat = Eigen::Matrix<float_t,dimobs,dimobs>;
71 
73  using siMat = Eigen::Matrix<float_t,dimstate,diminput>;
74 
76  using oiMat = Eigen::Matrix<float_t,dimobs,diminput>;
77 
79  using obsStateSizeMat = Eigen::Matrix<float_t,dimobs,dimstate>;
80 
82  using stateObsSizeMat = Eigen::Matrix<float_t,dimstate,dimobs>;
83 
84 
86 
89  kalman();
90 
91 
93 
96  kalman(const ssv &initStateMean, const ssMat &initStateVar);
97 
98 
102  virtual ~kalman();
103 
104 
109  float_t getLogCondLike() const;
110 
111 
116  ssv getFiltMean() const;
117 
118 
123  ssMat getFiltVar() const;
124 
125 
130  osv getPredYMean(const ssMat &stateTrans,
131  const obsStateSizeMat &obsMat,
132  const siMat &stateInptAffector,
133  const oiMat &obsInptAffector,
134  const isv &inputData) const;
135 
136 
141  osMat getPredYVar(const ssMat &stateTrans,
142  const ssMat &cholStateVar,
143  const obsStateSizeMat &obsMat,
144  const osMat &cholObsVar) const;
145 
146 
148 
158  void update(const osv &yt,
159  const ssMat &stateTrans,
160  const ssMat &cholStateVar,
161  const siMat &stateInptAffector,
162  const isv &inputData,
163  const obsStateSizeMat &obsMat,
164  const oiMat &obsInptAffector,
165  const osMat &cholObsVar);
166 
167 private:
168 
171 
174 
177 
180 
183 
185  bool m_fresh;
186 
188  const float_t m_pi;
189 
201  void updatePrior(const ssMat &stateTransMat,
202  const ssMat &cholStateVar,
203  const siMat &stateInptAffector,
204  const isv &inputData);
205 
206 
215  void updatePosterior(const osv &yt,
216  const obsStateSizeMat &obsMat,
217  const oiMat &obsInptAffector,
218  const isv &inputData,
219  const osMat &cholObsVar);
220 };
221 
222 
223 template<size_t dimstate, size_t dimobs, size_t diminput, typename float_t>
225  : cf_filter<dimstate,dimobs,float_t>()
226  , m_predMean(ssv::Zero())
227  , m_predVar(ssMat::Zero())
228  , m_fresh(true)
229  , m_pi(3.14159265358979)
230 {
231 }
232 
233 
234 template<size_t dimstate, size_t dimobs, size_t diminput, typename float_t>
235 kalman<dimstate,dimobs,diminput,float_t>::kalman(const ssv &initStateMean, const ssMat &initStateVar)
236  : cf_filter<dimstate,dimobs,float_t>()
237  , m_predMean(initStateMean)
238  , m_predVar(initStateVar)
239  , m_fresh(true)
240  , m_pi(3.14159265358979)
241 {
242 }
243 
244 
245 template<size_t dimstate, size_t dimobs, size_t diminput, typename float_t>
247 
248 
249 template<size_t dimstate, size_t dimobs, size_t diminput, typename float_t>
251  const ssMat &cholStateVar,
252  const siMat &stateInptAffector,
253  const isv &inputData)
254 {
255  ssMat Q = cholStateVar.transpose() * cholStateVar;
256  m_predMean = stateTransMat * m_filtMean + stateInptAffector * inputData;
257  m_predVar = stateTransMat * m_filtVar * stateTransMat.transpose() + Q;
258 }
259 
260 
261 template<size_t dimstate, size_t dimobs, size_t diminput, typename float_t>
263  const obsStateSizeMat &obsMat,
264  const oiMat &obsInptAffector,
265  const isv &inputData,
266  const osMat &cholObsVar)
267 {
268  osMat R = cholObsVar.transpose() * cholObsVar; //obs
269  osMat sigma = obsMat * m_predVar * obsMat.transpose() + R; // pred or APA' + R
270  osMat symSigma = (sigma.transpose() + sigma )/2.0; // ensure symmetric
271  osMat siginv = symSigma.inverse();
272  stateObsSizeMat K = m_predVar * obsMat.transpose() * siginv;
273  osv obsPred = obsMat * m_predMean + obsInptAffector * inputData;
274  osv innov = yt - obsPred;
275  m_filtMean = m_predMean + K*innov;
276  m_filtVar = m_predVar - K * obsMat * m_predVar;
277 
278  // conditional likelihood stuff
279  osMat quadForm = innov.transpose() * siginv * innov;
280  osMat cholSig ( sigma.llt().matrixL() );
281  float_t logDet = 2.0*cholSig.diagonal().array().log().sum();
282  m_lastLogCondLike = -.5*innov.rows()*log(2*m_pi) - .5*logDet - .5*quadForm(0,0);
283 }
284 
285 
286 template<size_t dimstate, size_t dimobs, size_t diminput, typename float_t>
288 {
289  return m_lastLogCondLike;
290 }
291 
292 
293 template<size_t dimstate, size_t dimobs, size_t diminput, typename float_t>
295 {
296  return m_filtMean;
297 }
298 
299 
300 template<size_t dimstate, size_t dimobs, size_t diminput, typename float_t>
302 {
303  return m_filtVar;
304 }
305 
306 
307 template<size_t dimstate, size_t dimobs, size_t diminput, typename float_t>
309  const ssMat &stateTrans,
310  const ssMat &cholStateVar,
311  const siMat &stateInptAffector,
312  const isv &inData,
313  const obsStateSizeMat &obsMat,
314  const oiMat &obsInptAffector,
315  const osMat &cholObsVar)
316 {
317  // this assumes that we have latent states x_{1:...} and y_{1:...} (NOT x_{0:...})
318  // for that reason, we don't have to run updatePrior() on the first iteration
319  if (m_fresh == true)
320  {
321  this->updatePosterior(yt, obsMat, obsInptAffector, inData, cholObsVar);
322  m_fresh = false;
323  }else
324  {
325  this->updatePrior(stateTrans, cholStateVar, stateInptAffector, inData);
326  this->updatePosterior(yt, obsMat, obsInptAffector, inData, cholObsVar);
327  }
328 }
329 
330 
331 template<size_t dimstate, size_t dimobs, size_t diminput, typename float_t>
333  const ssMat &stateTrans,
334  const obsStateSizeMat &obsMat,
335  const siMat &stateInptAffector,
336  const oiMat &obsInptAffector,
337  const isv &futureInputData) const -> osv
338 {
339  return obsMat * (stateTrans * m_filtMean + stateInptAffector * futureInputData) + obsInptAffector * futureInputData;
340 }
341 
342 
343 template<size_t dimstate, size_t dimobs, size_t diminput, typename float_t>
345  const ssMat &stateTrans,
346  const ssMat &cholStateVar,
347  const obsStateSizeMat &obsMat,
348  const osMat &cholObsVar) const -> osMat
349 {
350  return obsMat * (stateTrans * m_filtVar * stateTrans.transpose() + cholStateVar.transpose()*cholStateVar) * obsMat.transpose() + cholObsVar.transpose() * cholObsVar;
351 }
352 
353 
355 
361 template<size_t dimstate, size_t dimobs, typename float_t>
362 class hmm : public cf_filter<dimstate,dimobs,float_t>
363 {
364 
365 public:
366 
368  using ssv = Eigen::Matrix<float_t,dimstate,1>;
369 
371  using osv = Eigen::Matrix<float_t,dimobs,1>;
372 
374  using ssMat = Eigen::Matrix<float_t,dimstate,dimstate>;
375 
376 
378 
381  hmm();
382 
383 
385 
390  hmm(const ssv &initStateDistr, const ssMat &transMat);
391 
392 
396  virtual ~hmm();
397 
398 
400 
403  float_t getLogCondLike() const;
404 
405 
407 
411  ssv getFilterVec() const;
412 
413 
415 
419  void update(const ssv &condDensVec);
420 
421 
422 private:
423 
426 
429 
431  float_t m_lastCondLike;
432 
434  bool m_fresh;
435 
436 };
437 
438 
439 template<size_t dimstate, size_t dimobs, typename float_t>
441  : cf_filter<dimstate,dimobs,float_t>::cf_filter()
442  , m_filtVec(ssv::Zero())
443  , m_transMatTranspose(ssMat::Zero())
444  , m_lastCondLike(0.0)
445  , m_fresh(true)
446 {
447 }
448 
449 
450 template<size_t dimstate, size_t dimobs, typename float_t>
451 hmm<dimstate,dimobs,float_t>::hmm(const ssv &initStateDistr, const ssMat &transMat)
452  : cf_filter<dimstate,dimobs,float_t>()
453  , m_filtVec(initStateDistr)
454  , m_transMatTranspose(transMat.transpose())
455  , m_lastCondLike(0.0)
456  , m_fresh(true)
457 {
458 }
459 
460 
461 template<size_t dimstate, size_t dimobs, typename float_t>
463 
464 
465 template<size_t dimstate, size_t dimobs, typename float_t>
467 {
468  return std::log(m_lastCondLike);
469 }
470 
471 
472 template<size_t dimstate, size_t dimobs, typename float_t>
474 {
475  return m_filtVec;
476 }
477 
478 
479 template<size_t dimstate, size_t dimobs, typename float_t>
481 {
482  if (m_fresh) // hasn't seen data before and so filtVec is just time 1 state prior
483  {
484  m_filtVec = m_filtVec.cwiseProduct( condDensVec ); // now it's p(x_1, y_1)
485  m_lastCondLike = m_filtVec.sum();
487  m_fresh = false;
488 
489  } else { // has seen data before
490  m_filtVec = m_transMatTranspose * m_filtVec; // now p(x_t |y_{1:t-1})
491  m_filtVec = m_filtVec.cwiseProduct( condDensVec ); // now p(y_t,x_t|y_{1:t-1})
492  m_lastCondLike = m_filtVec.sum();
493  m_filtVec /= m_lastCondLike; // now p(x_t|y_{1:t})
494  }
495 }
496 
497 
498 
499 
500 
501 
502 
503 
505 
511 template<size_t dim_pred, typename float_t>
512 class gamFilter : public cf_filter<1,1,float_t>
513 {
514 
515 public:
516 
518  using psv = Eigen::Matrix<float_t,dim_pred,1>;
519 
521  using tsv = Eigen::Matrix<float_t,2,1>;
522 
523 
525 
528  //gamFilter();
529 
530 
532 
537  gamFilter(const float_t &nOneTilde, const float_t &dOneTilde);
538 
539 
543  virtual ~gamFilter();
544 
545 
547 
550  float_t getLogCondLike() const;
551 
552 
554 
558  tsv getFilterVec() const;
559 
560 
562 
570  void update(const float_t& yt, const psv &xt, const psv& beta, const float_t& sigmaSquared, const float_t& delta);
571 
572 
573 private:
574 
577 
580 
582  bool m_fresh;
583 
584 };
585 
586 
587 //template<typename dim_pred, typename float_t>
588 //gamFilter<dim_pred,float_t>::gamFilter()
589 // : cf_filter<1,1,float_t>::cf_filter()
590 // , m_filtVec(tsv::Zero())
591 // , m_lastCondLike(0.0)
592 // , m_fresh(true)
593 //{
594 //}
595 
596 
597 template<size_t dim_pred, typename float_t>
598 gamFilter<dim_pred,float_t>::gamFilter(const float_t &nOneTilde, const float_t &dOneTilde)
599  : cf_filter<1,1,float_t>()
600  , m_lastLogCondLike(0.0)
601  , m_fresh(true)
602 {
603  m_filtVec(0) = nOneTilde;
604  m_filtVec(1) = dOneTilde;
605 }
606 
607 
608 template<size_t dim_pred, typename float_t>
610 
611 
612 template<size_t dim_pred, typename float_t>
614 {
615  return m_lastLogCondLike;
616 }
617 
618 
619 template<size_t dim_pred, typename float_t>
621 {
622  return m_filtVec;
623 }
624 
625 
626 
627 template<size_t dim_pred, typename float_t>
628 void gamFilter<dim_pred,float_t>::update(const float_t& yt, const psv &xt, const psv& beta, const float_t& sigmaSquared, const float_t& delta)
629 {
630 
631  if(sigmaSquared <= 0 || delta <= 0)
632  throw std::invalid_argument("ME: both sigma squared and delta have to be positive.\n");
633 
634  if (m_fresh) // hasn't seen data before and so filtVec is just time 1 state prior
635  {
636  float_t tmpScale = std::sqrt(sigmaSquared*m_filtVec(1)/m_filtVec(0));
637  m_lastLogCondLike = rveval::evalScaledT<float_t>(yt, xt.dot(beta), tmpScale, m_filtVec(0), true);
638  m_filtVec(0) += 1;
639  m_filtVec(1) += (yt - xt.dot(beta))*(yt - xt.dot(beta))/sigmaSquared;
640  m_fresh = false;
641 
642  } else { // has seen data before
643 
644  m_filtVec(0) *= delta;
645  m_filtVec(1) *= delta;
646  float_t tmpScale = std::sqrt(sigmaSquared*m_filtVec(1)/m_filtVec(0));
647  m_lastLogCondLike = rveval::evalScaledT<float_t>(yt, xt.dot(beta), tmpScale, m_filtVec(0), true);
648  m_filtVec(0) += 1;
649  m_filtVec(1) += (yt - xt.dot(beta))*(yt - xt.dot(beta))/sigmaSquared;
650  }
651 }
652 
653 
655 // it's for a multivariate response.
662 template<size_t dim_obs, size_t dim_pred, typename float_t>
663 class multivGamFilter : public cf_filter<1,dim_obs,float_t>
664 {
665 
666 public:
667 
669  using psv = Eigen::Matrix<float_t,dim_pred,1>;
670 
672  using bsm = Eigen::Matrix<float_t,dim_obs, dim_pred>;
673 
675  using tsv = Eigen::Matrix<float_t,2,1>;
676 
678  using osv = Eigen::Matrix<float_t, dim_obs, 1>;
679 
681  using osm = Eigen::Matrix<float_t, dim_obs, dim_obs>;
682 
683 
685 
690  multivGamFilter(const float_t &nOneTilde, const float_t &dOneTilde);
691 
692 
696  virtual ~multivGamFilter();
697 
698 
700 
703  float_t getLogCondLike() const;
704 
705 
707 
711  tsv getFilterVec() const;
712 
713 
715 
723  void update(const osv& yt, const psv &xt, const bsm& B, const osm& Sigma, const float_t& delta);
724 
725 
727 
735  osv getFcastMean(const psv &xtp1, const bsm& B, const osm& Sigma, const float_t& delta);
736 
737 
739 
747  osm getFcastCov(const psv &xtp1, const bsm& B, const osm& Sigma, const float_t& delta);
748 
749 
750 private:
751 
754 
757 
759  bool m_fresh;
760 
761 };
762 
763 
764 template<size_t dim_obs, size_t dim_pred, typename float_t>
765 multivGamFilter<dim_obs,dim_pred,float_t>::multivGamFilter(const float_t &nOneTilde, const float_t &dOneTilde)
766  : cf_filter<1,dim_obs,float_t>()
767  , m_lastLogCondLike(0.0)
768  , m_fresh(true)
769 {
770  m_filtVec(0) = nOneTilde;
771  m_filtVec(1) = dOneTilde;
772 }
773 
774 
775 template<size_t dim_obs, size_t dim_pred, typename float_t>
777 
778 
779 template<size_t dim_obs, size_t dim_pred, typename float_t>
781 {
782  return m_lastLogCondLike;
783 }
784 
785 
786 template<size_t dim_obs, size_t dim_pred, typename float_t>
788 {
789  return m_filtVec;
790 }
791 
792 
793 
794 template<size_t dim_obs, size_t dim_pred, typename float_t>
795 void multivGamFilter<dim_obs,dim_pred,float_t>::update(const osv& yt, const psv &xt, const bsm& B, const osm& Sigma, const float_t& delta)
796 {
797 
798  // TODO: doesn't check that Sigma is positive definite or symmetric!
799  if(delta <= 0)
800  throw std::invalid_argument("ME: delta has to be positive (you're not even checking Sigma).\n");
801 
802  if (m_fresh) // hasn't seen data before and so filtVec is just time 1 state prior
803  {
804  osm scaleMat = Sigma * m_filtVec(1)/m_filtVec(0);
805  osv modeVec = B*xt;
806  m_lastLogCondLike = rveval::evalMultivT<dim_obs,float_t>(yt, modeVec, scaleMat, m_filtVec(0), true);
807  m_filtVec(0) += 1;
808  m_filtVec(1) += (Sigma.ldlt().solve(yt - modeVec)).squaredNorm();
809  m_fresh = false;
810 
811  } else { // has seen data before
812 
813  m_filtVec(0) *= delta;
814  m_filtVec(1) *= delta;
815 
816  osm scaleMat = Sigma * m_filtVec(1)/m_filtVec(0);
817  osv modeVec = B*xt;
818  m_lastLogCondLike = rveval::evalMultivT<dim_obs,float_t>(yt, modeVec, scaleMat, m_filtVec(0), true);
819 
820  m_filtVec(0) += 1;
821  m_filtVec(1) += (Sigma.ldlt().solve(yt - modeVec)).squaredNorm();
822  }
823 }
824 
825 
826 template<size_t dim_obs, size_t dim_pred, typename float_t>
827 auto multivGamFilter<dim_obs,dim_pred,float_t>::getFcastMean(const psv &xtp1, const bsm& B, const osm& Sigma, const float_t& delta) -> osv
828 {
829  if(delta*m_filtVec(0) > 1.0)
830  return B*xtp1;
831 }
832 
833 
834 template<size_t dim_obs, size_t dim_pred, typename float_t>
835 auto multivGamFilter<dim_obs,dim_pred,float_t>::getFcastCov(const psv &xtp1, const bsm& B, const osm& Sigma, const float_t& delta) -> osm
836 {
837  if(delta * m_filtVec(0) > 2.0)
838  return Sigma * delta * m_filtVec(1) / (delta * m_filtVec(0) - 2.0);
839 }
840 
841 
842 
843 
844 #endif //CF_FILTERS_H
Eigen::Matrix< float_t, 2, 1 > tsv
"two by 1 vector"
Definition: cf_filters.h:521
void update(const ssv &condDensVec)
Perform a HMM filter update.
Definition: cf_filters.h:480
tsv m_filtVec
filter vector (shape and rate)
Definition: cf_filters.h:576
Eigen::Matrix< float_t, dimstate, 1 > osv
Definition: cf_filters.h:25
Eigen::Matrix< float_t, dimstate, dimstate > ssMat
"state size matrix"
Definition: cf_filters.h:374
ssv m_filtMean
filter mean
Definition: cf_filters.h:173
ssv getFiltMean() const
Get the current filter mean.
Definition: cf_filters.h:294
Eigen::Matrix< float_t, dim_obs, 1 > osv
"observation size vector"
Definition: cf_filters.h:678
Eigen::Matrix< float_t, dim_obs, dim_obs > osm
"observation size matrix"
Definition: cf_filters.h:681
float_t m_lastLogCondLike
last log of the conditional likelihood
Definition: cf_filters.h:579
float_t m_lastCondLike
last conditional likelihood
Definition: cf_filters.h:431
float_t getLogCondLike() const
Get the latest conditional likelihood.
Definition: cf_filters.h:780
Eigen::Matrix< float_t, dim_obs, dim_pred > bsm
"beta size matrix"
Definition: cf_filters.h:672
Eigen::Matrix< float_t, diminput, 1 > isv
Definition: cf_filters.h:64
Another class template for Gamma filtering, but this time.
Definition: cf_filters.h:663
multivGamFilter(const float_t &nOneTilde, const float_t &dOneTilde)
Constructor.
Definition: cf_filters.h:765
Eigen::Matrix< float_t, dim_pred, 1 > psv
"predictor size vector"
Definition: cf_filters.h:669
ssMat m_filtVar
filter var matrix
Definition: cf_filters.h:179
virtual ~hmm()
The (virtual) desuctor.
Definition: cf_filters.h:462
Eigen::Matrix< float_t, dimobs, dimstate > obsStateSizeMat
Definition: cf_filters.h:79
A class template for Gamma filtering.
Definition: cf_filters.h:512
Eigen::Matrix< float_t, dimstate, 1 > ssv
Definition: cf_filters.h:23
void update(const osv &yt, const psv &xt, const bsm &B, const osm &Sigma, const float_t &delta)
Perform a filtering update.
Definition: cf_filters.h:795
ssMat getFiltVar() const
Get the current filter variance-covariance matrix.
Definition: cf_filters.h:301
Eigen::Matrix< float_t, 2, 1 > tsv
"two by 1 vector" to store size and shapes of gamma distributions
Definition: cf_filters.h:675
Eigen::Matrix< float_t, dimstate, dimstate > ssMat
Definition: cf_filters.h:67
void update(const osv &yt, const ssMat &stateTrans, const ssMat &cholStateVar, const siMat &stateInptAffector, const isv &inputData, const obsStateSizeMat &obsMat, const oiMat &obsInptAffector, const osMat &cholObsVar)
Perform a Kalman filter predict-and-update.
Definition: cf_filters.h:308
A class template for Kalman filtering.
Definition: cf_filters.h:53
Eigen::Matrix< float_t, dimobs, diminput > oiMat
Definition: cf_filters.h:76
osm getFcastCov(const psv &xtp1, const bsm &B, const osm &Sigma, const float_t &delta)
Get the forecast covariance matrix (assuming filtering has been performed already) ...
Definition: cf_filters.h:835
ssv m_filtVec
filter vector
Definition: cf_filters.h:425
virtual ~kalman()
The (virtual) destructor.
Definition: cf_filters.h:246
ssv m_predMean
predictive state mean
Definition: cf_filters.h:170
bool m_fresh
has data been observed?
Definition: cf_filters.h:759
ssv getFilterVec() const
Get the current filter vector.
Definition: cf_filters.h:473
osv getFcastMean(const psv &xtp1, const bsm &B, const osm &Sigma, const float_t &delta)
Get the forecast mean (assuming filtering has been performed already)
Definition: cf_filters.h:827
tsv m_filtVec
filter vector (shape and rate)
Definition: cf_filters.h:753
osMat getPredYVar(const ssMat &stateTrans, const ssMat &cholStateVar, const obsStateSizeMat &obsMat, const osMat &cholObsVar) const
get the one-step-ahead forecast variance
Definition: cf_filters.h:344
Eigen::Matrix< float_t, dim_pred, 1 > psv
"predictor size vector"
Definition: cf_filters.h:518
bool m_fresh
has data been observed?
Definition: cf_filters.h:434
bool m_fresh
has data been observed?
Definition: cf_filters.h:185
float_t m_lastLogCondLike
latest log conditional likelihood
Definition: cf_filters.h:182
virtual float_t getLogCondLike() const =0
returns the log of the most recent conditional likelihood
kalman()
Default constructor.
Definition: cf_filters.h:224
A class template for HMM filtering.
Definition: cf_filters.h:362
Eigen::Matrix< float_t, dimstate, diminput > siMat
Definition: cf_filters.h:73
float_t getLogCondLike() const
returns the log of the latest conditional likelihood.
Definition: cf_filters.h:287
float_t getLogCondLike() const
Get the latest conditional likelihood.
Definition: cf_filters.h:466
gamFilter(const float_t &nOneTilde, const float_t &dOneTilde)
Default constructor.
Definition: cf_filters.h:598
void update(const float_t &yt, const psv &xt, const psv &beta, const float_t &sigmaSquared, const float_t &delta)
Perform a filtering update.
Definition: cf_filters.h:628
ssMat m_transMatTranspose
transition matrix
Definition: cf_filters.h:428
const float_t m_pi
pi
Definition: cf_filters.h:188
void updatePosterior(const osv &yt, const obsStateSizeMat &obsMat, const oiMat &obsInptAffector, const isv &inputData, const osMat &cholObsVar)
Turns prediction into new filtering distribution.
Definition: cf_filters.h:262
osv getPredYMean(const ssMat &stateTrans, const obsStateSizeMat &obsMat, const siMat &stateInptAffector, const oiMat &obsInptAffector, const isv &inputData) const
get the one-step-ahead point forecast for y
Definition: cf_filters.h:332
virtual ~multivGamFilter()
The (virtual) desuctor.
Definition: cf_filters.h:776
float_t m_lastLogCondLike
last log of the conditional likelihood
Definition: cf_filters.h:756
Abstract Base Class for all closed-form filters.
Definition: cf_filters.h:18
tsv getFilterVec() const
Get the current filter vector.
Definition: cf_filters.h:787
hmm()
Default constructor.
Definition: cf_filters.h:440
virtual ~cf_filter()
The (virtual) destructor.
Definition: cf_filters.h:42
Eigen::Matrix< float_t, dimobs, dimobs > osMat
Definition: cf_filters.h:70
ssMat m_predVar
predictive var matrix
Definition: cf_filters.h:176
virtual ~gamFilter()
The (virtual) desuctor.
Definition: cf_filters.h:609
tsv getFilterVec() const
Get the current filter vector.
Definition: cf_filters.h:620
bool m_fresh
has data been observed?
Definition: cf_filters.h:582
Eigen::Matrix< float_t, dimstate, dimobs > stateObsSizeMat
Definition: cf_filters.h:82
void updatePrior(const ssMat &stateTransMat, const ssMat &cholStateVar, const siMat &stateInptAffector, const isv &inputData)
Predicts the next state.
Definition: cf_filters.h:250
float_t getLogCondLike() const
Get the latest conditional likelihood.
Definition: cf_filters.h:613