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

Содержание

 

Введение

 

Аппроксимация, или приближение —математический метод, состоящий в замене одних математических объектов другими, в том или ином смысле близкими к исходным, но более простыми. В нашем случае, аппроксимация заключается в том, что используя имеющуюся информацию по f(x) можно рассмотреть другую функцию φ(y) близкую в некотором смысле к f(x), позволяющую выполнить над ней соответствующие операции и получить оценку  погрешность такой замены.

       Аппроксимация позволяет исследовать числовые характеристики и качественные свойства объекта, сводя задачу к изучению более простых или более удобных объектов (например, таких, характеристики которых легко вычисляются, или свойства которых уже известны).

        Как известно, между величинами может существовать точная (функциональная) связь, когда одному значению аргумента соответствует одно определенное значение, и менее точная (корреляционная) связь, когда одному конкретному значению аргумента соответствует приближенное значение или некоторое множество значений функции, в той или иной степени близких друг к другу.

      При   выборе  аппроксимации  следует   исходить  из  конкретной   задачи  исследования. Обычно,  чем  более  простое   уравнение   используется   для аппроксимации, тем более  приблизительно  получаемое  описание  зависимости. Поэтому важно считывать, насколько существенны и чем обусловлены отклонения конкретных  значений  от  получаемого  тренда. Также очень важно уловить общую закономерность,  которая в данном случае  наиболее  логично и с приемлемой  точностью   выражается   именно двухпараметрическим уравнением степенной функции. 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

1. Постановка задачи

 

В результате эксперимента получена зависимость, представленная в таблице:

x	1	1.5	2	2.5	3	3.5	4	4.5	5
Y	6.5	20.38	46.4	88.63	151.1	237.9	535	500.3	648.5

 

где величина Y является функцией аргумента x.

Ax³+Cx+B

Необходимо:

-функцию, заданную  таблично, аппроксимировать, используя метод наименьших квадратов;

- построить линию тренда;

- вычислить числовые  характеристики функции Y(x);

- написать программу  на языке Паскаль;

-сравнить результаты  работы программы на языке  Паскаль и их соответствии  контрольному просчету.

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

2. Описание математической модели  решения задачи

При анализе эмпирических данных возникает необходимость  найти в явном виде функциональную зависимость между величинами х  и Y, которые получены в результате измерений.

Метод наименьших квадратов - метод  оценки параметров модели на основании  экспериментальных данных, содержащих случайные ошибки.

В основе метода лежат следующие  рассуждения: при замене точного (неизвестного) параметра модели приблизительным  значением необходимо минимизировать разницу между экспериментальными данными и теоретическими. Это позволяет рассчитать параметры модели с помощью метода наименьших квадратов с минимальной погрешностью.

Мерой разницы в методе наименьших квадратов служит сумма квадратов отклонений действительных (экспериментальных) значений от теоретических. Выбираются такие значения параметров модели, при которых сумма квадратов разностей будет наименьшей.

 

Необходимо определить параметры  функции Y = A*x³+C*x+B.

Составим функцию: ∑ (Yp -Yt ) = ∑ (Yp - A*x ³-C* x -B)    =>min

 

F' =-2∑ (Yp - A* x -C* x -B)* x =-2∑ (Yp * x -A* x -C* x -B* x ) = 0;

F'c=-3∑ (Yp - A* x ³-C* x - B)* x³ =-3∑ (Yp * x ³-A* x ⁶-C* x ⁴-B* x ³) = 0;

F'b=-3∑ (Yp - A* x ³-C* x -B)* x=-3∑ (Yp *x - A* x ⁴-C* x ²-B) = 0;

Запишем систему уравнений:

        n*B+C*∑ x +A*∑ x ³= Yp ;

        B*∑ x +C*∑ x +A*∑ x ⁴=∑ Yp * x ;

        B*∑ x ³+C*∑ x ⁴+A*∑ x ⁶=∑ Yp * x ³.

