Stan Math Library  2.20.0
reverse mode automatic differentiation
gp_exp_quad_cov.hpp
Go to the documentation of this file.
1 #ifndef STAN_MATH_PRIM_MAT_FUN_GP_EXP_QUAD_COV_HPP
2 #define STAN_MATH_PRIM_MAT_FUN_GP_EXP_QUAD_COV_HPP
3 
14 #include <cmath>
15 #include <vector>
16 
17 namespace stan {
18 namespace math {
34 template <typename T_x, typename T_sigma, typename T_l>
35 inline
36  typename Eigen::Matrix<typename stan::return_type<T_x, T_sigma, T_l>::type,
37  Eigen::Dynamic, Eigen::Dynamic>
38  gp_exp_quad_cov(const std::vector<T_x> &x, const T_sigma &sigma,
39  const T_l &length_scale) {
40  using std::exp;
41  check_positive("gp_exp_quad_cov", "magnitude", sigma);
42  check_positive("gp_exp_quad_cov", "length scale", length_scale);
43 
44  size_t x_size = x.size();
45  Eigen::Matrix<typename stan::return_type<T_x, T_sigma, T_l>::type,
46  Eigen::Dynamic, Eigen::Dynamic>
47  cov(x_size, x_size);
48 
49  if (x_size == 0)
50  return cov;
51 
52  for (size_t n = 0; n < x.size(); ++n)
53  check_not_nan("gp_exp_quad_cov", "x", x[n]);
54 
55  T_sigma sigma_sq = square(sigma);
56  T_l neg_half_inv_l_sq = -0.5 / square(length_scale);
57 
58  for (size_t j = 0; j < x_size; ++j) {
59  cov(j, j) = sigma_sq;
60  for (size_t i = j + 1; i < x_size; ++i) {
61  cov(i, j)
62  = sigma_sq * exp(squared_distance(x[i], x[j]) * neg_half_inv_l_sq);
63  cov(j, i) = cov(i, j);
64  }
65  }
66  return cov;
67 }
68 
83 template <typename T_x, typename T_sigma, typename T_l>
84 inline
85  typename Eigen::Matrix<typename stan::return_type<T_x, T_sigma, T_l>::type,
86  Eigen::Dynamic, Eigen::Dynamic>
87  gp_exp_quad_cov(const std::vector<Eigen::Matrix<T_x, Eigen::Dynamic, 1>> &x,
88  const T_sigma &sigma,
89  const std::vector<T_l> &length_scale) {
90  using std::exp;
91 
92  check_positive_finite("gp_exp_quad_cov", "magnitude", sigma);
93  check_positive_finite("gp_exp_quad_cov", "length scale", length_scale);
94 
95  size_t x_size = x.size();
96  Eigen::Matrix<typename stan::return_type<T_x, T_sigma, T_l>::type,
97  Eigen::Dynamic, Eigen::Dynamic>
98  cov(x_size, x_size);
99  if (x_size == 0)
100  return cov;
101 
102  size_t l_size = length_scale.size();
103  check_size_match("gp_exp_quad_cov", "x dimension", x[0].size(),
104  "number of length scales", l_size);
105 
106  T_sigma sigma_sq = square(sigma);
107  std::vector<
108  Eigen::Matrix<typename return_type<T_x, T_l>::type, Eigen::Dynamic, 1>>
109  x_new = divide_columns(x, length_scale);
110 
111  for (size_t j = 0; j < x_size; ++j) {
112  cov(j, j) = sigma_sq;
113  for (size_t i = j + 1; i < x_size; ++i) {
114  cov(i, j) = sigma_sq * exp(-0.5 * squared_distance(x_new[i], x_new[j]));
115  cov(j, i) = cov(i, j);
116  }
117  }
118  return cov;
119 }
120 
141 template <typename T_x1, typename T_x2, typename T_sigma, typename T_l>
142 inline typename Eigen::Matrix<
143  typename stan::return_type<T_x1, T_x2, T_sigma, T_l>::type, Eigen::Dynamic,
144  Eigen::Dynamic>
145 gp_exp_quad_cov(const std::vector<T_x1> &x1, const std::vector<T_x2> &x2,
146  const T_sigma &sigma, const T_l &length_scale) {
147  using std::exp;
148 
149  const char *function_name = "gp_exp_quad_cov";
150  check_positive(function_name, "magnitude", sigma);
151  check_positive(function_name, "length scale", length_scale);
152 
153  size_t x1_size = x1.size();
154  size_t x2_size = x2.size();
155  Eigen::Matrix<typename stan::return_type<T_x1, T_x2, T_sigma, T_l>::type,
156  Eigen::Dynamic, Eigen::Dynamic>
157  cov(x1_size, x2_size);
158  if (x1_size == 0 || x2_size == 0)
159  return cov;
160 
161  for (size_t i = 0; i < x1_size; ++i)
162  check_not_nan(function_name, "x1", x1[i]);
163  for (size_t i = 0; i < x2_size; ++i)
164  check_not_nan(function_name, "x2", x2[i]);
165 
166  T_sigma sigma_sq = square(sigma);
167  T_l neg_half_inv_l_sq = -0.5 / square(length_scale);
168 
169  for (size_t i = 0; i < x1.size(); ++i) {
170  for (size_t j = 0; j < x2.size(); ++j) {
171  cov(i, j)
172  = sigma_sq * exp(squared_distance(x1[i], x2[j]) * neg_half_inv_l_sq);
173  }
174  }
175  return cov;
176 }
177 
197 template <typename T_x1, typename T_x2, typename T_s, typename T_l>
198 inline typename Eigen::Matrix<
199  typename stan::return_type<T_x1, T_x2, T_s, T_l>::type, Eigen::Dynamic,
200  Eigen::Dynamic>
201 gp_exp_quad_cov(const std::vector<Eigen::Matrix<T_x1, Eigen::Dynamic, 1>> &x1,
202  const std::vector<Eigen::Matrix<T_x2, Eigen::Dynamic, 1>> &x2,
203  const T_s &sigma, const std::vector<T_l> &length_scale) {
204  using std::exp;
205  size_t x1_size = x1.size();
206  size_t x2_size = x2.size();
207  size_t l_size = length_scale.size();
208 
209  Eigen::Matrix<typename stan::return_type<T_x1, T_x2, T_s, T_l>::type,
210  Eigen::Dynamic, Eigen::Dynamic>
211  cov(x1_size, x2_size);
212  if (x1_size == 0 || x2_size == 0)
213  return cov;
214 
215  const char *function_name = "gp_exp_quad_cov";
216  for (size_t i = 0; i < x1_size; ++i)
217  check_not_nan(function_name, "x1", x1[i]);
218  for (size_t i = 0; i < x2_size; ++i)
219  check_not_nan(function_name, "x2", x2[i]);
220  check_positive_finite(function_name, "magnitude", sigma);
221  check_positive_finite(function_name, "length scale", length_scale);
222  check_size_match(function_name, "x dimension", x1[0].size(),
223  "number of length scales", l_size);
224  check_size_match(function_name, "x dimension", x2[0].size(),
225  "number of length scales", l_size);
226 
227  T_s sigma_sq = square(sigma);
228 
229  std::vector<Eigen::Matrix<typename return_type<T_x1, T_l, T_s>::type,
230  Eigen::Dynamic, 1>>
231  x1_new = divide_columns(x1, length_scale);
232  std::vector<Eigen::Matrix<typename return_type<T_x2, T_l, T_s>::type,
233  Eigen::Dynamic, 1>>
234  x2_new = divide_columns(x2, length_scale);
235 
236  for (size_t i = 0; i < x1_size; ++i) {
237  for (size_t j = 0; j < x2_size; ++j) {
238  cov(i, j) = sigma_sq * exp(-0.5 * squared_distance(x1_new[i], x2_new[j]));
239  }
240  }
241  return cov;
242 }
243 } // namespace math
244 } // namespace stan
245 #endif
Eigen::Matrix< typename stan::return_type< T_x, T_sigma, T_l >::type, Eigen::Dynamic, Eigen::Dynamic > gp_exp_quad_cov(const std::vector< T_x > &x, const T_sigma &sigma, const T_l &length_scale)
Returns a squared exponential kernel.
void check_size_match(const char *function, const char *name_i, T_size1 i, const char *name_j, T_size2 j)
Check if the provided sizes match.
fvar< T > square(const fvar< T > &x)
Definition: square.hpp:12
void check_positive_finite(const char *function, const char *name, const T_y &y)
Check if y is positive and finite.
std::vector< Eigen::Matrix< typename return_type< T_x, T_v, double >::type, Eigen::Dynamic, 1 > > divide_columns(const std::vector< Eigen::Matrix< T_x, Eigen::Dynamic, 1 >> &x, const std::vector< T_v > &vec)
Takes Stan data type vector[n] x[D] and divides column vector in x element-wise by the values in vec...
fvar< T > squared_distance(const Eigen::Matrix< fvar< T >, R, C > &v1, const Eigen::Matrix< double, R, C > &v2)
Returns the squared distance between the specified vectors of the same dimensions.
boost::math::tools::promote_args< double, typename scalar_type< T >::type, typename return_type< Types_pack... >::type >::type type
Definition: return_type.hpp:36
fvar< T > exp(const fvar< T > &x)
Definition: exp.hpp:11
void check_not_nan(const char *function, const char *name, const T_y &y)
Check if y is not NaN.
int size(const std::vector< T > &x)
Return the size of the specified standard vector.
Definition: size.hpp:17
void check_positive(const char *function, const char *name, const T_y &y)
Check if y is positive.

     [ Stan Home Page ] © 2011–2018, Stan Development Team.