Kovarianz, Korrelation, (lineare) Regression

Die Kovarianz und die Korrelation sind zwei (verwandte) Maße für die Stärke der Beziehung zwischen 2 numerischen Variablen.

head(epg)

##        F1   F2       COG        DP SUM1278 V
## 1 391.213 1640 1.0555556 0.3750000       8 a
## 2 760.000 1820 1.3823529 0.6666667      12 a
## 3 660.000 1720 1.2857143 0.5416667      11 a
## 4 487.437 1740 0.7857143 0.2916667       6 a
## 5 300.036 2280 1.9736842 0.5833333      16 E
## 6 311.557 2280 1.9000000 0.6250000      16 E

epg enthält F1- und F2-Werte zum Vokaloffset für V/k/-Folgen (mit V = /a ɛ ɪ i ɔ ʊ/), gesprochen von einem deutschen Muttersprachler.
COG (“centre of gravity”) ist ein gemessener Parameter aus der Bereich der Elektropalatographie. Er gibt ein Maß für den Schwerpunkt der Kontakte auf einem 8 x 8 Kontaktflächen großen künstlichem Gaumen.
SUM1278 ist die Summe der aktivierten Kontaktflächen auf den jeweils zwei äußersten Kontaktflächenspalten des künstlichen Gaumens (da es acht Spalten sind, handelt es sich um die Spalten 1, 2, 7, und 8).

1. Kovarianz

Die Kovarianz ist ein Maß für die Assoziation zwischen zwei Variablen. Je mehr die Kovarianz von 0 (Null) abweicht, umso deutlicher ist die lineare Beziehung zwischen den beiden Variablen.

# mit ggplot() und gridExtra
library(gridExtra)

p1 = ggplot(epg) + 
  aes(y = F2, x = COG) + 
  geom_point()
p2 = ggplot(epg) + 
  aes(y = F1, x = COG) + 
  geom_point()

temp = with(epg, V %in% c("i", "I", "E", "a"))
p3 = ggplot(epg[temp, ]) + 
  aes(y = F1, x = SUM1278) + 
  geom_point()

# um die Bilder in 3 Spalten und einer Reihe zu arrangieren
grid.arrange(p1, p2, p3,  ncol = 3, nrow = 1)

Links: Kovarianz hoch und positiv (509,6908, s.u.); Mitte: Kovarianz nah an 0 (-24,26598); Rechts: Kovarianz mittelstark und negativ (-289,516)

Die Kovarianz berechnet sich als die Produkt-Summe der Abweichungen vom Mittelwert

y = epg$F2
x = epg$COG
n = length(y)

# Mittelwert:
mx = mean(x)
my = mean(y)

# Abweichungen vom Mittelwert
dx = x - mean(x)
dy = y - mean(y)

# Kovarianz = Produkt-Summe der Abweichungen dividiert durch n-1
covxy = sum(dx*dy)/(n-1)
covxy

## [1] 509.6908

# all das vorgenannte als Funktion in R
cov(x, y)

## [1] 509.6908

Einige Merkmale der Kovarianz

cov(x,y) == cov(y,x)
cov(x,x) == var(x) ( == sd(x)^2 )
var(x+y) == var(x)+var(y) + 2 * cov(x,y)

Daher: wenn es keine lineare Beziehung zwischen x und y gibt, ist cov(x,y) 0 (Null), sodass

var(x+y) == var(x) + var(y)

2. Kovarianz und Korrelation

Die Korrelation (Pearson’s product-moment correlation), r, ist dasselbe wie die Kovarianz, aber sie normalisiert für die Größe von x und y

r ist die Kovarianz von x, y, dividiert durch deren Standardabweichungen
r variiert zwischen -1 und +1

cov(x,y)

## [1] 509.6908

r = cov(x,y)/(sd(x) * sd(y))
r

## [1] 0.8917474

cor(x,y)

## [1] 0.8917474

cov(x,y)

## [1] 509.6908

xgross = x*1000
cov(xgross,y)

## [1] 509690.8

cor(xgross,y)

## [1] 0.8917474

3. Regression

y-auf-x Regression: y (die abhängige Variable) soll durch x (die unabhängige Variable) modelliert werden, also y durch die Werte von x eingeschätzt werden.

Die Regressionslinie ist jene gerade Linie durch die Verteilung, sodass der Abstand der Stichprobe zu der Linie minimiert wird.

Diese Regressionslinie durchschneidet immer (mx, my), also den Mittelwert (X) der Verteilung

Die Regressionslinie wird berechnet über

\[ ŷ = bx + k \]

\(b\) ist die Steigung (Engl.: ‘slope’):

b = r * sd(y)/sd(x)
b

## [1] 670.267

# oder
b = cov(x,y)/var(x)
b

## [1] 670.267

\(k\) ist das sog. Intercept (y-Achsenabschnitt):

k = my - b*mx
k

## [1] 610.6845

\(ŷ\) sind die eingeschätzten Werte, die auf der Regressionslinie liegen:

yhut = b * x + k
yhut

##  [1] 1318.1885 1537.2299 1472.4562 1137.3228 1933.5798 1884.1917 1911.0023
##  [8] 1792.4709 2329.5949 2453.9186 2428.5297 2183.2338 2183.2338 1951.2184
## [15] 1884.1917 1137.3228 1094.7661 1137.3228 1029.6013 1169.2402 1169.2402
## [22] 1094.7661 1328.8276 1113.3847 1392.6626 1213.9247 2051.7584 2109.9658
## [29] 1988.4554 1988.4554 2459.3240 2733.1965 2751.8150 2495.8103 2727.3169
## [36] 2315.4939 2484.9494 1113.3847  945.8179  945.8179  945.8179 1169.2402
## [43] 1137.3228 1169.2402 1169.2402

plot(y ~ x)
abline(k, b)
points(x, yhut, col = "red")

Error und SSE

Der error (auch residuals) ist der Unterschied zwischen den tatsächlichen und den eingeschätzten Werten.

error = y - yhut
error

##  [1]  321.81154  282.77005  247.54375  602.67723  346.42025  395.80834
##  [7]  448.99766  299.54908  130.40514   26.08143  -88.52967   96.76616
## [13]  376.76616  448.78164  215.80834 -136.81277 -214.76614 -149.47677
## [19] -109.60130 -209.24025 -409.24025 -134.76614  211.17238  146.61533
## [25]  -12.66256   46.07529  185.00160 -127.80579  190.99459  595.86459
## [31]  -99.32395 -438.03647 -713.17499 -195.81026 -268.88693 -177.93387
## [37] -304.94945 -158.40367 -205.81793 -145.81793 -145.81793 -149.24025
## [43] -277.32277 -353.23225 -389.24025

