Polynomial regression: Difference between revisions

=={{header|Modula-2}}==

<lang modula2>MODULE PolynomialRegression;

FROM FormatString IMPORT FormatString;

FROM RealStr IMPORT RealToStr;

FROM Terminal IMPORT WriteString,WriteLn,ReadChar;

PROCEDURE Eval(a,b,c,x : REAL) : REAL;

BEGIN

RETURN a + b*x + c*x*x;

END Eval;

PROCEDURE Regression(x,y : ARRAY OF INTEGER);

VAR

n,i : INTEGER;

xm,x2m,x3m,x4m : REAL;

ym : REAL;

xym,x2ym : REAL;

sxx,sxy,sxx2,sx2x2,sx2y : REAL;

a,b,c : REAL;

buf : ARRAY[0..63] OF CHAR;

BEGIN

n := SIZE(x)/SIZE(INTEGER);

xm := 0.0;

ym := 0.0;

x2m := 0.0;

x3m := 0.0;

x4m := 0.0;

xym := 0.0;

x2ym := 0.0;

FOR i:=0 TO n-1 DO

xm := xm + FLOAT(x[i]);

ym := ym + FLOAT(y[i]);

x2m := x2m + FLOAT(x[i]) * FLOAT(x[i]);

x3m := x3m + FLOAT(x[i]) * FLOAT(x[i]) * FLOAT(x[i]);

x4m := x4m + FLOAT(x[i]) * FLOAT(x[i]) * FLOAT(x[i]) * FLOAT(x[i]);

xym := xym + FLOAT(x[i]) * FLOAT(y[i]);

x2ym := x2ym + FLOAT(x[i]) * FLOAT(x[i]) * FLOAT(y[i]);

END;

xm := xm / FLOAT(n);

ym := ym / FLOAT(n);

x2m := x2m / FLOAT(n);

x3m := x3m / FLOAT(n);

x4m := x4m / FLOAT(n);

xym := xym / FLOAT(n);

x2ym := x2ym / FLOAT(n);

sxx := x2m - xm * xm;

sxy := xym - xm * ym;

sxx2 := x3m - xm * x2m;

sx2x2 := x4m - x2m * x2m;

sx2y := x2ym - x2m * ym;

b := (sxy * sx2x2 - sx2y * sxx2) / (sxx * sx2x2 - sxx2 * sxx2);

c := (sx2y * sxx - sxy * sxx2) / (sxx * sx2x2 - sxx2 * sxx2);

a := ym - b * xm - c * x2m;

WriteString("y = ");

RealToStr(a, buf);

WriteString(buf);

WriteString(" + ");

RealToStr(b, buf);

WriteString(buf);

WriteString("x + ");

RealToStr(c, buf);

WriteString(buf);

WriteString("x^2");

WriteLn;

FOR i:=0 TO n-1 DO

FormatString("%2i %3i ", buf, x[i], y[i]);

WriteString(buf);

RealToStr(Eval(a,b,c,FLOAT(x[i])), buf);

WriteString(buf);

WriteLn;

END;

END Regression;

TYPE R = ARRAY[0..10] OF INTEGER;

VAR

x,y : R;

BEGIN

x := R{0,1,2,3,4,5,6,7,8,9,10};

y := R{1,6,17,34,57,86,121,162,209,262,321};

Regression(x,y);

ReadChar;

END PolynomialRegression.</lang>