SHOGUN  4.2.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages
GaussianLikelihood.cpp
Go to the documentation of this file.
1 /*
2  * Copyright (c) The Shogun Machine Learning Toolbox
3  * Written (W) 2013 Roman Votyakov
4  * Written (W) 2012 Jacob Walker
5  * All rights reserved.
6  *
7  * Redistribution and use in source and binary forms, with or without
8  * modification, are permitted provided that the following conditions are met:
9  *
10  * 1. Redistributions of source code must retain the above copyright notice, this
11  * list of conditions and the following disclaimer.
12  * 2. Redistributions in binary form must reproduce the above copyright notice,
13  * this list of conditions and the following disclaimer in the documentation
14  * and/or other materials provided with the distribution.
15  *
16  * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND
17  * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
18  * WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
19  * DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR
20  * ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
21  * (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
22  * LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND
23  * ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
24  * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
25  * SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
26  *
27  * The views and conclusions contained in the software and documentation are those
28  * of the authors and should not be interpreted as representing official policies,
29  * either expressed or implied, of the Shogun Development Team.
30  *
31  */
33 
34 
37 
38 using namespace shogun;
39 using namespace Eigen;
40 
42 {
43  init();
44 }
45 
47 {
48  init();
49  set_sigma(sigma);
50 }
51 
52 void CGaussianLikelihood::init()
53 {
54  m_log_sigma=0.0;
55  SG_ADD(&m_log_sigma, "log_sigma", "Observation noise in log domain", MS_AVAILABLE, GRADIENT_AVAILABLE);
56 }
57 
59 {
60 }
61 
63  CLikelihoodModel* lik)
64 {
65  ASSERT(lik!=NULL);
66 
67  if (lik->get_model_type()!=LT_GAUSSIAN)
68  SG_SERROR("Provided likelihood is not of type CGaussianLikelihood!\n")
69 
70  SG_REF(lik);
71  return (CGaussianLikelihood*)lik;
72 }
73 
75  SGVector<float64_t> mu, SGVector<float64_t> s2, const CLabels* lab) const
76 {
77  return SGVector<float64_t>(mu);
78 }
79 
81  SGVector<float64_t> mu, SGVector<float64_t> s2, const CLabels* lab) const
82 {
83  SGVector<float64_t> result(s2);
84  Map<VectorXd> eigen_result(result.vector, result.vlen);
85 
86  eigen_result=eigen_result.array()+CMath::exp(m_log_sigma*2.0);
87 
88  return result;
89 }
90 
92  SGVector<float64_t> func) const
93 {
94  // check the parameters
95  REQUIRE(lab, "Labels are required (lab should not be NULL)\n")
97  "Labels must be type of CRegressionLabels\n")
98  REQUIRE(lab->get_num_labels()==func.vlen, "Number of labels must match "
99  "length of the function vector\n")
100 
101  Map<VectorXd> eigen_f(func.vector, func.vlen);
102 
103  SGVector<float64_t> result(func.vlen);
104  Map<VectorXd> eigen_result(result.vector, result.vlen);
105 
106  SGVector<float64_t> y=((CRegressionLabels*)lab)->get_labels();
107  Map<VectorXd> eigen_y(y.vector, y.vlen);
108 
109  // compute log probability: lp=-(y-f).^2./sigma^2/2-log(2*pi*sigma^2)/2
110  eigen_result=eigen_y-eigen_f;
111  eigen_result=-eigen_result.cwiseProduct(eigen_result)/(2.0*CMath::exp(m_log_sigma*2.0))-
112  VectorXd::Ones(result.vlen)*log(2.0*CMath::PI*CMath::exp(m_log_sigma*2.0))/2.0;
113 
114  return result;
115 }
116 
118  const CLabels* lab, SGVector<float64_t> func, index_t i) const
119 {
120  // check the parameters
121  REQUIRE(lab, "Labels are required (lab should not be NULL)\n")
123  "Labels must be type of CRegressionLabels\n")
124  REQUIRE(lab->get_num_labels()==func.vlen, "Number of labels must match "
125  "length of the function vector\n")
126  REQUIRE(i>=1 && i<=3, "Index for derivative should be 1, 2 or 3\n")
127 
128  Map<VectorXd> eigen_f(func.vector, func.vlen);
129 
130  SGVector<float64_t> result(func.vlen);
131  Map<VectorXd> eigen_result(result.vector, result.vlen);
132 
133  SGVector<float64_t> y=((CRegressionLabels*)lab)->get_labels();
134  Map<VectorXd> eigen_y(y.vector, y.vlen);
135 
136  // set result=y-f
137  eigen_result=eigen_y-eigen_f;
138 
139  // compute derivatives of log probability wrt f
140  if (i == 1)
141  eigen_result/=CMath::exp(m_log_sigma*2.0);
142  else if (i == 2)
143  eigen_result=-VectorXd::Ones(result.vlen)/CMath::exp(m_log_sigma*2.0);
144  else if (i == 3)
145  eigen_result=VectorXd::Zero(result.vlen);
146 
147  return result;
148 }
149 
151  SGVector<float64_t> func, const TParameter* param) const
152 {
153  REQUIRE(lab, "Labels are required (lab should not be NULL)\n")
155  "Labels must be type of CRegressionLabels\n")
156  REQUIRE(lab->get_num_labels()==func.vlen, "Number of labels must match "
157  "length of the function vector\n")
158 
159  Map<VectorXd> eigen_f(func.vector, func.vlen);
160 
161  SGVector<float64_t> result(func.vlen);
162  Map<VectorXd> eigen_result(result.vector, result.vlen);
163 
164  if (strcmp(param->m_name, "log_sigma"))
165  return SGVector<float64_t>();
166 
167  SGVector<float64_t> y=((CRegressionLabels*)lab)->get_labels();
168  Map<VectorXd> eigen_y(y.vector, y.vlen);
169 
170  // compute derivative of log probability wrt log_sigma:
171  // dlp_dlogsigma
172  // lp_dsigma=(y-f).^2/sigma^2-1
173  eigen_result=eigen_y-eigen_f;
174  eigen_result=eigen_result.cwiseProduct(eigen_result)/CMath::exp(m_log_sigma*2.0);
175  eigen_result-=VectorXd::Ones(result.vlen);
176 
177  return result;
178 }
179 
181  SGVector<float64_t> func, const TParameter* param) const
182 {
183  REQUIRE(lab, "Labels are required (lab should not be NULL)\n")
185  "Labels must be type of CRegressionLabels\n")
186  REQUIRE(lab->get_num_labels()==func.vlen, "Number of labels must match "
187  "length of the function vector\n")
188 
189  if (strcmp(param->m_name, "log_sigma"))
190  return SGVector<float64_t>();
191 
192  Map<VectorXd> eigen_f(func.vector, func.vlen);
193 
194  SGVector<float64_t> result(func.vlen);
195  Map<VectorXd> eigen_result(result.vector, result.vlen);
196 
197  SGVector<float64_t> y=((CRegressionLabels*)lab)->get_labels();
198  Map<VectorXd> eigen_y(y.vector, y.vlen);
199 
200  // compute derivative of (the first log_sigma derivative of log probability) wrt f:
201  // d2lp_dlogsigma_df == d2lp_df_dlogsigma
202  // dlp_dsigma=2*(f-y)/sigma^2
203  eigen_result=2.0*(eigen_f-eigen_y)/CMath::exp(m_log_sigma*2.0);
204 
205  return result;
206 }
207 
209  SGVector<float64_t> func, const TParameter* param) const
210 {
211  REQUIRE(lab, "Labels are required (lab should not be NULL)\n")
213  "Labels must be type of CRegressionLabels\n")
214  REQUIRE(lab->get_num_labels()==func.vlen, "Number of labels must match "
215  "length of the function vector\n")
216 
217  if (strcmp(param->m_name, "log_sigma"))
218  return SGVector<float64_t>();
219 
220  Map<VectorXd> eigen_f(func.vector, func.vlen);
221 
222  SGVector<float64_t> result(func.vlen);
223  Map<VectorXd> eigen_result(result.vector, result.vlen);
224 
225  // compute derivative of (the derivative of the first log_sigma derivative of log probability) wrt f:
226  // d3lp_dlogsigma_df_df == d3lp_df_df_dlogsigma
227  // d2lp_dsigma=2/sigma^2
228  eigen_result=2.0*VectorXd::Ones(result.vlen)/CMath::exp(m_log_sigma*2.0);
229 
230  return result;
231 }
232 
234  SGVector<float64_t> mu, SGVector<float64_t> s2, const CLabels *lab) const
235 {
237 
238  if (lab)
239  {
240  REQUIRE((mu.vlen==s2.vlen) && (mu.vlen==lab->get_num_labels()),
241  "Length of the vector of means (%d), length of the vector of "
242  "variances (%d) and number of labels (%d) should be the same\n",
243  mu.vlen, s2.vlen, lab->get_num_labels())
245  "Labels must be type of CRegressionLabels\n")
246 
247  y=((CRegressionLabels*)lab)->get_labels();
248  }
249  else
250  {
251  REQUIRE(mu.vlen==s2.vlen, "Length of the vector of means (%d) and "
252  "length of the vector of variances (%d) should be the same\n",
253  mu.vlen, s2.vlen)
254 
256  y.set_const(1.0);
257  }
258 
259  // create eigen representation of y, mu and s2
260  Map<VectorXd> eigen_mu(mu.vector, mu.vlen);
261  Map<VectorXd> eigen_s2(s2.vector, s2.vlen);
262  Map<VectorXd> eigen_y(y.vector, y.vlen);
263 
264  SGVector<float64_t> result(mu.vlen);
265  Map<VectorXd> eigen_result(result.vector, result.vlen);
266 
267  // compule lZ=-(y-mu).^2./(sn2+s2)/2-log(2*pi*(sn2+s2))/2
268  eigen_s2=eigen_s2.array()+CMath::exp(m_log_sigma*2.0);
269  eigen_result=-(eigen_y-eigen_mu).array().square()/(2.0*eigen_s2.array())-
270  (2.0*CMath::PI*eigen_s2.array()).log()/2.0;
271 
272  return result;
273 }
274 
276  SGVector<float64_t> s2, const CLabels *lab, index_t i) const
277 {
278  // check the parameters
279  REQUIRE(lab, "Labels are required (lab should not be NULL)\n")
280  REQUIRE((mu.vlen==s2.vlen) && (mu.vlen==lab->get_num_labels()),
281  "Length of the vector of means (%d), length of the vector of "
282  "variances (%d) and number of labels (%d) should be the same\n",
283  mu.vlen, s2.vlen, lab->get_num_labels())
284  REQUIRE(i>=0 && i<=mu.vlen, "Index (%d) out of bounds!\n", i)
286  "Labels must be type of CRegressionLabels\n")
287 
288  SGVector<float64_t> y=((CRegressionLabels*)lab)->get_labels();
289 
290  // compute 1st moment
291  float64_t Ex=mu[i]+s2[i]*(y[i]-mu[i])/(CMath::exp(m_log_sigma*2.0)+s2[i]);
292 
293  return Ex;
294 }
295 
297  SGVector<float64_t> s2, const CLabels *lab, index_t i) const
298 {
299  // check the parameters
300  REQUIRE(lab, "Labels are required (lab should not be NULL)\n")
301  REQUIRE((mu.vlen==s2.vlen) && (mu.vlen==lab->get_num_labels()),
302  "Length of the vector of means (%d), length of the vector of "
303  "variances (%d) and number of labels (%d) should be the same\n",
304  mu.vlen, s2.vlen, lab->get_num_labels())
305  REQUIRE(i>=0 && i<=mu.vlen, "Index (%d) out of bounds!\n", i)
307  "Labels must be type of CRegressionLabels\n")
308 
309  // compute 2nd moment
310  float64_t Var=s2[i]-CMath::sq(s2[i])/(CMath::exp(m_log_sigma*2.0)+s2[i]);
311 
312  return Var;
313 }
314 
virtual ELabelType get_label_type() const =0
Class that models Gaussian likelihood.
Real Labels are real-valued labels.
virtual SGVector< float64_t > get_log_probability_f(const CLabels *lab, SGVector< float64_t > func) const
int32_t index_t
Definition: common.h:62
The class Labels models labels, i.e. class assignments of objects.
Definition: Labels.h:43
virtual int32_t get_num_labels() const =0
real valued labels (e.g. for regression, classifier outputs)
Definition: LabelTypes.h:22
static T sq(T x)
Definition: Math.h:450
Definition: SGMatrix.h:20
virtual ELikelihoodModelType get_model_type() const
parameter struct
#define REQUIRE(x,...)
Definition: SGIO.h:206
void set_sigma(float64_t sigma)
#define SG_REF(x)
Definition: SGObject.h:54
virtual float64_t get_second_moment(SGVector< float64_t > mu, SGVector< float64_t > s2, const CLabels *lab, index_t i) const
virtual SGVector< float64_t > get_first_derivative(const CLabels *lab, SGVector< float64_t > func, const TParameter *param) const
index_t vlen
Definition: SGVector.h:492
virtual SGVector< float64_t > get_predictive_variances(SGVector< float64_t > mu, SGVector< float64_t > s2, const CLabels *lab=NULL) const
#define ASSERT(x)
Definition: SGIO.h:201
virtual SGVector< float64_t > get_predictive_means(SGVector< float64_t > mu, SGVector< float64_t > s2, const CLabels *lab=NULL) const
double float64_t
Definition: common.h:50
virtual float64_t get_first_moment(SGVector< float64_t > mu, SGVector< float64_t > s2, const CLabels *lab, index_t i) const
static CGaussianLikelihood * obtain_from_generic(CLikelihoodModel *lik)
virtual SGVector< float64_t > get_second_derivative(const CLabels *lab, SGVector< float64_t > func, const TParameter *param) const
virtual SGVector< float64_t > get_log_zeroth_moments(SGVector< float64_t > mu, SGVector< float64_t > s2, const CLabels *lab) const
all of classes and functions are contained in the shogun namespace
Definition: class_list.h:18
#define SG_SERROR(...)
Definition: SGIO.h:179
static float64_t exp(float64_t x)
Definition: Math.h:621
virtual SGVector< float64_t > get_log_probability_derivative_f(const CLabels *lab, SGVector< float64_t > func, index_t i) const
#define SG_ADD(...)
Definition: SGObject.h:84
virtual SGVector< float64_t > get_third_derivative(const CLabels *lab, SGVector< float64_t > func, const TParameter *param) const
The Likelihood model base class.
void set_const(T const_elem)
Definition: SGVector.cpp:150
static const float64_t PI
Definition: Math.h:2055

SHOGUN Machine Learning Toolbox - Documentation