1 Cieľ a dáta

Zaujíma nás, čo vysvetľuje počet bodov, ktoré hráč NBA dá v jednom zápase. Analyzujeme 28 716 individuálnych výkonov 582 hráčov v období od 21.10.2025 do 13.06.2026. Jedno pozorovanie je jeden hráč v jednom zápase.

1.1 Ekonometrický model

Body hráča vysvetľujeme jeho odohratými minútami a počtom pokusov o skórovanie:

\[ PTS_i = \beta_0 + \beta_1\, MP_i + \beta_2\, FGA_i + \beta_3\, FTA_i + u_i \]

kde \(MP\) sú odohrané minúty, \(FGA\) pokusy o strelu z hry a \(FTA\) pokusy z trestných hodov. Očakávame kladné koeficienty pri všetkých troch premenných, pretože viac času a viac pokusov znamená viac príležitostí skórovať.

1.2 Popisné štatistiky

Zobraziť kód

vars <- c(PTS = "Body (PTS)", MP = "Minúty (MP)", FGA = "Pokusy z hry (FGA)",
          FTA = "Pokusy z trestných hodov (FTA)", X3PA = "Pokusy za 3 body (X3PA)")

desc <- data.frame(
  `Premenná`        = unname(vars),
  `Priemer`         = sapply(nba_model[names(vars)], mean),
  `Medián`          = sapply(nba_model[names(vars)], median),
  `Smer. odchýlka`  = sapply(nba_model[names(vars)], sd),
  `Min`             = sapply(nba_model[names(vars)], min),
  `Max`             = sapply(nba_model[names(vars)], max),
  check.names = FALSE
)
tab(desc, caption = "Popisné štatistiky použitých premenných", digits = 2,
    row.names = FALSE, align = "lrrrrr")

Popisné štatistiky použitých premenných
Premenná	Priemer	Medián	Smer. odchýlka	Max
Body (PTS)	10,60	9,0	8,68	83,0
Minúty (MP)	22,22	23,3	10,53	52,3
Pokusy z hry (FGA)	8,18	7,0	5,86	45,0
Pokusy z trestných hodov (FTA)	2,17	2,0	2,86	43,0
Pokusy za 3 body (X3PA)	3,39	3,0	3,02	22,0

Zobraziť kód

ggplot(nba_model, aes(x = FGA, y = PTS)) +
  geom_jitter(alpha = 0.12, size = 0.7, width = 0.25, height = 0.25, colour = "#1d428a") +
  geom_smooth(method = "lm", colour = "#c8102e", se = FALSE, linewidth = 1.2) +
  labs(x = "Pokusy o strelu z hry (FGA)", y = "Body (PTS)") +
  theme_minimal(base_size = 13)

Body a pokusy o strelu z hry: vzťah je silný a približne lineárny

2 Model 1: základný model

Zobraziť kód

model <- lm(PTS ~ MP + FGA + FTA, data = nba_model)

2.1 Odhadnuté koeficienty

Zobraziť kód

coef_tab(model, "Odhady regresných koeficientov – model 1")

Odhady regresných koeficientov – model 1
Premenná	Odhad	Štand. chyba	t-hodnota	p-hodnota	Signif.
Absolútny člen	-0,5608	0,0486	-11,5512	< 0,0001	***
MP – minúty	0,0324	0,0032	10,1929	< 0,0001	***
FGA – pokusy z hry	1,0640	0,0059	181,4972	< 0,0001	***
FTA – pokusy z trestných hodov	0,7998	0,0085	94,3086	< 0,0001	***

Signif.: *** p < 0,001; ** p < 0,01; * p < 0,05; . p < 0,1. Pre každý koeficient testujeme \(H_0:\beta_j=0\) proti \(H_1:\beta_j\neq 0\) a \(H_0\) zamietame, ak je p-hodnota menšia ako 0,05.

Čo to znamená ekonometricky

