A base-class for Auxiliary Particle Filtering. Filtering only, no smoothing.
More...
#include <auxiliary_pf.h>
|
using | ssv = Eigen::Matrix< double, dimx, 1 > |
|
using | osv = Eigen::Matrix< double, dimy, 1 > |
|
using | Mat = Eigen::Matrix< double, dimx, dimx > |
|
using | arrayDouble = std::array< double, nparts > |
|
using | arrayVec = std::array< ssv, nparts > |
|
using | arrayUInt = std::array< unsigned int, nparts > |
|
|
| APF (const unsigned int &rs=1) |
| The constructor. More...
|
|
double | getLogCondLike () const |
| Get the latest log conditional likelihood. More...
|
|
std::vector< Mat > | getExpectations () const |
| return all stored expectations (taken with respect to $p(x_t|y_{1:t})$ More...
|
|
void | filter (const osv &data, const std::vector< std::function< const Mat(const ssv &)> > &fs=std::vector< std::function< const Mat(const ssv &)> >()) |
| Use a new datapoint to update the filtering distribution (or smoothing if pathLength > 0). More...
|
|
virtual double | logMuEv (const ssv &x1)=0 |
| Evaluates the log of mu. More...
|
|
virtual ssv | propMu (const ssv &xtm1)=0 |
| Evaluates the proposal distribution taking a Eigen::Matrix<double,dimx,1> from the previous time's state, and returning a state for the current time. More...
|
|
virtual ssv | q1Samp (const osv &y1)=0 |
| Samples from q1. More...
|
|
virtual ssv | fSamp (const ssv &xtm1)=0 |
| Samples from f. More...
|
|
virtual double | logQ1Ev (const ssv &x1, const osv &y1)=0 |
| Evaluates the log of q1. More...
|
|
virtual double | logGEv (const osv &yt, const ssv &xt)=0 |
| Evaluates the log of g. More...
|
|
template<size_t nparts, size_t dimx, size_t dimy, typename resampT>
class pf::APF< nparts, dimx, dimy, resampT >
A base-class for Auxiliary Particle Filtering. Filtering only, no smoothing.
- Author
- taylor
template<size_t nparts, size_t dimx, size_t dimy, typename resampT >
type alias for array of doubles
template<size_t nparts, size_t dimx, size_t dimy, typename resampT >
using pf::APF< nparts, dimx, dimy, resampT >::arrayUInt = std::array<unsigned int, nparts> |
type alias for array of unsigned ints
template<size_t nparts, size_t dimx, size_t dimy, typename resampT >
type alias for array of state vectors
template<size_t nparts, size_t dimx, size_t dimy, typename resampT >
using pf::APF< nparts, dimx, dimy, resampT >::Mat = Eigen::Matrix<double,dimx,dimx> |
type alias for linear algebra stuff (dimension of the state ^2)
template<size_t nparts, size_t dimx, size_t dimy, typename resampT >
using pf::APF< nparts, dimx, dimy, resampT >::osv = Eigen::Matrix<double,dimy,1> |
"observation size vector" type alias for linear algebra stuff
template<size_t nparts, size_t dimx, size_t dimy, typename resampT >
using pf::APF< nparts, dimx, dimy, resampT >::ssv = Eigen::Matrix<double,dimx,1> |
"state size vector" type alias for linear algebra stuff
template<size_t nparts, size_t dimx, size_t dimy, typename resampT >
pf::APF< nparts, dimx, dimy, resampT >::APF |
( |
const unsigned int & |
rs = 1 | ) |
|
The constructor.
- Parameters
-
rs | resampling schedule (e.g. resample every rs time points). |
template<size_t nparts, size_t dimx, size_t dimy, typename resampT >
void pf::APF< nparts, dimx, dimy, resampT >::filter |
( |
const osv & |
data, |
|
|
const std::vector< std::function< const Mat(const ssv &)> > & |
fs = std::vector<std::function<const Mat(const ssv&)> >() |
|
) |
| |
Use a new datapoint to update the filtering distribution (or smoothing if pathLength > 0).
- Parameters
-
data | a Eigen::Matrix<double,dimy,1> representing the data |
fs | a std::vector of callback functions that are used to calculate expectations with respect to the filtering distribution. |
template<size_t nparts, size_t dimx, size_t dimy, typename resampT >
virtual ssv pf::APF< nparts, dimx, dimy, resampT >::fSamp |
( |
const ssv & |
xtm1 | ) |
|
|
pure virtual |
Samples from f.
- Parameters
-
xtm1 | a Eigen::Matrix<double,dimx,1> representing the previous time's state. |
- Returns
- a Eigen::Matrix<double,dimx,1> state sample for the current time.
template<size_t nparts, size_t dimx, size_t dimy, typename resampT >
auto pf::APF< nparts, dimx, dimy, resampT >::getExpectations |
( |
| ) |
const |
return all stored expectations (taken with respect to $p(x_t|y_{1:t})$
- Returns
- return a std::vector<Mat> of expectations. How many depends on how many callbacks you gave to
template<size_t nparts, size_t dimx, size_t dimy, typename resampT >
double pf::APF< nparts, dimx, dimy, resampT >::getLogCondLike |
( |
| ) |
const |
Get the latest log conditional likelihood.
- Returns
- a double of the most recent conditional likelihood.
template<size_t nparts, size_t dimx, size_t dimy, typename resampT >
virtual double pf::APF< nparts, dimx, dimy, resampT >::logGEv |
( |
const osv & |
yt, |
|
|
const ssv & |
xt |
|
) |
| |
|
pure virtual |
Evaluates the log of g.
- Parameters
-
yt | a Eigen::Matrix<double,dimy,1> representing time t's data observation. |
xt | a Eigen::Matrix<double,dimx,1> representing time t's state. |
- Returns
- a double evaluation.
template<size_t nparts, size_t dimx, size_t dimy, typename resampT >
virtual double pf::APF< nparts, dimx, dimy, resampT >::logMuEv |
( |
const ssv & |
x1 | ) |
|
|
pure virtual |
Evaluates the log of mu.
- Parameters
-
x1 | a Eigen::Matrix<double,dimx,1> representing time 1's state. |
- Returns
- a double evaluation.
template<size_t nparts, size_t dimx, size_t dimy, typename resampT >
virtual double pf::APF< nparts, dimx, dimy, resampT >::logQ1Ev |
( |
const ssv & |
x1, |
|
|
const osv & |
y1 |
|
) |
| |
|
pure virtual |
Evaluates the log of q1.
- Parameters
-
x1 | a Eigen::Matrix<double,dimx,1> representing time 1's state. |
y1 | a Eigen::Matrix<double,dimy,1> representing time 1's data observation. |
- Returns
- a double evaluation.
template<size_t nparts, size_t dimx, size_t dimy, typename resampT >
virtual ssv pf::APF< nparts, dimx, dimy, resampT >::propMu |
( |
const ssv & |
xtm1 | ) |
|
|
pure virtual |
Evaluates the proposal distribution taking a Eigen::Matrix<double,dimx,1> from the previous time's state, and returning a state for the current time.
- Parameters
-
xtm1 | a Eigen::Matrix<double,dimx,1> representing the previous time's state. |
- Returns
- a Eigen::Matrix<double,dimx,1> representing a likely current time state, to be used by the observation density.
template<size_t nparts, size_t dimx, size_t dimy, typename resampT >
virtual ssv pf::APF< nparts, dimx, dimy, resampT >::q1Samp |
( |
const osv & |
y1 | ) |
|
|
pure virtual |
Samples from q1.
- Parameters
-
y1 | a Eigen::Matrix<double,dimy,1> representing time 1's data point. |
- Returns
- a Eigen::Matrix<double,dimx,1> sample for time 1's state.
The documentation for this class was generated from the following file: