apaam/mlpack
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

Browse Source
This commit is contained in:
houyang
2008-01-18 17:46:31 +00:00
parent 6ebd4223e4
commit 24fedb504c
4 changed files with 267 additions and 693 deletions
Download Patch File Download Diff File Expand all files Collapse all files
fastlib/u/garryb/svm/smo.h
+90 -106
View File
@@ -4,14 +4,13 @@
#include "fastlib/fastlib.h"
/* TODO: I don't actually want these to be public */
/* but sometimes we should provide freedoms for our advanced users */
const double SMO_ZERO = 1.0e-8;
const double SMO_EPS = 1.0e-4;
const double SMO_TOLERANCE = 1.0e-4;
template<typename TKernel>
class SMO {
FORBID_ACCIDENTAL_COPIES(SMO);
FORBID_COPY(SMO);
public:
typedef TKernel Kernel;
@@ -20,11 +19,9 @@ class SMO {
Matrix kernel_cache_sign_;
Kernel kernel_;
const Dataset *dataset_;
index_t n_data_; // number of data samples
index_t n_classes_; // number of classes
Matrix matrix_; // alias for the data matrix
Vector alpha_; // the alphas, to be optimized
index_t n_sv_; // number of support vectors
index_t n_data_;
Matrix matrix_;
Vector alpha_;
Vector error_;
double thresh_;
double c_;
@@ -40,16 +37,26 @@ class SMO {
*
* You must initialize separately the kernel.
*/
void Init(int n_classes_in, double c_in, int budget_in) {
void Init(const Dataset* dataset_in, double c_in, int budget_in) {
c_ = c_in;
n_classes_ = n_classes_in;
budget_ = budget_in;
thresh_ = 0.0;
n_sv_ = 0;
dataset_ = dataset_in;
matrix_.Alias(dataset_->matrix());
n_data_ = matrix_.n_cols();
budget_ = min(budget_in, n_data_);
alpha_.Init(n_data_);
alpha_.SetZero();
sum_alpha_ = 0;
error_.Init(n_data_);
error_.SetZero();
thresh_ = 0;
}
void Train(const Dataset* dataset_in);
void Train();
const Kernel& kernel() const {
return kernel_;
@@ -63,15 +70,11 @@ class SMO {
return thresh_;
}
//index_t num_sv() const {
// return n_sv_;
//}
void GetSVM(ArrayList<index_t> &dataset_index, ArrayList<double> &coef, ArrayList<bool> &sv_indicator);
void GetSVM(Matrix *support_vectors, Vector *support_alpha) const;
private:
index_t TrainIteration_(bool examine_all);
bool TryChange_(index_t j);
bool TakeStep_(index_t i, index_t j, double error_j);
@@ -89,7 +92,6 @@ class SMO {
return alpha <= 0 || alpha >= c_;
}
// labels: the last row of the data matrix, 0 or 1
int GetLabelSign_(index_t i) const {
return matrix_.get(matrix_.n_rows()-1, i) != 0 ? 1 : -1;
}
@@ -102,7 +104,7 @@ class SMO {
double val;
if (!IsBound_(alpha_[i])) {
val = error_[i];
VERBOSE_MSG(0, "error values %f and %f", error_[i], Evaluate_(i) - GetLabelSign_(i));
DEBUG_MSG(0, "error values %f and %f", error_[i], Evaluate_(i) - GetLabelSign_(i));
} else {
val = CalculateError_(i);
}
@@ -129,6 +131,7 @@ class SMO {
Vector v_j;
GetVector_(j, &v_j);
double k = kernel_.Eval(v_i, v_j);
kernel_cache_sign_.set(j, i, k * GetLabelSign_(i) * GetLabelSign_(j));
}
}
@@ -136,50 +139,72 @@ class SMO {
}
};
// Budget SMO training for 2-classes
template<typename TKernel>
void SMO<TKernel>::Train(const Dataset* dataset_in) {
void SMO<TKernel>::GetSVM(Matrix *support_vectors, Vector *support_alpha) const {
index_t n_support = 0;
index_t i_support = 0;
for (index_t i = 0; i < n_data_; i++) {
if (unlikely(alpha_[i] != 0)) {
n_support++;
}
}
support_vectors->Init(matrix_.n_rows() - 1, n_support);
support_alpha->Init(n_support);
for (index_t i = 0; i < n_data_; i++) {
if (unlikely(alpha_[i] != 0)) {
Vector source;
Vector dest;
GetVector_(i, &source);
support_vectors->MakeColumnVector(i_support, &dest);
dest.CopyValues(source);
(*support_alpha)[i_support] = alpha_[i] * GetLabelSign_(i);
i_support++;
}
}
}
template<typename TKernel>
double SMO<TKernel>::Evaluate_(index_t i) const {
// TODO: This only handles linear
Vector kernel_values;
double summation = 0;
kernel_cache_sign_.MakeColumnVector(i, &kernel_values);
summation = la::Dot(alpha_, kernel_values) * GetLabelSign_(i);
return (summation - thresh_);
}
template<typename TKernel>
void SMO<TKernel>::Train() {
bool examine_all = true;
index_t num_changed = 0;
int n_iter = 0;
// data-dependent initialization
dataset_ = dataset_in;
matrix_.Alias(dataset_->matrix());
n_data_ = matrix_.n_cols();
budget_ = min(budget_, n_data_);
alpha_.Init(n_data_);
alpha_.SetZero();
sum_alpha_ = 0;
error_.Init(n_data_);
error_.SetZero();
// calculate kernel_cache_sign_: [k_ij* y_i* y_j]
CalcKernels_();
while (num_changed > 0 || examine_all) { // TODO: other stopping criteria
VERBOSE_GOT_HERE(0);
// SMO iterations
while ((num_changed > 0 || examine_all)) {
DEBUG_GOT_HERE(0);
num_changed = TrainIteration_(examine_all);
if (examine_all) {
examine_all = false;
} else if (num_changed == 0) {
} else if (num_changed == 0) {
examine_all = true;
}
// if exceed the maximum number of iterations, finished
// 200...TODO
if (++n_iter == 200) {
fprintf(stderr, "Max iterations %f!!!!!!!!!!!!!!!!!!!!!!!!!!\n",
sum_alpha_);
break;
}
// for every 100 iterations, do budget shrink
if (n_iter % 100 == 0 && budget_ < n_data_) {
MinHeap<double, int> alphas;
@@ -202,16 +227,14 @@ void SMO<TKernel>::Train(const Dataset* dataset_in) {
if (!IsBound_(alpha_[i])) {
error_[i] = CalculateError_(i);
} else {
error_[i] = 0;
error_[i] = 0;
}
sum_alpha_ += alpha_[i];
}
} // if
} // while
}
}
}
// SMO training iterations
template<typename TKernel>
index_t SMO<TKernel>::TrainIteration_(bool examine_all) {
index_t num_changed = 0;
@@ -225,16 +248,15 @@ index_t SMO<TKernel>::TrainIteration_(bool examine_all) {
return num_changed;
}
// try to find working set (maximal violating pair)
template<typename TKernel>
bool SMO<TKernel>::TryChange_(index_t j) {
double error_j = Error_(j); // -y_j g_j^* -thresh
double error_j = Error_(j);
double rj = error_j * GetLabelSign_(j);
VERBOSE_GOT_HERE(0);
DEBUG_GOT_HERE(0);
if (!( (rj < -SMO_TOLERANCE && alpha_[j] < c_)
||(rj > SMO_TOLERANCE && alpha_[j] > 0) )) {
if (!((rj < -SMO_TOLERANCE && alpha_[j] < c_)
|| (rj > SMO_TOLERANCE && alpha_[j] > 0))) {
return false; // nothing to change
}
@@ -245,11 +267,10 @@ bool SMO<TKernel>::TryChange_(index_t j) {
double diff_max = 0;
// find the max(abs(y_i g_i^*))
for (index_t k = 0; k < n_data_; k++) {
if (!IsBound_(alpha_[k])) { // if 0 < alpha_[k] < c_
if (!IsBound_(alpha_[k])) {
double error_k = error_[k];
double diff_k = fabs(error_k - error_j); // abs(y_j g_j^* - y_k g_k^*)
double diff_k = fabs(error_k - error_j);
if (unlikely(diff_k > diff_max)) {
diff_max = diff_k;
i = k;
@@ -261,7 +282,7 @@ bool SMO<TKernel>::TryChange_(index_t j) {
}
}
VERBOSE_GOT_HERE(0);
DEBUG_GOT_HERE(0);
// try searching through non-bound examples
index_t start_i = rand() % n_data_;
index_t i = start_i;
@@ -274,7 +295,7 @@ bool SMO<TKernel>::TryChange_(index_t j) {
i = (i + 1) % n_data_;
} while (i != start_i);
VERBOSE_GOT_HERE(0);
DEBUG_GOT_HERE(0);
// try searching through all examples
start_i = rand() % n_data_;
i = start_i;
@@ -289,11 +310,10 @@ bool SMO<TKernel>::TryChange_(index_t j) {
return false;
}
// search direction, update gradient
template<typename TKernel>
bool SMO<TKernel>::TakeStep_(index_t i, index_t j, double error_j) {
if (i == j) {
VERBOSE_GOT_HERE(0);
DEBUG_GOT_HERE(0);
return false;
}
@@ -325,8 +345,8 @@ bool SMO<TKernel>::TakeStep_(index_t i, index_t j, double error_j) {
if (l >= u - SMO_TOLERANCE) {
// TODO: might put in some tolerance
VERBOSE_MSG(0, "l=%f, u=%f, r=%f, c_=%f, s=%f", l, u, r, c_, s);
VERBOSE_GOT_HERE(0);
DEBUG_MSG(0, "l=%f, u=%f, r=%f, c_=%f, s=%f", l, u, r, c_, s);
DEBUG_GOT_HERE(0);
return false;
}
@@ -334,18 +354,17 @@ bool SMO<TKernel>::TakeStep_(index_t i, index_t j, double error_j) {
double kii = EvalKernel_(i, i);
double kij = EvalKernel_(i, j);
double kjj = EvalKernel_(j, j);
// second derivative of the objective function
// second derivative of objective function
double eta = +2*kij - kii - kjj;
VERBOSE_MSG(0, "kij=%f, kii=%f, kjj=%f", kij, kii, kjj);
DEBUG_MSG(0, "kij=%f, kii=%f, kjj=%f", kij, kii, kjj);
// update alpha_j
if (likely(eta < 0)) {
VERBOSE_MSG(0, "Common case");
DEBUG_MSG(0, "Common case");
alpha_j = alpha_[j] - yj * (error_i - error_j) / eta;
alpha_j = FixAlpha_(math::ClampRange(alpha_j, l, u));
} else {
VERBOSE_MSG(0, "Uncommon case");
DEBUG_MSG(0, "Uncommon case");
//abort();
double c1 = eta/2;
double c2 = yj * (error_i - error_j) - eta * alpha_j;
@@ -368,11 +387,10 @@ bool SMO<TKernel>::TakeStep_(index_t i, index_t j, double error_j) {
// check if there is progress
if (fabs(delta_alpha_j) < SMO_EPS*(alpha_j + alpha_[j] + SMO_EPS)) {
VERBOSE_GOT_HERE(0);
DEBUG_GOT_HERE(0);
return false;
}
// update alpha_i
alpha_i = alpha_i - (s)*(delta_alpha_j);
if (alpha_i < SMO_ZERO) {
alpha_j += s * alpha_i;
@@ -403,7 +421,6 @@ bool SMO<TKernel>::TakeStep_(index_t i, index_t j, double error_j) {
kernel_cache_sign_.MakeColumnVector(i, &kernel_i);
kernel_cache_sign_.MakeColumnVector(j, &kernel_j);
// update gradient
for (index_t k = 0; k < n_data_; k++) {
if (likely(k != i) && likely(k != j) && !IsBound_(alpha_[k])) {
error_[k] += (delta_alpha_i*kernel_i[k] + delta_alpha_j*kernel_j[k]) * GetLabelSign_(k) - delta_thresh;
@@ -421,41 +438,8 @@ bool SMO<TKernel>::TakeStep_(index_t i, index_t j, double error_j) {
error_[i] = 0;
error_[j] = 0;
VERBOSE_GOT_HERE(0);
DEBUG_GOT_HERE(0);
return true;
}
template<typename TKernel>
double SMO<TKernel>::Evaluate_(index_t i) const {
// TODO: This only handles linear
Vector kernel_values;
double summation = 0;
kernel_cache_sign_.MakeColumnVector(i, &kernel_values);
summation = la::Dot(alpha_, kernel_values) * GetLabelSign_(i);
return (summation - thresh_);
}
// Get SVM results:coefficients, number and indecies of SVs
template<typename TKernel>
void SMO<TKernel>::GetSVM(ArrayList<index_t> &dataset_index, ArrayList<double> &coef, ArrayList<bool> &sv_indicator) {
coef.Init(n_data_);
for (index_t i = 0; i < n_data_; i++) {
if (alpha_[i] != 0) { // support vectors
coef[i] = alpha_[i] * GetLabelSign_(i);
sv_indicator[dataset_index[i]] = true;
n_sv_++;
//sv_index[i_sv]= dataset_bi_index[i];
//GetVector_(i, &source);
//support_vectors->MakeColumnVector(i_support, &dest);
//dest.CopyValues(source);
}
else {
coef[i] = 0;
}
}
}
#endif
fastlib/u/garryb/svm/svm.cc
+59 -155
View File
@@ -26,6 +26,7 @@ void DoSvmNormalize(Dataset* dataset) {
Matrix cov;
la::MulTransBInit(m, m, &cov);
Vector d;
@@ -33,7 +34,7 @@ void DoSvmNormalize(Dataset* dataset) {
Matrix ui; // the inverse of eigenvectors
//cov.PrintDebug("cov");
PASSED(la::EigenvectorsInit(cov, &d, &u));
la::EigenvectorsInit(cov, &d, &u);
la::TransposeInit(u, &ui);
for (index_t i = 0; i < d.length(); i++) {
@@ -57,170 +58,73 @@ void DoSvmNormalize(Dataset* dataset) {
}
//dataset->matrix().PrintDebug("m");
if (fx_param_bool(NULL, "save", 0)) {
fx_default_param(NULL, "kfold/save", "1");
dataset->WriteCsv("m_normalized.csv");
}
}
void CreateArtificialDataset(Dataset* dataset){
Matrix m;
index_t n = fx_param_int(NULL, "n", 30);
double offset = fx_param_double(NULL, "offset", 0.0);
double range = fx_param_double(NULL, "range", 1.0);
double slope = fx_param_double(NULL, "slope", 1.0);
double margin = fx_param_double(NULL, "margin", 1.0);
double var = fx_param_double(NULL, "var", 1.0);
double intercept = fx_param_double(NULL, "intercept", 0.0);
// 2 dimensional dataset, size n, 3 classes
m.Init(3, n);
for (index_t i = 0; i < n; i += 3) {
double x;
double y;
x = (rand() * range / RAND_MAX) + offset;
y = margin / 2 + (rand() * var / RAND_MAX);
m.set(0, i, x);
m.set(1, i, x*slope + y + intercept);
m.set(2, i, 0); // labels
x = (rand() * range / RAND_MAX) + offset;
y = margin / 2 + (rand() * var / RAND_MAX);
m.set(0, i+1, 10*x);
m.set(1, i+1, x*slope + y + intercept);
m.set(2, i+1, 1); // labels
x = (rand() * range / RAND_MAX) + offset;
y = margin / 2 + (rand() * var / RAND_MAX);
m.set(0, i+2, 20*x);
m.set(1, i+2, x*slope + y + intercept);
m.set(2, i+2, 2); // labels
}
data::Save("m.csv", m);
dataset->OwnMatrix(&m);
}
int LoadData(Dataset* dataset, String datafilename){
if (fx_param_exists(NULL, datafilename)) {
// when a data file is specified, use it.
if ( !PASSED(dataset->InitFromFile(fx_param_str_req(NULL, datafilename))) ) {
fprintf(stderr, "Couldn't open the data file.\n");
return 0;
}
}
else {
fprintf(stderr, "No data file exist. Generating artificial dataset.\n");
// otherwise, create an artificial dataset and save it to "m.csv"
CreateArtificialDataset(dataset);
}
if (fx_param_bool(NULL, "normalize", 1)) {
fprintf(stderr, "Normalizing\n");
DoSvmNormalize(dataset);
} else {
fprintf(stderr, "Skipping normalize\n");
}
return 1;
}
int main(int argc, char *argv[]) {
fx_init(argc, argv);
srand(time(NULL));
String mode = fx_param_str_req(NULL, "mode");
String kernel = fx_param_str_req(NULL, "kernel");
// TODO: more kernels to be supported
// Cross Validation Mode, need cross validation data
if(mode == "cv") {
fprintf(stderr, "SVM Cross Validation... \n");
// Load cross validation data
Dataset cvset;
if (LoadData(&cvset, "cv_data") == 0)
return 1;
if (kernel == "linear") {
SimpleCrossValidator< SVM<SVMLinearKernel> > cross_validator;
// Initialize n_folds_, confusion_matrix_; k_cv: number of cross-validation folds, need k_cv>1
cross_validator.Init(&cvset,cvset.n_labels(),fx_param_int_req(NULL,"k_cv"), fx_root, "svm");
// k_cv folds cross validation; (true): do training set permutation
cross_validator.Run(true);
cross_validator.confusion_matrix().PrintDebug("confusion matrix");
}
else if (kernel == "gaussian") {
SimpleCrossValidator< SVM<SVMRBFKernel> > cross_validator;
// Initialize n_folds_, confusion_matrix_; k_cv: number of cross-validation folds
cross_validator.Init(&cvset,cvset.n_labels(),fx_param_int_req(NULL,"k_cv"), fx_root, "svm");
// k_cv folds cross validation; (true): do training set permutation
cross_validator.Run(true);
cross_validator.confusion_matrix().PrintDebug("confusion matrix");
}
}
// Training Mode, need training data | Training + Testing(online) Mode, need training data + testing data
else if (mode=="train" || mode=="train_test"){
fprintf(stderr, "SVM Training... \n");
// Load training data
Dataset trainset;
if (LoadData(&trainset, "train_data") == 0) // TODO:param_req
Dataset dataset;
if (fx_param_exists(NULL, "data")) {
// if a data file is specified, use it.
if (!PASSED(dataset.InitFromFile(fx_param_str_req(NULL, "data")))) {
fprintf(stderr, "Couldn't open the data file.\n");
return 1;
}
} else {
// create an artificial dataset and save it to "m.csv"
// Begin SVM Training | Training and Testing
datanode *svm_module = fx_submodule(fx_root, NULL, "svm");
if (kernel == "linear") {
SVM<SVMLinearKernel> svm;
svm.InitTrain(trainset, trainset.n_labels(), svm_module);
if (mode=="train_test"){ // training and testing, thus no need to load model from file
fprintf(stderr, "SVM Classifying... \n");
// Load testing data
Dataset testset;
if (LoadData(&testset, "test_data") == 0) // TODO:param_req
return 1;
svm.BatchClassify(&testset, "test_labels");
}
}
else if (kernel == "gaussian") {
SVM<SVMRBFKernel> svm;
svm.InitTrain(trainset, trainset.n_labels(), svm_module);
if (mode=="train_test"){ // training and testing, thus no need to load model from file
fprintf(stderr, "SVM Classifying... \n");
// Load testing data
Dataset testset;
if (LoadData(&testset, "test_data") == 0) // TODO:param_req
return 1;
svm.BatchClassify(&testset, "test_labels"); // TODO:param_req
}
Matrix m;
index_t n = fx_param_int(NULL, "n", 30);
double offset = fx_param_double(NULL, "offset", 0.0);
double range = fx_param_double(NULL, "range", 1.0);
double slope = fx_param_double(NULL, "slope", 1.0);
double margin = fx_param_double(NULL, "margin", 1.0);
double var = fx_param_double(NULL, "var", 1.0);
double intercept = fx_param_double(NULL, "intercept", 0.0);
// 3 dimensional dataset, size n
m.Init(3, n);
for (index_t i = 0; i < n; i += 2) {
double x;
double y;
x = (rand() * range / RAND_MAX) + offset;
y = margin / 2 + (rand() * var / RAND_MAX);
m.set(0, i, x);
m.set(1, i, x*slope + y + intercept);
m.set(2, i, 0);
x = (rand() * range / RAND_MAX) + offset;
y = margin / 2 + (rand() * var / RAND_MAX);
m.set(0, i+1, x);
m.set(1, i+1, x*slope - y + intercept);
m.set(2, i+1, 1);
}
data::Save("m.csv", m);
dataset.OwnMatrix(&m);
}
// Testing(offline) Mode, need loading model file and testing data
else if (mode=="test") {
fprintf(stderr, "SVM Classifying... \n");
// Load testing data
Dataset testset;
if (LoadData(&testset, "test_data") == 0) // TODO:param_req
return 1;
// Begin Classification
datanode *svm_module = fx_submodule(fx_root, NULL, "svm");
if (kernel == "linear") {
SVM<SVMLinearKernel> svm;
svm.Init(testset, testset.n_labels(), svm_module); // TODO:n_labels() -> num_classes_
svm.LoadModelBatchClassify(&testset, "svm_model", "test_labels"); // TODO:param_req
}
else if (kernel == "gaussian") {
SVM<SVMRBFKernel> svm;
svm.Init(testset, testset.n_labels(), svm_module); // TODO:n_labels() -> num_classes_
svm.LoadModelBatchClassify(&testset, "svm_model", "test_labels"); // TODO:param_req
}
if (fx_param_bool(NULL, "normalize", 1)) {
fprintf(stderr, "Normalizing\n");
DoSvmNormalize(&dataset);
} else {
fprintf(stderr, "Skipping normalize\n");
}
if (fx_param_bool(NULL, "save", 0)) {
fx_default_param(NULL, "kfold/save", "1");
dataset.WriteCsv("normalized.csv");
}
SimpleCrossValidator< SVM<SVMRBFKernel> > cross_validator;
// k_cv: number of cross-validation folds
cross_validator.Init(&dataset, 2,fx_param_int_req(NULL,"k_cv"), fx_root, "svm");
cross_validator.Run(true);
cross_validator.confusion_matrix().PrintDebug("confusion matrix");
fx_done();
}
fastlib/u/garryb/svm/svm.h
+70 -413
View File
@@ -7,450 +7,107 @@
#include <typeinfo>
#define ID_LINEAR 0
#define ID_GAUSSIAN 1
struct SVMLinearKernel {
ArrayList<double> kpara_; // kernel parameters
void Init(datanode *node) { //TODO:node
kpara_.Init();
}
void GetName(String* kname) {
kname->Copy("linear");
}
int GetTypeId() {
return ID_LINEAR;
}
void Init(datanode *node) {}
void Copy(const SVMLinearKernel& other) {}
double Eval(const Vector& a, const Vector& b) const {
return la::Dot(a, b);
}
void SaveParam(FILE* fp) {
}
};
class SVMRBFKernel {
private:
double sigma_;
double gamma_;
public:
ArrayList<double> kpara_; // kernel parameters
void Init(datanode *node) { //TODO: node..
kpara_.Init(2);
kpara_[0] = fx_param_double_req(NULL, "sigma"); //sigma
kpara_[1] = -1.0 / (2 * math::Sqr(kpara_[0])); //gamma
void Init(datanode *node) {
sigma_ = fx_param_double_req(NULL, "sigma");
gamma_ = -1.0 / (2 * math::Sqr(sigma_));
}
void GetName(String* kname) {
kname->Copy("gaussian");
}
int GetTypeId() {
return ID_GAUSSIAN;
void Copy(const SVMRBFKernel& other) {
sigma_ = other.sigma_;
gamma_ = other.gamma_;
}
double Eval(const Vector& a, const Vector& b) const {
double distance_squared = la::DistanceSqEuclidean(a, b);
return exp(kpara_[1] * distance_squared);
return exp(gamma_ * distance_squared);
}
void SaveParam(FILE* fp) {
fprintf(fp, "sigma %g\n", kpara_[0]);
fprintf(fp, "gamma %g\n", kpara_[1]);
double sigma() const {
return sigma_;
}
};
template<typename TKernel>
class SVM {
private:
// models for the binary classifiers
struct SVM_MODELS {
double thresh_;
ArrayList<double> bi_coef_; // all coefficients of the binary dataset, not necessarily thoes of SVs
};
ArrayList<SVM_MODELS> models_;
ArrayList<int> train_labels_list_; // set of labels, need to be integers
// total set of support vectors and their coefficients
Matrix sv_;
Matrix sv_coef_;
index_t total_num_sv_;
ArrayList<bool> sv_indicator_;
ArrayList<index_t> sv_index_;
ArrayList<index_t> sv_list_startpos_;
ArrayList<index_t> sv_list_ct_;
struct SVM_PARAMETERS {
TKernel kernel_;
String kernelname_;
int kerneltypeid_;
double c_;
int b_;
};
SVM_PARAMETERS param_; // same for every binary model
int num_classes_;
int num_models_;
int num_features_;
public:
typedef TKernel Kernel;
void Init(const Dataset& dataset, int n_classes, datanode *module);
private:
TKernel kernel_;
double c_;
int b_;
Matrix support_vectors_;
double thresh_;
Vector alpha_;
public:
void InitTrain(const Dataset& dataset, int n_classes, datanode *module);
void SaveModel(String modelfilename);
void LoadModel(Dataset* testset, String modelfilename);
int Classify(const Vector& vector);
void BatchClassify(Dataset* testset, String testlabelfilename);
void LoadModelBatchClassify(Dataset* testset, String modelfilename, String testlabelfilename);
};
template<typename TKernel>
void SVM<TKernel>::Init(const Dataset& dataset, int n_classes, datanode *module){
models_.Init();
sv_indicator_.Init(dataset.n_points());
sv_index_.Init();
void SVM<TKernel>::InitTrain(
const Dataset& dataset, int n_classes, datanode *module) {
DEBUG_ASSERT_MSG(n_classes == 2, "SVM is only a binary classifier");
fx_set_param(module, "kernel_type", typeid(TKernel).name());
kernel_.Init(fx_submodule(module, "kernel", "kernel"));
c_ = fx_param_double(module, "c", 1.0);
b_ = fx_param_int(module, "b", dataset.n_points());
SMO<Kernel> smo;
smo.Init(&dataset, c_, b_);
smo.kernel().Copy(kernel_);
smo.Train();
thresh_ = smo.threshold();
smo.GetSVM(&support_vectors_, &alpha_);
DEBUG_ASSERT(alpha_.length() != 0);
DEBUG_ASSERT(alpha_.length() == support_vectors_.n_cols());
DEBUG_ONLY(fprintf(stderr, "----------------------\n"));
DEBUG_ONLY(support_vectors_.PrintDebug("support vectors"));
DEBUG_ONLY(alpha_.PrintDebug("support vector weights"));
DEBUG_ONLY(fprintf(stderr, "-- THRESHOLD: %f\n", thresh_));
param_.kernel_.Init(fx_submodule(module, "kernel", "kernel"));
param_.kernel_.GetName(&param_.kernelname_);
param_.kerneltypeid_ = param_.kernel_.GetTypeId();
param_.c_ = fx_param_double_req(NULL, "c");
// budget parameter, contorls # of support vectors; default: # of data samples (use all)
param_.b_ = fx_param_int(module, "b", dataset.n_points()); // TODO: param_req
num_classes_ = n_classes;
num_models_ = num_classes_ * (num_classes_-1) / 2;
num_features_ = 0;
sv_list_startpos_.Init(n_classes);
sv_list_ct_.Init(n_classes);
for (index_t i=0; i<dataset.n_points(); i++)
sv_indicator_[i] = false;
total_num_sv_ = 0;
fx_format_result(module, "n_support", "%"LI"d", alpha_.length());
}
// Multiclass SVM Classifier. Initilization and Training
// use One-vs-One, or called All-vs-All method
template<typename TKernel>
void SVM<TKernel>::InitTrain(const Dataset& dataset, int n_classes, datanode *module) {
Init(dataset, n_classes, module);
num_features_ = dataset.n_features()-1; // last column in dataset is for labels
// Group labels, split the training dataset for training bi-class SVM classifiers
ArrayList<index_t> train_labels_index; // e.g. [c1[0,5,6,7,10,13,17],c2[1,2,4,8,9],c3[...]]
ArrayList<index_t> train_labels_ct; // counter, e.g. [7,5,8]
ArrayList<index_t> train_labels_startpos; // start position, e.g. [0,7,12]
dataset.GetLabels(train_labels_list_, train_labels_index, train_labels_ct, train_labels_startpos);
// Train n_classes*(n_classes-1)/2 bi-class(labels: 0, 1) models using SMO
index_t ct = 0;
index_t i; index_t j;
for (i = 0; i < n_classes; i++) {
for (j = i+1; j < n_classes; j++) {
models_.AddBack();
// Initialize parameters c_, budget_, alpha_, error_, thresh_
SMO<Kernel> smo;
smo.Init(n_classes, param_.c_, param_.b_);
smo.kernel().Init(fx_submodule(module, "kernel", "kernel"));
// Construct dataset consists of two classes i and j (reassign labels 0 and 1)
Dataset dataset_bi;
dataset_bi.InitBlank();
dataset_bi.info().Init();
dataset_bi.matrix().Init(num_features_+1, train_labels_ct[i]+train_labels_ct[j]);
ArrayList<index_t> dataset_bi_index;
dataset_bi_index.Init(train_labels_ct[i]+train_labels_ct[j]);
for (index_t m = 0; m < train_labels_ct[i]; m++) {
Vector source, dest;
dataset_bi.matrix().MakeColumnVector(m, &dest);
dataset.matrix().MakeColumnVector(train_labels_index[train_labels_startpos[i]+m], &source);
dest.CopyValues(source);
dataset_bi.matrix().set(num_features_, m, 0); // last row for labels 0
dataset_bi_index[m] = train_labels_index[train_labels_startpos[i]+m];
}
for (index_t n = 0; n < train_labels_ct[j]; n++) {
Vector source, dest;
dataset_bi.matrix().MakeColumnVector(n+train_labels_ct[i], &dest);
dataset.matrix().MakeColumnVector(train_labels_index[train_labels_startpos[j]+n], &source);
dest.CopyValues(source);
dataset_bi.matrix().set(num_features_, n+train_labels_ct[i], 1); // last row for labels 1
dataset_bi_index[n+train_labels_ct[i]] = train_labels_index[train_labels_startpos[j]+n];
}
// 2-classes SVM-SMO training
smo.Train(&dataset_bi);
// Get the trained bi-class model
models_[ct].thresh_ = smo.threshold();
smo.GetSVM(dataset_bi_index, models_[ct].bi_coef_, sv_indicator_);
//smo.GetSVM(dataset_bi_index, models_[ct].bi_coef_, sv_indicator_, &(models_[ct].num_sv_bi_));
//DEBUG_ASSERT(models_[ct].alpha_.length() != 0);
//DEBUG_ASSERT(models_[ct].alpha_.length() == models_[ct].support_vectors_.n_cols());
//DEBUG_ONLY(fprintf(stderr, "----------------------\n"));
//DEBUG_ONLY(models_[ct].support_vectors_.PrintDebug("support vectors"));
//DEBUG_ONLY(models_[ct].alpha_.PrintDebug("support vector weights"));
//DEBUG_ONLY(fprintf(stderr, "-- THRESHOLD: %f\n", models_[ct].thresh_));
//fx_format_result(module, "n_support", "%"LI"d", models_[ct].alpha_.length());
//fx_format_result(module, "n_support", "%"LI"d", models_[ct].sv_bi_coef_.size());
ct++;
}
}
// Get total set of SVs;
index_t k;
sv_list_startpos_[0] = 0;
total_num_sv_ = 0;
for (i = 0; i < n_classes; i++) {
ct = 0;
for (j = 0; j < train_labels_ct[i]; j++) {
if (sv_indicator_[train_labels_startpos[i] + j]) {
*sv_index_.AddBack() = train_labels_index[train_labels_startpos[i] + j];
total_num_sv_++;
ct++;
}
}
sv_list_ct_[i] = ct;
if (i >= 1)
sv_list_startpos_[i] = sv_list_startpos_[i-1] + sv_list_ct_[i-1];
}
sv_.Init(num_features_, total_num_sv_);
for (i = 0; i < total_num_sv_; i++) {
Vector source, dest;
sv_.MakeColumnVector(i, &dest);
dataset.matrix().MakeColumnSubvector(sv_index_[i], 0, num_features_, &source); // last row of dataset is for labels
dest.CopyValues(source);
}
// Get coefficients for the total set of SVs, i.e. models_[x].bi_coef_ -> sv_coef_
index_t p;
index_t ct_model = 0;
sv_coef_.Init(n_classes-1, total_num_sv_);
sv_coef_.SetZero();
for (i = 0; i < n_classes; i++) {
for (j = i+1; j < n_classes; j++) {
p = sv_list_startpos_[i];
for (k = 0; k < train_labels_ct[i]; k++) {
if (sv_indicator_[train_labels_startpos[i]+k]) {
sv_coef_.set(j-1, p++, models_[ct_model].bi_coef_[k]);
}
}
p = sv_list_startpos_[j];
for (k = 0; k < train_labels_ct[j]; k++) {
if (sv_indicator_[train_labels_startpos[j]+k]) {
sv_coef_.set(i, p++, models_[ct_model].bi_coef_[k]);
}
}
ct_model++;
}
}
// Save models
SaveModel("svm_model"); // TODO: param_req, and for CV mode
}
// Save SVM model to file, TODO: use XML
template<typename TKernel>
void SVM<TKernel>::SaveModel(String modelfilename) {
FILE *fp = fopen(modelfilename,"w");
if (fp==NULL)
fprintf(stderr, "Cannot save trained model to file!");
index_t i; index_t j;
fprintf(fp, "svm_type svm_c\n"); // TODO: svm-mu, svm-regression...
fprintf(fp, "num_classes %d\n", num_classes_); // TODO: only for svm_c
fprintf(fp, "kernel_name %s\n", param_.kernelname_.c_str());
fprintf(fp, "kernel_typeid %d\n", param_.kerneltypeid_);
// save kernel parameters
param_.kernel_.SaveParam(fp);
fprintf(fp, "total_num_sv %d\n", total_num_sv_);
fprintf(fp, "labels ");
for (i = 0; i < num_classes_; i++)
fprintf(fp, "%d ", train_labels_list_[i]);
fprintf(fp, "\n");
// save models
fprintf(fp, "thresholds ");
for (i = 0; i < num_models_; i++)
fprintf(fp, "%f ", models_[i].thresh_);
fprintf(fp, "\n");
fprintf(fp, "sv_list_startpos ");
for (i =0; i < num_classes_; i++)
fprintf(fp, "%d ", sv_list_startpos_[i]);
fprintf(fp, "\n");
fprintf(fp, "sv_list_ct ");
for (i =0; i < num_classes_; i++)
fprintf(fp, "%d ", sv_list_ct_[i]);
fprintf(fp, "\n");
// save coefficients and support vectors
fprintf(fp, "SV_coefs\n");
for (i = 0; i < total_num_sv_; i++) {
for (j = 0; j < num_classes_-1; j++) {
fprintf(fp, "%f ", sv_coef_.get(j,i));
}
fprintf(fp, "\n");
}
fprintf(fp, "SVs\n");
for (i = 0; i < total_num_sv_; i++) {
for (j = 0; j < num_features_; j++) { // n_rows-1
fprintf(fp, "%f ", sv_.get(j,i));
}
fprintf(fp, "\n");
}
fclose(fp);
}
// Load SVM model file, TODO: use XML
template<typename TKernel>
void SVM<TKernel>::LoadModel(Dataset* testset, String modelfilename) {
//Init
train_labels_list_.Init(num_classes_);
num_features_ = testset->n_features();
//load model file
FILE *fp = fopen(modelfilename, "r");
char cmd[80];
int i, j; int temp_d; double temp_f;
for (i = 0; i < num_models_; i++) {
models_.AddBack();
}
while (1) {
fscanf(fp,"%80s",cmd);
if(strcmp(cmd,"svm_type")==0) {
fscanf(fp,"%80s",cmd);
if(strcmp(cmd,"svm_c")==0) {
fprintf(stderr, "SVM_C\n");
}
}
else if (strcmp(cmd, "num_classes")==0) {
fscanf(fp,"%d",&num_classes_);
}
else if (strcmp(cmd, "kernel_name")==0) {
fscanf(fp,"%80s",param_.kernelname_.c_str());
}
else if (strcmp(cmd, "kernel_typeid")==0) {
fscanf(fp,"%d",&param_.kerneltypeid_);
}
else if (strcmp(cmd, "sigma")==0) {
fscanf(fp,"%lf",&param_.kernel_.kpara_[0]); // for gaussian kernels
}
else if (strcmp(cmd, "gamma")==0) {
fscanf(fp,"%lf",&param_.kernel_.kpara_[1]); // for gaussian kernels
}
else if (strcmp(cmd, "total_num_sv")==0) {
fscanf(fp,"%d",&total_num_sv_);
}
else if (strcmp(cmd, "labels")==0) {
for (i=0; i<num_classes_; i++) {
fscanf(fp,"%d",&temp_d);
train_labels_list_[i] = temp_d;
}
}
else if (strcmp(cmd, "thresholds")==0) {
for ( i= 0; i < num_models_; i++) {
fscanf(fp,"%d",&temp_d);
models_[i].thresh_= temp_d;
}
}
else if (strcmp(cmd, "sv_list_startpos")==0) {
for ( i= 0; i < num_models_; i++) {
fscanf(fp,"%d",&temp_d);
sv_list_startpos_[i]= temp_d;
}
}
else if (strcmp(cmd, "sv_list_ct")==0) {
for ( i= 0; i < num_models_; i++) {
fscanf(fp,"%d",&temp_d);
sv_list_ct_[i]= temp_d;
}
break;
}
}
sv_coef_.Init(num_classes_-1, total_num_sv_);
sv_coef_.SetZero();
sv_.Init(num_features_, total_num_sv_);
while (1) {
fscanf(fp,"%80s",cmd);
if (strcmp(cmd, "SV_coefs")==0) {
for (i = 0; i < total_num_sv_; i++) {
for (j = 0; j < num_classes_-1; j++) {
fscanf(fp,"%lf",&temp_f);
sv_coef_.set(j, i, temp_f);
}
}
}
else if (strcmp(cmd, "SVs")==0) {
for (i = 0; i < total_num_sv_; i++) {
for (j = 0; j < num_features_; j++) {
fscanf(fp,"%lf",&temp_f);
sv_.set(j, i, temp_f);
}
}
break;
}
}
}
// Multiclass SVM Classifier. Testing for one sample
template<typename TKernel>
int SVM<TKernel>::Classify(const Vector& datum) {
index_t i; index_t j; index_t k;
ArrayList<double> keval;
keval.Init(total_num_sv_);
for (i = 0; i < total_num_sv_; i++) {
Vector support_vector_i;
sv_.MakeColumnVector(i, &support_vector_i);
keval[i] = param_.kernel_.Eval(datum, support_vector_i);
double summation = 0;
for (index_t i = 0; i < alpha_.length(); i++) {
Vector support_vector;
support_vectors_.MakeColumnVector(i, &support_vector);
double term = alpha_[i] * kernel_.Eval(datum, support_vector);
DEBUG_MSG(0, "alpha %f, term %f", alpha_[i], term);
summation += term;
}
ArrayList<double> values;
values.Init(num_models_);
index_t ct = 0;
for (i = 0; i < num_classes_; i++) {
for (j = i+1; j < num_classes_; j++) {
double sum = 0;
for(k = 0; k < sv_list_ct_[i]; k++) {
sum += sv_coef_.get(j-1, sv_list_startpos_[i]+k) * keval[sv_list_startpos_[i]+k];
}
for(k = 0; k < sv_list_ct_[j]; k++)
sum += sv_coef_.get(i, sv_list_startpos_[j]+k) * keval[sv_list_startpos_[j]+k];
sum -= models_[ct].thresh_;
values[ct] = sum;
fprintf(stderr, "%f\n", values[ct]);
ct++;
}
}
ArrayList<index_t> vote;
vote.Init(num_classes_);
ct = 0;
for (i = 0; i < num_classes_; i++) {
for (j = i+1; j < num_classes_; j++) {
if(values[ct] > 0)
++vote[i];
else
++vote[j];
ct++;
}
}
index_t vote_max_idx = 0;
for (i = 1; i < num_classes_; i++) {
if (vote[i] > vote[vote_max_idx])
vote_max_idx = i;
}
return train_labels_list_[vote_max_idx];
}
// Batch classification for multiple testing vectors, no need to load model file
// Note: for test set, if no true test labels provided, just put some dummy labels (e.g. all -1) in the last row of testset
template<typename TKernel>
void SVM<TKernel>::BatchClassify(Dataset* testset, String testlablefilename) {
FILE *fp = fopen(testlablefilename,"w");
if (fp==NULL)
fprintf(stderr, "Cannot save test labels to file!");
num_features_ = testset->n_features()-1;
for (index_t i=0; i<testset->n_points(); i++) {
Vector testvec;
testset->matrix().MakeColumnSubvector(i, 0, num_features_, &testvec);
int testlabel = Classify(testvec);
fprintf(fp, "%d\n", testlabel);
}
fclose(fp);
}
// Load models from a file, and perform batch classification for multiple testing vectors
template<typename TKernel>
void SVM<TKernel>::LoadModelBatchClassify(Dataset* testset, String modelfilename, String testlabelfilename) {
LoadModel(testset, modelfilename);
BatchClassify(testset, testlabelfilename);
DEBUG_MSG(0, "summation=%f, thresh_=%f", summation, thresh_);
return (summation - thresh_ > 0.0) ? 1 : 0;
}
#endif
fastlib/u/houyang/svm/svm.h
+48 -19
View File
@@ -108,7 +108,7 @@ void SVM<TKernel>::Init(const Dataset& dataset, int n_classes, datanode *module)
param_.kernel_.Init(fx_submodule(module, "kernel", "kernel"));
param_.kernel_.GetName(&param_.kernelname_);
param_.kerneltypeid_ = param_.kernel_.GetTypeId();
// c; default:1
param_.c_ = fx_param_double_req(NULL, "c");
// budget parameter, contorls # of support vectors; default: # of data samples (use all)
param_.b_ = fx_param_int(module, "b", dataset.n_points()); // TODO: param_req
@@ -255,25 +255,33 @@ void SVM<TKernel>::SaveModel(String modelfilename) {
param_.kernel_.SaveParam(fp);
fprintf(fp, "total_num_sv %d\n", total_num_sv_);
fprintf(fp, "labels ");
for (i=0; i<num_classes_; i++)
for (i = 0; i < num_classes_; i++)
fprintf(fp, "%d ", train_labels_list_[i]);
fprintf(fp, "\n");
// save models
fprintf(fp, "thresholds ");
for (i=0; i< num_models_; i++)
for (i = 0; i < num_models_; i++)
fprintf(fp, "%f ", models_[i].thresh_);
fprintf(fp, "\n");
fprintf(fp, "sv_list_startpos ");
for (i =0; i < num_classes_; i++)
fprintf(fp, "%d ", sv_list_startpos_[i]);
fprintf(fp, "\n");
fprintf(fp, "sv_list_ct ");
for (i =0; i < num_classes_; i++)
fprintf(fp, "%d ", sv_list_ct_[i]);
fprintf(fp, "\n");
// save coefficients and support vectors
fprintf(fp, "SV_coefs\n");
for (i=0; i<total_num_sv_; i++) {
for (j=0; j<num_classes_-1; j++) {
for (i = 0; i < total_num_sv_; i++) {
for (j = 0; j < num_classes_-1; j++) {
fprintf(fp, "%f ", sv_coef_.get(j,i));
}
fprintf(fp, "\n");
}
fprintf(fp, "SVs\n");
for (i=0; i<total_num_sv_; i++) {
for (j =0; j<num_features_; j++) { // n_rows-1
for (i = 0; i < total_num_sv_; i++) {
for (j = 0; j < num_features_; j++) { // n_rows-1
fprintf(fp, "%f ", sv_.get(j,i));
}
fprintf(fp, "\n");
@@ -287,11 +295,12 @@ void SVM<TKernel>::LoadModel(Dataset* testset, String modelfilename) {
//Init
train_labels_list_.Init(num_classes_);
num_features_ = testset->n_features();
//load model file
FILE *fp = fopen(modelfilename, "r");
char cmd[80];
int i, j; int temp_d; double temp_f;
for (i=0; i<num_models_; i++) {
for (i = 0; i < num_models_; i++) {
models_.AddBack();
}
while (1) {
@@ -299,7 +308,7 @@ void SVM<TKernel>::LoadModel(Dataset* testset, String modelfilename) {
if(strcmp(cmd,"svm_type")==0) {
fscanf(fp,"%80s",cmd);
if(strcmp(cmd,"svm_c")==0) {
fprintf(stderr, "SVM Classification");
fprintf(stderr, "SVM_C\n");
}
}
else if (strcmp(cmd, "num_classes")==0) {
@@ -327,29 +336,48 @@ void SVM<TKernel>::LoadModel(Dataset* testset, String modelfilename) {
}
}
else if (strcmp(cmd, "thresholds")==0) {
for (i=0; i<num_models_; i++) {
for ( i= 0; i < num_models_; i++) {
fscanf(fp,"%d",&temp_d);
models_[i].thresh_= temp_d;
}
}
else if (strcmp(cmd, "SV_coef")==0) {
for (i=0; i<total_num_sv_; i++) {
for (j=0; j<num_classes_-1; j++) {
fscanf(fp,"%f",&temp_f);
else if (strcmp(cmd, "sv_list_startpos")==0) {
for ( i= 0; i < num_models_; i++) {
fscanf(fp,"%d",&temp_d);
sv_list_startpos_[i]= temp_d;
}
}
else if (strcmp(cmd, "sv_list_ct")==0) {
for ( i= 0; i < num_models_; i++) {
fscanf(fp,"%d",&temp_d);
sv_list_ct_[i]= temp_d;
}
break;
}
}
sv_coef_.Init(num_classes_-1, total_num_sv_);
sv_coef_.SetZero();
sv_.Init(num_features_, total_num_sv_);
while (1) {
fscanf(fp,"%80s",cmd);
if (strcmp(cmd, "SV_coefs")==0) {
for (i = 0; i < total_num_sv_; i++) {
for (j = 0; j < num_classes_-1; j++) {
fscanf(fp,"%lf",&temp_f);
sv_coef_.set(j, i, temp_f);
}
}
}
else if (strcmp(cmd, "SVs")==0) {
for (i=0; i<total_num_sv_; i++) {
for (j=0; j<num_features_; j++) {
fscanf(fp,"%f",&temp_f);
for (i = 0; i < total_num_sv_; i++) {
for (j = 0; j < num_features_; j++) {
fscanf(fp,"%lf",&temp_f);
sv_.set(j, i, temp_f);
}
}
break;
}
}// while
}
}
// Multiclass SVM Classifier. Testing for one sample
@@ -369,8 +397,9 @@ int SVM<TKernel>::Classify(const Vector& datum) {
for (i = 0; i < num_classes_; i++) {
for (j = i+1; j < num_classes_; j++) {
double sum = 0;
for(k = 0; k < sv_list_ct_[i]; k++)
for(k = 0; k < sv_list_ct_[i]; k++) {
sum += sv_coef_.get(j-1, sv_list_startpos_[i]+k) * keval[sv_list_startpos_[i]+k];
}
for(k = 0; k < sv_list_ct_[j]; k++)
sum += sv_coef_.get(i, sv_list_startpos_[j]+k) * keval[sv_list_startpos_[j]+k];
sum -= models_[ct].thresh_;