Retrieving lower_ci, upper_ci (p or years) for return periods when providing return values

[rosapietroiusti]

Dear all,

I am wondering if there is a way to retrieve a lower_ci and upper_ci for return periods, when these are calculated based on return values (e.g., based on our model.extremes).

That is, if there is something similar to model.get_summary() but that would provide return periods and their CI (in years or as probabilities) as a function of return values instead of the other way around, and while being sure that modelling parameters applied are the same as in our model object, whereas if I understand correctly get_return_periods gives plotting position return period and not best estimate and CI based on model fit.

Many thanks !

Example output from summary = model.get_summary() I would like to invert, by providing the return values instead of the return period:

summary
Out[428]: 
               return value  lower ci  upper ci
return period                                  
2.0                0.280462  0.215774  0.346861
5.0                0.555262  0.449663  0.658352
10.0               0.737203  0.595743  0.869465
25.0               0.967087  0.776496  1.134802
50.0               1.137628  0.909908  1.334365
100.0              1.306909  1.044441  1.530192
250.0              1.529796  1.221623  1.786949
500.0              1.698092  1.355618  1.983863
1000.0             1.866268  1.490158  2.180612

Example output from get_return_periods with plotting positions outputted:

return_periods = get_return_periods(
		ts=data,
		extremes=model.extremes,
		extremes_method=model.extremes_method,
		extremes_type="high",
		block_size="365.2425D",
		return_period_size="365.2425D",
		plotting_position="weibull",
	)

return_periods
Out[429]: 
            extreme  exceedance probability  return period
date                                                       
1949-01-03    -0.065                0.945946       1.057143
1950-12-17     0.120                0.729730       1.370370
1952-01-01     0.580                0.175676       5.692308
1952-09-29     0.815                0.067568      14.800000
1953-01-08     0.100                0.770270       1.298246
             ...                     ...            ...
2017-12-05    -0.160                0.972973       1.027778
2018-06-12     0.417                0.324324       3.083333
2020-01-01     0.488                0.229730       4.352941
2020-05-16     1.392                0.027027      37.000000
2021-01-01     0.608                0.141892       7.047619

[73 rows x 3 columns]

[georgebv]

That's an interesting question. Should be absolutely doable. We would need to implement a method of the model instances - it would take return value and return upper tail probability using survival function of a given distribution (see this).

Confidence bounds would be much trickier to get - if you use direct inverse method (SF) then those would be different bounds than what you get normally. I don't know how that would work.

I am currently low on time I can spend on this package development, mainly because it already goes above and beyond what I personally needed from it. I would gladly accept contributions and could help with code review, though, if you (or anyone else reading this) is willing to take a shot at this problem.

[rosapietroiusti]

Hi, thanks for the reply! Alright, I see. I will keep following the package to be up to date if this gets implemented by someone, and if I am working on this in the future I will let you know.
Best!