Pokusy z hry (FGA): ďalší pokus o strelu z hry zvyšuje očakávaný počet bodov o 1,06 bodu pri nezmenených minútach a trestných hodoch (štatisticky významný, p < 0,001). Koeficient možno čítať ako očakávanú „výnosnosť“ jedného pokusu.
Trestné hody (FTA): jeden ďalší pokus z trestných hodov zvyšuje očakávané body o 0,80 bodu (štatisticky významný, p < 0,001). Pre porovnanie, liga trafí 78,2 % trestných hodov, takže jeden pokus prináša priamo najviac 1 bod. Koeficient môže zachytávať aj iné javy, napríklad rolu hráča, ktorý často vyvoláva fauly.
Minúty (MP): každá ďalšia minúta zvyšuje očakávané body o 0,032 bodu (štatisticky významný, p < 0,001). Minúty pôsobia na body hlavne nepriamo, cez to, koľko pokusov hráč dostane, a koeficient preto vyjadruje len ich zostatkový, priamy efekt.
Absolútny člen: -0,56 bodu je predikcia pre hráča s nulovými minútami a nulovými pokusmi. Má len technický význam, pretože takéto pozorovanie v dátach nie je.

2.2 Kvalita modelu

Zobraziť kód

s <- summary(model); fst <- s$fstatistic
kv <- data.frame(
  `Ukazovateľ` = c("R²", "Upravené R²", "Reziduálna smer. odchýlka σ̂ (body)",
                   "F-štatistika", "p-hodnota F-testu", "AIC", "BIC", "Počet pozorovaní"),
  `Hodnota` = c(fn(s$r.squared, 4), fn(s$adj.r.squared, 4), fn(s$sigma, 2),
                fn(fst[1], 1), ifelse(fpval(model) < 0.0001, "< 0,0001", fn(fpval(model), 4)),
                fn(AIC(model), 1), fn(BIC(model), 1), fn(nobs(model), 0)),
  `Čo vyjadruje` = c("podiel variability bodov vysvetlený modelom",
                     "R² očistené o počet premenných",
                     "typická veľkosť chyby predikcie",
                     "spoločná významnosť všetkých smernicových koeficientov",
                     "pravdepodobnosť takejto F-štatistiky pri platnosti H₀",
                     "informačné kritérium (nižšie = lepšie)",
                     "informačné kritérium s prísnejšou penalizáciou",
                     "použité pozorovania"),
  check.names = FALSE)
tab(kv, caption = "Súhrnné charakteristiky modelu 1", row.names = FALSE, align = "llr")

Súhrnné charakteristiky modelu 1
Ukazovateľ	Hodnota	Čo vyjadruje
R²	0,8421	podiel variability bodov vysvetlený modelom
Upravené R²	0,8421	R² očistené o počet premenných
Reziduálna smer. odchýlka σ̂ (body)	3,45	typická veľkosť chyby predikcie
F-štatistika	51 047,1	spoločná významnosť všetkých smernicových koeficientov
p-hodnota F-testu	< 0,0001	pravdepodobnosť takejto F-štatistiky pri platnosti H₀
AIC	152 576,8	informačné kritérium (nižšie = lepšie)
BIC	152 618,1	informačné kritérium s prísnejšou penalizáciou
Počet pozorovaní	28 716	použité pozorovania

Čo to znamená ekonometricky

Vysvetľovacia schopnosť: model vysvetľuje 84,2 % variability bodov. Zvyšných 15,8 % pripadá na faktory mimo modelu, napríklad kvalitu obrany súpera, tempo zápasu, rolu hráča alebo jeho formu v daný deň. Upravené \(\bar R^2\) (0,8421) je takmer rovnaké, lebo máme veľa pozorovaní a len tri premenné.
Presnosť predikcie: skutočné body sa od predikcie typicky líšia o 3,45 bodu, pričom priemerný hráč dá v zápase 10,6 bodu.
F-test: testujeme \(H_0:\beta_1=\beta_2=\beta_3=0\). P-hodnota je p < 0,001, takže H₀ zamietame: model ako celok je štatisticky významný a aspoň jedna premenná body skutočne vysvetľuje.

3 Diagnostika modelu 1

Klasické t-testy a F-test sa opierajú o predpoklady lineárneho modelu: správnu špecifikáciu, konštantný rozptyl chýb a (pre presnú inferenciu) normalitu. Overujeme ich pomocou grafov a testov.

3.1 Grafická analýza rezíduí

Zobraziť kód

par(mfrow = c(2, 2))
plot(model, cex = 0.4, pch = 16, col = adjustcolor("#1d428a", alpha.f = 0.35))

Zobraziť kód

par(mfrow = c(1, 1))

