hops
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
DegenerateMultivariateGaussian.hpp
Go to the documentation of this file.
1 #ifndef HOPS_DEGENERATEMULTIVARIATEGAUSSIAN_HPP
2 #define HOPS_DEGENERATEMULTIVARIATEGAUSSIAN_HPP
3 
4 #include <Eigen/Cholesky>
5 #include <Eigen/Core>
6 #include <Eigen/LU>
7 
8 #define _USE_MATH_DEFINES
9 #include <math.h> // Using deprecated math for windows
10 #include <utility>
11 
12 #include "MultivariateGaussian.hpp"
15 
16 
17 namespace hops {
19  public:
20 
22  std::vector<long> inactive = std::vector<long>(0));
23 
24  [[nodiscard]] MatrixType::Scalar computeNegativeLogLikelihood(const VectorType &x) const override;
25 
26  [[nodiscard]] std::optional<VectorType> computeLogLikelihoodGradient(const VectorType &x) const override;
27 
28  [[nodiscard]] std::optional<MatrixType> computeExpectedFisherInformation(const VectorType &) const override;
29 
30  private:
31  std::optional<MultivariateGaussian> gaussian;
32  std::vector<long> inactive;
33 
34  void removeRow(Eigen::MatrixXd &matrix, unsigned int rowToRemove) const {
35  unsigned int numRows = matrix.rows() - 1;
36  unsigned int numCols = matrix.cols();
37 
38  if (rowToRemove < numRows) {
39  matrix.block(rowToRemove, 0, numRows - rowToRemove, numCols) = matrix.bottomRows(numRows - rowToRemove);
40  }
41 
42  matrix.conservativeResize(numRows, numCols);
43  }
44 
45  void removeColumn(Eigen::MatrixXd &matrix, unsigned int colToRemove) const {
46  unsigned int numRows = matrix.rows();
47  unsigned int numCols = matrix.cols() - 1;
48 
49  if (colToRemove < numCols) {
50  matrix.block(0, colToRemove, numRows, numCols - colToRemove) = matrix.rightCols(numCols - colToRemove);
51  }
52 
53  matrix.conservativeResize(numRows, numCols);
54  }
55 
56  void removeRow(Eigen::VectorXd &vector, unsigned int rowToRemove) const {
57  unsigned int numRows = vector.rows() - 1;
58 
59  if (rowToRemove < numRows) {
60  vector.segment(rowToRemove, numRows - rowToRemove) = vector.tail(numRows - rowToRemove);
61  }
62 
63  vector.conservativeResize(numRows);
64  }
65 
66  void stripInactive(Eigen::MatrixXd &matrix) const {
67  for (auto &i : inactive) {
68  removeRow(matrix, i);
69  removeColumn(matrix, i);
70  }
71  }
72 
73  void stripInactive(Eigen::VectorXd &vector) const {
74  for (auto &i : inactive) {
75  removeRow(vector, i);
76  }
77  }
78  };
79 
81  MatrixType covariance,
82  std::vector<long> inactive) :
83  inactive(std::move(inactive)) {
84  stripInactive(mean);
85  stripInactive(covariance);
86  gaussian = MultivariateGaussian(mean, covariance);
87  }
88 
89  MatrixType::Scalar
91  VectorType _x = x;
92  stripInactive(_x);
93  return gaussian.value().computeNegativeLogLikelihood(_x);
94  }
95 
96  std::optional<MatrixType>
98  // Saves performance and skips stripping x here, because the FIM is constant anyways.
99  return gaussian.value().computeExpectedFisherInformation(x);
100  }
101 
102  std::optional<VectorType>
104  VectorType _x = x;
105  stripInactive(_x);
106  return gaussian.value().computeLogLikelihoodGradient(_x);
107  }
108 }
109 
110 #endif //HOPS_DEGENERATEMULTIVARIATEGAUSSIAN_HPP
hops::MatrixType
Eigen::MatrixXd MatrixType
Definition: MatrixType.hpp:7
hops::DegenerateMultivariateGaussian::computeExpectedFisherInformation
std::optional< MatrixType > computeExpectedFisherInformation(const VectorType &) const override
Definition: DegenerateMultivariateGaussian.hpp:97
hops::DegenerateMultivariateGaussian
Definition: DegenerateMultivariateGaussian.hpp:18
hops::DegenerateMultivariateGaussian::DegenerateMultivariateGaussian
DegenerateMultivariateGaussian(VectorType mean, MatrixType covariance, std::vector< long > inactive=std::vector< long >(0))
Definition: DegenerateMultivariateGaussian.hpp:80
hops::MultivariateGaussian
Definition: MultivariateGaussian.hpp:13
MultivariateGaussian.hpp
hops
Definition: CsvReader.hpp:8
hops::DegenerateMultivariateGaussian::computeLogLikelihoodGradient
std::optional< VectorType > computeLogLikelihoodGradient(const VectorType &x) const override
Definition: DegenerateMultivariateGaussian.hpp:103
VectorType.hpp
MatrixType.hpp
hops::Model
Definition: Model.hpp:12
hops::VectorType
Eigen::VectorXd VectorType
Definition: VectorType.hpp:7
hops::DegenerateMultivariateGaussian::computeNegativeLogLikelihood
MatrixType::Scalar computeNegativeLogLikelihood(const VectorType &x) const override
Evaluates the negative log likelihood for input x.
Definition: DegenerateMultivariateGaussian.hpp:90