svm_c.cpp

支持向量机(SVM)的VC源代码

      cout<<"Precision : "<<the_result.precision<<" ("<<the_result.pred_precision<<")"<<endl;
      cout<<"Recall    : "<<the_result.recall<<" ("<<the_result.pred_recall<<")"<<endl;
      if(parameters->verbosity>= 2){
	// nice printout ;-)
	int rows = (int)(1+log10((SVMFLOAT)(total_pos+total_neg)));
	int now_digits = rows+2;
	int i,j;
	cout<<endl;
	cout<<"Predicted values:"<<endl;
	cout<<"   |";
	for(i=0;i<rows;i++){ cout<<" "; };
	cout<<"+  |";
	for(j=0;j<rows;j++){ cout<<" "; };
	cout<<"-"<<endl;

	cout<<"---+";
	for(i=0;i<now_digits;i++){ cout<<"-"; };
	cout<<"-+-";
	for(i=0;i<now_digits;i++){ cout<<"-"; };
	cout<<endl;

	cout<<" + |  ";
	now_digits=rows-(int)(1+log10((SVMFLOAT)correct_pos))-1;
	for(i=0;i<now_digits;i++){ cout<<" "; };
	cout<<correct_pos<<"  |  ";
	now_digits=rows-(int)(1+log10((SVMFLOAT)(total_pos-correct_pos)))-1;
	for(i=0;i<now_digits;i++){ cout<<" "; };
	cout<<total_pos-correct_pos<<"    (true pos)"<<endl;

	cout<<" - |  ";
	now_digits=rows-(int)(1+log10((SVMFLOAT)(total_neg-correct_neg)))-1;
	for(i=0;i<now_digits;i++){ cout<<" "; };
	cout<<(total_neg-correct_neg)<<"  |  ";
	now_digits=rows-(int)(1+log10((SVMFLOAT)correct_neg))-1;
	for(i=0;i<now_digits;i++){ cout<<" "; };
	cout<<correct_neg<<"    (true neg)"<<endl;
	cout<<endl;
      };
    };
  };

  SVMINT dim = examples->get_dim();
  if((dim<30) && (parameters->verbosity>= 2)){ // && (parameters->is_linear != 0)){
    // print hyperplane
    SVMINT j;
    svm_example example;
    SVMFLOAT* w = new SVMFLOAT[dim];
    SVMFLOAT b = examples->get_b();
    for(j=0;j<dim;j++) w[j] = 0;
    for(i=0;i<examples_total;i++){
      example = examples->get_example(i);
      alpha = examples->get_alpha(i);
      for(j=0;j<example.length;j++){
	w[((example.example)[j]).index] += alpha*((example.example)[j]).att;
      };
    };
    if(examples->initialised_scale()){
      SVMFLOAT* exp = examples->get_exp();
      SVMFLOAT* var = examples->get_var();
      for(j=0;j<dim;j++){
	if(var[j] != 0){
	  w[j] /= var[j];
	};
	b -= w[j]*exp[j];
      };
    };
    for(j=0;j<dim;j++){
      cout << "w["<<j<<"] = " << w[j] << endl;
    };
    cout << "b = "<<b<<endl;
    if(dim==1){
      cout<<"y = "<<w[0]<<"*x+"<<b<<endl;
    };
    if((dim==2) && (is_pattern)){
      cout<<"x1 = "<<-w[0]/w[1]<<"*x0+"<<-b/w[1]<<endl;
    };
    delete []w;
  };

  if(parameters->verbosity>= 2){
    cout<<"Time for learning:"<<endl
	<<"init        : "<<(time_init/100)<<"s"<<endl
	<<"optimizer   : "<<(time_optimize/100)<<"s"<<endl
	<<"convergence : "<<(time_convergence/100)<<"s"<<endl
	<<"update ws   : "<<(time_update/100)<<"s"<<endl
	<<"calc ws     : "<<(time_calc/100)<<"s"<<endl
	<<"============="<<endl
	<<"all         : "<<(time_all/100)<<"s"<<endl;
  }
  else if(parameters->verbosity>=2){
    cout<<"Time for learning: "<<(time_all/100)<<"s"<<endl;
  };


  //////////////////////////////////////

