8

我需要获得类似于以下链接中 Excel 函数的实现方式的线性回归斜率：

http://office.microsoft.com/en-gb/excel-help/slope-function-HP010342903.aspx

是否有 C++ 库或有人创建的简单编码解决方案可以做到这一点？

我已经根据这个公式实现了代码，但是它并不总是给我正确的结果（取自这里http://easycalculation.com/statistics/learn-regression.php）......

Slope(b) = (NΣXY - (ΣX)(ΣY)) / (NΣX2 - (ΣX)2)
         = ((5)*(1159.7)-(311)*(18.6))/((5)*(19359)-(311)2)
         = (5798.5 - 5784.6)/(96795 - 96721)
         = 13.9/74
         = 0.19

如果我对以下向量进行尝试，我会得到错误的结果（我应该期待 0.305556）： x = 6,5,11,7,5,4,4 y = 2,3,9,1,8,7 ,5

提前致谢。

4

4 回答 4

25

这是一个 C++11 实现：

#include <algorithm>
#include <iostream>
#include <numeric>
#include <vector>

double slope(const std::vector<double>& x, const std::vector<double>& y) {
    const auto n    = x.size();
    const auto s_x  = std::accumulate(x.begin(), x.end(), 0.0);
    const auto s_y  = std::accumulate(y.begin(), y.end(), 0.0);
    const auto s_xx = std::inner_product(x.begin(), x.end(), x.begin(), 0.0);
    const auto s_xy = std::inner_product(x.begin(), x.end(), y.begin(), 0.0);
    const auto a    = (n * s_xy - s_x * s_y) / (n * s_xx - s_x * s_x);
    return a;
}

int main() {
    std::vector<double> x{6, 5, 11, 7, 5, 4, 4};
    std::vector<double> y{2, 3, 9, 1, 8, 7, 5};
    std::cout << slope(x, y) << '\n';  // outputs 0.305556
}

您可以为数学要求添加一个测试（x.size() == y.size()并且x不是恒定的），或者如上面的代码，假设用户会处理这个问题。

于 2013-09-26T22:02:07.363 回答
18

你为什么不写一个这样的简单代码（当然不是最好的解决方案，只是一个基于帮助文章的例子）：

double slope(const vector<double>& x, const vector<double>& y){
    if(x.size() != y.size()){
        throw exception("...");
    }
    size_t n = x.size();

    double avgX = accumulate(x.begin(), x.end(), 0.0) / n;
    double avgY = accumulate(y.begin(), y.end(), 0.0) / n;

    double numerator = 0.0;
    double denominator = 0.0;

    for(size_t i=0; i<n; ++i){
        numerator += (x[i] - avgX) * (y[i] - avgY);
        denominator += (x[i] - avgX) * (x[i] - avgX);
    }

    if(denominator == 0.0){
        throw exception("...");
    }

    return numerator / denominator;
}

注意，accumulate 函数的第三个参数必须是 0.0 而不是 0，否则编译器会扣除它的类型 as int，accumulate 调用的结果很有可能是错误的（实际上在传递 0 时使用 MSVC2010 和 mingw-w64 是错误的作为第三个参数）。

于 2013-09-24T06:11:00.713 回答
6

以下是我用于线性回归（拟合）的模板化函数。它需要 std::vector 的数据

template <typename T>
std::vector<T> GetLinearFit(const std::vector<T>& data)
{
    T xSum = 0, ySum = 0, xxSum = 0, xySum = 0, slope, intercept;
    std::vector<T> xData;
    for (long i = 0; i < data.size(); i++)
    {
        xData.push_back(static_cast<T>(i));
    }
    for (long i = 0; i < data.size(); i++)
    {
        xSum += xData[i];
        ySum += data[i];
        xxSum += xData[i] * xData[i];
        xySum += xData[i] * data[i];
    }
    slope = (data.size() * xySum - xSum * ySum) / (data.size() * xxSum - xSum * xSum);
    intercept = (ySum - slope * xSum) / data.size();
    std::vector<T> res;
    res.push_back(slope);
    res.push_back(intercept);
    return res;
}

该函数返回一个向量，第一个元素是斜率，第二个元素是线性回归的截距。

使用它的示例：

std::vector<double> myData;
myData.push_back(1);
myData.push_back(3);
myData.push_back(4);
myData.push_back(2);
myData.push_back(5);

std::vector<double> linearReg = GetLinearFit(myData);
double slope = linearReg[0];
double intercept = linearReg[1];

请注意，该函数假定您的 x 轴有一系列数字（这是我需要的）。如果您愿意，您可以在函数中更改它。

于 2013-09-24T18:37:52.653 回答
0

我必须创建一个类似的函数，但我需要它来处理一堆近乎垂直的斜坡。我从 Cassio Neri 的代码开始，然后对其进行了修改，以在围绕线 x=y 镜像每个点后重新计算看起来比 1 更陡的斜率（这可以通过切换 x 和 y 值轻松完成）。然后它会将其镜像回来并返回更准确的斜率。

#include <algorithm>
#include <iostream>
#include <numeric>
#include <vector>

double slope(const std::vector<double>& x, const std::vector<double>& y) {

    const double n     = x.size();
    const double s_x   = std::accumulate(x.begin(), x.end(), 0.0);
    const double s_y   = std::accumulate(y.begin(), y.end(), 0.0);
    const double s_xx  = std::inner_product(x.begin(), x.end(), x.begin(), 0.0);
    const double s_xy  = std::inner_product(x.begin(), x.end(), y.begin(), 0.0);
    const double numer = n * s_xy - s_x * s_y;  // The same regardless of inversion (both terms here are commutative)
    const double denom = n * s_xx - s_x * s_x;  // Will change if inverted; use this for now
    double       a;

    if (denom == 0) a = 2;  // If slope is vertical, force variable inversion calculation
    else a = numer / denom;

    if (std::abs(a) > 1) {  // Redo with variable inversion if slope is steeper than 1
        const double s_yy = std::inner_product(y.begin(), y.end(), y.begin(), 0.0);
        const double new_denom = n * s_yy - s_y * s_y;
        a = new_denom / numer;  // Invert the fraction because we've mirrored it around x=y
    }
    return a;
}

int main() {
    std::vector<double> x{6, 5, 11, 7, 5, 4, 4};
    std::vector<double> y{2, 3, 9, 1, 8, 7, 5};
    std::cout << slope(x, y) << '\n';
}
于 2014-01-25T00:00:44.810 回答