Graf	Čo overuje	Varovný signál
Residuals vs Fitted	linearitu a strednú hodnotu chýb	zakrivená červená čiara
Q-Q plot	normalitu rezíduí	body mimo priamky, hlavne na koncoch
Scale-Location	konštantný rozptyl	rastúca červená čiara, tvar lievika
Residuals vs Leverage	vplyvné pozorovania	body za kontúrou Cookovej vzdialenosti 0,5 resp. 1

3.2 Konštantnosť rozptylu

Testujeme \(H_0\): rozptyl chýb je konštantný (homoskedasticita), pomocou car::ncvTest().

Hypotézu o konštantnom rozptyle zamietame (p < 0,001): rezíduá sú heteroskedastické. Odhady koeficientov zostávajú nevychýlené, no klasické štandardné chyby, t-testy a F-test nie sú spoľahlivé. Pri štatistickej inferencii by sme preto mali použiť heteroskedasticite-robustné štandardné chyby.

3.3 Normalita rezíduí: Jarqueov-Berov test

Test porovnáva šikmosť \(S\) a špicatosť \(K\) rezíduí s normálnym rozdelením (\(S=0\), \(K=3\)):

\[ JB=\frac{n}{6}\left(S^2+\frac{(K-3)^2}{4}\right)\ \overset{a}{\sim}\ \chi^2(2) \]

Normalitu rezíduí zamietame (JB = 2 042,8, p < 0,001). Šikmosť je 0,15 – kladná (pravostranne šikmé rozdelenie: občas sa objaví mimoriadne vysoké skóre, ktoré model nepredpovedal) – a špicatosť 4,27 (pri normálnom rozdelení 3), takže rezíduá majú ťažšie chvosty ako normálne rozdelenie. Pri takom veľkom počte pozorovaní (n = 28 716) však t-testy a F-test z hľadiska normality zostávajú približne platné vďaka asymptotickým vlastnostiam OLS. Normalita by bola kľúčová skôr pre presné predikčné intervaly.

3.4 Odľahlé a vplyvné pozorovania

Funkcia car::outlierTest() testuje \(H_0\), že pozorovanie nie je odľahlé, pomocou externe študentizovaných rezíduí a Bonferroniho korekcie.

Zobraziť kód

tab(out_tab, caption = "Najextrémnejšie pozorovania podľa outlierTest (model 1)",
    row.names = FALSE, align = "llrrrr")

Najextrémnejšie pozorovania podľa outlierTest (model 1)
Hráč	Dátum	Skutočné body	Predikované body	Študentizované rezíduum	Bonferroniho p
Aaron Gordon	23.10.2025	50	27,8	6,43	< 0,0001
Jamal Murray	03.12.2025	52	31,2	6,04	< 0,0001
Sam Merrill	11.02.2026	32	13,8	5,28	0,0038
Jamal Murray	25.03.2026	53	35,4	5,12	0,0089
Nikola Jokić	12.11.2025	55	37,8	4,99	0,0170
Cade Cunningham	10.11.2025	46	63,2	-4,99	0,0177
Jalen Green	21.02.2026	16	33,1	-4,96	0,0201
Immanuel Quickley	20.01.2026	40	23,2	4,89	0,0294
Stephen Curry	14.12.2025	48	31,5	4,80	0,0451

Čo to znamená ekonometricky

Po Bonferroniho korekcii je štatisticky odľahlých 9 z 28 716 pozorovaní (0,03 %). Najextrémnejším pozorovaním je Aaron Gordon (23.10.2025): dal 50 bodov, zatiaľ čo model predpovedal 27,8 bodu.
Odľahlé pozorovanie neznamená chybné pozorovanie. Ide o skutočné výnimočné zápasy, takže ich nemažeme, a skôr nás upozorňujú, že model nezachytáva všetko, čo body ovplyvňuje.
Najväčšia Cookova vzdialenosť je 0,016, teda hlboko pod hranicou 1: žiadne pozorovanie neprimerane neovplyvňuje odhady koeficientov.

4 Hat matica, leverage a vplyvné pozorovania

Graf Residuals vs Leverage z predchádzajúcej kapitoly vychádza z tzv. hat matice. Ukážeme, ako vzniká, a použijeme ju na zistenie, ktoré zápasy najviac ovplyvňujú odhadnutý model.

4.1 Hat matica a leverage

