McNemar-próba

A McNemar-próba két összetartozó bináris változó eloszlását vizsgálja egy csoportban. A McNemar-próba lényegében egy \(\chi^2\) homogenitásvizsgálatot végez egy \(2 \times 2\)-es kontingencia táblázaton, ahol a sorok és oszlopok a két bináris változó egy-egy értékeinek felelnek meg. A próba célja megállapítani, hogy az összetartozó (matched/paired) mintán a két változónak ugyanolyan-e az eloszlása.

Legyen \((X, Y)\) kétdimenziós bináris háttérváltozók együttes eloszlása. A következő kontingenciatáblázat foglalja össze az \(n\) megfigyelés eredményét.

X \ Y	Y = 0	Y = 1	Sorösszeg
X = 0	a	b	a + b
X = 1	c	d	c + d
Oszlopösszeg	a + c	b + d	n

Hipotézisek

A nullhipotézis szerint a két marginális eloszlás (valószínűségek) megegyeznek mindkét kimenetelre, azaz \[p_a + p_b = p_a + p_c \quad \text{ és } \quad p_c + p_d = p_b + p_d\] azaz, \[ p_b = p_c \quad \text{ és } \quad p_c = p_b\] Itt a \(p_a\), \(p_b\), \(p_c\) és \(p_d\) az \(a\), \(b\), \(c\) és \(d\) cellákhoz tartozó elméleti/várható gyakoriságokat jelölik. Tehát a hipotézisek a következők:

\(H_0\): \(p_b = p_c\)
\(H_1\): \(p_b \neq p_c\)

A McNemar-próba próbastatisztikája a következő:

\[\chi^2 = \frac{(b-c)^2}{b+c} \overset{H_0}{\sim} \text{aszimp.}\ \chi^2\text{-eloszlás 1 szabadságfokkal}\]

Alkalmazások

A McNemar-próbát általában olyankor érdemes használni, ha egy kísérleti beavatkozás hatására vagyunk kíváncsiak: ugyanazokkal a személyekkel veszünk fel adatokat a beavatkozás előtt és után, majd megvizsgáljuk, hogy jelentősen változott-e a személyek állapota/válasza a beavatkozás hatására (vagyis különbözik-e az eloszlás). Másik tipikus helyzet, ha egy időpontban 2 vagy több adatot rögzítünk (pl. orvosi tesztet végzünk) ugyanazon személyeknél (kategorizált változók, ugyanazokkal a lehetséges értékekkel pl. pozitív és negatív), és ezek eloszlása érdekel bennünket.

R példa

A McNemar-próba elvégzéséhez a mcnemar.test függvényt fogjuk használni.

A következő táblázat két kérdőív eredményét tartalmazza, amiben az amerikai elnök hivatalban töltött teljesítményéről kérdezték a szavazókat 1 hónap eltéréssel. Lehetséges válaszok: helyeselik (approve) vagy helytelenítik (disapprove) az elnök tevékenységét.

## Agresti (1990), p. 350.
## Presidential Approval Ratings.
##  Approval of the President's performance in office in two surveys,
##  one month apart, for a random sample of 1600 voting-age Americans.
performance <-matrix(c(794, 86, 150, 570),
                     nrow = 2,
                     dimnames = list("1st Survey" = c("Approve", "Disapprove"),
                                     "2nd Survey" = c("Approve", "Disapprove")))
performance

##             2nd Survey
## 1st Survey   Approve Disapprove
##   Approve        794        150
##   Disapprove      86        570

mcnemar.test(performance)

## 
##  McNemar's Chi-squared test with continuity correction
## 
## data:  performance
## McNemar's chi-squared = 16.818, df = 1, p-value = 4.115e-05

A \(p\)-érték alapján látjuk hogy szignifikáns változás történt az értékelésekben. A táblázatból pedig az is kiolvasható hogy csökkent az elnök népszerűsége.

Hipotézisvizsgálat összefoglalás

Az alábbi táblázatban összegyűjtöttük és rendszereztük a tanult próbákat. A próbákhoz tartozó próbastatisztikák ezen a linken érhetőek el.

