1
Umweltsystemmodellierung (Subject) / Termin 4 Hausaufgaben (Lesson)
- 1. Univariate Sensitivitätsanalyse Vorbereitung: Die Datei „Zweispeicher2.R“ enthält den Code für das Speichermodell mit zwei Speichern aus der Hausaufgabe der vergangenen Woche. Die in dieser Datei enthaltene Funktion zweispeicher2 unterscheidet sich von zweispeicher dadurch, dass sie die Werte aller Anfangsbedingungen und Parameter (sa0, sb0, ka, kb, u) als einen gemeinsamen Vektor x einliest. Diese Eigenschaft erleichtert es, in einer Schleife die Werte der einzelnen Parameter nacheinander zu verändern. Lesen Sie die Funktion zweispeicher2 (über source) und den Datensatz „Carnation_Creek_Daten.csv“ ein. Extrahieren Sie aus den Daten einen Vektor „niedvec“, der die Niederschlagswerte für den Zeitraum vom 1. 8. 2008 bis 31. 7. 2009 enthält. a) Führen Sie eine Sensitivitätsanalyse für den Output „mittlerer Gesamt-Abfluss“ durch. Benutzen Sie für die Sensitivitätsanalyse folgende Referenzwerte der Parameter und Anfangsbedingungen: sa0=5, sb0=10, ka=0.9, kb=0.5, u=0.9. Berechnen Sie den Mittelwert des Gesamt-Abflusses für diese Referenzwerte und den Niederschlagsvektor „niedvec“. Ändern Sie nacheinander jeweils einen dieser Werte um 5% nach oben und berechnen Sie den mittleren Gesamt-Abfluss mit den geänderten Werten. Berechnen Sie jeweils den Sensitivitätsindex des mittleren Abflusses für die Änderungen nach der folgenden Formel: S=Aa−Ar/Ar//Pa−Pr/Pr S: Sensitivitätsindex A: Mittelwert des Gesamt-Abflusses P: Parameterwert Index r: Referenzwert Index a: veränderter Wert („altered“) b) Berechnen Sie auch die Sensitivitätsindizes für eine 5%ige Änderung nach unten. c) Auf welchen Parameter reagiert der mittlere Abfluss am stärksten? Unterscheidet sich die Sensitivität zwischen Änderung nach oben und Änderung nach unten? zweispeicher2source(zweispeicher2)carnation zeit<-as.Date(carnation$zeit, format="%d.%m.%Y")begin<-which(zeit=="2008-08-01")end<-which(zeit=="2009-07-31")niedervec<-carnation$nied[begin:end]niedervec x_ref<-c(sa0=5, sb0=10, ka=0.9, kb=0.5, u=0.9)output_ref<-mean(zweispeicher2(niedervec, x_ref)$Abfluss_ges)S_index_oben<-numeric(length(x_ref)) for (i in 1:length(x_ref)) {x_a<-x_refx_a[i]<-x_ref[i]*1.05output_a<-mean(zweispeicher2(niedervec,x_a)$Abfluss_ges)S_index_oben[i]<-((output_a-output_ref)/output_ref)/((x_a[i]-x_ref[i])/x_ref[i])}S_index_oben#S_index_untenS_index_unten<-numeric(length(x_ref))for (i in 1:length(x_ref)) {x_a<-x_refx_a[i]<-x_ref[i]*0.95output_a<-mean(zweispeicher2(niedervec,x_a)$Abfluss_ges)S_index_unten[i]<-((output_a-output_ref)/output_ref)/((x_a[i]-x_ref[i])/x_ref[i])}S_index_unten # Modifizierte Funktion des Zwei-Speicher-Modells einlesen:# (Diese Funktion liest die Parameterwerte und Anfangsbedingungen nicht als einzelne Werte,# sondern gemeinsam als Vektor x ein. Das erleichtert es, für die lokale # Sensitivitätsanalyse eine Schleife über die Parameter laufen zu lassen.)source("Zweispeicher2.R") # Carnation Creek - Daten einlesen:carnation <- read.csv("Carnation_Creek_Daten.csv")tage <- as.Date(carnation$zeit, format = "%d.%m.%Y") # Sage R, welches Format die Zeitangaben haben # Vektor mit Niederschlagsdaten:begin <- which(tage=="2008-08-01")end <- which(tage=="2009-07-31")niedvec <- carnation$nied[begin:end] # a) Ãndere nacheinander die Parameterwerte um 5% nach oben,# berechne den Sensitivitätsindex für den mittleren Gesamt-Abfluss:x_ref <- c(sa0=5, sb0=10, ka=0.9, kb=0.5, u=0.9) # Referenzwerte der Parameter und Anfangsbedingungenoutput_ref <- mean(zweispeicher2(niedvec, x_ref)$Abfluss_ges) # Mittlerer Abfluss mit den Referenz-ParameterwertenS_index_oben <- numeric(length(x_ref)) # Vektor, der die Werte des Sensitivitäts-Index für die Ãnderung der einzelnen Parameter aufnehmen sollnames(S_index_oben) <- names(x_ref)for (i in 1:length(x_ref)){ # Schleife über alle Parameterwertex_a <- x_refx_a[i] <- x_ref[i]*1.05output_a <- mean(zweispeicher2(niedvec, x_a)$Abfluss_ges) # Mittlerer Abfluss mit geändertem Wert von Parameter x_a[i]S_index_oben[i] <- ((output_a-output_ref)/output_ref)/((x_a[i]-x_ref[i])/x_ref[i]) # Sensitivitätsindex} # b) Dasselbe für eine 5%ige Ãnderung nach unten:S_index_unten <- numeric(length(x_ref)) # Vektor, der die Werte des Sensitivitäts-Index für die Ãnderung der einzelnen Parameter aufnehmen sollnames(S_index_unten) <- names(x_ref)for (i in 1:length(x_ref)){x_a <- x_refx_a[i] <- x_ref[i]*0.95output_a <- mean(zweispeicher2(input=niedvec, x=x_a)$Abfluss_ges) # Maximaler Abfluss mit geändertem Wert von Parameter x_a[i]S_index_unten[i] <- ((output_a-output_ref)/output_ref)/((x_a[i]-x_ref[i])/x_ref[i])} S_index_oben# sa0 sb0 ka kb u # 1.834917e-03 3.669835e-03 3.751378e-04 2.931989e-05 1.904116e-04 S_index_unten# sa0 sb0 ka kb u # 1.834917e-03 3.669835e-03 3.858391e-04 2.976862e-05 1.904116e-04 # c) Schlussfolgerung: Der mittlere Abfluss reagiert am stärksten auf Ãnderungen von sb0,# gefolgt von sa0. Die Sensitivität für Ãnderung nach oben oder unten unterscheidet sich nur geringfügig.. # Allerdings gelten diese Aussagen nur für die gewählten Referenzwerte der Parameter.