Odhadnuté hodnoty v lineárnom modeli získame z pozorovaných hodnôt vynásobením maticou \(\mathbf H\):

\[ \widehat{\mathbf y}=\mathbf X\widehat{\boldsymbol\beta} =\underbrace{\mathbf X(\mathbf X^\mathsf T\mathbf X)^{-1}\mathbf X^\mathsf T}_{\mathbf H}\,\mathbf y=\mathbf H\mathbf y \]

Matica \(\mathbf H\) je symetrická a idempotentná (\(\mathbf H^2=\mathbf H\)). Jej diagonálne prvky \(h_{ii}\) sa nazývajú leverage (pákový účinok) a platí

\[ 0\le h_{ii}\le 1,\qquad \sum_{i=1}^n h_{ii}=p,\qquad \bar h=\frac{p}{n}, \]

kde \(p\) je počet parametrov vrátane konštanty. Leverage závisí iba od matice \(\mathbf X\), nie od hodnôt \(y\). Vysoká leverage teda znamená neobvyklú kombináciu vysvetľujúcich premenných a sama osebe neznamená, že pozorovanie je chybné alebo vplyvné.

Pri \(n=\) 28 716 pozorovaniach by mala celá matica \(\mathbf H\) približne 6,6 GB, preto ju nepočítame celú. Diagonálu vieme získať priamo zo vzorca \(h_{ii}=\mathbf x_i^\mathsf T(\mathbf X^\mathsf T\mathbf X)^{-1}\mathbf x_i\) a overíme ju proti hatvalues().

Zobraziť kód

X <- model.matrix(model)
h_rucne <- rowSums((X %*% solve(crossprod(X))) * X)
h <- hatvalues(model)
stopifnot(isTRUE(all.equal(unname(h_rucne), unname(h))))

k_par <- length(coef(model)); n_obs <- nobs(model)
n_h2 <- sum(h > 2 * k_par / n_obs)
n_h3 <- sum(h > 3 * k_par / n_obs)

lev_tab <- data.frame(
  `Ukazovateľ` = c("Počet parametrov p", "Súčet leverage (stopa H)", "Priemerná leverage p/n",
                   "Minimálna leverage", "Maximálna leverage",
                   "Orientačný prah 2p/n", "Pozorovania nad 2p/n",
                   "Orientačný prah 3p/n", "Pozorovania nad 3p/n"),
  `Hodnota` = c(k_par, fn(sum(h), 2), fn(k_par / n_obs, 5), fn(min(h), 5), fn(max(h), 4),
                fn(2 * k_par / n_obs, 5), n_h2, fn(3 * k_par / n_obs, 5), n_h3),
  check.names = FALSE)
tab(lev_tab, caption = "Leverage v modeli 1", row.names = FALSE, align = "lr")

Leverage v modeli 1
Ukazovateľ	Hodnota
Počet parametrov p	4
Súčet leverage (stopa H)	4,00
Priemerná leverage p/n	0,00014
Minimálna leverage	3e-05
Maximálna leverage	0,0088
Orientačný prah 2p/n	0,00028
Pozorovania nad 2p/n	2164
Orientačný prah 3p/n	0,00042
Pozorovania nad 3p/n	914

Vlastnosti samotnej matice \(\mathbf H\) si overíme na náhodnej vzorke 200 pozorovaní:

Zobraziť kód

set.seed(123)
vz <- sample(nrow(nba_model), 200)
Xs <- model.matrix(~ MP + FGA + FTA, data = nba_model[vz, ])
Hs <- Xs %*% solve(crossprod(Xs)) %*% t(Xs)

vl <- data.frame(
  `Vlastnosť` = c("Symetria: H = Hᵀ", "Idempotentnosť: H² = H", "Stopa H (= počet parametrov p)"),
  `Výsledok`  = c(ifelse(isTRUE(all.equal(Hs, t(Hs))), "splnená", "nesplnená"),
                  ifelse(isTRUE(all.equal(Hs %*% Hs, Hs)), "splnená", "nesplnená"),
                  fn(sum(diag(Hs)), 2)),
  check.names = FALSE)
tab(vl, caption = "Vlastnosti hat matice (vzorka 200 pozorovaní)", row.names = FALSE, align = "lr")

