C++实现逻辑回归：从数学推导到高性能实现

1次阅读

没有评论

共计 2975 个字符，预计需要花费 8 分钟才能阅读完成。

在实际项目中，我们经常会遇到需要快速迭代模型或对性能有极高要求的场景。虽然 Python 的 scikit-learn 等库提供了现成的逻辑回归实现，但在以下情况下自己实现更有优势：

需要深度定制算法细节（如特殊正则化策略）
对预测延迟有严苛要求（C++ 通常比 Python 快 10-100 倍）
需要将模型集成到现有 C ++ 工程中
希望完全掌控内存管理和计算优化

逻辑回归虽然名字里有 ” 回归 ”，实际上是解决二分类问题的经典方法。其核心是 sigmoid 函数：

$$
\sigma(z) = \frac{1}{1+e^{-z}}
$$

它将线性组合 $z = w^Tx + b$ 映射到 (0,1) 区间，可以解释为样本属于正类的概率。我们使用交叉熵作为损失函数：

$$
J(w,b) = -\frac{1}{m}\sum_{i=1}^m [y^{(i)}\log h(x^{(i)}) + (1-y^{(i)})\log(1-h(x^{(i)}))]
$$

良好的数据预处理能显著提升模型性能。我们主要做两件事：

特征标准化：将各特征缩放到相近范围，加速梯度下降收敛
添加偏置项：将 b 并入 w 向量，简化计算

// 使用 Eigen 进行标准化
MatrixXd standardize(const MatrixXd& X) {RowVectorXd mean = X.colwise().mean();
    RowVectorXd std = ((X.rowwise() - mean).array().square().colwise().sum() 
                      / X.rows()).sqrt();
    return (X.rowwise() - mean).array().rowwise() / std.array();}

批量梯度下降的核心是迭代更新权重：

$$
w := w – \alpha \frac{\partial J}{\partial w}
$$

其中梯度计算公式为：

$$
\frac{\partial J}{\partial w} = \frac{1}{m}X^T(\sigma(Xw) – y)
$$

C++ 实现要点：

使用 Eigen 的矩阵运算避免显式循环
实现学习率衰减策略
添加早停机制

void LogisticRegression::fit(const MatrixXd& X, const VectorXd& y, 
                            int max_iter, double lr, double tol) {for (int i = 0; i < max_iter; ++i) {
        VectorXd z = X * weights;
        VectorXd predictions = sigmoid(z);
        VectorXd error = predictions - y;
        VectorXd gradient = X.transpose() * error / X.rows();

        weights -= lr * gradient;

        // 计算当前损失
        double loss = compute_loss(predictions, y);
        if (i % 100 == 0) {std::cout << "Iter" << i << ", loss:" << loss << std::endl;}

        // 早停检查
        if (gradient.norm() < tol) break;
    }
}

为防止过拟合，我们添加 L2 正则化（Ridge）：

$$
J_{reg}(w) = J(w) + \frac{\lambda}{2m} ||w||^2_2
$$

对应的梯度更新变为：

$$
w := w(1-\alpha\frac{\lambda}{m}) – \alpha \frac{\partial J}{\partial w}
$$

建议的项目结构：

logistic_regression/
├── include/
│   └── logistic_regression.hpp  # 类声明
├── src/
│   ├── logistic_regression.cpp  # 实现
│   └── main.cpp                # 示例用法
├── test/
│   └── test_lr.cpp             # 单元测试
└── CMakeLists.txt

关键类的声明：

class LogisticRegression {
public:
    LogisticRegression(double lambda = 0.0);

    void fit(const MatrixXd& X, const VectorXd& y, 
            int max_iter = 1000, double lr = 0.01, double tol = 1e-4);

    VectorXd predict_proba(const MatrixXd& X) const;
    VectorXi predict(const MatrixXd& X, double threshold = 0.5) const;

    void save(const std::string& filename) const;
    void load(const std::string& filename);

private:
    VectorXd weights;
    double lambda;  // 正则化系数

    VectorXd sigmoid(const VectorXd& z) const;
    double compute_loss(const VectorXd& preds, const VectorXd& y) const;
};

矩阵运算优化：
使用 Eigen 的 Map 类避免数据拷贝
利用表达式模板延迟计算
并行计算：
使用 OpenMP 并行化预测阶段
对于大型数据集，实现 mini-batch 梯度下降

#pragma omp parallel for
for (int i = 0; i < X.rows(); ++i) {predictions[i] = sigmoid(X.row(i).dot(weights));
}

内存优化：
重用中间矩阵
对稀疏数据使用 SparseMatrix

数值溢出问题：
在计算 sigmoid 时，对极大 / 极小值做截断

VectorXd sigmoid(const VectorXd& z) {return z.unaryExpr([](double x) {if (x > 20) return 1.0;
        if (x < -20) return 0.0;
        return 1.0 / (1.0 + exp(-x));
    });
}

学习率选择：
实现学习率衰减：lr = initial_lr / (1 + decay_rate * iter)
使用线搜索自动确定最佳学习率
特征相关性：
添加 L1 正则化（LASSO）自动特征选择
使用 PCA 降维

模型序列化：
使用 Eigen 的原生格式或 Protocol Buffers

void LogisticRegression::save(const std::string& filename) const {std::ofstream ofs(filename, std::ios::binary);
    Eigen::MatrixXd w = weights;
    Eigen::write_binary(ofs, w);
}