<html>
<head>
<meta http-equiv="Content-Type" content="text/html; charset=UTF-8">
</head>
<body>
<p>Hello everyone,<br>
</p>
<p><font face="Helvetica, Arial, sans-serif">I was hoping someone
would be able to help my with my script for extracting the 99th
percentiles (peak-over-threshold method) from my data, fit to a
GEV distribution using the maximum likelihood method. and than
from this distribution i want to obtain the uncertainty
estimates of the extreme value distribution (i.e the 5% and the
95% confidence bounds).</font></p>
<p><font face="Helvetica, Arial, sans-serif">Can someone please
advice me?</font></p>
<p><font face="Helvetica, Arial, sans-serif">Really appreciate this!</font></p>
<p><font face="Helvetica, Arial, sans-serif">Sri</font></p>
<p><font face="Helvetica, Arial, sans-serif"><br>
</font> </p>
<p>
<style type="text/css">p { margin-bottom: 0.1in; line-height: 120%; }a:link { }</style></p>
<div class="moz-cite-prefix">On 11/26/19 3:40 PM, Sri.durgesh
Nandini-Weiss via ncl-talk wrote:<br>
</div>
<blockquote type="cite"
cite="mid:001c0893-cef2-76c0-2ad3-caf63b82f86d@uni-hamburg.de">
<meta http-equiv="content-type" content="text/html; charset=UTF-8">
<p> </p>
<p style="margin-bottom: 0.2in; line-height: 100%"><font
face="Accanthis ADF Std"><font style="font-size: 13pt"
size="3">Hello everyone,</font></font></p>
<p style="margin-bottom: 0.2in; line-height: 100%"><font
face="Accanthis ADF Std"><font style="font-size: 13pt"
size="3"><br>
Im trying to calculate extreme value analysis and plot the
CDF of the long term changes in extreme sea level by 2100
using the 99<sup>th</sup> percentiles and the GEV
distribution method.</font></font></p>
<p style="margin-bottom: 0.2in; line-height: 100%"> <font
face="Accanthis ADF Std"><font style="font-size: 13pt"
size="3">According to previous studies (<a
href="https://www.frontiersin.org/articles/10.3389/feart.2018.00216/full#B33"
moz-do-not-send="true">Vogel et al., 1993</a>; <a
href="https://www.frontiersin.org/articles/10.3389/feart.2018.00216/full#B16"
moz-do-not-send="true">Huang et al., 2008</a>; <a
href="https://www.frontiersin.org/articles/10.3389/feart.2018.00216/full#B10"
moz-do-not-send="true">Feng and Jiang, 2015</a>) the GEV
distribution was used to get the return levels of the
extreme sea levels. Percentile analysis method is widely
used to assess the extreme sea level changes (<a
href="https://www.frontiersin.org/articles/10.3389/feart.2018.00216/full#B28"
moz-do-not-send="true">Menéndez and Woodworth, 2010</a>; <a
href="https://www.frontiersin.org/articles/10.3389/feart.2018.00216/full#B11"
moz-do-not-send="true">Feng et al., 2015</a>; <a
href="https://www.frontiersin.org/articles/10.3389/feart.2018.00216/full#B26"
moz-do-not-send="true">Marcos and Woodworth, 2017</a>).<br>
<br>
</font></font><span style="display: inline-block; border:
none; padding: 0in"><font style="font-size: 13pt" size="3"><font
face="Accanthis ADF Std">I modified NCL example extvalv_6
</font></font><font style="font-size: 13pt" size="3"><font
face="Accanthis ADF Std">but it uses the Annual maxima
series and I have used the highest percentile method
(threshold method): see my script attached below. My steps
are:<br>
</font></font></span></p>
<p style="margin-bottom: 0.2in; line-height: 100%"><span
style="display: inline-block; border: none; padding: 0in"><font
style="font-size: 13pt" size="3"><font face="Accanthis ADF
Std">1) Calculate and extract the 99th percentile value of
the anomaly data using stat_dispersion</font></font></span><span
style="display: inline-block; border: none; padding: 0in"><font
style="font-size: 13pt" size="3"><font face="Accanthis ADF
Std">2) Basic statistics of the original sample</font></font></span></p>
<p style="margin-bottom: 0.2in; line-height: 100%"><span
style="display: inline-block; border: none; padding: 0in"></span><font
style="font-size: 13pt" size="3"><font face="Accanthis ADF
Std">3<span style="display: inline-block; border: none;
padding: 0in">) Using extval_mlegev and extval_gev</span></font></font><font
face="Accanthis ADF Std"><font style="font-size: 13pt"
size="3"><br>
</font></font></p>
<p style="margin-bottom: 0.2in; line-height: 100%"><font
face="Accanthis ADF Std"><font style="font-size: 13pt"
size="3">4) Estimate GEV distribution parameters</font></font><font
face="Accanthis ADF Std"><font style="font-size: 13pt"
size="3">: Can I use the GEV distribution to get the
return levels of the extreme sea levels or use this?</font></font>
</p>
<pre class="western" style="margin-bottom: 0.2in"> <font face="Accanthis ADF Std"><font style="font-size: 13pt" size="3">Pr = 0.99 ; user specified </font></font>
<font face="Accanthis ADF Std"><font style="font-size: 13pt" size="3">Nr = extval_return_period(anom1_1d, Pr) ; Nr=?? years at 0.99 probability</font></font>
<font face="Accanthis ADF Std"><font style="font-size: 13pt" size="3">print(Nr) </font></font>
<font face="Accanthis ADF Std"><font style="font-size: 13pt" size="3">The risks associated with extreme sea levels can be assessed from the estimates of return levels and return periods. The return period of extreme sea level is defined as the sea level statistically expected to be equaled or exceeded every specific year. </font></font>
<font face="Accanthis ADF Std"><font style="font-size: 13pt" size="3">5) Plot the CDF and the return value with uncertainty estimates</font></font></pre>
<p>At the moment i have errors when calculating my return value:
"(23998) extval_return_period: Tr must be > 0: Tr=4.20953"
and so on.</p>
<p>Can someone please help me with this, or use my script on a
sample temperature data?</p>
<p>f = addfile ("anom_hist_zo_slr.nc", "r")<br>
hist_anom = f->hist_anom
;(time|240,ens|100,lat|90,lon|45)<br>
printVarSummary(hist_anom)<br>
<br>
f1 = addfile ("anom_rcp85_zo_slr.nc", "r")<br>
rcp45_anom = f1->rcp85_anom
;(time|240,ens|100,lat|90,lon|45)<br>
printVarSummary(rcp45_anom)<br>
<br>
</p>
<p><br>
;==================================================================<br>
print ("Data read in - start Extreme value calculations")<br>
;==================================================================<br>
; Region domain lat and lon<br>
;==================================================================<br>
<br>
hist_anom1=hist_anom(:,:,{18},{178}) ;(ntim,
ens) <br>
printVarSummary(hist_anom1)<br>
<br>
rcp45_anom1=rcp45_anom(:,:,{18},{178}) ;(ntim,
ens)<br>
printVarSummary(rcp45_anom1)<br>
<br>
;=================================================================<br>
<br>
anom1_1d = ndtooned(hist_anom1) <br>
printVarSummary(anom1_1d)<br>
anom2_1d = ndtooned(rcp45_anom1)<br>
printVarSummary(anom2_1d) <br>
<br>
;=================================================================<br>
; Use stat_dispersion over each region: to calculate 99
percentiles <br>
; (a) For each time period<br>
; (b) For each statistic,<br>
; (c) For each gridpoint? or for all gridpoints? [45*90] sort
into ascending order and pick the highest 1%<br>
; STATX(22 to 29) added for v6.1.0<br>
; STATX(22 to 27) require
more than 1000 observations<br>
; ; Otherwise, they are set
to _FillValue<br>
;<br>
; Lower 0.1%="+STATX(22) ; lower 0.1% => 0.1th
percentile<br>
; Lower 1.0%="+STATX(23) ; lower 1.0% => 1th<br>
; Lower 5.0%="+STATX(24) ; lower 5.0% => 5th<br>
; Upper 5.0%="+STATX(25) ; upper 5.0% => 95th
percentile<br>
; Upper 1.0%="+STATX(26) ; upper 1.0% => 99th <br>
; Upper 0.1%="+STATX(27) ; Upper 0.1% => 99.9th
<br>
<br>
;=================================================================<br>
<br>
hist_perc = stat_dispersion(anom1_1d, True)<br>
printVarSummary(hist_perc) <br>
hist_99= hist_perc(26)<br>
<br>
rcp45_perc = stat_dispersion(anom2_1d, True)<br>
printVarSummary(rcp45_perc) <br>
rcp45_99= hist_perc(26)</p>
<p>; In extreme analysis, a decision has to be made regarding the
definition of the extreme events,either in the form of the
number of extremes (when a GEV distribution is used) the method
of maximum likelihood is used to estimate GEV parameters
following Kharin and Zwiers [2005].</p>
<p> <br>
vals = extval_mlegev ( anom1_1d, 0, False) ; dims=0<br>
print(vals) ; vals(6)<br>
<br>
center = vals(0) ; extract
for clarity<br>
scale = vals(1)<br>
shape = vals(2)</p>
<p><br>
;================================================================================<br>
; extval_return_period:Calculates the period of an event (eg,
flood, heat wave, drought) occurring given an average event
recurrence interval and specified probability level. <br>
; The risks associated with extreme sea levels can be assessed
from the estimates of return levels and return periods. The
return period of extreme sea level is defined as the sea level
statistically expected to be equaled or exceeded every specific
year. <br>
;================================================================================<br>
<br>
Pr = 0.99 ; user specified <br>
Nr = extval_return_period(anom1_1d, Pr) ; Nr=?? years at 0.99
probability<br>
print(Nr) <br>
</p>
<p>;================================================================================
pdfcdf = extval_gev( anom1_1d, shape, scale, center,0)<br>
pdf = pdfcdf[0] ; extract from
list for convenience<br>
cdf = pdfcdf[1]<br>
delete(pdfcdf)<br>
<br>
printVarSummary(pdf)<br>
printMinMax(pdf,1)<br>
<br>
printVarSummary(cdf)<br>
printMinMax(cdf,1)<br>
<br>
;================================================================================<br>
; Plot the PDF and CDF<br>
;================================================================================<br>
<br>
nscl = dimsizes( scale )<br>
labels = new( nscl , "string")<br>
plot = new( 2 , "graphic")<br>
<br>
wks = gsn_open_wks ("png","extval") <br>
<br>
res = True <br>
res@gsnDraw = False<br>
res@gsnFrame = False<br>
<br>
res@xyLineThicknessF = 2.0<br>
res@xyMonoDashPattern = True<br>
res@xyLineColors = (/"green"/)<br>
res@trYMinF =-0.02 ; min value on y-axis<br>
res@trYMaxF = 0.5 ; max value on y-axis<br>
<br>
res@pmLegendDisplayMode = "Always" ; turn on legend<br>
res@pmLegendSide = "Top" ; Change location of
<br>
res@pmLegendParallelPosF = 0.675 ; move units right<br>
res@pmLegendOrthogonalPosF = -0.30 ; more neg = down<br>
res@pmLegendWidthF = 0.12 ; Change width and<br>
res@pmLegendHeightF = 0.150 ; height of legend.<br>
res@lgLabelFontHeightF = 0.0175 ; change font height<br>
res@lgPerimOn = True <br>
res@lgItemOrder = ispan(nscl-1,0,1)<br>
res@xyExplicitLegendLabels = " shp="+sprintf("%3.1f", shape)+
\<br>
"; scl="+sprintf("%3.1f",
scale)+"; ctr="+sprintf("%3.1f", center) <br>
it = dim_pqsort_n(anom1_1d, 1, 0) ; for plot
reasons make 'x' is ascending<br>
printVarSummary(it)<br>
<br>
res@tiYAxisString = "PDF: extval_gev"<br>
plot(0)= gsn_csm_xy (wks,anom1_1d(it),pdf(it),res) ; create
plot<br>
<br>
res@trYMinF = 0.0 <br>
res@trYMaxF = 1.0 <br>
res@pmLegendOrthogonalPosF = -0.95 ; more neg = down<br>
res@tiYAxisString = "CDF: extval_gev"<br>
plot(1)= gsn_csm_xy (wks,anom1_1d(it),cdf(it),res) ; create
plot<br>
<br>
resP = True ; modify
the panel plot<br>
resP@gsnMaximize = True ; ps, eps,
pdf<br>
resP@gsnPanelMainString = "GEV: PDF and CDF"<br>
gsn_panel(wks,plot,(/1,2/),resP) <br>
<br>
<br>
</p>
<br>
<p>
<style type="text/css">pre.cjk { font-family: "Courier New", monospace; }p { margin-bottom: 0.1in; line-height: 120%; }a:link { }</style></p>
<pre class="moz-signature" cols="72">Dr. Sri Nandini-Weiß
Research Scientist
Universität Hamburg
Center for Earth System Research and Sustainability (CEN)
Cluster of Excellence 'Climate, Climatic Change, and Society' (CLICCS)
Bundesstrasse 53, 20146 Hamburg
Tel: +49 (0) 40 42838 7472</pre>
<br>
<fieldset class="mimeAttachmentHeader"></fieldset>
<pre class="moz-quote-pre" wrap="">_______________________________________________
ncl-talk mailing list
<a class="moz-txt-link-abbreviated" href="mailto:ncl-talk@ucar.edu">ncl-talk@ucar.edu</a>
List instructions, subscriber options, unsubscribe:
<a class="moz-txt-link-freetext" href="http://mailman.ucar.edu/mailman/listinfo/ncl-talk">http://mailman.ucar.edu/mailman/listinfo/ncl-talk</a></pre>
</blockquote>
<pre class="moz-signature" cols="72">--
Dr. Sri Nandini-Weiß
Research Scientist
Universität Hamburg
Center for Earth System Research and Sustainability (CEN)
Cluster of Excellence 'Climate, Climatic Change, and Society' (CLICCS)
Bundesstrasse 53, 20146 Hamburg
Tel: +49 (0) 40 42838 7472</pre>
</body>
</html>