Vlastnosti hat matice (vzorka 200 pozorovaní)
Vlastnosť	Výsledok
Symetria: H = Hᵀ	splnená
Idempotentnosť: H² = H	splnená
Stopa H (= počet parametrov p)	4,00

Čo to znamená ekonometricky

Priemerná leverage je \(p/n=\) 0,00014. Nad orientačným prahom \(2p/n\) je 2164 pozorovaní (7,54 %) a nad \(3p/n\) 914. Prahy sú len pomôcka, nie kritické hodnoty.
Najvyššiu leverage (0,0088) má Bam Adebayo (10.03.2026): odohral 41,9 minút a mal 43 pokusov z hry. Ide o neobvyklú kombináciu vysvetľujúcich premenných voči ostatným výkonom, a nie o chybu.

4.2 Štandardizované rezíduá

Rezíduá nemajú rovnaký rozptyl, platí \(\operatorname{Var}(e_i\mid\mathbf X)=\sigma^2(1-h_{ii})\). Pri porovnávaní preto používame interne štandardizované rezíduá:

\[ r_i=\frac{e_i}{\widehat\sigma\sqrt{1-h_{ii}}} \]

Pri normálnych chybách by \(|r_i|>2\) malo nastať u približne 4,6 % pozorovaní a \(|r_i|>3\) u približne 0,3 %.

Zobraziť kód

rs <- rstandard(model)
sh2 <- mean(abs(rs) > 2); sh3 <- mean(abs(rs) > 3)

sr <- data.frame(
  `Hranica` = c("|r| > 2", "|r| > 3"),
  `Počet pozorovaní` = c(sum(abs(rs) > 2), sum(abs(rs) > 3)),
  `Podiel v dátach` = paste0(fn(100 * c(sh2, sh3), 2), " %"),
  `Pri normálnom rozdelení` = paste0(fn(100 * 2 * pnorm(-c(2, 3)), 2), " %"),
  check.names = FALSE)
tab(sr, caption = "Štandardizované rezíduá: porovnanie s normálnym rozdelením",
    row.names = FALSE, align = "lrrr")

Štandardizované rezíduá: porovnanie s normálnym rozdelením
Hranica	Počet pozorovaní	Podiel v dátach	Pri normálnom rozdelení
\|r\| > 2	1 587	5,53 %	4,55 %
\|r\| > 3	245	0,85 %	0,27 %

Čo to znamená ekonometricky

Mimo hranice \(\pm 3\) je 0,85 % pozorovaní, pričom pri normálnom rozdelení by ich bolo 0,27 %. Extrémne rezíduá sú teda častejšie, než pripúšťa normalita, čo zodpovedá ťažkým chvostom zisteným Jarqueovým-Berovým testom: model občas výrazne podstrelí skóre mimoriadne vydareného zápasu. Rezíduum mimo \(\pm 2\) ani \(\pm 3\) nie je dôvodom pozorovanie vyradiť.

4.3 Cookova vzdialenosť

Cookova vzdialenosť spája veľkosť rezídua a leverage a meria, ako výrazne by sa zmenili odhady po vynechaní pozorovania:

\[ D_i=\frac{e_i^2}{p\,\mathrm{MSE}}\cdot\frac{h_{ii}}{(1-h_{ii})^2}=\frac{r_i^2}{p}\cdot\frac{h_{ii}}{1-h_{ii}} \]

Orientačne sa používajú hranice \(4/n\), 0,5 a 1. Nejde o pevné pravidlá na vyraďovanie údajov.

Zobraziť kód

mse <- sum(residuals(model)^2) / df.residual(model)
D_rucne <- residuals(model)^2 / (k_par * mse) * h / (1 - h)^2
cd <- cooks.distance(model)
stopifnot(isTRUE(all.equal(unname(D_rucne), unname(cd))))

top <- order(cd, decreasing = TRUE)[1:5]
top_tab <- data.frame(
  `Hráč` = nba_model$Player[top],
  `Dátum` = format(nba_model$Date[top], "%d.%m.%Y"),
  `MP` = fn(nba_model$MP[top], 1), `FGA` = nba_model$FGA[top],
  `FTA` = nba_model$FTA[top], `PTS` = nba_model$PTS[top],
  `Leverage` = fn(h[top], 4),
  `Štand. rezíduum` = fn(rs[top], 2),
  `Cookova vzdialenosť` = fn(cd[top], 4),
  check.names = FALSE)
