Labordatenanalyse mit R

Niels Schwab | Ronja Gottschalk
RLab - Skriptbasierte modulare Umweltstatistik (Universitätskolleg 2.0)
Universität Hamburg
CC BY-SA 4.0| 2020

Inhaltsverzeichnis

Ziel

In diesem Digitalen Skript erfährst Du in Ergänzung beziehungsweise Fortsetzung zu Vorbereitung einer komplexen Tabelle für den Import in R, wie bei Analysen von zum Beispiel Bodenproben oder atmosphärischer Deposition im Labor des Instituts für Geographie der Universität Hamburg gewonnene Daten mit R ausgewertet werden können. Dazu ist R-Code erforderlich, der auf den folgenden Seiten präsentiert wird.

Es werden R-Skripte für die Analyse folgender Labordaten bereit gestellt (Du kannst die Links nutzen, um schnell zum Code zu gelangen, verpasst dann aber wichtige Infos auf der nächsten Seite):

Normalerweise kann die Bedienung von R mit RLab-Kursen direkt in R gelernt werden. Die Komplexität der Ausgangsdaten und des für die Datenanalyse erforderlichen R-Codes geht bei der Labordatenanalyse über das für einen RLab-Kurs sinnvolle Maß hinaus. R-Lab-Kurse bestehen aus Schritt-für-Schritt-Anleitungen zum Aufbau eher kurzer R-Funktionen. Sie bieten sich daher eher für grundlagenorientierte, einführende Inhalte an, da das Frustrationspotential bei einem sich lange hinziehenden Schritt-für-Schritt-Aufbau groß ist.

Daher wird der R-Code hier zusammen mit umfangreichen, erklärenden Kommentaren anschaulich zu präsentieren. Der Code kann als Ausgangspunkt für die Analyse, Visualisierung und Präsentation der eigenen Daten verwendet und entsprechend angepasst werden. Damit das gut klappt, sind ein paar Voraussetzungen erforderlich, die auf der folgenden Seite aufgelistet sind.

Voraussetzungen

Das Verständnis der Inhalte dieses Digitalen Skripts setzt einige Vorkenntnisse voraus:

Du hast R und RStudio installiert und RStudio mit R verbunden.
Du kennst Dich mit dem Programm RStudio und seinen vier Fenstern aus.
Du weißt, wie eine R-Funktion aufgebaut ist und welche Tipps und Tricks bei der Nutzung von R zu beachten sind.
Gut ist auch, wenn Du etwas praktische Erfahrung mit der Nutzung von R hast, zum Beispiel durch die Bearbeitung von RLab-Kursen.
Du besitzt grundlegende Kenntnisse zum Import von Daten in R. Diese Kenntnisse kannst Du Dir aneignen, indem Du das Digitale Skript Import von Daten in R und den RLab-Kurs Daten einlesen und kennenlernen durcharbeitest.
Beim Labor-Datenmaster handelt es sich um eine recht komplexe Excel-Tabelle, die eine etwas aufwändigere Vorbereitung, als im Digitalen Skript Import von Daten in R beschrieben, erforderlich macht. Um den auf den folgenden Seiten präsentierten R-Code mit Daten aus einem Datenmaster zu nutzen, solltest Du zunächst die Tabelle anhand des entsprechenden Digitalen Skripts aufbereiten.

Schaue gegebenenfalls in den oben verlinkten Ressourcen nach und kehre dann hierher zurück! Im R-Code, der für verschiedene Anwendungsfälle auf den nächsten Seiten zur Verfügung gestellt wird, geht es immer mit dem Import des aufbereiteten Datenmasters los.

Am besten öffnest Du spätestens jetzt parallel RStudio, da im Folgenden R-Code vorgestellt wird, den Du per Copy&Paste übernehmen und dann direkt ausprobieren kannst.

Für Anregungen und Kommentare zur Verbesserung dieses Digitalen Skripts ist das RLab-Team immer dankbar! Du kannst auch Fragen zu den Inhalten stellen! Nutze für all das gerne die Kommentar-Funktion!

R-Code zur Analyse der Parameter der Bodenlösung

Hier findest Du R-Code, mit dem Du im ,,Datenmaster Bodenlösung" enthaltene Werte von verschiedenen Bodenlösungsparametern analysieren und Ergebnisse visualisieren kannst.

Der Code wirkt sehr umfangreich, was vor allem daran liegt, dass er sehr ausführlich kommentiert wurde, damit Du ganz genau nachvollziehen kannst, wie die einzelnen R-Funktionen ausgeführt werden, welche Argumente sie benötigen und was es sonst noch zu beachten gibt. Wenn Du den Code in einer Kurzversion, ohne jegliche Kommentare haben möchtest, kannst Du ihn Dir hier als R-Skript runterladen und verwenden.

Zum Ausführen des Codes erforderliche Beispieldaten

