pf
rbpf_hmm< nparts, dimnss, dimss, dimy, resamp_t, float_t > Class Template Referenceabstract

Rao-Blackwellized/Marginal Particle Filter with inner HMMs. More...

#include <rbpf.h>

Inheritance diagram for rbpf_hmm< nparts, dimnss, dimss, dimy, resamp_t, float_t >:
Collaboration diagram for rbpf_hmm< nparts, dimnss, dimss, dimy, resamp_t, float_t >:

Public Types

using sssv = Eigen::Matrix< float_t, dimss, 1 >
 
using nsssv = Eigen::Matrix< float_t, dimnss, 1 >
 
using osv = Eigen::Matrix< float_t, dimy, 1 >
 
using nsssMat = Eigen::Matrix< float_t, dimnss, dimnss >
 
using Mat = Eigen::Matrix< float_t, Eigen::Dynamic, Eigen::Dynamic >
 
using arrayMod = std::array< hmm< dimnss, dimy, float_t >, nparts >
 
using arrayVec = std::array< sssv, nparts >
 
using arrayfloat_t = std::array< float_t, nparts >
 
- Public Types inherited from rbpf_base< float_t, dimss, dimnss, dimy >
using osv = Eigen::Matrix< float_t, dimobs, 1 >
 
using sssv = Eigen::Matrix< float_t, dim_s_state, 1 >
 
using nsssv = Eigen::Matrix< float_t, dim_ns_state, 1 >
 
using Mat = Eigen::Matrix< float_t, Eigen::Dynamic, Eigen::Dynamic >
 
using func = std::function< const Mat(const nsssv &, const sssv &)>
 
using funcs = std::vector< func >
 

Public Member Functions

 rbpf_hmm (const unsigned int &resamp_sched=1)
 The constructor. More...
 
virtual ~rbpf_hmm ()
 The (virtual) destructor.
 
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. More...
 
float_t getLogCondLike () const
 Get the latest conditional likelihood. More...
 
std::vector< MatgetExpectations () const
 Get vector of expectations. More...
 
virtual float_t logMuEv (const sssv &x21)=0
 Evaluates the first time state density. More...
 
virtual sssv q1Samp (const osv &y1)=0
 Sample from the first sampler. More...
 
virtual nsssv initHMMProbVec (const sssv &x21)=0
 Provides the initial mean vector for each HMM filter object. More...
 
virtual nsssMat initHMMTransMat (const sssv &x21)=0
 Provides the transition matrix for each HMM filter object. More...
 
virtual sssv qSamp (const sssv &x2tm1, const osv &yt)=0
 Samples the time t second component. More...
 
virtual float_t logQ1Ev (const sssv &x21, const osv &y1)=0
 Evaluates the proposal density of the second state component at time 1. More...
 
virtual float_t logFEv (const sssv &x2t, const sssv &x2tm1)=0
 Evaluates the state transition density for the second state component. More...
 
virtual float_t logQEv (const sssv &x2t, const sssv &x2tm1, const osv &yt)=0
 Evaluates the proposal density at time t > 1. More...
 
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. More...
 
- Public Member Functions inherited from rbpf_base< float_t, dimss, dimnss, dimy >
virtual void filter (const osv &data, const funcs &fs=funcs())=0
 the filtering function that must be defined
 

Private Attributes

unsigned int m_now
 
float_t m_lastLogCondLike
 
unsigned int m_rs
 
arrayMod m_p_innerMods
 
arrayVec m_p_samps
 
arrayfloat_t m_logUnNormWeights
 
resamp_t m_resampler
 
std::vector< Matm_expectations
 

Additional Inherited Members

- Static Public Attributes inherited from rbpf_base< float_t, dimss, dimnss, dimy >
static constexpr unsigned int dim_sampled_state
 
static constexpr unsigned int dim_not_sampled_state
 
static constexpr unsigned int dim_obs
 

Detailed Description

template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t, typename float_t>
class rbpf_hmm< nparts, dimnss, dimss, dimy, resamp_t, float_t >

Rao-Blackwellized/Marginal Particle Filter with inner HMMs.

Author
t

Member Typedef Documentation

◆ arrayfloat_t