tab(top_tab, caption = "Päť pozorovaní s najvyššou Cookovou vzdialenosťou",
    row.names = FALSE, align = "llrrrrrrr")

Päť pozorovaní s najvyššou Cookovou vzdialenosťou
Hráč	Dátum	MP	FGA	FTA	PTS	Leverage	Štand. rezíduum	Cookova vzdialenosť
Cade Cunningham	10.11.2025	44,6	45	18	46	0,0025	-4,98	0,0158
Nikola Jokić	25.12.2025	42,6	21	23	56	0,0020	4,19	0,0090
Nikola Jokić	12.11.2025	33,6	23	16	55	0,0010	4,99	0,0063
Luka Dončić	19.03.2026	37,6	30	19	60	0,0016	3,55	0,0050
James Harden	22.11.2025	35,2	26	14	55	0,0009	4,51	0,0044

Zobraziť kód

par(mfrow = c(1, 2))
plot(h, type = "h", col = adjustcolor("#1d428a", alpha.f = 0.5),
     xlab = "Pozorovanie", ylab = "Leverage", main = "Pákový účinok")
abline(h = 2 * k_par / n_obs, lty = 2, col = "#c8102e", lwd = 2)
plot(cd, type = "h", col = adjustcolor("#1d428a", alpha.f = 0.5),
     xlab = "Pozorovanie", ylab = "Cookova vzdialenosť", main = "Vplyv pozorovaní")
abline(h = 4 / n_obs, lty = 2, col = "#c8102e", lwd = 2)

Leverage a Cookova vzdialenosť jednotlivých pozorovaní (čiarkovaná čiara je orientačný prah)

Zobraziť kód

par(mfrow = c(1, 1))

4.3.1 Je vysoká leverage to isté ako vplyvné pozorovanie?

Nie. Leverage závisí len od \(\mathbf X\), rezíduum hovorí o odchýlke \(y_i\) od modelu a Cookova vzdialenosť zohľadňuje obe. Rozdiel najlepšie ukáže priame porovnanie: model odhadneme znova bez 10 pozorovaní s najvyššou Cookovou vzdialenosťou.

Zobraziť kód

vyrad <- order(cd, decreasing = TRUE)[1:10]
model_bez <- lm(PTS ~ MP + FGA + FTA, data = nba_model[-vyrad, ])
zm <- 100 * (coef(model_bez) / coef(model) - 1)

cmp <- data.frame(
  `Parameter` = unname(lab[names(coef(model))]),
  `Všetky pozorovania` = fn(coef(model), 4),
  `Bez 10 najvplyvnejších` = fn(coef(model_bez), 4),
  `Zmena v %` = fn(zm, 2),
  check.names = FALSE)
tab(cmp, caption = "Koeficienty s najvplyvnejšími pozorovaniami a bez nich",
    row.names = FALSE, align = "lrrr")

Koeficienty s najvplyvnejšími pozorovaniami a bez nich
Parameter	Všetky pozorovania	Bez 10 najvplyvnejších	Zmena v %
Absolútny člen	-0,5608	-0,5600	-0,15
MP – minúty	0,0324	0,0330	1,67
FGA – pokusy z hry	1,0640	1,0633	-0,07
FTA – pokusy z trestných hodov	0,7998	0,7953	-0,56

Čo to znamená ekonometricky

Najväčšia Cookova vzdialenosť je 0,016, teda hlboko pod hranicou 1: žiadne pozorovanie neprimerane neovplyvňuje odhady koeficientov. Pozorovaní nad hranicou \(4/n\) je 1859 (6,47 %).
Pozorovanie s najvyššou Cookovou vzdialenosťou (Cade Cunningham, 10.11.2025) má leverage 0,0025 a štandardizované rezíduum -4,98: vplyv vzniká kombináciou nezvyčajných hodnôt vysvetľujúcich premenných a veľkej chyby predikcie.
Po vynechaní 10 najvplyvnejších pozorovaní sa smernicové koeficienty zmenili najviac o 1,67 % (pri premennej MP – minúty). Odhady sú teda stabilné: žiadne jednotlivé pozorovanie neurčuje výsledky modelu.

5 Model 2: rozšírenie a porovnanie modelov

