Метод наименьших квадратов

Размещено на http://www.allbest.ru/

Московский Авиационный Институт

(государственный технический университет)

Курсовая работа

по теории вероятностей и математической статистике на тему

Метод наименьших квадратов

Выполнила:

Студентка гр. 05-212

Зятькова Н.Е.

Проверил:

Серонин А.Н.

Москва, 2003

Содержание

Исходные данные

Этап 1. Моделирование измерений

Этап 2. Оценка МНК

Этап 3. Доверительные интервалы

Список литературы

Исходные данные

квадрат матричный интервал доверительный

Вариант: К=4;

N=41;

=K*10^-3=0.004;

?=K*10^-2=0.04;

b1=[K/2]=2;

b2=[K/3]=4/3;

b3=[K/5]=4/5;

=0.1

Этап 1. Моделирование измерений

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

МНК допускает простую геометрическую интерпретацию, так как напрямую связан с проектированием в конечном евклидовом пространстве на некоторое его подпространство. евклидово пространств - это n-мерное векторное пространство, в котором определено скалярное произведение, элементы которого u=(u₁,…,u_n)^T и v+(v₁,…,v_n)^T складываются и умножаются на действительные числа по обычным законам, естественным правилам, а скалярное произведение задается соотношением (u,v)=u₁v₂+u₂v₂+…+u_nv_n. Скалярное произведение векторов a и b - число (скаляр) (a,b), равное произведению длин этих векторов на косинус угла между ними, т. е. (a,b) = a*b*cos .

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

Модель измерений. Предположим, что переменная y является функцией переменных х₁,х₂,…,х_n: y=f(x₁,x₂,…,x_n).

Данное соотношение часто сравнивается с некоторым техническим объектом, для которого переменные x₁,x₂,…,x_n (в данной курсовой работе - 1, t, t²) являются входными параметрами, y - выходной параметр, f - неизвестная характеристика объекта. Но можно осуществить некоторое количество опытов (41), позволяющих по заданным значениям входных параметров (1, t, t²) получить соответствующее значение выходного параметра. рассматривается полиномиальная аппроксимация (т.е. замена одних объектов другими, более простыми), т. е. Функция вида y = b₁+b₂*t+b₃*t².

В данной части курсовой работы функция f известна - y_i = b₁+b₂*t_i+b₃*t_i, и есть независимые ошибки _i. Другими словами, в первом этапе проведение опыта сводится к расчету y_i, который представляет собой следующую стохастическую (случайную, вероятную) функцию: y_i = b₁+b₂*t_i+b₃*t_i²+_i.

Рассчитываем и заполняем таблицу 1.

Таблица 1. Расчет параметра y_i

