Als de gegevens representatief zijn voor de populatie kan men in de meervoudige lineaire regressiecontext eveneens aantonen dat de kleinste kwadraten schatters voor het intercept en de hellingen onvertekend zijn, m.a.w
\[E[\hat \beta_j]=\beta_j,\quad j=0,\ldots,p-1.\]Het feit dat de schatters gemiddeld (over een groot aantal vergelijkbare studies) niet afwijken van de waarden in de populatie, impliceert niet dat ze niet rond die waarde variëren. Om inzicht te krijgen hoe dicht we de parameterschatters bij het werkelijke intercept \(\beta_0\) en de werkelijke hellingen \(\beta_j\) mogen verwachten, wensen we bijgevolg ook haar variabiliteit te kennen.
Net zoals in Hoofdstuk 6 is het op basis van de puntschatters voor de hellingen niet duidelijk of de verbanden werkelijk voorkomen in de populatie of indien we de verbanden door toeval hebben geobserveerd in de dataset. De schatting van de hellingen is immers onnauwkeurig en zal variëren van steekproef tot steekproef. Zoals steeds is het resultaat van een data-analyse dus niet interpreteerbaar zonder die variabiliteit in kaart te brengen.
Om de resultaten uit de steekproef te kunnen veralgemenen naar de populatie zullen we in deze context eveneens inzicht nodig hebben op de verdeling van de parameterschatters. Om op basis van slechts één steekproef te kunnen voorspellen hoe de parameterschatters variëren van steekproef tot steekproef zullen we naast de onderstelling van
Lineariteit
bijkomende aannames moeten maken over de verdeling van de gegevens, met name
Onafhankelijkheid: de metingen \((X_{11},\dots, X_{1p-1}, Y_1), ..., (X_{n1},\ldots,X_{np-1},Y_n)\) werden gemaakt bij n onafhankelijke subjecten/observationele eenheden
Homoscedasticiteit of gelijkheid van variantie: de observaties variëren met een gelijke variantie rond het regressievlak. De residuen \(\epsilon_i\) hebben dus een gelijke variantie \(\sigma^2\) voor elk covariaat patroon \((X_1=x_1, ..., X_{p-1}=x_{p-1})\). Dat impliceert ook dat de conditionele variantie van \(Y\) gegeven \(X_1,\ldots,X_{p-1}\), \(\text{var}(Y\vert X_1,\ldots,X_{p-1})\) dus gelijk is, met name \(\text{var}(Y\vert X_1,\ldots,X_{p-1}) = \sigma^2\) voor elk covariaat patroon \((X_1=x_1, ..., X_{p-1}=x_{p-1})\). De constante \(\sigma\) wordt opnieuw de residuele standaarddeviatie genoemd.
Normaliteit: de residuen \(\epsilon_i\) zijn normaal verdeeld.
Uit aannames 2, 3 en 4 volgt dus dat de residuen \(\epsilon_i\) onafhankelijk zijn en dat ze allen eenzelfde Normale verdeling volgen
\[\epsilon_i \sim N(0,\sigma^2).\]Als we ook steunen op de veronderstelling van lineariteit weten we dat de originele observaties conditioneel op \(X_1,\ldots,X_{p-1}\) eveneens Normaal verdeeld zijn
\[Y_i\sim N(\beta_0+\beta_1 X_{i1}+\ldots+\beta_{p-1} X_{ip-1},\sigma^2),\]met een gemiddelde dat varieert in functie van de waarde van de onafhankelijke variabelen \(X_{i1},\ldots,X_{ip-1}\).
Merk op dat de onzekerheid op de hellingen af zal nemen wanneer er meer observaties zijn en/of wanneer de observaties meer gespreid zijn. Voor het opzetten van een experiment kan dit belangrijke informatie zijn. Uiteraard wordt de precisie ook beïnvloed door de grootte van de variabiliteit van de observaties rond het regressievlak, \(\sigma^2\), maar dat heeft een onderzoeker typisch niet in de hand.
De conditionele variantie (\(\sigma^2\)) is echter niet gekend en is noodzakelijk voor de berekening van de variantie op de parameterschatters. We kunnen \(\sigma^2\) echter opnieuw schatten op basis van de mean squared error (MSE):
\[\hat\sigma^2=MSE=\frac{\sum\limits_{i=1}^n \left(y_i-\hat\beta_0-\hat\beta_1 X_{i1}-\ldots-\hat\beta_{p-1} X_{ip-1}\right)^2}{n-p}=\frac{\sum\limits_{i=1}^n e^2_i}{n-p}.\]Analoog als in Hoofdstuk 6 kunnen we opnieuw toetsen en betrouwbaarheidsintervallen construeren op basis van de teststatistieken
\[T_k=\frac{\hat{\beta}_k-\beta_k}{SE(\hat{\beta}_k)} \text{ met } k=0, \ldots, p-1.\]Als aan alle aannames is voldaan dan volgen deze statistieken \(T_k\) een t-verdeling met \(n-p\) vrijheidsgraden. Wanneer niet is voldaan aan de veronderstelling van normaliteit maar wel aan lineariteit, onafhankelijkheid en homoscedasticiteit dan kunnen we voor inferentie opnieuw beroep doen op de centrale limietstelling die zegt dat de statistiek \(T_k\) bij benadering een standaard Normale verdeling zal volgen wanneer het aantal observaties voldoende groot is.
Voor het prostaatkanker voorbeeld kunnen we de effecten in de steekproef opnieuw veralgemenen naar de populatie toe door betrouwbaarheidsintervallen te bouwen voor de hellingen:
\[[\hat\beta_j - t_{n-p,\alpha/2} \text{SE}_{\hat\beta_j},\hat\beta_j + t_{n-p,\alpha/2} \text{SE}_{\hat\beta_j}]\]confint(lmVWS)
## 2.5 % 97.5 %
## (Intercept) -1.3473509 0.8112061
## lcavol 0.4033628 0.6999144
## lweight 0.2103288 0.8067430
## sviinvasion 0.2495824 1.0827342
Gezien nul niet in de intervallen ligt weten we eveneens dat de associaties tussen lpsa \(\leftrightarrow\) lcavol, lpsa \(\leftrightarrow\) lweight, lpsa \(\leftrightarrow\) svi, statistisch significant zijn op het 5% significantieniveau.
Anderzijds kunnen we ook formele hypothesetoetsen uitvoeren. Onder de nulhypothese veronderstellen we dat er geen associatie is tussen lpsa en de predictor \(X_j\):
\[H_0: \beta_j=0\]en onder de alternatieve hypothese is er een associatie tussen response en predictor \(X_j\):
\[H_1: \beta_j\neq0\]Met de test statistiek
\[T=\frac{\hat{\beta}_j-0}{SE(\hat{\beta}_j)}\]kunnen we de nulhypothese falsifiëren. Onder \(H_0\) volgt de statistiek een t-verdeling met \(n-p\) vrijheidsgraden, waarbij p het aantal model parameters is van het regressiemodel inclusief het intercept.
Deze tweezijdige testen zijn standaard geïmplementeerd in de standaard output van R.
summary(lmVWS)
##
## Call:
## lm(formula = lpsa ~ lcavol + lweight + svi, data = prostate)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1.72966 -0.45767 0.02814 0.46404 1.57012
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.26807 0.54350 -0.493 0.62301
## lcavol 0.55164 0.07467 7.388 6.3e-11 ***
## lweight 0.50854 0.15017 3.386 0.00104 **
## sviinvasion 0.66616 0.20978 3.176 0.00203 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.7168 on 93 degrees of freedom
## Multiple R-squared: 0.6264, Adjusted R-squared: 0.6144
## F-statistic: 51.99 on 3 and 93 DF, p-value: < 2.2e-16
De testen geven weer dat de associaties tussen lpsa\(\leftrightarrow\)lcavol, lpsa\(\leftrightarrow\)lweight en lpsa\(\leftrightarrow\)svi, respectievelijk extreem significant (\(p << 0.001\)) en sterk significant (\(p=0.001\) en \(p=0.002\)) zijn.