//   cout<<"############################################"<<endl;

//   cout<<"statistics for SVs:"<<endl;
//   cout<<"alpha\tf(x)-y\tlambda"<<endl;
//   SVMFLOAT sum_diff;
//   SVMFLOAT* kernel_row;
//   for(SVMINT i=0;i<examples_total;i++){
//     if(all_alphas[i] != 0){
//       sum_diff = sum[i];
//       kernel_row = kernel->get_row(i);
//       for(SVMINT j=0;j<examples_total;j++){
// 	sum_diff -= all_alphas[j]*kernel_row[j];
//       };
//       cout<<all_alphas[i]<<"\t"<<predict(i)-all_ys[i]<<"\t"<<lambda(i)<<"\t"<<sum_diff<<endl;
//     };
//   };


//   cout<<"statistics for non-SVs:"<<endl;
//   cout<<"alpha\tf(x)-y\tlambda+\tlambda-"<<endl;
//   for(SVMINT i=0;i<examples_total;i++){
//     if(all_alphas[i] == 0){
//       sum_diff = sum[i];
//       kernel_row = kernel->get_row(i);
//       for(SVMINT j=0;j<examples_total;j++){
// 	sum_diff -= all_alphas[j]*kernel_row[j];
//       };
//       cout<<all_alphas[i]<<"\t"<<predict(i)-all_ys[i]<<"\t"<<sum[i] - all_ys[i] + epsilon_neg + lambda_eq<<"\t"<<-sum[i] + all_ys[i] + epsilon_pos-lambda_eq<<"\t"<<sum_diff<<endl;
//     };
//   };


//   cout<<"############################################"<<endl;

  ///////////////////////////////////////


  return the_result;
};


/**
 *
 * Pattern SVM
 *
 */

int svm_pattern_c::is_alpha_neg(const SVMINT i){
  // variable i is alpha*

  int result;
  if(all_ys[i] > 0){
    result = 1;
  }
  else{
    result = -1;
  };
  return result;
};


SVMFLOAT svm_pattern_c::nabla(const SVMINT i){
  // = is_alpha_neg(i) * sum[i] - 1
  if(all_ys[i]>0){
    return(sum[i]-1);
  }
  else{
    return(-sum[i]-1);
  };
  //  return(all_ys[i]*sum[i] - 1);
};


SVMFLOAT svm_pattern_c::lambda(const SVMINT i){
  // size lagrangian multiplier of the active constraint

  SVMFLOAT alpha;
  SVMFLOAT result=0;

  alpha=all_alphas[i];

  if(alpha>is_zero){
    // alpha*
    if(alpha-Cneg >= - is_zero){
      // upper bound active
      result = -lambda_eq - sum[i] + 1;
    };
  }
  else if(alpha >= -is_zero){
    // lower bound active
    if(all_ys[i]>0){
      result = (sum[i]+lambda_eq) - 1;
    }
    else{
      result = -(sum[i]+lambda_eq) - 1;
    };
  }
  else if(alpha+Cpos <= is_zero){
    // upper bound active
    result = lambda_eq + sum[i] + 1;
  }
  else{
    if(all_ys[i]>0){
      result = -abs(sum[i]+lambda_eq - 1);
    }
    else{
      result = -abs(-sum[i]-lambda_eq - 1);
    };
  };

  return result;
};