I	ti	оi=у*оi*	b1+b2ti+b3ti2	yi=b1+b2ti+b3ti2+оi	оi*
1	-0,8	-0,0012	1,4453333	1,444132405	-0,30023
2	-0,76	-0,00511	1,4487467	1,443635934	-1,27768
3	-0,72	0,00098	1,45472	1,455697029	0,244257
4	-0,68	0,00511	1,4632533	1,468359227	1,276474
5	-0,64	0,00479	1,4743467	1,479140068	1,19835
6	-0,6	0,00693	1,488	1,494932532	1,733133
7	-0,56	-0,00873	1,5042133	1,495478983	-2,18359
8	-0,52	-0,00094	1,5229867	1,522049942	-0,23418
9	-0,48	0,00438	1,54432	1,54870009	1,095023
10	-0,44	-0,00435	1,5682133	1,563866531	-1,0867
11	-0,4	-0,00276	1,5946667	1,59190585	-0,6902
12	-0,36	-0,00676	1,62368	1,616918271	-1,69043
13	-0,32	-0,00739	1,6552533	1,64786569	-1,84691
14	-0,28	-0,00391	1,6893867	1,685476149	-0,97763
15	-0,24	-0,00309	1,72608	1,722985972	-0,77351
16	-0,2	-0,00847	1,7653333	1,756861608	-2,11793
17	-0,16	-0,00227	1,8071467	1,804874967	-0,56792
18	-0,12	-0,00162	1,85152	1,84990381	-0,40405
19	-0,08	0,00054	1,8984533	1,898992746	0,134853
20	-0,04	-0,00146	1,9479467	1,946484695	-0,36549
21	0	-0,00131	2	1,998692037	-0,32699
22	0,04	-0,00148	2,0546133	2,053132371	-0,37024
23	0,08	0,00537	2,1117867	2,117157233	1,342642
24	0,12	-0,00034	2,17152	2,171178862	-0,08528
25	0,16	-0,00074	2,2338133	2,233068703	-0,18616
26	0,2	-0,00205	2,2986667	2,296613837	-0,51321
27	0,24	0,00789	2,36608	2,373968848	1,972212
28	0,28	0,00346	2,4360533	2,439516025	0,865673
29	0,32	0,0095	2,5085867	2,518089286	2,375655
30	0,36	-0,00262	2,58368	2,581060373	-0,65491
31	0,4	0,00665	2,6613333	2,667979157	1,661456
32	0,44	-0,00645	2,7415467	2,735097076	-1,6124
33	0,48	0,00216	2,82432	2,826475793	0,538948
34	0,52	0,00361	2,9096533	2,913262099	0,902191
35	0,56	0,00768	2,9975467	3,005222329	1,918916
36	0,6	-0,00034	3,088	3,087661932	-0,08452
37	0,64	-0,0021	3,1810133	3,178918153	-0,5238
38	0,68	0,0027	3,2765867	3,27928722	0,675138
39	0,72	-0,00153	3,37472	3,373194705	-0,38132
40	0,76	0,00303	3,4754133	3,478443779	0,757611
41	0,8	-0,00578	3,5786667	3,57288992	-1,44419

Этап 2. Оценка МНК

Статистика - любая функция результатов опытов, которая не зависит от неизвестных статистических характеристик. Оценкой статистической характеристики называется статистика, реализация которой, полученная в результате опытов, принимается за неизвестное истинное значение параметра. Если математическое ожидание оценки равно характеристике , то оценка несмещенная. Разность М^- - смещение оценки. Оценка статистической характеристики называется состоятельной, если она сходится по вероятности к при неограниченном увеличении опытов.

Проведя серию «опытов» в первом этапе, получили значение выходного параметра у и значение входных параметров 1, t, t². Теперь по этим данным найдем значения коэффициентов b^ в аппроксимации y_i = b₁^+b₂^*t_i + b₃^ *t²_i, то есть оценку b^ вектора параметров b.

По известным входным и выходным параметрам найдем коэффициенты b₁^, b₂^, b₃^ при входных параметрах 1, t, t² соответственно, которые будут составлять вектор оценок вектора b и будут точечными оценками параметров b₁, b₂, b₃ в исходном уравнении y_i = b₁+b₂*t_i+b₃*t_i²+_i.. Так как значение y_i, t_i, t_i², _i переменных в I-ом по счету опытах, то уравнение y_i = b₁+b₂*t_i+b₃*t_i²+_i можно записать в виде системы уравнений:

y₁ = b₁+b₂*t₁+b₃*t₁²+₁

y₂ = b₁+b₂*t₂+b₃*t₂²+₂

…

y_n = b₁+b₂*t_n+b₃*t_n²+_n

Или их можно записать в векторной записи: y_i = X^Tb + , где ~N(0,²I_N) - нормально распределенный случайный вектор ошибок измерений такой, что E = 0, K = E(^T) = ²I_N, I_NR^N*N - единичная матрица, параметр b=(b₁, b₂, b₃)^T R₃, X^T = X_ij R_3XN.

