Matrixalgebra

Bisher haben wir gelernt, dass viele statistische Größen mithilfe von Summen, Mittelwerten, Quadraten, Abweichungen und weiteren recht einfachen Rechenoperationen bestimmt werden können. Diese können oft noch einfacher mithilfe von Matrizen und Vektoren dargestellt werden. Vielleicht sind Ihnen Matrizen aus der Schule bekannt. Falls nicht, ist das kein Problem, denn dieser Beitrag wird Vektoren und Matrizen ausführlich erklären. Matrixalgebra wird bspw. auch in Eid et al. (2017) im Anhang B: Matrixalgebra ab Seite 1051 behandelt.

Vektoren kennen wir bereits aus den ersten Datensätzen, die wir kennengelernt haben. Bspw. enthält ein Variablenvektor einer Person einfach nur die Einträge der Variablen aus dem Datensatz dieser spezifischen Person. Matrizen haben wir im Grunde auch schon kennengelernt. Ein Datensatz ist eng mit einer Matrix verwandt. Genauso wie ein Datensatz, besteht auch eine Matrix aus Zeilen und Spalten. Der Hauptunterschied ist, dass bei einer Matrix nur numerische Inhalte, also Zahlen, erlaubt sind. Wenn wir Daten in R verarbeiten wollen, wird der Datensatz oft in Matrizen umtransformiert (falls er vorher nicht-numerische Inhalte enthielt) und dann mit geeigneten Operationen, sogenannten Matrixoperationen, verarbeitet.

R ist eine vektorbasierte Programmiersprache, was bedeutet, dass möglichst viel mit Vektor- oder Matrixoperationen durchgeführt werden soll, da diese besonders optimiert (und damit besonders schnell) sind. Um davon Gebrauch zu machen, müssen wir uns mit diesen Operationen vertraut machen:

Vektoren

Vektoren werden häufig (aber nicht immer, Ausnahmen bestätigen die Regel) in Kleinbuchstaben dargestellt. Seien x und y zwei Vektoren, die wir mit dem Zuordnungspfeil <- und mit der Vektorfunktion c() erstellen:

x <- c(4, 2, 3)
y <- c(10, 8, 6)

Vektoren werden meistens als sogenannte Spaltenvektoren (dazu später mehr) dargestellt. x und y sehen also so aus

Die Elemente werden mit der jeweiligen Position im Vektor nummeriert. Das Element $x_1$ ist also das 1. Element von $x$, nämlich 4. Dies kann auch so dargestellt werden:

y[2]

## [1] 8

y[2:3]

## [1] 8 6

y[c(1,3)]

## [1] 10  6

Addition und Subtraktion bei Vektoren

Die Addition von Vektoren funktioniert elementenweise. Das bedeutet, dass das 1. Element des 1. Vektors und das 1. Element des 2. Vektors miteinander addiert werden und das 2. Element des 1. Vektors mit dem 2. Element des 2. Vektors miteinander addiert werden, etc.

Elementeweise Additionen funktionieren super simpel, indem wir x und y einfach mit + verknüpfen.

x + y  # Addition

## [1] 14 10  9

Wenn x und y nicht dieselbe Länge haben, ist es in R oft so, dass die Vektoren künstlich verlängert werden, um verrechnet zu werden. Dies sollten wir immer im Hinterkopf behalten.

z <- c(1:6) # Zahlen 1 bis 6
z + y

## [1] 11 10  9 14 13 12

z ist hier doppelt so lang wie y, sodass in der Addition y einfach zweimal hintereinander geschrieben wird, damit die Addition möglich ist, denn eine Addition bei Vektoren (und auch Matrizen) funktioniert nur, wenn die beiden Elemente das gleiche Format haben! Im Kontext der Matrixalgebra ist die künstliche Verlängerung eines Vektors zur Addition kein zulässiger Rechenschritt, sondern eine Besonderheit von R, was hin und wieder zu Problemen oder Fehlern führen kann. Die Länge eines Vektors entspricht seiner Anzahl an Elementen. Diese können wir erhalten über

length(x)

## [1] 3

