slope 0.31.0
Loading...
Searching...
No Matches
slope_path.h
Go to the documentation of this file.
1
11#pragma once
12
13#include "slope_fit.h"
14#include <Eigen/Dense>
15#include <Eigen/SparseCore>
16
17namespace slope {
18
31{
32private:
33 std::vector<SlopeFit> fits;
34
35public:
39 SlopePath() = default;
40
47 SlopePath(const std::vector<SlopeFit>& fits)
48 : fits{ fits }
49 {
50 }
51
58 const SlopeFit& operator()(const size_t step) const
59 {
60 assert(step < fits.size());
61
62 return fits[step];
63 }
64
74 std::vector<Eigen::VectorXd> getIntercepts(
75 const bool original_scale = true) const
76 {
77 std::vector<Eigen::VectorXd> intercepts;
78
79 for (const auto& fit : fits) {
80 intercepts.emplace_back(fit.getIntercepts(original_scale));
81 }
82
83 return intercepts;
84 }
85
96 std::vector<Eigen::SparseMatrix<double>> getCoefs(
97 const bool original_scale = true) const
98 {
99 std::vector<Eigen::SparseMatrix<double>> coefs;
100
101 for (const auto& fit : fits) {
102 coefs.emplace_back(fit.getCoefs(original_scale));
103 }
104
105 return coefs;
106 }
107
113 std::vector<Clusters> getClusters() const
114 {
115 std::vector<Clusters> clusters;
116
117 for (const auto& fit : fits) {
118 clusters.emplace_back(fit.getClusters());
119 }
120
121 return clusters;
122 }
123
127 Eigen::ArrayXd getAlpha() const
128 {
129
130 Eigen::ArrayXd alphas(fits.size());
131
132 for (size_t i = 0; i < fits.size(); i++) {
133 alphas(i) = fits[i].getAlpha();
134 }
135
136 return alphas;
137 }
138
142 const Eigen::ArrayXd& getLambda() const { return fits.front().getLambda(); }
143
147 std::vector<double> getDeviance() const
148 {
149 std::vector<double> deviances;
150
151 for (const auto& fit : fits) {
152 deviances.emplace_back(fit.getDeviance());
153 }
154
155 return deviances;
156 }
157
161 double getNullDeviance() const { return fits.front().getNullDeviance(); }
162
166 std::vector<std::vector<double>> getPrimals() const
167 {
168 std::vector<std::vector<double>> primals;
169
170 for (const auto& fit : fits) {
171 primals.emplace_back(fit.getPrimals());
172 }
173
174 return primals;
175 }
176
180 std::vector<std::vector<double>> getDuals() const
181 {
182 std::vector<std::vector<double>> duals;
183
184 for (const auto& fit : fits) {
185 duals.emplace_back(fit.getDuals());
186 }
187
188 return duals;
189 }
190
194 std::vector<std::vector<double>> getTime() const
195 {
196 std::vector<std::vector<double>> time;
197
198 for (const auto& fit : fits) {
199 time.emplace_back(fit.getTime());
200 }
201
202 return time;
203 }
204
208 std::vector<int> getPasses() const
209 {
210 std::vector<int> passes;
211
212 for (const auto& fit : fits) {
213 passes.emplace_back(fit.getPasses());
214 }
215
216 return passes;
217 }
218
224 const std::vector<double> getDevianceRatios() const
225 {
226 std::vector<double> ratios;
227
228 for (const auto& fit : fits) {
229 ratios.emplace_back(fit.getDevianceRatio());
230 }
231
232 return ratios;
233 }
234
240 std::vector<std::vector<double>> getGaps() const
241 {
242 std::vector<std::vector<double>> gaps;
243
244 for (const auto& fit : fits) {
245 gaps.emplace_back(fit.getGaps());
246 }
247
248 return gaps;
249 }
250
255 std::size_t size() const { return fits.size(); }
256};
257
258} // namespace slope
A class representing the results of SLOPE (Sorted L1 Penalized Estimation) fitting.
Definition slope_fit.h:27
Container class for SLOPE regression solution paths.
Definition slope_path.h:31
std::vector< Eigen::VectorXd > getIntercepts(const bool original_scale=true) const
Returns the vector of intercept terms for each solution in the path.
Definition slope_path.h:74
std::vector< double > getDeviance() const
Gets the deviance values for each solution.
Definition slope_path.h:147
const SlopeFit & operator()(const size_t step) const
Get one of the coefficients along the path.
Definition slope_path.h:58
std::vector< std::vector< double > > getGaps() const
Computes the duality gaps (primal - dual objectives) for each solution.
Definition slope_path.h:240
std::vector< std::vector< double > > getDuals() const
Gets the dual objective values during optimization.
Definition slope_path.h:180
std::vector< int > getPasses() const
Gets the number of iterations for each solution.
Definition slope_path.h:208
std::vector< std::vector< double > > getPrimals() const
Gets the primal objective values during optimization.
Definition slope_path.h:166
std::vector< std::vector< double > > getTime() const
Gets the computation times for each solution.
Definition slope_path.h:194
const std::vector< double > getDevianceRatios() const
Computes the deviance ratio (explained deviance) for each solution.
Definition slope_path.h:224
SlopePath()=default
Constructs an empty SlopePath object.
const Eigen::ArrayXd & getLambda() const
Gets the lambda (regularization) weights.
Definition slope_path.h:142
Eigen::ArrayXd getAlpha() const
Gets the alpha parameter sequence.
Definition slope_path.h:127
std::size_t size() const
Gets the number of solutions in the path.
Definition slope_path.h:255
std::vector< Clusters > getClusters() const
Returns the clusters for each solution in the path.
Definition slope_path.h:113
std::vector< Eigen::SparseMatrix< double > > getCoefs(const bool original_scale=true) const
Returns the vector of coefficient matrices for each solution in the path.
Definition slope_path.h:96
SlopePath(const std::vector< SlopeFit > &fits)
Constructs a SlopePath object containing SLOPE regression solution path data.
Definition slope_path.h:47
double getNullDeviance() const
Gets the null model deviance.
Definition slope_path.h:161
Namespace containing SLOPE regression implementation.
Definition clusters.cpp:5
SLOPE (Sorted L-One Penalized Estimation) fitting results.