Коэффициенты b₁^, b₂^, b₃^ должны быть такими, чтобы вектор y - y^ = y-b₁^+b₂^*t+b₃^*t² , рассматриваемый как элемент евклидового пространства. Решение этой задачи известно по курсу линейной алгебры. Искомый вектор y^= b₁^+b₂^*t+b₃^*t² , является ортогональной проекцией вектора у на подпространство, порожденное векторами 1, t, t² (x₁, x₂, …,x_n), так что скалярное произведение разности у - y^ с каждым из векторов 1, t, t² (x₁, x₂, …,x_n) равно нулю: (x_i,y - y^)=0 (i=1,2,…,m).

Данное соотношение выражает необходимое и достаточное условие минимума длины вектора у - y^. Перепишем их в виде системы относительно искомых уравнений:

(x₁,x₁)b^₁+(x₁,x₂)b^₂+…+(x₁,x_m)b^_m=(x₁,y)

(x₂,x₁)b^₁+(x₂,x₂)b^₂+…+(x₂,x_m)b^_m=(x₂,y)

…

(x_m,x₁)b^₁+(x_m,x₂)b^₁+…(x_m,x_m)b^_m=(x_m,y)

Данная система всегда совместна и однозначно определяет вектор y^= b₁^х₁+b₂^*х₂+…+b_m^*x_m, хотя сами коэффициенты b₁^, b₂^, b₃^… b_m^ могут находиться неоднозначно - если векторы x₁, x₂, …,x_m линейно не зависимы.

Матрица

(Х₁,Х₁) (Х₁,Х₂) … (Х₁,Х_m)

А = (Х₂,Х₁) (Х₂,Х₂) … (Х₂,Х_m)

…

(Х_m,Х₁) (Х_m,Х₂) …(Х_m,Х_m)

называется матрицей Грамма для системы векторов x₁, x₂, …,x_m. Эта матрица является невырожденной в том и только том случае, если x₁, x₂, …,x_m - линейно независимая система. Для удобства вводим матрицу:

Х₁₁ Х₁₂ … Х_1n

Х = Х₂₁ Х₂₂ … Х_2n

…

Х_m1 Х_m2 … Х_mn

Тогда систему уравнений можно переписать в виде так называемой нормальной системы уравнений:

ХХ^Tb^=Xy или

Ab^=Xy,

где

A=XX^T, так как предположили, что вектора линейно не зависимы, а значит det A 0, получаем:

b^=A^-1Xy.

Расчет оценки b^ вектора параметров b.

Уравнение y_i = b₁^+b₂^*t_i+b₃^*t²_i можно записать в матричном виде у= Х^Tb^, где параметр b^=( b₁^, b₂^, b₃^)^Т R_N и X^T = X_ij R_3XN - получены в этапе 1 и представлены в таблице 2.

Решением данного уравнения будет уравнение в матричной форме b^=A^-1Xy, где A=XX^T. Матрица = (1/detA)A^*, где detA - определитель матрица, А^* - присоединенная матрица.

Для расчета матрицы A^-1 воспользуемся таблицей 3. Перемножим сначала матрицы Х и У, а затем полученную матрицу умножим на матрицу A^-1.

Матрица А:

	41	0	9,184
A=	0	9,184	0
	9,184	0	3,70004992

Матрица A^-1:

	0,05493	-0	-0,13634975
A-1=	0	0,108885	0
	-0,1363	0	0,608704257

Определитель матрицы А:

det A=
618,5992556

312,843646

Решение системы: Ху = 124,346784

266,283407

Расчет оценки вектора параметра b:

b^ =	0,05493	-0	-0,13634975	*	312,843646	=	1,636941
	0	0,108885	0		124,346784		1,508154
	-0,1363	0	0,608704257		266,283407		1,675816

Так как у= ХТ + и b^ = A^-1Xy, то b^ = A^-1X (ХТ + ) = b+ A^-1X, b^ - является смещенной оценкой вектора.

Оценка b^ вектора параметров b: b₁^ =1,636941; b₂^ =1,508154; b₃^ =1,675816.

Или в матричной форме: b^ = (1,636941;1,508154;1,675816)^Т.