int svm_pattern_c::feasible(const SVMINT i){
  // is direction i feasible to minimize the target function
  // (includes which_alpha==0)

  if(at_bound[i] >= shrink_const){ return 0; };

  SVMFLOAT alpha;
  SVMFLOAT result;

  alpha=all_alphas[i];

  if(alpha-Cneg >= - is_zero){
    // alpha* at upper bound
    if(all_ys[i]>0){
      result = -lambda_eq - sum[i] + 1;
    }
    else{
      result = -lambda_eq + sum[i] + 1;
    };
    if(result>=-feasible_epsilon){
      return 0; 
    };
  }
  else if((alpha<=is_zero) && (alpha >= -is_zero)){
    // lower bound active
    if(all_ys[i]>0){
      result = sum[i]+lambda_eq - 1;
    }
    else{
      result = -sum[i]-lambda_eq - 1;
    };
    if(result>=-feasible_epsilon){
      return 0; 
    };
  }
  else if(alpha+Cpos <= is_zero){
    // alpha at upper bound
    if(all_ys[i]>0){
      result = lambda_eq - sum[i] + 1;
    }
    else{
      result = lambda_eq + sum[i] + 1;
    };
    if(result>=-feasible_epsilon){
      return 0; 
    };
  }
  else{
    // not at bound
    if(all_ys[i]>0){
      result= abs(sum[i]-1+lambda_eq);
    }
    else{
      result= abs(sum[i]+1+lambda_eq);
    };
    if(result<=feasible_epsilon){
      return 0; 
    };
  };
  return 1;
};


int svm_pattern_c::feasible(const SVMINT i, SVMFLOAT* the_nabla, SVMFLOAT* the_lambda, int* atbound){
  // is direction i feasible to minimize the target function
  // (includes which_alpha==0)
  int is_feasible=1;

  if(at_bound[i] >= shrink_const){ is_feasible = 0; };

  SVMFLOAT alpha;

  alpha=all_alphas[i];
  //  feasible_epsilon=-1;
  if(alpha == Cneg){ //alpha-Cneg >= - is_zero){
    // alpha* at upper bound
    *atbound = 1;
    *the_nabla = sum[i] - 1;
    *the_lambda = -lambda_eq - *the_nabla; //sum[i] + 1;
    if(*the_lambda >= lambda_threshold){
      at_bound[i]++;
      if(at_bound[i] == shrink_const) to_shrink++;
    }
    else{
      at_bound[i] = 0;
    };
  }
  else if(alpha == 0){
    // lower bound active
    *atbound = -1;
    if(all_ys[i]>0){
      *the_nabla = sum[i] - 1;
      *the_lambda = lambda_eq + *the_nabla; //sum[i]+lambda_eq - 1;
    }
    else{
      *the_nabla = -sum[i] - 1;
      *the_lambda = -lambda_eq + *the_nabla; //-sum[i]-lambda_eq - 1;
    };
    if(*the_lambda >= lambda_threshold){
      at_bound[i]++;
      if(at_bound[i] == shrink_const) to_shrink++;
    }
    else{
      at_bound[i] = 0;
    };
  }
  else if(alpha == -Cpos){ //alpha+Cpos <= is_zero){
    *atbound = 1;
    *the_nabla = -sum[i] - 1;
    *the_lambda = lambda_eq - *the_nabla; //sum[i] - 1;
    if(*the_lambda >= lambda_threshold){
      at_bound[i]++;
      if(at_bound[i] == shrink_const) to_shrink++;
    }
    else{
      at_bound[i] = 0;
    };
  }
  else{
    // not at bound
    *atbound = 0;
    if(all_ys[i]>0){
      *the_nabla = sum[i] - 1;
      *the_lambda = -abs(*the_nabla+lambda_eq);
    }
    else{
      *the_nabla = -sum[i] - 1;
      *the_lambda = -abs(lambda_eq - *the_nabla);
    };
    at_bound[i] = 0;
  };
  if(*the_lambda >= -feasible_epsilon){
    is_feasible = 0; 
  };

  return is_feasible;
};


void svm_pattern_c::init_optimizer(){
  // Cs are dived by examples_total in init_optimizer
  svm_c::init_optimizer();
  for(SVMINT i=0;i<working_set_size;i++){
    qp.l[i] = 0;
  };
};


/**
 *
 * Regression SVM
 *
 */

