Аппроксимация функции, метод наименьших квадратов

uses graphabc;
const M=4;
      N=11;
var  Tochki, PrTochki, AprTochki: array[1..N] of real;
     i,j,k:integer;
     MK:array[1..n,1..n] of real;
    a:array[1..M,1..M] of real;
    b,x:array[1..M] of real;
    h:real;
    F:array[1..N,1..M] of real;
    Ft, Fu:array[1..M,1..N] of real;
 
function y(x:real):real;
begin
    y:=cos(sqr(x)+1)+exp(0.3*x);
end;
 
function fi(x:real; i:integer):real;
var j:integer;
    p:real;
begin
     p:=0;
     if i<1 then
     fi:=x
     else begin
      for j:=2 to i do
           p:=p+x*x;
      fi:=p;
     end;
end;
 
procedure NaytiTochki;
var i,j:integer;
    v:real;
begin
     for i:=1 to N do  begin
         Tochki[i]:=y(i);
         randomize;
         v:=0;
         for j:=1 to 6 do
             v:=v+Random-3;
         PrTochki[i]:=Tochki[i] + (v*sqrt(2)*Tochki[i]*0.5*0.01);
     end;
end;
 
procedure Gaus;
begin
{for i:=1 to n do begin
   writeln('Элементы матрицы:',i);
   for j:=1 to n do begin
     write('Введите A(',i,',',j,'):');
     readln(a[i,j]);
    end;
   write('Введите B(',i,'):');
   readln(b[i]);
  end;}
 
 for i:=1 to M do begin     //Прямой метол
   for j:=i+1 to M do begin
     a[j,i]:=-a[j,i]/a[i,i];
     for k:=i+1 to M do
       a[j,k]:=a[j,k]+a[j,i]*a[i,k];
     b[j]:=b[j]+a[j,i]*b[i];
    end;
  end;
 
 x[M]:=b[M]/a[M,M];      //Обратный
 for i:=M-1 downto 1 do begin
   h:=b[i];
   for j:=I+1 to M do
    h:=h-x[j]*a[i,j];
   x[i]:=h/a[i,i];
  end;
  
   {writeln(' Переменные:');        // Вывод переменных
 for i:=1 to M do
   writeln('X(',i,')=',x[i]:3:5); }
end;
 
procedure Proizvedenie;
var i,j,z:integer;
    sig:real;
begin
sig:=0.5;
for i:=1 to N do             //Заполняем матрицу k
    for j:=1 to N do
        if i = j then
           MK[i,j]:=1/sqr(sig);
 
for i:=1 to N do             //Заполняем матрицу F и Ft
    for j:=1 to M do  begin
        F[i,j]:=0;
        z:=0;
        while z<j do begin
              if z = 0 then
                 F[i,j]:=F[i,j] + Tochki[i]
              else  F[i,j]:=F[i,j]*Tochki[i];
              inc(z);
        end;
        Ft[j,i]:=F[i,j];
    end;
 
for i:=1 to M do       //Ft * K
    for j:=1 to N do begin
        Fu[i,j]:=0;
        for z:=1 to N do
            Fu[i,j]:=Fu[i,j] + Ft[i,z]*MK[z,j];
    end;
    
for i:=1 to M do       //Fu * F
    for j:=1 to M do begin
        a[i,j]:=0;
        for z:=1 to N do
            a[i,j]:=a[i,j] + Fu[i,z]*F[z,j];
    end;
    
for i:=1 to M do begin       //Fu * y
    b[i]:=0;
    for z:=1 to N do
        b[i]:=b[i] + Fu[i,z]*PrTochki[z];
    end;
end;
 
procedure Grafik;
var x, x1, y, y1,i,r : integer;
    Xn: real;
begin
Line(100,600,1250,600); //OX
Line(100,600,100,0); //OY
//Деления
Line(100,590,100,610);
TextOut(100,615, '0');
Line(200,590,200,610);
TextOut(200,615, '1');
Line(300,590,300,610);
TextOut(300,615, '2');
Line(400,590,400,610);
TextOut(400,615, '3');
Line(500,590,500,610);
TextOut(500,615, '4');
Line(600,590,600,610);
TextOut(600,615, '5');
Line(700,590,700,610);
TextOut(700,615, '6');
Line(800,590,800,610);
TextOut(800,615, '7');
Line(900,590,900,610);
TextOut(900,615, '8');
Line(1000,590,1000,610);
TextOut(1000,615, '9');
Line(1100,590,1100,610);
TextOut(1100,615, '10');
Line(1200,590,1200,610);
TextOut(1200,615, '11');
 
Line(90,100,110,100);
TextOut(75,95, '25');
Line(90,200,110,200);
TextOut(75,195, '20');
Line(90,300,110,300);
TextOut(75,295, '15');
Line(90,400,110,400);
TextOut(75,395, '10');
Line(90,500,110,500);
TextOut(75,495, '5');
 
//Точные точки
SetBrushColor(1);
SetBrushStyle(bsSolid);
x:=200;
r:=4;
for i:=1 to N do begin
    y:=600 - round(20*Tochki[i]);
    Circle(x, y, r);
    x:=x+100;
end;
 
//Зашумленные точки
SetBrushColor(clBlue);
SetBrushStyle(bsSolid);
x:=200;
for i:=1 to N do begin
    y:=600 - round(20*PrTochki[i]);
    Circle(x, y, r);
    x:=x+100;
end;
 
//Апроксимирующие точки
SetBrushColor(clRed);
SetBrushStyle(bsSolid);
x:=200;
for i:=1 to N do begin
    y:=600 - round(20*AprTochki[i]);
    Circle(x, y, r);
    x:=x+100;
end;
 
end;
 
begin
NaytiTochki; //Находим точные и зашумленные точки
Proizvedenie;
Gaus;
for i:=1 to N do begin   // Апроксимирующие точки
    AprTochki[i]:=0;
    for j:=1 to M do
        AprTochki[i]:=AprTochki[i] + x[j]*F[i,j];
end;
Grafik;
end.