quantitative-methods-hs22/report.rmd

---
title: "Quantitative Methods HS22"
author: "Marc Gauch"
date: "`r Sys.Date()`"
output:
  html_document:
    toc: yes
    toc_depth: 4
    toc_float: yes
    number_sections: yes
  pdf_document:
    toc: yes
    toc_depth: '2'
---

```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE)
```

# Preparation

```{r, message=FALSE}
if (!require(tidyverse)) {
  install.packages("tidyverse")
  library(tidyverse)
}
if (!require(moments)) {
  install.packages("moments")
  library(moments)
}
```

# Own functions

## Frequency Table ordered from wish.com
```{r}
freq <- function(data, rounded_digits = 2) {
  # counts
  total_count <- length(data)
  na_count <- length(data[is.na(data)])
  valid_count <- total_count - na_count
  
  frequency <- table(data)
  p <- prop.table(frequency)
  valid_percent <- round(p * 100, digits = rounded_digits)
  na_percent <- round(na_count / length(data) * 100, digits = rounded_digits)
  
  percent_raw <- frequency/total_count*100
  percent <- round(percent_raw, digits = rounded_digits)
  cumulative_percent <- round(cumsum(p) * 100, digits = rounded_digits)
  freq_table <- cbind(frequency, percent, percent_raw, valid_percent, cumulative_percent)

  valid_percent_sum <- round(sum(as.data.frame(freq_table)$percent_raw), digits = rounded_digits)
  Valid_Total <- c(valid_count, valid_percent_sum, 100, NaN)
  
  freq_table <- subset(freq_table, select = -c(percent_raw))
  
  NAs <- c(na_count, na_percent, NaN, NaN)
  Total <- c(total_count, 100, NaN, NaN)
  print(rbind(freq_table, Valid_Total, NAs, Total))
}
```
*Source:  https://tellmi.psy.lmu.de/tutorials/deskriptive-statistiken-und-grafiken.html#haeufigkeiten-diskret and adapted*

## Modus
```{r}
getmode <- function(v) {
  uniqv <- unique(v)
  x <- tabulate(match(v, uniqv))
  uniqv[which(x == max(x))]
}
```

# Load Data
## Load from CSV
``` {r loadData}
litdata <- read_csv("DataLit_R.csv", show_col_types = FALSE)
litdata <- as_tibble(litdata)
```

## First inspection of data {.tabset}
### Summary
``` {r}
summary(litdata)
```
### Glimpse
``` {r}
glimpse(litdata)
```
### Print
``` {r}
print(litdata)
```
### Head
``` {r}
head(litdata)
```

# Data cleaning
## Converting Strings to numbers and *Keine Antwort* zu *NaN*
``` {r}
litdata <- litdata %>%
  mutate_all(~ replace(., . == "Stimme voll zu5", 5)) %>%
  mutate_all(~ replace(., . == "Stimme überhaupt nicht zu1", 1)) %>%
  mutate_all(~ replace(., . == "Keine Antwort-", NaN))
```

## Make it numeric
The following code will **NOT** be run. The Idea is to show a way to automatically edit all columns. It works but some columns are NOT numeric.
```{r, eval=FALSE}
# All colnames that exist
litdataColnames <- colnames(litdata)
# the ones we don't want to change
litdataNonNumericCols <- c("submitdate", "startlanguage", "startdate", "datestamp", "lastpage", "seed")
# the colnames that should be changed
litdataColsToMakeNumeric <- litdataColnames[!(litdataColnames %in% litdataNonNumericCols)]
print(litdataColsToMakeNumeric)
litdataColsToMakeNumeric <- c("R1")
for (col in litdataColsToMakeNumeric) {
  litdata[[col]] <- as.numeric(litdata[[col]])
}
```

First we rename all the columns
```{r}
litdata <- litdata %>% rename(
  "A1" = "W001",
  "A2" = "W002",
  "A3" = "W003",
  "A4" = "W004",
  "A5" = "W005",
  "A6" = "W006",
  "A7" = "W007",
  "A8" = "W008",
  "A9" = "W009",
  "B1" = "K001",
  "B2" = "K002",
  "B3" = "K003",
  "B4" = "K004",
  "B5" = "K005",
  "B6" = "K006",
  "B7" = "K007",
  "B8" = "K008",
  "B9" = "K009",
  "C1_1" = "TK001_01",
  "C1_2" = "TK001_02",
  "C1_3" = "TK001_03",
  "C1_4" = "TK001_04",
  "C2_1" = "TK002_01",
  "C2_2" = "TK002_02",
  "C2_3" = "TK002_03",
  "C2_4" = "TK002_04",
  "C3_1" = "TK003_01",
  "C3_2" = "TK003_02",
  "C3_3" = "TK003_03",
  "C3_4" = "TK003_04",
  "C4_1" = "TK004_01",
  "C4_2" = "TK004_02",
  "C4_3" = "TK004_03",
  "C4_4" = "TK004_04",
  "C5_1" = "TK005_01",
  "C5_2" = "TK005_02",
  "C5_3" = "TK005_03",
  "C5_4" = "TK005_04",
  "C6_1" = "TK006_01",
  "C6_2" = "TK006_02",
  "C6_3" = "TK006_03",
  "C6_4" = "TK006_04",
  "D1_1" = "H001_001",
  "D1_2" = "H001_002",
  "D1_3" = "H001_003",
  "D1_4" = "H001_004",
  "D1_5" = "H001_005",
  "D1_6" = "H001_006",
  "D1_7" = "H001_007",
  "D2" = "H002",
  "D3" = "H003",
  "D4" = "H004",
  "D4_comment" = "H004_other",
  "D5" = "H005",
  "D5_comment" = "H005_other",
  "D6" = "H006",
  "D7" = "H007",
  "D8" = "H008",
  "E1" = "R1"
)
```
Then we change the datatype and fix the values
```{r}
litdata$A1 <- as.numeric(litdata$A1)
litdata$A2 <- as.numeric(litdata$A2)
litdata$A3 <- as.numeric(litdata$A3)
litdata$A4 <- as.numeric(litdata$A4)
litdata$A5 <- as.numeric(litdata$A5)
litdata$A6 <- as.numeric(litdata$A6)
litdata$A7 <- as.numeric(litdata$A7)
litdata$A8 <- as.numeric(litdata$A8)
litdata$A9 <- as.numeric(litdata$A9)

litdata$B1 <- as.numeric(litdata$B1)
litdata$B2 <- as.numeric(litdata$B2)
litdata$B3 <- as.numeric(litdata$B3)
litdata$B4 <- as.numeric(litdata$B4)
litdata$B5 <- as.numeric(litdata$B5)
litdata$B6 <- as.numeric(litdata$B6)
litdata$B7 <- as.numeric(litdata$B7)
litdata$B8 <- as.numeric(litdata$B8)
litdata$B9 <- as.numeric(litdata$B9)

litdata$C1_1 <- as.numeric(litdata$C1_1)
litdata$C1_2 <- as.numeric(litdata$C1_2)
litdata$C1_3 <- as.numeric(litdata$C1_3)
litdata$C1_4 <- as.numeric(litdata$C1_4)
litdata$C2_1 <- as.numeric(litdata$C2_1)
litdata$C2_2 <- as.numeric(litdata$C2_2)
litdata$C2_3 <- as.numeric(litdata$C2_3)
litdata$C2_4 <- as.numeric(litdata$C2_4)
litdata$C3_1 <- as.numeric(litdata$C3_1)
litdata$C3_2 <- as.numeric(litdata$C3_2)
litdata$C3_3 <- as.numeric(litdata$C3_3)
litdata$C3_4 <- as.numeric(litdata$C3_4)
litdata$C4_1 <- as.numeric(litdata$C4_1)
litdata$C4_2 <- as.numeric(litdata$C4_2)
litdata$C4_3 <- as.numeric(litdata$C4_3)
litdata$C4_4 <- as.numeric(litdata$C4_4)
litdata$C5_1 <- as.numeric(litdata$C5_1)
litdata$C5_2 <- as.numeric(litdata$C5_2)
litdata$C5_3 <- as.numeric(litdata$C5_3)
litdata$C5_4 <- as.numeric(litdata$C5_4)
litdata$C6_1 <- as.numeric(litdata$C6_1)
litdata$C6_2 <- as.numeric(litdata$C6_2)
litdata$C6_3 <- as.numeric(litdata$C6_3)
litdata$C6_4 <- as.numeric(litdata$C6_4)

litdata <- litdata %>% mutate(D1_1 = ifelse(D1_1 == "Ja", TRUE, ifelse(D1_1 == "Nicht Gewählt", FALSE, D1_1)))
litdata$D1_1 <- as.logical(litdata$D1_1)

litdata <- litdata %>% mutate(D1_2 = ifelse(D1_2 == "Ja", TRUE, ifelse(D1_2 == "Nicht Gewählt", FALSE, D1_2)))
litdata$D1_2 <- as.logical(litdata$D1_2)

litdata <- litdata %>% mutate(D1_3 = ifelse(D1_3 == "Ja", TRUE, ifelse(D1_3 == "Nicht Gewählt", FALSE, D1_3)))
litdata$D1_3 <- as.logical(litdata$D1_3)


litdata <- litdata %>% mutate(D1_4 = ifelse(D1_4 == "Ja", TRUE, ifelse(D1_4 == "Nicht Gewählt", FALSE, D1_4)))
litdata$D1_4 <- as.logical(litdata$D1_4)

litdata <- litdata %>% mutate(D1_5 = ifelse(D1_5 == "Ja", TRUE, ifelse(D1_5 == "Nicht Gewählt", FALSE, D1_5)))
litdata$D1_5 <- as.logical(litdata$D1_5)

litdata <- litdata %>% mutate(D1_6 = ifelse(D1_6 == "Ja", TRUE, ifelse(D1_6 == "Nicht Gewählt", FALSE, D1_6)))
litdata$D1_6 <- as.logical(litdata$D1_6)

litdata <- litdata %>% mutate(D1_7 = ifelse(D1_7 == "Ja", TRUE, ifelse(D1_7 == "Nicht Gewählt", FALSE, D1_7)))
litdata$D1_7 <- as.logical(litdata$D1_7)

litdata <- litdata %>% mutate(D2 = ifelse(D2 == "Ja", TRUE, ifelse(D2 == "Nein", FALSE, D2)))
litdata$D2 <- as.logical(litdata$D2)

# skipping D3 because it's just a free text

litdata$D4 <- as.factor(litdata$D4)

# skipping D4_comment because it's a free text

litdata$D5 <- as.factor(litdata$D5)

# skipping D5_comment because it's a free text

# can't be a number as there is a 2010 or earlier option.
litdata$D6 <- as.factor(litdata$D6)

litdata$D7 <- as.numeric(litdata$D7)

litdata$D8 <- as.factor(litdata$D8)

# skipping E1 because it's a free text
```


## Second inspection of data {.tabset}
### Summary
``` {r}
summary(litdata)
```
### Glimpse
``` {r}
glimpse(litdata)
```
### Print
``` {r}
print(litdata)
```
### Head
``` {r}
head(litdata)
```

# Selbststudium 1
*Berechnen Sie die Häufigkeiten für die Variablen W003, K003, H001_001, H005, H007 und H008.*

## Data {.tabset}
```{r}
displayFunction1 <- function(table, column) {
  tmp <- table[column]
  tmp <- rename(tmp, value = all_of(column))
  tmp <- tmp %>%
    count(value) %>%
    mutate(percentage = prop.table(n) * 100)
  print(tmp, n = 100)
  ggplot(
    tmp,
    aes(x = value, y = n)
  ) +
    geom_bar(stat = "identity") +
    theme(axis.text.x = element_text(angle = 45, hjust = 1))
}
```
### A3 (W003) {-}
```{r}
displayFunction1(litdata, "A3")
```

### B3 (K003) {-}
```{r}
displayFunction1(litdata, "B3")
```

### D1_1 (H001_001) {-}
```{r}
displayFunction1(litdata, "D1_1")
```

### D5 (H005) {-}
```{r}
displayFunction1(litdata, "D5")
```

### D7 (H007) {-}
```{r}
displayFunction1(litdata, "D7")
```

Die Warnung resultiert daraus, dass es sehr viele *NA* gibt.

### D8 (H008) {-}
```{r}
displayFunction1(litdata, "D8")
```

# Selbststudium 2.1
We have the year 2021
## Preparation
```{r}
birthyears <- litdata$D7
# remove NAs
birthyears <- birthyears[!is.na(birthyears)]
age <- 2021 - birthyears
```

## Frequency
```{r}
freq(age)
# with NA
freq(2021 - litdata$D7)
```

## Selbststudium 2
https://www.beratung-statistik.de/statistik-beratung-infos/r-tutorial/deskriptive-statistik-r/

### Modalwert
```{r}
# own method
getmode(age)
```

### Median
```{r}
median(age)
```

### Arithmetischer Mittelwert
```{r}
mean(age)
```

### Spannweite
```{r}
max(age) - min(age)
```

### Quartilsabstand
```{r}
IQR(age)
# just for fun
summary(age)
```

### Varianz
```{r}
var(age)
```

### Standardabweichung
```{r}
sd(age)
```

### Schiefe
```{r}
library(moments)
skewness(age)
```
> Die Kennzahl Schiefe ist wird Null bei einer perfekt symmetrischen Verteilung, größer als Null bei einer rechtsschiefen und kleiner als Null bei einer linksschiefen Verteilung. 
https://www.beratung-statistik.de/statistik-beratung-infos/r-tutorial/deskriptive-statistik-r/