Analog zur Addition erfolgt auch die Subtraktion von Vektoren elementenweise. Wollen wir $y - x$ berechnen, so würden wir entsprechend zur Addition das 1. Element von $y$ verrechnen mit dem 1. Element von $x$ usw. und erhalten folgendes Resultat:

y-x

## [1] 6 6 3

Multiplikation und Division bei Vektoren

Bezüglich der Multiplikation und Division bei Vektoren wird unterschieden zwischen zwei Rechenoperationen: (1) Die Multiplikation bzw. Division eines Vektors mit einer Zahl, also einem Skalar, und (2) die Multiplikation zweier Vektoren. Zunächst befassen wir uns mit der ersten Variante. Wenn wir einen Vektor mit einer Zahl bzw. einem Skalar multiplizieren, so bewirkt dies eine elementenweise Multiplikation mit dieser Zahl. Wenn wir beispielsweise $k$ mit $x$ multiplizieren wollen, so erhalten wir:

Für $k=3$ erhalten wir bspw.

In R sieht das so aus

3*x

## [1] 12  6  9

Genauso können wir auch jedes Element durch 2 teilen, indem wir den Vektor $x$ mit dem Skalar $\frac{1}{2}$ also 0.5 multiplizieren.

1/2*x

## [1] 2.0 1.0 1.5

Im Kontrast ist die zweite Rechenoperation, nämlich die Multiplikation zweier Vektoren, komplexer. An dieser Stelle werden wir das Mysterium noch nicht lüften, weil es sich anbietet, diese besondere Form der Multiplikation später bei den Matrizen im Detail zu behandeln. Wir möchten jedoch darauf verweisen, dass es in R möglich ist, zwei Vektoren der gleichen Länge miteinander elementenweise zu multiplizieren über folgenden Ausdruck:

x*y 

## [1] 40 16 18

Hierbei handelt es sich jedoch nicht um eine matrixalgebraische Multiplikation, sondern eine R-spezifische Art der Multiplikation.

Matrizen

Ein Vektor ist eine eindimensionale Sammlung von Zahlen. Die Elemente werden einfach durchnummeriert. Eine Matrix ist ein zweidimensionales Objekt, welche aus einer Vielzahl von Vektoren gleicher Länge besteht, die aneinander “geklebt” werden. Matrizen werden oft in Großbuchstaben beschrieben. Elemente von Matrizen hingegen in Kleinbuchstaben. Auch hier ist das nicht wirklich einheitlich geregelt.

Würden wir die beiden Vektoren $x$ und $y$ aneinander “kleben”, erhalten wir eine Matrix $A$, welche 3 Zeilen und 2 Spalten enthält.

Analog zu Vektoren ist die Position in einer Matrix auch über Indizes definiert, wobei nun zwei Laufnummern vergeben werden für jedes Element $a_{ij}$. Der erste Index $i = 1,2,…,n$ bezieht sich auf die Zeile, in welchem sich das Element befindet und der zweite Index $j=1,…,m$ die Spalte. Demnach beschreibt $a_{32}$ das Matrixelement in der 3. Zeile (1. Index) und 2. Spalte (2. Index). Für die Matrix $A$ entspricht dies dem Wert 6. Das allgemeine Format einer Matrix kann über ihre Typangaben beschrieben werden: Auch hier wird zuerst die Anzahl der Zeilen der Matrix angegeben und dann ihre Spalten, sprich in der Form $z \times s$. Bei Matrix $A$ handelt es sich also um eine Matrix vom Typ $3 \times 2$ (3 Zeilen und 2 Spalten). Im Folgenden sind ein paar weitere Beispiele aufgelistet:

Typ $4\times2$: $$\begin{pmatrix} 2 & 4 \\ 5 & 1 \\ -2 & 3\\ 1 & 2 \end{pmatrix}$$

Typ $2\times3$: $$\begin{pmatrix} 5 & 7 & 9 \\ 1 & 9 & 4\end{pmatrix}$$

Typ $3\times1$ (sog. Spaltenvektor): $$\begin{pmatrix} 4 \\ 2 \\ 3 \end{pmatrix}$$

Typ $1\times3$ (sog. Zeilenvektor): $$\begin{pmatrix} 2 & 5 & 7 \end{pmatrix}$$