Cél	Próba neve	Paraméteres/nemparam	Leírás
Illeszkedésvizsgálat	1-mintás \(u\)-próba	paraméteres	Megállapítja, hogy egy normális eloszlású \(\sigma\) szórású minta várhatóértéke \(\mu\)-e
	1-mintás \(t\)-próba	paraméteres	Megállapítja, hogy egy normális eloszlású \(\sigma\) szórású minta várhatóértéke \(\mu\)-e ismeretlen (emiatt becsült) \(\sigma\)-ra
	Khí-négyzet próba (illeszkedésvizsgálatra)	nemparaméteres	Megállapítja, hogy az adathalmaz illeszkedik-e a megadott diszkrét eloszlásra
	Kolmogorov–Szmirnov próba	nemparaméteres	Folytonos háttérváltozó eloszlásának illeszkedését vizsgálja
	Shapiro–Wilk teszt	nemparaméteres	Azt vizsgálja, hogy a minta normális eloszlásból származik-e
	Binomiális teszt	paraméteres	Megállapítja, hogy egy Bernoulli változóra vett minta adott paraméterrel következik-e be

Homogenitásvizsgálat	2-mintás \(u\)-próba	paraméteres	Összehasonlítja két fgtlen normális eloszlású \(\sigma_1\), \(\sigma_2\) szórású minták átlagát
	2-mintás független \(t\)-próba	paraméteres	Két független normális eloszlású minta átlagát hasonlítja össze
	Összetartozó 2-mintás \(t\)-próba	paraméteres	Két összetartozó minta átlagát hasonlítja össze. Különbség normális eloszlású
	F-próba	paraméteres	Két független normális eloszlású minta szórását hasonlítja össze
	Mann–Whitney–Wilcoxon próba (fgtlen \(t\)-próba alternatíva)	nemparaméteres	Megállapítja, hogy két független minta azonos eloszlásból származik-e
	Wilcoxon (előjeles rang) próba (összet. \(t\)-próba alternatíva)	nemparaméteres	Megállapítja, hogy két összefüggő minta azonos eloszlásból származik-e
	Khí-négyzet próba (homogenitásvizsgálatra)	paraméteres	Több csoport egy közös kategorikus tulajdonságának eloszlását hasonlítja össze
	Binomiális teszt	nemparaméteres	Két független bináris mintára megállapítja, hogy a Bernoulli változók egyenlő valószínűséggel következnek-e be
	McNemar-próba	nemparaméteres	Összehasonlítja két összetartozó bináris változó eloszlását

Függetlenségvizsgálat	Khí-négyzet próba (függetlenségvizsgálatra)	paraméteres	Megállapítja, hogy két változó megfigyelései függetlenek-e egymástól

Lineáris regresszió és lineáris modell

Az elméleti modell

Az \(Y\) függő változót szeretnénk az \(X_1,\dots,X_p\) változók függvényével (\(f\)) közelíteni úgy, hogy a valódi \(Y\) és a becslés \(\hat{Y}=f(X_1,\dots,X_p)\) közötti várható négyzetes eltérés minimális legyen. A probléma elméletileg megoldott: a feltételes várható érték, mint függvény (változók a feltételben) minimalizálja a várható négyzetes eltérést: \[ \underset{f}{\operatorname{arg\,min\,}} \mathbb{E}\left(\big\lVert Y-f(X_1,\dots,X_p)\big\rVert_2\right) = \mathbb{E}\left(Y\mid X_1,\dots,X_p\right). \]

A lineáris regresszió annyiban különbözik az általános regressziótól, hogy csak lineáris függvények körében keressük az előbbi minimumot (azaz \(f(X_1,\dots,X_p) = a_0 + \sum_{i=1}^p a_i X_i\), és az \(a_i\) együtthatókat keressük). Elméletileg ez a probléma is megoldott, a korrelációk és várható értékek segítségével könnyen ki lehet számolni a lineáris függvény együtthatóit és konstansát: \[ \min_{f:f\ \text{lineáris fv.}} \mathbb{E}\left(\big\lVert Y-f(X_1,\dots,X_p)\big\rVert_2\right) = \min_{a_0,a_1,\dots,a_p} \mathbb{E}\left(\big\lVert Y-\left(a_0+a_1X_1+\dots+a_pX_p\right)\big\rVert_2\right), \] ahol \[ \underset{a_0,a_1,\dots,a_p}{\operatorname{arg\,min\,}} \left\lVert Y-\left(a_0+a_1X_1+\dots+a_pX_p\right)\right\rVert_2 = \begin{cases} \left[ \begin{matrix} a_1^*\\ \vdots\\ a_p^* \end{matrix} \right]= \left[ \begin{matrix} \vdots \\ \cdots \left\{Cov(X_j, X_k)\right\}_{j,k=1}^{p,p}\cdots\\ \vdots \end{matrix} \right]^{-1} \left[ \begin{matrix} Cov(X_1, Y)\\ \vdots\\ Cov(X_p, Y) \end{matrix} \right]\\ \quad a_0^*=\mathbb{E}Y- \left[a_1^*\mathbb{E}X_1+\dots +a_p^*\mathbb{E}X_p\right] \end{cases} \]

