93 lines
2.0 KiB
C
93 lines
2.0 KiB
C
#include "common.h"
|
|
#include "percep.h"
|
|
#include <math.h>
|
|
#include <gsl/gsl_matrix.h>
|
|
#include <gsl/gsl_blas.h>
|
|
|
|
/* `iterate(data, w, b, d, r)` performs one
|
|
* iteration of the perceptron training.
|
|
*
|
|
* For each point xi compute:
|
|
*
|
|
* 1. yi = θ((w, xi) + b)
|
|
*
|
|
* 2. Δi = r⋅(d - yi)
|
|
*
|
|
* 3. w = w + Δi⋅xi, b = b + Δi
|
|
*
|
|
* The results are computed in-place.
|
|
*/
|
|
void iterate(
|
|
gsl_matrix *data,
|
|
gsl_vector *weight,
|
|
double *bias,
|
|
int expected,
|
|
double rate) {
|
|
|
|
double proj, delta;
|
|
gsl_vector *x = gsl_vector_alloc(2);
|
|
|
|
for (size_t i = 0; i < data->size1; i++) {
|
|
/* Get a vector view of the
|
|
* current row in the data matrix.
|
|
*/
|
|
gsl_vector row = gsl_matrix_const_row(data, i).vector;
|
|
gsl_vector_memcpy(x, &row);
|
|
|
|
/* Project x onto the weight vector. */
|
|
gsl_blas_ddot(weight, x, &proj);
|
|
|
|
/* Calculate Δ
|
|
* Note: the step function θ(x) is computed
|
|
* by negating the sign bit of the floating
|
|
* point.
|
|
*/
|
|
delta = rate * (expected - !signbit(proj + *bias));
|
|
|
|
/* Update weight and bias. */
|
|
*bias += delta;
|
|
gsl_vector_scale(x, delta);
|
|
gsl_vector_add(weight, x);
|
|
}
|
|
|
|
gsl_vector_free(x);
|
|
}
|
|
|
|
/* `percep_train(sig, noise, iter, &cut)`
|
|
* train a single perceptron to discriminate `sig`
|
|
* from `noise`.
|
|
*
|
|
* The weigths are adjusted `iter` times and
|
|
* returned, the bias/threshold is stored in the
|
|
* `cut` argument.
|
|
*/
|
|
gsl_vector *percep_train(
|
|
sample_t *signal, sample_t *noise,
|
|
int iter, double *cut) {
|
|
|
|
/* Initially set weights/bias to zero
|
|
* and a high learning rate.
|
|
*/
|
|
gsl_vector *w = gsl_vector_calloc(2);
|
|
double bias = 0;
|
|
double rate = 0.8;
|
|
|
|
/* Go trough the sample `iter` times
|
|
* and recalculate the weights.
|
|
*/
|
|
for (int i = 0; i < iter; i++) {
|
|
iterate( noise->data, w, &bias, 0, rate);
|
|
iterate(signal->data, w, &bias, 1, rate);
|
|
}
|
|
|
|
/* Normalise the weights and force
|
|
* positiveness for easier
|
|
* comparison with other methods.
|
|
*/
|
|
double norm = gsl_blas_dnrm2(w);
|
|
gsl_vector_scale(w, 1/norm);
|
|
*cut = -bias/norm;
|
|
|
|
return w;
|
|
}
|