13#include <Eigen/SparseCore>
29 Eigen::VectorXd intercepts;
30 Eigen::SparseMatrix<double> coefs;
33 Eigen::ArrayXd lambda;
40 std::vector<double> time;
42 std::string centering_type;
43 std::string scaling_type;
44 std::string loss_type;
46 Eigen::VectorXd x_centers;
47 Eigen::VectorXd x_scales;
74 const Eigen::SparseMatrix<double>& coefs,
77 const Eigen::ArrayXd& lambda,
78 const double deviance,
79 const double null_deviance,
80 const std::vector<double>& primals,
81 const std::vector<double>& duals,
82 const std::vector<double>& time,
84 const std::string& centering_type,
85 const std::string& scaling_type,
86 const bool has_intercept,
87 const Eigen::VectorXd& x_centers,
88 const Eigen::VectorXd& x_scales)
89 : intercepts{ intercepts }
91 , clusters{ clusters }
94 , deviance{ deviance }
95 , null_deviance{ null_deviance }
100 , centering_type{ centering_type }
101 , scaling_type{ scaling_type }
102 , has_intercept{ has_intercept }
103 , x_centers{ x_centers }
104 , x_scales{ x_scales }
116 if (original_scale) {
121 this->has_intercept);
133 Eigen::SparseMatrix<double>
getCoefs(
const bool original_scale =
true)
const
135 if (original_scale) {
141 this->has_intercept);
142 return beta_out.sparseView();
156 const Eigen::ArrayXd&
getLambda()
const {
return lambda; }
181 const std::vector<double>&
getPrimals()
const {
return primals; }
186 const std::vector<double>&
getDuals()
const {
return duals; }
191 const std::vector<double>&
getTime()
const {
return time; }
214 std::vector<double> gaps(primals.size());
215 for (
size_t i = 0; i < primals.size(); i++) {
216 gaps[i] = primals[i] - duals[i];
237 Eigen::MatrixXd
predict(T& x,
const std::string& type =
"response")
const
241 Eigen::MatrixXd eta = x *
getCoefs();
247 if (type ==
"linear") {
252 std::unique_ptr<Loss> loss =
setupLoss(this->loss_type);
254 return loss->predict(eta);
Representation of the clusters in SLOPE.
A class representing the results of SLOPE (Sorted L1 Penalized Estimation) fitting.
const Clusters & getClusters() const
Gets the clusters.
bool hasIntercept() const
Checks if model has intercept.
Eigen::VectorXd getIntercepts(const bool original_scale=true) const
Gets the intercept terms for this SLOPE fit.
const Eigen::ArrayXd & getLambda() const
Gets the lambda (regularization) parameter used.
SlopeFit(const Eigen::VectorXd &intercepts, const Eigen::SparseMatrix< double > &coefs, const Clusters &clusters, const double alpha, const Eigen::ArrayXd &lambda, const double deviance, const double null_deviance, const std::vector< double > &primals, const std::vector< double > &duals, const std::vector< double > &time, const int passes, const std::string ¢ering_type, const std::string &scaling_type, const bool has_intercept, const Eigen::VectorXd &x_centers, const Eigen::VectorXd &x_scales)
Construct a new Slope Fit object.
double getAlpha() const
Gets the alpha (mixing) parameter used.
Eigen::SparseMatrix< double > getCoefs(const bool original_scale=true) const
Gets the sparse coefficient matrix for this fit.
const std::vector< double > & getPrimals() const
Gets the sequence of primal objective values during optimization.
const std::vector< double > & getDuals() const
Gets the sequence of dual objective values during optimization.
const std::vector< double > & getTime() const
Gets the sequence of computation times during optimization.
std::vector< double > getGaps() const
Calculate the duality gaps during optimization.
Eigen::MatrixXd predict(T &x, const std::string &type="response") const
Predict the response for a given input matrix.
const std::string & getLossType() const
Gets the model's loss type.
double getDevianceRatio() const
Calculate the deviance ratio (1 - deviance/null_deviance)
double getNullDeviance() const
Gets the null model deviance.
double getDeviance() const
Gets the model deviance.
int getPasses() const
Gets the total number of optimization iterations.
The declaration of the Clusters class.
Namespace containing SLOPE regression implementation.
std::unique_ptr< Loss > setupLoss(const std::string &loss)
Factory function to create the appropriate loss function based on the distribution family.
void validateOption(const std::string &value, const std::set< std::string > &valid_options, const std::string ¶meter_name)
Validates if a given value exists in a set of valid options.
std::tuple< Eigen::VectorXd, Eigen::MatrixXd > rescaleCoefficients(const Eigen::VectorXd &beta0, const Eigen::SparseMatrix< double > &beta, const Eigen::VectorXd &x_centers, const Eigen::VectorXd &x_scales, const bool intercept)
Rescales the coefficients using the given parameters.
Functions to normalize the design matrix and rescale coefficients in case the design was normalized.
Factory function to create the appropriate loss function based on.
Various utility functions.