Az \(a_0\) az ú.n. konstans tag, vagy intercept/tengelymetszet. Nulla várható értékű változók esetén erre nincs szükség.

Ha \(p=1\), akkor az egyetlen együttható (\(a_1\)) egyszerűbb: \(a_1^*=\frac{Cov(X_1, Y)}{Var(X_1)}=r_{X_1Y}\frac{\mathbb{SD}(Y)}{\mathbb{SD}(X_1)}\)

Természetesen más a lineáristól eltérő függvényosztályokon is lehet vizsgálni a minimalizáló függvény kérdését.

Fontos tény, hogy ha a vizsgált változók (\(X_1,\dots,X_p,Y\)) együttes eloszlása normális, akkor a legjobb közelítő függvény lineáris, vagyis ugyanaz az általános és lineáris regressziós probléma megoldása.

Megj.: az elméleti regresszió sem ad mindig jó közelítést. Például ha \(Y\) független az \((X_1,\dots,X_p)\) változó vektortól, akkor a legjobb négyzetes közelítés az \(\mathbb{E}(Y)\).

A fenti képletben szereplő általában ismeretlen paramétereket (kovarianciák és várható értékek) valahogy becsülni kell, ha az elméleti modellt használjuk.

Ezt az elméleti modellt a Bolla-Krámli könyv 6.3. része (251.o.) járja körül.

A lineáris modell determinisztikus és random változata

A gyakorlatban nem ismerjük az elméleti eloszlást, csak az \((Y, X_1,\dots,X_p)\) vektorváltozóról vannak megfigyeléseink (\(n\) darab), és felteszzük, hogy a változók közötti összefüggés lineárisan modellezhető: \[ Y^{(i)}=a_0+a_1X_1^{(i)}+\dots+a_pX_p^{(i)}+\epsilon^{(i)},\qquad\forall i=1,\dots,n; \] ahol az \(\epsilon^{(i)}\) az i-dik megfigyelés hibája. A változók viszonya alapján két fő esetet különböztetünk most csak meg:

Determinisztikus esetben az \(X_1,\dots,X_p\) változók konkrét értékei ismertek, nem randomak. Ez például valamilyen ipari-mérnöki környezetben lehetséges, ahol a mérés paramétereit (\(X\)-ek értékeit) a mérnök beállítja, majd a beállított paraméterekkel megméri az \(Y\) célváltozót. Erről szól a Bolla-Krámli könyv 6.4. része (257.o.). Itt csak az \(Y\)-ok és a hibák randomak.
Random esetben esetben az \(X_1,\dots,X_p\) valószínűségi változók, és ezek valamint \(Y\) \(n\) db független realizációja ismert, egy random mintavételezés. Erről a változatról a Bolla-Krámli könyv konkrétan nem beszél. Sok más irodalom eléggé összemossa a két esetet. Egy jó könyv, ami jól elválasztja a két modellt pl. Greene Econometric Analysis c. könyve.

A megkülöböztetés azért fontos, mert nem minden eredmény valid mind a két alaphelyzetben, vagy csak más feltételek mellett.

