slope 0.29.0
Loading...
Searching...
No Matches
loss.h
Go to the documentation of this file.
1
7#pragma once
8
9#include <Eigen/Core>
10
11namespace slope {
12
20class Loss
21{
22public:
26 virtual ~Loss() = default;
27
38 virtual double loss(const Eigen::MatrixXd& eta, const Eigen::MatrixXd& y) = 0;
39
51 virtual double dual(const Eigen::MatrixXd& theta,
52 const Eigen::MatrixXd& y,
53 const Eigen::VectorXd& w) = 0;
54
65 virtual Eigen::MatrixXd residual(const Eigen::MatrixXd& eta,
66 const Eigen::MatrixXd& y) = 0;
67
79 virtual void updateWeightsAndWorkingResponse(Eigen::VectorXd& w,
80 Eigen::VectorXd& z,
81 const Eigen::VectorXd& eta,
82 const Eigen::VectorXd& y) = 0;
83
89 virtual Eigen::MatrixXd preprocessResponse(const Eigen::MatrixXd& y) = 0;
90
99 virtual void updateIntercept(Eigen::VectorXd& beta0,
100 const Eigen::MatrixXd& eta,
101 const Eigen::MatrixXd& y)
102 {
103 Eigen::MatrixXd residual = this->residual(eta, y);
104 beta0 -= residual.colwise().mean() / this->lipschitz_constant;
105 };
106
113 virtual Eigen::MatrixXd link(const Eigen::MatrixXd& mu) = 0;
114
121 virtual Eigen::MatrixXd inverseLink(const Eigen::MatrixXd& eta) = 0;
122
128 virtual Eigen::MatrixXd predict(const Eigen::MatrixXd& eta) = 0;
129
137 virtual double deviance(const Eigen::MatrixXd& eta, const Eigen::MatrixXd& y)
138 {
139 return 2.0 * (loss(eta, y) - loss(link(y), y));
140 }
141
142protected:
150 explicit Loss(double lipschitz_constant)
151 : lipschitz_constant(lipschitz_constant)
152 {
153 }
154
155private:
156 const double lipschitz_constant;
157};
158
159} // namespace slope
Loss function interface.
Definition loss.h:21
Loss(double lipschitz_constant)
Constructs an loss function with specified Lipschitz constant.
Definition loss.h:150
virtual double loss(const Eigen::MatrixXd &eta, const Eigen::MatrixXd &y)=0
Calculates the loss function.
virtual Eigen::MatrixXd preprocessResponse(const Eigen::MatrixXd &y)=0
Preprocess response.
virtual double deviance(const Eigen::MatrixXd &eta, const Eigen::MatrixXd &y)
Computes deviance, which is 2 times the difference between the loglikelihood of the model and the log...
Definition loss.h:137
virtual void updateIntercept(Eigen::VectorXd &beta0, const Eigen::MatrixXd &eta, const Eigen::MatrixXd &y)
Updates the intercept with a gradient descent update. Also updates the linear predictor (but not the ...
Definition loss.h:99
virtual void updateWeightsAndWorkingResponse(Eigen::VectorXd &w, Eigen::VectorXd &z, const Eigen::VectorXd &eta, const Eigen::VectorXd &y)=0
Updates the weights and working response.
virtual Eigen::MatrixXd inverseLink(const Eigen::MatrixXd &eta)=0
The inverse link function, also known as the mean function.
virtual Eigen::MatrixXd predict(const Eigen::MatrixXd &eta)=0
Return predicted response.
virtual double dual(const Eigen::MatrixXd &theta, const Eigen::MatrixXd &y, const Eigen::VectorXd &w)=0
Calculates the dual loss.
virtual Eigen::MatrixXd residual(const Eigen::MatrixXd &eta, const Eigen::MatrixXd &y)=0
Calculates the residual.
virtual ~Loss()=default
Destructor for the Loss class.
virtual Eigen::MatrixXd link(const Eigen::MatrixXd &mu)=0
The link function.
Namespace containing SLOPE regression implementation.
Definition clusters.cpp:5