Line data Source code
1 : /**
2 : * @file gensvm_optimize.c
3 : * @author G.J.J. van den Burg
4 : * @date 2013-08-09
5 : * @brief Main functions for training the GenSVM solution.
6 : *
7 : * @details
8 : * Contains update and loss functions used to actually find
9 : * the optimal V.
10 : *
11 : * @copyright
12 : Copyright 2016, G.J.J. van den Burg.
13 :
14 : This file is part of GenSVM.
15 :
16 : GenSVM is free software: you can redistribute it and/or modify
17 : it under the terms of the GNU General Public License as published by
18 : the Free Software Foundation, either version 3 of the License, or
19 : (at your option) any later version.
20 :
21 : GenSVM is distributed in the hope that it will be useful,
22 : but WITHOUT ANY WARRANTY; without even the implied warranty of
23 : MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
24 : GNU General Public License for more details.
25 :
26 : You should have received a copy of the GNU General Public License
27 : along with GenSVM. If not, see <http://www.gnu.org/licenses/>.
28 :
29 : */
30 :
31 : #include "gensvm_optimize.h"
32 :
33 : /**
34 : * Iteration frequency with which to print to stdout
35 : */
36 : #ifndef GENSVM_PRINT_ITER
37 : #define GENSVM_PRINT_ITER 100
38 : #endif
39 :
40 : /**
41 : * @brief The main training loop for GenSVM
42 : *
43 : * @details
44 : * This function is the main training function. This function
45 : * handles the optimization of the model with the given model parameters, with
46 : * the data given. On return the matrix GenModel::V contains the optimal
47 : * weight matrix.
48 : *
49 : * In this function, step doubling is used in the majorization algorithm after
50 : * a burn-in of 50 iterations.
51 : *
52 : * @param[in,out] model the GenModel to be trained. Contains optimal
53 : * V on exit.
54 : * @param[in] data the GenData to train the model with.
55 : */
56 3 : void gensvm_optimize(struct GenModel *model, struct GenData *data)
57 : {
58 3 : long it = 0;
59 : double L, Lbar;
60 :
61 3 : long n = model->n;
62 3 : long m = model->m;
63 3 : long K = model->K;
64 :
65 : // initialize the workspace
66 3 : struct GenWork *work = gensvm_init_work(model);
67 :
68 : // print some info on the dataset and model configuration
69 3 : note("Starting main loop.\n");
70 3 : note("Dataset:\n");
71 3 : note("\tn = %i\n", n);
72 3 : note("\tm = %i\n", m);
73 3 : note("\tK = %i\n", K);
74 3 : note("Parameters:\n");
75 3 : note("\tkappa = %f\n", model->kappa);
76 3 : note("\tp = %f\n", model->p);
77 3 : note("\tlambda = %15.16f\n", model->lambda);
78 3 : note("\tepsilon = %g\n", model->epsilon);
79 3 : note("\n");
80 :
81 : // compute necessary simplex vectors
82 3 : gensvm_simplex(model);
83 3 : gensvm_simplex_diff(model);
84 :
85 : // get initial loss
86 3 : L = gensvm_get_loss(model, data, work);
87 3 : Lbar = L + 2.0*model->epsilon*L;
88 :
89 : // run main loop
90 470 : while ((it < model->max_iter) && (Lbar - L)/L > model->epsilon)
91 : {
92 : // ensures V contains newest V and Vbar contains V from
93 : // previous
94 464 : gensvm_get_update(model, data, work);
95 464 : if (it > 50)
96 311 : gensvm_step_doubling(model);
97 :
98 464 : Lbar = L;
99 464 : L = gensvm_get_loss(model, data, work);
100 :
101 464 : if (it % GENSVM_PRINT_ITER == 0)
102 6 : note("iter = %li, L = %15.16f, Lbar = %15.16f, "
103 6 : "reldiff = %15.16f\n", it, L, Lbar, (Lbar - L)/L);
104 464 : it++;
105 : }
106 :
107 : // status == 0 means training was successful
108 3 : model->status = 0;
109 :
110 : // print warnings if necessary
111 3 : if (L > Lbar) {
112 0 : err("[GenSVM Warning]: Negative step occurred in "
113 : "majorization.\n");
114 0 : model->status = 1;
115 : }
116 :
117 3 : if (it >= model->max_iter) {
118 0 : err("[GenSVM Warning]: maximum number of iterations "
119 : "reached.\n");
120 0 : model->status = 2;
121 : }
122 :
123 : // print final iteration count and loss
124 3 : note("Optimization finished, iter = %li, loss = %15.16f, "
125 : "rel. diff. = %15.16f\n", it-1, L,
126 3 : (Lbar - L)/L);
127 :
128 : // compute and print the number of SVs in the model
129 3 : note("Number of support vectors: %li\n", gensvm_num_sv(model));
130 :
131 : // store the training error in the model
132 3 : model->training_error = (Lbar - L)/L;
133 :
134 : // store the iteration count in the model
135 3 : model->elapsed_iter = it - 1;
136 :
137 : // free the workspace
138 3 : gensvm_free_work(work);
139 3 : }
140 :
141 : /**
142 : * @brief Calculate the current value of the loss function
143 : *
144 : * @details
145 : * The current loss function value is calculated based on the matrix V in the
146 : * given model. Note that the matrix ZV is passed explicitly to avoid having
147 : * to reallocate memory at every step.
148 : *
149 : * @param[in] model GenModel structure which holds the current
150 : * estimate V
151 : * @param[in] data GenData structure
152 : * @param[in] work allocated workspace with the ZV matrix to use
153 : * @returns the current value of the loss function
154 : */
155 469 : double gensvm_get_loss(struct GenModel *model, struct GenData *data,
156 : struct GenWork *work)
157 : {
158 : long i, j;
159 469 : long n = model->n;
160 469 : long K = model->K;
161 469 : long m = model->m;
162 :
163 469 : double value, rowvalue, loss = 0.0;
164 :
165 469 : gensvm_calculate_errors(model, data, work->ZV);
166 469 : gensvm_calculate_huber(model);
167 :
168 4891 : for (i=0; i<n; i++) {
169 4422 : rowvalue = 0;
170 4422 : value = 0;
171 22094 : for (j=0; j<K; j++) {
172 17672 : if (j == (data->y[i]-1))
173 4422 : continue;
174 13250 : value = matrix_get(model->H, K, i, j);
175 13250 : value = pow(value, model->p);
176 13250 : rowvalue += value;
177 : }
178 4422 : rowvalue = pow(rowvalue, 1.0/(model->p));
179 4422 : rowvalue *= model->rho[i];
180 4422 : loss += rowvalue;
181 : }
182 469 : loss /= ((double) n);
183 :
184 469 : value = 0;
185 2317 : for (i=1; i<m+1; i++) {
186 7386 : for (j=0; j<K-1; j++) {
187 5538 : value += pow(matrix_get(model->V, K-1, i, j), 2.0);
188 : }
189 : }
190 469 : loss += model->lambda * value;
191 :
192 469 : return loss;
193 : }
194 :
195 : /**
196 : * @brief Use step doubling
197 : *
198 : * @details
199 : * Step doubling can be used to speed up the maorization algorithm. Instead of
200 : * using the value at the minimimum of the majorization function, the value
201 : * ``opposite'' the majorization point is used. This can essentially cut the
202 : * number of iterations necessary to reach the minimum in half.
203 : *
204 : * @param[in] model GenModel containing the augmented parameters
205 : */
206 312 : void gensvm_step_doubling(struct GenModel *model)
207 : {
208 : long i, j;
209 : double value;
210 :
211 312 : long m = model->m;
212 312 : long K = model->K;
213 :
214 1845 : for (i=0; i<m+1; i++) {
215 6128 : for (j=0; j<K-1; j++) {
216 4595 : matrix_mul(model->V, K-1, i, j, 2.0);
217 4595 : value = - matrix_get(model->Vbar, K-1, i, j);
218 4595 : matrix_add(model->V, K-1, i, j, value);
219 : }
220 : }
221 312 : }
222 :
223 : /**
224 : * @brief Calculate the Huber hinge errors
225 : *
226 : * @details
227 : * For each of the scalar errors in Q the Huber hinge errors are
228 : * calculated. The Huber hinge is here defined as
229 : * @f[
230 : * h(q) =
231 : * \begin{dcases}
232 : * 1 - q - \frac{\kappa + 1}{2} & \text{if } q \leq
233 : * -\kappa \\
234 : * \frac{1}{2(\kappa + 1)} ( 1 - q)^2 & \text{if } q \in
235 : * (-\kappa, 1] \\
236 : * 0 & \text{if } q > 1
237 : * \end{dcases}
238 : * @f]
239 : *
240 : * @param[in,out] model the corresponding GenModel
241 : */
242 472 : void gensvm_calculate_huber(struct GenModel *model)
243 : {
244 : long i, j;
245 : double q, value;
246 :
247 4915 : for (i=0; i<model->n; i++) {
248 22178 : for (j=0; j<model->K; j++) {
249 17735 : q = matrix_get(model->Q, model->K, i, j);
250 17735 : value = 0.0;
251 17735 : if (q <= -model->kappa) {
252 59 : value = 1.0 - q - (model->kappa+1.0)/2.0;
253 17676 : } else if (q <= 1.0) {
254 15363 : value = 1.0/(2.0*model->kappa+2.0)*pow(1.0 - q,
255 : 2.0);
256 : }
257 17735 : matrix_set(model->H, model->K, i, j, value);
258 : }
259 : }
260 472 : }
261 :
262 : /**
263 : * @brief Calculate the scalar errors
264 : *
265 : * @details
266 : * Calculate the scalar errors q based on the current estimate of V, and
267 : * store these in Q. It is assumed that the memory for Q has already been
268 : * allocated. In addition, the matrix ZV is calculated here. It is assigned
269 : * to a pre-allocated block of memory, which is passed to this function.
270 : *
271 : * @param[in,out] model the corresponding GenModel
272 : * @param[in] data the corresponding GenData
273 : * @param[in,out] ZV a pointer to a memory block for ZV. On exit
274 : * this block is updated with the new ZV matrix
275 : * calculated with GenModel::V
276 : */
277 472 : void gensvm_calculate_errors(struct GenModel *model, struct GenData *data,
278 : double *ZV)
279 : {
280 : long i, j;
281 472 : double q, *uu_row = NULL;
282 :
283 472 : long n = model->n;
284 472 : long K = model->K;
285 :
286 472 : gensvm_calculate_ZV(model, data, ZV);
287 :
288 4918 : for (i=0; i<n; i++) {
289 22190 : for (j=0; j<K; j++) {
290 17744 : if (j == (data->y[i]-1))
291 4446 : continue;
292 13298 : uu_row = &model->UU[((data->y[i]-1)*K+j)*(K-1)];
293 13298 : q = cblas_ddot(K-1, &ZV[i*(K-1)], 1, uu_row, 1);
294 13298 : matrix_set(model->Q, K, i, j, q);
295 : }
296 : }
297 472 : }
298 :
|
