17.1 OLS a lineáris algebra és statisztika találkozásánál

17.1.1 Cél és konkrét helyzet

A Northstar rendelési adatból azt szeretnénk tudni, hogyan függ össze a kosár értéke néhány, már a rendelés pillanatában ismert jellemzővel. A megfigyelési egység most egy rendelés: minden sor egy order_id-hez tartozó kosár, a célváltozó pedig a basket_minor, azaz a kosár értéke egész eurócentben. Bemeneti jellemzőnek a tételdarabszámot (quantity_sum) és a prioritásjelzőt (priority_flag, értéke 0 vagy 1) választjuk; a távolság distance_km a feladatoknál csatlakozik. Célunk kettős: ugyanazt a lineáris illesztést előbb a NumPy lineáris algebrai megoldásával, majd a statsmodels formula API-jával elvégezni, és a két útból származó együtthatót, reziduumot és tanítóhibát egymással ellenőrizni.

Induló becslésként fogalmazzuk meg, mit várunk: a quantity_sum együtthatója várhatóan pozitív, mert több tétel általában nagyobb kosárértékkel jár; a priority_flag együtthatója a prioritásos és nem prioritásos rendelések feltételes különbségét méri azonos darabszám mellett. Azt nem tudjuk előre, mekkora ez a különbség, ezért a modell feladata éppen ennek becslése. Az alábbi mini-tábla tisztán oktatási adat, amelyben a 0 darabos sor az intercept értelmezését teszi kézzel láthatóvá; a valódi Northstar orders.csv csak pozitív darabszámokat tartalmaz.

17.1.2 Fogalmi magyarázat: a modell és a mátrixok

A lineáris modell minden $i=1,\dots,n$ rendelésre egy egyenletet ír fel:

$y_i=\beta_0+\beta_1 x_{1i}+\beta_2 x_{2i}+e_i,$

ahol $y_i$ a basket_minor eurócentben, $x_{1i}$ a quantity_sum darabszámban, $x_{2i}$ a priority_flag (0 vagy 1), $\beta_0,\beta_1,\beta_2$ ismeretlen együtthatók, $e_i$ pedig a reziduum. A reziduum ugyanabban a mértékegységben van, mint a cél: eurócentben. Az egyenlet azt mondja, hogy a megfigyelt kosárérték egy lineáris, szisztematikus rész és egy zaj összege.

Az egyenleteket mátrix alakba rendezve kapjuk a dizájnmátrixot (design matrix):

$\mathbf{y}=X\boldsymbol\beta+\mathbf{e},$

ahol $\mathbf{y}$ alakja $(n,)$ , $X$ alakja $(n,p)$ , $\boldsymbol\beta$ alakja $(p,)$ , és $\mathbf{e}$ alakja $(n,)$ . A $p$ a paraméterek száma; három paraméterünk van, tehát $p=3$ . A dizájnmátrix első oszlopa csupa egyes:

$X=\begin{pmatrix}1&x_{11}&x_{21}\\1&x_{12}&x_{22}\\\vdots&\vdots&\vdots\\1&x_{1n}&x_{2n}\end{pmatrix}.$

Ez az egyesekből álló oszlop adja az interceptet (tengelymetszet). A modellben $\beta_0$ tehát nem valami különleges mennyiség, hanem egy további együttható, amelyhez a konstans magyarázó „változó” tartozik. Értelmezése: a basket_minor várható értéke, ha quantity_sum = 0 és priority_flag = 0. Ha a 0 darabos cella kívül esik a megfigyelt tartományon, az intercept extrapoláció, ezért óvatosan kell értelmezni.

A coefficient (együttható) jelentése parciális meredekség: $\beta_1$ azt adja meg, mennyivel változik a modell szerinti várható kosárérték, ha a darabszám egy egységgel nő, miközben a másik prediktor rögzített. Mértékegysége eurócent/darab. $\beta_2$ az a feltételes különbség eurócentben, amely a prioritásos rendeléshez tartozik azonos darabszám mellett. Ez asszociációs értelmezés, nem oksági bizonyíték.

A reziduum (residual) a megfigyelt és a modell által becsült érték különbsége:

$\hat y_i=x_i\hat\beta,\qquad e_i=y_i-\hat y_i.$

