最小二乗法でy=ax+bの回帰直線を求める ~Delphiソースコード集
最小二乗法のユニット
ファイル名「ULeastSquares.pas」で保存してください。
unit ULeastSquares;
{**************************************************************
最小2乗法 ユニット
y = ax + b
で、x,yの配列(TLeastSquaresValues)を与えると、
傾きTLeastSquaresAns.aと切片TLeastSquaresAns.bを返す
**************************************************************}
interface
type
TLeastSquaresValue=record
x,y:Extended;
end;
TLeastSquaresValues=TArray<TLeastSquaresValue>;
TLeastSquaresAns=record
a,b:Extended;
end;
TLeastSquares=class
public
class function Calc(XYArr:TLeastSquaresValues):TLeastSquaresAns;
private
class function SumXY(XYArr:TLeastSquaresValues):Extended;
class function SumX(XYArr:TLeastSquaresValues):Extended;
class function SumY(XYArr:TLeastSquaresValues):Extended;
class function SumX2(XYArr:TLeastSquaresValues):Extended;
end;
implementation
class function TLeastSquares.Calc(XYArr: TLeastSquaresValues): TLeastSquaresAns;
var x,y,xy,x2:Extended;
e:Extended;
n:Extended;
begin
n:=Length(XYArr);//xyデータの数
x:=SumX(XYArr); //Σx
y:=SumY(XYArr); //Σy
xy:=SumXY(XYArr);//Σx*y
x2:=SumX2(XYArr);//Σx*x
e:=(n*x2-x*x);
if e<>0 then
begin
result.a:=(n*xy-x*y)/e;
result.b:=(x2*y-xy*x)/e;
end
else
begin
result.a:=0;
result.b:=0;
end;
end;
class function TLeastSquares.SumX(XYArr: TLeastSquaresValues): Extended;
var i:Integer;
begin
//Σx
result:=0;
for i := Low(XYArr) to High(XYArr) do
result:=result+XYArr[i].x;
end;
class function TLeastSquares.SumY(XYArr: TLeastSquaresValues): Extended;
var i:Integer;
begin
//Σy
result:=0;
for i := Low(XYArr) to High(XYArr) do
result:=result+XYArr[i].y;
end;
class function TLeastSquares.SumXY(XYArr: TLeastSquaresValues): Extended;
var i:Integer;
begin
//Σx*y
result:=0;
for i := Low(XYArr) to High(XYArr) do
result:=result+XYArr[i].x*XYArr[i].y;
end;
class function TLeastSquares.SumX2(XYArr: TLeastSquaresValues): Extended;
var i:Integer;
begin
//Σx*x
result:=0;
for i := Low(XYArr) to High(XYArr) do
result:=result+XYArr[i].x*XYArr[i].x;
end;
end.
使用例
ファイル⇒新規作成⇒VCLフォームアプリケーションをクリックします(新規プロジェクトが作成されます)。フォームにTButton×1個、TMemo×1個をドラッグ&ドロップします。
Button1をダウブルクリックして以下ソースコードを記述します。
unit Unit1;
interface
uses
Winapi.Windows, Winapi.Messages, System.SysUtils, System.Variants, System.Classes,
Vcl.Graphics, Vcl.Controls, Vcl.Forms, Vcl.Dialogs, Vcl.StdCtrls;
type
TForm1 = class(TForm)
Button1: TButton;
Memo1: TMemo;
procedure Button1Click(Sender: TObject);
private
{ Private 宣言 }
public
{ Public 宣言 }
end;
var
Form1: TForm1;
implementation
{$R *.dfm}
uses ULeastSquares;
procedure TForm1.Button1Click(Sender: TObject);
var xy:TLeastSquaresValues;
ans:TLeastSquaresAns;
begin
SetLength(xy,10);
xy[0].x:= 0.6094; xy[0].y:= 0.8393;
xy[1].x:=-0.2344; xy[1].y:= 0.0354;
xy[2].x:=-0.6029; xy[2].y:=-0.1950;
xy[3].x:= 0.1168; xy[3].y:= 0.2681;
xy[4].x:=-0.7346; xy[4].y:=-0.8405;
xy[5].x:=-0.8232; xy[5].y:= 0.4612;
xy[6].x:= 0.6606; xy[6].y:=-0.3923;
xy[7].x:= 0.926 ; xy[7].y:= 0.9928;
xy[8].x:=-0.5316; xy[8].y:= 0.5900;
xy[9].x:=-0.8263; xy[9].y:=-0.0484;
ans:=TLeastSquares.Calc(xy);
Memo1.Lines.Add(
Format('y = %8.4fx+%8.4f',[ans.a,ans.b])
);
end;
end.
実行する
実行して、Button1を押すと、最小二乗法で傾きと切片が計算されます。
参考:エクセルでの近似直線結果
