多元 Huber Regression 带截距和预测的C++带类的完整实现及可以得到正确结果的案例

以下是多元 Huber Regression 的 C++ 类的实现及一个简单的案例，可以用于带截距和预测的情况。这个类使用了 Armadillo 线性代数库来实现矩阵运算。

#include <iostream>
#include <armadillo>

class HuberRegression {
private:
    arma::mat X; // 输入矩阵
    arma::colvec y; // 输出向量
    arma::colvec beta; // 参数向量
    double tau; // 分位数
    double epsilon; // 步长
    int max_iter; // 最大迭代次数
public:
    HuberRegression(const arma::mat& X, const arma::colvec& y, double tau, double epsilon, int max_iter) {
        this->X = arma::join_horiz(arma::ones<arma::mat>(X.n_rows, 1), X); // 添加截距项
        this->y = y;
        this->beta = arma::zeros<arma::colvec>(X.n_cols + 1);
        this->tau = tau;
        this->epsilon = epsilon;
        this->max_iter = max_iter;
    }

    void fit() {
        for (int i = 0; i < max_iter; i++) {
            arma::colvec beta_new = beta;
            for (int j = 0; j < beta.n_rows; j++) {
                arma::colvec r = y - X * beta_new; // 残差向量
                double w = 1.0;
                if (arma::median(arma::abs(r)) > 0) {
                    w = std::min(tau, arma::median(arma::abs(r)) / 0.6745);
                }
                double z = arma::as_scalar(X.col(j).t() * r) / w; // z 统计量
                beta_new(j, 0) = soft_threshold(z, epsilon) / (X.col(j).t() * X.col(j)) * X.col(j).t() * y; // 软阈值回归更新
            }
            if (arma::norm(beta_new - beta) < epsilon) { // 收敛条件
                break;
            }
            beta = beta_new;
        }
    }

    arma::colvec predict(const arma::mat& X_test) {
        arma::mat X_test_new = arma::join_horiz(arma::ones<arma::mat>(X_test.n_rows, 1), X_test);
        return X_test_new * beta;
    }

private:
    double soft_threshold(double z, double epsilon) {
        if (z > epsilon) {
            return z - epsilon;
        }
        else if (z < -epsilon) {
            return z + epsilon;
        }
        else {
            return 0;
        }
    }
};

接下来是一个简单的例子：

int main() {
    // 构造训练数据
    arma::mat X = arma::randn<arma::mat>(100, 3);
    arma::colvec y = X * arma::colvec({ 1, 2, 3 }) + arma::randn<arma::colvec>(100) * 0.5;

    // 构造测试数据
    arma::mat X_test = arma::randn<arma::mat>(10, 3);

    // 构造 Huber Regression 模型
    HuberRegression model(X, y, 1.345, 1e-6, 1000);

    // 拟合模型
    model.fit();

    // 预测并输出结果
    arma::colvec y_pred = model.predict(X_test);
    std::cout << "Predicted values:" << std::endl;
    std::cout << y_pred << std::endl;

    return 0;
}

这个例子中，我们生成了一个随机的输入矩阵 X 和一个真实的输出向量 y（由 X 的三列线性组合和一些噪声生成）。然后，我们构造了一个 Huber Regression 模型并拟合它，最后用测试数据进行预测。