Az $x_i$ itt az $X$ mátrix $i$ -edik sora, alakja $(p,)$ . A reziduumvektor hossza tehát $n$ , nem $p$ . A legkisebb négyzetek (least squares) elve az az együtthatóvektort választja, amely a reziduális négyzetösszeget, az $RSS=\sum_i e_i^2$ mennyiséget minimalizálja. Az RSS mértékegysége eurócent a négyzeten.

A gradiensmódszer logikájához visszanyúlva Optimalizáló választás és konvergenciaellenőrzés az RSS gradiense akkor nulla, ha teljesül a normálegyenlet:

$X'X\hat\beta=X'\mathbf{y},$

amelynek megoldása teljes oszloprangú $X$ esetén

$\hat\beta=(X'X)^{-1}X'\mathbf{y}.$

Itt $X'X$ alakja $(p,p)$ , $X'\mathbf{y}$ alakja $(p,)$ , tehát egy kis lineáris egyenletrendszert oldunk meg. A megoldás geometriailag vetítés Sajátérték, SVD és dimenziócsökkentés előképe: $H=X(X'X)^{-1}X'$ a kalapmátrix, $\hat{\mathbf{y}}=H\mathbf{y}$ az $\mathbf{y}$ vetülete az $X$ oszlopterére, a reziduum pedig az erre merőleges komponens. Az ortogonalitásból következik a hasznos ellenőrző tulajdonság:

$X'\mathbf{e}=\mathbf{0}.$

Mivel az $X$ első oszlopa csupa egyes, ebből különösen $\sum_i e_i=0$ is adódik: a reziduumok összege elméletileg pontosan nulla, ha van intercept a modellben.

A likelihood-kapcsolat Paraméter, becslő és likelihood szerint, ha a hibák független, azonos $N(0,\sigma^2)$ eloszlásúak, akkor a legkisebb négyzetes $\hat\beta$ egyben a maximum likelihood becslés is. A bizonytalanság alapja ugyanebben a modellben:

$\operatorname{Cov}(\hat\beta)=\sigma^2(X'X)^{-1}.$

A $\sigma^2$ torzítatlan becslése a tanítóhiba négyzete:

$s^2=\frac{RSS}{n-p},$

ahol $n-p$ a reziduális szabadsági fok. Ezt nevezzük mse_resid-nek a statsmodels-ban. A belőle vont $s$ négyzetgyök a tanítóhiba jellemző nagysága eurócentben; ez a minta belső illeszkedési hibája, nem az új rendelésre vonatkozó predikciós hiba Predikciós, paraméter- és Monte Carlo intervallum. Az együtthatók standard hibája:

$SE(\hat\beta_j)=s\sqrt{\bigl[(X'X)^{-1}\bigr]_{jj}},\qquad j=0,1,2,$

a $t$ statisztika pedig $\hat\beta_j/SE(\hat\beta_j)$ . A summary() táblázat éppen ezeket mutatja: coef, std err, t, P>|t| és a 95%-os konfidenciaintervallum oszlopai. A P>|t| a nullhipotézis ( $\beta_j=0$ ) alatti feltételes szélsőség Nullmodell és p-érték, nem annak valószínűsége, hogy az együttható nulla. A klasszikus summary standard hibái feltételezik a hibák függetlenségét, azonos varianciáját és normalitását; a Northstar rendeléseknél az ügyfélen belüli ismételt rendelések miatt a függetlenség sérülhet, ezért a summary standard hibáit később robusztus vagy klaszterezett változattal is ellenőrizni kell.

17.1.3 Teljes minimális munkapélda: NumPy és statsmodels ugyanazon az adaton

A példaadatok négy megfigyelést tartalmaznak:

basket_minor	quantity_sum	priority_flag
100	0	0
110	1	0
130	1	1
150	2	1

A statsmodels importjai és a formula-alapú illesztés:

import numpy as np
import pandas as pd
import statsmodels.formula.api as smf

data = pd.DataFrame({
    "basket_minor": [100.0, 110.0, 130.0, 150.0],
    "quantity_sum": [0.0, 1.0, 1.0, 2.0],
    "priority_flag": [0.0, 0.0, 1.0, 1.0],
})

fit_smf = smf.ols(
    "basket_minor ~ quantity_sum + priority_flag",
    data=data,
).fit()
print(fit_smf.params)

Ugyanez NumPy-val, a normálegyenletből:

y = data["basket_minor"].to_numpy(dtype=float)
X = np.column_stack((
    np.ones(len(data)),
    data["quantity_sum"].to_numpy(dtype=float),
    data["priority_flag"].to_numpy(dtype=float),
))
beta = np.linalg.solve(X.T @ X, X.T @ y)
yhat = X @ beta
resid = y - yhat
rss = resid @ resid
sigma2 = rss / (len(y) - X.shape[1])

print("beta:", beta)
print("resid:", resid)
print("sum(resid):", resid.sum())
print("X.T @ resid:", X.T @ resid)
print("rss:", rss, "sigma2:", sigma2, "sqrt:", np.sqrt(sigma2))

assert np.allclose(beta, fit_smf.params.to_numpy(), atol=1e-12)
assert np.allclose(resid, fit_smf.resid.to_numpy(), atol=1e-12)
assert abs(resid.sum()) < 1e-12
assert np.allclose(X.T @ resid, 0.0, atol=1e-12)
assert np.allclose(yhat + resid, y, atol=1e-12)

A beta elméletileg várt értéke a normálegyenlet megoldásából:

$X'X=\begin{pmatrix}4&4&2\\4&6&3\\2&3&2\end{pmatrix},\qquad X'\mathbf{y}=\begin{pmatrix}490\\540\\280\end{pmatrix},\qquad \hat\beta=\begin{pmatrix}97{,}5\\15\\20\end{pmatrix}.$

A kód futtatása nélkül is ellenőrizhető: a reziduumok $[2{,}5,-2{,}5,-2{,}5,2{,}5]$ , összegük $0$ , $X'\mathbf{e}=\mathbf{0}$ , $RSS=25$ , $n-p=4-3=1$ , tehát $s^2=25$ , azaz $s=5$ eurócent. A statsmodels summary() kivonata ezen a kis adaton elméletileg a következő szerkezetet adná (nem futtatott kimenet, a számok a fenti kézi levezetésből):

                         OLS Regression Results
==============================================================================
Dep. Variable:           basket_minor   R-squared:                       0.983
Model:                            OLS   Adj. R-squared:                  0.949
No. Observations:                   4   Df Residuals:                       1
Df Model:                           2
==============================================================================
                 coef    std err          t      P>|t|      [0.025      0.975]
------------------------------------------------------------------------------
Intercept     97.5000      4.330     22.516      0.028      42.497     152.503
quantity_sum  15.0000      5.000      3.000      0.205     -48.531      78.531
priority_flag 20.0000      7.071      2.828      0.216     -69.843     109.843
==============================================================================

A std err a fenti inverz diagonális elemeiből adódik; a meglehetősen széles konfidenciaintervallumok jól mutatják, hogy négy megfigyelés és egy reziduális szabadsági fok mellett a paraméter-bizonytalanság nagy. Részletes diagnosztikával később foglalkozunk; itt a kis minta eleve nem alkalmas például a normalitásteszt értelmezésére.

17.1.4 Eredmény értelmezése

A quantity_sum együtthatója 15 eurócent/darab: azonos prioritási státusz mellett egy további darab a modell szerint átlagosan 15 centtel magasabb kosárértékhez kapcsolódik. A priority_flag együtthatója 20 eurócent: azonos darabszám mellett a prioritásos rendelés várható kosara 20 centtel nagyobb. Az intercept 97,5 eurócent lenne a 0 darabos, nem prioritásos cellában; ezt a konkrét Northstar-adatokon extrapolációként, nem üzleti tényként kell kezelni.

A reziduumok összege nulla, és a vektor merőleges a dizájnmátrix minden oszlopára. A tanítóhiba $s=5$ eurócent: a modell a négy tanítórendelésen átlagosan ekkora tipikus eltéréssel illeszkedik. Az $R^2=0{,}983$ azt mondja, hogy a kosárérték varianciájának 98,3%-át magyarázza a két változó ezen a négy ponton. Az alacsony szabadsági fok miatt a p-értékek és a konfidenciaintervallumok inkább a mechanika bemutatására valók, mint érdemi következtetésre.

Numerikus stabilitás: oktatásban a normálegyenlet levezethető, de gyakorlati illesztéshez előnyös a np.linalg.lstsq(X, y, rcond=None) QR/SVD-alapú megoldása; az X.T @ X explicit invertálása rosszul kondicionált adatoknál rontja a pontosságot.

17.1.5 Feladatok

Értelmezés. A fenti mini-táblán a becsült modell $\hat y=97{,}5+15\cdot\text{quantity\_sum}+20\cdot\text{priority\_flag}$ . Fogalmazd meg saját szavakkal, mit jelent a 20-as együttható, és miért nem jelenti azt, hogy a prioritásos ügyfelek minden rendelése 20 centtel drágább. Melyik megfigyelt cellában érvényes az intercept?
Megvalósítás. Írj tiszta NumPy függvényt ols_by_hand(data, y_col, x_cols) néven, amely a kapott pandas DataFrame-ből felépíti a dizájnmátrixot intercept oszloppal, kiszámítja a normálegyenlet megoldását, és visszaadja a (beta, resid, sigma2) hármast. A függvény ne használja a statsmodels-t, és a reziduumot a y - X @ beta művelettel adja. Teszteld a mini-táblára, és ellenőrizd a sum(resid) és X.T @ resid invariánsokat.
Hibadiagnózis. Egy hallgató így akarja reprodukálni a formula API eredményét:
```
import statsmodels.api as sm
y = data["basket_minor"].to_numpy()
X = data[["quantity_sum", "priority_flag"]].to_numpy()
fit_bad = sm.OLS(y, X).fit()
```
Mi a különbség ez és a smf.ols("basket_minor ~ quantity_sum + priority_flag", data).fit() között? Milyen ellenőrzéssel vennéd észre a hibát a reziduumokon? Hogyan javítanád a sm.OLS verziót?
Megváltozott feltétel: Northstar-adat. Töltsd be a labor orders.csv állományát, és illeszd a
```
basket_minor ~ quantity_sum + distance_km + priority_flag
```
modellt a formula API-val. A CSV pontos elérési útja a saját projektgyökeredhez képest data/orders.csv; ha a labor nálad más almappában van, cseréld az útvonalat. A generált orders.csv 2400 sort tartalmaz, nincs benne hiányzó érték. Ellenőrizd, hogy a reziduális szabadsági fok $2400-4=2396$ , a reziduum hossza 2400, és az összege közel nulla. Értelmezd a distance_km együttható mértékegységét és asszociációs jellegét. A tényleges együttható a generált adattól függ; ezt a saját futtatásod adja, itt nem mechanikusan várt szám.

17.1.6 Megoldások és indoklás

1. feladat megoldása. A 20-as együttható a feltételes különbség ugyanazon quantity_sum mellett: ha két rendelés darabszáma egyenlő, de az egyik prioritásos, a másik nem, a modell szerint az előbbi várható kosara 20 eurócenttel nagyobb. Nem minden prioritásos rendelés 20 centtel drágább, mert a darabszám is hat a kosárra, és a reziduum egyedi zajt is tartalmaz. Az intercept a quantity_sum = 0, priority_flag = 0 cellához tartozik; a mini-táblában ez a 100-as megfigyelés cellája, a Northstar-adatokon viszont extrapolált cella.

2. feladat megoldása. A függvény placeholderek nélküli lehetséges változata:

def ols_by_hand(data, y_col, x_cols):
    y = data[y_col].to_numpy(dtype=float)
    X = np.column_stack(
        [np.ones(len(data))]
        + [data[c].to_numpy(dtype=float) for c in x_cols]
    )
    beta = np.linalg.solve(X.T @ X, X.T @ y)
    resid = y - X @ beta
    rss = resid @ resid
    sigma2 = rss / (len(y) - X.shape[1])
    return beta, resid, sigma2

Az ellenőrzés:

beta, resid, sigma2 = ols_by_hand(
    data, "basket_minor", ["quantity_sum", "priority_flag"]
)
assert abs(resid.sum()) < 1e-12
assert np.allclose(X.T @ resid, 0.0, atol=1e-12)  # X a függvényből kinyerhető
assert np.isclose(sigma2, 25.0)

A np.linalg.solve kis, jól kondicionált mátrixra alkalmas; nagyobb vagy majdnem kollineáris dizájnmátrixra a np.linalg.lstsq stabilabb választás.

3. feladat megoldása. Az sm.OLS osztály nem ad hozzá interceptet automatikusan, a formula API viszont igen. Ezért a hibás kód az origón átmenő síkot illeszt, a becsült együtthatók és reziduumok is eltérnek. Diagnózisnak kiváló a fit_bad.resid.sum(): intercept nélkül a reziduumok összege általában nem nulla, mert a dizájnmátrixból hiányzik az egyesekből álló oszlop, így az $X'\mathbf{e}=\mathbf{0}$ feltétel nem tartalmazza a $\sum e_i=0$ egyenletet. Javítás:

X = sm.add_constant(
    data[["quantity_sum", "priority_flag"]].to_numpy(dtype=float)
)
fit_good = sm.OLS(y, X).fit()

Ezzel a fit_good.params és a fit_good.resid a formula API eredményével numerikusan egyezik.

4. feladat megoldása. A betöltés és illesztés:

orders = pd.read_csv(
    "data/orders.csv",
    parse_dates=["order_ts"],
)
fit_northstar = smf.ols(
    "basket_minor ~ quantity_sum + distance_km + priority_flag",
    data=orders,
).fit()

assert len(orders) == 2400
assert fit_northstar.df_resid == 2396
assert fit_northstar.resid.shape == (2400,)
assert abs(fit_northstar.resid.sum()) < 1e-6  # eurócent, lebegőpontos tűrés
print(fit_northstar.params)

A df_resid azért 2396, mert négy paraméter van: intercept, quantity_sum, distance_km, priority_flag. A reziduumösszeg az intercept megléte miatt elméletileg nulla; a gyakorlatban a lebegőpontos kerekítések miatt csak közel nulla. A distance_km együtthatójának mértékegysége eurócent/kilométer: azt mutatja, mennyivel más várható kosárérték tartozik egy kilométerrel nagyobb távolsághoz, ha a darabszám és a prioritásjelző rögzített. Ez leíró asszociáció, nem annak bizonyítéka, hogy a távolság növelése okozza a kosárérték változását. A 180 napos megfigyelt rendelési ablak Sokaság, minta és adatfelvételi folyamat korlátja is érvényben marad.

17.1.7 Gyakori tévedés és ellenőrzés

A leggyakoribb hiba az API-választásból fakad: a felhasználó azt hiszi, hogy az intercept mindig „benne van” a modellben, pedig az sm.OLS(y, X) verzió nem ad hozzá. Az ellenkező irányú hiba is előfordul: ha pd.get_dummies-szal hozunk létre dummyváltozókat, alapértelmezésben minden kategória kap egy oszlopot, és ezeket egy intercept oszloppal együtt tesszük a mátrixba, akkor tökéletes kollinearitást kapunk. Például három régió esetén a három dummy összege minden sorban 1, ami megegyezik az intercept oszloppal; az $X'X$ mátrix szinguláris, a np.linalg.solve hibát dob, a statsmodels formulaillesztése viszont általában pszeudoinverzzel kezelheti a ranghiányos mátrixot, és nem feltétlenül jelez PatsyError hibát. A paraméterek ilyenkor nem azonosíthatók egyértelműen; ellenőrizni kell a rangot és a kondíciószámot. A javítás a referencia-kategória elhagyása vagy a formula C(region) jelölése. Önellenőrzés: a reziduum vektor hossza $n$ , nem $p$ ; intercept esetén az összege nulla; minden becsült együtthatóhoz tartozik mértékegység; a summary p-értéke a nullhipotézis alatti feltételes szélsőség, nem az együttható igazságának valószínűsége.

17.1.8 Rövid összefoglaló

Az OLS a reziduális négyzetösszeget minimalizálja, és a normálegyenlet $X'X\hat\beta=X'\mathbf{y}$ megoldását adja. A dizájnmátrix első oszlopa csupa egyes; így az intercept közönséges együttható. A statsmodels formula API automatikusan hozzáadja, a tiszta NumPy vagy sm.OLS megoldásnál explicit oszlop kell. A reziduum merőleges a dizájnmátrix oszlopaira, ezért az összege nulla. A summary a coef, standard hiba, t, p-érték és konfidenciaintervallum oszlopait az $(X'X)^{-1}$ és a $RSS/(n-p)$ tanítóhiba alapján számolja; a klasszikus standard hibák független, homoszkedasztikus, normális hibákat feltételeznek. A következő alfejezetben formula, interakció és kategorikus változó segítségével bővítjük ezt a modellt.