A cél az elméleti együtthatók \((a_0,a_1,\dots,a_p)\) mintából való becslése, illetve az elméleti együtthatókhoz kötődő hipotézisvizsgálat. Az egyik fő módszer az együtthatók mintából való becslésére a legkisebb négyzetek módszerével határozza meg a minimalizálandó függvényt, amiből az alábbi képlet adódik egyszerű parciális deriválással: \[ \underset{a_0,a_1,\dots,a_p}{\operatorname{arg\,min\,}} \sum_{i=1}^n \left[ Y^{(i)}-\left(a_0+a_1X_1^{(i)}+\dots+a_pX_p^{(i)}\right)\right]^2 = \left[ \begin{matrix} \hat{a_0}\\ \hat{a_1}\\ \vdots\\ \hat{a_p} \end{matrix} \right]= \left[ \left[1\mid \mathbf{X}\right]^T \left[1\mid \mathbf{X}\right] \right]^{-1} \left[1\mid \mathbf{X}\right]^T \left[\mathbf{Y}\right], \] ahol: \[ \left[1\mid \mathbf{X}\right]:= \left[ \begin{matrix} 1\ &X_1^{(1)}\dots &X_p^{(1)}\\ \vdots &\vdots &\vdots\\ 1\ &X_1^{(n)}\dots &X_p^{(n)} \end{matrix} \right] \qquad\text{és}\qquad \left[\mathbf{Y}\right]:= \left[ \begin{matrix} Y^{(1)}\\ \vdots\\ Y^{(n)} \end{matrix} \right]. \]

A \(\left[1\mid \mathbf{X}\right]\) mátrixot szokták design mátrixnak is nevezni angol szakirodalomban a determinisztikus esetben. Szokás még szimplán \(\mathbf{X}\)-nek jelölni, de egyelőre a fenti jelölés beszédesebb. A csupa 1-es első oszlopra a konstans tag miatt van szükség. Ha konstans nélkül becslünk, erre az oszlopra sincs szükség.

A fenti ú.n. lineáris modellben az \(i\)-dik adatpont hibája: \[ \epsilon^{(i)}:=Y^{(i)}-\left(a_0+a_1X_1^{(i)}+\dots+a_pX_p^{(i)}\right) \] önmagában is egy valószínűségi változó. Minden randomitás benne van, amit az ismert változóink nem magyaráznak, és minden \(i\)-re más és más. Egy adott minta mellett az \(\mathbf{\epsilon}=(\epsilon_1,\dots,\epsilon_n)\) vektorváltozóra vonatkozó feltételeket szokás tenni a lineáris modellben. Determinisztikus esetben (ahol pl. \(\mathbf{X}\) konkrét beállított értékeket tartalmaz) a feltételek közvetlenül \(\mathbf{\epsilon}\)-ra vonatkoznak, míg random esetben (amikor pl. \(\mathbf{X}\) random mintavételezés eredménye) az \(\mathbf{\epsilon}|\mathbf{X}\) feltételes változóra vonatkoznak a feltevések:

a várható értéke 0, \(\mathbb{E}_{\mathbf{X}}(\epsilon_i)=0\quad \forall i\);
a varianciája koordinátánként azonos, \(Var_{\mathbf{X}}(\epsilon_i)=\sigma\quad \forall i\);
páronként feltételesen korrelálatlanok, \(Cov_{\mathbf{X}}(\epsilon_i,\epsilon_j)=0\quad \forall i,j\).

Szokás még feltenni, hogy a hibák legyenek együttesen normálisan eloszlásúak az iménti paraméterekkel, de ez a modellnek önmagában nem, csak a kiszámolt együtthatók becsléseire vonatkozó bizonyos hipotézisteszteknek a feltétele.

A mintára vonatkozó általános feltétel, hogy a \(\left[1\mid \mathbf{X}\right]\) mátrix legyen teljes rangú, így invertálható. Ez a determinisztikus esetben egyértelmű feltételezés, random esetben majdnem biztosan értendő. Ha ez nem teljesül egy adott mintára, azt a gyakorlatban általánosított inverz használatával meg szokták néha kerülni, ha a cél csak a predikció, de az együtthatókról levont következtetések többsége érvényét veszti ilyenkor.

A random változatban bizonyos eredmények feltétele lehet még a minta f.a.e tulajdonsága.

