ČESKÉ VYSOKÉ UČENÍ TECHNICKÉ V PRAZE FAKULTA STAVEBNÍ, OBOR GEODÉZIE A KARTOGRAFIE KATEDRA VYŠŠÍ GEODÉZIE název předmětu
TEORIE CHYB A VYROVNÁVACÍ POČET 2 č. úlohy název úlohy
6 zadání
68
TESTOVÁNÍ STATISTICKÝCH HYPOTÉZ studijní skupina
G-263
zpracoval
SETNIČKA M ARTIN
e-mail
[email protected]
datum
31.5. 2007
ZADÁNÍ
TECHNICKÁ ZPRÁVA Numerické řešení bylo provedeno pomocí skriptu programu Matlab. Byl použit statistický toolbox. Následuje teoretický postup. Číselné hodnoty jsou sestaveny přehledně v závěru. Ve všech případech byly prováděny oboustranné testy. Příklad 1 Část a. - Pro výpočet zde i v části d. byly použity standardní funkce 'mean' a 'std'. lA = mean( vektor_mereni ) % vyrovnaná délka sgmA = std( vektor_mereni ) % aposteriorni smerodatna odchylka jednotkova Část b. lA: Výpočet hranic intervalu pomocí inverzní distribuční funkce studentova rozdělení 'tinv': lA +/- tinv( 1-alfa/2 , m-1 ) *sgmA /sqrt( m ) kde 'alfa' je hladina významnosti, 'm' počet pozorování a 'm-1' počet nadbytečných měření σA: Výpočet hranic intervalu pomocí inverzní distribuční funkce chí-kvadrát rozdělení 'chi2inv': sqrt( m *sgmA^2 /chi2inv( alfa/2 , m-1 ) ) % dolní hranice sqrt( m *sgmA^2 /chi2inv( 1-alfa/2 , m-1 ) ) % horní hranice Část c. t_test = abs( lA-l0 ) *sqrt( m ) /sgmA % hodnota testovaciho kriteria tinv( 1-alfa/2 , m-1 ) % kriticka hodnota ( 1 -tcdf( t_test , m-1 ) )*2 % p-hodnota pomoci distribucni fce Část d. - Hodnota testovacího kritéria bylo volena < 1 sgmA^2 /sgmB^2; % hodnota testovaciho kriteria .. a do inverzní distribuční funkce Fisherova rozdělení byla použita pravděpodobnost odpovídající dolní mezi. finv( alfa/2, m-1, n-1 ) % kriticka hodnota fcdf( f_test, m-1, n-1 )*2 % p-hodnota pomoci distribucni fce Příklad 2 Část a. 1) Rovnice oprav ( Jde o opravy provozní doby 'T', malé 't' je teplota, 'b' jsou odhadované koeficienty ) v = b0 + b1 .*t + b2 .*t.^2 - T 2) Matice plánu: A = [ ones(t), t, t.^2 ] 3) Úloha je lineární, proto volíme redukovaná měření jsou 'T' a výsledkem jsou přímo koeficienty 'b' b = inv( A'*A ) * A'*T Část b. 1) v = A*b - T % opravy 2) m02 = v'*v /(n-3) % aposteriorni smerodatna odchylka jednotkova ve druhe mocnine 3) Qb = m02 *inv( A'*A ) % kovarianční matice vyrovnaných koeficientů 4) p = 2:3; % pořadí testovaných koeficientů ve vektru 'b' a matici 'Qb' q = length(p); % počet testovaných koeficientů Hodnota testovacího kritéria bylo > 1 test = 1/q/m02 *b(p)'*Qb(p,p)^-1*b(p) % hodnota testovaciho kriteria .. a do inverzní distribuční funkce Fisherova rozdělení byla použita pravděpodobnost odpovídající horní mezi. kh = finv(1-alfa/2,q,n-3) % kriticka hodnota kde 'n' je počet dvojic měřen a 'n-3' nadbytečný počet měření Část c. ( viz příklad 1 c. ) test = abs( b(3)-0 ) / Qb(3,3)^0.5 * sqrt(n) % hodnota testovaciho kriteria kh = tinv( 1-alfa/2, n-3 ) % kriticka hodnota Část d. 1) Výpočet jednotlivých směrodatných odchylek 'T' : sgmT = diag( A*Qb*A' ).^0.5 1) Pro výpošet mezí viz příklad 1 b. T_aprox +/- tinv( 1-alfa/2 , n-3 ) .* sgmA kde 'T_aprox' jsou provozní doby spočtené pomocí vyrovnaných koeficientů
Příklad 3 Část a. - Pro přimky vyrovnané minimalizací oprav 'x' a 'y' byla použita funkce 'polyfit' ( 3. argument je stupen regresního polynomu) Označení: x ~ výška, y ~ váha, a, b ~ vyrovnané koeficienty regresních přímek 1) Minimalizace oprav 'y'
ay = ( AT A) −1 . A. y b y
; kde A je matice plánu Ay = [ones(n,1),x]
2) Minimalizace oprav 'x'
ax = ( AT A) −1 . A. x bx
; kde A je matice plánu Ay = [ones(n,1),y]
3) Minimalizace oprav 'x' a 'y' ( minimalizace vzálenosti od regresní přímky ) xt = mean(x) % x-souradnice teziste yt = mean(y) % y-souradnice teziste vx = xt-x % redukce k tezisti vy = yt-y % redukce k tezisti alfa = atan2(2*sum(vx.*vy),((sum(vx.^2)-sum(vy.^2))) )/2; axy(1) = yt - tan(alfa)*xt bxy(2) = tan(alfa) Část b. % q .. sumy oprav k tezisti qxy = vx' * vy qxx = vx' * vx qyy = vy' * vy r1 = qxy / sqrt( qxx*qyy ) % vypocet primo (z oprav) r2 = sqrt( ay/ax ) % druhy vypocet ze smernic Část c. test = r1 *sqrt(n-2) /sqrt(1-r1^2) % hodnota testovaciho kriteria kh = tinv( 1-alfa/2,n-2 ) % kriticka hodnota
Z Á V Ě R: V Ý S L E D K Y: alfa pro všechny výpočty: α = 0.05 Příklad 1 Část a. lA = 100.0016 m σA = 0.0028 m Část b. P { 99.9999 m < lA < 100.0033 m } = 95% P { 0.0021 m< σA < 0.0048 m } = 95% Část c. H0: l0 = lA testovací kritérium ... 2.0474 ktitická hodnota ... 2.1788 p-hodnota ... 0.0632 > alfa Část d. H0: σA = σB testovací kritérium ... 0.3301 ktitická hodnota ... 0.3540 p-hodnota ... 0.0368 < alfa Příklad 2 Část a.
=> Nezamítáme nulovou hypotézu.
=> Zamítáme nulovou hypotézu. Přijímáme hypotézu H1: σA ≠ σB
Část b.
β0 = 12.4874 σβ0 = 0.3226 β1 = 0.1855 σβ1 = 0.0318 β2 = -0.0027 σβ2 = 0.00074 H0: β2 = 0 ∧ β1 = 0 testovací kritérium ... 123.329 ktitická hodnota ... 5.71471 p-hodnota ... 2.88e-007 < alfa => Zamítáme nulovou hypotézu. Přijímáme hypotézu H1: β2 ≠ 0 ∨ β1 ≠ 0
Část c.
H0: β2 = 0
testovací kritérium ... 3.7065 ktitická hodnota ... 2.2622 p-hodnota ... 0.0049 < alfa
=> Zamítáme nulovou hypotézu. Přijímáme hypotézu H1: β2 ≠ 0
Část d.
Příklad 3 Část a.
Část b. Část c.
y=ax+b 1) Minimalizace oprav 'y': a = 1.1022 b = -109.6608 2) Minimalizace oprav 'x': a = 2.6554 b = -388.0151 3) Minimalizace oprav 'x' a 'y': a = 2.2022 b = -306.7968 Výsledky obou postupů se dle předpokladu neliší r = 0.6443 H0: r = 0 ktitická hodnota ... P (-1,98 < r < 1,98) = 95% p-hodnota ... 0.2604 > alfa => Nezamítáme nulovou hypotézu.