[ncl-talk] GEV assistance needed

Tue Nov 23 13:07:32 MST 2021

As far as I can tell, the NCL-computed GEV parameters based on MLE do not match prior studies' results.
I say that with the huge caveat that I am not a GEV expert.

For example, based on the annual max values for 1950-2014 in Table 1 of this paper:
https://www.researchgate.net/publication/288991817_Predicted_Return_Level_of_Annual_Maximum_Earthquake_in_Northern_California/link/56881a6208aebccc4e155539/download
the location, scale, and shape parameters are 5.62, 0.54, and -0.098, respectively, in Table 2.

When I perform this task below with those data values,
yearmax=(/5.3,5.8,7.2,6.1,6.3,5.8,5.8,5.8,5.5,5.5,5.7,6.4,5.5,5.5,5.1,6.0,5.0, \
          5.3,6.0,6.2,5.6,5.8,6.0,5.3,7.2,5.5,6.0,5.4,4.9,6.3,6.2,5.6,5.6,5.4, \
          5.5,5.0,6.0,5.6,5.9,5.7,5.4,5.2,5.2,5.0,5.2,5.7,6.3,5.0,5.8,5.8,6.9, \
          5.9,5.5,6.7,6.6,5.6,6.6,6.8,6.0,6.6,6.3,6.6,7.7,6.2,6.8,6.3,5.6,5.6, \
          5.5,7.3,6.1,6.6,5.9,6.5,6.0,7.0,5.0,5.5,5.4,4.9,6.5,4.7,5.6,5.7,6.8/)
param=extval_mlegev(yearmax,0,False)
print(param)
I get the following, which does not match the results from the study at all.
(0) 3.7994          location
(1) 2.233068      scale
(2) 0.4374145    shape

Is there anyone with experience with GEV and NCL that can offer some insights?

Thanks,
Michael

Michael Notaro
Associate Director
Nelson Institute Center for Climatic Research
University of Wisconsin-Madison
Phone: (608) 261-1503
Email: mnotaro at wisc.edu

________________________________
From: Dennis Shea <shea at ucar.edu>
Sent: Wednesday, November 17, 2021 9:58 AM
To: Michael Notaro <mnotaro at wisc.edu>
Cc: ncl-talk at ucar.edu <ncl-talk at ucar.edu>
Subject: Re: [ncl-talk] GEV assistance needed

Hi Michael,

Not sure I can help on this. I am not very knowledgeable on MLE or GEV

As noted in the documentation, NCL invokes a August 1994 version of:

Algorithm 215: Maximum-Likelihood Estimation of the Parameters of the Generalized Extreme-Value Distribution

J. R. M. Hosking
  Journal of the Royal Statistical Society. Series C (Applied Statistics)
  Vol. 34, No. 3 (1985), pp. 301-310
  URL:  http://www.jstor.org/stable/2347483
  Code: http://ftp.uni-bayreuth.de/math/statlib/apstat/215

My recollection is that we tested this on a few 'book' examples and compared the output to results generated by R's gev.fit function.

====

extval_mlegev<https://www.ncl.ucar.edu/Document/Functions/Built-in/extval_mlegev.shtml>: Example 1 illustrates a comparison.  I believe you have used R.  Perhaps trying your data using R
will provide another (NCL, MATLAB) viewpoint.

====

Again, I am not very knowledgeable on either MLE or GEV so I have no idea what the -1 does within the algorithm.

Regards

D

On Tue, Nov 16, 2021 at 3:24 PM Michael Notaro via ncl-talk <ncl-talk at mailman.ucar.edu<mailto:ncl-talk at mailman.ucar.edu>> wrote:
In the example below, I have 20 values of precipitation in inches.
Each value is the highest daily precipitation for one of 20 years.
I want to estimate the 1-day/100-year GEV precipitation level.
For some reason, my older script included multiplying the shape
parameter by -1.  That gave a GEV value of 65757", definitely wrong.
If I get rid of the *-1 part, the GEV value becomes more reasonable,
1.73", but too low compared to the 20 annual max values.  So I have
two questions.

  1.  Is there something wrong with the script below?
  2.  Is there ever a need to multiply any of the parameters by -1?  I
  3.  did it once before but have no recollection of why I did it.

Thanks, Michael

yearmax=(/1.83,1.69,3.24,1.99,1.69,1.79,2.78,1.70,4.85,2.35,3.37,1.72,3.96,2.35,2.29,3.33,3.97,3.53,3.11,1.69/)
param=extval_mlegev(yearmax,0,False)
param_mle_location=param(0)
param_mle_scale=param(1)
param_mle_shape=param(2)*-1.
coeff=-1.*log(1.-1./100.)
gev=param_mle_scale/param_mle_shape*[coeff^(-1.*param_mle_shape)-1.]+param_mle_location
print("param_mle_location="+param_mle_location+"  param_mle_scale="+param_mle_scale+"  param_mle_shape="+param_mle_shape+"  coeff="+coeff+"  gev= "+gev)

param_mle_location=1.65715  param_mle_scale=0.216018  param_mle_shape=2.98224  coeff=0.0100503  gev= 65757.5

earmax=(/1.83,1.69,3.24,1.99,1.69,1.79,2.78,1.70,4.85,2.35,3.37,1.72,3.96,2.35,2.29,3.33,3.97,3.53,3.11,1.69/)
param=extval_mlegev(yearmax,0,False)
param_mle_location=param(0)
param_mle_scale=param(1)
param_mle_shape=param(2)*-1.
coeff=-1.*log(1.-1./100.)
gev=param_mle_scale/param_mle_shape*[coeff^(-1.*param_mle_shape)-1.]+param_mle_location
print("param_mle_location="+param_mle_location+"  param_mle_scale="+param_mle_scale+"  param_mle_shape="+param_mle_shape+"  coeff="+coeff+"  gev= "+gev)

param_mle_location=1.65715  param_mle_scale=0.216018  param_mle_shape=-2.98224  coeff=0.0100503  gev= 1.72958

Michael Notaro
Associate Director
Nelson Institute Center for Climatic Research
University of Wisconsin-Madison
Phone: (608) 261-1503
Email: mnotaro at wisc.edu<mailto:mnotaro at wisc.edu>
_______________________________________________
ncl-talk mailing list
ncl-talk at mailman.ucar.edu<mailto:ncl-talk at mailman.ucar.edu>
List instructions, subscriber options, unsubscribe:
https://mailman.ucar.edu/mailman/listinfo/ncl-talk
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <https://mailman.ucar.edu/pipermail/ncl-talk/attachments/20211123/5fc43e25/attachment.html>