Do modelu pridáme pokusy za tri body (\(X3PA\)) a v modeli 3 aj kvadratický člen minút (\(MP^2\)), ktorý pripúšťa, že dodatočná minúta má iný efekt pri krátkom a dlhom čase na palubovke:

\[ PTS_i = \beta_0 + \beta_1 MP_i + \beta_2 FGA_i + \beta_3 FTA_i + \beta_4 X3PA_i + u_i \]

Zobraziť kód

model2 <- lm(PTS ~ MP + FGA + FTA + X3PA, data = nba_model)
model3 <- lm(PTS ~ MP + I(MP^2) + FGA + FTA + X3PA, data = nba_model)

Zobraziť kód

coef_tab(model2, "Odhady regresných koeficientov – model 2")

Odhady regresných koeficientov – model 2
Premenná	Odhad	Štand. chyba	t-hodnota	p-hodnota	Signif.
Absolútny člen	-0,5802	0,0484	-11,9781	< 0,0001	***
MP – minúty	0,0298	0,0032	9,3864	< 0,0001	***
FGA – pokusy z hry	1,0144	0,0070	144,4276	< 0,0001	***
FTA – pokusy z trestných hodov	0,8249	0,0087	94,9959	< 0,0001	***
X3PA – pokusy za 3 body	0,1263	0,0099	12,7411	< 0,0001	***

Čo to znamená ekonometricky

Pokusy za tri body (X3PA): koeficient 0,13 (štatisticky významný, p < 0,001) vyjadruje dodatočný efekt trojkového pokusu oproti dvojkovému, keďže FGA už trojky zahŕňa. Koeficient je kladný: pri danom celkovom počte pokusov je trojkový pokus v priemere výnosnejší než dvojkový.
Zmena ostatných odhadov: koeficient pri FGA sa po pridaní X3PA zmenil z 1,06 na 1,01. Pri korelovaných vysvetľujúcich premenných sa odhady posúvajú, lebo sa medzi ne inak rozdelí spoločný vplyv.
Multikolinearita: najvyšší VIF je 4,11 (FGA), teda multikolinearita nie je problém.
Štatistická vs. ekonomická významnosť: pri 28 716 pozorovaniach vyjde významným aj veľmi malý efekt. Pri hodnotení preto posudzujeme nielen p-hodnotu, ale aj veľkosť koeficientu v bodoch.

5.1 Porovnanie modelov

Zobraziť kód

m_row <- function(m, nm, d) {
  s <- summary(m)
  data.frame(Model = nm,
             `R²` = fn(s$r.squared, 4), `Upravené R²` = fn(s$adj.r.squared, 4),
             `σ̂ (body)` = fn(s$sigma, 3), AIC = fn(AIC(m), 1), BIC = fn(BIC(m), 1),
             `JB p-hodnota` = fpv(d$jb_p), `Heterosk. p-hodnota` = fpv(d$ncv_p),
             `Odľahlé (Bonf.)` = d$n_out, check.names = FALSE)
}
porov <- rbind(
  m_row(model,  "Model 1: MP + FGA + FTA", ds1),
  m_row(model2, "Model 2: + X3PA", ds2),
  m_row(model3, "Model 3: + X3PA + MP²", ds3))
tab(porov, caption = "Porovnanie modelov 1 až 3", row.names = FALSE, align = "lrrrrrrrr")

Porovnanie modelov 1 až 3
Model	R²	Upravené R²	σ̂ (body)	AIC	BIC	JB p-hodnota	Heterosk. p-hodnota	Odľahlé (Bonf.)
Model 1: MP + FGA + FTA	0,8421	0,8421	3,447	152 576,8	152 618,1	< 0,001	< 0,001	9
Model 2: + X3PA	0,8430	0,8430	3,438	152 416,9	152 466,5	< 0,001	< 0,001	7
Model 3: + X3PA + MP²	0,8430	0,8430	3,437	152 410,2	152 468,0	< 0,001	< 0,001	8

Výber modelu

Podľa AIC je najlepší Model 3 a podľa BIC Model 2. Rozdiel AIC oproti modelu 1 je 166,6 (rozdiel nad 10 sa považuje za silný dôkaz v prospech modelu s nižším AIC). Kvadratický člen minút je štatisticky významný (p = 0,003), takže vplyv minút na body nie je striktne lineárny. Heteroskedasticita (p-hodnota testu < 0,001) a odchýlky od normality pretrvávajú vo všetkých troch modeloch, takže žiadna zo zmien špecifikácie ich neodstráni.