Die gerade behandelten Vektoren können wir ganz leicht zu einer Matrix machen, indem wir den Befehl as.matrix bspw. auf x anwenden. Dieser Befehl erzeugt eine $3\times1$ Matrix - also aus mathematischer Sicht einen Spaltenvektor.

as.matrix(x)

##      [,1]
## [1,]    4
## [2,]    2
## [3,]    3

Wir können die beiden Vektoren auch zu einer Matrix kombinieren, indem wir sie mit dem Befehl cbind (was für column binding steht) zusammenfügen - genauso geht dies auch mit rbind (was für row binding steht):

A <- cbind(x, y)
A

##      x  y
## [1,] 4 10
## [2,] 2  8
## [3,] 3  6

B <- rbind(x, y)
B

##   [,1] [,2] [,3]
## x    4    2    3
## y   10    8    6

Man erkennt, dass mit cbind die beiden Vektoren als Spaltenvektoren nebeneinander verknüpft werden, während bei rbind die Vektoren als Zeilenvektoren untereinander verknüpft wurden.

B[1, 2] 

## x 
## 2

B[1, ] # 1. Zeile

## [1] 4 2 3

B[, 2] # 2. Spalte

## x y 
## 2 8

Transposition, Addition und Subtraktion

So wie A und B erzeugt wurden, ist ersichtlich, dass die Spalten von A den Zeilen von B entspricht. Wir können Zeilen und Spalten auch vertauschen, indem wir die Matrix mit dem Befehl t() transponieren:

A

##      x  y
## [1,] 4 10
## [2,] 2  8
## [3,] 3  6

t(A)

##   [,1] [,2] [,3]
## x    4    2    3
## y   10    8    6

B

##   [,1] [,2] [,3]
## x    4    2    3
## y   10    8    6

Wir erkennen, dass die Matrix B gerade die Transponierte von A ist! Im Prinzip wurden die Indizes von jedem Element der Matrix getauscht, sprich durch die Transposition von $A$ wurde jedes Element $a_{ij}$ zu $a_{ji}$.

Die Matrixaddition und Matrixsubtraktion funktioniert genauso wie die von Vektoren. Sie wird elementenweise durchgeführt. Allerdings müssen dafür die Matrizen vom selben Typ sein, also gleich viele Zeilen und Spalten besitzen.

Die beiden Matrizen A und B lassen sich beispielsweise nicht addieren, da sie nicht das richtige Format haben:

A + B

## Error in A + B : non-conformable arrays

Eine $3\times 2$ Matrix lässt sich nicht mit einer $2\times 3$ Matrix addieren. Wird eine Matrix transponiert, so ändert sich damit auch ihr Typ. Die Zeilen werden zu Spalten und umgekehrt. Die Transponierte von $A$ ist also eine $2\times 3$ Matrix und lässt sich damit mit $B$ addieren:

t(A) + B

##   [,1] [,2] [,3]
## x    8    4    6
## y   20   16   12

Da die beiden Matrizen t(A) und B identisch sind, würden wir zum selben Ergebnis kommen, wenn wir jeden Eintrag von $B$ mit dem Skalar 2 multiplizieren.

B * 2 # skalare Multiplikation

##   [,1] [,2] [,3]
## x    8    4    6
## y   20   16   12

Multiplikation zweier Matrizen

Wie bereits im Abschnitt zur [Multiplikation und Division bei Vektoren] angedeutet, funktioniert die Multiplikation von Vektoren, und als Erweiterung die Multiplikation von zwei Matrizen, nicht elementenweise. Stellen wir uns vor, wir möchten die beiden Matrizen $A$ und $B$ miteinander multiplizieren, also $AB$ rechnen. Man kann zwei Matrizen nur dann miteinander multiplizieren, wenn ihre inneren Typangaben miteinander korrespondieren, d.h. wenn die ,,erste’’ Matrix (also $A$) genauso viele Spalten wie die ,,zweite’’ Matrix (also $B$) Zeilen hat: $A B = (p \times q) (q \times r) = (3 \times 2) (2 \times 3)$. Da die Matrix $A$ zwei Spalten und die Matrix $B$ zwei Zeilen hat, können die beiden Matrizen miteinander multipliziert werden.