Vor der der Verwendung vorhandener Rohdaten (,,Datenmaster") in R müssen diese im Rahmen einer Vorbereitung bereinigt werden. Dies ist bei den folgenden Beispieldaten schon für Dich erledigt worden:

bereinigte Daten der Bodenlösung, wie sie im R-Code (siehe unten) verwendet werden (zum Download als .csv-Datei ggf. Rechtsklick und "Speichern unter...")

Falls Du die Vorbereitung selbst nachvollziehen möchtest, steht hier auch der Rohdatensatz zur Verfügung:

Rohdaten der Bodenlösung, die vor den R-Berechnungen mit z.B. MS-Excel bereinigt werden müssen (zum Download als .xlsx-Datei ggf. Rechtsklick und "Speichern unter...")

Code

Kurz-Version des Bodenlösungs-Skripts, ohne Kommentare (zum Download als .R-Datei ggf. Rechtsklick und "Speichern unter...")

Auch den kommentierten Code gibt es als R-Skript zum Download:

Kommentierter Code des Bodenlösungs-Skripts (zum Download als .R-Datei ggf. Rechtsklick und "Speichern unter...")

Alternativ kannst Du den kommentierten Code hier direkt aus dem Digitalen Skript übernehmen (Copy&Paste):


 # Course: RLab - Labor
 # Lesson: Bodenloesung
 # Author: RLab-Team
 # Co-Authors: Ronja Gottschalk (programming, design), Niels Schwab (design),
 #             Elke Fischer (design)
 # Organization: Universität Hamburg, project 'RLab - Skriptbasierte modulare
 #               Umweltstatistik' Funded within 'Lehrlabor' of
 #               Universitätskolleg 2.0 by BMBF (01PL17033).

 #####################################################################

 # Herzlich Willkommen im R-Skript zur Analyse der Bodenlösung.
 # Bevor Du mit der Arbeit beginnst, schließe bitte noch einen
 # Vorschritt ab, indem Du auf File --> 'Reopen with Encoding'
 # klickst. Wähle in dem dann erscheinenden Dialogfenster
 # bitte 'UTF-8' aus.
 # Damit wird das Skript erneut in der UTF-8 Kodierung geladen,
 # sodass Umlaute, etc. korrekt dargestellt werden.
 # Um dieses Skript und die darin enthaltenen Kommentare gut zu
 # verstehen, ist ein grundlegendes Verständnis für die Nutzung von
 # R und RStudio erforderlich. Solltest Du Dir diese aneignen oder
 # sie auffrischen wollen, ist die Bearbeitung der Grundlagen-Inhalte  
 # des RLab-Angebots (https://rlab.blogs.uni-hamburg.de/grundlagen/)
 # und ggf. weiterer RLab-Kurse
 # (https://rlab.blogs.uni-hamburg.de/rlab-kurse/) zu empfehlen.

 #####################################################################
 #### I. Vorarbeit ####

 # Workspace aufräumen
 rm(list = ls())

 # Kontrolle des Workspace
 ls()

 # Abfrage des aktuellen 'Working Directory'
 getwd()

 # neues Working Directory festlegen
 setwd('F:/RLabPortable/R-Labor/Bodenlösung')
 # Beachte hier bitte, dass Du Deinen Dateipfad entsprechend Deiner
 # Ordnerstruktur anpassen musst.
 # Gute Hilfen zum Einlesen von Daten erhältst Du im RLab-Kurs
 # 'Daten einlesen'.

 # erneute Abfrage des 'Working Directory' zur Kontrolle, ob die
 # Änderungen übernommen wurden
 getwd()

 # Datensatz einladen
 BL <- read.table('Bodenloesung2.csv',
                  skip = 1,
                  na.strings=c('',' ','NA'),
                  header = T,
                  sep = ';',
                  dec = ',',
                  stringsAsFactors = F)

 # Auch hier sollte der Dateipfad entsprechend Deines Verzeichnisses
 # angepasst werden, sofern Deine Datei in einem anderen Ordner liegt,
 # als Dein Working Directory eingestellt ist.
 # skip = 1 bewirkt, dass die erste Spalte wird nicht mit eingelesen
 # wird.
 # na.strings definierte leere Zellen (''), Blanks (' ') und NA-Werte
 # ('NA') für R bereits beim Einlesen als 'NA' übergeben, andernfalls
 # werden diese bloß als Zeichenkette erkannt.
 # stringsAsFactors = F verhindert, dass Variablen automatisch als
 # Faktoren übernommen werden.
 # Ggf. sep = und dec = anpassen, weitere Infos dazu im Hinweis-
 # Feld am Ende der hier verlinkten Seite:
 # https://rlab.blogs.uni-hamburg.de/dig-skripte/Vorbereitung_Datenmaster/index.html?s=Schritt%202:%20Doppelte%20Kopfzeile%20entfernen

 #####################################################################
 #### II. Überblick über den Datensatz verschaffen #####

 # Datensatz in R als Tabelle ansehen
 View(BL)
 # es gibt einige NAs und leere Zellen. Ziel: Auffüllen der leeren
 # Zellen mit NA# Struktur der einzelnen Spalten ansehen

 # Struktur des Datensatzes ansehen
 str(BL)
 # die Mehrheit der Variablen wird als 'chr', also character
 # angegeben.
 # Gewollt sind hier für die Spalten mit Zahlenwerten jedoch
 # (wie im weiteren Verlauf noch deutlich wird) numerische Variablen.
 # Dies wird im weiteren Verlauf entsprechend zugewiesen.

 #####################################################################
 #### III. Überprüfung des Datensatzes ####
 # Inhalt: 3.1 Test auf fehlende Werte (NA) - 3.2 Variablen als
 # numerisch ausweisen - 3.3 Datumsformat einstellen

 #####################################################################
 #### 3.1 'NA'- WERTE

 # Abfrage, ob für die numerischen Variablen NAs erkannt werden
 table(is.na(BL[,c(6:9,11:27)])) # False 5765 True 1102
 # Es werden Variablen als NAs erkannt.
 # Zusätzlich ist es immer ratsam eine Sichtkontrolle im Datenblatt
 # 'BL' durchzuführen, um die Daten auf eventuell leere Zellen zu
 # kontrollieren. Das würde dann darauf hindeuten, dass für eine
 # bestimmte Zeichenfolge für R nicht erkennbar ist, dass es sich
 # dabei um 'NA' handelt. In diesem Fall sollte der Datenmaster dort,
 # wo keine Zahlenwerte vorliegen, vollständig mit 'NA' befüllt sein.

 # Nachfolgend sollen dennoch kurz einige Möglichkeiten vorgestellt
 # werden, wie man bei einer großen Datenmenge spezifischer überprüfen
 # kann, ob NA-Werte vorliegen oder nicht. Dies kann ggf. zur
 # Identifikation noch nicht als 'NA' definierter Zeichenketten
 # hilfreich sein.

 # Beispielhaft soll nun anhand einer Variable mit Datenlücken
 # überprüft werden, ob auch die Datenlücken als 'NA' erkannt werden.
 is.na(BL$Aluminium.Al..mg.l.)
 # Die Logikabfrage gibt für einige Zellen ein 'True' für die Abfrage
 # nach NAs aus. NAs werden also als solche erkannt.

 # Weiterführend soll nun noch überprüft werden, ob die zu Beginn beim
 # Daten einlesen als 'NA' definierten Zeichenketten noch vorliegen
 # oder bereits als 'NA' definiert sind.

 # Abfrage der Länge der Variable
 length(BL$Aluminium.Al..mg.l.)
 # 327 Werte liegen in dieser Spalte vor

 # Überprüfen, ob leere Zellen vorliegen
 sum(BL$Aluminium.Al..mg.l. == '')
 # Die Abfrage ergibt 'NA'. Das bedeutet in diesem Fall, dass keine
 # leeren Zellen in der Spalte 'Aluminium' vorliegen.

 # Überprüfen der Anzahl an 'NA' für die Spalte der Aluminiumwerte
 sum(is.na(BL$Aluminium.Al..mg.l.))
 # Die Abfrage ergibt 31 'NA'-Werte. Anhand der Sichtkontrolle der
 # Datentabelle lässt sich diese Zahl als plausibel im Verhältnis zum
 # Gesamtumfang der Variable von 327 bewerten.

 # Vorangehend wurden Dir somit eine Möglichkeit gezeigt, wie Du bei
 # großen Datensätzen ein Gefühl dafür bekommen kannst, ob Deine
 # Daten hinsichtlich der fehlenden Werte plausibel sind. Dies ist
 # deshalb wichtig, weil Du mit den Daten nicht vernünftig rechnen
 # kannst, wenn undefinierte Zellen vorliegen oder Fehlwerte wie
 # z.B. '-999' fälschlicherweise als gültige Werte behandelt werden
 # und Dir Deine statistischen Parameter verzerren.
 # Generell lässt sich festhalten, dass es hierfür keine
 # 'Pauschal-Funktion' in R gibt, da R case-sensitiv ist
 # (also z.B. Groß- und Kleinschreibung in R bereits Unterschiede
 # machen; z.B. ist 'Test' für R ungleich 'test')
 # und insofern unendlich viele Zeichenketten korrekterweise
 # oder fälschlicherweise als 'NA' erkennen kann oder eben nicht.
 # Von daher ist es erforderlich, dass Du Dich im Vorwege mit Deinen
 # Daten auseinandersetzt, um beurteilen zu können, welche Werte
 # möglicherweise fehlerhaft von R verstanden werden.

 #####################################################################
 #### 3.2 VARIABLEN ALS NUMERISCH AUSWEISEN -
 # Umwandlung von 'chr' zu 'num'

 # Erneute Abfrage der Struktur von 'BL'
 str(BL)
 # Die Abfrage ergibt, dass nahezu alle Spalten das Format 'chr'
 # (character) haben. Das bedeutet das R die Werte als Buchstaben
 # behandelt.
 # Da hier Zahlenwerte vorliegen, ist es für weitere Operationen in
 # R erforderlich, dass Spalten welche ausschließlich Zahlenwerte
 # enthalten, auch als numerisches Format von R erkannt werden.  

 # Daher soll mittels der Funktion as.numeric an R übergeben werden,
 # welche Spalten numerisch sind.
 # In diesem Fall sind das die Spalten 6-9 und 11-27.
 BL[,c(6:9,11:27)]=as.data.frame(sapply(BL[,c(6:9,11:27)], as.numeric))
 # Die Spalte 4 wird hierbei ausgelassen, da bei den Entnahmetiefen
 # auch der Buchstabe "H" vorhanden ist, weshalb die Spalte 4 nicht
 # vollständig numerisch ist.
 # Die Warning Message: 'NAs introduced by coercion' kann ignoriert
 # werden. Sie dient lediglich dem Hinweis, dass es bei der Auswahl
 # Zellen gibt, welche NA-Werte enthalten.

 # Erneutes testen der Struktur
 str(BL)
 # Die Abfrage ergibt nun, dass die Spalten 6-9 und 11-27 das Format
 # 'num' (numeric) haben. R erkennt die numerischen Spalten nun
 # also wie gewünscht in dem Datenformat.

 #####################################################################
 #### 3.3 DATUMSFORMAT EINSTELLEN

 # Bezüglich der Datenformate in R ist nun noch ein weiterer Schritt
 # durchzuführen. Dir ist sicherlich aufgefallen, dass es eine Spalte
 # namens 'Entnahme.datum' gibt.

 # Abfrage der Struktur von 'Entnahme.datum'
 str(BL$Entnahme.datum)
 # Die Abfrage ergibt, dass die Spalte als character definiert wird.
 # Hier ist das Datumsformat gewünscht.

 # 'Entnahmedatum' als Datumsformat umschreiben und damit die
 # Datumsspalte überschreiben.
 BL$Entnahme.datum=as.Date(BL$Entnahme.datum, format='%d.%m.%Y')

 # Finales Überprüfen der Struktur
 str(BL$Entnahme.datum)
 # Die Strukturabfrage der variable ergibt, dass sie nun als
 # Datumsspalte geführt wird.

 #####################################################################
 #### IV. Grafiken generieren ####
 # Inhalt: 4.1 Vorbereitung: Aufteilen des Datensatzes - 4.2 Boxplots

 # Hier werden Boxplots erstellt, die beispielhaft Kaliumwerte der
 # unterschiedlichen Flächen 'Buche' und 'Douglasie' darstellen.
 # Die besondere Herausforderung dabei ist, dass die Daten in
 # zwei „Ebenen“ aufgeteilt, nebeneinander dargestellt werden sollen:
 # Zum einen werden die Kaliumwerte geplottet und zum anderen
 # gleichzeitig nach Entnahmetiefe differenziert.

 #####################################################################
 # 4.1 VORBEREITUNG: AUFTEILEN DES DATENSATZES

 ####
 ## Datensatz aufteilen nach 'Standort' -
 ## VARIANTE 1 mit split.dataframe()
 ####

 # Diese Variante ist gut geeignet, wenn die Daten in einer Spalte
 # angeordnet sind und es einen gemeinsamen 'Identifier' gibt,
 # nach dem aufgeteilt werden kann. In diesem Fall ist es die
 # Variable 'Standort'.
 Standortart<-split.data.frame(BL, BL$Standort)

 # Struktur des neu generierten Objekts prüfen, der Datensatz liegt
 # nun in 'gesplitteter' Form vor
 str(Standortart)

 # Die Objekte Buche und Douglasie gesondert generieren
 Buche <- Standortart$Buche
 Douglasie <- Standortart$Douglasie

 # Kurze Strukturprüfung und statistischer Überblick über die beiden
 # Objekte
 str(Buche)
 summary(Buche)
 str(Douglasie)
 summary(Douglasie)

 ####
 ## Datensatz aufteilen nach "Standort" - VARIANTE 2 "Filtern" /
 ## Subsetten mit []
 ####

 # Diese Variante eignet sich, wenn Daten extrahiert werden sollen,
 # auf die eine bestimmte Bedingung zutrifft.
 # In diesem Fall sollen die Subsets die Values für die Zeilen
 # rausschreiben, für die gilt Standort =='Buche' oder == 'Douglasie'.
 # Die Stärke dieser Funktion liegt allerdings in ihrer
 # Adaptierbarkeit, so können auch Bedingungen wie > oder < in die
 # Abfrage eingebaut werden.

 # In das Objekt Buche 2 werden nun alle Zeilen aus dem Datensatz 'BL'
 # herausgeschrieben, für die innerhalb der Spalte 'Standort'
 # zutrifft, dass sie den Standort 'Buche' haben (=='Buche')-
 Buche2 <- BL[BL[,"Standort"]=="Buche",]
 # Strukturen überprüfen
 str(Buche2)

 # Äquivalent werden in das Objekt 'Douglasie 2' alle Zeilen
 # herausgeschrieben, für die innerhalb der Spalte 'Standort'
 # zutrifft, dass sie den Standort 'Douglasie' haben (=='Douglasie')  
 Douglasie2=BL[BL[,"Standort"]=="Douglasie",]

 # Da im weiteren Verlauf mit den Objekten 'Buche' und 'Douglasie'
 # gearbeitet werden soll, können diese mit den soeben erstellten
 # Objekten 'Buche2' und 'Douglasie2' gleichgesetzt werden, da sie
 # in diesem Fall dasselbe enthalten.
 # Dies kann über die Funktion all.equal() erreicht werden.

 # Die Objekte Buche und Buche2 gleichsetzen
 all.equal(Buche,Buche2)
 all.equal(Douglasie,Douglasie2)

 #####################################################################
 #### 4.2 BOXPLOTS

 ####
 ## VARIANTE 1 - Funktion boxplot()
 ####

 # Zunächst soll mittels der bereits in R implementierten Funktion
 # boxplot() gearbeitet werden.
 # Diese bietet eine vergleichsweise einfache Befehlsstruktur, ist
 # aber in Abhängigkeit der Komplexität der Daten nur begrenzt
 # einsetzbar. Wichtig ist, dass sie nur für das Plotten vollständig
 # numerischer Variablen geeignet ist.
 # Aus diesem Grund kann die nachfolgende Funktion nur für den
 # Douglasienstandort durchgeführt werden.
 # Mit der nicht vollständig numerischen Variable 'Entnahmetiefe'
 # (H, 30, 60, 105) beim Buchenstandort ist diese z.B. nicht
 # zum Plotten der Buchendaten geeignet.

 # Zunächst werden die Farben für die unterschiedlichen zu plottenden
 # Variablenausprägungen eingestellt.
 # In diesem Fall sollen die Kaliumwerte der Douglasienfläche
 # unterteilt in die verschiedenen Entnahmetiefen geplottet werden.

 # Dazu kannst Du Dir die Struktur der Entnahmetiefe mit str()
 # ins Gedächtnis rufen.
 str(Douglasie$Entnahme.tiefe)
 # Die Entnahmetiefen am Douglasienstandort untergliedern sich
 # in '30', '60' und '105'

 # Daher sollen nun 3 Farbwerte für die nachfolgenden Plots festgelegt
 # werden.
 # Dies geschieht mit der Funktion 'rainbow'. Sie ist eine Funktion,
 # die Farbpaletten an R übergibt.
 # Diese wird dem Objekt 'cols' zugewiesen.
 # In () folgt dann das Argument für die Anzahl an Farben, in
 # diesem Fall also drei, da drei unterschiedliche Entnahmetiefen
 # gegenüber gestellt werden sollen.
 cols <- rainbow(3, s = 0.5)

 # Es folgt das Plotten der Kaliumwerte der Douglasienfläche
 # unterteilt nach Entnahmetiefen
 # Die Funktion benötigt als erstes Argument innerhalb der Klammern
 # die Werte, die auf der y-Achse abgetragen werden sollen.
 # Die Tilde ~ gibt dann an, nach welchen nachfolgenden Faktoren
 # die Plots unterteilt werden sollen.
 # In diesem Fall soll nach Entnahmetiefen unterschieden werden.
 # Da die jeweiligen Entnahmetiefen jedoch mehrfach vorkommen,
 # wird ein zusätzlicher Identifier benötigt, damit zwischen den
 # Werten für dieselbe Entnahmetiefe unterschieden werden kann.
 # Im vorliegenden Beispiel kommt jede Entnahmetiefe je Datum
 # für den Douglasienstandort nur einmal vor.
 # Daher können alle Kaliumwerte differenziert werden, wenn sie
 # sowohl nach 'Entnahmetiefe' als auch zusätzlich nach 'Datum'
 # unterteilt dargestellt werden.
 # Diese Verknüpfung sowohl nach 'Entnahmetiefe' als auch zusätzlich
 # nach 'Datum' zu unterteilen, geschieht mit +
 # Das Argument data= benötigt den Datensatz, auf den zugegriffen
 # werden soll.
 # col= ist das Argument für color, deren Farbwerte zuvor mittels
 # Rainbow festgelegt und zugewiesen wurden.
 boxplot(Douglasie$Kalium.K..mg.l. ~ Douglasie$Entnahme.tiefe +
           Douglasie$Entnahme.datum, data = Douglasie, col = cols)

 # Die Funktion legend() fügt dem Plot eine Legende hinzu  
 # 'topright' gibt an, dass die Legende in die obere rechte Ecke
 # geplottet werden soll
 # fill= gibt an, welche Farben dargestellt werden sollen und legend=
 # welche Argumente in der Legende geplottet werden sollen.
 legend('topright', fill = cols, legend = c(30,60,105), horiz = F)

 # Wie Du sehen kannst ist der Plot so noch sehr rudimentär.
 # Das Datumsformat passt noch nicht richtig.

 # Daher soll nun ein etwas komplexerer Plot erstellt werden.
 # Dieser soll direkt für die Buchenfläche erstellt werden, welche
 # die Entnahmetiefe um 'H' also einen Buchstaben erweitert.
 # Die Variable kann somit nicht numerisch sein.
 # Zunächst kannst Du Dir die Struktur des Buchendatensatzes noch
 # einmal ansehen.
 str(Buche)

 # Wie Du in der Console sehen kannst, ist das Format der Spalte
 # 'Entnahmetiefe' momentan 'chr' also Character. Damit nach der
 # Entnahmetiefe differenziert werden kann, muss das Format der
 # Spalte auf Factor geändert werden.
 # Dabei werden die Daten für den Standort 'Buche' direkt in das
 # gleichnamige Objekt geschrieben.
 Buche <- within(Buche, Entnahme.tiefe <-
                   factor(Entnahme.tiefe,
                          levels=c('H','30','60','105')))
 str(Buche$Entnahme.tiefe)

 # Ziel ist nun Spalten für jedes Datum zu erhalten. Dies geschieht
 # mit der Funktion spread() aus dem Package 'tidyr'
 # ggf. das Package installieren, ansonsten emfielht es sich die
 # nachfolgende Zeile mit einem '#' auszukommentieren

 #install.packages('tidyr')
 # Package 'tidyr' aktivieren
 library(tidyr)

 # Den Datensatz mittels der Funktion spread() aufteilen.
 Bu_Ka=spread(key=Entnahme.datum, value= Kalium.K..mg.l.,data=Buche
              [,c('Labor_Nummer','Entnahme.tiefe','Entnahme.datum',
                  'Kalium.K..mg.l.')])
 # Die Labornummer wird hier erstmals verwendet, da ein einzigartiger
 # Identifier benötigt wird, die Werte aller anderen Spalten kommen
 # über mehrere Zeilen hinweg mehrfach vor (spread() funktioniert
 # nicht, wenn kein einzigartiger Identifier vorliegt)

 # äquivalent erfolgt das Splitten nun für den Buchenstandort
 Bu_Ka=spread(key=Entnahme.datum, value= Kalium.K..mg.l.,
              data=BL[BL[,'Standort']=='Buche',c('Labor_Nummer',
                                                 'Entnahme.tiefe',
                                                 'Entnahme.datum',
                                                 'Kalium.K..mg.l.')])

 # Größe des Grafikfensters mittels par() für allgemeine
 # Grafikfenstereinstellungen und mar() für die Optionen zum 'Margin'
 # festlegen

 dev.off()
 par(mar=c(10.5,4.5,1,1))

 # Nun wird das Boxplot-Fenster, gewissermaßen der Rahmen erstellt
 # und die Textoptionen festgelegt
 boxplot(Bu_Ka[,-c(1:2)], xlim = c(0.5, ncol(Bu_Ka[,-c(1:2)])+0.5),
         boxfill=rgb(1, 1, 1, alpha=1), border=rgb(1, 1, 1, alpha=1),
         las=3,
         ylab='Buche Kalium [mg/l]')
 mtext(side = 1, 'Entnahmedatum', line = 6)
 # Jetzt werden nacheinander jeweils die Boxplots der einzelnen
 # Entnahmetiefen hinzugefügt
 # Dazu war es wie zuvor ausgeführt notwendig, den Datensatz mittels
 # spread() aufzuteilen, sodass nun auf die einzelnen Entnahmetiefen
 # zugegriffen werden kann
 # Dabei gibt das Argument add = T an, dass die Box jeweils der
 # bestehenden Grafik hinzugefügt werden soll
 # Das Argument 'boxfill' dient dazu die Farbe der einzelnen Boxen
 # festzulegen
 boxplot(Bu_Ka[which(Bu_Ka$Entnahme.tiefe=='H'), -c(1:2)],
         xaxt = 'n', add = TRUE, boxfill='purple', boxwex=0.13,
         at = 1:ncol(Bu_Ka[,-c(1:2)]) - 0.3)
 boxplot(Bu_Ka[which(Bu_Ka$Entnahme.tiefe=='30'), -c(1:2)],
         xaxt = 'n', add = TRUE, boxfill='red', boxwex=0.13,
         at = 1:ncol(Bu_Ka[,-c(1:2)]) - 0.1)
 boxplot(Bu_Ka[which(Bu_Ka$Entnahme.tiefe=='60'), -c(1:2)],
         xaxt = 'n', add = TRUE, boxfill='blue', boxwex=0.13,
         at = 1:ncol(Bu_Ka[,-c(1:2)]) + 0.1)
 boxplot(Bu_Ka[which(Bu_Ka$Entnahme.tiefe=='105'), -c(1:2)],
         xaxt = 'n', add = TRUE, boxfill='green', boxwex=0.13,
         at = 1:ncol(Bu_Ka[,-c(1:2)]) + 0.3)

 # Nun soll noch eine Legende hinzugefügt werden
 # Das Argument xpd = T legt fest, dass die Legende an die
 # 'figure-Region' angefügt werden soll und nicht in den Plot
 par(xpd=T)

 # Die Legende wird mittels 'legend()' hinzugefügt
 # Das Argument 'bottom' legt die Position der Legende am unteren
 # Grafikrahmen fest, dieses kann nach Bedarf auch in 'top' etc.
 # verändert werden
 # fill legt die Farben fest, welche in der Legende dargestellt
 # werden sollen und der Vektor in legend die Beschriftung dieser
 # Beide Abfolgen, sollten der Reihenfolge entsprechen, in der Du
 # Deine einzelnen Entnahmetiefen geplottet hast
 legend('bottom', inset=c(0.2,-0.5),fill = c('purple','red','blue'
                                             ,'green'),
        legend = c('H','30 cm','60 cm','105 cm'), horiz = T,bty='o',
        title='Entnahmetiefe')

 #####################################################################
 #### 4.2 BOXPLOTS
 ## VARIANTE 2 - Gesamtplots mit ggplot

 # Nun konntest Du bereits verschiedene Boxplots erstellen
 # Nachfolgend soll nun mit dem Package ggplot2 gearbeitet werden, um
 # Alternativen der vorherigen Plots zu erzeugen
 # ggplot2 ist ein Package für grafische Darstellungen in R, das von
 # der Funktionsweise her allerdings etwas anders aufgebaut ist,
 # als die sonstigen Funktionen für Plots in R.

 # ggf. Package ggplot2 installieren
 #install.packages('ggplot2')
 # Package ggplot2 aktivieren
 library(ggplot2)

 # erneute Abfrage der Datenstruktur vom Datensatz Buche
 str(Buche)
 head(Buche)

 # Zunächst soll eine gruppierter Boxplot erstellt werden
 # Dazu wird mittels ggplot das Package ggplot angesteuert. geom_ legt
 # dann fest, welche Art des Plots folgen soll, in unserem Fall also
 # geom_boxplot
 # aes legt dann die 'aesthetics' des Plots fest. Hier übergibst Du
 # die Daten an R, welche geplottet werden sollen und legst fest,
 # welche Daten auf welcher Achse abgetragen werden sollen.

 ggplot(Buche, aes(x=Entnahme.datum, y=Kalium.K..mg.l.,
                   fill=Entnahme.tiefe)) +
   geom_boxplot()
 # Dies ist nun ein relativ rudimentärer Plot, der das Entnahmedatum
 # noch etwas ungünstig darstellt.

 # Damit das Datum für ggplot vernünftig übernommen wird, sollte es
 # zunächst als 'factor' formatiert werden und nicht wie momentan als
 # 'chr'.
 # Dabei kannst Du die umgewandelte Variable direkt einem neuen
 # Objekt zuweisen
 Buche2 <- within(Buche, Entnahme.datum <- factor(Entnahme.datum))
 str(Buche2)

 # Nun kannst Du die obige Funktion erneut abschicken, dieses Mal
 # allerdings für das Objekt Buche2, in dem die Datumsspalte das
 # benötigte Format hat.
 ggplot(Buche2, aes(x=Entnahme.datum, y=Kalium.K..mg.l.,
                    fill=Entnahme.tiefe)) +
   geom_boxplot()
 # Nun werden die einzelnen Tage beim Datum bereits korrekt
 # abgetragen.
 # Es ergeben sich hier allerdings unterschiedliche Boxbreiten,
 # da es zum Teil vollständig fehlende Werte für einzelne
 # Tiefen / Datum gibt
 # In ggplot werden die Boxen dann bei fehlenden nachbarwerten
 # entsprechend größer dargestellt
 # Um das zu verhindern und äquivalent breite Boxen zu erhalten,
 # kannst Du die Funktion geom_boxplot um das Argument 'position'
 # erweitern
 # Indem innerhalb der Position preserve = 'single' festgelegt wird,
 # wird das Strecken der einzelnen Boxen verhindert.
 ggplot(Buche2, aes(x=Entnahme.datum, y=Kalium.K..mg.l.,
                    fill=Entnahme.tiefe)) +
   geom_boxplot(position = position_dodge(preserve = "single"))

 # Nun sieht der Plot schon sehr viel stimmiger aus.
 # Was nun noch verbessert werden kann ist die Achsenbeschriftung
 # Dazu wird die Funktion wiederum erweitert, sodass nun die Labels
 # festgelegt werden
 # Diese Erweiterung geschieht mittels '+' und wird nicht,
 # wie Du es sonst kennst innerhalb der Klammern der Funktion
 # eingefügt
 ggplot(Buche2, aes(x=Entnahme.datum, y=Kalium.K..mg.l.,
                    fill=Entnahme.tiefe)) +
   geom_boxplot(position = position_dodge(preserve = "single")) +
   xlab('Datum') + ylab('Kalium [mg/l]')

 # Nun ist die Achsenbenennung bereits angepasst, zudem kann nun
 # noch die Beschriftung verbessert dargestellt werden, indem sie
 # etwas rotiert wird
 # Dazu wird die Funktion um das Argument theme(axis.text) erweitert.
 # Somit kann die Beschriftung in R für die entsprechende Achse
 #(hier axis.x) rotiert werden.
 # Auch dieses Argument wird mittels '+' hinzugefügt.
 ggplot(Buche2, aes(x=Entnahme.datum, y=Kalium.K..mg.l.,
                    fill=Entnahme.tiefe)) +
   geom_boxplot(position = position_dodge(preserve = "single")) +
   xlab('Datum') + ylab('Kalium [mg/l]')+
   theme(axis.text.x=element_text(angle=60, hjust=1))

 # Nun fehlt der Grafik der Vollständigkeit halber noch ein Titel
 # Auch dieser wird mittels '+' hinzugefügt und über das Argument
 # 'ggtitle' zum Plot hinzugefügt
 ggplot(Buche2, aes(x=Entnahme.datum, y=Kalium.K..mg.l.,
                    fill=Entnahme.tiefe)) +
   geom_boxplot(position = position_dodge(preserve = "single")) +
   xlab('Datum') + ylab('Kalium [mg/l]')+
   theme(axis.text.x=element_text(angle=60, hjust=1))+
   ggtitle('Kaliumwerte des Buchenstandorts')

 # Nun bildet der Plt die wichtigsten Informationen ab. Weitere
 # Beschriftungen und Formatierungen können vorgenommen werden.
 # Hierzu können wird Dir empfehlen, zum einen den RLab-Kurs zu
 # ggplots 'Daten Visualisieren mit ggplot2' zu machen und
 # zum anderen für gezielte Fragestellungen zu einzelnen
 # ggplot-Funktionen das Internet zu befragen.

 # In nachfolgenden Schritt soll es nun noch darum gehen, zwei
 # ggplots in einem übergeordneten Frame zu arrangieren
 # Dazu soll zunächst noch ein äquivalenter Plot für die
 # Douglasienfläche erstellt werden. Dabei sollen sowohl der
 # Douglasienplot, als auch der Buchenplot in Objekte geschrieben
 # werden.
 # Der Plot erscheint dann allerdings nur wenn man ihn als
 # namentliches Objekt ausführt.
 str(Douglasie)
 Douglasie2 <- within(Douglasie, Entnahme.datum <-
                        factor(Entnahme.datum))
 str(Douglasie2)

 KaliumDouglasie <-
   ggplot(Douglasie2, aes(x=Entnahme.datum, y=Kalium.K..mg.l.,
                          fill=Entnahme.tiefe)) +
   geom_boxplot(position = position_dodge(preserve = "single")) +
   xlab('Datum') + ylab('Kalium [mg/l]')+
   theme(axis.text.x=element_text(angle=60, hjust=1))+
   ggtitle('Kaliumwerte des Douglasienstandorts')
 KaliumDouglasie

 KaliumBuche <-
   ggplot(Buche2, aes(x=Entnahme.datum, y=Kalium.K..mg.l.,
                      fill=Entnahme.tiefe)) +
   geom_boxplot(position = position_dodge(preserve = "single")) +
   xlab('Datum') + ylab('Kalium [mg/l]')+
   theme(axis.text.x=element_text(angle=60, hjust=1))+
   ggtitle('Kaliumwerte des Buchenstandorts')
 KaliumBuche

 # Um nun beide Plots in einem Fenster gemeinsam darzustellen,
 # werden das Packages grid & grid. extra benötigt
 # ggf. Packages installieren

 #install.packages('gridExtra')
 #install.packages('grid')

 # Packages aktivieren
 library(gridExtra)
 library(grid)

 # Nun können die beiden Plotobjekte mittels der Funktion grid.arrange
 # in einen Frame geschrieben werden
 # Hierbei kann eine Gesamtüberschrift mittels des Arguments 'top'
 # eingefügt werden.
 grid.arrange(arrangeGrob(KaliumDouglasie, KaliumBuche,
                          top=textGrob("Kaliumwerte
                                       unterteilt nach Standort",
                                       gp=gpar(fontsize=20,
                                               font = 3))))

 # Nun wurden die Kaliumwerte unterteilt nach Standort und innerhalb
 # dessen unterteilt nach Entnahmetiefe dargestellt.

 # Für einige Fragestellungen kann es interessant sein, die
 # Entnahmetiefen untereinander zu vergleichen und somit einzeln
 # größer darzustellen
 # Dafür kann es sich anbieten, einen Plot je Entnahmetiefe zu
 # erstellen. Dies lässt sich mit der folgenden ggplot-Funktion sehr
 # gut realisieren.
 # Das einzige, was in der nachfolgenden Funktion verändert wird, ist
 # die Ergänzung der bisherigen ggplot-Funktion um das Argument
 # 'facet_wrap'
 ggplot(Buche2, aes(x=Entnahme.datum, y=Kalium.K..mg.l.,
                    fill=Entnahme.tiefe)) +
   geom_boxplot()+
   facet_wrap(~Entnahme.tiefe) +
   theme(axis.text.x=element_text(angle=60, hjust=1))
 # Damit wäre je Entnahmetiefe ein Plot generiert. Alternativ kann
 # somit auch ein Plot je Datum etc. erzeugt werden.
 # Dann muss entsprechend das Argument angepasst werden.  

 # Noch ein Hinweis zum Erstellen von Grafiken in R: Wenn Du die
 # Abbildungen aus R exportieren möchtest, findest Du im RLab eine
 # gute Anleitung dazu.
 # https://rlab.blogs.uni-hamburg.de/dig-skripte/Tipps_und_Tricks/index.html

 #####################################################################
 #### 5. Datenverteilung und Voraussetzungen für weitere
 # Rechenverfahren
 # 1. Normalverteilung -
 # 2. Varianzhomogenität -
 # 3. Test auf Unterschiede (Kruscal-Test)

 #### 5.1 Normalverteilung

 # Für weitere Berechnungen aber auch für einen Eindruck über die
 # Verteilung der Daten ist es notwendig die Daten auf
 # Normalverteilung zu testen.
 # Zum einen kann dies visuell geschehen, alternativ stehen
 # Testverfahren wie der Lilliefors-Test zur Verfügung
 # Nachfolgend soll exemplarisch der Lilliefors-Test angewendet
 # werden und die Daten anschließend in Abhängigkeit des Vorliegens
 # von Normalverteilung oder nicht näher analysiert werden

 # Lillie-Test
 # ggf. Package installieren
 #install.packages("nortest")
 # Package aktivieren
 library(nortest)

 # Exemplarisch sollen nun die Kaliumwerte der beiden unterteilten
 #Datensätze und die Kaliumwerte für den Gesamtdatensatz 'BL' auf
 # Normalverteilung getestet werden.
 lillie.test(Buche$Kalium.K..mg.l.)
 # nicht normalverteilt, da p < 0.05 (p = p-value < 2.2e-16)

 lillie.test(Douglasie$Kalium.K..mg.l.)
 # nicht normalverteilt, da p < 0.05 (p-value = 0.0001524)

 lillie.test(BL$Kalium.K..mg.l.)
 # nicht normalverteilt, da p < 0.05 (p = p-value < 2.2e-16)

 #### 5.2 Test auf Varianzhomogenität

 bartlett.test(BL$Kalium.K..mg.l., g = BL$Standort)
 #Bartlett's K-squared = 267.75, df = 1, p-value < 2.2e-16,
 # also < 0.05, ergo varianzheterogen

 bartlett.test(BL$Kalium.K..mg.l., g = BL$Entnahme.tiefe)

 # Die Daten sind demnach varianzheterogen und nicht normalverteilt.
 # Demnach sind die verwendbaren Verfahren auf solche beschränkt,
 # die diese Voraussetzungen nicht benötigen. Daher wird nun der
 # Kruscal- Wallis Test angewendet.

 #### 5.3 Test auf signifikante Unterschiede zwischen den Standorten.

 # Bevor der Kruscal-Test gerechnet wird, sieh Dir noch einmal die
 # Struktur der BL- Gesamtdatensatzes, insbesondere hinsichtlich
 # des Standortes an.
 str(BL)
 # Wie Du sehen kannst, ist die Variable 'Standort' im
 # Gesamtdatensatz im Format 'chr'. Für alle weiteren Testverfahren
 # wird die Variable im Format factor benötigt.
 # Daher soll die Variable zunächst ins Faktorformat umgeformt
 # werden und anschließend die Struktur erneut geprüft werden.
 BL$Standort <- as.factor(BL$Standort)
 BL$Entnahme.tiefe <- as.factor(BL$Entnahme.tiefe)
 str(BL)

 kruskal.test(BL$Kalium.K..mg.l., g = BL$Standort)
 # Kruskal-Wallis chi-squared = 30.885, df = 1, p-value = 2.738e-08
 # das heißt p < 0.05. Somit gibt es signifikante Unterschiede
 # zwischen den zwei Standorten hinsichtlich der Kaliumwerte.
 # hierzu hättest Du jedoch keinen Kruscal-Test benötigt, da hier nur
 # zwei Variablen gegeneinander getestet werden
 # Wenn jedoch mehrere nicht normalverteilte Variablen auf
 # Unterschiede getestet werden sollen, kommt der Kruscal-Wallis
 # Test ins Spiel

 # Nun kannst Du Dir die Unterschiede der Entnahmetiefen
 # beispielhaft für den Buchenstandort ansehen.

 kruskal.test(BL$Kalium.K..mg.l., g = BL$Entnahme.tiefe)
 # p-value = 1.601e-13 d.h. < 0.05.
 # Demnach gibt es signifikante Unterschied hinsichtlich der
 # Kaliumwerte zwischen den unterschiedlichen Entnahmetiefen

 # Um zu identifizieren, wo diese Unterschiede liegen, kannst Du
 # nun den Post-Hoc Test des Kruscal-Wallis Tests machen.
 # Dazu wird zunächst das Package 'pgirmes' benötigt

 # ggf. Installieren des Packages pgirmess
 # install.packages('pgirmess')

 #Aktivieren des Packages pgirmes
 library(pgirmess)

 kruskalmc(BL$Kalium.K..mg.l. ~ BL$Entnahme.tiefe)
 #für diesen Posthoc-Test Package 'pgirmess' installieren

 # Multiple comparison test after Kruskal-Wallis
 # p.value: 0.05
 # Comparisons
 # obs.dif critical.dif difference
 # 105-30  27.935787     34.50762      FALSE
 # 105-60   2.055538     34.92945      FALSE
 # 105-H  121.040241     45.22238       TRUE
 # 30-60   29.991325     32.56268      FALSE
 # 30-H    93.104455     43.42034       TRUE
 # 60-H   123.095779     43.75633       TRUE

 # Die logische Ausgabe 'False' oder 'True' gibt nun jeweils an, ob es
 # signifikante Unterschiede hinsichtlich der Kaliumwerte für das
 # jeweilige Paar gibt.
 # So gibt es beispielsweise signifikante Unterschiede in den
 # Kaliumwerten zwischen den Entnahmetiefen 'H' und '105'.

 # Im Regelfall wird Du Unterschiede für nicht-normalverteilte
 # Variablen untersuchen. Sollten jedoch normalverteilte und
 # varianzhomogene Daten vorliegen, so kannst Du hier mittels ANOVA
 # auf Unterschiede testen. Die Funktion dafür wäre wie folgt:

 # Variable_aov <- aov (Datensatz$Variable ~
 # Datensatz$VariableunterschiedlicheGruppen)

 # summary(Variable_aov)

 # hier muss man den Umweg gehen, das Ergebnis von aov() als Objekt
 # abzuspeichern und sich hinterher die summary davon anzuschauen
 # Bei Signifikanz wäre der zugehörige Post-Hoc Test
 # TukeyHSD(aov(Datensatz$Variable ~
 # Datensatz$VariableunterschiedlicheGruppen))

 #####################################################################
 #### 6. Zusammenhänge untersuchen - Korrelationen

 # Abhängig von der eigenen Fragestellung kann es ebenso wichtig sein,
 # Variablen auf ihre statistischen Zusammenhänge zu untersuchen.

 # Dies soll wieder exemplarisch anhand der Variable der Kaliumwerte
 # durchgeführt werden.
 # Da die Kaliumwerte wie bereits zuvor getestet nicht normalverteilt
 # sind, wird hier die Korrelation nach Spearman zugrunde gelegt.
 # Bei normalverteilten Werten dürfte hier mit der Pearson-Korrelation
 # gerechnet werden.
 # Zunächst soll eine simple Korrelation zwischen Kalium und Mangan
 # ermittelt werden.
 # Da die Kaliumwerte nicht normalverteilt sind, wird die
 # Spearman-Methode verwendet.

 cor.test(BL$Kalium.K..mg.l., BL$Mangan_.Mn..mg.l.)
 # p-value < 2.2e-16, cor 0.8261355

 # Nun soll noch ein Korrelationsplot erstellt werden. Dazu wird das
 # Package 'PerformanceAnalytics' benötigt.
 # Der Einfachheit halber, werden nun die numerischen Variablen in
 # den Spalten 8,9 und 10 bis 20 ausgewertet.
 # Achte darauf, nur numerische Variablen in den Plot mit
 # einzubeziehen.

 # Ggf. Package 'PerformanceAnalytics' installieren
 # install.packages('PerformanceAnalytics', dependencies = T)

 # Package aktivieren
 library(PerformanceAnalytics)

 chart.Correlation(BL[,c(8,9,11:20)], histogram = T, pch = 19,
                   method = 'spearman')

 # Die Chart bildet alle wichtigen Informationen in einem ab.
 # Signifikante Korrelationen werden mit * bis *** gekennzeichnet.
 # Die Schriftgröße zeigt die Stärke der Zusammenhänge an.
 # Zudem werden die Histogramme zur jeweiligen Variable dargestellt.

 # Beispielhaft kann beschrieben werden, dass ein
 # höchstsignifikanter, positiver und starker statistischer
 # Zusammenhang (*** p = .001, r = 0.83) zwischen Kalium und Mangan
 # besteht.

 #####################################################################
 #### Hinweise

 # Damit sind nun die wichtigsten Analysen, Datenüberprüfungen und
 # Grafiken generiert.
 # Bei Fragen ziehe gerne das Internet zurate oder wende Dich im
 # RLab-Blog mit der Kommentarfunktion an andere User oder das
 # RLab-Team (https://rlab.blogs.uni-hamburg.de/)

 #####################################################################
 #### ENDE

Du hast Fragen zu diesem Code? Nutze gerne die Kommentar-Funktion!

Trennzeichen in .csv-Dateien - Ursache für fehlerhaften Datenimport
.csv-Dateien besitzen (leider) nicht immer identische Zeichen für die Trennung von Spalten und Dezimalstellen. Häufig werden Spalten durch Semikola getrennt, üblich sind aber auch Kommata oder andere Zeichen. Die Dezimalstellen werden meistens mit Punkt oder Komma abgetrennt, aber auch hier sind prinzipiell andere Zeichen möglich.

Bei der Anzeige in MS-Excel ist die Identifikation der Trennzeichen nicht eindeutig möglich. Wir empfehlen, eine .csv-Datei dazu mit einem Texteditor (z.B. Notepad oder Editor) zu öffnen. In der Anzeige dieser Programme sind die Zeichen eindeutig ablesbar und können per copy&paste in den R-Code übernommen werden: Beim Einlesen einer .csv-Datei in R werden die Trennzeichen in der Funktion read.table() mit den Argumenten sep = und dec = angegeben. Für die Beispieldateien sind in den Code-Beispielen die passenden Trennzeichen angegeben. Bei anderen Dateien ist es ggf. erforderlich, sie anzupassen.

Ausführliche Infos zu .csv-Dateien findest Du im Digitalen Skript ,,Import von Daten in R” auf der Seite csv-Dateien im Detail.

R-Code zur Analyse der bodenchemischen Parameter eines Leitprofils

Hier findest Du R-Code, mit dem Du im ,,Datenmaster Leitprofile" enthaltene Werte von Leitprofilen analysieren und Ergebnisse visualisieren kannst.

Zum Ausführen des Codes erforderliche Beispieldaten

bereinigte Daten der Leitprofile, wie sie im R-Code (siehe unten) verwendet werden (zum Download als .csv-Datei ggf. Rechtsklick und "Speichern unter...")

Falls Du die Vorbereitung selbst nachvollziehen möchtest, steht hier auch der Rohdatensatz zur Verfügung:

Rohdaten Bodenchemie der Leitprofile, die vor den R-Berechnungen mit z.B. MS-Excel bereinigt werden müssen (zum Download als .xlsx-Datei ggf. Rechtsklick und "Speichern unter...")

Code

Kurz-Version des Bodenchemie-Skripts, ohne Kommentare (zum Download als .R-Datei ggf. Rechtsklick und "Speichern unter...")

Auch den kommentierten Code gibt es als R-Skript zum Download:

Kommentierter Code des Bodenchemie-Skripts (zum Download als .R-Datei ggf. Rechtsklick und "Speichern unter...")

Alternativ kannst Du den kommentierten Code hier direkt aus dem Digitalen Skript übernehmen (Copy&Paste):


 # Course: RLab - Labor
 # Lesson: Leitprofile Bodenchemische Parameter
 # Author: RLab-Team
 # Co-Authors: Ronja Gottschalk (programming, design), Niels Schwab (design),
 #             Elke Fischer (design)
 # Organization: Universität Hamburg, project 'RLab - Skriptbasierte modulare
 #               Umweltstatistik' Funded within 'Lehrlabor' of
 #               Universitätskolleg 2.0 by BMBF (01PL17033).

 #####################################################################

 # Herzlich Willkommen im R-Skript zur Analyse des Leitprofils.
 # Bevor Du mit der Arbeit beginnst, schließe bitte noch einen
 # Vorschritt ab, indem Du auf File --> 'Reopen with Encoding' klickst.
 # Wähle in dem dann erscheinenden Dialogfenster bitte 'UTF-8' aus.
 # Damit wird das Skript erneut in der UTF-8 Kodierung geladen,
 # sodass Umlaute, etc. korrekt dargestellt werden. Um dieses Skript
 # und die darin enthaltenen Kommentare gut zu verstehen, ist ein
 # grundlegendes Verständnis für die Nutzung von R und RStudio
 # erforderlich. Solltest Du Dir diese aneignen oder sie auffrischen
 # wollen, ist die Bearbeitung der Grundlagen-Inhalte des RLab-Angebots
 # (https://rlab.blogs.uni-hamburg.de/grundlagen/) und ggf. weiterer
 # RLab-Kurse (https://rlab.blogs.uni-hamburg.de/rlab-kurse/) zu empfehlen.


 # #### Vorarbeit ----------------------------------------------------------

 # Workspace aufräumen
 rm(list = ls())

 # Kontrolle des Workspace
 ls()

 # Abfrage des aktuellen 'Working Directory'
 getwd()

 # neues Working Directory festlegen
 setwd('E:/RLabPortable/R-Labor/Bodenphysik/')
 # Beachte hier bitte, dass Du Deinen Dateipfad entsprechend Deiner
 # Ordnerstruktur anpassen musst.
 # Gute Hilfen zum Einlesen von Daten erhältst Du im RLab-Kurs
 # 'Daten einlesen'.

 # Erneute Abfrage des 'Working Directory' zur Kontrolle, ob die
 # Änderungen übernommen wurden
 getwd()

 # Datensatz einladen
 LP <- read.table('LP_Datenmaster_final.csv',
                  na.strings=c('',' ','NA'),
                  header = T,
                  sep = ';',
                  dec = ',',
                  stringsAsFactors = F)

 # Auch hier sollte der Dateipfad entsprechend Deines Verzeichnisses
 # angepasst werden, sofern Deine Datei in einem anderen Ordner liegt,
 # als Dein Working Directory eingestellt ist.
 # na.strings definierte leere Zellen (''), Blanks (' ') und NA-Werte
 # ('NA') für R bereits beim Einlesen als 'NA' übergeben, andernfalls
 # werden diese bloß als Zeichenkette erkannt.
 # stringsAsFactors = F verhindert, dass Variablen automatisch als
 # Faktoren übernommen werden.
 # Ggf. sep = und dec = anpassen, weitere Infos dazu im Hinweis-
 # Feld am Ende der hier verlinkten Seite:
 # https://rlab.blogs.uni-hamburg.de/dig-skripte/Vorbereitung_Datenmaster/index.html?s=Schritt%202:%20Doppelte%20Kopfzeile%20entfernen
 #
 # Eine Anleitung zum Aufbereiten Deines Datensatzes für R findest Du
 # im RLab in den digitalen Skripten.
 # Hier ist allerdings darüber hinausgehend zu beachten, dass es sinnvoll ist,
 # bereits im Vorfeld in einem geeigneten Tabellenkalkulationsprogramm ein
 # Subset zu erstellen, das nur die Spalten enthält, die Du tatsächlich in
 # deiner Analyse in R benötigst.
 # Im hiesigen Skript arbeitest Du bereits mit einer gekürzten und optimierten
 # Variante des Datenmasters. Dennoch wird der Datensatz auch hier in R weiter
 # gesubsettet werden. Eine Vorabreduktion ist jedoch sinnvoll um den Überblick
 # zu behalten. Das ist gerade bei einem Datenmaster wie dem zu den Bodenparametern
 # wichtig, da dieser sehr umfangreich ist.

 # Überblick über den Datensatz  -------------------------------------------

 # Datensatz in R als Tabelle ansehen
 View(LP)

 # Struktur des Datensatzes ansehen
 str(LP)
 # die Mehrheit der Variablen wird als 'chr', also character
 # angegeben.
 # Gewollt sind hier für die Spalten mit Zahlenwerten jedoch
 # (wie im weiteren Verlauf noch deutlich wird) numerische Variablen.
 # Dies wird im weiteren Verlauf entsprechend zugewiesen.

 # Nachfolgend werden zunächst Subsets generiert. Das dient dazu den hiesigen
 # Ausgangsdatenmaster weiter zu reduzieren, damit er übersichtlicher wird.

 # Subset des Datensatzes - Bodenchemie ------------------------------------

 # Für die Analyse und Darstellung der bodenchemikalischen Parameter werden
 # allgemeine Daten wie beispielsweise der Standort benötigt und die Spalten zu den Laboranalysen
 # der einelnen Elemente.
 # Diese werden im nachfolgenden Subset aus dem Dataframe extrahiert und als neues Objekt
 # 'LP_chem' abgespeichert.
 LP_chem <- LP[,c(1:9,40:57)]

 # Überprüfung des Datensatzes ---------------------------------------------

 # Für den Datensatz werden nun zunächstg Sichtkontrollen und Kontrollen der
 # Datenstruktur sowie die deskriptiven Statistiken erstellt.

 View(LP_chem)

 str(LP_chem)

 summary(LP_chem)

 # Wie Du sehen kannst, ist der Datensatz nun dennoch sehr groß und unhandlich.
 # Es empfielt sich daher sich zu überlegen, welche Parameter tatsächlich im
 # Zuge der Analysen und Visualisierungen benötigt werden.
 # Nachfolgend soll der Datenmaster daher weiter gesubsettet werden,
 # sodass lediglich die interessierenden Parameter im Datensatz enthalten
 # bleiben.

 # Für die bodenchemischen Parameter sind vor allem einzelne chemische Parameter von Interesse.
 # Je nachdem, welche Parameter Du verarbeiten möchtest, werden
 # natürlich die jeweiligen Spalten benötigt.
 # Damit der Datensatz in diesem Skript nachvollziehbar bleibt, wird nachfolgend mit den
 # Profilkenngrößen (Spalten 1:6) und den Parametern Kationenaustauschkapazität ('KAK', Spalte 26),
 # dem pH-Wert ('pH', Spalte 11) und der Basensättigung ('BS', Spalte 27) gearbeitet.
 # Zunächst wird der Datensatz auf die zuvor genannten Spalten gesubsettet.
 LP_chem1 <- LP_chem[,c(1:6,11,26,27)]

 # Kurze Strukturprüfung des Subsets
 str(LP_chem1)

 # Wie Du in der Ausgabe des str() sehen kannst, sind die drei letzten Parameter noch als 'chr'
 # formatiert. Diese sollen nun in einem Schritt zu numerischen Variablen werden.
 # Mittels sapply wird die as.numeric Funktion gleichzeitig auf alle angegebenen
 # Spalten angewendet.
 LP_chem1[,c(7:9)]=as.data.frame(sapply(LP_chem1[,c(7:9)], as.numeric))

 # Im Verlauf Deiner Arbeit, kann es hilfreich sein die Spaltennamen intern in R umzubenennen.
 # Dies geschieht mit dem Argument colnames(). Hier werden exemplarisch die Namen angepasst.
 colnames(LP_chem1) <- c('Labornr','Probe','Profilnr','Flaeche','HorizontOG',
                         'HorizontUG','pH','KAKeff','BS')

 # Nachfolgend sollen einige Plots erstellt werden. Dabei ist
 # zunächst das Arbeiten mit maximal drei Variablen, die gleichzeitig in einem Plot
 # dargestellt werden, vorgesehen.Hierzu müssen die in den Varibalen enthaltenen Werte über eine
 # 'Iditifier'-Spalte eindeutig einer Gruppe zuzuordnen sein.
 # Würde man nun aber beispielsweise die Horizontuntergrenze nehmen, würde diese die Werte der
 # Proben a, b und c enthalten. Daher soll zunächst eine Spalte erzeugt werden, die den
 # Probenbuchstaben enthält, sodass später ein Subset je Probenbuchstabe erstellt werden kann.
 str(LP_chem1)
 LP_chem1$Flaeche <- as.factor(LP_chem1$Flaeche)

 # Die Spalte 'Probe' wird mit dem Objekt x gleichgesetzt.
 x <- LP_chem1$Probe

 # Der Buchstabe wird aus den Zeichenketten in der Spalte 'Probe' extrahiert
 Buchstaben_extrahieren <- function(x, n){
   substr(x, nchar(x)-n+1, nchar(x))
 }

 x2 <- Buchstaben_extrahieren(x, 1)

 # In dem Dataframe wird der Dataframe LP_chem1 mit der Spalte des extrahierten
 # Buchstabens verbunden.
 LP_chem_id <- cbind(LP_chem1, x2)

 # Der spaltenname wird für die neue Spalte entsprechend angepasst
 colnames(LP_chem_id)[10] <- 'Probenbuchstabe'

 # Es wird erneut ein Subset erstellt. Hier sollen nur noch die Proben 'a' in den Datensatz
 # eingehen.
 # Nachfoglend wird dann mit dem Subset weitergearbeitet.
 LP_chem_a <- LP_chem_id[LP_chem_id[,'Probenbuchstabe']=='a',]


 # Visualisierung ----------------------------------------------------------

 # Nun sollen Plots erstellt werden.
 # Dazu wird das Package 'ggplot2' benötigt.
 # install.packages('ggplot2')
 library(ggplot2)
 summary(LP_chem_id)
 # Nun soll ein Plot der Verteilung der chemischen Parameter in die Tiefe erstellt werden.
 # Im weitesten Sinne ist das ein Liniendiagramm, dessen X-Achse nach oben versetzt wird.
 # Da der pH-Wert üblicherweise zwischen 2-10 liegt, wird der Achsenbereich zunächst auf
 # 0-12 skaliert.
 ggplot(LP_chem_id, aes(x = pH, y = HorizontOG)) + # x-Achse = pH-Wert; y-Achse = Horizontobergrenze
     geom_path(aes(color= Flaeche)) + # farbliche Unterscheidung Der Gruppen anhand der 'Fläche'
   ylab('Profiltiefe [cm]') + xlab('pH-Wert') + # Achsenbeschriftungen
   theme(axis.text.x=element_text(angle = 45, vjust = 0, hjust = 0)) + # rotierte x-Achsenbeschriftung
   ggtitle('pH-Werte nach Standort ') + # Titel des Plots
   theme(plot.title = element_text(hjust = 0.5)) +
   geom_path(aes(color= Flaeche)) +
   scale_y_reverse() + # y-Achse invertiert darstellen
   scale_x_continuous(position = 'top', breaks = c(0,2,4,6,8,10,12), # Skalierung der x-Achse festlegen
   labels= c(0,2,4,6,8,10,12)) +
   coord_cartesian(xlim=c(0,12), ylim=c(0,110)) # Grenzen der Achsen festlegen

 # Wertebereich im Detail
 ggplot(LP_chem_id, aes(x = pH, y = HorizontOG)) +
   geom_path(aes(color= Flaeche)) +
   ylab('Profiltiefe [cm]') + xlab('pH-Wert') +
   theme(axis.text.x=element_text(angle = 45, vjust = 0, hjust = 0)) +
   ggtitle('pH-Werte nach Standort ') +
   theme(plot.title = element_text(hjust = 0.5)) +
   geom_path(aes(color= Flaeche)) +
   scale_y_reverse() +
   scale_x_continuous(position = 'top', breaks = c(2,4,6),
                      labels= c(2,4,6)) +
   coord_cartesian(xlim=c(2,6), ylim=c(0,110))


 # Damit konntest Du die Bodenchemischen Parameter über die Tiefe visualisieren.
 # Noch ein Hinweis zum Erstellen von Grafiken in R: Wenn Du die
 # Abbildungen aus R exportieren möchtest, findest Du im RLab eine
 # gute Anleitung dazu.
 # https://rlab.blogs.uni-hamburg.de/dig-skripte/Tipps_und_Tricks/index.html

 # Bei Fragen ziehe gerne das Internet zurate oder wende Dich im
 # RLab-Blog mit der Kommentarfunktion an andere User oder das
 # RLab-Team (https://rlab.blogs.uni-hamburg.de/)

 #####################################################################
 #### ENDE

Du hast Fragen zu diesem Code? Nutze gerne die Kommentar-Funktion!

Trennzeichen in .csv-Dateien - Ursache für fehlerhaften Datenimport
.csv-Dateien besitzen (leider) nicht immer identische Zeichen für die Trennung von Spalten und Dezimalstellen. Häufig werden Spalten durch Semikola getrennt, üblich sind aber auch Kommata oder andere Zeichen. Die Dezimalstellen werden meistens mit Punkt oder Komma abgetrennt, aber auch hier sind prinzipiell andere Zeichen möglich.

Bei der Anzeige in MS-Excel ist die Identifikation der Trennzeichen nicht eindeutig möglich. Wir empfehlen, eine .csv-Datei dazu mit einem Texteditor (z.B. Notepad oder Editor) zu öffnen. In der Anzeige dieser Programme sind die Zeichen eindeutig ablesbar und können per copy&paste in den R-Code übernommen werden: Beim Einlesen einer .csv-Datei in R werden die Trennzeichen in der Funktion read.table() mit den Argumenten sep = und dec = angegeben. Für die Beispieldateien sind in den Code-Beispielen die passenden Trennzeichen angegeben. Bei anderen Dateien ist es ggf. erforderlich, sie anzupassen.

Ausführliche Infos zu .csv-Dateien findest Du im Digitalen Skript ,,Import von Daten in R” auf der Seite csv-Dateien im Detail.

R-Code zur Analyse der bodenphysikalischen Parameter eines Leitprofils

Hier findest Du R-Code, mit dem Du im ,,Datenmaster Leitprofile" enthaltene Werte von Leitprofilen analysieren und Ergebnisse visualisieren kannst.

Zum Ausführen des Codes erforderliche Beispieldaten

bereinigte Daten der Leitprofile, wie sie im R-Code (siehe unten) verwendet werden (zum Download als .csv-Datei ggf. Rechtsklick und "Speichern unter...")

Falls Du die Vorbereitung selbst nachvollziehen möchtest, steht hier auch der Rohdatensatz zur Verfügung:

Rohdaten der Bodenphysik Leitprofile, die vor den R-Berechnungen mit z.B. MS-Excel bereinigt werden müssen (zum Download als .xlsx-Datei ggf. Rechtsklick und "Speichern unter...")

Code

Kurz-Version des Bodenphysik-Skripts, ohne Kommentare (zum Download als .R-Datei ggf. Rechtsklick und "Speichern unter...")

Auch den kommentierten Code gibt es als R-Skript zum Download:

Kommentierter Code des Bodenphysik-Skripts (zum Download als .R-Datei ggf. Rechtsklick und "Speichern unter...")

Alternativ kannst Du den kommentierten Code hier direkt aus dem Digitalen Skript übernehmen (Copy&Paste):


 # Course: RLab - Labor
 # Lesson: Leitprofile Bodenphysikalische Parameter
 # Author: RLab-Team
 # Co-Authors: Ronja Gottschalk (programming, design), Niels Schwab (design),
 #             Elke Fischer (design)
 # Organization: Universität Hamburg, project 'RLab - Skriptbasierte modulare
 #               Umweltstatistik' Funded within 'Lehrlabor' of
 #               Universitätskolleg 2.0 by BMBF (01PL17033).

 #####################################################################

 # Herzlich Willkommen im R-Skript zur Analyse des Leitprofils.
 # Bevor Du mit der Arbeit beginnst, schließe bitte noch einen
 # Vorschritt ab, indem Du auf File --> 'Reopen with Encoding' klickst.
 # Wähle in dem dann erscheinenden Dialogfenster bitte 'UTF-8' aus.
 # Damit wird das Skript erneut in der UTF-8 Kodierung geladen,
 # sodass Umlaute, etc. korrekt dargestellt werden. Um dieses Skript
 # und die darin enthaltenen Kommentare gut zu verstehen, ist ein
 # grundlegendes Verständnis für die Nutzung von R und RStudio
 # erforderlich. Solltest Du Dir diese aneignen oder sie auffrischen
 # wollen, ist die Bearbeitung der Grundlagen-Inhalte des RLab-Angebots
 # (https://rlab.blogs.uni-hamburg.de/grundlagen/) und ggf. weiterer
 # RLab-Kurse (https://rlab.blogs.uni-hamburg.de/rlab-kurse/) zu empfehlen.


 # #### Vorarbeit ----------------------------------------------------------

 # Workspace aufräumen
 rm(list = ls())

 # Kontrolle des Workspace
 ls()

 # Abfrage des aktuellen 'Working Directory'
 getwd()

 # neues Working Directory festlegen
 setwd('E:/RLabPortable/R-Labor/Bodenphysik/')
 # Beachte hier bitte, dass Du Deinen Dateipfad entsprechend Deiner
 # Ordnerstruktur anpassen musst.
 # Gute Hilfen zum Einlesen von Daten erhältst Du im RLab-Kurs
 # 'Daten einlesen'.

 # Erneute Abfrage des 'Working Directory' zur Kontrolle, ob die
 # Änderungen übernommen wurden
 getwd()

 # Datensatz einladen
 LP <- read.table('LP_Datenmaster_final.csv',
                  na.strings=c('',' ','NA'),
                  header = T,
                  sep = ';',
                  dec = ',',
                  stringsAsFactors = F)

 # Auch hier sollte der Dateipfad entsprechend Deines Verzeichnisses
 # angepasst werden, sofern Deine Datei in einem anderen Ordner liegt,
 # als Dein Working Directory eingestellt ist.
 # na.strings definierte leere Zellen (''), Blanks (' ') und NA-Werte
 # ('NA') für R bereits beim Einlesen als 'NA' übergeben, andernfalls
 # werden diese bloß als Zeichenkette erkannt.
 # stringsAsFactors = F verhindert, dass Variablen automatisch als
 # Faktoren übernommen werden.
 # Ggf. sep = und dec = anpassen, weitere Infos dazu im Hinweis-
 # Feld am Ende der hier verlinkten Seite:
 # https://rlab.blogs.uni-hamburg.de/dig-skripte/Vorbereitung_Datenmaster/index.html?s=Schritt%202:%20Doppelte%20Kopfzeile%20entfernen
 #
 # Eine Anleitung zum Aufbereiten Deines Datensatzes für R findest Du
 # im RLab in den digitalen Skripten.
 # Hier ist allerdings darüber hinausgehend zu beachten, dass es sinnvoll ist,
 # bereits im Vorfeld in einem geeigneten Tabellenkalkulationsprogramm ein
 # Subset zu erstellen, das nur die Spalten enthält, die Du tatsächlich in
 # deiner Analyse in R benötigst.
 # Im hiesigen Skript, arbeitest Du bereits mit einer gekürzten und optimierten
 # Variante des Datenmasters. Dennoch wird der Datensatz auch hier in R weiter
 # gesubsettet werden. Eine Vorabreduktion ist jedoch sinnvoll um den Überblick
 # zu behalten. Das ist gerade bei einem Datenmaster wie dem zu den Bodenparametern
 # wichtig, da dieser sehr umfangreich ist.

 # Überblick über den Datensatz  -------------------------------------------

 # Datensatz in R als Tabelle ansehen
 View(LP)

 # Struktur des Datensatzes ansehen
 str(LP)
 # die Mehrheit der Variablen wird als 'chr', also character
 # angegeben.
 # Gewollt sind hier für die Spalten mit Zahlenwerten jedoch
 # (wie im weiteren Verlauf noch deutlich wird) numerische Variablen.
 # Dies wird im weiteren Verlauf entsprechend zugewiesen.

 # Nachfolgend werden zunächst Subsets generiert. Das dient dazu den hiesigen
 # Ausgangsdatenmaster weiter zu reduzieren, damit er übersichtlicher wird.

 # Subset des Datensatzes - Bodenphysik ------------------------------------

 # Für die Analyse und Darstellung der bodenphysikalischen Parameter werden
 # allgemeine Daten wie bespw. der Standort benötigt und die Spalten zur Konrgrößenverteilung.
 # Diese werden im nachfolgenden Subset aus dem Dataframe extrahiert und als neues Objekt
 # 'Boden_phys' abgespeichert.
 LP_phys <- LP[,c(1:40)]

 # Überprüfung des Datensatzes ---------------------------------------------

 # Für den Datensatz werden nun zunächstg Sichtkontrollen und Kontrollen der
 # Datenstruktur sowie die deskriptiven Statistiken erstellt.

 View(LP_phys)

 str(LP_phys)

 summary(LP_phys)

 # Wie Du sehen kannst, ist der Datensatz nun dennoch sehr groß und unhandlich.
 # Es empfielt sich daher sich zu überlegen, welche Parameter tatsächlich im
 # Zuge der Analysen und Visualisierungen benötigt werden.
 # Nachfolgend soll der Datenmaster daher weiter gesubsettet werden,
 # sodass lediglich die interessierenden Parameter im Datensatz enthalten
 # bleiben.
 # Zunächst soll die Analytik und die Visualisierung für den physischen Datensatz
 # erfolgen und anschließend das selbe Prozedere für den bodenchemischen Teil.
 # Je nachdem für welchen Part des Datensatzdes Du Dich inhaltlich interessierst,
 # folge daher bitte weiter dem Skript oder springe weiter nach unten vor.

 # Für den bodenphysikalischen Part des Leitprofils ist vor allem die
 # Korngrößenzusammensetzung von Interesse. Daher wird nachfolgend ein Untersubset
 # erstellt, welches nur die Kerndaten und den prozentualen Anteil der
 # Bodenart am Gesamtboden betrachtet. Hiermit wird nachfolgend weitergearbeitet.

 LP_phys1 <- LP_phys[,c(1:6, 18,20,22,26,29,32,34)]

 # erneute Strukturprüfung
 str(LP_phys1)

 # Wie Du in der Ausgabe des str() sehen kannst, sind die meisten Parameter noch als 'chr'
 # formatiert. Diese sollen nun in einem Schritt zu numerischen Variablen werden.
 # Mittels sapply wird die as.numeric Funktion gleichzeitig auf alle angegebenen
 # Spalten angewendet.
 LP_phys1[,c(7:13)]=as.data.frame(sapply(LP_phys1[,c(7:13)], as.numeric))

 # Im Verlauf Deiner Arbeit, kann es hilfreich sein die Spaltennamen intern in R umzubenennen.
 # Dies geschieht mit dem Argument colnames(). Hier werden exemplarisch die Namen angepasst.
 colnames(LP_phys1) <- c('Labornr','Probe','Profilnr','Flaeche','HorizontOG',
                         'HorizontUG','fS','mS','gS','gU','mU','fU','T')

 # Momentan sind die einzelnen Korngrößenfraktionen noch in Spalten formatiert.
 # Für die gleichzeitige Darstellung in einem ggplot ist die Umstrukurierung dieser Daten in
 # eine gemeinsame Spalte erfoderlich. Das geschieht mit der Funktion gather aus dem Package (tidyr)
 #install.packages('tidyr')
 library(tidyr)
 # Umstrukturieren der Daten von breit auf lang
 LP_phys_reshaped <- gather(LP_phys1, Korngroessenfraktion, Anteil_prozentual,
                               'fS','mS','gS','gU','mU','fU','T', factor_key=TRUE)

 # Strukturprüfung der Daten zur KOntrolle der Funktion
 str(LP_phys_reshaped)

 # Darüber hinaus fällt Dir vielleicht auf, dass die Spalte 'Fläche' als 'chr' formatiert ist.
 # Der Datensatz soll für spätere Arbeiten anhand der verschiedenen Flächen in
 # Einzeldatensätze aufgeteilt werden.
 # für die Funktion 'split.data.frame' ist für die variable, anhand derer aufgeteilt werden soll,
 # ein Faktorformat erfoderlich.
 LP_phys_reshaped$Flaeche <- as.factor(LP_phys_reshaped$Flaeche)

 # jetzt kann der Datensatz anhand der Flächen in Einzeldatensätze aufgeteilt werden.
 Korngroesse <- split.data.frame(LP_phys_reshaped, LP_phys_reshaped$Flaeche)
 Buche <- Korngroesse$Eichen_Buchen_Flaeche
 Douglasie <- Korngroesse$Douglasien_Flaeche
 Loess <- Korngroesse$Loess_Profil
 Podsol <- Korngroesse$Podsol_Profil

 # Exemplarische Strukturprüfung für den Datensatz 'Buche'
 str(Buche)

 # Nachfolgend sollen einige Plots erstellt werden. Dabei ist
 # zunächst das Arbeiten mit maximal drei Variablen, die gleichzeitig in einem Plot
 # dargestellt werden, vorgesehen. Hierzu müssen die in den Varibalen enthaltenen Werte über eine
 # 'Iditifier'-Spalte eindeutig einer Gruppe zuzuordnen sein.
 # Würde man nun aber beispielsweise die Horizontuntergrenze nehmen, würde diese die Werte der
 # Proben a, b und c enthalten. Daher soll zunächst eine Spalte erzeugt werden, die den
 # Probenbuchstaben enthält, sodass später ein Subset je Probenbuchstabe erstellt werden kann.

 # Die Spalte 'Probe' wird mit dem Objekt x gleichgesetzt.
 x <- Buche$Probe

 # Funktion zur Extraktion schreiben
 Buchstaben_extrahieren <- function(x, n){
   substr(x, nchar(x)-n+1, nchar(x))
 }

 # Der Buchstabe wird aus den Zeichenketten in der Spalte 'Probe' extrahiert,
 # indem die Funktion angewendet wird.
 x2 <- Buchstaben_extrahieren(x, 1)

 # In dem Dataframe Buche2 wird der Dataframe LP_Buche mit der Spalte des extrahierten
 # Buchstabens verbunden.
 Buche2 <- cbind(Buche, x2)

 # Anpassung des Spaltennamens
 colnames(Buche2)[9] <- 'Probenbuchstabe'

 # Es wird erneut ein Subset erstellt. Hier sollen nur noch die Proben 'a' in den
 # Datenasatz eingehen.
 # Nachfoglend wird dann mit dem Subset weitergearbeitet.
 Buche_a <- Buche2[Buche2[,'Probenbuchstabe']=='a',]


 # Visualisierung ----------------------------------------------------------

 # Nun sollen Plots erstellt werden.
 # Dazu wird das Package 'ggplot2' benötigt.
 # install.packages('ggplot2')
 library(ggplot2)
 # Nun soll ein Plot der Verteilung der verschieden Korngrößen in die Tiefe erstellt werden.
 # Im weitesten Sinne ist das ein Liniendiagramm, dessen X-Achse nach oben versetzt wird.
 ggplot(Buche_a, aes(x = Buche_a$Anteil_prozentual, y = HorizontOG)) + # x-Achse = pH-Wert; y-Achse = Horizontobergrenze
   geom_path(aes(color= Korngroessenfraktion)) + # farbliche Unterscheidung der Gruppen anhand der 'Korngroessenfraktion'
   ylab('Profiltiefe [cm]') + xlab('Anteil [%]') + # Achsenbeschriftungen
   theme(axis.text.x=element_text(angle = 45, vjust = 0, hjust = 0)) +
   ggtitle('Korngroessenzusammensetzung des Profils der Buchenflaeche') + # Titel des Plots
   theme(plot.title = element_text(hjust = 0.5)) +
   geom_path(aes(color= Korngroessenfraktion)) +
   scale_y_reverse() + # y-Achse invertiert darstellen
   scale_x_continuous(position = 'top', breaks = c(-10,0,25,50,75,100), # Skalierung der x-Achse festlegen
   labels= c(-10,0,25,50,75,100)) +
   coord_cartesian(xlim=c(-10,100), ylim=c(0,95)) # Wertegrenzen der Achsen festlegen

 # Wie du sehen kannst, bewegt sich einer der Graphen im negativen Bereich.
 # Dies ist durch die Datengrundlage zu erklären, bei der rechnerisch über die Subtraktion der Anteile
 # der einzelnen Fraktionen ermittelt werden.
 # Da sich rechnerisch auch negative Werte ergeben können, gibt es in diesem Datensatz negative Werte,
 # obwohl dies inhaltlich nicht plausibel ist.

 # Damit konntest Du die Korngrößenfraktion über die Tiefe visualisieren.
 # Alternativ könnte man hier auch andere psychikalische Kenngrößen, wie die Lagerungsdichte,
 # den Anteil organischer Substanz oder ähnliches darstellen.
 # Je nachdem muss der Plot entsprechend angepasst werden.
 # Wenn die Daten bereits in einer gemeinsamen Spalte formatiert sind und es zusätzlich eine
 # Unterscheidungsvariable gibt, können die Daten direkt geplottet werden.
 # Alternativ müssen die Daten wie oben umstrukturiert werden.

 # Exemplarisch sollen daher nun noch für die übrigen Flächen der
 # Gesamtanteil des Sandes in die Tiefe dargestellt werden.
 # Dazu wird nun der ursprüngliche Datensatz LP_phys verwendet.
 # Auf der X-Achse soll die Spalte 10 geplottet werden.

 # Die Spalte 'Fläche' ist ein Beispiel, für das die Daten nicht neu strukturiert werden müssen.
 str(LP_phys)
 LP_phys$Sand.gesamt_.g. <- as.numeric(LP_phys$Sand.gesamt_.g.)

 # Damit der Schritt mit dem Extrahieren nach Probenbuchstabe entfällt,
 # wird für dieses Beispiel die Fläche 'Eiche-Buche' ausgeschlossen.
 LP_phys2 <- LP_phys[19:32,]

 ggplot(LP_phys2, aes(x = LP_phys2$Sand.gesamt_.g., y = LP_phys2$Horizontobergrenze.cm.)) + # x-Achse = pH-Wert; y-Achse = Horizontobergrenze
   geom_path(aes(color= Flaeche)) + # farbliche Unterscheidung der Gruppen anhand der 'Korngroessenfraktion'
   ylab('Profiltiefe [cm]') + xlab('Anteil [%]') + # Achsenbeschriftungen
   theme(axis.text.x=element_text(angle = 45, vjust = 0, hjust = 0)) +
   ggtitle('Korngroessenzusammensetzung des Profils der Buchenflaeche') + # Titel des Plots
   theme(plot.title = element_text(hjust = 0.5)) +
   geom_path(aes(color= Flaeche)) +
   scale_y_reverse() + # y-Achse invertiert darstellen
   scale_x_continuous(position = 'top', breaks = c(-10,0,25,50,75,100), # Skalierung der x-Achse festlegen
                      labels= c(-10,0,25,50,75,100)) +
   coord_cartesian(xlim=c(-10,100), ylim=c(0,95)) # Wertegrenzen der Achsen festlegen

 # Damit hast du zwei unterschiedliche Datenstrukturen aufbereitet und visualisiert.

 # Noch ein Hinweis zum Erstellen von Grafiken in R: Wenn Du die
 # Abbildungen aus R exportieren möchtest, findest Du im RLab eine
 # gute Anleitung dazu.
 # https://rlab.blogs.uni-hamburg.de/dig-skripte/Tipps_und_Tricks/index.html

 # Bei Fragen ziehe gerne das Internet zurate oder wende Dich im
 # RLab-Blog mit der Kommentarfunktion an andere User oder das
 # RLab-Team (https://rlab.blogs.uni-hamburg.de/)

 #####################################################################
 #### ENDE

Du hast Fragen zu diesem Code? Nutze gerne die Kommentar-Funktion!

Trennzeichen in .csv-Dateien - Ursache für fehlerhaften Datenimport
.csv-Dateien besitzen (leider) nicht immer identische Zeichen für die Trennung von Spalten und Dezimalstellen. Häufig werden Spalten durch Semikola getrennt, üblich sind aber auch Kommata oder andere Zeichen. Die Dezimalstellen werden meistens mit Punkt oder Komma abgetrennt, aber auch hier sind prinzipiell andere Zeichen möglich.

Bei der Anzeige in MS-Excel ist die Identifikation der Trennzeichen nicht eindeutig möglich. Wir empfehlen, eine .csv-Datei dazu mit einem Texteditor (z.B. Notepad oder Editor) zu öffnen. In der Anzeige dieser Programme sind die Zeichen eindeutig ablesbar und können per copy&paste in den R-Code übernommen werden: Beim Einlesen einer .csv-Datei in R werden die Trennzeichen in der Funktion read.table() mit den Argumenten sep = und dec = angegeben. Für die Beispieldateien sind in den Code-Beispielen die passenden Trennzeichen angegeben. Bei anderen Dateien ist es ggf. erforderlich, sie anzupassen.

Ausführliche Infos zu .csv-Dateien findest Du im Digitalen Skript ,,Import von Daten in R” auf der Seite csv-Dateien im Detail.

R-Code zur Darstellung der Korngrößenzusammensetzung eines Leitprofils

Hier findest Du R-Code, mit dem Du im ,,Datenmaster Bodenphysik" enthaltene Korngrößenwerte von Leitprofilen analysieren und Ergebnisse visualisieren kannst.

Zum Ausführen des Codes erforderliche Beispieldaten

bereinigte Daten der Leitprofile, wie sie im R-Code (siehe unten) verwendet werden (zum Download als .csv-Datei ggf. Rechtsklick und "Speichern unter...")

Falls Du die Vorbereitung selbst nachvollziehen möchtest, steht hier auch der Rohdatensatz zur Verfügung:

Rohdaten Bodenphysik der Leitprofile, die vor den R-Berechnungen mit z.B. MS-Excel bereinigt werden müssen (zum Download als .xlsx-Datei ggf. Rechtsklick und "Speichern unter...")

Code

Kurz-Version des Korngrößen-Skripts, ohne Kommentare (zum Download als .R-Datei ggf. Rechtsklick und "Speichern unter...")

Auch den kommentierten Code gibt es als R-Skript zum Download:

Kommentierter Code des Korngrößen-Skripts (zum Download als .R-Datei ggf. Rechtsklick und "Speichern unter...")

Alternativ kannst Du den kommentierten Code hier direkt aus dem Digitalen Skript übernehmen (Copy&Paste):


 # Course: RLab - Labor
 # Lesson: Leitprofile Bodenphysikalische Parameter Korngroessendreieck
 # Author: RLab-Team
 # Co-Authors: Ronja Gottschalk (programming, design), Niels Schwab (design),
 #             Elke Fischer (design)
 # Organization: Universität Hamburg, project 'RLab - Skriptbasierte modulare
 #               Umweltstatistik' Funded within 'Lehrlabor' of
 #               Universitätskolleg 2.0 by BMBF (01PL17033).

 #####################################################################

 # Herzlich Willkommen im R-Skript zur Analyse des Leitprofils.
 # Bevor Du mit der Arbeit beginnst, schließe bitte noch einen
 # Vorschritt ab, indem Du auf File --> 'Reopen with Encoding' klickst.
 # Wähle in dem dann erscheinenden Dialogfenster bitte 'UTF-8' aus.
 # Damit wird das Skript erneut in der UTF-8 Kodierung geladen,
 # sodass Umlaute, etc. korrekt dargestellt werden. Um dieses Skript
 # und die darin enthaltenen Kommentare gut zu verstehen, ist ein
 # grundlegendes Verständnis für die Nutzung von R und RStudio
 # erforderlich. Solltest Du Dir diese aneignen oder sie auffrischen
 # wollen, ist die Bearbeitung der Grundlagen-Inhalte des RLab-Angebots
 # (https://rlab.blogs.uni-hamburg.de/grundlagen/) und ggf. weiterer
 # RLab-Kurse (https://rlab.blogs.uni-hamburg.de/rlab-kurse/) zu empfehlen.


 # #### Vorarbeit ----------------------------------------------------------

 # Workspace aufräumen
 rm(list = ls())

 # Kontrolle des Workspace
 ls()

 # Abfrage des aktuellen 'Working Directory'
 getwd()

 # neues Working Directory festlegen
 setwd('E:/RLabPortable/R-Labor/Bodenphysik/')
 # Beachte hier bitte, dass Du Deinen Dateipfad entsprechend Deiner
 # Ordnerstruktur anpassen musst.
 # Gute Hilfen zum Einlesen von Daten erhältst Du im RLab-Kurs
 # 'Daten einlesen'.

 # Erneute Abfrage des 'Working Directory' zur Kontrolle, ob die
 # Änderungen übernommen wurden
 getwd()

 # Datensatz einladen
 LP <- read.table('LP_Datenmaster_final.csv',
                  na.strings=c('',' ','NA'),
                  header = T,
                  sep = ';',
                  dec = ',',
                  stringsAsFactors = F)

 # Auch hier sollte der Dateipfad entsprechend Deines Verzeichnisses
 # angepasst werden, sofern Deine Datei in einem anderen Ordner liegt,
 # als Dein Working Directory eingestellt ist.
 # na.strings definierte leere Zellen (''), Blanks (' ') und NA-Werte
 # ('NA') für R bereits beim Einlesen als 'NA' übergeben, andernfalls
 # werden diese bloß als Zeichenkette erkannt.
 # stringsAsFactors = F verhindert, dass Variablen automatisch als
 # Faktoren übernommen werden.
 # Ggf. sep = und dec = anpassen, weitere Infos dazu im Hinweis-
 # Feld am Ende der hier verlinkten Seite:
 # https://rlab.blogs.uni-hamburg.de/dig-skripte/Vorbereitung_Datenmaster/index.html?s=Schritt%202:%20Doppelte%20Kopfzeile%20entfernen
 #
 # Eine Anleitung zum Aufbereiten Deines Datensatzes für R findest Du
 # im RLab in den digitalen Skripten.
 # Hier ist allerdings darüber hinausgehend zu beachten, dass es sinnvoll ist,
 # bereits im Vorfeld in einem geeigneten Tabellenkalkulationsprogramm ein
 # Subset zu erstellen, das nur die Spalten enthält, die Du tatsächlich in
 # deiner Analyse in R benötigst.
 # Im hiesigen Skript, arbeitest Du bereits mit einer gekürzten und optimierten
 # Variante des Datenmasters. Dennoch wird der Datensatz auch hier in R weiter
 # gesubsettet werden. Eine Vorabreduktion ist jedoch sinnvoll um den Überblick
 # zu behalten. Das ist gerade bei einem Datenmaster wie dem zu den Bodenparametern
 # wichtig, da dieser sehr umfangreich ist.

 # Überblick über den Datensatz  -------------------------------------------

 # Datensatz in R als Tabelle ansehen
 View(LP)

 # Struktur des Datensatzes ansehen
 str(LP)
 # die Mehrheit der Variablen wird als 'chr', also character
 # angegeben.
 # Gewollt sind hier für die Spalten mit Zahlenwerten jedoch
 # (wie im weiteren Verlauf noch deutlich wird) numerische Variablen.
 # Dies wird im weiteren Verlauf entsprechend zugewiesen.

 # Nachfolgend werden zunächst Subsets generiert. Das dient dazu den hiesigen
 # Ausgangsdatenmaster weiter zu reduzieren, damit er übersichtlicher wird.

 # Subset des Datensatzes - Bodenphysik ------------------------------------

 # Für die Analyse und Darstellung der bodenphysikalischen Parameter werden
 # allgemeine Daten wie bespw. der Standort benötigt und die Spalten zur Konrgrößenverteilung.
 # Diese werden im nachfolgenden Subset aus dem Dataframe extrahiert und als neues Objekt
 # 'Boden_phys' abgespeichert.
 LP_phys <- LP[,c(1:40)]

 # Überprüfung des Datensatzes ---------------------------------------------

 # Für den Datensatz werden nun zunächstg Sichtkontrollen und Kontrollen der
 # Datenstruktur sowie die deskriptiven Statistiken erstellt.

 View(LP_phys)

 str(LP_phys)

 summary(LP_phys)

 # Wie Du sehen kannst, ist der Datensatz nun dennoch sehr groß und unhandlich.
 # Es empfielt sich daher sich zu überlegen, welche Parameter tatsächlich im
 # Zuge der Analysen und Visualisierungen benötigt werden.
 # Nachfolgend soll der Datenmaster daher weiter gesubsettet werden,
 # sodass lediglich die interessierenden Parameter im Datensatz enthalten
 # bleiben.
 # Zunächst soll die Analytik und die Visualisierung für den physischen Datensatz
 # erfolgen und anschließend das selbe Prozedere für den bodenchemischen Teil.
 # Je nachdem für welchen Part des Datensatzdes Du Dich inhaltlich interessierst,
 # folge daher bitte weiter dem Skript oder springe weiter nach unten vor.

 # Für den bodenphysikalischen Part des Leitprofils ist vor allem die
 # Korngrößenzusammensetzung von Interesse. Daher wird nachfolgend ein Untersubset
 # erstellt, welches nur die Kerndaten und den prozentualen Anteil der
 # Bodenart am Gesamtboden betrachtet. Hiermit wird nachfolgend weitergearbeitet.

 LP_phys1 <- LP_phys[,c(1:6, 18,20,22,26,29,32,34)]

 # Bevor mit den eigentlichen Daten weitergearbeitet wird, soll zunächst die Korngrößendreieckfunktion
 # thematisiert werden.
 # Hierzu wird das Beispiel aus dem Handbuch des packages übernommen, um zunächst zu verstehen,
 # wie die Daten aufgebaut sein sollten.
 # Das Handbuch findet sich unter folgendem Link:
 # ftp://ftp.rediris.org/mirror/CRAN/web/packages/soiltexture/vignettes/soiltexture_vignette.pdf

 # Zunächst muss das erforderliche Package 'soiltexture' geladen werden.
 install.packages("soiltexture")

 # Package aktivieren
 library(soiltexture)

 # Beispieldataframe zum besseren Verständnis der erforderlichen Datenstruktur.
 # Entnommen aus: Handbuch 'soiltexture' [verändert]
 # Link: https://cran.r-project.org/web/packages/soiltexture/soiltexture.pdf
 KGS <- data.frame("CLAY"  = c(05,60,15,05,25),
                       "SILT"  = c(05,08,15,25,55),
                       "SAND"  = c(90,32,70,70,20),
                       "Org.Car."    = c(20,14,15,05,12),
                       'Horizontnummer'= c(1,2,3,4,5))

 # Ansehen des Dataframes
 View(KGS)

 # strukturprüfung des Dataframes
 str(KGS)

 # Das Package erfordert die Eingabe der drei Korngrößen zeilenweise sortiert nach Profiltiefe und
 # Spaltenweise sortiert nach Korngröße

 TT.plot(class.sys = "DE.BK94.TT", # Argument für Konrgrößendreieck nach KA5
         class.p.bg.col  = TRUE,
         tri.data = KGS, # Dataframe der dargestellt werden soll
         main = "Soil texture data", # Überschrift
         z.name = "Horizontnummer") # Dritte Variable anhand derer unterschieden wird, hier die einzelnen Horizonte

 # Die Punkte im Konrgrößendreieck visualisieren nach Größe aufsteigend, d.h. je tiefer der Horizont, desto,
 # größer der Punkt.

 # Nun soll das auf die eigentlichen Daten übertragen werden. Dazu werden die drei Hauptkorngrößen benötigt.
 # Da die Korngrößen im LP_phys1 Datenbeispiel teilweise (Ton nicht) stärker differenziert sind,
 # müssen sie zunächst zusammengefasst werden.

 str(LP_phys1)
 # Wie Du in der Ausgabe des str() sehen kannst, sind die meisten Parameter noch als 'chr'
 # formatiert. Diese sollen nun in einem Schritt zu numerischen Variablen werden.
 # Mittels sapply wird die as.numeric Funktion gleichzeitig auf alle angegebenen
 # Spalten angewendet.
 LP_phys1[,c(7:13)]=as.data.frame(sapply(LP_phys1[,c(7:13)], as.numeric))

 # Dazu werden die Grob- Fein und MIttel- Anteile der einzelnen Fraktionen addiert.
 LP_phys1$UProzent <- LP_phys1$gU_korrigiertProzent + LP_phys1$mU_corrProzent + LP_phys1$fU_corrProzent
 LP_phys1$SProzent <- LP_phys1$gS_Prozent.1 + LP_phys1$mS_Prozent.1 + LP_phys1$fS_Prozent.1

 KGS2 <- LP_phys1[,c(4,5,13:15)]
 KGS2 <- KGS2[24:28,]

 TT.plot(class.sys = "DE.BK94.TT", # Argument für Konrgrößendreieck nach KA5
         class.p.bg.col  = TRUE,
         tri.data = KGS2, # Dataframe der dargestellt werden soll
         main = "Korngrößendreieck", # Überschrift
         z.name = "Horizontobergrenze.cm.") # Dritte Variable anhand derer unterschieden wird, hier die einzelnen Horizonte


 # Wie du sehen kannst funktioniert der Plot in dieser Form noch nicht. Erforderlich ist zudem,
 # dass die Spalten eindeutig mit den englischen Fraktionsnamen 'SAND', 'SILT', und 'CLAY' benannt werden.

 # Dies geschieht nachfolgend.
 str(KGS2)
 colnames(KGS2) <- c('Flaeche','Horizontobergrenze','CLAY','SILT','SAND')
 str(KGS2)

 # Nun kann die Plotfunktion erneut abgeschickt werden
 TT.plot(class.sys = "DE.BK94.TT", # Argument für Konrgrößendreieck nach KA5
         class.p.bg.col  = TRUE,
         tri.data = KGS2, # Dataframe der dargestellt werden soll
         main = "Korngrößendreieck", # Überschrift
         z.name = "Horizontobergrenze.cm.") # Dritte Variable anhand derer unterschieden wird, hier die einzelnen Horizonte

 # Der nächste Fallstrick ist die folgende Fehlermeldung:
 # Error in TT.data.test(tri.data = tri.data, css.names = css.names, text.sum = text.sum,  :
 # The sum of the 3 plotted variables should be around 100: check the data, or change 'text.tol' parameter.

 # Der Hintergrund ist, dass die Werte im Labor rechnerisch über Differenzen bestimmt werden.
 # Da neben Siebung auch die Pipettmethode zum Einsatz kommt und dabei Korngrößenanteile
 # unter- bzw. überschätz werden, wird die folgende Fehlermeldung generiert.
 # Das Package erlaubt jedoch nur Werte die addiert über die drei Fraktionen annähernd 100 ergeben.
 # Entweder müssen die Daten korrigiert werden oder andernfalls muss die default-Einstellung dieser
 # Beschränkung deaktiviert werden.
 # Dies geschieht mit dem Argument: tri.sum.tst = F

 TT.plot(class.sys = "DE.BK94.TT", # Argument für Konrgrößendreieck nach KA5
         class.p.bg.col  = TRUE,
         tri.data = KGS2, # Dataframe der dargestellt werden soll
         tri.sum.tst=F,
         main = "Korngrößendreieck", # Überschrift
         z.name = "Horizontobergrenze") # Dritte Variable anhand derer unterschieden wird, hier die einzelnen Horizonte

 # Nun konnte das Konrgrößendreieck für die hiesigen Daten erstellt werden.

 # Noch ein Hinweis zum Erstellen von Grafiken in R: Wenn Du die
 # Abbildungen aus R exportieren möchtest, findest Du im RLab eine
 # gute Anleitung dazu.
 # https://rlab.blogs.uni-hamburg.de/dig-skripte/Tipps_und_Tricks/index.html

 # Bei Fragen ziehe gerne das Internet zurate oder wende Dich im
 # RLab-Blog mit der Kommentarfunktion an andere User oder das
 # RLab-Team (https://rlab.blogs.uni-hamburg.de/)

 #####################################################################
 #### ENDE

Du hast Fragen zu diesem Code? Nutze gerne die Kommentar-Funktion!

Trennzeichen in .csv-Dateien - Ursache für fehlerhaften Datenimport
.csv-Dateien besitzen (leider) nicht immer identische Zeichen für die Trennung von Spalten und Dezimalstellen. Häufig werden Spalten durch Semikola getrennt, üblich sind aber auch Kommata oder andere Zeichen. Die Dezimalstellen werden meistens mit Punkt oder Komma abgetrennt, aber auch hier sind prinzipiell andere Zeichen möglich.

Bei der Anzeige in MS-Excel ist die Identifikation der Trennzeichen nicht eindeutig möglich. Wir empfehlen, eine .csv-Datei dazu mit einem Texteditor (z.B. Notepad oder Editor) zu öffnen. In der Anzeige dieser Programme sind die Zeichen eindeutig ablesbar und können per copy&paste in den R-Code übernommen werden: Beim Einlesen einer .csv-Datei in R werden die Trennzeichen in der Funktion read.table() mit den Argumenten sep = und dec = angegeben. Für die Beispieldateien sind in den Code-Beispielen die passenden Trennzeichen angegeben. Bei anderen Dateien ist es ggf. erforderlich, sie anzupassen.

Ausführliche Infos zu .csv-Dateien findest Du im Digitalen Skript ,,Import von Daten in R” auf der Seite csv-Dateien im Detail.

R-Code zur Analyse der atmosphärischen Deposition

Hier findest Du R-Code, mit dem Du im ,,Datenmaster Depo" enthaltene Werte von Leitprofilen analysieren und Ergebnisse visualisieren kannst.

Zum Ausführen des Codes erforderliche Beispieldaten

bereinigte Daten der atmosphärischen Deposition, wie sie im R-Code (siehe unten) verwendet werden (zum Download als .csv-Datei ggf. Rechtsklick und "Speichern unter...")

Falls Du die Vorbereitung selbst nachvollziehen möchtest, steht hier auch der Rohdatensatz zur Verfügung:

Rohdaten der atmosphärischen Deposition, die vor den R-Berechnungen mit z.B. MS-Excel bereinigt werden müssen (zum Download als .xlsx-Datei ggf. Rechtsklick und "Speichern unter...")

Code

Kurz-Version des Depo-Skripts, ohne Kommentare (zum Download als .R-Datei ggf. Rechtsklick und "Speichern unter...")

Auch den kommentierten Code gibt es als R-Skript zum Download:

Kommentierter Code des Depo-Skripts (zum Download als .R-Datei ggf. Rechtsklick und "Speichern unter...")

Alternativ kannst Du den kommentierten Code hier direkt aus dem Digitalen Skript übernehmen (Copy&Paste):


 # Course: RLab - Labor
 # Lesson: Atmosphaerische Deposition
 # Author: RLab-Team
 # Co-Authors: Ronja Gottschalk (programming, design), Niels Schwab (design),
 #             Elke Fischer (design)
 # Organization: Universität Hamburg, project 'RLab - Skriptbasierte modulare
 #               Umweltstatistik' Funded within 'Lehrlabor' of
 #               Universitätskolleg 2.0 by BMBF (01PL17033).

 #####################################################################

 # Herzlich Willkommen im R-Skript zur Analyse der Atmosphaerischen Deposition.
 # Bevor Du mit der Arbeit beginnst, schließe bitte noch einen
 # Vorschritt ab, indem Du auf File --> 'Reopen with Encoding' klickst.
 # Wähle in dem dann erscheinenden Dialogfenster bitte 'UTF-8' aus.
 # Damit wird das Skript erneut in der UTF-8 Kodierung geladen,
 # sodass Umlaute, etc. korrekt dargestellt werden. Um dieses Skript
 # und die darin enthaltenen Kommentare gut zu verstehen, ist ein
 # grundlegendes Verständnis für die Nutzung von R und RStudio
 # erforderlich. Solltest Du Dir diese aneignen oder sie auffrischen
 # wollen, ist die Bearbeitung der Grundlagen-Inhalte des RLab-Angebots
 # (https://rlab.blogs.uni-hamburg.de/grundlagen/) und ggf. weiterer
 # RLab-Kurse (https://rlab.blogs.uni-hamburg.de/rlab-kurse/) zu empfehlen.


 # #### Vorarbeit ----------------------------------------------------------

 # Workspace aufräumen
 rm(list = ls())

 # Kontrolle des Workspace
 ls()

 # Abfrage des aktuellen 'Working Directory'
 getwd()

 # neues Working Directory festlegen
 setwd('F:/RLabPortable/R-Labor/Klima/AD/')
 # Beachte hier bitte, dass Du Deinen Dateipfad entsprechend Deiner
 # Ordnerstruktur anpassen musst.
 # Gute Hilfen zum Einlesen von Daten erhältst Du im RLab-Kurs
 # 'Daten einlesen'.

 # erneute Abfrage des 'Working Directory' zur Kontrolle, ob die
 # Änderungen übernommen wurden
 getwd()

 # Datensatz einladen
 AD <- read.table('2016_Depo_Roptimiert.csv',
                  na.strings=c('',' ','NA'),
                  header = T,
                  sep = ';',
                  dec = ',',
                  stringsAsFactors = F)

 # Auch hier sollte der Dateipfad entsprechend Deines Verzeichnisses
 # angepasst werden, sofern Deine Datei in einem anderen Ordner liegt,
 # als Dein Working Directory eingestellt ist.
 # na.strings definierte leere Zellen (''), Blanks (' ') und NA-Werte
 # ('NA') für R bereits beim Einlesen als 'NA' übergeben, andernfalls
 # werden diese bloß als Zeichenkette erkannt.
 # stringsAsFactors = F verhindert, dass Variablen automatisch als
 # Faktoren übernommen werden.
 # Ggf. sep = und dec = anpassen, weitere Infos dazu im Hinweis-
 # Feld am Ende der hier verlinkten Seite:
 # https://rlab.blogs.uni-hamburg.de/dig-skripte/Vorbereitung_Datenmaster/index.html?s=Schritt%202:%20Doppelte%20Kopfzeile%20entfernen
 #
 # Eine Anleitung zum Aufbereiten Deines Datensatzes für R findest Du
 # im RLab in den digitalen Skripten.


 # Überblick über den Datensatz  -------------------------------------------

 # Datensatz in R als Tabelle ansehen
 View(AD)

 # Struktur des Datensatzes ansehen
 str(AD)
 # die Mehrheit der Variablen wird als 'chr', also character
 # angegeben.
 # Gewollt sind hier für die Spalten mit Zahlenwerten jedoch
 # (wie im weiteren Verlauf noch deutlich wird) numerische Variablen.
 # Dies wird im weiteren Verlauf entsprechend zugewiesen.


 # Überprüfung des Datensatzes ---------------------------------------------

 # Inhalt:
 # Test auf fehlende Werte (NA) - Variablen als numerisch ausweisen und
 # Datumsformat einstellen

 #### 'NA'- WERTE

 # Abfrage, ob für die numerischen Variablen NAs erkannt werden
 table(is.na(AD[,c(6:9,11:27)])) # False  6479 - True 3958
 # Es werden Variablen als NAs erkannt.
 # Zusätzlich ist es immer ratsam eine Sichtkontrolle im Datenblatt
 # 'AD' durchzuführen, um die Daten auf eventuell leere Zellen zu
 # kontrollieren. Das würde dann darauf hindeuten, dass für eine
 # bestimmte Zeichenfolge für R nicht erkennbar ist, dass es sich
 # dabei um 'NA' handelt. In diesem Fall sollte der Datenmaster dort,
 # wo keine Zahlenwerte vorliegen, vollständig mit 'NA' befüllt sein.

 # Nachfolgend sollen dennoch kurz einige Möglichkeiten vorgestellt
 # werden, wie man bei einer großen Datenmenge spezifischer überprüfen
 # kann, ob NA-Werte vorliegen oder nicht. Dies kann ggf. zur
 # Identifikation noch nicht als 'NA' definierter Zeichenketten
 # hilfreich sein.

 # Beispielhaft soll nun anhand einer Variable mit Datenlücken
 # überprüft werden, ob auch die Datenlücken als 'NA' erkannt werden.
 is.na(AD$Aluminium_Al_.mg.l.)
 # Die Logikabfrage gibt für einige Zellen ein 'True' für die Abfrage
 # nach NAs aus. NAs werden also als solche erkannt.

 # Vorangehend wurden Dir somit eine Möglichkeit gezeigt, wie Du bei
 # großen Datensätzen ein Gefühl dafür bekommen kannst, ob Deine
 # Daten hinsichtlich der fehlenden Werte plausibel sind. Dies ist
 # deshalb wichtig, weil Du mit den Daten nicht vernünftig rechnen
 # kannst, wenn undefinierte Zellen vorliegen oder Fehlwerte wie
 # z.B. '-999' fälschlicherweise als gültige Werte behandelt werden
 # und Dir Deine statistischen Parameter verzerren.
 # Generell lässt sich festhalten, dass es hierfür keine
 # 'Pauschal-Funktion' in R gibt, da R case-sensitiv ist
 # (also z.B. Groß- und Kleinschreibung in R bereits Unterschiede
 # machen; z.B. ist 'Test' für R ungleich 'test')
 # und insofern unendlich viele Zeichenketten korrekterweise
 # oder fälschlicherweise als 'NA' erkennen kann oder eben nicht.
 # Von daher ist es erforderlich, dass Du Dich im Vorwege mit Deinen
 # Daten auseinandersetzt, um beurteilen zu können, welche Werte
 # möglicherweise fehlerhaft von R verstanden werden.


 #### Variablen als numerisch ausweisen
 #### Umwandlung von 'chr' zu 'num'

 # Erneute Abfrage der Struktur von 'AD'
 str(AD)
 # Die Abfrage ergibt, dass nahezu alle Spalten das Format 'chr'
 # (character) haben. Das bedeutet das R die Werte als Buchstaben
 # behandelt.
 # Da hier Zahlenwerte vorliegen, ist es für weitere Operationen in
 # R erforderlich, dass Spalten, welche ausschließlich Zahlenwerte
 # enthalten, auch als numerisches Format von R erkannt werden.  

 # Daher soll mittels der Funktion as.numeric an R übergeben werden,
 # welche Spalten numerisch sind.
 # In diesem Fall sind das die Spalten 6-9 und 11-27.
 AD[,c(6:9,11:27)]=as.data.frame(sapply(AD[,c(6:9,11:27)], as.numeric))
 # Die Warning Message: 'NAs introduced by coercion' kann ignoriert
 # werden. Sie dient lediglich dem Hinweis, dass es bei der Auswahl
 # Zellen gibt, welche NA-Werte enthalten.

 # Erneutes Testen der Struktur
 str(AD)
 # Die Abfrage ergibt nun, dass die Spalten 6-9 und 11-27 das Format
 # 'num' (numeric) haben. R erkennt die numerischen Spalten nun
 # also wie gewünscht in dem Datenformat.

 #### Datumsformat einstellen

 # Bezüglich der Datenformate in R ist nun noch ein weiterer Schritt
 # durchzuführen. Dir ist sicherlich aufgefallen, dass es eine Spalte
 # namens 'Entnahmedatum' gibt.

 # Abfrage der Struktur von 'Entnahmedatum'
 str(AD$Entnahmedatum)
 # Die Abfrage ergibt, dass die Spalte als 'character' definiert wird.
 # Hier ist das Datumsformat (date oder posixct) gewünscht.
 # Da in diesem Fall keine Zeitangaben vorliegen, sondern lediglich das
 # Datum bietet sich die Konvertierung in das 'date'- Fomrat an
 # Zudem unterstützt das Grafik-Package ggplot nur 'date'-Formate und keine
 # 'posixct'-Formate

 # 'Entnahmedatum' als Datumsformat umschreiben und damit die
 # Datumsspalte überschreiben.
 AD$Entnahmedatum=as.Date(AD$Entnahmedatum, format='%d.%m.%Y')

 # Finales Überprüfen der Struktur
 str(AD$Entnahmedatum)
 # Die Strukturabfrage der variable ergibt, dass sie nun als
 # Datumsspalte geführt wird.

 # Zudem noch ein Zusatz, mit dem Du komplizierte oder ungünstig
 # geschriebene Spaltennamen verändern kannst.
 # Beispielhaft soll nun der Spalte 6 'Vol_.ml.' ein neuer Name
 # zugewiesen werden
 colnames(AD)[6] <- 'Volumen[ml]'

 # Zur Kontrolle kannst Du Dir den Datensatz nun erneut ansehen.
 View(AD)


 # Grafiken generieren (ggplot) --------------------------------------------
 # Inhalt: Vorbereitung: Aufteilen des Datensatzes - Boxplots

 # Hier werden Boxplots erstellt, die beispielhaft die pH-Werte der
 # unterschiedlichen Standorte 'Buche', 'Douglasie' und 'Freifläche'
 # darstellen.
 # Die besondere Herausforderung dabei ist, dass die Daten in
 # zwei "Ebenen" aufgeteilt, nebeneinander dargestellt werden sollen:
 # Zum einen werden die pH-Werte geplottet und zum anderen
 # gleichzeitig nach Standort differenziert. In ggplot funktioniert das
 # für die Unterscheidung nach Standorten ohne das Datenset aufzuteilen.
 # Möchte man aber Werte z.B. nur für den Buchenstandort plotten, empfiehlt
 # es sich, den Dataframe vorher in drei Einzeldatensätze nach Standort zu
 # unterteilen. Daher werden hier beide Möglichkeiten aufgezeigt.

 #### Vorbereitung - Aufteilen des Datensatzes

 ## Datensatz aufteilen nach 'Standort' -  split.dataframe()

 # Diese Variante ist gut geeignet, wenn die Daten in einer Spalte
 # angeordnet sind und es somit einen gemeinsamen 'Identifier' gibt,
 # nach dem aufgeteilt werden kann. In diesem Fall sind die Daten günstig
 # strukturiertist und es kann direkt anhand der variable Standort
 # aufgesplittet werden.

 Standortart<-split.data.frame(AD, AD$Standort)

 # Struktur des neu generierten Objekts prüfen, der Datensatz liegt
 # nun in 'gesplitteter' Form vor
 str(Standortart)

 # Die Objekte je Standort gesondert generieren
 Buche <- Standortart$Buche
 Douglasie <- Standortart$Douglasie
 Freiflaeche <- Standortart$Freiflaeche

 # Kurze Strukturprüfung und statistischer Überblick über die drei
 # Objekte
 str(Buche)
 summary(Buche)
 str(Douglasie)
 summary(Douglasie)
 str(Freiflaeche)
 summary(Freiflaeche)

 ## Datensatz aufteilen nach "Standort" - VARIANTE 2 "Filtern" /
 ## Subsetten mit []

 # Diese Variante eignet sich, wenn Daten extrahiert werden sollen,
 # auf die eine bestimmte Bedingung zutrifft.
 # In diesem Fall sollen die Subsets die Values für die Zeilen
 # rausschreiben, für die gilt Standort =='Buche', bzw. == 'Douglasie'
 # oder == 'Freiflaeche'.
 # Die Stärke dieser Funktion liegt allerdings in ihrer
 # Adaptierbarkeit, so können auch Bedingungen wie > oder < in die
 # Abfrage eingebaut werden.

 # In das Objekt Buche 2 werden nun beispielhaft  alle Zeilen aus dem
 # Datensatz 'AD'herausgeschrieben, für die innerhalb der Spalte 'Standort'
 # zutrifft, dass sie den Standort 'Buche' haben (=='Buche')
 Buche2 <- AD[AD[,'Standort']=='Buche',]
 # Strukturen überprüfen
 str(Buche2)

 # Nachfolgend soll es um die Visualisierung der Daten und
 # deren statistische Überprüfung gehen.

 # Visualisierung - Plots ----------------------------------------------------------------
 #### Boxplots


 # Das Generieren von Boxplots mit den in R implementierten Grafikoptionen
 # wurde im Script 'Bodenloesung' eingehender thematisiert.
 # Nachfolgend soll nun ausschließlich mit dem Package 'ggplot2' gearbeitet
 # werden, um Plots zu erzeugen.
 # ggplot2 ist ein Package für grafische Darstellungen in R, das von
 # der Funktionsweise her allerdings etwas anders aufgebaut ist,
 # als die sonstigen Funktionen für Plots in R.

 # ggf. Package ggplot2 installieren
 #install.packages('ggplot2')
 # Package ggplot2 aktivieren
 library(ggplot2)
 str(AD)

 # Zunächst soll ein gruppierter Boxplot erstellt werden
 # Dazu wird mittels 'ggplot' das Package 'ggplot' angesteuert. geom_ legt
 # dann fest, welche Art des Plots folgen soll, in unserem Fall also
 # geom_boxplot
 # aes legt dann die 'aesthetics' des Plots fest. Hier übergibst Du
 # die Daten an R, welche geplottet werden sollen und legst fest,
 # welche Daten auf welcher Achse abgetragen werden sollen.

 ggplot(AD, aes(x='', y=pH,
                   fill=Standort)) +
   theme(axis.text.x=element_text(angle=60, hjust=1)) +
   geom_boxplot()

 # Alternativ bietet sich hier an die einezlen Standorte in einem Frame
 # voneinander abgesondert darzustellen. Dies geschiet mit dem Argument
 # facet-wrap

 ggplot(AD, aes(x='', y=pH,
                fill=Standort)) +
   facet_wrap(~Standort) +
   theme(axis.text.x=element_text(angle=60, hjust=1)) +
   geom_boxplot()

 # Auch das Plotten eines einzelnen Standorts ist möglich
 # So kannst Du hier exemplarisch die pH-Werte des Standorts 'Buche'
 # einzeln über die Zeit darstellen. Dazu kannst Du als Daten
 # einfach den Subset 'Buche' ansteuren

 ggplot(Buche, aes(x=Entnahmedatum, y=pH)) +
   theme(axis.text.x=element_text(angle=60, hjust=1)) +
   geom_boxplot()

 # Die eigentlich interessierende Darstellung ist jedoch das Plotten von Zeitreihen,
 # also der Darstellung des Datums auf der X-Achse
 # Dazu sollen Gruppierte Boxplots über die Zeit erstellt werden

 ggplot(AD, aes(x=Entnahmedatum, y=pH,
                fill=Standort)) +
   theme(axis.text.x=element_text(angle=60, hjust=1)) +
   scale_x_date(breaks = '1 month') +
   geom_boxplot()

 # Zunächst wird das Datenset mit einem Subset auf einen interessierenden zweiwöchigen
 # Zeitraum gekürzt


 str(AD)
 # Subset Zeitreihe 2 Wochen
 AD_sub <- AD[AD$Entnahmedatum >= '2016-01-08' &
                          AD$Entnahmedatum <= '2016-01-22',]
 # Damit das Datum auf der X-Achse korrekt dargestellt wird, wird die Variable
 # Entnahmedatum nun wieder als 'chr' formatiert.
 # In ggplot gibt es die Funktion Daten in verschiedenen Varianten darzustellen,
 # da hier zwei-Wöchige Zeitabstände gefragt sind, muss dieser Umweg genommen werden,
 # da die Darstellung
 AD_sub$Entnahmedatum_chr <- as.character(AD_sub$Entnahmedatum)
 str(AD_sub)

 # Das bietet den Vorteil, dass Du genau kontrollieren kannst, welche Daten in den Plot
 # eingehen. Zwar gibt es in ggplot implementierte Argumente zum 'aufbrechen' von Zeitreihen
 # Jedoch ist deren Funktionsweise gerade anfangs sehr undurchschaubar

 # Nun soll die gesubsetette Variable mittels Boxplots dargestellt werden
 ggplot(AD_sub, aes(x=Entnahmedatum_chr, y=pH,
                fill=Standort)) +
   theme(axis.text.x=element_text(angle=60, hjust=1)) +
   geom_boxplot()

 # Der Zeitrum könnte selbstverständlich auch beliebig größer gewählt werden
 # Dazu müsstest du lediglich das Subset anpassen

 # Nun soll das Entnahmedatum noch eine vernünftige Achselnbeschriftung erhalten.
 # Dies geschiet mittels 'xlab'
 ggplot(AD_sub, aes(x=Entnahmedatum_chr, y=pH,
                    fill=Standort)) +
   theme(axis.text.x=element_text(angle=60, hjust=1)) +
   geom_boxplot() +
   xlab('Entnahmedatum')

 # Nun kannst Du noch einen Titel für Deinen Plot hinzufügen.
 # Beachte, dass man jedachdem, wo die Grafik später eingebunden wird,
 # nicht zwangsläufig einen Titel braucht
 ggplot(AD_sub, aes(x=Entnahmedatum_chr, y=pH,
                    fill=Standort)) +
   theme(axis.text.x=element_text(angle=60, hjust=1)) +
   geom_boxplot() +
   xlab('Entnahmedatum') +
   ggtitle('pH-Werte differenziert nach Standort \n im Zeitraum vom 08.01. - 22.01.2016')

 # Damit konntest Du bereits Plots erstellen.
 # Noch ein Hinweis zum Erstellen von Grafiken in R: Wenn Du die
 # Abbildungen aus R exportieren möchtest, findest Du im RLab eine
 # gute Anleitung dazu.
 # https://rlab.blogs.uni-hamburg.de/dig-skripte/Tipps_und_Tricks/index.html

 # Abschließend soll nun noch darauf eingegangen werden, wie Du Variablen
 # auf signifikante Unterschied beispielsweise zwischen den Gruppen testen kannst.

 # Voraussetzungen für Test auf Mittelwertsunterschiede --------------------
 # 1. Normalverteilung - 2. Varianzhomogenität

 #### Normalverteilung

 # Für weitere Berechnungen aber auch für einen Eindruck über die Verteilung der
 # Daten ist es notwendig die Daten auf Normalverteilung zu testen.
 # Zum einen kann dies visuell geschehen, alternativ stehen Testverfahren wie der
 # Lillefors-Test zur Verfügung
 # Nachfolgend soll exemplarisch der Lillefors-Test angewendet werden und die
 # Daten anschließend in Abhängigkeit des Vorliegens von Normalverteilung näher
 # analysiert werden
 # Hierzu wird beispielhaft die Variable 'Aluminiumgehalt' verwendet

 # Lillie-Test
 # ggf. package installieren
 #install.packages("nortest")
 # Package aktivieren
 library(nortest)

 colnames(AD)[11] <- 'Al'
 str(AD)

 # Exemplarisch sollen nun die Aluminiumwerte für den Gesamtdatensatz 'AD' auf
 # Normalverteilung getestet werden.
 lillie.test(AD$Al) # p-value < 2.2e-16 --> nicht normalverteilt, da p < 0.05


 #### Test auf Varianzhomogenität

 bartlett.test(AD$pH, g = AD$Standort)
 #Bartlett's K-squared = 267.75, df = 1, p-value < 0.001443, also < 0.05, ergo varianzheterogen

 # Die Daten sind demnach varianzheterogen und nicht normalverteilt. Demanch sind
 # die verwendbaren Verfahren auf solche
 # beschränkt, die diese Voraussetzungen nicht benötigen.
 # Daher wird nachfolgend der Kruscal- Wallis Test angewendet.


 #### Test auf signifikante Unterschiede zwischen den Standorten.

 # Bevor der Kruscal Test gerechnet wird, sieh Dir noch einmal die Struktur der
 # AD- Gesamtdatensatzes,
 # insbesondere hinsichtlich des Standortes an.
 str(AD)
 # Wie Du sehen kannst, ist die Variable 'Standort' im Gesamtdatensatz im Format
 # 'chr'. Für alle weiteren Testverfahren wird die Variable im Format 'factor' benötigt.
 # Daher soll die Variable zunächst ins Faktorformat umgeformt werden und an-
 # schließend die Struktur erneut geprüft werden.
 AD$Standort <- as.factor(AD$Standort)
 str(AD)

 kruskal.test(AD$Al, g = AD$Standort)
 # Kruskal-Wallis chi-squared = 30.885, df = 1, p-value = 2.854e-11
 # das heißt p < 0.05. Somit gibt es signifikante Unterschiede zwischen den drei
 # Standorten hinsichtlich der Aluminiumwerte.
 # Wenn mehrere nicht normalverteilte Variablen auf Unterschiede getestet werden sollen,
 # bietet sich der Kruscal-Wallis Test an, um zwischen verschiedenen Gruppen auf
 # signifikante Unterschiede zu testen.

 # Um zu identifizieren, wo diese Unterschiede liegen, kannst Du nun den Post-Hoc
 # Test des Kruscal-Wallis Tests durchführen.
 # Dazu wird zunächst das Package 'pgirmes' benötigt.

 # ggf. Installieren des Packages pgirmess
 #install.packages('pgirmess')

 #Aktivieren des Packages pgirmes
 library(pgirmess)

 kruskalmc(AD$Al ~ AD$Standort) #für diesen Posthoctest package pgirmess installieren

 # Multiple comparison test after Kruskal-Wallis
 #p.value: 0.05
 #Comparisons
 #obs.dif critical.dif difference
 #Buche-Douglasie        28.34517     28.80798      FALSE
 #Buche-Freiflaeche      89.00336     40.39357       TRUE
 #Douglasie-Freiflaeche 117.34853     40.32850       TRUE

 # Die logische Ausgabe 'False' oder 'True' gibt nun jeweils an, ob es signifi-
 # kante Unterschiede hinsichtlich der Aluminumwerte für das jeweilige Paar gibt.
 # So gibt es beispielsweise signifikante Unterschiede in den Aliminiumwerten
 # zwischen den Standorten 'Buche' und 'Freifläche'.

 # Im Regelfall wirst Du Unterschiede für nicht-normalverteilte Variablen unter-
 # suchen. Sollten jedoch normalverteilte und varianzhomogene Daten vorliegen,
 # so kannst du hier mittels ANOVA auf Unterschiede testen.
 # Die Funktion dafür wäre wie folgt:

 # Variable_aov <- aov (Datensatz$Variable ~ Datensatz$VariableunterschiedlicheGruppen)
 # summary(Variable_aov)
 # hier muss man den Umweg gehen, das Ergebnis von aov() als Objekt abzuspeichern
 # und sich hinterher die summary davon anzuschauen
 # Bei Signifikanz des Tests wäre der zugehörige Post-Hoc Test:
 # TukeyHSD(aov(Datensatz$Variable ~ Datensatz$VariableunterschiedlicheGruppen))

 #####################################################################
 #### Hinweise

 # Damit sind nun die wichtigsten Analysen, Datenüberprüfungen und
 # Grafiken generiert.
 # Bei Fragen ziehe gerne das Internet zurate oder wende Dich im
 # RLab-Blog mit der Kommentarfunktion an andere User oder das
 # RLab-Team (https://rlab.blogs.uni-hamburg.de/)

 #####################################################################
 #### ENDE

Du hast Fragen zu diesem Code? Nutze gerne die Kommentar-Funktion!

Trennzeichen in .csv-Dateien - Ursache für fehlerhaften Datenimport
.csv-Dateien besitzen (leider) nicht immer identische Zeichen für die Trennung von Spalten und Dezimalstellen. Häufig werden Spalten durch Semikola getrennt, üblich sind aber auch Kommata oder andere Zeichen. Die Dezimalstellen werden meistens mit Punkt oder Komma abgetrennt, aber auch hier sind prinzipiell andere Zeichen möglich.

Bei der Anzeige in MS-Excel ist die Identifikation der Trennzeichen nicht eindeutig möglich. Wir empfehlen, eine .csv-Datei dazu mit einem Texteditor (z.B. Notepad oder Editor) zu öffnen. In der Anzeige dieser Programme sind die Zeichen eindeutig ablesbar und können per copy&paste in den R-Code übernommen werden: Beim Einlesen einer .csv-Datei in R werden die Trennzeichen in der Funktion read.table() mit den Argumenten sep = und dec = angegeben. Für die Beispieldateien sind in den Code-Beispielen die passenden Trennzeichen angegeben. Bei anderen Dateien ist es ggf. erforderlich, sie anzupassen.

Ausführliche Infos zu .csv-Dateien findest Du im Digitalen Skript ,,Import von Daten in R” auf der Seite csv-Dateien im Detail.

R-Code zur Analyse der Klimadaten von HOBO-Loggern

Hier findest Du R-Code, mit dem Du im ,,Datenmaster Klima" enthaltene Werte von HOBO-Loggern analysieren und Ergebnisse visualisieren kannst.

Zum Ausführen des Codes erforderliche Beispieldaten

bereinigte Klimadaten der HOBO-Logger, wie sie im R-Code (siehe unten) verwendet werden (zum Download der .zip-Datei ggf. Rechtsklick und "Speichern unter...", dann Datei entpacken)

Falls Du die Vorbereitung selbst nachvollziehen möchtest, steht hier auch der Rohdatensatz zur Verfügung:

Rohdaten der HOBO-Logger, die vor den R-Berechnungen mit z.B. MS-Excel bereinigt werden müssen (zum Download der .zip-Datei ggf. Rechtsklick und "Speichern unter...", dann Datei entpacken)

Code

Kurz-Version des HOBO-Skripts, ohne Kommentare (zum Download als .R-Datei ggf. Rechtsklick und "Speichern unter...")

Auch den kommentierten Code gibt es als R-Skript zum Download:

Kommentierter Code des HOBO-Skripts (zum Download als .R-Datei ggf. Rechtsklick und "Speichern unter...")

Alternativ kannst Du den kommentierten Code hier direkt aus dem Digitalen Skript übernehmen (Copy&Paste):


 # Course: RLab - Labor
 # Lesson: HOBO-Logger Auswertung
 # Author: RLab-Team
 # Co-Authors: Ronja Gottschalk (programming, design), Niels Schwab (design),
 #             Elke Fischer (design)
 # Organization: Universität Hamburg, project 'RLab - Skriptbasierte modulare
 #               Umweltstatistik' Funded within 'Lehrlabor' of Universitäts-
 #               kolleg 2.0 by BMBF (01PL17033).

 #####################################################################

 # Herzlich Willkommen im R-Skript zur Analyse der Hobo-Logger.
 # Bevor Du mit der Arbeit beginnst, schließe bitte noch einen
 # Vorschritt ab, indem Du auf File --> 'Reopen with Encoding' klickst.
 # Wähle in dem dann erscheinenden Dialogfenster bitte 'UTF-8' aus.
 # Damit wird das Skript erneut in der UTF-8 Kodierung geladen,
 # sodass Umlaute, etc. korrekt dargestellt werden. Um dieses Skript
 # und die darin enthaltenen Kommentare gut zu verstehen, ist ein
 # grundlegendes Verständnis für die Nutzung von R und RStudio
 # erforderlich. Solltest Du Dir diese aneignen oder sie auffrischen
 # wollen, ist die Bearbeitung der Grundlagen-Inhalte des RLab-Angebots
 # (https://rlab.blogs.uni-hamburg.de/grundlagen/) und ggf. weiterer
 # RLab-Kurse (https://rlab.blogs.uni-hamburg.de/rlab-kurse/) zu empfehlen.


 # #### Vorarbeit ----------------------------------------------------------

 # Workspace aufräumen
 rm(list = ls())

 # Kontrolle des Workspace
 ls()

 # Abfrage des aktuellen 'Working Directory'
 getwd()

 # neues Working Directory festlegen
 setwd('F:/RLabPortable/R-Labor/Klima/Hobos/')
 # Beachte hier bitte, dass Du Deinen Dateipfad entsprechend Deiner
 # Ordnerstruktur anpassen musst.
 # Gute Hilfen zum Einlesen von Daten erhältst Du im RLab-Kurs
 # 'Daten einlesen'.

 # erneute Abfrage des 'Working Directory' zur Kontrolle, ob die
 # Änderungen übernommen wurden
 getwd()

 # Datensatz einladen
 Hobos_FF <- read.table('Hobo_FF.csv',
                  na.strings=c('',' ','NA'),
                  header = T,
                  sep = ',',
                  dec = '.',
                  stringsAsFactors = F)

 Hobos_Buche <- read.table('Hobo_Buche.csv',
                        na.strings=c('',' ','NA'),
                        header = T,
                        sep = ',',
                        dec = '.',
                        stringsAsFactors = F)

 Hobos_Douglasie <- read.table('Hobo_Douglasie.csv',
                           na.strings=c('',' ','NA'),
                           header = T,
                           sep = ',',
                           dec = '.',
                           stringsAsFactors = F)

 # Auch hier sollte der Dateipfad entsprechend Deines Verzeichnisses
 # angepasst werden, sofern Deine Datei in einem anderen Ordner liegt
 # als Dein Working Directory eingestellt ist.
 # na.strings definierte leere Zellen (''), Blanks (' ') und NA-Werte
 # ('NA') für R bereits beim Einlesen als 'NA' übergeben, andernfalls
 # werden diese bloß als Zeichenkette erkannt.
 # stringsAsFactors = F verhindert, dass Variablen automatisch als
 # Faktoren übernommen werden.
 # !!!
 # !!! Bei HOBO-Daten unbedingt beachten:
 # Ggf. sep = und dec = anpassen, weitere Infos dazu im Hinweis-
 # Feld am Ende der hier verlinkten Seite:
 # https://rlab.blogs.uni-hamburg.de/dig-skripte/Vorbereitung_Datenmaster/index.html?s=Schritt%202:%20Doppelte%20Kopfzeile%20entfernen
 # !!!
 #
 # Eine Anleitung zum Aufbereiten Deines Datensatzes für R findest Du
 # im RLab in den digitalen Skripten.

 # Überblick über den Datensatz  -------------------------------------------

 # Datensätze in R als Tabelle ansehen
 View(Hobos_FF)
 View(Hobos_Douglasie)
 View(Hobos_Buche)

 # Strukturen der Datensätze ansehen
 str(Hobos_FF)
 str(Hobos_Douglasie)
 str(Hobos_Buche)

 # Teilweise werden Variablen als 'chr', also character
 # erkannt.
 # Gewollt sind hier für die Spalten mit Zahlenwerten jedoch
 # (wie im weiteren Verlauf noch deutlich wird) numerische Variablen.
 # Das kann je nach Deinen Default-Einstellungen, bereits der Fall sein.
 # Sicherhalbshalber kannst Du die Funktion as.numeric abschicken, damit die
 # numerischen Daten in jedem Fall im richtigen Format sind.

 Hobos_FF[,c(3,4)]= as.data.frame(sapply(Hobos_FF[,c(3,4)], as.numeric))
 Hobos_Buche[,c(3,4)]= as.data.frame(sapply(Hobos_Buche[,c(3,4)], as.numeric))
 Hobos_Douglasie[,c(3,4)]= as.data.frame(sapply(Hobos_Douglasie[,c(3,4)], as.numeric))

 # Teildatensätze anpassen und in einen Dataframe zusammenführen:  ---------

 # Ein neuer Dataframe 'Hobos_Buche2' wird aus 'Hobos_Buche' Spalten
 # 2 bis 4 generiert. Anschließend wird der ursprüngliche DF durch
 # den reduzierten DF ersetzt und die Überschriften der Spalten neu
 # geschrieben.
 Hobos_Buche2 <- Hobos_Buche[,c(2:4)]
 Hobos_Buche <- Hobos_Buche2
 colnames(Hobos_Buche) <- c('Datetime','Temp_B','RH_B')

 # Das Vorgehen entspricht hier dem bei 'Hobos_Buche'
 Hobos_Douglasie2 <- Hobos_Douglasie[,c(2:4)]
 Hobos_Douglasie <- Hobos_Douglasie2
 colnames(Hobos_Douglasie) <- c('Datetime','Temp_D','RH_D')

 Hobos_FF2 <- Hobos_FF[,c(2:4)]
 Hobos_FF <- Hobos_FF2
 colnames(Hobos_FF) <- c('Datetime','Temp_FF','RH_FF')

 # Die lokale Zeitzone, die sich hier am R ausführenden Betriebssystem orientiert,
 # wird gesetzt (=TRUE).
 Sys.timezone(location = TRUE)
 # Es folgt eine Kontrolle
 Sys.getlocale()

 # Es folgt die Zuweisung des Datumsformats.
 # Aufbau des Datumsformats nach %Month.%Day.%Year %Hour:%Minute:%Second
 # %p (steht für AM oder PM). %I in Kombination mit %p zeigt an,
 # dass es sich um das 12 stündige Datumsformat handelt.
 Hobos_FF$Datetime <- as.POSIXct(Hobos_FF$Datetime, format='%m.%d.%y %I:%M:%S %p')
 Hobos_Buche$Datetime <- as.POSIXct(Hobos_Buche$Datetime, format='%m.%d.%y %I:%M:%S %p')
 Hobos_Douglasie$Datetime <- as.POSIXct(Hobos_Douglasie$Datetime, format='%m.%d.%y %I:%M:%S %p')

 # Kontrolle der Datenstruktur auf korrekte Übernahme des Datumsformats
 str(Hobos_FF$Datetime)
 str(Hobos_Buche$Datetime)
 str(Hobos_Douglasie$Datetime)

 # Nun kannst Du Dir noch ansehen, wie die Zentral- und Streuungsmaße für
 # die Temperatur und Relative Luftfeuchte der einzelnen Standorte ausfallen
 summary(Hobos_FF)
 summary(Hobos_Buche)
 summary(Hobos_Douglasie)

 # Die Daten sehen zunächst plausibel aus. Würden dir hier unplausible Minima der Maxima
 # auffallen, wäre es hier angezeigt, die Werte als NA zu definieren.

 # In der Folge werden die drei DF (DF= Dataframe) zu einem Gesamt-DF 'Hobos' zusammengeführt.
 Hobos <- merge(Hobos_FF, Hobos_Buche, by = 'Datetime', all = T, fill = 'NA')
 Hobos <- merge(Hobos, Hobos_Douglasie, by = 'Datetime', all = T, fill = 'NA')


 # Aggregierung zu Stundenmitteln ------------------------------------------
 # Die Funktion 'aggregate' ist in der Lage eine Funktion auf Daten anzuwenden
 # und diese Stapelzuverarbeiten und dabei zusammenzufassen. Zunächst werden
 # hier mittels des Arguments 'hour' stündliche Mittel (mean) erstellt.
 # Durch das Setzen von '-1' wird die Zeile mit der Überschrift von der
 # Aggregierung ausgeschlossen.
 # Zusätzlich wird hier eine Liste angesteuert, damit mehrere Spalten
 # gleichzeitig aggregiert werden können.
 Hobos_hourly <- aggregate(list(Hobos$Temp_FF, Hobos$RH_FF,
                                Hobos$Temp_B, Hobos$RH_B,
                                Hobos$Temp_D, Hobos$RH_D),
                   list(hour=cut(as.POSIXct(Hobos$Datetime)-1, "hour")),
                   mean)

 # Aggregierung zu Tagesmitteln --------------------------------------------
 # Hier wird durch Anpassung des Arguments 'hour' entsprechnend des Vorhabens
 # das Argument zu'day' veränder, sodass nun  Tagesmittel gebildet werden.
 Hobos_daily <- aggregate(list(Hobos$Temp_FF, Hobos$RH_FF,
                               Hobos$Temp_B, Hobos$RH_B,
                               Hobos$Temp_D, Hobos$RH_D),
                          list(hour=cut(as.POSIXct(Hobos$Datetime)-1, "day")),
                          mean)

 # Mittels colnames() werden die Spaltennamen neu gesetzt, da diese durch die
 # Aggregierung verfremdet werden.  
 colnames(Hobos_daily) <- c('Datetime', 'Temp_FF','RH_FF', 'Temp_B','RH_B', 'Temp_D','RH_D')

 # install.packages('tidyr') - Wenn noch nicht geschehen, Installation des Packages 'tidyr'
 library(tidyr)

 # Im aktuellen DF 'Hobos_daily' oder 'Hobos_hourly' sind die unterschiedlichen Standorte
 # in Spalten nebeneinander angeordnet. ZU Visualisierung mit 'ggplot' ist es jedoch
 # notwendig, dass die Spalten  untereinander in einer Spalte angeordnet werden.
 # Dies geschieht nachfolgend mittels der Funktion 'gather' (aus dem Package 'tidyr')
 # Durch 'gather' wird ein neuer Identifier eingefügt, hier der 'Standort'.
 # Damit werden die ehemaligen Spaltenüberschriften entsprechend der Anordnung
 # untereinander in der neuen Identifiery Spalte angeordnet.

 # Die 'gather'-Funktion verbindet inmmer die angegebenen Spalten zu einer
 # gemeinsamen Spalte.
 # Daher werden hier zunächst die Temperaturspalten und dann separat die Spalten zur
 # relativen Luftfeuchte (RH) in zwei unterschiedlichen Objekten zusammengeführt.
 # Da die Objekte nachfolgend in einen Gesamt-DF zusammengeführt werden sollen,
 # werden sie direkt auf die relevanten Spalten gesubsettet.

 hobos_reshaped_temp <- gather(Hobos_daily, Standort, Temp,
                               'Temp_FF', 'Temp_B', 'Temp_D', factor_key=TRUE)
 hobos_reshaped_temp <- hobos_reshaped_temp[,c(1,5,6)]

 hobos_reshaped_RH <- gather(Hobos_daily, Standort, RH,
                               'RH_FF', 'RH_B', 'RH_D', factor_key=TRUE)
 hobos_reshaped_RH <- hobos_reshaped_RH[,c(1,5,6)]

 # Es folgt die erneute Zuweisung des Datumsformats - %Year-%Month-%Day
 hobos_reshaped_temp$Datetime <- as.POSIXct(hobos_reshaped_temp$Datetime, format='%Y-%m-%d')
 hobos_reshaped_RH$Datetime <- as.POSIXct(hobos_reshaped_RH$Datetime, format='%Y-%m-%d')

 # Kontrolle der Datenstruktur
 str(hobos_reshaped_temp)
 str(hobos_reshaped_RH)

 # Zusammenführen der beiden DF zu einem Gesamtobjekt, by = Datum, also der
 # Identifyer zum Zusammenführen
 hobos_daily_reshaped <- merge(hobos_reshaped_temp, hobos_reshaped_RH,
                               by='Datetime', all = T, fill = 'NA')

 # Kontrolle der Struktur der Daten
 str(hobos_daily_reshaped)

 # Erneutes Zuweisen der Spaltennamen
 colnames(hobos_daily_reshaped) <- c('Datetime', 'Standort','Temperatur', 'Standort2','RH')

 # Nun kannst Du Dir den neu erzeugten und umstrukturierten Datensatz ansehen und
 # die str() prüfen
 View(hobos_daily_reshaped)
 str(hobos_daily_reshaped)

 # Nun kannst Du Dir noch die deskpitive Statistik der Tagesmittel ausgeben lassen.
 summary(hobos_daily_reshaped)
 sd(hobos_daily_reshaped$Temperatur, na.rm = T)
 sd(hobos_daily_reshaped$RH, na.rm = T)

 # Visualisierung des DF mittels 'ggplot2'
 # Zur Anwendung schaue Dir gerne den RLab-Kurs 'Daten visualisieren mit
 # ggplo2' an. Diesen findest Du, wie alle RLab-Kurse im RLab-Blog.
 # Hier werdn exemplarisch die Tagesmitteltemperaturen als Linien in einem Diagramm
 # mit dem Datum auf der X-Achse und der Temperatur auf der Y-Achse visualisiert.
 # 'geom_line' ist dabei das Argument, welches Liniendiagramme in ggplot generiert.
 library(ggplot2)
 ggplot(hobos_daily_reshaped, aes(x = Datetime, y = Temperatur)) +
   geom_line(aes(color= Standort)) +
   ylab('Temperatur [°C]') + xlab('Datum') +
   theme(axis.text.x=element_text(angle=60, hjust=1)) +
   ggtitle('Tagesmitteltemperaturen (2016)') +
   theme(plot.title = element_text(hjust = 0.5))

 # Exemplarisch wird nun die relative Luftfeuchte im Tagesmittel dargestellt.
 # Merke, es haben sich lediglich die in 'y' und 'color' angesteuerten
 # Variablen geändert, um eine neue Variable darzustellen.
 library(ggplot2)
 ggplot(hobos_daily_reshaped, aes(x = Datetime, y = RH)) +
   geom_line(aes(color= Standort2)) +
   ylab('Relative Feuchte [%]') + xlab('Datum') +
   theme(axis.text.x=element_text(angle=60, hjust=1)) +
   ggtitle('Tagesmittel der relativen Luftfeuchte (2016)') +
   theme(plot.title = element_text(hjust = 0.5))

 # Abschließende Hinweise --------------------------------------------------

 # Damit sind nun die wichtigsten Schritte zum Datenhandling in R und
 # grundlegende Grafiken generiert.
 # Bei Fragen ziehe gerne das Internet zurate oder wende Dich im
 # RLab-Blog mit der Kommentarfunktion an andere User oder das
 # RLab-Team (https://rlab.blogs.uni-hamburg.de/)

 # ENDE

Du hast Fragen zu diesem Code? Nutze gerne die Kommentar-Funktion!

Trennzeichen in .csv-Dateien - Ursache für fehlerhaften Datenimport
.csv-Dateien besitzen (leider) nicht immer identische Zeichen für die Trennung von Spalten und Dezimalstellen. Häufig werden Spalten durch Semikola getrennt, üblich sind aber auch Kommata oder andere Zeichen. Die Dezimalstellen werden meistens mit Punkt oder Komma abgetrennt, aber auch hier sind prinzipiell andere Zeichen möglich.

Bei der Anzeige in MS-Excel ist die Identifikation der Trennzeichen nicht eindeutig möglich. Wir empfehlen, eine .csv-Datei dazu mit einem Texteditor (z.B. Notepad oder Editor) zu öffnen. In der Anzeige dieser Programme sind die Zeichen eindeutig ablesbar und können per copy&paste in den R-Code übernommen werden: Beim Einlesen einer .csv-Datei in R werden die Trennzeichen in der Funktion read.table() mit den Argumenten sep = und dec = angegeben. Für die Beispieldateien sind in den Code-Beispielen die passenden Trennzeichen angegeben. Bei anderen Dateien ist es ggf. erforderlich, sie anzupassen.

Ausführliche Infos zu .csv-Dateien findest Du im Digitalen Skript ,,Import von Daten in R” auf der Seite csv-Dateien im Detail.

Besonders häufig treten Probleme beim Einlesen von HOBO-Daten auf, da beim Export von Daten aus der HOBO-Software unterschiedliche Trennzeichen angegeben werden können. Zusätzlich können sich die Zeichen bei einer weiteren Bearbeitung in MS-Excel unbemerkt ändern. Daher: Am besten kurz per Texteditor checken und Trennzeichen im Code entsprechend anpassen.

Zusammenfassung

In diesem Digitalen Skript wurde ausführlich kommentierter R-Code zur Verfügung gestellt. Wenn Du Grundlagenwissen zur Nutzung von R erhalten möchtest, geht das mit vielen anderen Lernressourcen des RLab. Dazu zählen unter anderem die folgenden Digitalen Skripte:

Auch die meisten RLab-Kurse, mit denen R direkt in R gelernt werden kann, eignen sich zum Einstieg sehr gut:

RLab-Kurse

Labordatenanalyse mit R

elearn.js Template

Labordatenanalyse mit R

Inhaltsverzeichnis

Ziel

Voraussetzungen

R-Code zur Analyse der Parameter der Bodenlösung

Zum Ausführen des Codes erforderliche Beispieldaten

Code

R-Code zur Analyse der bodenchemischen Parameter eines Leitprofils

Zum Ausführen des Codes erforderliche Beispieldaten

Code

R-Code zur Analyse der bodenphysikalischen Parameter eines Leitprofils

Zum Ausführen des Codes erforderliche Beispieldaten

Code

R-Code zur Darstellung der Korngrößenzusammensetzung eines Leitprofils

Zum Ausführen des Codes erforderliche Beispieldaten

Code

R-Code zur Analyse der atmosphärischen Deposition

Zum Ausführen des Codes erforderliche Beispieldaten

Code

R-Code zur Analyse der Klimadaten von HOBO-Loggern

Zum Ausführen des Codes erforderliche Beispieldaten

Code

Zusammenfassung