SHOGUN  4.2.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages
GaussianARDKernel.cpp
Go to the documentation of this file.
1 /*
2  * This program is free software; you can redistribute it and/or modify
3  * it under the terms of the GNU General Public License as published by
4  * the Free Software Foundation; either version 3 of the License, or
5  * (at your option) any later version.
6  *
7  * Written (W) 2015 Wu Lin
8  * Written (W) 2012 Jacob Walker
9  *
10  * Adapted from WeightedDegreeRBFKernel.cpp
11  */
12 
15 
16 #ifdef HAVE_LINALG_LIB
18 #endif
19 
20 using namespace shogun;
21 
23 {
24  init();
25 }
26 
28 {
29 }
30 
31 void CGaussianARDKernel::init()
32 {
33 
34 #ifdef HAVE_LINALG_LIB
35  m_sq_lhs=SGVector<float64_t>();
36  m_sq_rhs=SGVector<float64_t>();
37  SG_ADD(&m_sq_lhs, "sq_lhs", "squared left-hand side", MS_NOT_AVAILABLE);
38  SG_ADD(&m_sq_rhs, "sq_rhs", "squared right-hand side", MS_NOT_AVAILABLE);
39 #endif
40 }
41 
42 float64_t CGaussianARDKernel::distance(int32_t idx_a, int32_t idx_b)
43 {
44  float64_t result=0.0;
45 #ifdef HAVE_LINALG_LIB
46  REQUIRE(lhs, "Left features (lhs) not set!\n")
47  REQUIRE(rhs, "Right features (rhs) not set!\n")
48 
49  if (lhs==rhs && idx_a==idx_b)
50  return result;
51 
52  if (m_ARD_type==KT_SCALAR)
53  {
54  result=(m_sq_lhs[idx_a]+m_sq_rhs[idx_b]-2.0*CDotKernel::compute(idx_a,idx_b));
55  result*=CMath::exp(2.0*m_log_weights[0]);
56  }
57  else
58  {
61  avec=linalg::add(avec, bvec, 1.0, -1.0);
62  result=compute_helper(avec, avec);
63  }
64 #endif /* HAVE_LINALG_LIB */
65  return result/2.0;
66 }
67 
68 #ifdef HAVE_LINALG_LIB
70  : CExponentialARDKernel(size)
71 {
72  init();
73 }
74 
76  CDotFeatures* r, int32_t size)
77  : CExponentialARDKernel(size)
78 {
79  init();
80 }
81 
82 bool CGaussianARDKernel::init(CFeatures* l, CFeatures* r)
83 {
84  bool status=CExponentialARDKernel::init(l,r);
85 
86  if (m_ARD_type==KT_SCALAR)
87  precompute_squared();
88 
89  return status;
90 }
91 
92 SGVector<float64_t> CGaussianARDKernel::precompute_squared_helper(CDotFeatures* df)
93 {
94  REQUIRE(df, "Features not set\n")
95  int32_t num_vec=df->get_num_vectors();
96  SGVector<float64_t> sq(num_vec);
97  for (int32_t i=0; i<num_vec; i++)
98  sq[i]=df->dot(i,df, i);
99  return sq;
100 }
101 
102 void CGaussianARDKernel::precompute_squared()
103 {
104  if (!lhs || !rhs)
105  return;
106  m_sq_lhs=precompute_squared_helper((CDotFeatures*) lhs);
107 
108  if (lhs==rhs)
109  m_sq_rhs=m_sq_lhs;
110  else
111  m_sq_rhs=precompute_squared_helper((CDotFeatures*) rhs);
112 }
113 
114 
116 {
117  if (kernel->get_kernel_type()!=K_GAUSSIANARD)
118  {
119  SG_SERROR("Provided kernel is not of type CGaussianARDKernel!\n");
120  }
121 
122  /* since an additional reference is returned */
123  SG_REF(kernel);
124  return (CGaussianARDKernel*)kernel;
125 }
126 
127 float64_t CGaussianARDKernel::compute_helper(SGVector<float64_t> avec, SGVector<float64_t>bvec)
128 {
129  SGMatrix<float64_t> left;
130  SGMatrix<float64_t> left_transpose;
131  float64_t scalar_weight=1.0;
132  if (m_ARD_type==KT_SCALAR)
133  {
134  left=SGMatrix<float64_t>(avec.vector,1,avec.vlen,false);
135  scalar_weight=CMath::exp(m_log_weights[0]);
136  }
137  else if(m_ARD_type==KT_FULL || m_ARD_type==KT_DIAG)
138  {
139  left_transpose=get_weighted_vector(avec);
140  left=SGMatrix<float64_t>(left_transpose.matrix,1,left_transpose.num_rows,false);
141  }
142  else
143  SG_ERROR("Unsupported ARD type\n");
144  SGMatrix<float64_t> right=compute_right_product(bvec, scalar_weight);
145  SGMatrix<float64_t> res=linalg::matrix_product(left, right);
146  return res[0]*scalar_weight;
147 }
148 
149 float64_t CGaussianARDKernel::compute_gradient_helper(SGVector<float64_t> avec,
151 {
152  float64_t result=0.0;
153 
154  if(m_ARD_type==KT_DIAG)
155  result=2.0*avec[index]*bvec[index]*CMath::exp(2.0*m_log_weights[index]);
156  else
157  {
159 
160  if (m_ARD_type==KT_SCALAR)
161  {
162  SGMatrix<float64_t> left(avec.vector,1,avec.vlen,false);
163  SGMatrix<float64_t> right(bvec.vector,bvec.vlen,1,false);
164  res=linalg::matrix_product(left, right);
165  result=2.0*res[0]*CMath::exp(2.0*m_log_weights[0]);
166  }
167  else if(m_ARD_type==KT_FULL)
168  {
169  int32_t row_index=0;
170  int32_t col_index=index;
171  int32_t offset=m_weights_rows;
172  int32_t total_offset=0;
173  while(col_index>=offset && offset>0)
174  {
175  col_index-=offset;
176  total_offset+=offset;
177  offset--;
178  row_index++;
179  }
180  col_index+=row_index;
181 
182  SGVector<float64_t> row_vec=SGVector<float64_t>(m_log_weights.vector+total_offset,m_weights_rows-row_index,false);
183  row_vec[0]=CMath::exp(row_vec[0]);
184 
185  SGMatrix<float64_t> row_vec_r(row_vec.vector,row_vec.vlen,1,false);
186  SGMatrix<float64_t> left(avec.vector+row_index,1,avec.vlen-row_index,false);
187 
188  res=linalg::matrix_product(left, row_vec_r);
189  result=res[0]*bvec[col_index];
190 
191  SGMatrix<float64_t> row_vec_l(row_vec.vector,1,row_vec.vlen,false);
192  SGMatrix<float64_t> right(bvec.vector+row_index,bvec.vlen-row_index,1,false);
193 
194  res=linalg::matrix_product(row_vec_l, right);
195  result+=res[0]*avec[col_index];
196 
197  if(row_index==col_index)
198  result*=row_vec[0];
199  row_vec[0]=CMath::log(row_vec[0]);
200  }
201  else
202  {
203  SG_ERROR("Unsupported ARD type\n");
204  }
205 
206  }
207  return result*scale;
208 }
209 
210 
212  const TParameter* param, index_t index)
213 {
214  REQUIRE(param, "Param not set\n");
215  REQUIRE(lhs , "Left features not set!\n");
216  REQUIRE(rhs, "Right features not set!\n");
217 
218  if (lhs==rhs)
219  {
220  if (!strcmp(param->m_name, "log_weights"))
221  {
222  SGVector<float64_t> derivative(num_lhs);
223  derivative.zero();
224  return derivative;
225  }
226  }
227  else
228  {
229  int32_t length=CMath::min(num_lhs, num_rhs);
230  SGVector<float64_t> derivative(length);
231  check_weight_gradient_index(index);
232  for (index_t j=0; j<length; j++)
233  {
234  if (!strcmp(param->m_name, "log_weights") )
235  {
236  if (m_ARD_type==KT_SCALAR)
237  {
238  float64_t dist=distance(j,j);
239  derivative[j]=CMath::exp(-dist)*(-dist*2.0);
240  }
241  else
242  {
245  derivative[j]=get_parameter_gradient_helper(param,index,j,j,avec,bvec);
246  }
247 
248  }
249  }
250  return derivative;
251  }
252 
253  SG_ERROR("Can't compute derivative wrt %s parameter\n", param->m_name);
254  return SGVector<float64_t>();
255 }
256 
257 
258 float64_t CGaussianARDKernel::get_parameter_gradient_helper(
259  const TParameter* param, index_t index, int32_t idx_a,
260  int32_t idx_b, SGVector<float64_t> avec, SGVector<float64_t> bvec)
261 {
262  REQUIRE(param, "Param not set\n");
263 
264  if (!strcmp(param->m_name, "log_weights"))
265  {
266  bvec=linalg::add(avec, bvec, 1.0, -1.0);
267  float64_t scale=-kernel(idx_a,idx_b)/2.0;
268  return compute_gradient_helper(bvec, bvec, scale, index);
269  }
270  else
271  {
272  SG_ERROR("Can't compute derivative wrt %s parameter\n", param->m_name);
273  return 0.0;
274  }
275 }
276 
278  const TParameter* param, index_t index)
279 {
280  REQUIRE(param, "Param not set\n");
281  REQUIRE(lhs , "Left features not set!\n");
282  REQUIRE(rhs, "Right features not set!\n");
283 
284  if (!strcmp(param->m_name, "log_weights"))
285  {
286  SGMatrix<float64_t> derivative(num_lhs, num_rhs);
287  check_weight_gradient_index(index);
288  for (index_t j=0; j<num_lhs; j++)
289  {
291  for (index_t k=0; k<num_rhs; k++)
292  {
293  if (m_ARD_type==KT_SCALAR)
294  {
295  float64_t dist=distance(j,k);
296  derivative(j,k)=CMath::exp(-dist)*(-dist*2.0);
297  }
298  else
299  {
301  derivative(j,k)=get_parameter_gradient_helper(param,index,j,k,avec,bvec);
302  }
303  }
304  }
305  return derivative;
306  }
307  else
308  {
309  SG_ERROR("Can't compute derivative wrt %s parameter\n", param->m_name);
310  return SGMatrix<float64_t>();
311  }
312 }
313 #endif /* HAVE_LINALG_LIB */
SGVector< float64_t > m_log_weights
int32_t index_t
Definition: common.h:62
int32_t num_rhs
number of feature vectors on right hand side
Definition: Kernel.h:1070
Vector::Scalar dot(Vector a, Vector b)
Definition: Redux.h:58
virtual float64_t distance(int32_t idx_a, int32_t idx_b)
parameter struct
#define SG_ERROR(...)
Definition: SGIO.h:129
#define REQUIRE(x,...)
Definition: SGIO.h:206
float64_t kernel(int32_t idx_a, int32_t idx_b)
Definition: Kernel.h:207
virtual float64_t compute(int32_t idx_a, int32_t idx_b)
Definition: DotKernel.h:134
Features that support dot products among other operations.
Definition: DotFeatures.h:44
#define SG_REF(x)
Definition: SGObject.h:54
index_t num_rows
Definition: SGMatrix.h:374
Gaussian Kernel with Automatic Relevance Detection computed on CDotFeatures.
index_t vlen
Definition: SGVector.h:492
void add(Matrix A, Matrix B, Matrix C, typename Matrix::Scalar alpha=1.0, typename Matrix::Scalar beta=1.0)
Definition: Core.h:66
double float64_t
Definition: common.h:50
virtual SGVector< float64_t > get_feature_vector(int32_t idx, CFeatures *hs)
int32_t num_lhs
number of feature vectors on left hand side
Definition: Kernel.h:1068
CFeatures * rhs
feature vectors to occur on right hand side
Definition: Kernel.h:1062
static CKernel * obtain_from_generic(CSGObject *kernel)
Definition: Kernel.cpp:897
all of classes and functions are contained in the shogun namespace
Definition: class_list.h:18
virtual EKernelType get_kernel_type()=0
CFeatures * lhs
feature vectors to occur on left hand side
Definition: Kernel.h:1060
The class Features is the base class of all feature objects.
Definition: Features.h:68
static T min(T a, T b)
Definition: Math.h:157
#define SG_SERROR(...)
Definition: SGIO.h:179
void scale(Matrix A, Matrix B, typename Matrix::Scalar alpha)
Definition: Core.h:94
static float64_t exp(float64_t x)
Definition: Math.h:621
virtual SGMatrix< float64_t > get_parameter_gradient(const TParameter *param, index_t index=-1)
Definition: Kernel.h:851
static float64_t log(float64_t v)
Definition: Math.h:922
The Kernel base class.
Definition: Kernel.h:159
#define SG_ADD(...)
Definition: SGObject.h:84
Exponential Kernel with Automatic Relevance Detection computed on CDotFeatures.
virtual SGVector< float64_t > get_parameter_gradient_diagonal(const TParameter *param, index_t index=-1)
Definition: Kernel.h:865

SHOGUN Machine Learning Toolbox - Documentation