template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
using rbpf_hmm< nparts, dimnss, dimss, dimy, resamp_t, float_t >::arrayfloat_t = std::array<float_t,nparts>

array of weights

◆ arrayMod

template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
using rbpf_hmm< nparts, dimnss, dimss, dimy, resamp_t, float_t >::arrayMod = std::array<hmm<dimnss,dimy,float_t>,nparts>

array of model objects

◆ arrayVec

template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
using rbpf_hmm< nparts, dimnss, dimss, dimy, resamp_t, float_t >::arrayVec = std::array<sssv,nparts>

array of samples

◆ Mat

template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
using rbpf_hmm< nparts, dimnss, dimss, dimy, resamp_t, float_t >::Mat = Eigen::Matrix<float_t,Eigen::Dynamic,Eigen::Dynamic>

Dynamic size matrix

◆ nsssMat

template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
using rbpf_hmm< nparts, dimnss, dimss, dimy, resamp_t, float_t >::nsssMat = Eigen::Matrix<float_t,dimnss,dimnss>

"not sampled state size matrix"

◆ nsssv

template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
using rbpf_hmm< nparts, dimnss, dimss, dimy, resamp_t, float_t >::nsssv = Eigen::Matrix<float_t,dimnss,1>

"not sampled state size vector"

◆ osv

template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
using rbpf_hmm< nparts, dimnss, dimss, dimy, resamp_t, float_t >::osv = Eigen::Matrix<float_t,dimy,1>

"observation size vector"

◆ sssv

template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
using rbpf_hmm< nparts, dimnss, dimss, dimy, resamp_t, float_t >::sssv = Eigen::Matrix<float_t,dimss,1>

"sampled state size vector"

Constructor & Destructor Documentation

◆ rbpf_hmm()

template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
rbpf_hmm< nparts, dimnss, dimss, dimy, resamp_t, float_t >::rbpf_hmm ( const unsigned int &  resamp_sched = 1)

The constructor.

constructor.

Parameters
resamp_schedhow often to resample (e.g. once every resamp_sched time periods)

Member Function Documentation

◆ filter()

template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
void rbpf_hmm< nparts, dimnss, dimss, dimy, resamp_t, float_t >::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.

filters everything based on a new data point.

Parameters
datathe most recent time series observation.
fsa vector of functions computing E[h(x_1t, x_2t^i)| x_2t^i,y_1:t] to be averaged to yield E[h(x_1t, x_2t)|,y_1:t]. Will access the probability vector of x_1t

◆ getExpectations()

template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
auto rbpf_hmm< nparts, dimnss, dimss, dimy, resamp_t, float_t >::getExpectations ( ) const

Get vector of expectations.

Returns
vector of expectations

◆ getLogCondLike()

template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
float_t rbpf_hmm< nparts, dimnss, dimss, dimy, resamp_t, float_t >::getLogCondLike ( ) const

Get the latest conditional likelihood.

Get the latest conditional likelihood.

Returns
the latest conditional likelihood.

◆ initHMMProbVec()

template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
virtual nsssv rbpf_hmm< nparts, dimnss, dimss, dimy, resamp_t, float_t >::initHMMProbVec ( const sssv x21)
pure virtual

Provides the initial mean vector for each HMM filter object.

provides the initial probability vector for each HMM filter object.

Parameters
x21the second state componenent at time 1.
Returns
a Vec representing the probability of each state element.

◆ initHMMTransMat()

template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
virtual nsssMat rbpf_hmm< nparts, dimnss, dimss, dimy, resamp_t, float_t >::initHMMTransMat ( const sssv x21)
pure virtual

Provides the transition matrix for each HMM filter object.

provides the transition matrix for each HMM filter object.

Parameters
x21the second state component at time 1.
Returns
a transition matrix where element (ij) is the probability of transitioning from state i to state j.

◆ logFEv()

template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
virtual float_t rbpf_hmm< nparts, dimnss, dimss, dimy, resamp_t, float_t >::logFEv ( const sssv x2t,
const sssv x2tm1 
)
pure virtual

Evaluates the state transition density for the second state component.

Evaluates the state transition density for the second state component.