Konkret werden bei der Matrixmultiplikation die Zeilen der ersten Matrix $A$ mit den Spalten der zweiten Matrix $B$ elementenweise multipliziert und diese Elemente werden anschließend addiert. Die resultierende Matrix hat dann so viele Zeilen wie die erste Matrix und so viele Spalten wie die zweite. Das resultierende Format der neuen Matrix ist also $p\times r$.

Für unseren Fall entsteht eine neue Matrix der Dimension $3\times 3$.

Drehen wir das Ganze um, erhalten wir

Hier entsteht eine Matrix vom Format $2\times 2$. Wir sehen deutlich, dass Matrixmultiplikation im Allgemeinen nicht kommutativ ist, also $AB \neq BA$.

Der Operator in R hierfür heißt %*% (verwenden wir stattdessen *, so wird eine elementenweise Multiplikation durchgeführt, was etwas komplett anderes ist!):

A %*% B # Matrixprodukt A*B

##      [,1] [,2] [,3]
## [1,]  116   88   72
## [2,]   88   68   54
## [3,]   72   54   45

B %*% A # Matrixprodukt B*A

##    x   y
## x 29  74
## y 74 200

An den Ergebnissen erkennen wir erneut, dass Matrixprodukte nicht kommutativ sind, also die Reihenfolge, in der multipliziert wird, wichtig ist. Analog zur Multiplikation zweier Matrizen kann man zwei Vektoren miteinander multiplizieren. Einen Exkurs hierfür finden Sie im Appendix.

Spezielle Matrizen

Eine quadratische Matrix ist eine Matrix mit gleich vielen Zeilen wie Spalten. Eine wichtige quadratische Matrix ist die Einheitsmatrix $I$, welche nur 1en auf der Diagonalen und sonst 0en hat. Diese ist gerade das Element, mit welchem wir getrost multiplizieren können (falls die Dimensionen stimmen), weil dann nichts passiert (wie die Multiplikation mit 1 bei Zahlen). Wir erhalten sie mit diag, was eigentlich eine (quadratische) Diagonalmatrix erzeugt mit beliebigen Elementen auf der Diagonalen:

diag(3) # Einheitsmatrix 3x3

##      [,1] [,2] [,3]
## [1,]    1    0    0
## [2,]    0    1    0
## [3,]    0    0    1

diag(1:3) # Diagonalmatrix mit Elementen 1,2,3 auf der Diagonalen

##      [,1] [,2] [,3]
## [1,]    1    0    0
## [2,]    0    2    0
## [3,]    0    0    3

Bisher haben wir Matrizen dadurch konstruiert, dass wir zwei Vektoren mit cbind() oder rbind() miteinander verknüpft haben. Wir können jedoch eine Matrix mit dem matrix() direkt definieren. Diesem übergeben wir einen Vektor und die Dimensionen der Matrix. Dem Argument data werden die Daten, die wir in die Matrix schreiben wollen, übergeben, nrow und ncol bestimmen die Anzahl der Zeilen und Spalten und mit byrow = T legen wir fest, dass wir die Matrix zeilenweise gefüllt bekommen möchten:

C <- matrix(data = c(1:9), # Daten/Inhalt des Vektors (hier die Zahlen 1 bis 9)
            nrow = 3,      # Zeilenanzahl der Matrix
            ncol = 3,      # Spaltenanzahl der Matrix
            byrow = T)     # Zeilen der resultierenden Matrix werden zuerst mit data gefüllt
C

##      [,1] [,2] [,3]
## [1,]    1    2    3
## [2,]    4    5    6
## [3,]    7    8    9

Wir können mit diag auch wieder die Diagonalelemente einer Matrix erfahren:

diag(C)

## [1] 1 5 9

Determinanten und Invertierung

Die Inverse, also jene Matrix, mit der wir eine Matrix $X$ multiplizieren müssen, um die Einheitsmatrix $I$ zu erhalten ($XX^{-1} = I$), lässt sich in R mit dem solve Befehl erhalten (dies geht nur bei quadratischen Matrizen):

