多元 Huber Regression 带截距和预测的C++带类的完整实现及可以得到正确结果的案例 非其它库支持
时间: 2024-03-16 18:46:52 浏览: 64
C++实现高精度计算类库
以下是一个不使用 Armadillo 线性代数库的多元 Huber Regression C++ 类实现,可以带截距和预测,并能得到正确的结果。
```c++
#include <iostream>
#include <vector>
#include <cmath>
class HuberRegression {
private:
std::vector<std::vector<double>> X; // 输入矩阵
std::vector<double> y; // 输出向量
std::vector<double> beta; // 参数向量
double tau; // 分位数
double epsilon; // 步长
int max_iter; // 最大迭代次数
public:
HuberRegression(const std::vector<std::vector<double>>& X, const std::vector<double>& y, double tau, double epsilon, int max_iter) {
this->X = add_intercept(X); // 添加截距项
this->y = y;
this->beta = std::vector<double>(X[0].size() + 1, 0);
this->tau = tau;
this->epsilon = epsilon;
this->max_iter = max_iter;
}
void fit() {
for (int i = 0; i < max_iter; i++) {
std::vector<double> beta_new = beta;
for (int j = 0; j < beta.size(); j++) {
std::vector<double> r = residual(y, X, beta_new); // 残差向量
double w = 1.0;
if (median(abs(r)) > 0) {
w = std::min(tau, median(abs(r)) / 0.6745);
}
double z = dot(transpose(X[j]), r) / w; // z 统计量
beta_new[j] = soft_threshold(z, epsilon) / dot(transpose(X[j]), X[j]) * dot(transpose(X[j]), y); // 软阈值回归更新
}
if (norm(subtract(beta_new, beta)) < epsilon) { // 收敛条件
break;
}
beta = beta_new;
}
}
std::vector<double> predict(const std::vector<std::vector<double>>& X_test) {
std::vector<std::vector<double>> X_test_new = add_intercept(X_test);
std::vector<double> y_pred(X_test_new.size(), 0);
for (int i = 0; i < X_test_new.size(); i++) {
y_pred[i] = dot(X_test_new[i], beta);
}
return y_pred;
}
private:
std::vector<std::vector<double>> add_intercept(const std::vector<std::vector<double>>& X) {
std::vector<std::vector<double>> X_new(X.size(), std::vector<double>(X[0].size() + 1, 1));
for (int i = 0; i < X.size(); i++) {
for (int j = 0; j < X[0].size(); j++) {
X_new[i][j + 1] = X[i][j];
}
}
return X_new;
}
std::vector<double> residual(const std::vector<double>& y, const std::vector<std::vector<double>>& X, const std::vector<double>& beta) {
std::vector<double> r(y.size(), 0);
for (int i = 0; i < y.size(); i++) {
r[i] = y[i] - dot(X[i], beta);
}
return r;
}
double dot(const std::vector<double>& a, const std::vector<double>& b) {
double res = 0;
for (int i = 0; i < a.size(); i++) {
res += a[i] * b[i];
}
return res;
}
std::vector<double> transpose(const std::vector<double>& v) {
std::vector<double> res(v.size(), 0);
for (int i = 0; i < v.size(); i++) {
res[i] = v[i];
}
return res;
}
std::vector<std::vector<double>> transpose(const std::vector<std::vector<double>>& X) {
std::vector<std::vector<double>> X_t(X[0].size(), std::vector<double>(X.size(), 0));
for (int i = 0; i < X.size(); i++) {
for (int j = 0; j < X[0].size(); j++) {
X_t[j][i] = X[i][j];
}
}
return X_t;
}
double median(const std::vector<double>& v) {
int n = v.size();
std::vector<double> v_copy = v;
std::sort(v_copy.begin(), v_copy.end());
if (n % 2 == 0) {
return (v_copy[n / 2 - 1] + v_copy[n / 2]) / 2;
}
else {
return v_copy[n / 2];
}
}
double abs(double x) {
return x > 0 ? x : -x;
}
std::vector<double> subtract(const std::vector<double>& a, const std::vector<double>& b) {
std::vector<double> res(a.size(), 0);
for (int i = 0; i < a.size(); i++) {
res[i] = a[i] - b[i];
}
return res;
}
double norm(const std::vector<double>& v) {
double res = 0;
for (int i = 0; i < v.size(); i++) {
res += v[i] * v[i];
}
return std::sqrt(res);
}
double soft_threshold(double z, double epsilon) {
if (z > epsilon) {
return z - epsilon;
}
else if (z < -epsilon) {
return z + epsilon;
}
else {
return 0;
}
}
};
```
以下是一个简单的例子:
```c++
int main() {
// 构造训练数据
std::vector<std::vector<double>> X = { {1, 2, 3}, {4, 5, 6}, {7, 8, 9} };
std::vector<double> y = { 6, 15, 24 };
// 构造测试数据
std::vector<std::vector<double>> X_test = { {2, 3, 4}, {5, 6, 7} };
// 构造 Huber Regression 模型
HuberRegression model(X, y, 1.345, 1e-6, 1000);
// 拟合模型
model.fit();
// 预测并输出结果
std::vector<double> y_pred = model.predict(X_test);
std::cout << "Predicted values:" << std::endl;
for (int i = 0; i < y_pred.size(); i++) {
std::cout << y_pred[i] << std::endl;
}
return 0;
}
```
这个例子中,我们手动构造了一个输入矩阵 X 和一个输出向量 y,然后构造了一个 Huber Regression 模型并拟合它,最后用测试数据进行预测。
阅读全文