Parameters
x2tthe current second state component.
x2tm1the previous second state component.
Returns
a float_t evaluation.

◆ logMuEv()

template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
virtual float_t rbpf_hmm< nparts, dimnss, dimss, dimy, resamp_t, float_t >::logMuEv ( const sssv x21)
pure virtual

Evaluates the first time state density.

evaluates mu.

Parameters
x21component two at time 1
Returns
a float_t evaluation

◆ logQ1Ev()

template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
virtual float_t rbpf_hmm< nparts, dimnss, dimss, dimy, resamp_t, float_t >::logQ1Ev ( const sssv x21,
const osv y1 
)
pure virtual

Evaluates the proposal density of the second state component at time 1.

Evaluates the proposal density of the second state component at time 1.

Parameters
x21the second state component at time 1 you sampled.
y1time 1 observation.
Returns
a float_t evaluation of the density.

◆ logQEv()

template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
virtual float_t rbpf_hmm< nparts, dimnss, dimss, dimy, resamp_t, float_t >::logQEv ( const sssv x2t,
const sssv x2tm1,
const osv yt 
)
pure virtual

Evaluates the proposal density at time t > 1.

Evaluates the proposal density at time t > 1.

Parameters
x2tthe current second state component.
x2tm1the previous second state component.
ytthe current time series observation.
Returns
a float_t evaluation.

◆ q1Samp()

template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
virtual sssv rbpf_hmm< nparts, dimnss, dimss, dimy, resamp_t, float_t >::q1Samp ( const osv y1)
pure virtual

Sample from the first sampler.

samples the second component of the state at time 1.

Parameters
y1most recent datum.
Returns
a sssv sample for x21.

◆ qSamp()

template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
virtual sssv rbpf_hmm< nparts, dimnss, dimss, dimy, resamp_t, float_t >::qSamp ( const sssv x2tm1,
const osv yt 
)
pure virtual

Samples the time t second component.

Samples the time t second component.

Parameters
x2tm1the previous time's second state component.
ytthe current observation.
Returns
a Vec sample of the second state component at the current time.

◆ updateHMM()

template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
virtual void rbpf_hmm< nparts, dimnss, dimss, dimy, resamp_t, float_t >::updateHMM ( hmm< dimnss, dimy, float_t > &  aModel,
const osv yt,
const sssv x2t 
)
pure virtual

How to update your inner HMM filter object at each time.

How to update your inner HMM filter object at each time.

Parameters
aModela HMM filter object describing the conditional closed-form model.
ytthe current time series observation.
x2tthe current second state component.

Member Data Documentation

◆ m_expectations

template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
std::vector<Mat> rbpf_hmm< nparts, dimnss, dimss, dimy, resamp_t, float_t >::m_expectations
private

the vector of expectations

◆ m_lastLogCondLike

template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
float_t rbpf_hmm< nparts, dimnss, dimss, dimy, resamp_t, float_t >::m_lastLogCondLike
private

last conditional likelihood

◆ m_logUnNormWeights

template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
arrayfloat_t rbpf_hmm< nparts, dimnss, dimss, dimy, resamp_t, float_t >::m_logUnNormWeights
private

the array of unnormalized log-weights

◆ m_now

template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
unsigned int rbpf_hmm< nparts, dimnss, dimss, dimy, resamp_t, float_t >::m_now
private

the current time period

◆ m_p_innerMods

template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
arrayMod rbpf_hmm< nparts, dimnss, dimss, dimy, resamp_t, float_t >::m_p_innerMods
private

the array of inner closed-form models

◆ m_p_samps

template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
arrayVec rbpf_hmm< nparts, dimnss, dimss, dimy, resamp_t, float_t >::m_p_samps
private

the array of samples for the second state portion

◆ m_resampler

template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
resamp_t rbpf_hmm< nparts, dimnss, dimss, dimy, resamp_t, float_t >::m_resampler
private

the resampler object

◆ m_rs

template<size_t nparts, size_t dimnss, size_t dimss, size_t dimy, typename resamp_t , typename float_t >
unsigned int rbpf_hmm< nparts, dimnss, dimss, dimy, resamp_t, float_t >::m_rs
private

resampling schedue


The documentation for this class was generated from the following file: