Musterlösung
# Speicherort sowie Dateiname anpassen
<- read_delim("./HIER RELATIVEN DATEIPFAD EINGEBEN", "HIER SEPERATOR EINGEBEN") depo
Die Projektstruktur steht. Nun können die Daten eingelesen und die nötigen Datentypen definiert werden.
Lädt die Daten zuerst von Moodle herunter:
Zähldaten zu eurem Standort (211_sihlwaldstrasse_2017_2024.csv, 502_sihluferweg_2016_2024.csv)
Meteodaten und Legende
Hinweis: Siehe [Einführung] für den Standort der Zähler 211 und 502.
Die Zähldaten des WPZ wurden vorgängig bereinigt (z.B. wurden Stundenwerte entfernt, an denen am Zähler Wartungsarbeiten stattgefunden haben). Das macht es für uns einfach, denn wir können die Daten ohne vorgängige Bereinigung einlesen. Behaltet aber im Hinterkopf, dass die Datenaufbereitung, die Datenbereinigung mit viel Aufwand verbunden ist.
# Speicherort sowie Dateiname anpassen
<- read_delim("./HIER RELATIVEN DATEIPFAD EINGEBEN", "HIER SEPERATOR EINGEBEN") depo
Hinweis: Im Stundenformat zeigen die Werte bei 11:00 die Zähldaten zwischen 11:00 bis 12:00 Uhr.
# lese die Daten ein
<- read_delim("datasets/fallstudie_s/211_sihlwaldstrasse.csv", ";")
depo
# erstes Sichten und anpassen der Datentypen
str(depo)
<- depo |>
depo mutate(
Datetime = as.POSIXct(DatumUhrzeit, format = "HIER STEHT DAS DATUMSFORMAT", tz = "CET"),
# nun schreiben wir uns das Datum in eine seperate Spalte
Datum = as.Date(Datetime)
)
# hier der code mit dem richtigen Format
<- depo |>
depo mutate(
Datetime = as.POSIXct(as.character(Datetime), format = "%Y%m%d%H", tz = "CET"),
Datum = as.Date(Datetime)
)
Ihr könnt selbst wählen, ob ihr Fussgänger:innen oder Fahrräder untersuchen wollt (je nachdem ob sie in eurem Datensatz vorhanden sind).
Hinweis: mit select() können Spalten gewählt werden, mit filter() Zeilen.
# In dieser Auswertung werden nur Personen zu Fuss betrachtet!
# it select werden spalten ausgewaehlt oder eben fallengelassen
<- depo |>
depo ::select(-c(Velo_IN, Velo_OUT)) dplyr
Tipp: Wenn man R sagt: “addiere mir Spalte x mit Spalte y”, dann macht R das für alle Zeilen in diesen zwei Spalten. Wenn man nun noch sagt: “speichere mir das Ergebnis dieser Addition in einer neuen Spalte namens Total”, dann hat man die Aufgabe bereits gelöst. Arbeitet mit mutate()).
# Berechnen des Totals, da dieses in den Daten nicht vorhanden ist
<- depo |>
depo mutate(Total = Fuss_IN + Fuss_OUT)
# Entferne die NA's in dem df.
<- na.omit(depo) depo
# Einlesen
<- read_delim("datasets/fallstudie_s/order_124839_data.txt", ";") meteo
Tipp: Das Datum wird als Integer erkannt. Zuerst muss es in Text umgewandelt werden aus dem dann das eigentliche Datum herausgelesen werden kann. Das ist mühsam - darum hier der Code.
<- mutate(meteo, time = as.Date(as.character(time), "%Y%m%d")) meteo
Hinweise:
Die Zeitangaben sind in UTC: 00:40 UTC = 02:40 Sommerzeit = 01:40 Winterzeit, Beispiel: 13 = beinhaltet Messperiode von 12:01 bis 13:00
Da wir mit Tageshöchstwerten oder -summen rechnen, können wir zum Glück ignorieren, dass das nicht mit den Daten der Zählstellen übereinstimmt. Learning: es ist zentral immer die Metadaten zu checken.
Was ist eigentlich Niederschlag:
Werden den anderen Spalten die richtigen Typen zugewiesen? Falls nicht, ändert die Datentypen.
Nun schneiden wir den Datensatz auf die Untersuchungsdauer zu.
|>
... filter(time >= depo_start, time <= depo_end)
Dann müssen auch hier alle nicht verfügbare Werte (NA’s) herausgefiltert werden. Macht das wieder mit na.omit()
Prüft nun, wie die Struktur des data.frame (df) aussieht und ob alle NA Werte entfernt wurden:
sum(is.na(df$Variable))
# Die eigentlichen Messwerte sind alle nummerisch
<- meteo |>
meteo mutate(
tre200nx = as.numeric(tre200nx),
tre200jx = as.numeric(tre200jx),
rre150n0 = as.numeric(rre150n0),
rre150j0 = as.numeric(rre150j0),
sremaxdv = as.numeric(sremaxdv)
|>
) filter(time >= depo_start, time <= depo_end) # schneide dann auf Untersuchungsdauer
# Was ist eigentlich Niederschlag:
# https://www.meteoschweiz.admin.ch/home/wetter/wetterbegriffe/niederschlag.html
# Filtere Werte mit NA
<- na.omit(meteo)
meteo # Pruefe ob alles funktioniert hat
str(meteo)
sum(is.na(meteo)) # zeigt die Anzahl NA's im data.frame an
Jetzt fügen wir viele Convenience Variablen hinzu. Wir brauchen:
Der Code dazu könnte so aussehen:
|>
...mutate(Wochenende = ifelse(Wochentag %in% c(6,7), "Wochenende", "Werktag")) |>
# 1 means Monday and 7 means Sunday (default)
mutate(Wochenende = as.factor(Wochenende)) |>
...
je als Faktor:
<- depo |>
depo # wday sortiert die Wochentage automatisch in der richtigen Reihenfolge
mutate(
Wochentag = wday(Datetime, week_start = 1),
Wochentag = factor(Wochentag),
# Werktag oder Wochenende hinzufuegen
Wochenende = ifelse(Wochentag %in% c(6, 7), "Wochenende", "Werktag"),
Wochenende = as.factor(Wochenende),
# Kalenderwoche hinzufuegen
KW = isoweek(Datetime),
KW = factor(KW),
# monat und Jahr
Monat = month(Datetime),
Monat = factor(Monat),
Jahr = year(Datetime),
Jahr = factor(Jahr))
Dies machen wir auch mit dem “meteo” Datensatz. Wiederum bitte Wochentag, Werktag oder Wochenende, Kalenderwoche, Monat und Jahr. Ebenfalls alles als Faktor speichern.
# Wir gruppieren die Meteodaten noch nach Kalenderwoche und Werktag / Wochenende
# Dafür brauchen wir zuerst diese als Convenience Variablen
<- meteo |>
meteo # wday sortiert die Wochentage automatisch in der richtigen Reihenfolge
mutate(
Wochentag = wday(time, week_start = 1),
Wochentag = factor(Wochentag),
# Werktag oder Wochenende hinzufuegen
Wochenende = ifelse(Wochentag %in% c(6, 7), "Wochenende", "Werktag"),
Wochenende = as.factor(Wochenende),
# Kalenderwoche hinzufuegen
KW = isoweek(time),
KW = factor(KW),
# monat und Jahr
Monat = month(time),
Monat = factor(Monat),
Jahr = year(time),
Jahr = factor(Jahr))
Wieder zurück zum depo-Datensazt.
Ich mache den folgenden Punkt nachgelagert zu den voerherigen Convenience Variablen, da zu viele Operationen in einem Schritt auch schon mal etwas durcheinander erzeugen können.
Phasen der Covid-Pandemie (Code untenstehend). Wir definieren 5 Phasen:
von Anfang Untersuchungsperiode bis vor Lockdown 1
Lockdown 1
zwischen den Lockdowns
Lockdown 2
Ende 2. Lockdown bis Ende Untersuchungsperiode
Wir packen alle Phasen in eine Spalte –> long-format ist praktischer für das plotten als wide-format.
Später im multivariaten Modell werden die Levels der Variablen per “default” alphabetisch geordnet und die Effektstärken der einzelnen Levels gegenüber dem ersten Level gerechnet. Das macht wenig Sinn, den die Levels sind nicht alphabetisch, sondern gemäss der Liste oben (später mehr dazu). Das passen wir ebenfalls an.
Hier der Code dazu:
<- depo |>
depo mutate(Phase = case_when(
< lock_1_start ~ "Pre",
Datetime >= lock_1_start & Datetime <= lock_1_end ~ "Lockdown_1",
Datetime > lock_1_end & Datetime < lock_2_start ~ "Inter",
Datetime >= lock_2_start & Datetime <= lock_2_end ~ "Lockdown_2",
Datetime > lock_2_end ~ "Post"
Datetime
))
# hat das gepklappt?!
unique(depo$Phase)
<- depo |>
depo # mit factor() koennen die levels direkt einfach selbst definiert werden.
# wichtig: speizfizieren, dass aus R base, ansonsten kommt es zu einem
# mix-up mit anderen packages
mutate(Phase = base::factor(Phase, levels = c("Pre", "Lockdown_1", "Inter", "Lockdown_2", "Post")))
str(depo)
Neben dem Lockdown können auch die Schulferien einen Einfluss auf die Besuchszahlen haben. Wir haben die Schulferien bereits als .csv eingelesen. Allerdings können wir die Schulferien nicht mit der case_when()-Funktion zuweisen, da diese mit dieser Funktion alle Vektoren im Datensatz “schulferien” verglichen werden, und nicht elementweise für jede Zeile im “depo”-Datensatz. Dies führt dazu, dass die Bedingungen nur einmal überprüft werden und dann auf den gesamten Vektor angewendet werden, anstatt Zeile für Zeile.
# schreibe nun eine Funktion zur zuweisung Ferien. WENN groesser als start UND kleiner als
# ende, DANN schreibe ein 1
for (i in 1:nrow(schulferien)) {
$Ferien[depo$Datum >= schulferien[i, "Start"] & depo$Datum <= schulferien[i, "Ende"]] <- 1
depo
}$Ferien[is.na(depo$Ferien)] <- 0
depo
# als faktor speichern
$Ferien <- factor(depo$Ferien) depo
# Fuer einige Auswertungen muss auf die Stunden als nummerischer Wert zurueckgegriffen werden
$Stunde <- hour(depo$Datetime)
depo# hour gibt uns den integer
typeof(depo$Stunde)
Die Daten wurden durch den WPZ kalibriert (Nachkommastellen). Unser späteres Modell kann nicht mit Nachkommastellen in der abhängigen Variable umgehen (später dazu mehr).
Rundet die Zähldaten in der Spalte “Total” auf 0 Nachkommastellen. Der Befehl lautet round()
Definiert sie sicherheitshalber als Integer (= Ganzzahl)
Macht das nun noch für IN und OUT.
$Total <- as.integer(round(depo$Total, digits = 0))
depo
$Fuss_IN <- as.integer(round(depo$Fuss_IN, digits = 0))
depo
$Fuss_OUT <- as.integer(round(depo$Fuss_OUT, digits = 0)) depo
Wir setzen den Fokus unserer Untersuchung auf die Veränderung der Besuchszahlen in der Abend- und Morgendämmerung sowie der Nacht. Dafür müssen wir diese tageszeitliche Einteilung der Daten erst machen. Da dies über den Umfang dieser Fallstudie hinaus geht, liefere ich euch hier den Code dazu.
Die wichtigsten Punkte:
Die Tageslänge wurde für den Standort Zürich (Zeitzone CET) mit dem Package “suncalc” berechnet. Dabei wurden Sommer- und Winterzeit berücksichtigt.
Die Einteilung der Tageszeit beruht auf dem Start und dem Ende der astronomischen Dämmerung sowie der Golden Hour. Der Morgen und der Abend wurden nach dieser Definition berechnet und um je eine Stunde Richtung Tag verlängert.
Untenstehenden Code könnt ihr einfach kopieren.
Beschreibt in eurem Bericht später, dass ihr die Einteilung der Tageszeit gemäss den Dämmerungszeiten in Zürich und gemäss meinem Code gemacht habt.
Hinweis: damit case_when() funktioniert, müsst ihr dplyr Version als 1.1.1 oder neuer haben. Das könnt ihr unter “Packages” (neben dem Reiter “Plots”, unten rechts) prüfen.
# Einteilung Standort Zuerich
<- 47.38598
Latitude <- 8.50806
Longitude
# Start und das Ende der Sommerzeit:
# https://www.schulferien.org/schweiz/zeit/zeitumstellung/
# Welche Zeitzone haben wir eigentlich?
# Switzerland uses Central European Time (CET) during the winter as standard time,
# which is one hour ahead of Coordinated Universal Time (UTC+01:00), and
# Central European Summer Time (CEST) during the summer as daylight saving time,
# which is two hours ahead of Coordinated Universal Time (UTC+02:00).
# https://en.wikipedia.org/wiki/Time_in_Switzerland
# Was sind Astronomische Dämmerung und Golden Hour ueberhaupt?
# https://sunrisesunset.de/sonne/schweiz/zurich-kreis-1-city/
# https://www.rdocumentation.org/packages/suncalc/versions/0.5.0/topics/getSunlightTimes
# Wir arbeiten mit folgenden Variablen:
# "nightEnd" : night ends (morning astronomical twilight starts)
# "goldenHourEnd" : morning golden hour (soft light, best time for photography) ends
# "goldenHour" : evening golden hour starts
# "night" : night starts (dark enough for astronomical observations)
<-
lumidata getSunlightTimes(
date = seq.Date(depo_start, depo_end, by = 1),
keep = c("nightEnd", "goldenHourEnd", "goldenHour", "night"),
lat = Latitude,
lon = Longitude,
tz = "CET"
|>
) as_tibble()
# jetzt haben wir alle noetigen Angaben zu Sonnenaufgang, Tageslaenge usw.
# diese Angaben koennen wir nun mit unseren Zaehldaten verbinden:
<- depo |>
depo left_join(lumidata, by = c(Datum = "date"))
<- depo |>
depo mutate(Tageszeit = case_when(
>= nightEnd & Datetime <= goldenHourEnd ~ "Morgen",
Datetime > goldenHourEnd & Datetime < goldenHour ~ "Tag",
Datetime >= goldenHour & Datetime <= night ~ "Abend",
Datetime .default = "Nacht"
|>
)) mutate(Tageszeit = factor(Tageszeit, levels = c("Morgen", "Tag", "Abend", "Nacht"), ordered = TRUE))
# behalte die relevanten Var
<- depo |> dplyr::select(-nightEnd, -goldenHourEnd, -goldenHour, -night, -lat, -lon)
depo
# Plotte zum pruefn ob das funktioniert hat
ggplot(depo, aes(y = Datetime, color = Tageszeit, x = Stunde)) +
geom_jitter() +
scale_color_manual(values = mycolors)
sum(is.na(depo))
# bei mir hat der Zusatz der Tageszeit noch zu einigen NA-Wertren gefueht.
# Diese loesche ich einfach:
<- na.omit(depo)
depo # hat das funktioniert?
sum(is.na(depo))
Unsere Daten liegen, wie ihr wisst, im Stundenformat vor. Für einige Auswertungen müssen wir aber auf ganze Tage zurückgreifen.
Hinweis: Wir gruppieren nur nach Datum, da ich mit den vielen weiteren Gruppierungen hier Probleme hatte, eine korrekte Summe zu erhalten.
<- depo |>
depo_d group_by(Datum) |> # Gruppieren nach der Variable Datum
summarise(Total = sum(Fuss_IN + Fuss_OUT),# Berechnen der gewünschten Werte
Fuss_IN = sum(Fuss_IN),
...
# hier werden also pro Nutzergruppe und Richtung die Stundenwerte pro Tag aufsummiert
<- depo |>
depo_d group_by(Datum) |>
summarise(
Total = sum(Fuss_IN + Fuss_OUT),
Fuss_IN = sum(Fuss_IN),
Fuss_OUT = sum(Fuss_OUT)
)
<- depo_d |>
depo_d mutate(Tage_bis_Neujahr = as.numeric(difftime(ymd(paste0(year(Datum), "-12-31")), Datum, units = "days")))
<- depo_d |>
depo_d mutate(
Wochentag = wday(Datum, week_start = 1),
Wochentag = factor(Wochentag),
# Werktag oder Wochenende hinzufuegen
Wochenende = ifelse(Wochentag %in% c(6, 7), "Wochenende", "Werktag"),
Wochenende = as.factor(Wochenende),
# Kalenderwoche hinzufuegen
KW = isoweek(Datum),
KW = factor(KW),
# monat und Jahr
Monat = month(Datum),
Monat = factor(Monat),
Jahr = year(Datum),
Jahr = factor(Jahr))
<- depo_d |>
depo_d mutate(Phase = case_when(
< lock_1_start ~ "Pre",
Datum >= lock_1_start & Datum <= lock_1_end ~ "Lockdown_1",
Datum > lock_1_end & Datum < lock_2_start ~ "Inter",
Datum >= lock_2_start & Datum <= lock_2_end ~ "Lockdown_2",
Datum > lock_2_end ~ "Post"
Datum
))
<- depo_d |>
depo_d mutate(Phase = base::factor(Phase, levels = c("Pre", "Lockdown_1", "Inter", "Lockdown_2", "Post")))
for (i in 1:nrow(schulferien)) {
$Ferien[depo_d$Datum >= schulferien[i, "Start"] & depo_d$Datum <= schulferien[i, "Ende"]] <- 1
depo_d
}$Ferien[is.na(depo_d$Ferien)] <- 0
depo_d
$Ferien <- factor(depo_d$Ferien)
depo_d
# pruefe das df
head(depo_d)
<- depo |>
depo_daytime group_by(Jahr, Monat, KW, Phase, Ferien, Wochenende, Tageszeit) |>
summarise(
Total = sum(Fuss_IN + Fuss_OUT),
Fuss_IN = sum(Fuss_IN),
Fuss_OUT = sum(Fuss_OUT))
<- depo_daytime |>
mean_phase_d group_by(Phase, Tageszeit) |>
summarise(
Total = mean(Total),
IN = mean(Fuss_IN),
OUT = mean(Fuss_OUT))
Tipp: Braucht wiederum group_by() und summarise().
<- depo |>
depo_m group_by(Jahr, Monat) |>
summarise(Total = sum(Total))
Hier der fertige Code dazu (da etwas umständlich):
<- depo_m |>
depo_m mutate(
Ym = paste(Jahr, Monat), # und mache eine neue Spalte, in der Jahr und Monat sind
Ym = lubridate::ym(Ym)
# formatiere als Datum )
# Gruppiere die Werte nach Monat und TAGESZEIT
<- depo |>
depo_m_daytime group_by(Jahr, Monat, Tageszeit) |>
summarise(Total = sum(Total))
# sortiere das df aufsteigend (nur das es sicher stimmt)
<- depo_m_daytime |>
depo_m_daytime mutate(
Ym = paste(Jahr, Monat), # und mache eine neue Spalte, in der Jahr und Monat sind
Ym = lubridate::ym(Ym)
# formatiere als Datum )
Macht euch mit den Daten vertraut. Plottet sie, seht euch die df’s an, versteht, was sie repräsentieren.
Z.B. sind folgende Befehle und Plots wichtig:
str()
summarize()
head()
Scatterplot, x = Datum, y = Anzahl pro Zeiteinheit
Histrogram
usw.
Hinweis: Geht noch nicht zu weit mit euren Plots. Die Idee ist, dass man sich einen Überblick über die Daten verschafft und noch keine “analysierenden” Plots erstellt.
Nachdem nun alle Daten vorbereitet sind folgt im nächsten Schritt die deskriptive Analyse.