---
title: "Faktorenanalyse Übung 3"
format: html
editor: visual
---

# Vorbereitungen {.unnumbered .unlisted visibility="hidden"}

Der folgende R-Chunk soll Ihnen bei der Installation von benötigten Paketen helfen und eine geeignete Dateistruktur anlegen, wenn sie nicht schon da ist. Speichern Sie die qmd-Datei in einem Ordner auf Ihrem Rechner ab (am besten in einem Ordner mit Cloudsynchronisation wie Dropbox, SWITCH oder oneDrive). Wenn Sie den folgenden Chunk ausführen, werden Ihnen in dem Ordner die Unterordner "files", "data", und "images" angelegt (wenn sie nicht schon da sind und im Fall unangetastet bleiben). Dann wird die \_common.R heruntergeladen, in der ein paar Optionen gespeichert sind. Wenn Sie daran Änderungen vornehmen und diese für später behalten möchten, müssen Sie die Zeile mit "download.file..." beginnt auskommentieren (# davor).

Versuchen Sie die folgenden Zeilen zu verstehen und führen Sie die Befehle aus. Wenn Sie etwas nicht verstehen, googlen Sie die Pakte (zB "r corrr") und schauen sich dort die Vignette des Paketes an.

Mit dem Die Installation von benötigten Paketen ist weiter vereinfacht.

```{r}
#| label: Installation
#| echo: true 
#| eval: false

# Prüfe, ob es in dem Ordner in der die Uebung_2.qmd gespeichert ist, folgende Unterordner gibt und wenn nicht, dann lege sie an
if(dir.exists("files")){} else {dir.create("files")}
if(dir.exists("data")){} else {dir.create("data")}
if(dir.exists("images")){} else {dir.create("images")}

# lade die _common.R und speichere sie im Ordner "files"
download.file("https://stat.ikmz.uzh.ch/Aufbau/Folien/Sitzung_06/files/_common.R", destfile= "_common.R")

# Wenn die folgenen Pakete jeweils nicht schon installiert sind, werden sie installiert

if (!require("devtools")) install.packages("devtools")

if(!require("corrr")) install.packages("corrr")
if(!require("kableExtra")) install.packages("kableExtra")
if(!require("psych")) install.package("psych")

# jetzt noch die Problembär-Pakete forcierter installieren. 
# Wenn die Nachricht kommt, das die Pakete seit der letzten Installation nicht geändert wurden und daher geskippt werden, is das perfekt in Ordnung

if(!require("sjmisc"))devtools::install_github("strengejacke/sjmisc")
```

# Übung 3

## 3.0 Laden Sie die Daten unserer Befragung {.unnumbered .unlisted .smaller}

(Für den Fall, dass es Probleme gibt, laden Sie die Daten [hier](https://stat.ikmz.uzh.ch/Aufbau/Folien/Sitzung_08/data/Stat_Aufbau_Befragung.RDS) herunter. Wenn das auch zu Problemen führt, melden Sie sich im Forum.)

```{r}
#| echo: true
#| code-summary: "R-Code anzeigen"
#| code-fold: false

# tidyverse laden 
suppressPackageStartupMessages(library("tidyverse"))

# Wenn wir die Dateien einmal geladen und aufbereitet haben (mit der Extradatei "Aufbereitung.qmd"), dann ist es besser, nicht immer wieder neu die Rohdaten zu laden, sondern die Aufbereiteten. 

# Prüfe, ob es in dem Ordner in der die Uebung_1_ab.qmd gespeichert ist, einen Unterordner "files" gibt und, wenn nicht, dann lege ihn an.
if(dir.exists("data")){} else {dir.create("data")}

download.file("https://stat.ikmz.uzh.ch/Aufbau/Folien/Sitzung_04/data/Stat_Aufbau_Befragung.RDS", "data/Stat_Aufbau_Befragung.RDS")

DATEN <- readRDS("data/Stat_Aufbau_Befragung.RDS")  |>
  haven::zap_formats()

# DATEN  |>  sjlabelled::get_label()

DATEN |> 
  saveRDS("data/Stat_Aufbau_Befragung.RDS")
```


## «psych» laden {.unnumbered .unlisted}

Das Paket psych brauchen wir im Folgenden mehrfach. Es kann für die Sitzung mit Hilfe von `library()`-funktion geladen werden. Danach braucht man eigentlich nicht mehr immer `psych::` vor jede Funktion schreiben, aber ich empfehle es und der Styleguide von R-Studio auch. Dann weiss man 1. immer, welches Paket gerade verwendet wird bei einer Funktion (für andere ist das auch deutlich leichter lesbar). 2. kann man sicher sein, dass es keine Konflikte gibt, weil zB mehrere Pakete die Funktion `filter` verwenden und 3. braucht man eben nicht mehr irgendwo am Anfang die ganzen Pakete mit `library` laden. Also, den Aufruf kann man machen, muss aber nicht.

```{r}
library(psych)
```

## Orientierung in den Daten {.unnumbered .unlisted}

Schauen Sie in den Fragebogen, welche Variablen zur RAQ-Skala gehören: [FB_STATA.pdf (uzh.ch)](https://stat.ikmz.uzh.ch/Aufbau/files/FB_STATA.pdf)

Dann rufen Sie immer vor einer höheren Analyse die Häufigkeitstabellen einmal auf und schauen sich die Variablen an. Was haben sie für einen Range? Gibt es Werte in den Variablen, die eigentlich nicht sein können (zB sind befragte teils über 2000 Jahre als)? Wie viele fehlende Werte (NA) haben die Variablen? Gibt es zB negative Variablen (-9) die als Wert mitverarbeitet werden, aber eigentlich für einen fehlenden Wert (Missing stehen)? Also:

```{r}
DATEN |> 
  sjmisc::frq(R101_01:R101_24)
```


## Aufgabenstellung: Führen Sie eine Faktorenanalyse für unsere RAQ aus mit Hilfe des Paketes `psych` {.unnumbered .unlisted .smaller}

Ziehen Sie dazu die Hilfe dieser Website heran: https://md.psych.bio.uni-goettingen.de/mv/unit/fa/efa.html

3.1 Schauen Sie sich die Korrelationsmatrix (mit `cor()`) an.

3.2 Führen Sie einen Bartlett's-Test aus

3.3 Führen Sie einen KMO-Test aus (`psych::KMO()`)

3.4. Machen Sie eine Parallel-Analyse (`psych::fa.parallel()`), um die Anzahl der Faktoren zu bestimmen.

3.5. Führen Sie eine Faktorenanalyse durch mit `psych::fa()` und rotieren Sie mit `rotate = "promax"`. 5.1. Schauen Sie sich die Faktorladungen an. 5.2. Wie viel Varianz erklären die Faktoren?

3.6. Erstellen Sie ein Diagramm für die Faktorenanalyse (`psych::fa.diagramm`).

3.7. Führen Sie eine Faktorenanalyse mit Rotation nach dem Kriterium `rotate = "oblimin"` durch.


## 3.1 Korrelationstabelle {.unnumbered .unlisted}

Erste Versuche, die nicht funktionieren:

```{r}
#| eval: false

# geht nicht:
DATEN |> 
  cor(R101_01:R101_24)

# geht auch nicht:
head(round(cor(DATEN$R101_01:DATEN$R101_24), 2))
```

So geht es:

Ich lege ein neues Datenset «items» an, das gebildet (\<-) wird indem\
1. die Daten genommen werden \|\> 2. die Variablen ausgewählt werden, die miteinander korreliert werden sollen \|\> 3. Fehlende Werte (NA) ausgelassen werden (omit) mit der Funktion na.omit()

Dann

```{r}
#|label: geht!

# erstmal die items extrahieren
items <- DATEN |> 
  select(R101_01:R101_24) |> 
  na.omit() # sorgt dafür, dass keine NAs mehr in den Daten sind

# jetzt funktioniert die Vorlage:
round(cor(items), 2)

# Schöner ist aber mit Pipes:

items |> 
  cor() |> 
  round(2)

# Und hier noch zur Freude aller eine grafische Darstellung der Korrelationen:
items |> 
  corrr::correlate() |> 
  corrr::rplot()

```

### 3.2 Bartlett-Test

Unterscheiden sich die Korrelationen aller Variablen mit allen anderen Variablen von 0? Etwas technischer: Hebt sich die tatsächliche Korrelationsmatrix (mit den 1 in der Diagonale und Werten daneben für die Kombis) von einer Matrix ab, die nur Nullen neben der Diagonalen hat?

```{r}
items |> 
  cor() |> 
  psych::cortest.bartlett(nrow(DATEN)) 
```

Volle Lotte!

### 3.3 KMO-Test

Der KMO fragt im Grunde, wie viel eine Variable mit irgendwelchen übrigen Variablen gemeinsam hat. Wenn der Wert nahe 1 ist, ist eine Variable für eine Faktorenanalyse mit den übrigen Variablen geeignet. Den MSA (Measure of Sampling Adequacy) gibt es für jedes item, aber auch für als Gesamtmass für alle Variablen. Unser Wert hier von .78 ist schon ganz gut.

```{r}
# aus der Vorlage:

psych::KMO(items)

# oder schöner mit Pipe statt Schachtel:

items |> 
  psych::KMO()

```

### 3.4 Parallelanalyse

Mit der Parallelanalyse kann festgestellt werden, wie viele Faktoren aus einem item-set extrahiert werden können bzw. sollten. Ein Orientierungspunkt ist der Eigenwert, also der Wert der angibt, ob ein Faktor mehr Varianz (\>1) oder weniger (\<1) als eine einzelne Variable auf sich vereint. Wenn der Eigenwert fast gleich 1 ist, kann man sagen, dass es schon noch sinnvoll ist den aufzunehmen, aber auch sagen, dass in dem Fall keine Dimensionsreduktion mehr stattfindet und der Faktor nicht mehr aufgenommen wird. Generell empfehle ich in Zweifelsfällen, die Analysen mal mit dem einen Faktor mehr und einmal mit dem Faktor weniger durchzuführen.

```{r}
items |> 
  psych::fa.parallel(# fa = "fa" # 
                     )
```

Wir gehen erstmal von 3 Hauptkomponten aus und probieren aber auch 4.

### 3.5 Die Faktorenanalyse

Direkt einfach nur die Faktorenanalyse.

```{r}
items |> 
  psych::fa(nfactors = 4, 
            rotate = "none") 
```

Erstmal die (unrotiert) Faktorenanalyse in m.ml.u speichern und dann ausgeben mit `print()`:

```{r}
m.ml.u <- items |> 
  psych::fa( 
    nfactors=4,
    n.obs = nrow(DATEN),
  SMC=TRUE,
  rotate="none",
  max.iter=100
  )

print(m.ml.u)
```

Hier mal nur die Faktorladungen (alle und nur für Ladungen über .4):

```{r}
print(m.ml.u$loadings)

print(m.ml.u$loadings, cutoff=.4)
```

### 3.6 Rotierte Faktorenlösung

Die initiale Faktorextraktion geht in der Regel schrittweise vor, indem der erste Faktor das Maximale an Varianz versucht zu erfassen. Der zweite Faktor kann noch den Rest (praktisch die Residuen der der ersten Faktorgeraden) abbilden. Der dritte Faktor den dann noch übrigen Rest und so weiter. Wenn das gemacht ist, kann man schauen, wie viele Faktoren noch ausreichend Varianz aufklären (Eigenwerte \>1). Dann haben wir aber das Problem, dass die aufgeklärte Varianz nicht gut auf die Faktoren verteilt ist. Aber wir können unsere extrahierten Faktoren (sagen wir 4) noch gut in dem 4-dimensionalen Raum drehen und anpassen, so dass jeder Faktor möglichst gleich viel Varianzproportionen abbildet (VARIMAX-Methode), während die Faktoren unkorreliert bleiben müssen (orthogonal). Etwas realistischer ist allerdings die Methode, bei der jeder Faktor maximale Varianzproportionen (promax) auf sich vereint und die Faktoren leicht miteinander korrelieren dürfen (oblimin):

```{r}

m.ml.r <- psych::fa(items, 
    nfactors=4,
    n.obs = nrow(items),
  SMC=TRUE,
  fm="ml",
    rotate="varmax",
    max.iter=100
  )

print(m.ml.r)
```

Jetzt schauen wir mal die Faktorenladungen an und unterdrücken in der Ausgabe alle Werte, die kleiner gleich .3 sind und lassen nach der Grösse der Faktorladungen sortieren:

```{r}
psych::print.psych(m.ml.r, cut=0.3, sort=T)
```

Man kann auch Diagramme rauslassen, die allerdings noch etwas Nacharbeit brauchen würden, bis sie gut lesbar sind. Dazu sollte man im Fall die Vignette von psych konsultieren.

```{r}
m.ml.r |> 
  psych::fa.diagram()
```

## 3.7 Führen Sie eine FA mit Oblimin-Rotation durch {.unnumbered .unlisted}

```{r}

m.ml.oblimin <- psych::fa(items, 
    nfactors=4,
    n.obs = nrow(items),
  SMC=TRUE,
  fm="ml",
    rotate="promax",
    max.iter=100
  )

psych::print.psych(m.ml.oblimin, cut=0.3, sort=T)
```

## Faktoren in die Daten speichern {.unnumbered}

Das Speichern ist nicht so einfach, weil die Faktoren der Faktoranalyse nur als Vectoren ausgegeben werden und dann nur für die Fälle, die keine fehlenden Werte in den Variablen der Faktorenanlyse haben. Darum muss man (wenn man fehlende Werte hat) erstmal einen Hilfsdatensatz erstellen mit den Fall-IDs (bei uns die CASES). Dann bindet man die neuen Faktor\$scores mit column bind `cbind()` zusammen und kann das dann mit `left_join` an den ursprünglichen Datensatz matchen.

```{r}

Cases <- DATEN |> 
  select(CASE, R101_01:R101_24) |> 
  na.omit() 
  
Faktoren <- items |> 
  psych::factor.scores(f = m.ml.r, method = "Thurstone")

Faktoren_der_Faelle <- cbind(Cases, Faktoren$scores)

DATEN <- DATEN |> 
  left_join(Faktoren_der_Faelle, by = "CASE")

DATEN |> saveRDS("data/Stat_A-Befragung_mit_FA.RDS")
```