Diagnostické grafy modelu 2 (kliknutím rozbaliť)

Zobraziť kód

par(mfrow = c(2, 2))
plot(model2, cex = 0.4, pch = 16, col = adjustcolor("#1d428a", alpha.f = 0.35))

Zobraziť kód

par(mfrow = c(1, 1))

6 Uloženie workspace

Aby sme nemuseli v ďalšom cvičení znovu importovať a čistiť dáta, uložíme celý workspace a vieme ho kedykoľvek obnoviť.

aaa <- 6; bbb <- "text"
rm(list = c("aaa", "bbb"))     # zmazanie vybraných objektov
save.image("pamat.Rdata")      # uloženie celého workspace
rm(list = ls())                # vyčistenie celej pamäte
load("pamat.Rdata")            # opätovné načítanie uloženého workspace
cat("Po načítaní je v pamäti", length(ls()), "objektov.\n")

Po načítaní je v pamäti 67 objektov.

7 Záver

Čo určuje body: model 1 vysvetľuje 84,2 % variability bodov v zápase. Každý ďalší pokus z hry zvyšuje očakávané body o 1,06 a každý pokus z trestných hodov o 0,80. Minúty pôsobia hlavne cez počet pokusov.
Spoľahlivosť inferencie: test odhalil heteroskedasticitu, preto by sa t-testy mali opierať o robustné štandardné chyby. Normalitu rezíduí zamietame, no pri veľkej vzorke to inferenciu zásadne neohrozuje.
Odľahlé hodnoty: 9 pozorovaní je po Bonferroniho korekcii odľahlých. Ide o výnimočné zápasy, nie o chyby v dátach, a preto ich v dátach ponechávame.
Najlepší model podľa AIC: Model 3.
Obmedzenie: dáta obsahujú opakované pozorovania tých istých hráčov, takže predpoklad nezávislosti chýb nemusí byť úplne splnený. Možným rozšírením sú hráčske fixné efekty, robustné štandardné chyby alebo premenné o súperovi a tempe hry.

8 Príloha: ako čítať normálny Q-Q graf

Štyri typické tvary rozdelenia (kliknutím rozbaliť)

Nasledujúce príklady sú simulované dáta a nesúvisia priamo s NBA. Ukazujú, aký obraz vznikne v Q-Q grafe pri rôznych tvaroch rozdelenia. Pri normálnom rozdelení by body približne sledovali červenú priamku.

Zobraziť kód

set.seed(123)
ukazky <- list(
  "Ľavostranná šikmosť"  = -rchisq(1000, df = 5),
  "Pravostranná šikmosť" = rchisq(1000, df = 5),
  "Ťažké chvosty"        = rt(1000, df = 2),
  "Ľahké chvosty"        = runif(1000, -3, 3))
farby <- c("skyblue", "orange", "lightgreen", "lightblue")

par(mfrow = c(2, 4), mar = c(4, 4, 3, 1))
for (i in seq_along(ukazky))
  hist(ukazky[[i]], breaks = 30, main = names(ukazky)[i], xlab = "Hodnota", col = farby[i])
for (i in seq_along(ukazky)) {
  qqnorm(ukazky[[i]], main = "Q-Q graf")
  qqline(ukazky[[i]], col = "red", lwd = 2)
}

Histogramy (hore) a normálne Q-Q grafy (dole) pre štyri tvary rozdelenia

Zobraziť kód

par(mfrow = c(1, 1))

Tvar rozdelenia	Obraz v Q-Q grafe
ľavostranná šikmosť	oblúk v tvare kopca, ľavý koniec výrazne pod priamkou
pravostranná šikmosť	oblúk v tvare úsmevu, pravý koniec výrazne nad priamkou
ťažké chvosty	tvar písmena S: ľavý koniec pod priamkou, pravý nad ňou
ľahké chvosty	obrátené S: ľavý koniec nad priamkou, pravý pod ňou

Rozdelenie \(t\) s dvomi stupňami voľnosti nemá konečnú kurtózu, preto tu hovoríme len o ťažkých chvostoch. Rezíduá z bodov v zápase bývajú pravostranne šikmé, lebo body sú nezáporné a občas padne mimoriadne vysoké skóre.