Der SSE (‘sum of the squares of the error’ ist die sogenannte Summe der Quadrate, d.h. die Summe der quadrierten Fehler:

SSE = sum(error^2)
SSE

## [1] 3871019

Die Regressionslinie wird auf eine solche Weise berechnet, dass SSE minimiert wird.

Der SSE wird übrigens manchmal auch RSS (‘residual sum of squares’) genannt.

Regression mit `lm()`

reg = lm(y ~ x)
plot(y ~ x)
abline(reg)

#Regressionskoeffizienten (Intercept und Steigung)
coef(reg)

## (Intercept)           x 
##    610.6845    670.2670

# Eingeschätzte Werte
yhut = b*x + k
yhut

##  [1] 1318.1885 1537.2299 1472.4562 1137.3228 1933.5798 1884.1917 1911.0023
##  [8] 1792.4709 2329.5949 2453.9186 2428.5297 2183.2338 2183.2338 1951.2184
## [15] 1884.1917 1137.3228 1094.7661 1137.3228 1029.6013 1169.2402 1169.2402
## [22] 1094.7661 1328.8276 1113.3847 1392.6626 1213.9247 2051.7584 2109.9658
## [29] 1988.4554 1988.4554 2459.3240 2733.1965 2751.8150 2495.8103 2727.3169
## [36] 2315.4939 2484.9494 1113.3847  945.8179  945.8179  945.8179 1169.2402
## [43] 1137.3228 1169.2402 1169.2402

yhut = predict(reg)
yhut

##         1         2         3         4         5         6         7 
## 1318.1885 1537.2299 1472.4562 1137.3228 1933.5798 1884.1917 1911.0023 
##         8         9        10        11        12        13        14 
## 1792.4709 2329.5949 2453.9186 2428.5297 2183.2338 2183.2338 1951.2184 
##        15        16        17        18        19        20        21 
## 1884.1917 1137.3228 1094.7661 1137.3228 1029.6013 1169.2402 1169.2402 
##        22        23        24        25        26        27        28 
## 1094.7661 1328.8276 1113.3847 1392.6626 1213.9247 2051.7584 2109.9658 
##        29        30        31        32        33        34        35 
## 1988.4554 1988.4554 2459.3240 2733.1965 2751.8150 2495.8103 2727.3169 
##        36        37        38        39        40        41        42 
## 2315.4939 2484.9494 1113.3847  945.8179  945.8179  945.8179 1169.2402 
##        43        44        45 
## 1137.3228 1169.2402 1169.2402

# Error
error = y - yhut
error

##          1          2          3          4          5          6 
##  321.81154  282.77005  247.54375  602.67723  346.42025  395.80834 
##          7          8          9         10         11         12 
##  448.99766  299.54908  130.40514   26.08143  -88.52967   96.76616 
##         13         14         15         16         17         18 
##  376.76616  448.78164  215.80834 -136.81277 -214.76614 -149.47677 
##         19         20         21         22         23         24 
## -109.60130 -209.24025 -409.24025 -134.76614  211.17238  146.61533 
##         25         26         27         28         29         30 
##  -12.66256   46.07529  185.00160 -127.80579  190.99459  595.86459 
##         31         32         33         34         35         36 
##  -99.32395 -438.03647 -713.17499 -195.81026 -268.88693 -177.93387 
##         37         38         39         40         41         42 
## -304.94945 -158.40367 -205.81793 -145.81793 -145.81793 -149.24025 
##         43         44         45 
## -277.32277 -353.23225 -389.24025

resid(reg)

##          1          2          3          4          5          6 
##  321.81154  282.77005  247.54375  602.67723  346.42025  395.80834 
##          7          8          9         10         11         12 
##  448.99766  299.54908  130.40514   26.08143  -88.52967   96.76616 
##         13         14         15         16         17         18 
##  376.76616  448.78164  215.80834 -136.81277 -214.76614 -149.47677 
##         19         20         21         22         23         24 
## -109.60130 -209.24025 -409.24025 -134.76614  211.17238  146.61533 
##         25         26         27         28         29         30 
##  -12.66256   46.07529  185.00160 -127.80579  190.99459  595.86459 
##         31         32         33         34         35         36 
##  -99.32395 -438.03647 -713.17499 -195.81026 -268.88693 -177.93387 
##         37         38         39         40         41         42 
## -304.94945 -158.40367 -205.81793 -145.81793 -145.81793 -149.24025 
##         43         44         45 
## -277.32277 -353.23225 -389.24025

# SSE:
sum(error^2)

## [1] 3871019

deviance(reg)

## [1] 3871019

Regression: drei wichtige Quantitäten

1. SSE (oder RSS) sum of the squared errors

SSE = sum(error^2)
SSE

## [1] 3871019

# oder
SSE = deviance(reg)
SSE

## [1] 3871019

1. SSY (oder SST): sum-of-the-squared deviations der tatsächlichen Werte

SSY = sum( (y - my)^2)
SSY

## [1] 18902691

1. SSR: sum of the squared-deviations von ŷ, also von den eingeschätzten Werten

SSR = sum((yhut - my)^2)
SSR

## [1] 15031672

\[SSY=SSR+SSE\]

SSY = SSR + SSE
SSY

## [1] 18902691

R-squared (\(R^2\))

\[SSY=SSR+SSE\] Je besser die Werte durch die Regressionslinie modelliert werden (also je geringer der Abstand zwischen y und ŷ) umso kleiner SSE, sodass im besten Fall SSE = 0 und SSY = SSR oder SSR/SSY = 1 (bedeutet: die tatsächlichen Werte sitzen auf der Linie).

\(R^2 = SSR/SSY\) beschreibt auch die Proportion der Varianz in \(y\) die durch die Regressionlinie erklärt werden kann.

\(R^2\) variiert zwischen 0 (keine ‘Erklärung’) und 1 (die Regressionslinie erklärt 100% der Varianz in y).

Diese Quantität \(SSR/SSY\) nennt man auch R-squared, weil sie denselben Wert hat wie der quadrierte Korrelationskoeffizient \(r\):

SSR/SSY

## [1] 0.7952134

cor(x, y)^2

## [1] 0.7952134

Da \(r\) zwischen -1 und +1 variiert, muss \(R^2\) zwischen 0 und 1 variieren.

Signifikanz-Test

Was ist die Wahrscheinlichkeit, dass ein lineares Verhältnis zwischen x und y besteht?

Dies kann mit einem t-test mit n-2 Freiheitsgraden berechnet werden: \[tstat = \frac{r}{\sqrt{\frac{1-r^2}{n-2}}}\]

rsb = sqrt( (1 - r^2)/(n-2))
tstat = r/rsb
tstat

## [1] 12.92187

2 * (1 - pt(tstat, n-2))

## [1] 2.220446e-16

fstat = tstat^2
1 - pf(fstat, 1, n-2)

## [1] 2.220446e-16

cor.test(x,y)

## 
##  Pearson's product-moment correlation
## 
## data:  x and y
## t = 12.922, df = 43, p-value < 2.2e-16
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  0.8103217 0.9393887
## sample estimates:
##       cor 
## 0.8917474

Die Wahrscheinlichkeit, dass die Variablen nicht miteinander linear assoziiert sind, ist fast 0. (Hoch signifikant, \(p < 0.001\)).

summary(reg)

## 
## Call:
## lm(formula = y ~ x)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -713.17 -195.81  -99.32  215.81  602.68 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)   610.68      94.65   6.452 8.03e-08 ***
## x             670.27      51.87  12.922  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 300 on 43 degrees of freedom
## Multiple R-squared:  0.7952, Adjusted R-squared:  0.7905 
## F-statistic:   167 on 1 and 43 DF,  p-value: < 2.2e-16

Es gibt eine signifikante lineare Beziehung zwischen COG und F2 (\(R^2 = 0.80\), \(F[1, 43] = 167\), \(p < 0.001\)).

Gültigkeit einer Regression¹

Die Residuen (‘residuals’) sollen:

(a) von einer Normalverteilung nicht signifikant abweichen:

shapiro.test(resid(reg))

## 
##  Shapiro-Wilk normality test
## 
## data:  resid(reg)
## W = 0.97037, p-value = 0.2987

(b) eine konstante Varianz aufweisen. Insbesondere sollten die Residuals nicht wesentlich größer am Anfang/Ende sein, sondern auf eine randomisierte Weise um die 0 Linie verteilt sein. Das ist hier nicht der Fall:

plot(resid(reg))
abline(h=0, lty=2)

(c) keine Autokorrelation aufweisen. Autokorrelation ist, wenn die Werte eines Signals mit sich selbst korreliert sind. Ein gutes Beispiel für autokorrelierte Daten wäre jene in einem stimmhaften, also periodischem, Sprachsignal.

acf(resid(reg))

Wenn die meisten ACF-Werte innerhalb der blauen Linien liegen, gibt es keine Autokorrelation. Insbesondere die Werte bei lag 1 und 2 beobachten: diese sollten innerhalb des Vertauensintervalls liegen (hier nicht der Fall).

Regression

Jonathan Harrington & Ulrich Reubold

Juni 2018

Kovarianz, Korrelation, (lineare) Regression

1. Kovarianz

Einige Merkmale der Kovarianz

2. Kovarianz und Korrelation

3. Regression

Error und SSE

Regression mit `lm()`

Regression: drei wichtige Quantitäten

R-squared (\(R^2\))

Signifikanz-Test

Gültigkeit einer Regression¹

Regression

Jonathan Harrington & Ulrich Reubold

Juni 2018

Kovarianz, Korrelation, (lineare) Regression

1. Kovarianz

Einige Merkmale der Kovarianz

2. Kovarianz und Korrelation

3. Regression

Error und SSE

Regression mit lm()

Regression: drei wichtige Quantitäten

R-squared (\(R^2\))

Signifikanz-Test

Gültigkeit einer Regression1

Regression mit `lm()`

Gültigkeit einer Regression¹