A hibákat sem fogjuk ismerni, ezekre is csak közelítéseket kapunk majd a becsült paraméterekkel: \[ \epsilon^{(i)}\approx Y^{(i)}-\widehat{Y}^{(i)}=Y^{(i)}-\left(\hat{a_0}+\hat{a_1}X_1^{(i)}+\dots+\hat{a_p}X_p^{(i)}\right). \]

\(R^2\) statisztika

Hogy az így kapott becslőmodell mennyire ad jó közelítést, azt pl. az \(R^2\) statisztika alapján dönthetjük el, amit az \(Y\) empirikus szórása, valamint az \(Y\) és a becsült értéke (\(\widehat{Y}\)) közötti eltérések segítségével számolunk ki: \[ R^2=1- \frac{\sum_{i=1}^n\left[Y^{(i)}-\widehat{Y}^{(i)}\right]^2} {\sum_{i=1}^n\left[Y^{(i)}-\overline{Y}\right]^2}. \] Ennek az értéke \(]-\infty, 1]\) közé esik:

1 az értéke, ha a becslőnk tökéletesen reprodukálja a mintát;
negatív, ha annyira rossz a becslőnk, hogy átlag is pontosabb becslést ad;
egyébként pedig 0 és 1 közötti értéket kapunk.

Minél nagyobb annál jobb a becslőnk a tesztelt mintán.

Lineáris regressziónál, abban az esetben, ha szerepel konstans a modellben, akkor az iménti \(R^2\) statisztika megegyezik az \(Y\) és a kapott legjobb lineáris közelítés (\(\widehat{Y}\)) empirikus korrelációjának a négyzetével (ez az ú.n. empirikus többszörös korrelációs együttható): \(R^2=\hat{r}^2_{Y\widehat{Y}}\).

Ha az elméleti modellt nézzük és egyetlen magyarázó változónk van akkor ez az \(R^2\) megegyezik az egyetlen magyarázó változó és a célváltozó empirikus korrelációjának négyzetével: \(R^2=\hat{r}^2_{YX}\). Mindez kvantitatíven alátámasztja azt, hogy a korreláció a lineáris összefüggőség mérőszáma.

R példa

Egy használt autók adatait tartalmazó adathalmazon fogjuk kipróbálni az eddigieket. Az ehhez tartalmazó Excel file megtalálható a labor anyagai között, itt a Kuiper data mappában, vagy Shonda Kuiper cikkének kiegészítő anyagai között.

Olvassuk be az adatot:

library(readxl)
carprices<-read_excel("kuiper.xlsx")

Az ebben található adatok random mintavételezés eredményének tekinthetőek a leírás alapján. Közelítsük a Price (ár) változót a

Mileage (megtett mérföldek száma)
Cylinder (hengerek száma)
Sound (fejlettebb hangrendszer)
Leather (bőrülés)

változók segítségével. Ezt az R lm() fv.-ével tudjuk megtenni, ami egy ú.n. formula objektumot vár bemenetként az alábbi módon. Természetesen alternatív megadás is lehetséges, de most maradjunk ennél.

Gyorsan nézzünk rá választott változóinkra:

head(carprices[,c('Price','Mileage','Cylinder','Sound','Leather')])

## # A tibble: 6 x 5
##    Price Mileage Cylinder Sound Leather
##    <dbl>   <dbl>    <dbl> <dbl>   <dbl>
## 1 17314.    8221        6     1       1
## 2 17542.    9135        6     1       0
## 3 16219.   13196        6     1       0
## 4 16337.   16342        6     0       0
## 5 16339.   19832        6     0       1
## 6 15709.   22236        6     1       0

A célváltozónk és a mérföldek száma folytonosnak tekinthető.

summary(carprices[,c('Price','Mileage','Cylinder','Sound','Leather')])

##      Price          Mileage         Cylinder         Sound       
##  Min.   : 8639   Min.   :  266   Min.   :4.000   Min.   :0.0000  
##  1st Qu.:14273   1st Qu.:14624   1st Qu.:4.000   1st Qu.:0.0000  
##  Median :18025   Median :20914   Median :6.000   Median :1.0000  
##  Mean   :21343   Mean   :19832   Mean   :5.269   Mean   :0.6791  
##  3rd Qu.:26717   3rd Qu.:25213   3rd Qu.:6.000   3rd Qu.:1.0000  
##  Max.   :70755   Max.   :50387   Max.   :8.000   Max.   :1.0000  
##     Leather      
##  Min.   :0.0000  
##  1st Qu.:0.0000  
##  Median :1.0000  
##  Mean   :0.7239  
##  3rd Qu.:1.0000  
##  Max.   :1.0000