Для нахождения параметров А, В, С необходимо решить эту систему линейных алгебраических уравнений.

Имеем две матрицы:

- матрица коэффициентов при  неизвестных

          n        ∑ x      ∑ x ³

А =    ∑ x      ∑ x    ∑ x ⁴

          ∑ x ³   ∑ x ⁴     ∑ x ⁶

• матрица свободных коэффициентов

                  ∑QY

        В =    ∑ Yp * x

                  ∑ Yp * x ³

 

Для нахождения коэффициентов  решим эту систему уравнений использую метод Крамера для решения матриц.

 

Составим главный определитель d1:

 

 

 

d=   x⁶     x⁴       x³

 x⁴       x²      x

      x³      x       n

Составим d1,d2,d3 используя матрицу В:

 

d1=   y x³  x⁴       x³

                 yx    x²      x

         y       x      n

 

d2=  x⁶     y x³   x³

        x⁴     yx    x

        x³        y        n

 

d3=  x⁶    x⁴         y x³   

        x⁴       x²         yx

        x³       x       y

 

Затем находим коэффициенты: формулам:

А=d1/d;

C=d2/d;

B=d3/d;

 

После чего можем записать уравнение: Y(x)=A*x³+C*x+B  с найденными коэффициентами.

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

3.Блок-схема

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

4.Описание алгоритма

 

Блок-схема реализует  циклический процесс вычисления функции и содержит следующие блоки:

1. Начало
2. Задание переменных S_x, S_xY, S_x2,S_x3, S_x4,S_x6, S_x3Y,S_Y значение 0;

3-4. Формирование цикла

5-6. Вычисление сумм S_x, S_xY, S_x2,S_x3, S_x4,S_x6, S_x3Y,S_Y

7-10.Вычисление матриц  Δ, Δ1,Δ2,Δ3

11-13. Вычисление коэффициентов A,C,B

14-15. Нахождение теоретического значения Yt

16.Вывод коэффициентов

17. Конец алгоритма

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

5. Характеристика данных и их  условные обозначения

№ п/п	Наименование данных	Обозначение в блок-схеме	Обозначение в програм-ме	Тип Переменной
1	Число аппроксимации, количество элементов в массиве	n	n	целочисленная
2	Счетчик цикла	i	i	целочисленная
3	Файловая переменная для исходного файла	ish	File_ish	текстовая
4	Файловая переменная для файла с результатами	rez	File_rez	текстовая
5	Аргумент функции	x	x[i]	matr
6	Ордината функции (практическое значение)	y	yp[i]	matr
7	Ордината функции (теоретическое  значение)	yti	yt[i]	matr
8	Суммы	S_x	S_x	действительная
9		S_y	S_yp	действительная
10		S_x2	S_x2	действительная
11		S_x3	S_x3	действительная
12		S_x4	S_x4	действительная
13		S_xy	S_xyp	действительная
14		S_x3y	S_x2yp	действительная
15	Элементы матрицы Δ	Δ1	d1	matr
16		Δ2	d2	matr
17		Δ3	d3	matr
18	Коэффициенты функции	А	At	действительная
19		В	Bt	действительная
20		С	Ct	действительная

 

 

 

 

 

 

 

 

 

 

 

6. Текст программы

unit Unit1;

interface

uses

  Windows, Messages, SysUtils, Classes, Graphics, Controls, Forms, Dialogs,

  StdCtrls;

type

  TForm1 = class(TForm)

    GroupBox1: TGroupBox;

    Label1: TLabel;

    Label2: TLabel;

    Memo_ish: TMemo;

    B_open: TButton;

    GroupBox2: TGroupBox;

    Memo_rez: TMemo;

    B_View: TButton;

    B_Save: TButton;

    OpenDialog1: TOpenDialog;

    SaveDialog1: TSaveDialog;

    OpenDialog2: TOpenDialog;

    B_close: TButton;

    B_clear: TButton;

 

 

    procedure B_openClick(Sender: TObject);

    procedure B_SaveClick(Sender: TObject);

    procedure B_viewClick(Sender: TObject);

    procedure B_closeClick(Sender: TObject);

    procedure B_clearClick(Sender: TObject);

    procedure B_calcClick(Sender: TObject);

  private

  public

  end;

 