Определение оценки ²^ дисперсии ошибок ²:

Оценка дисперсии ² рассчитывается по формуле ²^=S²/(n-m). Статистика S²/(n-m) говорит о том, что оценка ²^ несмещенная оценка дисперсии ². y^ рассчитывается так: y^ =X^T*b = X^TA^-1xy.

Расчет оценки ²^ дисперсии ²:

S² = (y-y^)^T(y-y^)=(y-y^)².

Рассчитываем при помощи таблицы 3.

n - число опытов;

m - число параметров;

²^=0.01514/(41-3) = 0.000398;

Таблица 2. Матрицы y, X^T, b,

	1,44413		1	-0,8	0,64	Оценка b^	1,636941	о =	-0,0012
	1,44364		1	-0,76	0,5776		1,508154		-0,00511
	1,4557		1	-0,72	0,5184		1,675816		0,000977
	1,46836		1	-0,68	0,4624			0,005106
	1,47914		1	-0,64	0,4096			0,004793
	1,49493		1	-0,6	0,36			0,006933
	1,49548		1	-0,56	0,3136			-0,00873
	1,52205		1	-0,52	0,2704			-0,00094
	1,5487		1	-0,48	0,2304			0,00438
	1,56387		1	-0,44	0,1936			-0,00435
	1,59191		1	-0,4	0,16			-0,00276
	1,61692		1	-0,36	0,1296			-0,00676
	1,64787		1	-0,32	0,1024			-0,00739
	1,68548		1	-0,28	0,0784			-0,00391
	1,72299		1	-0,24	0,0576			-0,00309
	1,75686		1	-0,2	0,04			-0,00847
	1,80487		1	-0,16	0,0256			-0,00227
	1,8499		1	-0,12	0,0144			-0,00162
	1,89899		1	-0,08	0,0064			0,000539
	1,94648		1	-0,04	0,0016			-0,00146
y =	1,99869	XT =	1	0	0			-0,00131
	2,05313		1	0,04	0,0016			-0,00148
	2,11716		1	0,08	0,0064			0,005371
	2,17118		1	0,12	0,0144			-0,00034
	2,23307		1	0,16	0,0256			-0,00074
	2,29661		1	0,2	0,04			-0,00205
	2,37397		1	0,24	0,0576			0,007889
	2,43952		1	0,28	0,0784			0,003463
	2,51809		1	0,32	0,1024			0,009503
	2,58106		1	0,36	0,1296			-0,00262
	2,66798		1	0,4	0,16			0,006646
	13,14		1	0,44	0,1936			-0,00645
	2,82648		1	0,48	0,2304			0,002156
	2,91326		1	0,52	0,2704			0,003609
	3,00522		1	0,56	0,3136			0,007676
	3,08766		1	0,6	0,36			-0,00034
	3,17892		1	0,64	0,4096			-0,0021
	3,27929		1	0,68	0,4624			0,002701
	3,37319		1	0,72	0,5184			-0,00153
	3,47844		1	0,76	0,5776			0,00303
	3,57289		1	0,8	0,64			-0,00578

Таблица 3. Вспомогательные расчеты для матрицы А^-1, произведения матриц Х и У, оценки ²^.