A fenti összegzés alapján gyanús, hogy a többi változónk diszkrét.

unique(carprices$Cylinder)

## [1] 6 8 4

unique(carprices$Sound)

## [1] 1 0

unique(carprices$Leather)

## [1] 1 0

A hengerek száma (Cylinder) ordinális változó, míg a Sound és Leather valamilyen bináris változók.

Rakjuk össze a lineáris modellünket:

model1<-lm(Price ~ Mileage + Cylinder + Sound + Leather, data = carprices)
summary(model1)

## 
## Call:
## lm(formula = Price ~ Mileage + Cylinder + Sound + Leather, data = carprices)
## 
## Residuals:
##    Min     1Q Median     3Q    Max 
## -11677  -4893  -2737   2654  35647 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  3.268e+03  1.434e+03   2.279 0.022900 *  
## Mileage     -1.565e-01  3.409e-02  -4.590 5.16e-06 ***
## Cylinder     3.891e+03  2.030e+02  19.172  < 2e-16 ***
## Sound       -2.131e+03  6.095e+02  -3.496 0.000498 ***
## Leather      2.933e+03  6.354e+02   4.616 4.56e-06 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 7910 on 799 degrees of freedom
## Multiple R-squared:  0.3628, Adjusted R-squared:  0.3596 
## F-statistic: 113.7 on 4 and 799 DF,  p-value: < 2.2e-16

A modell objektumra hívott summary() függvénnyel különböző statisztikákat kapunk a modellünkről és arról, hogyan viselkedünk a mintánkon.

Fenn látjuk a reziduálisokra (a hibák mintára számolt értékei) vonatkozó statisztikákat. Ez alapján az ember valamennyire tud következtetni, hogy nagyjából milyen feltételek teljesülnek ezekre, és az alapján pedig, hogy milyen következtetéseket von le a későbbiekben. A hibák komolyabb elemzésével az ú.n. “Residual Analysis” témakör foglalkozik.

A táblázat központi részén látjuk a becsült együtthatókat soronként, és annak a hipotézisvizsgálatnak a \(p\)-értékét (Pr(>|t|)), amely azt vizsgálja, hogy az egyes elméleti együtthatók nullák-e.

Az \((a_0,a_1,\dots,a_p)\) “valódi” együtthatók a lineáris modellünkben a mintától függő változók, és erre kaptunk egy-egy becslést a legkisebb négyzetek techinkájával: \((\hat{a_0},\hat{a_1},\dots,\hat{a_p})\), az Estimate oszlopban. Ha feltesszük a mintára vonatkozó rangfeltételt (\(\left[1\mid \mathbf{X}\right]\) teljes rangú) és a hibák (\(\epsilon\)) normalitását a fenn taglalt paraméterekkel (\(\mathbf{\epsilon}\sim\mathcal{N}(\underline{0},\sigma^2I_n)\)), akkor belátható, hogy \[ \left[ \begin{matrix} \hat{a_0}\\\vdots\\\hat{a_p} \end{matrix} \right] \sim\mathcal{N}_{p+1} \left( \left[ \begin{matrix} a_0\\\vdots\\a_p \end{matrix} \right], \sigma^2 \left[Q_{\mathbf{X}}\right]^{-1} \right), \] ahol \(\left[Q_{\mathbf{X}}\right]= \left[1\mid \mathbf{X}\right]^T\left[1\mid \mathbf{X}\right]\); vagyis az együtthatók becslései együttesen normális eloszlásúak.

Megj.: Ebben a formában az iménti eredmény a determinisztikus változatban igaz, majd látjátok előadáson is. Eközben a random esetben a becslésnek az \(\mathbf{X}\)-re vett feltételes eloszlása lesz ilyen eloszlású, ezzel a fenn írt Greene könyv foglalkozik részleteiben.