var

  Form1: TForm1;

 

 

implementation

 Const n=9;

 

    Var

 

       x,yp:array [1..10] of real;

        S_x3y,S_x3,S_x6,S_y,S_x,S_x2,S_x4,S_xy: real;

       d,d1,d2,d3:real;

        FName_ish:string;

        FName_rez:string;

        File_ish,File_rez:TextFile;

        A,C,B:real;

 

 

procedure TForm1.B_openClick(Sender: TObject);

 

    Var i:integer;

  begin

  if OpenDialog1.Execute then

     begin

     FName_ish:=OpenDialog1.FileName;

     Memo_ish.Lines.LoadFromFile(FName_ish);

     AssignFile(File_ish,FName_ish);

     Reset(File_ish);

      for i:=1 to n do

         read(File_ish,x[i],Yp[i]);

         readln(File_ish);

      end;

 

 

  end;

procedure TForm1.B_SaveClick(Sender: TObject);

   Var i:integer;

 begin

   if SaveDialog1.Execute then

      begin

       FName_rez:=SaveDialog1.FileName;

           AssignFile(File_rez,FName_rez);

           Rewrite(File_rez);

           writeln(File_rez,'Òåîðåòè÷åñêèå çíà÷åíèÿ y :');

           for i:=1 to n do

           writeln (File_rez,Yp[i]:5:3);

           writeln(File_rez);

           writeln (File_rez,'========================');

 

 

           writeln(File_rez);

           writeln(File_rez,'========================');

           writeln(File_rez,'Èñêîìîå óðàâíåíèå:');

           writeln(File_rez,'y='+FormatFloat('0.00', A)+' x*x*x +'

           +FormatFloat('0.00',C)+'x+'+FormatFloat('0.00',B));

           Memo_rez.Lines.Add('y='+FormatFloat('0.00',A)

          +' x*x*x+ '+FormatFloat('0.00',C)+'x+'+FormatFloat('0.00',B));

           CloseFile(File_rez);

        end;

  end;

 

 

procedure TForm1.B_closeClick(Sender: TObject);

begin

  close;

end;

 

procedure TForm1.B_clearClick(Sender: TObject);

  begin

  Memo_ish.clear;

  Memo_rez.clear;

  end;

procedure TForm1.B_ViewClick(Sender: TObject);

   begin

     if OpenDialog2.Execute then

       begin

         FName_ish:=OpenDialog2.FileName;

         Memo_rez.Lines.LoadFromFile(FName_ish);

        end;

    end;

 

procedure TForm1.B_calcClick(Sender: TObject);

var i: integer;

begin

 

S_x:=0;

S_x3y:=0;

S_x3:=0;

S_x6:=0;

S_y:=0;

S_x2:=0;

S_x4:=0;

S_xy:=0;

for i:=1 to n do

 

  begin

S_x3y:=S_x3y+x[i]*x[i]*x[i]*yp[i];

S_x3:=S_x3+x[i]*x[i]*x[i];

S_x6:=S_x6+x[i]*x[i]*x[i]*x[i]*x[i]*x[i];

S_y:=S_y+yp[i];

S_x:=S_x+x[i];

S_xy:=S_xy+x[i]*yp[i];

S_x2:=S_x2+x[i]*x[i];

S_x4:=S_x4+x[i]*x[i]*x[i]*x[i];

  end;

 

d:=S_x6*(S_x2*n-S_x*S_x)-S_x4*(S_x4*n-S_x3*S_x)+S_x3*(S_x4*S_x-S_x3*S_x2);

d1:=S_x3y*(S_x2*n-S_x*S_x)-S_x4*(S_xy*n-S_y*S_x)+S_x3*(S_xy*S_x-S_y*S_x2);

