1 #ifndef BOOTSTRAP_FILTER_H 2 #define BOOTSTRAP_FILTER_H 23 template<
size_t nparts,
size_t dimx,
size_t dimy,
size_t dimcov,
typename resamp_t,
typename float_t>
29 using ssv = Eigen::Matrix<float_t, dimx, 1>;
31 using osv = Eigen::Matrix<float_t, dimy, 1>;
33 using cvsv = Eigen::Matrix<float_t,dimcov,1>;
35 using Mat = Eigen::Matrix<float_t,Eigen::Dynamic,Eigen::Dynamic>;
41 using Funcs = std::vector<std::function<const Mat(const ssv&, const cvsv&)>>;
156 template<
size_t nparts,
size_t dimx,
size_t dimy,
size_t dimcov, typename resamp_t, typename float_t>
159 , m_logLastCondLike(0.0)
163 std::fill(m_logUnNormWeights.begin(), m_logUnNormWeights.end(), 0.0);
167 template<
size_t nparts,
size_t dimx,
size_t dimy,
size_t dimcov,
typename resamp_t,
typename float_t>
171 template<
size_t nparts,
size_t dimx,
size_t dimy,
size_t dimcov,
typename resamp_t,
typename float_t>
178 for(
size_t ii = 0; ii < nparts; ++ii)
190 for(
size_t i = 0; i < nparts; ++i){
202 unsigned int rows = testOutput.rows();
203 unsigned int cols = testOutput.cols();
204 Mat numer = Mat::Zero(rows,cols);
205 float_t weightNormConst (0.0);
206 for(
size_t prtcl = 0; prtcl < nparts; ++prtcl){
227 float_t maxOldLogUnNormWts(-std::numeric_limits<float_t>::infinity());
229 for(
size_t ii = 0; ii < nparts; ++ii)
245 float_t sumExp1(0.0);
246 float_t sumExp2(0.0);
247 for(
size_t i = 0; i < nparts; ++i){
249 sumExp2 += std::exp(oldLogUnNormWts[i] - maxOldLogUnNormWts);
251 m_logLastCondLike = maxNumer + std::log(sumExp1) - maxOldLogUnNormWts - std::log(sumExp2);
258 unsigned int rows = testOutput.rows();
259 unsigned int cols = testOutput.cols();
260 Mat numer = Mat::Zero(rows,cols);
261 float_t weightNormConst (0.0);
262 for(
size_t prtcl = 0; prtcl < nparts; ++prtcl){
280 template<
size_t nparts,
size_t dimx,
size_t dimy,
size_t dimcov,
typename resamp_t,
typename float_t>
287 template<
size_t nparts,
size_t dimx,
size_t dimy,
size_t dimcov,
typename resamp_t,
typename float_t>
295 #endif // BOOTSTRAP_FILTER_H virtual float_t logGEv(const osv &yt, const ssv &xt, const cvsv &zt)=0
Calculate gEv or logGEv.
std::vector< std::function< const Mat(const ssv &, const cvsv &)> > Funcs
Definition: bootstrap_filter_with_covariates.h:41
Eigen::Matrix< float_t, dimy, 1 > osv
Definition: bootstrap_filter_with_covariates.h:31
virtual ssv q1Samp(const osv &y1, const cvsv &z1)=0
Samples from time 1 proposal.
float_t m_logLastCondLike
log p(y_t|y_{1:t-1}) or log p(y1)
Definition: bootstrap_filter_with_covariates.h:138
A base class for the bootstrap particle filter with covariates.
Definition: bootstrap_filter_with_covariates.h:24
unsigned int m_resampSched
resampling schedule (e.g. resample every __ time points)
Definition: bootstrap_filter_with_covariates.h:147
virtual ssv fSamp(const ssv &xtm1, const cvsv &zt)=0
Sample from the state transition distribution.
Eigen::Matrix< float_t, dimx, 1 > ssv
Definition: bootstrap_filter_with_covariates.h:29
std::vector< Mat > m_expectations
expectations E[h(x_t) | y_{1:t}] for user defined "h"s
Definition: bootstrap_filter_with_covariates.h:144
std::array< float_t, nparts > arrayfloat_t
Definition: bootstrap_filter_with_covariates.h:39
Eigen::Matrix< float_t, dimcov, 1 > cvsv
Definition: bootstrap_filter_with_covariates.h:33
arrayfloat_t m_logUnNormWeights
particle unnormalized weights
Definition: bootstrap_filter_with_covariates.h:132
virtual float_t logQ1Ev(const ssv &x1, const osv &y1, const cvsv &z1)=0
Calculate q1Ev or log q1Ev.
virtual ~BSFilterWC()
The (virtual) destructor.
Definition: bootstrap_filter_with_covariates.h:168
float_t getLogCondLike() const
Returns the most recent (log-) conditiona likelihood.
Definition: bootstrap_filter_with_covariates.h:281
std::array< ssv, nparts > arrayStates
Definition: bootstrap_filter_with_covariates.h:37
resamp_t m_resampler
resampler object
Definition: bootstrap_filter_with_covariates.h:141
unsigned int m_now
time point
Definition: bootstrap_filter_with_covariates.h:135
All non Rao-Blackwellized particle filters inherit from this.
Eigen::Matrix< float_t, Eigen::Dynamic, Eigen::Dynamic > Mat
Definition: bootstrap_filter_with_covariates.h:35
void filter(const osv &ydata, const cvsv &covdata, const Funcs &fs=Funcs())
updates filtering distribution on a new datapoint. Optionally stores expectations of functionals...
Definition: bootstrap_filter_with_covariates.h:172
arrayStates m_particles
particle samples
Definition: bootstrap_filter_with_covariates.h:129
auto getExpectations() const -> std::vector< Mat >
return all stored expectations (taken with respect to $p(x_t|y_{1:t})$
Definition: bootstrap_filter_with_covariates.h:288
BSFilterWC(const unsigned int &rs=1)
The constructor.
Definition: bootstrap_filter_with_covariates.h:157
virtual float_t logMuEv(const ssv &x1, const cvsv &z1)=0
Calculate muEv or logmuEv.