直線回帰(単回帰)

 <Xデータ> TAB <Yデータ> 改行 形式のデータセットから直線回帰(単回帰)を求めるスクリプトです。
 Y = α * X + β 式の α、βを求めます。

 "kaiki.vvs":

//----------------------------------------------------------------------
//
//          File:           "kaiki.vvs"
//          Author:         Nobuhide Tsuda
//          Created:        11-Jun-2004
//          Description:    直線回帰を計算するスクリプト
//
//----------------------------------------------------------------------

//----------------------------------------------------------------------
//
//      使用方法:
//          <Xデータ> TAB <Yデータ> 改行
//          上記形式のデータ行を選択し、スクリプト実行
//
//----------------------------------------------------------------------

//----------------------------------------------------------------------
//
//      計算式:
//          α = Σ(X[i] - ave(X)) * (Y[i] - ave(Y)) / Σ(X(i) - ave(X))^2
//          β = ave(Y) - α * ave(X)
//
//----------------------------------------------------------------------

function main()
{
    var line1 = 1;
    var line2 = thisView.getLineCount();
    if( thisView.isSelected() ) {
        var range = thisView.getSelectedRange();
        //writeln(range);
        line1 = range.line1;
        line2 = range.line2;
        if( range.offset2 == 0 )
            line2 -= 1;
    }

    var nData = 0;
    var xVect, yVect;
    var xSum = 0.0;
    var ySum = 0.0;
    for(var line = line1; line <= line2; ++line) {
        var text = thisView.getLineString(line);
        if( text =~ /(\S+)\s+(\S+)/ ) {
            var x = double($1);
            var y = double($2);
            writeln("x = " + x + ", y = ", y);
            xVect[nData] = x;
            yVect[nData] = y;
            xSum += x;
            ySum += y;
            nData += 1;
        }
    }
    if( nData == 0 ) {
        writeln("No valid data.");
        return;
    }
    var aveX = xSum / nData;
    var aveY = ySum / nData;
    writeln(format("ave(X) = %f, ave(Y) = %f", aveX, aveY));
    var num = 0.0;  //  分子
    var den = 0.0;  //  分母
    for(var i = 0; i < nData; ++i) {
        num += (xVect[i] - aveX) * (yVect[i] - aveY);
        den += (xVect[i] - aveX) * (xVect[i] - aveX);
    }
    var alpha = num / den;
    var beta = aveY - alpha * aveX;
    writeln(format("alpha = %f, beta = %f", alpha, beta));
}

 データ:

0	56
25	69
115	119
140	136

 実行結果:

x = 0, y = 56
x = 25, y = 69
x = 115, y = 119
x = 140, y = 136
ave(X) = 70.000000, ave(Y) = 95.000000
alpha = 0.566787, beta = 55.324910