d2:=S_x6*(S_xy*n-S_y*S_x)-S_x3y*(S_x4*n-S_x3*S_x)+S_x3*(S_x4*S_y-S_x3*S_xy);

d3:=S_x6*(S_x2*S_y-S_x*S_xy)-S_x4*(S_x4*S_y-S_x3*S_xy)+S_x3y*(S_x4*S_x-S_x3*S_x2);

A:=d1/d;

C:=d2/d;

B:=d3/d;

 

 

  for i:=1 to n do

   begin

yp[i]:=A*x[i]*x[i]*x[i]+C*x[i]+B;

 

 

   end;

end;

 

end.

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

7.Контрольный пример

n	x	y	x^2	x^3	x^4	x^6	y*x	y*x^3
9	1	6,5	1	1	1	1	6,5	6,5
	1,5	20,38	2,25	3,375	5,0625	11,390625	30,57	68,7825
	2	46,4	4	8	16	64	92,8	371,2
	2,5	88,63	6,25	15,625	39,0625	244,140625	221,575	1384,84375
	3	151,1	9	27	81	729	453,3	4079,7
	3,5	237,9	12,25	42,875	150,063	1838,265625	832,65	10199,9625
	4	535	16	64	256	4096	2140	34240
	4,5	500,3	20,25	91,125	410,063	8303,765625	2251,35	45589,8375
	5	684,5	25	125	625	15625	3422,5	85562,5
сумма	27	2270,71	96	378	1583,25	30912,5625	9451,245	181503,326
Δ=	30912,5625	1583,25	378	213505,875		A=ΔA/Δ=	4,484287033
	1583,25	96	27
	378	27	9
ΔA=	181503,3263	1583,25	378	957421,6267		C=ΔC/Δ=	41,63632524
	9451,245	96	27
	2270,71	27	9
ΔC=	30912,56	181503,326	378	8889600,052		B=ΔB/Δ=	-60,94831882
	1583,25	9451,245	27
	378	2270,71	9
ΔB=	30912,56	1583,25	181503,326	-13012824,14
	1583,25	96	9451,245
	378	27	2270,71
Yэкспер
-14,82771
16,640638
58,198628
113,20948
185,03641
277,04263
392,59135
535,0458
707,76919

8.Анализ  результатов

Выполнив курсовую работу, мы нашли параметры функции Y= Ax³+Cx+B .

Искомое уравнение имеет  вид:

Y= 4.48x³+41.63x-60.94

Получили результаты расчетов программы на языке Паскаль , используя  программу Delphi, а также выполнили контрольный просчет с помощью программы Microsoft Excel. Результаты этих просчетов совпадают, что подтверждает правильность составления алгоритма и программы. Также мы можем сделать вывод о соответствии линии Тренда с экспериментальными данными, которые были получены с помощью метода наименьших квадратов.

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Заключение

В курсовой работе была выполнена поставленная задача: составлен алгоритм решения, приведены результаты аппроксимации заданной функции методом наименьших квадратов, полученные с использованием программы MicrosoftExcel,а также приведена программа решения, написанной на языке Delphi. Вышеприведенные способы аппроксимации таблично заданной функции являются взаимно проверяющими. Преимуществом использования MSExcel является простота применения этой программы. В процессе выполнения поставленной задачи, мы определили искомые коэффициенты:

А=4,484287033

 

В=-60,94831882

 

С=41,63632524

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

Приложения

 

            

 

 

 

           

 

 

 

 

 

 

Файл  с исходными значениями

 

1          6.5

1.5        20.38

2          46.4

2.5        88.63

3          151.1

3.5        237.9

4          535

4.5        500.3

5          684.5

 

 

 

Файл с результатами вычислений

 

 

Теоретические значения y :

-14.828

16.641

58.199

113.210

185.037

277.043

392.592

535.047

707.770

 

========================

 

========================

Искомое уравнение:

y=4,48 x*x*x +41,64x+-60,95