// 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; // ----------------------------------------------------------------------------------------voidgenerate_concentric_circles( std::vector<sample_type>& samples, std::vector<double>& labels, constintnum_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 !*/ // ----------------------------------------------------------------------------------------voidtest_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. !*/ // ----------------------------------------------------------------------------------------intmain(){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;}// ----------------------------------------------------------------------------------------voidtest_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 (unsignedlongi = 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 (unsignedlongi = 0; i < projected_samples.size(); ++i){if (labels[i] == 1) class1_center += projected_samples[i]; else class2_center += projected_samples[i];}constintpoints_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 (unsignedlongi = 0; i < projected_samples.size(); ++i){doubledistance_to_class1 = length(projected_samples[i] - class1_center);doubledistance_to_class2 = length(projected_samples[i] - class2_center);boolpredicted_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.doublebias = dot((class1_center+class2_center)/2, plane_normal_vector); // Now classify points by which side of the plane they are on. for (unsignedlongi = 0; i < projected_samples.size(); ++i){doubleside = dot(plane_normal_vector, projected_samples[i]) - bias;boolpredicted_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 (unsignedlongi = 0; i < samples.size(); ++i){doubleside = dec_funct(samples[i]); // And let's just check that the dec_funct really does compute the same thing as the previous equation.doubleside_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;boolpredicted_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;}}// ----------------------------------------------------------------------------------------voidgenerate_concentric_circles( std::vector<sample_type>& samples, std::vector<double>& labels, constintnum ){sample_type m(2,1); samples.clear(); labels.clear(); dlib::rand rnd; // make some samples near the origindoubleradius = 1.0; for (longi = 0; i < num; ++i){doublesign = 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 (longi = 0; i < num; ++i){doublesign = 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);}}// ----------------------------------------------------------------------------------------