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

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;
}

累積関数の 3 番目の引数は 0 ではなく 0.0 でなければならないことに注意してください。そうしないと、コンパイラはその型を as としてint推論し、累積呼び出しの結果が間違っている可能性が高くなります (0 を渡すときに MSVC2010 と mingw-w64 を使用すると実際には間違っています)。 3 番目のパラメーターとして)。

于 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;
}

この関数は、最初の要素が勾配で、2 番目の要素が線形回帰の切片であるベクトルを返します。

使用例:

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 に答える