mirror of https://github.com/davisking/dlib.git
356 lines
15 KiB
C++
356 lines
15 KiB
C++
// The contents of this file are in the public domain. See LICENSE_FOR_EXAMPLE_PROGRAMS.txt
|
|
/*
|
|
|
|
This is an example illustrating the use of the empirical_kernel_map
|
|
from the dlib C++ Library.
|
|
|
|
This example program assumes you are familiar with some general elements of
|
|
the library. In particular, you should have at least read the svm_ex.cpp
|
|
and matrix_ex.cpp examples.
|
|
|
|
|
|
Most of the machine learning algorithms in dlib are some flavor of "kernel machine".
|
|
This means they are all simple linear algorithms that have been formulated such
|
|
that the only way they look at the data given by a user is via dot products between
|
|
the data samples. These algorithms are made more useful via the application of the
|
|
so-called kernel trick. This trick is to replace the dot product with a user
|
|
supplied function which takes two samples and returns a real number. This function
|
|
is the kernel that is required by so many algorithms. The most basic kernel is the
|
|
linear_kernel which is simply a normal dot product. More interesting, however,
|
|
are kernels which first apply some nonlinear transformation to the user's data samples
|
|
and then compute a dot product. In this way, a simple algorithm that finds a linear
|
|
plane to separate data (e.g. the SVM algorithm) can be made to solve complex
|
|
nonlinear learning problems.
|
|
|
|
An important element of the kernel trick is that these kernel functions perform
|
|
the nonlinear transformation implicitly. That is, if you look at the implementations
|
|
of these kernel functions you won't see code that transforms two input vectors in
|
|
some way and then computes their dot products. Instead you will see a simple function
|
|
that takes two input vectors and just computes a single real number via some simple
|
|
process. You can basically think of this as an optimization. Imagine that originally
|
|
we wrote out the entire procedure to perform the nonlinear transformation and then
|
|
compute the dot product but then noticed we could cancel a few terms here and there
|
|
and simplify the whole thing down into a more compact and easily evaluated form.
|
|
The result is a nice function that computes what we want but we no longer get to see
|
|
what those nonlinearly transformed input vectors are.
|
|
|
|
The empirical_kernel_map is a tool that undoes this. It allows you to obtain these
|
|
nonlinearly transformed vectors. It does this by taking a set of data samples from
|
|
the user (referred to as basis samples), applying the nonlinear transformation to all
|
|
of them, and then constructing a set of orthonormal basis vectors which spans the space
|
|
occupied by those transformed input samples. Then if we wish to obtain the nonlinear
|
|
version of any data sample we can simply project it onto this orthonormal basis and
|
|
we obtain a regular vector of real numbers which represents the nonlinearly transformed
|
|
version of the data sample. The empirical_kernel_map has been formulated to use only
|
|
dot products between data samples so it is capable of performing this service for any
|
|
user supplied kernel function.
|
|
|
|
The empirical_kernel_map is useful because it is often difficult to formulate an
|
|
algorithm in a way that uses only dot products. So the empirical_kernel_map lets
|
|
us easily kernelize any algorithm we like by using this object during a preprocessing
|
|
step. However, it should be noted that the algorithm is only practical when used
|
|
with at most a few thousand basis samples. Fortunately, most datasets live in
|
|
subspaces that are relatively low dimensional. So for these datasets, using the
|
|
empirical_kernel_map is practical assuming an appropriate set of basis samples can be
|
|
selected by the user. To help with this dlib supplies the linearly_independent_subset_finder.
|
|
I also often find that just picking a random subset of the data as a basis works well.
|
|
|
|
|
|
|
|
In what follows, we walk through the process of creating an empirical_kernel_map,
|
|
projecting data to obtain the nonlinearly transformed vectors, and then doing a
|
|
few interesting things with the data.
|
|
*/
|
|
|
|
|
|
|
|
|
|
#include <dlib/svm.h>
|
|
#include <dlib/rand.h>
|
|
#include <iostream>
|
|
#include <vector>
|
|
|
|
|
|
using namespace std;
|
|
using namespace dlib;
|
|
|
|
// ----------------------------------------------------------------------------------------
|
|
|
|
// First let's make a typedef for the kind of samples we will be using.
|
|
typedef matrix<double, 0, 1> sample_type;
|
|
|
|
// We will be using the radial_basis_kernel in this example program.
|
|
typedef radial_basis_kernel<sample_type> kernel_type;
|
|
|
|
// ----------------------------------------------------------------------------------------
|
|
|
|
void generate_concentric_circles (
|
|
std::vector<sample_type>& samples,
|
|
std::vector<double>& labels,
|
|
const int num_points
|
|
);
|
|
/*!
|
|
requires
|
|
- num_points > 0
|
|
ensures
|
|
- generates two circles centered at the point (0,0), one of radius 1 and
|
|
the other of radius 5. These points are stored into samples. labels will
|
|
tell you if a given samples is from the smaller circle (its label will be 1)
|
|
or from the larger circle (its label will be 2).
|
|
- each circle will be made up of num_points
|
|
!*/
|
|
|
|
// ----------------------------------------------------------------------------------------
|
|
|
|
void test_empirical_kernel_map (
|
|
const std::vector<sample_type>& samples,
|
|
const std::vector<double>& labels,
|
|
const empirical_kernel_map<kernel_type>& ekm
|
|
);
|
|
/*!
|
|
This function computes various interesting things with the empirical_kernel_map.
|
|
See its implementation below for details.
|
|
!*/
|
|
|
|
// ----------------------------------------------------------------------------------------
|
|
|
|
int main()
|
|
{
|
|
std::vector<sample_type> samples;
|
|
std::vector<double> labels;
|
|
|
|
// Declare an instance of the kernel we will be using.
|
|
const kernel_type kern(0.1);
|
|
|
|
// create a dataset with two concentric circles. There will be 100 points on each circle.
|
|
generate_concentric_circles(samples, labels, 100);
|
|
|
|
empirical_kernel_map<kernel_type> ekm;
|
|
|
|
|
|
// Here we create an empirical_kernel_map using all of our data samples as basis samples.
|
|
cout << "\n\nBuilding an empirical_kernel_map with " << samples.size() << " basis samples." << endl;
|
|
ekm.load(kern, samples);
|
|
cout << "Test the empirical_kernel_map when loaded with every sample." << endl;
|
|
test_empirical_kernel_map(samples, labels, ekm);
|
|
|
|
|
|
|
|
|
|
|
|
|
|
// create a new dataset with two concentric circles. There will be 1000 points on each circle.
|
|
generate_concentric_circles(samples, labels, 1000);
|
|
// Rather than use all 2000 samples as basis samples we are going to use the
|
|
// linearly_independent_subset_finder to pick out a good basis set. The idea behind this
|
|
// object is to try and find the 40 or so samples that best spans the subspace containing all the
|
|
// data.
|
|
linearly_independent_subset_finder<kernel_type> lisf(kern, 40);
|
|
// populate lisf with samples. We have configured it to allow at most 40 samples but this function
|
|
// may determine that fewer samples are necessary to form a good basis. In this example program
|
|
// it will select only 26.
|
|
fill_lisf(lisf, samples);
|
|
|
|
// Now reload the empirical_kernel_map but this time using only our small basis
|
|
// selected using the linearly_independent_subset_finder.
|
|
cout << "\n\nBuilding an empirical_kernel_map with " << lisf.size() << " basis samples." << endl;
|
|
ekm.load(lisf);
|
|
cout << "Test the empirical_kernel_map when loaded with samples from the lisf object." << endl;
|
|
test_empirical_kernel_map(samples, labels, ekm);
|
|
|
|
|
|
cout << endl;
|
|
}
|
|
|
|
// ----------------------------------------------------------------------------------------
|
|
|
|
void test_empirical_kernel_map (
|
|
const std::vector<sample_type>& samples,
|
|
const std::vector<double>& labels,
|
|
const empirical_kernel_map<kernel_type>& ekm
|
|
)
|
|
{
|
|
|
|
std::vector<sample_type> projected_samples;
|
|
|
|
// The first thing we do is compute the nonlinearly projected vectors using the
|
|
// empirical_kernel_map.
|
|
for (unsigned long i = 0; i < samples.size(); ++i)
|
|
{
|
|
projected_samples.push_back(ekm.project(samples[i]));
|
|
}
|
|
|
|
// Note that a kernel matrix is just a matrix M such that M(i,j) == kernel(samples[i],samples[j]).
|
|
// So below we are computing the normal kernel matrix as given by the radial_basis_kernel and the
|
|
// input samples. We also compute the kernel matrix for all the projected_samples as given by the
|
|
// linear_kernel. Note that the linear_kernel just computes normal dot products. So what we want to
|
|
// see is that the dot products between all the projected_samples samples are the same as the outputs
|
|
// of the kernel function for their respective untransformed input samples. If they match then
|
|
// we know that the empirical_kernel_map is working properly.
|
|
const matrix<double> normal_kernel_matrix = kernel_matrix(ekm.get_kernel(), samples);
|
|
const matrix<double> new_kernel_matrix = kernel_matrix(linear_kernel<sample_type>(), projected_samples);
|
|
|
|
cout << "Max kernel matrix error: " << max(abs(normal_kernel_matrix - new_kernel_matrix)) << endl;
|
|
cout << "Mean kernel matrix error: " << mean(abs(normal_kernel_matrix - new_kernel_matrix)) << endl;
|
|
/*
|
|
Example outputs from these cout statements.
|
|
For the case where we use all samples as basis samples:
|
|
Max kernel matrix error: 7.32747e-15
|
|
Mean kernel matrix error: 7.47789e-16
|
|
|
|
For the case where we use only 26 samples as basis samples:
|
|
Max kernel matrix error: 0.000953573
|
|
Mean kernel matrix error: 2.26008e-05
|
|
|
|
|
|
Note that if we use enough basis samples we can perfectly span the space of input samples.
|
|
In that case we get errors that are essentially just rounding noise (Moreover, using all the
|
|
samples is always enough since they are always within their own span). Once we start
|
|
to use fewer basis samples we may begin to get approximation error. In the second case we
|
|
used 26 and we can see that the data doesn't really lay exactly in a 26 dimensional subspace.
|
|
But it is pretty close.
|
|
*/
|
|
|
|
|
|
|
|
// Now let's do something more interesting. The following loop finds the centroids
|
|
// of the two classes of data.
|
|
sample_type class1_center;
|
|
sample_type class2_center;
|
|
for (unsigned long i = 0; i < projected_samples.size(); ++i)
|
|
{
|
|
if (labels[i] == 1)
|
|
class1_center += projected_samples[i];
|
|
else
|
|
class2_center += projected_samples[i];
|
|
}
|
|
|
|
const int points_per_class = samples.size()/2;
|
|
class1_center /= points_per_class;
|
|
class2_center /= points_per_class;
|
|
|
|
|
|
// Now classify points by which center they are nearest. Recall that the data
|
|
// is made up of two concentric circles. Normally you can't separate two concentric
|
|
// circles by checking which points are nearest to each center since they have the same
|
|
// centers. However, the kernel trick makes the data separable and the loop below will
|
|
// perfectly classify each data point.
|
|
for (unsigned long i = 0; i < projected_samples.size(); ++i)
|
|
{
|
|
double distance_to_class1 = length(projected_samples[i] - class1_center);
|
|
double distance_to_class2 = length(projected_samples[i] - class2_center);
|
|
|
|
bool predicted_as_class_1 = (distance_to_class1 < distance_to_class2);
|
|
|
|
// Now print a message for any misclassified points.
|
|
if (predicted_as_class_1 == true && labels[i] != 1)
|
|
cout << "A point was misclassified" << endl;
|
|
|
|
if (predicted_as_class_1 == false && labels[i] != 2)
|
|
cout << "A point was misclassified" << endl;
|
|
}
|
|
|
|
|
|
|
|
// Next, note that classifying a point based on its distance between two other
|
|
// points is the same thing as using the plane that lies between those two points
|
|
// as a decision boundary. So let's compute that decision plane and use it to classify
|
|
// all the points.
|
|
|
|
sample_type plane_normal_vector = class1_center - class2_center;
|
|
// The point right in the center of our two classes should be on the deciding plane, not
|
|
// on one side or the other. This consideration brings us to the formula for the bias.
|
|
double bias = dot((class1_center+class2_center)/2, plane_normal_vector);
|
|
|
|
// Now classify points by which side of the plane they are on.
|
|
for (unsigned long i = 0; i < projected_samples.size(); ++i)
|
|
{
|
|
double side = dot(plane_normal_vector, projected_samples[i]) - bias;
|
|
|
|
bool predicted_as_class_1 = (side > 0);
|
|
|
|
// Now print a message for any misclassified points.
|
|
if (predicted_as_class_1 == true && labels[i] != 1)
|
|
cout << "A point was misclassified" << endl;
|
|
|
|
if (predicted_as_class_1 == false && labels[i] != 2)
|
|
cout << "A point was misclassified" << endl;
|
|
}
|
|
|
|
|
|
// It would be nice to convert this decision rule into a normal decision_function object and
|
|
// dispense with the empirical_kernel_map. Happily, it is possible to do so. Consider the
|
|
// following example code:
|
|
decision_function<kernel_type> dec_funct = ekm.convert_to_decision_function(plane_normal_vector);
|
|
// The dec_funct now computes dot products between plane_normal_vector and the projection
|
|
// of any sample point given to it. All that remains is to account for the bias.
|
|
dec_funct.b = bias;
|
|
|
|
// now classify points by which side of the plane they are on.
|
|
for (unsigned long i = 0; i < samples.size(); ++i)
|
|
{
|
|
double side = dec_funct(samples[i]);
|
|
|
|
// And let's just check that the dec_funct really does compute the same thing as the previous equation.
|
|
double side_alternate_equation = dot(plane_normal_vector, projected_samples[i]) - bias;
|
|
if (abs(side-side_alternate_equation) > 1e-14)
|
|
cout << "dec_funct error: " << abs(side-side_alternate_equation) << endl;
|
|
|
|
bool predicted_as_class_1 = (side > 0);
|
|
|
|
// Now print a message for any misclassified points.
|
|
if (predicted_as_class_1 == true && labels[i] != 1)
|
|
cout << "A point was misclassified" << endl;
|
|
|
|
if (predicted_as_class_1 == false && labels[i] != 2)
|
|
cout << "A point was misclassified" << endl;
|
|
}
|
|
|
|
}
|
|
|
|
// ----------------------------------------------------------------------------------------
|
|
|
|
void generate_concentric_circles (
|
|
std::vector<sample_type>& samples,
|
|
std::vector<double>& labels,
|
|
const int num
|
|
)
|
|
{
|
|
sample_type m(2,1);
|
|
samples.clear();
|
|
labels.clear();
|
|
|
|
dlib::rand rnd;
|
|
|
|
// make some samples near the origin
|
|
double radius = 1.0;
|
|
for (long i = 0; i < num; ++i)
|
|
{
|
|
double sign = 1;
|
|
if (rnd.get_random_double() < 0.5)
|
|
sign = -1;
|
|
m(0) = 2*radius*rnd.get_random_double()-radius;
|
|
m(1) = sign*sqrt(radius*radius - m(0)*m(0));
|
|
|
|
samples.push_back(m);
|
|
labels.push_back(1);
|
|
}
|
|
|
|
// make some samples in a circle around the origin but far away
|
|
radius = 5.0;
|
|
for (long i = 0; i < num; ++i)
|
|
{
|
|
double sign = 1;
|
|
if (rnd.get_random_double() < 0.5)
|
|
sign = -1;
|
|
m(0) = 2*radius*rnd.get_random_double()-radius;
|
|
m(1) = sign*sqrt(radius*radius - m(0)*m(0));
|
|
|
|
samples.push_back(m);
|
|
labels.push_back(2);
|
|
}
|
|
}
|
|
|
|
// ----------------------------------------------------------------------------------------
|
|
|