slope 0.29.0
Loading...
Searching...
No Matches
screening.h
Go to the documentation of this file.
1
9#pragma once
10
11#include "jit_normalization.h"
12#include <Eigen/SparseCore>
13#include <memory>
14#include <vector>
15
16namespace slope {
17
24std::vector<int>
25activeSet(const Eigen::VectorXd& beta);
26
35std::vector<int>
36strongSet(const Eigen::VectorXd& gradient_prev,
37 const Eigen::ArrayXd& lambda,
38 const Eigen::ArrayXd& lambda_prev);
39
48{
49public:
53 virtual ~ScreeningRule() = default;
54
63 virtual std::vector<int> initialize(const std::vector<int>& full_set,
64 int alpha_max_ind) = 0;
65
76 virtual std::vector<int> screen(Eigen::VectorXd& gradient,
77 const Eigen::ArrayXd& lambda_curr,
78 const Eigen::ArrayXd& lambda_prev,
79 const Eigen::VectorXd& beta,
80 const std::vector<int>& full_set) = 0;
81
98 virtual bool checkKktViolations(Eigen::VectorXd& gradient,
99 const Eigen::VectorXd& beta,
100 const Eigen::ArrayXd& lambda_curr,
101 std::vector<int>& working_set,
102 const Eigen::MatrixXd& x,
103 const Eigen::MatrixXd& residual,
104 const Eigen::VectorXd& x_centers,
105 const Eigen::VectorXd& x_scales,
106 JitNormalization jit_normalization,
107 const std::vector<int>& full_set) = 0;
123 virtual bool checkKktViolations(Eigen::VectorXd& gradient,
124 const Eigen::VectorXd& beta,
125 const Eigen::ArrayXd& lambda_curr,
126 std::vector<int>& working_set,
127 const Eigen::SparseMatrix<double>& x,
128 const Eigen::MatrixXd& residual,
129 const Eigen::VectorXd& x_centers,
130 const Eigen::VectorXd& x_scales,
131 JitNormalization jit_normalization,
132 const std::vector<int>& full_set) = 0;
133
138 virtual std::string toString() const = 0;
139
140protected:
142 std::vector<int> strong_set;
143};
144
150{
151public:
152 std::vector<int> initialize(const std::vector<int>& full_set,
153 int alpha_max_ind) override;
154
155 std::vector<int> screen(Eigen::VectorXd& gradient,
156 const Eigen::ArrayXd& lambda_curr,
157 const Eigen::ArrayXd& lambda_prev,
158 const Eigen::VectorXd& beta,
159 const std::vector<int>& full_set) override;
160
161 bool checkKktViolations(Eigen::VectorXd& gradient,
162 const Eigen::VectorXd& beta,
163 const Eigen::ArrayXd& lambda_curr,
164 std::vector<int>& working_set,
165 const Eigen::MatrixXd& x,
166 const Eigen::MatrixXd& residual,
167 const Eigen::VectorXd& x_centers,
168 const Eigen::VectorXd& x_scales,
169 JitNormalization jit_normalization,
170 const std::vector<int>& full_set) override;
171
172 bool checkKktViolations(Eigen::VectorXd& gradient,
173 const Eigen::VectorXd& beta,
174 const Eigen::ArrayXd& lambda_curr,
175 std::vector<int>& working_set,
176 const Eigen::SparseMatrix<double>& x,
177 const Eigen::MatrixXd& residual,
178 const Eigen::VectorXd& x_centers,
179 const Eigen::VectorXd& x_scales,
180 JitNormalization jit_normalization,
181 const std::vector<int>& full_set) override;
182
183 std::string toString() const override;
184};
185
191{
192public:
193 std::vector<int> initialize(const std::vector<int>& full_set,
194 int alpha_max_ind) override;
195
196 std::vector<int> screen(Eigen::VectorXd& gradient,
197 const Eigen::ArrayXd& lambda_curr,
198 const Eigen::ArrayXd& lambda_prev,
199 const Eigen::VectorXd& beta,
200 const std::vector<int>& full_set) override;
201
202 bool checkKktViolations(Eigen::VectorXd& gradient,
203 const Eigen::VectorXd& beta,
204 const Eigen::ArrayXd& lambda_curr,
205 std::vector<int>& working_set,
206 const Eigen::MatrixXd& x,
207 const Eigen::MatrixXd& residual,
208 const Eigen::VectorXd& x_centers,
209 const Eigen::VectorXd& x_scales,
210 JitNormalization jit_normalization,
211 const std::vector<int>& full_set) override;
212
213 bool checkKktViolations(Eigen::VectorXd& gradient,
214 const Eigen::VectorXd& beta,
215 const Eigen::ArrayXd& lambda_curr,
216 std::vector<int>& working_set,
217 const Eigen::SparseMatrix<double>& x,
218 const Eigen::MatrixXd& residual,
219 const Eigen::VectorXd& x_centers,
220 const Eigen::VectorXd& x_scales,
221 JitNormalization jit_normalization,
222 const std::vector<int>& full_set) override;
223
224 std::string toString() const override;
225
226private:
227 template<typename MatrixType>
228 bool checkKktViolationsImpl(Eigen::VectorXd& gradient,
229 const Eigen::VectorXd& beta,
230 const Eigen::ArrayXd& lambda_curr,
231 std::vector<int>& working_set,
232 const MatrixType& x,
233 const Eigen::MatrixXd& residual,
234 const Eigen::VectorXd& x_centers,
235 const Eigen::VectorXd& x_scales,
236 JitNormalization jit_normalization,
237 const std::vector<int>& full_set);
238};
239
247std::unique_ptr<ScreeningRule>
248createScreeningRule(const std::string& screening_type);
249
250} // namespace slope
No screening rule - uses all variables.
Definition screening.h:150
std::vector< int > initialize(const std::vector< int > &full_set, int alpha_max_ind) override
Initialize the screening rule at the start of the path algorithm.
Definition screening.cpp:77
std::vector< int > screen(Eigen::VectorXd &gradient, const Eigen::ArrayXd &lambda_curr, const Eigen::ArrayXd &lambda_prev, const Eigen::VectorXd &beta, const std::vector< int > &full_set) override
Screen for the next path step.
Definition screening.cpp:83
std::string toString() const override
Get string representation of the screening rule.
bool checkKktViolations(Eigen::VectorXd &gradient, const Eigen::VectorXd &beta, const Eigen::ArrayXd &lambda_curr, std::vector< int > &working_set, const Eigen::MatrixXd &x, const Eigen::MatrixXd &residual, const Eigen::VectorXd &x_centers, const Eigen::VectorXd &x_scales, JitNormalization jit_normalization, const std::vector< int > &full_set) override
Check for KKT violations and update working set if necessary.
Definition screening.cpp:94
Base class for screening rules in SLOPE.
Definition screening.h:48
virtual bool checkKktViolations(Eigen::VectorXd &gradient, const Eigen::VectorXd &beta, const Eigen::ArrayXd &lambda_curr, std::vector< int > &working_set, const Eigen::MatrixXd &x, const Eigen::MatrixXd &residual, const Eigen::VectorXd &x_centers, const Eigen::VectorXd &x_scales, JitNormalization jit_normalization, const std::vector< int > &full_set)=0
Check for KKT violations and update working set if necessary.
virtual bool checkKktViolations(Eigen::VectorXd &gradient, const Eigen::VectorXd &beta, const Eigen::ArrayXd &lambda_curr, std::vector< int > &working_set, const Eigen::SparseMatrix< double > &x, const Eigen::MatrixXd &residual, const Eigen::VectorXd &x_centers, const Eigen::VectorXd &x_scales, JitNormalization jit_normalization, const std::vector< int > &full_set)=0
Check for KKT violations with sparse matrix input.
std::vector< int > strong_set
Strong set of variables.
Definition screening.h:142
virtual ~ScreeningRule()=default
Virtual destructor.
virtual std::vector< int > initialize(const std::vector< int > &full_set, int alpha_max_ind)=0
Initialize the screening rule at the start of the path algorithm.
virtual std::string toString() const =0
Get string representation of the screening rule.
virtual std::vector< int > screen(Eigen::VectorXd &gradient, const Eigen::ArrayXd &lambda_curr, const Eigen::ArrayXd &lambda_prev, const Eigen::VectorXd &beta, const std::vector< int > &full_set)=0
Screen for the next path step.
Implements strong screening rules for SLOPE.
Definition screening.h:191
std::vector< int > screen(Eigen::VectorXd &gradient, const Eigen::ArrayXd &lambda_curr, const Eigen::ArrayXd &lambda_prev, const Eigen::VectorXd &beta, const std::vector< int > &full_set) override
Screen for the next path step.
std::vector< int > initialize(const std::vector< int > &full_set, int alpha_max_ind) override
Initialize the screening rule at the start of the path algorithm.
std::string toString() const override
Get string representation of the screening rule.
bool checkKktViolations(Eigen::VectorXd &gradient, const Eigen::VectorXd &beta, const Eigen::ArrayXd &lambda_curr, std::vector< int > &working_set, const Eigen::MatrixXd &x, const Eigen::MatrixXd &residual, const Eigen::VectorXd &x_centers, const Eigen::VectorXd &x_scales, JitNormalization jit_normalization, const std::vector< int > &full_set) override
Check for KKT violations and update working set if necessary.
Enums to control predictor standardization behavior.
Namespace containing SLOPE regression implementation.
Definition clusters.cpp:5
std::unique_ptr< ScreeningRule > createScreeningRule(const std::string &screening_type)
Creates a screening rule based on the provided type.
JitNormalization
Enums to control predictor standardization behavior.
std::vector< int > strongSet(const Eigen::VectorXd &gradient_prev, const Eigen::ArrayXd &lambda, const Eigen::ArrayXd &lambda_prev)
Determines the strong set using sequential strong rules.
Definition screening.cpp:27
std::vector< int > activeSet(const Eigen::VectorXd &beta)
Identifies previously active variables.
Definition screening.cpp:19