solve(C)

## Error in solve.default(C) : 
##   system is computationally singular: reciprocal condition number = 2.59052e-18

Obwohl $C$ eine quadratische Matrix ist, lässt sie sich nicht invertieren, da sie singulär und damit nicht invertierbar ist. Dies bedeutet, dass es lineare Abhängigkeiten der Zeilen bzw. Spalten gibt. Wir können dies explizit prüfen, indem wir die Determinante bestimmen mit det:

det(C)

## [1] 6.661338e-16

round(det(C), 14)

## [1] 0

Mit round runden wir das Ergebnis auf die 14. Nachkommastelle. Eine Matrix ist genau dann invertierbar (also regulär im Vergleich zu singulär), wenn die Determinante dieser (quadratischen) Matrix nicht Null ist. Lineare Abhängigkeit bedeutet, dass die Zeilen oder Spalten durch Addition, Subtraktion und skalare Multiplikationen auseinander hervorgehen. Die lineare Abhängigkeit zwischen den Spalten wird ersichtlich, wenn wir von der 2. Spalte die 1. Spalte abziehen und das Ergebnis zur 3. Spalte addieren - also de facto $2*2.Spalte - 1. Spalte$ rechnen:

2*C[, 2] - C[, 1]     # 2*2.Spalte - 1. Spalte rechnen ist gleich

## [1] 3 6 9

C[, 3]               # 3. Spalte

## [1] 3 6 9

Hätten wir C^-1 gerechnet, so hätten wir eine elementeweise Invertierung durchgeführt:

C^-1

##           [,1]  [,2]      [,3]
## [1,] 1.0000000 0.500 0.3333333
## [2,] 0.2500000 0.200 0.1666667
## [3,] 0.1428571 0.125 0.1111111

C^-1 %*% C # ist nicht die Einheitsmatrix

##          [,1]     [,2]     [,3]
## [1,] 5.333333 7.166667 9.000000
## [2,] 2.216667 2.833333 3.450000
## [3,] 1.420635 1.799603 2.178571

C^-1 * C   # elementeweise ergibt überall 1 - ist immer noch nicht die Einheitsmatrix!

##      [,1] [,2] [,3]
## [1,]    1    1    1
## [2,]    1    1    1
## [3,]    1    1    1

Dies bedeutet, dass C^-1 in R nicht die Invertierung betitelt, sondern solve!

Betrachten wir nun eine invertierbare Matrix D:

D <- matrix(c(1, 0, 0,
              1, 1, 1,
              2, 4, 5), 3, 3, byrow = T)
det(D)

## [1] 1

Die Determinante von D ist 1. Somit können wir D invertieren. Das Produkt aus D mit seiner Inversen ergibt gerade die 3x3 Einheitsmatrix:

solve(D)

##      [,1] [,2] [,3]
## [1,]    1    0    0
## [2,]   -3    5   -1
## [3,]    2   -4    1

D %*% solve(D)

##      [,1] [,2] [,3]
## [1,]    1    0    0
## [2,]    0    1    0
## [3,]    0    0    1

solve(D) %*% D

##      [,1] [,2] [,3]
## [1,]    1    0    0
## [2,]    0    1    0
## [3,]    0    0    1

Das Produkt von $D$ und $D^{-1}$ ist (ausnahmsweise) kommutativ: $DD^{-1}=D^{-1}D=I$. Im Allgemeinen ist die Bestimmung der Inversen einer Matrix komplex. Für eine $2 \times 2$-Matrix gibt es jedoch eine einfache Lösung. Eine $2\times 2$ Matrix $$M=\begin{pmatrix}a&b\\ c&d \end{pmatrix}$$

lässt sich genau dann invertieren, wenn die Determinante nicht 0 ist, also wenn keine lineare Abhängigkeit innerhalb der Zeilen oder der Spalten besteht. Die Determinante von $M$ bestimmen wir als $$\text{det}[M] = ad-bc.$$ Ist diese nicht 0, so erhalten wir die Inverse von $M$ durch $$M^{-1}=\frac{1}{\text{det}[M]}\begin{pmatrix}d&-b\\ -c&a\end{pmatrix}=\frac{1}{ad-bc}\begin{pmatrix}d&-b\\ -c&a\end{pmatrix}.$$ Wir müssen also die Diagonalelemente von $M$ vertauschen und die Nebendiagonalelemente mit einem Minus versehen. Anschließend müssen wir jeden Eintrag durch die Determinante teilen (oder wir nehmen einfach solve).

