Modely

smysl statistických modelů je

redukovat (komprimovat) data - popsat výsledek měření několika parametry
činit předpovědi - pokud závěry učiněné na daném vzorku extrapolujeme na zbytek populace
lze parametrizovat i simulovaná data (ušetření výpočetního času)

paradigma: získání spolehlivých údajů z nedokonalých dat

mohou být

parametrické (tušíme funkční závislost)
neparametrické (vyhlazování / kubický spline)

$Y=r(X) + e$

platí $E(e)=0$, i když obecně rozdělení $e$ závisí na $X$

z daného vzorku máme jen odhad $\hat{r}(X)$, který se liší od skutečné $r(X)$ vlivem nejistot parametrů

$$V(Y-\hat{r}) = \sigma^2 + E((\hat{r}-r)^2)= \sigma^2 + (E(\hat{r})-r)^2 + V(\hat{r})$$

bias - variance decomposition - větší počet parametrů snižuje rezidua, ale zvyšuje nejistoty modelu (v důsledku korelací parametrů)

mnohorozměrná regrese je analogií vzorce pro proložení přímkou ($\sigma_{xy}/\sigma_{xx}$)

$$\beta= v^{-1} D(\vec{X},Y)$$

transformace pomocí $v^{-1}$ odstraňuje korelace mezi komponentami modelu

vektor reziduí $Y-\beta X$ musí být dekorelován (ortogonální) s vektorem $\vec{X}$

konstrukce modelu

máme-li diskrétní hodnoty (kategorie) - může stačit pruměrování (očekávaná hodnota) a rozptyl v rámci kategorie

v praxi spíš hledáme pro model spojitou funkci (interpolace mezi měřenými hodnotami)

hledáme funkce z nějaké množiny "rozumných"

otázka "jak přesně sledujeme data" vede k "bias-variance" trade-off...

lineární model $r = A H^{-1} A^T Y$ patří do skupiny lineárního vyhlazování (linear smoother) daného obecnější formulí $$ \hat{r}(x) = \sum_i{y_i w(x_i,x)} $$

pro lin. model $w(x_i,x)=(x_i/n s^2_x) x$

k-nejbližších sousedů

$w(x_i,x)=1/k$ pokud jde o jednoho z k sousedů bodu $x$, jinak 0

vyhlazování s kernelem

$w(x_i,x)=K(x-x_i)/\sum_j K(x-x_j)$, kde nezáporný kernel splňuje podmínky $\int x K(x) dx=0$, $\int x^2 K(x) dx< \infty$

pro rovnoměrně rozdělená data vypadá jako konvoluce

parametrizace

náš odhad parametrů $\theta$ minimalizuje "loss function" $L(\bf{z}_n;\theta)$ (může to být např. záporný logaritmus věrohodnosti) v prostoru parametrů; přitom tato funkce se liší od její "střední hodnoty" určované nad celou populací dat (a která by minimalizací dala skutečné hodnoty parametrů)

$$L(\bf{z}_n;\theta)=E(L(\bf{Z};\theta)) + \eta_n(\theta)$$

skryté a přebytečné proměnné

rozdíl mezi věrohodností třídy $\pi$ modelů obecnějších a $\rho$ modelů s omezenými (fixovanými) parametry je úměrný $\chi^2_{p-q}$ (p-q je počet fixovaných parametrů)

rezidua jsou z definice nekorelována s modelem (nezávislými parametry), nicméně by měly splňovat také podmínky pro bílý šum

normální rozdělení
stejné rozdělení (zejména rozptyl) pro různá x
navzájem nekorelované

testování dalších parametrů

t-test parametrů s nejistotami
F-test tříd modelů

kovarianční matice klesá s 1/N (objemem měřených dat), stále více parametrů může být statisticky významných, ale skutečné fyzikální mechanismy na velikosti vzorku nezávisí

přidání dalšího parametru do modelu musí mít fyzikální opodstatnění
pouhá korelace neznamená příčinnou souvislost
pokud není splněna podmínka normálnosti rozdělení parametrů, standardní testy nepomůžou (lze aplikovat bootstraping)



In [28]:

    
%pylab inline

x=r_[-3:3:20j]
plot(x,7*x**2-0.5*x,'k')
x=uniform(-3,3,20)
tres=[7,-.5,0]
ytrue=polyval(tres,x)
y=ytrue+normal(size=x.shape)
plot(x,y,'*')
ords=arange(1,10)
res=[polyfit(x,y,i,cov=True) for i in ords]
[[round(p,3) for p in r[0][::-1]] for r in res]









    



Populating the interactive namespace from numpy and matplotlib






    Out[28]:





[[19.445, -0.848],
 [0.022, -0.381, 7.031],
 [-0.013, -0.663, 7.036, 0.051],
 [-0.091, -0.563, 7.167, 0.034, -0.018],
 [-0.044, -0.196, 6.99, -0.208, 0.011, 0.027],
 [-0.133, 0.009, 7.279, -0.344, -0.109, 0.044, 0.011],
 [-0.138, -0.413, 7.549, 0.199, -0.248, -0.118, 0.027, 0.013],
 [-0.307, -0.118, 8.46, -0.133, -0.902, -0.02, 0.171, 0.005, -0.01],
 [-0.312, -0.361, 8.702, 0.444, -1.099, -0.337, 0.217, 0.065, -0.013, -0.004]]

validace spočívá v aplikaci modelů na nový vzorek dat

pokud nemáme nový vzorek, můžeme dělat podvýběry ze stávajícího

rozdělení na poloviny, obecněji k-dílů, jeden díl se použije na testování
jackknife (viz redukce biasu)
bootstrap



In [38]:

    
chi2=r_[[((y-polyval(res[i-1][0],x))**2).sum()/(len(x)) for i in ords]]
semilogy(ords,chi2,'s')
grid()
semilogy(ords,chi2/(len(x)-ords-1)*len(x),'d')
ynew=ytrue+normal(size=x.shape)
gme=r_[[((ytrue-polyval(res[i-1][0],x))**2).sum()/(len(x)) for i in ords]]
semilogy(ords,gme,'o')#,fillcolor=None)
valme=r_[[((ynew-polyval(res[i-1][0],x))**2).sum()/(len(x)-i-1) for i in ords]]
semilogy(ords,valme,'v')
legend(['mse','chi2','gme','valid.'])
ylim(0.02,30)









    Out[38]:





(0.02, 30)

oblíbené omyly v souvislosti s lineární regresí

proměnná s (významně) nenulovým koeficientem ovlivňuje výsledek
proměnná s nulovým koeficientem nemůže ovlivňovat výsledek
změna proměnných ovlivní výsledek podle regresního vzorce

korelace 2 proměnných je často důsledkem závislosti na třetím faktoru problematika významnosti korelace viz zde