26 template<
size_t nparts,
size_t dimnss,
size_t dimss,
size_t dimy,
typename resamp_t,
typename float_t >
32 using sssv = Eigen::Matrix<float_t,dimss,1>;
34 using nsssv = Eigen::Matrix<float_t,dimnss,1>;
36 using osv = Eigen::Matrix<float_t,dimy,1>;
38 using nsssMat = Eigen::Matrix<float_t,dimnss,dimnss>;
40 using Mat = Eigen::Matrix<float_t,Eigen::Dynamic,Eigen::Dynamic>;
42 using arrayMod = std::array<hmm<dimnss,dimy,float_t>,nparts>;
54 rbpf_hmm(
const unsigned int &resamp_sched=1);
70 const std::vector<std::function<
const Mat(
const nsssv &x1tProbs,
const sssv &x2t)> >& fs
71 = std::vector<std::function<
const Mat(
const nsssv&,
const sssv&)> >());
195 template<
size_t nparts,
size_t dimnss,
size_t dimss,
size_t dimy,
typename resamp_t,
typename float_t>
205 template<
size_t nparts,
size_t dimnss,
size_t dimss,
size_t dimy,
typename resamp_t,
typename float_t>
209 template<
size_t nparts,
size_t dimnss,
size_t dimss,
size_t dimy,
typename resamp_t,
typename float_t>
218 float_t sumexpdenom(0.0);
219 float_t m1(-std::numeric_limits<float_t>::infinity());
221 for(
size_t ii = 0; ii < nparts; ++ii){
239 float_t sumexpnumer(0.0);
240 for(
size_t p = 0; p < nparts; ++p)
250 unsigned int rows = testOutput.rows();
251 unsigned int cols = testOutput.cols();
252 Mat numer = Mat::Zero(rows,cols);
254 for(
size_t prtcl = 0; prtcl < nparts; ++prtcl){
275 float_t m1(-std::numeric_limits<float_t>::infinity());
276 for(
size_t ii = 0; ii < nparts; ++ii){
292 for(
size_t p = 0; p < nparts; ++p)
294 m_lastLogCondLike = m1 + std::log(sumexp) - std::log(static_cast<float_t>(nparts));
303 unsigned int rows = testOutput.rows();
304 unsigned int cols = testOutput.cols();
305 Mat numer = Mat::Zero(rows,cols);
307 for(
size_t prtcl = 0; prtcl < nparts; ++prtcl){
326 template<
size_t nparts,
size_t dimnss,
size_t dimss,
size_t dimy,
typename resamp_t,
typename float_t>
333 template<
size_t nparts,
size_t dimnss,
size_t dimss,
size_t dimy,
typename resamp_t,
typename float_t>
352 template<
size_t nparts,
size_t dimnss,
size_t dimss,
size_t dimy,
typename resamp_t,
typename float_t>
358 using sssv = Eigen::Matrix<float_t,dimss,1>;
360 using nsssv = Eigen::Matrix<float_t,dimnss,1>;
362 using osv = Eigen::Matrix<float_t,dimy,1>;
364 using nsssMat = Eigen::Matrix<float_t,dimnss,dimnss>;
366 using Mat = Eigen::Matrix<float_t,Eigen::Dynamic,Eigen::Dynamic>;
368 using arrayMod = std::array<hmm<dimnss,dimy,float_t>,nparts>;
396 const std::vector<std::function<
const Mat(
const nsssv &x1tProbs,
const sssv &x2t)> >& fs
397 = std::vector<std::function<
const Mat(
const nsssv&,
const sssv&)> >());
420 virtual sssv muSamp() = 0;
446 virtual sssv fSamp(
const sssv &x2tm1) = 0;
480 template<
size_t nparts,
size_t dimnss,
size_t dimss,
size_t dimy,
typename resamp_t,
typename float_t>
490 template<
size_t nparts,
size_t dimnss,
size_t dimss,
size_t dimy,
typename resamp_t,
typename float_t>
494 template<
size_t nparts,
size_t dimnss,
size_t dimss,
size_t dimy,
typename resamp_t,
typename float_t>
502 float_t sumexpdenom(0.0);
503 float_t m1(-std::numeric_limits<float_t>::infinity());
505 for(
size_t ii = 0; ii < nparts; ++ii){
521 float_t sumexpnumer(0.0);
522 for(
size_t p = 0; p < nparts; ++p)
532 unsigned int rows = testOutput.rows();
533 unsigned int cols = testOutput.cols();
534 Mat numer = Mat::Zero(rows,cols);
536 for(
size_t prtcl = 0; prtcl < nparts; ++prtcl){
556 float_t m1(-std::numeric_limits<float_t>::infinity());
557 for(
size_t ii = 0; ii < nparts; ++ii){
573 for(
size_t p = 0; p < nparts; ++p)
575 m_lastLogCondLike = m1 + std::log(sumexp) - std::log(static_cast<float_t>(nparts));
584 unsigned int rows = testOutput.rows();
585 unsigned int cols = testOutput.cols();
586 Mat numer = Mat::Zero(rows,cols);
588 for(
size_t prtcl = 0; prtcl < nparts; ++prtcl){
607 template<
size_t nparts,
size_t dimnss,
size_t dimss,
size_t dimy,
typename resamp_t,
typename float_t>
614 template<
size_t nparts,
size_t dimnss,
size_t dimss,
size_t dimy,
typename resamp_t,
typename float_t>
633 template<
size_t nparts,
size_t dimnss,
size_t dimss,
size_t dimy,
typename resamp_t,
typename float_t>
640 using sssv = Eigen::Matrix<float_t,dimss,1>;
642 using nsssv = Eigen::Matrix<float_t,dimnss,1>;
644 using osv = Eigen::Matrix<float_t,dimy,1>;
646 using Mat = Eigen::Matrix<float_t,Eigen::Dynamic,Eigen::Dynamic>;
648 using nsssMat = Eigen::Matrix<float_t,dimnss,dimnss>;
650 using arrayMod = std::array<kalman<dimnss,dimy,0,float_t>,nparts>;
675 void filter(
const osv &data,
const std::vector<std::function<
const Mat(
const nsssv &x1t,
const sssv &x2t)> >& fs
676 = std::vector<std::function<
const Mat(
const nsssv &x1t,
const sssv &x2t)> >() );
699 virtual float_t logMuEv(
const sssv &x21) = 0;
708 virtual sssv q1Samp(
const osv &y1) = 0;
717 virtual nsssv initKalmanMean(
const sssv &x21) = 0;
726 virtual nsssMat initKalmanVar(
const sssv &x21) = 0;
736 virtual sssv qSamp(
const sssv &x2tm1,
const osv &yt) = 0;
746 virtual float_t logQ1Ev(
const sssv &x21,
const osv &y1) = 0;
756 virtual float_t logFEv(
const sssv &x2t,
const sssv &x2tm1) = 0;
767 virtual float_t logQEv(
const sssv &x2t,
const sssv &x2tm1,
const osv &yt) = 0;
800 template<
size_t nparts,
size_t dimnss,
size_t dimss,
size_t dimy,
typename resamp_t,
typename float_t>
810 template<
size_t nparts,
size_t dimnss,
size_t dimss,
size_t dimy,
typename resamp_t,
typename float_t>
814 template<
size_t nparts,
size_t dimnss,
size_t dimss,
size_t dimy,
typename resamp_t,
typename float_t>
823 float_t m1(-std::numeric_limits<float_t>::infinity());
825 float_t sumexpdenom(0.0);
826 for(
size_t ii = 0; ii < nparts; ++ii){
844 float_t sumexpnumer(0.0);
845 for(
size_t p = 0; p < nparts; ++p)
855 unsigned int rows = testOutput.rows();
856 unsigned int cols = testOutput.cols();
857 Mat numer = Mat::Zero(rows,cols);
859 for(
size_t prtcl = 0; prtcl < nparts; ++prtcl){
879 float_t m1(-std::numeric_limits<float_t>::infinity());
880 for(
size_t ii = 0; ii < nparts; ++ii){
896 for(
size_t p = 0; p < nparts; ++p)
898 m_lastLogCondLike = m1 + std::log(sumexp) - std::log(static_cast<float_t>(nparts));
907 unsigned int rows = testOutput.rows();
908 unsigned int cols = testOutput.cols();
909 Mat numer = Mat::Zero(rows,cols);
911 for(
size_t prtcl = 0; prtcl < nparts; ++prtcl){
930 template<
size_t nparts,
size_t dimnss,
size_t dimss,
size_t dimy,
typename resamp_t,
typename float_t>
937 template<
size_t nparts,
size_t dimnss,
size_t dimss,
size_t dimy,
typename resamp_t,
typename float_t>
957 template<
size_t nparts,
size_t dimnss,
size_t dimss,
size_t dimy,
typename resamp_t,
typename float_t>
964 using sssv = Eigen::Matrix<float_t,dimss,1>;
966 using nsssv = Eigen::Matrix<float_t,dimnss,1>;
968 using osv = Eigen::Matrix<float_t,dimy,1>;
970 using Mat = Eigen::Matrix<float_t,Eigen::Dynamic,Eigen::Dynamic>;
972 using nsssMat = Eigen::Matrix<float_t,dimnss,dimnss>;
974 using arrayMod = std::array<kalman<dimnss,dimy,0,float_t>,nparts>;
999 void filter(
const osv &data,
const std::vector<std::function<
const Mat(
const nsssv &x1t,
const sssv &x2t)> >& fs
1000 = std::vector<std::function<
const Mat(
const nsssv &x1t,
const sssv &x2t)> >() );
1022 virtual sssv muSamp() = 0;
1049 virtual sssv fSamp(
const sssv &x2tm1) = 0;
1082 template<
size_t nparts,
size_t dimnss,
size_t dimss,
size_t dimy,
typename resamp_t,
typename float_t>
1086 ,
m_rs(resamp_sched)
1092 template<
size_t nparts,
size_t dimnss,
size_t dimss,
size_t dimy,
typename resamp_t,
typename float_t>
1096 template<
size_t nparts,
size_t dimnss,
size_t dimss,
size_t dimy,
typename resamp_t,
typename float_t>
1105 float_t m1(-std::numeric_limits<float_t>::infinity());
1107 float_t sumexpdenom(0.0);
1108 for(
size_t ii = 0; ii < nparts; ++ii){
1127 float_t sumexpnumer(0.0);
1128 for(
size_t p = 0; p < nparts; ++p)
1133 unsigned int fId(0);
1138 unsigned int rows = testOutput.rows();
1139 unsigned int cols = testOutput.cols();
1140 Mat numer = Mat::Zero(rows,cols);
1142 for(
size_t prtcl = 0; prtcl < nparts; ++prtcl){
1162 float_t m1(-std::numeric_limits<float_t>::infinity());
1163 for(
size_t ii = 0; ii < nparts; ++ii){
1178 float_t sumexp(0.0);
1179 for(
size_t p = 0; p < nparts; ++p)
1181 m_lastLogCondLike = m1 + std::log(sumexp) - std::log(static_cast<float_t>(nparts));
1185 unsigned int fId(0);
1190 unsigned int rows = testOutput.rows();
1191 unsigned int cols = testOutput.cols();
1192 Mat numer = Mat::Zero(rows,cols);
1194 for(
size_t prtcl = 0; prtcl < nparts; ++prtcl){
1213 template<
size_t nparts,
size_t dimnss,
size_t dimss,
size_t dimy,
typename resamp_t,
typename float_t>
1220 template<
size_t nparts,
size_t dimnss,
size_t dimss,
size_t dimy,
typename resamp_t,
typename float_t>
resamp_t m_resampler
Definition: rbpf.h:794
virtual void updateHMM(hmm< dimnss, dimy, float_t > &aModel, const osv &yt, const sssv &x2t)=0
How to update your inner HMM filter object at each time.
std::array< kalman< dimnss, dimy, 0, float_t >, nparts > arrayMod
Definition: rbpf.h:974
virtual nsssv initHMMProbVec(const sssv &x21)=0
Provides the initial mean vector for each HMM filter object.
virtual void updateHMM(hmm< dimnss, dimy, float_t > &aModel, const osv &yt, const sssv &x2t)=0
How to update your inner HMM filter object at each time.
Eigen::Matrix< float_t, Eigen::Dynamic, Eigen::Dynamic > Mat
Definition: rbpf.h:970
arrayfloat_t m_logUnNormWeights
Definition: rbpf.h:788
arrayMod m_p_innerMods
Definition: rbpf.h:1066
Eigen::Matrix< float_t, dimnss, 1 > nsssv
Definition: rbpf.h:34
virtual sssv fSamp(const sssv &x2tm1)=0
Samples the time t second component.
resamp_t m_resampler
Definition: rbpf.h:1076
Eigen::Matrix< float_t, dimnss, 1 > nsssv
Definition: rbpf.h:642
std::array< hmm< dimnss, dimy, float_t >, nparts > arrayMod
Definition: rbpf.h:42
Eigen::Matrix< float_t, Eigen::Dynamic, Eigen::Dynamic > Mat
Definition: rbpf.h:366
arrayMod m_p_innerMods
Definition: rbpf.h:784
unsigned int m_now
Definition: rbpf.h:461
virtual nsssMat initHMMTransMat(const sssv &x21)=0
Provides the transition matrix for each HMM filter object.
float_t getLogCondLike() const
Get the latest log conditional likelihood.
Definition: rbpf.h:1214
Eigen::Matrix< float_t, dimnss, 1 > nsssv
Definition: rbpf.h:966
std::array< float_t, nparts > arrayfloat_t
Definition: rbpf.h:46
std::vector< Mat > m_expectations
Definition: rbpf.h:475
virtual sssv qSamp(const sssv &x2tm1, const osv &yt)=0
Samples the time t second component.
Eigen::Matrix< float_t, Eigen::Dynamic, Eigen::Dynamic > Mat
Definition: rbpf.h:646
Eigen::Matrix< float_t, dimnss, dimnss > nsssMat
Definition: rbpf.h:38
virtual ~rbpf_hmm()
The (virtual) destructor.
Definition: rbpf.h:206
std::array< float_t, nparts > arrayfloat_t
Definition: rbpf.h:654
std::array< sssv, nparts > arrayVec
Definition: rbpf.h:44
virtual nsssMat initHMMTransMat(const sssv &x21)=0
Provides the transition matrix for each HMM filter object.
rbpf_kalman(const unsigned int &resamp_sched=1)
The constructor.
Definition: rbpf.h:801
resamp_t m_resampler
Definition: rbpf.h:473
rbpf_kalman_bs(const unsigned int &resamp_sched=1)
The constructor.
Definition: rbpf.h:1083
std::vector< Mat > m_expectations
Definition: rbpf.h:796
Rao-Blackwellized/Marginal Bootstrap Filter with inner Kalman Filter objectss.
Definition: rbpf.h:958
float_t m_lastLogCondLike
Definition: rbpf.h:792
Eigen::Matrix< float_t, dimy, 1 > osv
Definition: rbpf.h:968
float_t m_lastLogCondLike
Definition: rbpf.h:178
Eigen::Matrix< float_t, dimy, 1 > osv
Definition: rbpf.h:644
std::vector< Mat > m_expectations
Definition: rbpf.h:1078
float_t getLogCondLike() const
Get the latest conditional likelihood.
Definition: rbpf.h:608
resamp_t m_resampler
Definition: rbpf.h:188
float_t m_lastLogCondLike
Definition: rbpf.h:463
A class template for Kalman filtering.
Definition: cf_filters.h:53
virtual void updateKalman(kalman< dimnss, dimy, 0, float_t > &kMod, const osv &yt, const sssv &x2t)=0
How to update your inner Kalman filter object at each time.
float_t m_lastLogCondLike
Definition: rbpf.h:1074
std::array< sssv, nparts > arrayVec
Definition: rbpf.h:370
Eigen::Matrix< float_t, dimss, 1 > sssv
Definition: rbpf.h:358
Eigen::Matrix< float_t, dimnss, dimnss > nsssMat
Definition: rbpf.h:972
float_t getLogCondLike() const
Get the latest log conditional likelihood.
Definition: rbpf.h:931
virtual float_t logFEv(const sssv &x2t, const sssv &x2tm1)=0
Evaluates the state transition density for the second state component.
unsigned int m_now
Definition: rbpf.h:790
arrayfloat_t m_logUnNormWeights
Definition: rbpf.h:471
virtual ~rbpf_hmm_bs()
The (virtual) destructor.
Definition: rbpf.h:491
Eigen::Matrix< float_t, dimy, 1 > osv
Definition: rbpf.h:362
Eigen::Matrix< float_t, dimnss, 1 > nsssv
Definition: rbpf.h:360
virtual nsssMat initKalmanVar(const sssv &x21)=0
Provides the initial covariance matrix for each Kalman filter object.
void filter(const osv &data, const std::vector< std::function< const Mat(const nsssv &x1tProbs, const sssv &x2t)> > &fs=std::vector< std::function< const Mat(const nsssv &, const sssv &)> >())
Filter.
Definition: rbpf.h:495
arrayVec m_p_samps
Definition: rbpf.h:1068
virtual ~rbpf_kalman_bs()
The (virtual) destructor.
Definition: rbpf.h:1093
unsigned int m_now
Definition: rbpf.h:1072
std::vector< Mat > m_expectations
Definition: rbpf.h:190
Rao-Blackwellized/Marginal Particle Filter with inner Kalman Filter objectss.
Definition: rbpf.h:634
std::vector< Mat > getExpectations() const
Get the latest filtered expectation E[h(x_1t, x_2t) | y_{1:t}].
Definition: rbpf.h:1221
virtual float_t logFEv(const sssv &x2t, const sssv &x2tm1)=0
Evaluates the state transition density for the second state component.
unsigned int m_rs
Definition: rbpf.h:1064
arrayVec m_p_samps
Definition: rbpf.h:469
float_t getLogCondLike() const
Get the latest conditional likelihood.
Definition: rbpf.h:327
void filter(const osv &data, const std::vector< std::function< const Mat(const nsssv &x1tProbs, const sssv &x2t)> > &fs=std::vector< std::function< const Mat(const nsssv &, const sssv &)> >())
Filter.
Definition: rbpf.h:210
virtual sssv muSamp()=0
Sample from the first sampler.
std::array< sssv, nparts > arrayVec
Definition: rbpf.h:652
std::array< float_t, nparts > arrayfloat_t
Definition: rbpf.h:372
virtual nsssv initHMMProbVec(const sssv &x21)=0
Provides the initial mean vector for each HMM filter object.
virtual sssv q1Samp(const osv &y1)=0
Sample from the first sampler.
A class template for HMM filtering.
Definition: cf_filters.h:362
Eigen::Matrix< float_t, dimss, 1 > sssv
Definition: rbpf.h:32
arrayMod m_p_innerMods
Definition: rbpf.h:182
All non Rao-Blackwellized particle filters inherit from this.
std::array< kalman< dimnss, dimy, 0, float_t >, nparts > arrayMod
Definition: rbpf.h:650
virtual float_t logMuEv(const sssv &x21)=0
Evaluates the first time state density.
std::vector< Mat > getExpectations() const
Get vector of expectations.
Definition: rbpf.h:334
virtual nsssMat initKalmanVar(const sssv &x21)=0
Provides the initial covariance matrix for each Kalman filter object.
unsigned int m_rs
Definition: rbpf.h:180
void filter(const osv &data, const std::vector< std::function< const Mat(const nsssv &x1t, const sssv &x2t)> > &fs=std::vector< std::function< const Mat(const nsssv &x1t, const sssv &x2t)> >())
Filter!
Definition: rbpf.h:815
rbpf_hmm(const unsigned int &resamp_sched=1)
The constructor.
Definition: rbpf.h:196
arrayfloat_t m_logUnNormWeights
Definition: rbpf.h:186
virtual float_t logQ1Ev(const sssv &x21, const osv &y1)=0
Evaluates the proposal density of the second state component at time 1.
Eigen::Matrix< float_t, dimnss, dimnss > nsssMat
Definition: rbpf.h:648
virtual sssv qSamp(const sssv &x2tm1, const osv &yt)=0
Samples the time t second component.
virtual sssv q1Samp(const osv &y1)=0
Sample from the first time's proposal distribution.
arrayVec m_p_samps
Definition: rbpf.h:786
Eigen::Matrix< float_t, dimss, 1 > sssv
Definition: rbpf.h:964
rbpf_hmm_bs(const unsigned int &resamp_sched=1)
The constructor.
Definition: rbpf.h:481
unsigned int m_rs
Definition: rbpf.h:465
virtual float_t logQEv(const sssv &x2t, const sssv &x2tm1, const osv &yt)=0
Evaluates the proposal density at time t > 1.
Eigen::Matrix< float_t, Eigen::Dynamic, Eigen::Dynamic > Mat
Definition: rbpf.h:40
virtual float_t logMuEv(const sssv &x21)=0
Evaluates the first time state density.
virtual void updateKalman(kalman< dimnss, dimy, 0, float_t > &kMod, const osv &yt, const sssv &x2t)=0
How to update your inner Kalman filter object at each time.
Rao-Blackwellized/Marginal Particle Filter with inner HMMs.
Definition: rbpf.h:27
std::array< float_t, nparts > arrayfloat_t
Definition: rbpf.h:978
std::vector< Mat > getExpectations() const
Get vector of expectations.
Definition: rbpf.h:615
std::array< hmm< dimnss, dimy, float_t >, nparts > arrayMod
Definition: rbpf.h:368
virtual nsssv initKalmanMean(const sssv &x21)=0
Provides the initial mean vector for each Kalman filter object.
Eigen::Matrix< float_t, dimss, 1 > sssv
Definition: rbpf.h:640
virtual nsssv initKalmanMean(const sssv &x21)=0
Provides the initial mean vector for each Kalman filter object.
unsigned int m_now
Definition: rbpf.h:176
forces structure on the closed-form filters.
std::array< sssv, nparts > arrayVec
Definition: rbpf.h:976
void filter(const osv &data, const std::vector< std::function< const Mat(const nsssv &x1t, const sssv &x2t)> > &fs=std::vector< std::function< const Mat(const nsssv &x1t, const sssv &x2t)> >())
Filter!
Definition: rbpf.h:1097
virtual float_t logQEv(const sssv &x2t, const sssv &x2tm1, const osv &yt)=0
Evaluates the proposal density at time t > 1.
Eigen::Matrix< float_t, dimnss, dimnss > nsssMat
Definition: rbpf.h:364
Rao-Blackwellized/Marginal Bootstrap Filter with inner HMMs.
Definition: rbpf.h:353
virtual float_t logQ1Ev(const sssv &x21, const osv &y1)=0
Evaluates the proposal density of the second state component at time 1.
virtual sssv muSamp()=0
Sample from the first time's proposal distribution.
Eigen::Matrix< float_t, dimy, 1 > osv
Definition: rbpf.h:36
arrayVec m_p_samps
Definition: rbpf.h:184
std::vector< Mat > getExpectations() const
Get the latest filtered expectation E[h(x_1t, x_2t) | y_{1:t}].
Definition: rbpf.h:938
arrayMod m_p_innerMods
Definition: rbpf.h:467
unsigned int m_rs
Definition: rbpf.h:782
arrayfloat_t m_logUnNormWeights
Definition: rbpf.h:1070
virtual sssv fSamp(const sssv &x2tm1)=0
Samples the time t second component.