Tehát a már tanult u/z-próbát tudnánk elvégezni minden koordinátára/együtthatóra, ha ismernénk \(\sigma\) értékét. Ezt többnyire csak becsülni tudjuk. Itt jön a képbe a táblázatban is szereplő Std. Error, ami a becslésünk ú.n. standard hibája változónként és ebben már \(\sigma\) becslése szerepel: \[ s_{\hat{a_j}}= \hat{\sigma}\sqrt{\left[Q_{\mathbf{X}}\right]^{-1}_{jj}}= \sqrt{ \frac{1}{n-(p+1)+2} \sum_{i=1}^n\left[Y^{(i)}-\widehat{Y}^{(i)}\right]^2} \sqrt{\left[Q_{\mathbf{X}}\right]^{-1}_{jj}} \] Ezt a közelítést használva t-próba használható, mivel belátható, hogy: \[ t=\frac{\hat{a_j}-a_j}{ s_{\hat{a_j}}}\sim t_{n-(p+1)}\text{-eloszlás}. \] Ez az eredmény a determinisztikus és a random változatban is él függetlenül a mintától. Tehát azt a (kétoldali) nullhipotézist vizsgáljuk, hogy \(a_j=0\), vagyis a számolt t value: \[ t=\frac{\hat{a_j}}{ s_{\hat{a_j}}}\underset{H_0}{\sim} t_{n-(p+1)}\text{-eloszlás}, \] és erre vonatkozik a sorvégi \(p\)-érték.

Lenn még több statisztikát kiír a függvény. Például láthatjuk a fenn taglalt \(R^2\) statisztikát, és ennek korrigált változatát is: \[ R_*^2=1-\frac{n-1}{n-(p+1)-1}\left(1-R^2\right), \] ami mindig kisebb a sima \(R^2\)-nél és azt “bünteti”, ha túl sok magyarázó változó alapján akarunk becsülni.

A Residual standard error a \(\sigma\) becslése tehát a fenti képletből a \(\hat{\sigma}\), és a szabadsági foka arra utal, hogy ez a becslő egy \(n-(p+1)\) szabadságú fokú \(\chi^2\)-eloszlással arányos.

Az alul található F-próba azt a nullhipotézist vizsgálja, hogy az összes együttható (a konstans kivéve) nulla egyszerre (\(a_1=\dots=a_p=0\)). Ennek részletei a determinisztikus esetet feltételezve a Bolla-Krámli könyv 263. oldala körül megtalálhatóak, de igaz a random esetben is.

Megj.: Nem a regressziós együtthatók nagysága a döntő abban, hogy melyik változó mennyire fontos, hiszen az a mértékegység függvénye. Ennek eldöntésében vagy a \(p\)-érték a mérvadó, vagy a regressziót standardizált változókra futtatjuk le, így bizonyos értelemben mértékegységtől függetlenül hasonlíthatjuk össze az egyes változók szerepét.

Ezt könnyen megtehetjük R-ben a scale() fv. segítségével, ami alapértelmezetten az átlaggal és a korr. emp. szórással standardizál. Ezt most alkalmazzuk a Mileage és a Cylinder változóra, mivel a másik két változó bináris (0-1 értékű).

model2<-lm(Price ~ scale(Mileage) + scale(Cylinder) + Sound + Leather, data = carprices)
summary(model2)

## 
## Call:
## lm(formula = Price ~ scale(Mileage) + scale(Cylinder) + Sound + 
##     Leather, data = carprices)
## 
## Residuals:
##    Min     1Q Median     3Q    Max 
## -11677  -4893  -2737   2654  35647 
## 
## Coefficients:
##                 Estimate Std. Error t value Pr(>|t|)    
## (Intercept)      20667.2      628.2  32.898  < 2e-16 ***
## scale(Mileage)   -1282.3      279.4  -4.590 5.16e-06 ***
## scale(Cylinder)   5399.2      281.6  19.172  < 2e-16 ***
## Sound            -2131.1      609.5  -3.496 0.000498 ***
## Leather           2933.0      635.4   4.616 4.56e-06 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 7910 on 799 degrees of freedom
## Multiple R-squared:  0.3628, Adjusted R-squared:  0.3596 
## F-statistic: 113.7 on 4 and 799 DF,  p-value: < 2.2e-16

A modell alapvetően nem változott meg, a változók szignifikanciájáét ez az átalakítás nem befolyásolta. A konstans tag vált még jobban szignifikánssá.