y	ti2	ti3	ti4	y*t	y*t2	y^	y-y^	(y-y^)2
1,44413	0,64	-0,512	0,4096	-1,15531	0,924245	1,34244	0,101692	0,010341
1,44364	0,5776	-0,43898	0,333622	-1,09716	0,833844	1,23258	0,211056	0,044545
1,4557	0,5184	-0,37325	0,268739	-1,0481	0,754633	1,218767	0,23693	0,056136
1,46836	0,4624	-0,31443	0,213814	-0,99848	0,678969	1,39876	0,069599	0,004844
1,47914	0,4096	-0,26214	0,167772	-0,94665	0,605856	1,43998	0,03916	0,001534
1,49493	0,36	-0,216	0,1296	-0,89696	0,538176	1,50763	-0,0127	0,000161
1,49548	0,3136	-0,17562	0,098345	-0,83747	0,468982	1,519507	-0,02403	0,000577
1,52205	0,2704	-0,14061	0,073116	-0,79147	0,411562	1,566026	-0,04398	0,001934
1,5487	0,2304	-0,11059	0,053084	-0,74338	0,356821	1,612544	-0,06384	0,004076
1,56387	0,1936	-0,08518	0,037481	-0,6881	0,302765	1,659062	-0,0952	0,009062
1,59191	0,16	-0,064	0,0256	-0,63676	0,254705	1,705581	-0,11367	0,012922
1,61692	0,1296	-0,04666	0,016796	-0,58209	0,209553	1,752099	-0,13518	0,018274
1,64787	0,1024	-0,03277	0,010486	-0,52732	0,168741	1,798617	-0,15075	0,022726
1,68548	0,0784	-0,02195	0,006147	-0,47193	0,132141	1,845136	-0,15966	0,025491
1,72299	0,0576	-0,01382	0,003318	-0,41352	0,099244	1,891654	-0,16867	0,028449
1,75686	0,04	-0,008	0,0016	-0,35137	0,070274	1,938172	-0,18131	0,032874
1,80487	0,0256	-0,0041	0,000655	-0,28878	0,046205	1,984691	-0,17982	0,032334
1,8499	0,0144	-0,00173	0,000207	-0,22199	0,026639	2,031209	-0,18131	0,032872
1,89899	0,0064	-0,00051	4,1E-05	-0,15192	0,012154	2,077728	-0,17873	0,031946
1,94648	0,0016	-6,4E-05	2,56E-06	-0,07786	0,003114	2,124246	-0,17776	0,031599
1,99869	0	0	0	0	0	2,003643	-0,00495	2,45E-05
2,05313	0,0016	0,000064	2,56E-06	0,082125	0,003285	2,217283	-0,16415	0,026945
2,11716	0,0064	0,000512	4,1E-05	0,169373	0,01355	2,263801	-0,14664	0,021504
2,17118	0,0144	0,001728	0,000207	0,260541	0,031265	2,310319	-0,13914	0,01936
2,23307	0,0256	0,004096	0,000655	0,357291	0,057167	2,356838	-0,12377	0,015319
2,29661	0,04	0,008	0,0016	0,459323	0,091865	2,403356	-0,10674	0,011394
2,37397	0,0576	0,013824	0,003318	0,569753	0,136741	2,449875	-0,07591	0,005762
2,43952	0,0784	0,021952	0,006147	0,683064	0,191258	2,496393	-0,05688	0,003235
2,51809	0,1024	0,032768	0,010486	0,805789	0,257852	2,542911	-0,02482	0,000616
2,58106	0,1296	0,046656	0,016796	0,929182	0,334505	2,58943	-0,00837	7E-05
2,66798	0,16	0,064	0,0256	1,067192	0,426877	2,635948	0,032031	0,001026
13,14	0,1936	0,085184	0,037481	5,7816	2,543904	2,682466	10,45753	109,36
2,82648	0,2304	0,110592	0,053084	1,356708	0,65122	2,728985	0,097491	0,009504
2,91326	0,2704	0,140608	0,073116	1,514896	0,787746	2,775503	0,137759	0,018978
3,00522	0,3136	0,175616	0,098345	1,682925	0,942438	2,822022	0,183201	0,033563
3,08766	0,36	0,216	0,1296	1,852597	1,111558	2,86854	0,219122	0,048014
3,17892	0,4096	0,262144	0,167772	2,034508	1,302085	2,915058	0,26386	0,069622
3,27929	0,4624	0,314432	0,213814	2,229915	1,516342	2,961577	0,317711	0,10094
3,37319	0,5184	0,373248	0,268739	2,4287	1,748664	3,008095	0,3651	0,133298
3,47844	0,5776	0,438976	0,333622	2,643617	2,009149	3,054613	0,42383	0,179632
3,57289	0,64	0,512	0,4096	2,858312	2,28665	3,101132	0,471758	0,222556