int svm_regression_c::is_alpha_neg(const SVMINT i){
  // variable i is alpha*
  int result;
  if(all_alphas[i] > 0){
    result = 1;
  }
  else if(all_alphas[i] == 0){
    if(sum[i] - all_ys[i] + lambda_eq>0){
      result = -1;
    }
    else{
      result = 1;
    };
    //    result = 2*((i+at_bound[i])%2)-1;
  }
  else{
    result = -1;
  };
  return result;
};


SVMFLOAT svm_regression_c::nabla(const SVMINT i){
  if(all_alphas[i] > 0){
    return( sum[i] - all_ys[i] + epsilon_neg);
  }
  else if(all_alphas[i] == 0){
    if(is_alpha_neg(i)>0){
      return( sum[i] - all_ys[i] + epsilon_neg);
    }
    else{
      return(-sum[i] + all_ys[i] + epsilon_pos);
    };
  }
  else{
    return(-sum[i] + all_ys[i] + epsilon_pos);
  };
};


SVMFLOAT svm_regression_c::lambda(const SVMINT i){
  // size lagrangian multiplier of the active constraint

  SVMFLOAT alpha;
  SVMFLOAT result = -abs(nabla(i)+is_alpha_neg(i)*lambda_eq);
    //= -infinity; // default = not at bound

  alpha=all_alphas[i];

  if(alpha>is_zero){
    // alpha*
    if(alpha-Cneg >= - is_zero){
      // upper bound active
      result = -lambda_eq-sum[i] + all_ys[i] - epsilon_neg;
    };
  }
  else if(alpha >= -is_zero){
    // lower bound active
    if(all_alphas[i] > 0){
      result = sum[i] - all_ys[i] + epsilon_neg + lambda_eq;
    }
    else if(all_alphas[i] == 0){
      if(is_alpha_neg(i)>0){
	result = sum[i] - all_ys[i] + epsilon_neg + lambda_eq;
      }
      else{
	result = -sum[i] + all_ys[i] + epsilon_pos-lambda_eq;
      };
    }
    else{
      result = -sum[i] + all_ys[i] + epsilon_pos-lambda_eq;
    };
  }
  else if(alpha+Cpos <= is_zero){
    // upper bound active
    result = lambda_eq + sum[i] - all_ys[i] - epsilon_pos;
  };

  return result;
};


int svm_regression_c::feasible(const SVMINT i){
  // is direction i feasible to minimize the target function
  // (includes which_alpha==0)

  if(at_bound[i] >= shrink_const){ return 0; };

  SVMFLOAT alpha;
  SVMFLOAT result;

  alpha=all_alphas[i];
  //  feasible_epsilon=-1;
  if(alpha-Cneg >= - is_zero){
    // alpha* at upper bound
    result = -lambda_eq-sum[i] + all_ys[i] - epsilon_neg;
    if(result>=-feasible_epsilon){
      return 0; 
    };
  }
  else if((alpha<=is_zero) && (alpha >= -is_zero)){
    // lower bound active
    if(all_alphas[i] > 0){
      result = sum[i] - all_ys[i] + epsilon_neg + lambda_eq;
    }
    else if(all_alphas[i] == 0){
      if(is_alpha_neg(i)>0){
	result = sum[i] - all_ys[i] + epsilon_neg + lambda_eq;
      }
      else{
	result = -sum[i] + all_ys[i] + epsilon_pos - lambda_eq;
      };
    }
    else{
      result = -sum[i] + all_ys[i] + epsilon_pos - lambda_eq;
    };
    if(result>=-feasible_epsilon){
      return 0; 
    };
  }
  else if(alpha+Cpos <= is_zero){
    // alpha at upper bound
    result = lambda_eq + sum[i] - all_ys[i] - epsilon_pos;
    if(result>=-feasible_epsilon){
      return 0; 
    };
  }
  else{
    // not at bound
    result = sum[i] - all_ys[i] + lambda_eq;
    if(is_alpha_neg(i)>0){
      result = abs(result + epsilon_neg);
    }
    else{
      result = abs(result - epsilon_pos);
    };
    if(result<=feasible_epsilon){
      return 0; 
    };
  };
  return 1;
};


int svm_regression_c::feasi