A Cylinder változót faktorrá is alakíthatjuk az as.factor() fv.-nyel. Így a Cylinder változó faktor lett, a lineáris regressziót végző parancs kicsit mást ad végeredményül. Technikailag az történt, hogy létrejött két új bináris változó (a Cylinder lehetséges értékei - 1 darab változó): Cylinder6 és Cylinder8. Ezek akkor 1-ek ha a Cylinder értéke 6 illetve 8 volt, egyébként 0-k. A regresszióba az eredeti Cylinder változó helyett ez a két új bináris változó került be. Azt mondhatjuk erre, hogy ez az átalakítás megengedi, hogy az eltérő Cylinder értékek esetén más legyen a lineáris regresszió konstans tagja.

model3<-lm(Price ~ scale(Mileage) + as.factor(Cylinder) + Sound + Leather, data = carprices)
summary(model3)

## 
## Call:
## lm(formula = Price ~ scale(Mileage) + as.factor(Cylinder) + Sound + 
##     Leather, data = carprices)
## 
## Residuals:
##    Min     1Q Median     3Q    Max 
## -16879  -4382  -2487   2405  29228 
## 
## Coefficients:
##                      Estimate Std. Error t value Pr(>|t|)    
## (Intercept)           18175.1      638.1  28.485  < 2e-16 ***
## scale(Mileage)        -1318.1      251.9  -5.234 2.13e-07 ***
## as.factor(Cylinder)6   2224.5      548.1   4.058 5.43e-05 ***
## as.factor(Cylinder)8  20683.0      822.7  25.142  < 2e-16 ***
## Sound                 -1030.4      555.3  -1.856   0.0639 .  
## Leather                 604.5      597.7   1.011   0.3121    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 7130 on 798 degrees of freedom
## Multiple R-squared:  0.4829, Adjusted R-squared:  0.4797 
## F-statistic: 149.1 on 5 and 798 DF,  p-value: < 2.2e-16

Láthatjuk, hogy ez \(R^2\)-ben sokkal jobb modellt eredményezett. Viszont a két eredeti bináris változónk elvesztette szignifikanciáját. Most ha elhagyjuk őket a modellből, tehát együtthatóikat kvázi nullára állítjuk:

model4<-lm(Price ~ scale(Mileage) + as.factor(Cylinder), data = carprices)
summary(model4)

## 
## Call:
## lm(formula = Price ~ scale(Mileage) + as.factor(Cylinder), data = carprices)
## 
## Residuals:
##    Min     1Q Median     3Q    Max 
## -16843  -4460  -2559   2555  28763 
## 
## Coefficients:
##                      Estimate Std. Error t value Pr(>|t|)    
## (Intercept)           17906.6      359.7  49.778  < 2e-16 ***
## scale(Mileage)        -1305.1      252.0  -5.178 2.84e-07 ***
## as.factor(Cylinder)6   2132.1      542.2   3.932 9.14e-05 ***
## as.factor(Cylinder)8  21020.6      799.5  26.293  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 7138 on 800 degrees of freedom
## Multiple R-squared:  0.4804, Adjusted R-squared:  0.4785 
## F-statistic: 246.6 on 3 and 800 DF,  p-value: < 2.2e-16

Nem kaptunk sokkal másabb vagy rosszabb modellt.

A következő heti anyagban még tovább fogjuk fejleszteni az elkezdett modellt több irányból is, de addig is otthon Ti is próbálkozzatok meg ezzel. Pl. a változókat lehet máshogy is transzformálni, ha pl. tudjuk, hogy az egyik változótól négyzetesen függ a cél.

Ábra

Az alábbi ábrán látjuk az \(x\)-tengelyen a Mileage változót az \(y\)-tengelyen a Price-t. A kék pontok az adatok, a piros vonal pedig az illesztett egyenes. Itt csak a Mileage alapján becsülünk.

plot(carprices$Mileage, carprices$Price, col = "blue",
     main = "Price plotted against Mileage", xlab = "Mileage", ylab = "Price")

abline(lm(Price ~ Mileage, data = carprices), col="red")

Statisztika 1 - Lab 11

Baranyi Máté & Nagy Marcell

2020 04 28