17 template<
size_t dimstate,
size_t dimobs,
typename float_t>
23 using ssv = Eigen::Matrix<float_t,dimstate,1>;
25 using osv = Eigen::Matrix<float_t,dimstate,1>;
41 template<
size_t dimstate,
size_t dimobs,
typename float_t>
52 template<
size_t dimstate,
size_t dimobs,
size_t diminput,
typename float_t>
58 using ssv = Eigen::Matrix<float_t,dimstate,1>;
61 using osv = Eigen::Matrix<float_t,dimobs,1>;
64 using isv = Eigen::Matrix<float_t,diminput,1>;
67 using ssMat = Eigen::Matrix<float_t,dimstate,dimstate>;
70 using osMat = Eigen::Matrix<float_t,dimobs,dimobs>;
73 using siMat = Eigen::Matrix<float_t,dimstate,diminput>;
76 using oiMat = Eigen::Matrix<float_t,dimobs,diminput>;
116 ssv getFiltMean()
const;
123 ssMat getFiltVar()
const;
130 osv getPredYMean(
const ssMat &stateTrans,
132 const siMat &stateInptAffector,
133 const oiMat &obsInptAffector,
134 const isv &inputData)
const;
142 const ssMat &cholStateVar,
144 const osMat &cholObsVar)
const;
158 void update(
const osv &yt,
159 const ssMat &stateTrans,
160 const ssMat &cholStateVar,
161 const siMat &stateInptAffector,
162 const isv &inputData,
164 const oiMat &obsInptAffector,
165 const osMat &cholObsVar);
201 void updatePrior(
const ssMat &stateTransMat,
202 const ssMat &cholStateVar,
203 const siMat &stateInptAffector,
204 const isv &inputData);
215 void updatePosterior(
const osv &yt,
217 const oiMat &obsInptAffector,
218 const isv &inputData,
219 const osMat &cholObsVar);
223 template<
size_t dimstate,
size_t dimobs,
size_t diminput,
typename float_t>
226 , m_predMean(
ssv::Zero())
227 , m_predVar(
ssMat::Zero())
229 , m_pi(3.14159265358979)
234 template<
size_t dimstate,
size_t dimobs,
size_t diminput,
typename float_t>
240 ,
m_pi(3.14159265358979)
245 template<
size_t dimstate,
size_t dimobs,
size_t diminput,
typename float_t>
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)
255 ssMat Q = cholStateVar.transpose() * cholStateVar;
261 template<
size_t dimstate,
size_t dimobs,
size_t diminput,
typename float_t>
264 const oiMat &obsInptAffector,
265 const isv &inputData,
266 const osMat &cholObsVar)
268 osMat R = cholObsVar.transpose() * cholObsVar;
270 osMat symSigma = (sigma.transpose() + sigma )/2.0;
271 osMat siginv = symSigma.inverse();
273 osv obsPred = obsMat *
m_predMean + obsInptAffector * inputData;
274 osv innov = yt - obsPred;
279 osMat quadForm = innov.transpose() * siginv * innov;
280 osMat cholSig ( sigma.llt().matrixL() );
281 float_t logDet = 2.0*cholSig.diagonal().array().log().sum();
286 template<
size_t dimstate,
size_t dimobs,
size_t diminput,
typename float_t>
293 template<
size_t dimstate,
size_t dimobs,
size_t diminput,
typename float_t>
300 template<
size_t dimstate,
size_t dimobs,
size_t diminput,
typename float_t>
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,
314 const oiMat &obsInptAffector,
315 const osMat &cholObsVar)
321 this->
updatePosterior(yt, obsMat, obsInptAffector, inData, cholObsVar);
325 this->
updatePrior(stateTrans, cholStateVar, stateInptAffector, inData);
326 this->
updatePosterior(yt, obsMat, obsInptAffector, inData, cholObsVar);
331 template<
size_t dimstate,
size_t dimobs,
size_t diminput,
typename float_t>
333 const ssMat &stateTrans,
335 const siMat &stateInptAffector,
336 const oiMat &obsInptAffector,
337 const isv &futureInputData)
const ->
osv 339 return obsMat * (stateTrans *
m_filtMean + stateInptAffector * futureInputData) + obsInptAffector * futureInputData;
343 template<
size_t dimstate,
size_t dimobs,
size_t diminput,
typename float_t>
345 const ssMat &stateTrans,
346 const ssMat &cholStateVar,
350 return obsMat * (stateTrans *
m_filtVar * stateTrans.transpose() + cholStateVar.transpose()*cholStateVar) * obsMat.transpose() + cholObsVar.transpose() * cholObsVar;
361 template<
size_t dimstate,
size_t dimobs,
typename float_t>
368 using ssv = Eigen::Matrix<float_t,dimstate,1>;
371 using osv = Eigen::Matrix<float_t,dimobs,1>;
374 using ssMat = Eigen::Matrix<float_t,dimstate,dimstate>;
390 hmm(
const ssv &initStateDistr,
const ssMat &transMat);
411 ssv getFilterVec()
const;
439 template<
size_t dimstate,
size_t dimobs,
typename float_t>
442 , m_filtVec(
ssv::Zero())
443 , m_transMatTranspose(
ssMat::Zero())
444 , m_lastCondLike(0.0)
450 template<
size_t dimstate,
size_t dimobs,
typename float_t>
461 template<
size_t dimstate,
size_t dimobs,
typename float_t>
465 template<
size_t dimstate,
size_t dimobs,
typename float_t>
472 template<
size_t dimstate,
size_t dimobs,
typename float_t>
479 template<
size_t dimstate,
size_t dimobs,
typename float_t>
491 m_filtVec = m_filtVec.cwiseProduct( condDensVec );
511 template<
size_t dim_pred,
typename float_t>
518 using psv = Eigen::Matrix<float_t,dim_pred,1>;
521 using tsv = Eigen::Matrix<float_t,2,1>;
537 gamFilter(
const float_t &nOneTilde,
const float_t &dOneTilde);
570 void update(
const float_t& yt,
const psv &xt,
const psv& beta,
const float_t& sigmaSquared,
const float_t& delta);
597 template<
size_t dim_pred,
typename float_t>
600 , m_lastLogCondLike(0.0)
608 template<
size_t dim_pred,
typename float_t>
612 template<
size_t dim_pred,
typename float_t>
619 template<
size_t dim_pred,
typename float_t>
627 template<
size_t dim_pred,
typename float_t>
631 if(sigmaSquared <= 0 || delta <= 0)
632 throw std::invalid_argument(
"ME: both sigma squared and delta have to be positive.\n");
639 m_filtVec(1) += (yt - xt.dot(beta))*(yt - xt.dot(beta))/sigmaSquared;
649 m_filtVec(1) += (yt - xt.dot(beta))*(yt - xt.dot(beta))/sigmaSquared;
662 template<
size_t dim_obs,
size_t dim_pred,
typename float_t>
669 using psv = Eigen::Matrix<float_t,dim_pred,1>;
672 using bsm = Eigen::Matrix<float_t,dim_obs, dim_pred>;
675 using tsv = Eigen::Matrix<float_t,2,1>;
678 using osv = Eigen::Matrix<float_t, dim_obs, 1>;
681 using osm = Eigen::Matrix<float_t, dim_obs, dim_obs>;
723 void update(
const osv& yt,
const psv &xt,
const bsm& B,
const osm& Sigma,
const float_t& delta);
735 osv getFcastMean(
const psv &xtp1,
const bsm& B,
const osm& Sigma,
const float_t& delta);
747 osm getFcastCov(
const psv &xtp1,
const bsm& B,
const osm& Sigma,
const float_t& delta);
764 template<
size_t dim_obs,
size_t dim_pred,
typename float_t>
775 template<
size_t dim_obs,
size_t dim_pred,
typename float_t>
779 template<
size_t dim_obs,
size_t dim_pred,
typename float_t>
786 template<
size_t dim_obs,
size_t dim_pred,
typename float_t>
794 template<
size_t dim_obs,
size_t dim_pred,
typename float_t>
800 throw std::invalid_argument(
"ME: delta has to be positive (you're not even checking Sigma).\n");
808 m_filtVec(1) += (Sigma.ldlt().solve(yt - modeVec)).squaredNorm();
821 m_filtVec(1) += (Sigma.ldlt().solve(yt - modeVec)).squaredNorm();
826 template<
size_t dim_obs,
size_t dim_pred,
typename float_t>
834 template<
size_t dim_obs,
size_t dim_pred,
typename float_t>
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