Этап 3. Доверительные интервалы

Случайный интервал, полностью определяемый результатами опытов и не зависящий от неизвестных характеристик, который с заданной вероятностью накрывает неизвестную скалярную статистическую характеристику ,называется доверительным интервалом для этой характеристики, соответствующим интервалом доверия .

Так как ~N(0,² I_N), то можно найти не только точечные оценки, но и интервальные оценки неизвестных параметров. Так как выполняются условия:

у= Х^Т* b+

Е=0; К=Е() =²I_N

Аb^=Ху

det А0, то S²(b)= (y- Х^Т b)^T(y- Х^Т b);

S²(b^)= (y- Х^Т b^)^T(y-Х^Т b^);

Т²= S²(b)- S²(b^), тогда

b^~ N(b ;² А^-1);

S²(b^)/ ²~(n-m);

Т²/ ²~(m).

При этом b^ и S²(b^), а также S²(b) и Т² независимы.

(r) - хи-квадрат распределение с r степенями свободы. То есть распределение случайной величины (- независимые и одинаково распределенные по закону N(0,1) случайные величины). Отсюда вытекает, что:

Р(S²(b^)/_1-/2,n-m<2<S2(b^)/_/2,n-m)=1-,

где 2n-квантиль уровня для 2- распределения с n- степенями свободы.

Далее,

(b_i^-b_i/c_ii)/( S²(b^)/)=*( b_i^- b_i)/c_ii S(b^)~t(n-m),

где с_ii- обозначает (i,i)- й элемент матрицы А^-1, а символ t(n-m)- распределение Стьюдента с n-m степенями свободы, распределение случайной величины. Отсюда:

Р(b_i^-t_1-/2,n-m*c_iiS(b^)/<b_i<b_i^+t_1-/2,n-m*c_iiS(b^)/)=1-,

где t- квантиль уровня для распределения Стьюдента с n степенями свободы.

Расчет доверительного интервала ² .

По формуле:

Р(S²(b^)/_1-/2,n-m<²<S2(b^)/_/2,n-m)=1-, где

S²(b^)- значение, которое мы берем из этапа 2(Таблица № 3 ( y- y ^)²) _1-/2,n-m; _/2,n-m- квантиль уровня 0,95 (0,05) для распределения с 38 степенями свободы и равно 53,38351 (24,88389), соответственно.

Р(S²(b^)/_0,95,38<²<S2(b^)/_0,05,38)=0,9,

Р( 123Е-7 <²< 170Е-6 )=0,9;

Доверительный интервал для ² - 81Е-6 <²< 111Е-6

Определяем доверительный интервал для b₁^, b₂^, b₃^.

По формуле:

Р(b_i^-t_1-/2,n-m*c_iiS(b^)/<b_i<b_i^+t_1-/2,n-m*c_iiS(b^)/)=1-;

где: t_1-- распределение Стьюдента с n-m степенями свободы.

t_0,95,38 = 2,04227

S(b^)- значение из этапа 2, S(b^)== 0,123044707

С_ii- соответствующий элемент матрицы А^-1, значение берем из этапа 2.

Р(b_i^- 2,04227 *c_ii 0,123044707 / <b_i<b_i^+2,04227 *c_ii *0,123044707 /)=0,9

Доверительный интервал для b₁^:

3,996503273 < b₁^< 4,00619438

Доверительный интервал для b₂^:

2,999021874 < b₂^< 3,015350856

Доверительный интервал для b₃^:

0,988198305 < b₃^< 1,00450879

Список литературы

Кочетков Е. С. Метод наименьших квадратов. Учебное пособие. М.: Издательство МАИ, 1993

Чистяков В. П. Курс теории вероятностей. М.: Наука, 1982.

Ивченко Г. И., Медведев Ю. Н. Математическая статистика. М.: Высшая школа, 1984.

Размещено на Allbest.ru