#include <ElasticNetModel.h>

Inheritance diagram for WireCell::ElasticNetModel:

Public Member Functions
	ElasticNetModel (double lambda=1., double alpha=1., int max_iter=100000, double TOL=1e-3, bool non_negtive=true)

	~ElasticNetModel ()

void	SetLambdaWeight (Eigen::VectorXd w)

void	SetLambdaWeight (int i, double weight)

void	SetX (Eigen::MatrixXd X)

virtual void	Fit ()

Public Member Functions inherited from WireCell::LinearModel
	LinearModel ()

virtual	~LinearModel ()

Eigen::VectorXd &	Gety ()

Eigen::MatrixXd &	GetX ()

Eigen::VectorXd &	Getbeta ()

virtual void	SetData (Eigen::MatrixXd X, Eigen::VectorXd y)

virtual void	Sety (Eigen::VectorXd y)

virtual void	Setbeta (Eigen::VectorXd beta)

Eigen::VectorXd	Predict ()

double	chi2_base ()

double	MeanResidual ()

Public Attributes
double	lambda

double	alpha

int	max_iter

double	TOL

bool	non_negtive

Eigen::VectorXd	lambda_weight

Public Attributes inherited from WireCell::LinearModel
std::string	name

Protected Member Functions
double	_soft_thresholding (double x, double lambda_)

Protected Attributes
std::vector< bool >	_active_beta

Protected Attributes inherited from WireCell::LinearModel
Eigen::VectorXd	_y

Eigen::MatrixXd	_X

Eigen::VectorXd	_beta

Detailed Description

Definition at line 10 of file ElasticNetModel.h.

Constructor & Destructor Documentation

WireCell::ElasticNetModel::ElasticNetModel	(	double	lambda = `1.`,
		double	alpha = `1.`,
		int	max_iter = `100000`,
		double	TOL = `1e-3`,
		bool	non_negtive = `true`
	)

Definition at line 17 of file ElasticNetModel.cxx.

 : lambda(lambda), alpha(alpha), max_iter(max_iter), TOL(TOL), non_negtive(non_negtive)
 {
     name = "Elastic net";
 }

WireCell::ElasticNetModel::~ElasticNetModel ( )

Definition at line 23 of file ElasticNetModel.cxx.

24 {}

Member Function Documentation

double WireCell::ElasticNetModel::_soft_thresholding	(	double	x,
		double	lambda_
	)

protected

Definition at line 98 of file ElasticNetModel.cxx.

 {
 
     if (delta > lambda_) {
         return delta - lambda_;
     }
     else {
         if (non_negtive) {
             return 0;
         }
         else {
             if (delta < -lambda_) {
                 return delta + lambda_;
             }
             else {
                 return 0;
             }
         }
     }
 }

void WireCell::ElasticNetModel::Fit ( )

virtual

Reimplemented from WireCell::LinearModel.

Reimplemented in WireCell::LassoModel.

Definition at line 26 of file ElasticNetModel.cxx.

 {
     // initialize solution to zero unless user set beta already
     Eigen::VectorXd beta = _beta;
     if (0 == beta.size()) {
         beta = VectorXd::Zero(_X.cols());
     }
 
     // initialize active_beta to true
     int nbeta = beta.size();
     _active_beta = vector<bool>(nbeta, true);
 
     // use alias for easy notation
     Eigen::VectorXd y = Gety();
     Eigen::MatrixXd X = GetX();
 
     // cooridate decsent
 
     //int N = y.size();
     VectorXd norm(nbeta);
     for (int j=0; j<nbeta; j++) {
         norm(j) = X.col(j).squaredNorm();
         if (norm(j) < 1e-6) {
             cerr << "warning: the " << j << "th variable is not used, please consider removing it." << endl;
             norm(j) = 1;
         }
     }
     double tol2 = TOL*TOL*nbeta;
 
     int double_check = 0;
     for (int i=0; i<max_iter; i++) {
         VectorXd betalast = beta;
         for (int j=0; j<nbeta; j++) {
             if (!_active_beta[j]) {continue;}
             VectorXd X_j = X.col(j);
             VectorXd beta_tmp = beta;
             beta_tmp(j) = 0;
             VectorXd r_j = (y - X * beta_tmp);
             double delta_j = X_j.dot(r_j);
             //            beta(j) = _soft_thresholding(delta_j, N*lambda*alpha*lambda_weight(j)) / (1+lambda*(1-alpha)) / norm(j);
             beta(j) = _soft_thresholding(delta_j/norm(j), lambda*alpha*lambda_weight(j)) / (1+lambda*(1-alpha));
             
             //cout << i << " " << j << " " << beta(j) << std::endl;
             if(fabs(beta(j)) < 1e-6) { _active_beta[j] = false; }
             // else { cout << beta(j) << endl;}
             // beta(j) = _soft_thresholding(delta_j, N*lambda*alpha, j) / (1+lambda*(1-alpha)) / norm(j);
             // if (j==0) cout << beta(j) << ", " << arg1 << endl;
         }
         double_check++;
         // cout << endl;
         VectorXd diff = beta - betalast;
 
         //std::cout << i << " " << diff.squaredNorm() << " " << tol2 << std::endl;
         if (diff.squaredNorm()<tol2) {
             if (double_check!=1) {
                 double_check = 0;
                 for (int k=0; k<nbeta; k++) {
                     _active_beta[k] = true;
                 }
             }
             else {
               //                cout << "found minimum at iteration: " << i << endl;
                 break;
             }
 
         }
     }
 
     // save results in the model
     Setbeta(beta);
 }