M <- matrix(c(2, 2, 3, 4), ncol = 2, nrow = 2, byrow = TRUE)
M

##      [,1] [,2]
## [1,]    2    2
## [2,]    3    4

det(M)

## [1] 2

M[1,1]*M[2,2] - M[1,2]*M[2,1] # Determinante mit Hand

## [1] 2

K <- matrix(c(M[2,2], -M[1,2], -M[2,1], M[1,1]), byrow = TRUE, ncol = 2, nrow = 2)
K # Kofaktorenmatrix

##      [,1] [,2]
## [1,]    4   -2
## [2,]   -3    2

1/(M[1,1]*M[2,2] - M[1,2]*M[2,1])*K # Inverse von M

##      [,1] [,2]
## [1,]  2.0   -1
## [2,] -1.5    1

solve(M) # Inverse von M

##      [,1] [,2]
## [1,]  2.0   -1
## [2,] -1.5    1

Wofür sind Matrixoperationen in der Statistik wichtig? Im letzten Abschnitt schauen wir uns für fb25 an, wie man einfache Statistiken mit Hilfe der Matrixalgebra bestimmen kann.

Statistiken mit Matrixoperationen bestimmen

# Daten laden
load(url('https://pandar.netlify.app/daten/fb25.rda'))

Summen lassen sich sehr leicht auch durch ein Matrixprodukt darstellen. Nehmen wir beispielsweise die Big Five Variablen des fb25 Datensatzes her und nennen diese X. Wenn wir einen 1en-Zeilenvektor der Dimension $1\times n$ bilden, wobei $n$ die Anzahl an Zeilen von X ist, dann können wir die Summe der Elemente pro Spalte von X durch ein Matrixprodukt ausdrücken:

X <- as.matrix(fb25[, c("extra", "vertr", "gewis", "neuro", "offen")])
X <- na.omit(X)  # fehlende Werte wurden entfernt
n <- nrow(X)
z <- t(rep(1, n)) # 1en-Zeilenvektor der Länge n
z %*% X

##      extra vertr gewis neuro offen
## [1,] 709.5   736 755.5 669.5   824

colSums(X)

## extra vertr gewis neuro offen 
## 709.5 736.0 755.5 669.5 824.0

In der Vorlesung haben Sie auch folgende Formel zur Berechnung der Kovarianzmatrix für eine beliebige Anzahl an Variablen kennengelernt

$$\Sigma = \frac{1}{n-1} (X-\overline{X})’(X-\overline{X}),$$ wobei $X$ eine Datenmatrix, $\overline{X}$ eine Matrix mit den Mittelwerten der Variablen und $n$ die Stichprobengröße ist. Für unsere Matrix $X$ müssten Sie, würden Sie die Formel nicht anwenden, individuelle Varianzen und Kovarianzen berechnen, was bei fünf Variablen insgesamt 15 Befehlen in R entsprechen würde. Glücklicherweise kann man in R mit cov die Kovarianzmatrix erhalten, eine manuelle Durchführung ist demnach nicht erforderlich:

cov(X) 

##             extra       vertr      gewis       neuro       offen
## extra  0.94810321  0.04483937 0.01549897 -0.39451470  0.11271360
## vertr  0.04483937  0.67222602 0.09283436 -0.01046936  0.08050809
## gewis  0.01549897  0.09283436 0.63037708  0.09877535  0.01823878
## neuro -0.39451470 -0.01046936 0.09877535  0.98957052 -0.01430850
## offen  0.11271360  0.08050809 0.01823878 -0.01430850  0.79321030