### Kurtosis
```{r}
kurtosis(age)
```
>  Eine weitere bekannte Kennzahl ist die Kurtosis. Um eine Vorstellung von der Bedeutung der Kurtosis zu erhalten, betrachten Sie nachfolgende Graphik.  
In dieser Graphik sind eine Normalverteilung, sowie eine steilgipflige (aka leptokurtisch) und eine flachgipflige (aka platykurtisch) dargestellt.
Die steilgipflige Verteilung ist in der Mitte spitzer als die Normalverteilung und an den Rändern breiter. Bei der flachgipligen Verteilung ist es anders herum. Die Kurtosis ist nun eine Kennzahl, mit der untersucht wird, ob eine Verteilung im Vergleich zur Normalverteilung flachgipflig oder steilgipflig ist:  
- Für eine Normalverteilung nimmt die Kurtosis genau den Wert 3 an.  
- Eine steilgipflige Verteilung hat eine Kurtosis, die größer als 3 ist.  
- Für eine flachgipflige Verteilung ist die Kurtosis kleiner als 3.  
- Beachten Sie: Anstatt der Kurtosis wird häufig auch der sogenannte Exzess verwendet. Dies ist eine weitere Kennzahl, die definiert ist durch die Formel: Exzess = Kurtosis - 3.   
- Der Exzess ist somit größer als Null, wenn die Verteilung steilgipflig ist, und kleiner als Null bei einer flachgipfligen Verteilung.  
![Abbildung Kurtosis](images/kurtosis.png)  
Frech kopiert von: https://www.beratung-statistik.de/statistik-beratung-infos/r-tutorial/deskriptive-statistik-r/

### QQ-Plot
```{r}
qqnorm(age)
qqline(age)
```

### Historam for age {.tabset}

#### Frequency {-}
```{r}
hist(age, freq = F)
lines(density(age), lwd = 2, col = "black")
```

#### Auto Breaks {-}
```{r}
hist(age)
```

#### 3 Breaks {-}
```{r}
hist(age, breaks = 3)
```

#### 5 Breaks {-}
```{r}
hist(age, breaks = 5)
```

#### 7 Breaks {-}
```{r}
hist(age, breaks = 7)
```

#### 10 Breaks {-}
```{r}
hist(age, breaks = 10)
```

#### 15 Breaks {-}
```{r}
hist(age, breaks = 15)
```

#### 20 Breaks {-}
```{r}
hist(age, breaks = 20)
```

#### 30 Breaks {-}
```{r}
hist(age, breaks = 30)
```

# Selbststudium 2.2

Auf der Grundlage von Daten aus einer Schweizer Schüllererhebung wird aus verschiedenen Variablen (z.B. Angaben zum Beruf der Eltern, zur Elternausbildung sowie zur Anzahl von Bücchern zu Hause) ein Index zur sozialen Herkunft erstellt. Dieser Index erscheint in einer neu gebildeten numerischen Variable im Datensatz, gibt also für jeden Fall in diesem Datensatz einen Skalenwert zur sozialen Herkunft an. Die neu gebildete Skala läuft von 0 (Wert mit der geringsten Ausprägung) bis 10 (Wert mit der höchsten Ausprägung).

Für die gesamte Schweiz liegt der arithmetische Mittelwert auf dieser Skala bei 5.6 und die Standardabweichung beträgt 1.8 (Zahlen sind von mir frei erfunden!). Die Verteilung entspricht einer Normalverteilung.

Der Mittelwert der Verteilung der Bündner Schüler liegt etwas tiefer als in der Gesamtschweiz, nämlich bei 5.1 mit einer Standardabweichung von 2.

## Aufgabenstellung 1

Wo in der Verteilung der Schweiz liegt der Bündner Mittelwert, bzw. wie viele Schweizer Schüler haben bzgl. der sozialen Herkunft einen tieferen Wert als der typische Bündner Schüler?

```{r}
pnorm(5.1, mean = 5.6, sd = 1.8)
```

## Aufgabenstellung 2

In einem Bündner Ort beträgt der Mittelwert auf der Skala zur sozialen Herkunft 5.6, er ist also genauso hoch wie in der Gesamtschweiz. Was ist zur Lage dieses Wertes bezogen auf die Verteilung in Graubünden zu sagen? Also: Wie viele Schüler in GR liegen mit ihrem Wert darunter?

```{r}
pnorm(5.6, mean = 5.1, sd = 2.0)
```