1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
|
/**
* @file msvmmaj_kernel.c
* @author Gertjan van den Burg
* @date October 18, 2013
* @brief Defines main functions for use of kernels in MSVMMaj.
*
* @details
* Functions for constructing different kernels using user-supplied
* parameters. Also contains the functions for decomposing the
* kernel matrix using several decomposition methods.
*
*/
#include <math.h>
#include "msvmmaj.h"
#include "msvmmaj_kernel.h"
#include "msvmmaj_lapack.h"
#include "msvmmaj_matrix.h"
#include "util.h"
/**
* @brief Create the kernel matrix
*
* Create a kernel matrix based on the specified kerneltype. Kernel parameters
* are assumed to be specified in the model.
*
* @param[in] model MajModel specifying the parameters
* @param[in] data MajData specifying the data.
*
*/
void msvmmaj_make_kernel(struct MajModel *model, struct MajData *data)
{
if (model->kerneltype == K_LINEAR)
return;
long i, j;
long n = model->n;
double value;
double *x1, *x2;
double *K = Calloc(double, n*n*sizeof(double));
for (i=0; i<n; i++) {
for (j=i; j<n; j++) {
x1 = &data->Z[i*(data->m+1)+1];
x2 = &data->Z[j*(data->m+1)+1];
if (model->kerneltype == K_POLY)
value = msvmmaj_compute_poly(x1, x2,
model->kernelparam, data->m);
else if (model->kerneltype == K_RBF)
value = msvmmaj_compute_rbf(x1, x2,
model->kernelparam, data->m);
else if (model->kerneltype == K_SIGMOID)
value = msvmmaj_compute_rbf(x1, x2,
model->kernelparam, data->m);
else {
fprintf(stderr, "Unknown kernel type in "
"msvmmaj_make_kernel\n");
exit(1);
}
matrix_set(K, n, i, j, value);
matrix_set(K, n, j, i, value);
}
}
// copy kernel to data
data->Z = realloc(data->Z, n*(n+1)*(sizeof(double)));
for (i=0; i<n; i++) {
for (j=0; j<n; j++)
matrix_set(data->Z, n+1, i, j+1,
matrix_get(K, n, i, j));
matrix_set(data->Z, n+1, i, 0, 1.0);
}
data->m = n;
// let data know what it's made of
data->kerneltype = model->kerneltype;
free(data->kernelparam);
switch (model->kerneltype) {
case K_LINEAR:
break;
case K_POLY:
data->kernelparam = Calloc(double, 3);
data->kernelparam[0] = model->kernelparam[0];
data->kernelparam[1] = model->kernelparam[1];
data->kernelparam[2] = model->kernelparam[2];
break;
case K_RBF:
data->kernelparam = Calloc(double, 1);
data->kernelparam[0] = model->kernelparam[0];
break;
case K_SIGMOID:
data->kernelparam = Calloc(double, 2);
data->kernelparam[0] = model->kernelparam[0];
data->kernelparam[1] = model->kernelparam[1];
}
model->m = n;
free(K);
}
/**
* @brief Compute the RBF kernel between two vectors
*
* @details
* The RBF kernel is computed between two vectors. This kernel is defined as
* @f[
* k(x_1, x_2) = \exp( -\gamma \| x_1 - x_2 \|^2 )
* @f]
* where @f$ \gamma @f$ is a kernel parameter specified.
*
* @param[in] x1 first vector
* @param[in] x2 second vector
* @param[in] kernelparam array of kernel parameters (gamma is first
* element)
* @param[in] n length of the vectors x1 and x2
* @returns kernel evaluation
*/
double msvmmaj_compute_rbf(double *x1, double *x2, double *kernelparam, long n)
{
long i;
double value = 0.0;
for (i=0; i<n; i++)
value += (x1[i] - x2[i]) * (x1[i] - x2[i]);
value *= -kernelparam[0];
return exp(value);
}
/**
* @brief Compute the polynomial kernel between two vectors
*
* @details
* The polynomial kernel is computed between two vectors. This kernel is
* defined as
* @f[
* k(x_1, x_2) = ( \gamma \langle x_1, x_2 \rangle + c)^d
* @f]
* where @f$ \gamma @f$, @f$ c @f$ and @f$ d @f$ are kernel parameters.
*
* @param[in] x1 first vector
* @param[in] x2 second vector
* @param[in] kernelparam array of kernel parameters (gamma, c, d)
* @param[in] n length of the vectors x1 and x2
* @returns kernel evaluation
*/
double msvmmaj_compute_poly(double *x1, double *x2, double *kernelparam, long n)
{
long i;
double value = 0.0;
for (i=0; i<n; i++)
value += x1[i]*x2[i];
value *= kernelparam[0];
value += kernelparam[1];
return pow(value, ((int) kernelparam[2]));
}
/**
* @brief Compute the sigmoid kernel between two vectors
*
* @details
* The sigmoid kernel is computed between two vectors. This kernel is defined
* as
* @f[
* k(x_1, x_2) = \tanh( \gamma \langle x_1 , x_2 \rangle + c)
* @f]
* where @f$ \gamma @f$ and @f$ c @f$ are kernel parameters.
*
* @param[in] x1 first vector
* @param[in] x2 second vector
* @param[in] kernelparam array of kernel parameters (gamma, c)
* @param[in] n length of the vectors x1 and x2
* @returns kernel evaluation
*/
double msvmmaj_compute_sigmoid(double *x1, double *x2, double *kernelparam, long n)
{
long i;
double value = 0.0;
for (i=0; i<n; i++)
value += x1[i]*x2[i];
value *= kernelparam[0];
value += kernelparam[1];
return tanh(value);
}
|