Der Vollständigkeit halber werden wir aber beispielhaft die Kovarianzmatrix mit Hilfe der obigen Formel manuell berechnen, um auch zu prüfen, ob cov tatsächlich richtige Lösungen ausgibt. Im ersten Schritt benötigen wir die Mittelwertmatrix $\overline{X}$, welche wir mit folgenden Formeln, die Sie auch in der Vorlesung kennengelernt haben, erhalten

$$\overline{x}’ = \frac{1}{n}u’X, \qquad \overline{X} = u\overline{x}’,$$ wobei $\overline{x}’$ ein Vektor der Mittelwerte und $u’$ ein Vektor bestehend aus den Werten 1 ist. Für unsere Matrix $X$ erhalten wir die Größen in R folgendermaßen:

u <- rep(1, nrow(X)) # u enthält so viele Einsen als Elemente, wie die Matrix X an Beobachtungen hat
x_mean <- 1/nrow(X) * t(u) %*% X
x_mean # Mittelwertevektor

##         extra    vertr    gewis    neuro   offen
## [1,] 3.378571 3.504762 3.597619 3.188095 3.92381

X_mean <- u %*% x_mean # Mittelwertematrix
head(X_mean) # um die ersten sechs Zeilen der Matrix anzuzeigen

##         extra    vertr    gewis    neuro   offen
## [1,] 3.378571 3.504762 3.597619 3.188095 3.92381
## [2,] 3.378571 3.504762 3.597619 3.188095 3.92381
## [3,] 3.378571 3.504762 3.597619 3.188095 3.92381
## [4,] 3.378571 3.504762 3.597619 3.188095 3.92381
## [5,] 3.378571 3.504762 3.597619 3.188095 3.92381
## [6,] 3.378571 3.504762 3.597619 3.188095 3.92381

Nun können wir die Differenz zwischen der Datenmatrix $X$ und seiner Mittelwertmatrix $\overline{X}$ berechnen ($X-\overline{X}$). Wir führen also eine Zentrierung durch, da wir lediglich den Mittelwert von jedem beobachteten Wert abziehen.

X_zentriert <- X - X_mean
head(X_zentriert)

##           extra      vertr       gewis      neuro      offen
## [1,]  0.1214286  0.4952381 -1.09761905 -1.1880952 -1.4238095
## [2,]  0.6214286 -0.5047619  0.90238095 -1.1880952  0.5761905
## [3,] -0.8785714  0.4952381 -0.09761905 -0.1880952  0.5761905
## [4,]  0.6214286 -0.5047619  0.90238095  0.8119048  1.0761905
## [5,]  1.1214286 -1.5047619 -0.59761905  0.8119048  0.5761905
## [6,] -0.3785714  0.9952381 -0.09761905  1.8119048  0.5761905

Zuletzt können wir nun das Kreuzprodukt $(X-\overline{X})’(X-\overline{X})$ berechnen und mit $1/(n-1)$ multiplizieren:

kreuzprodukt <- t(X_zentriert) %*% X_zentriert
X_cov <- 1/(nrow(X)-1) * kreuzprodukt
X_cov

##             extra       vertr      gewis       neuro       offen
## extra  0.94810321  0.04483937 0.01549897 -0.39451470  0.11271360
## vertr  0.04483937  0.67222602 0.09283436 -0.01046936  0.08050809
## gewis  0.01549897  0.09283436 0.63037708  0.09877535  0.01823878
## neuro -0.39451470 -0.01046936 0.09877535  0.98957052 -0.01430850
## offen  0.11271360  0.08050809 0.01823878 -0.01430850  0.79321030

Nicht überraschend kommen wir zu demselben Ergebnis wie cov. Mit Hilfe der Matrixalgebra besteht also die Möglichkeit, Statistiken wie die Kovarianzen und Varianzen über eine Formel direkt zu berechnen, statt sie individuell zu berechnen, was insbesondere bei einer hohen Anzahl an Variablen aufwendig wäre.

x <- matrix(data = c(3, -2, 1), nrow = 1) # Zeilenvektor x
y <- matrix(data = c(5, 1, 7), ncol = 1) # Spaltenvektor y

x %*% y

##      [,1]
## [1,]   20

y %*% x

##      [,1] [,2] [,3]
## [1,]   15  -10    5
## [2,]    3   -2    1
## [3,]   21  -14    7