How can I compute the probability $P(X > value)$?

[nucflash]

Hi, thank you for the great effort to bring Extreme Value Analysis to Python! It is extremely useful to have something like this in Python as it plays so well with all other data science tools that are around.

I'm following Fawcett's tutorial/lectures, and in one of his example websites he shows an example of "Choose a Wave height to find the probability of exceeding it every year", this translates to calculating $P(X > value)$, where $value$ is the height of a wave in feet.

I'm trying to figure out how to achieve the same in pyextremes. I looked at the tutorials and the code but I can't figure out a way, but I'm sure that one must exist. Can you please show me how?

Thanks again for the great work!

[nucflash]

I think I found it, I think it is:

1 - model.model.cdf(value)

CDF returns $P(X < value)$ so if I want the opposite I ask for $1 - P(X < value) = P(X > value)$

My next question is:

Is it possible to calculate this probability within a specific return period? The question I want to answer is "What is the probability of achieving value in the first year, the second year, etc.".

Thanks!

[georgebv]

You are correct. Also don't forget to factor in rate of extremes when converting return period into probability:

pyextremes/src/pyextremes/eva.py

Lines 1230 to 1259 in 3706ff8

	# Calculate rate of extreme events
	# as number of extreme events per `return_period_size`
	if self.extremes_method == "BM":
	extremes_rate = return_period_size / self.extremes_kwargs["block_size"]
	elif self.extremes_method == "POT":
	n_periods = (self.data.index[-1] - self.data.index[0]) / return_period_size
	extremes_rate = len(self.extremes) / n_periods
	else:
	raise AssertionError
	# Convert 'return_period' to ndarray
	return_period = np.asarray(a=return_period, dtype=np.float64).copy()
	if return_period.ndim == 0:
	return_period = return_period[np.newaxis]
	if return_period.ndim != 1:
	raise ValueError(
	f"invalid shape in {return_period.shape} "
	f"for the 'return_period' argument, must be 1D array"
	)
	# Calculate exceedance probability
	exceedance_probability = 1 / return_period / extremes_rate
	# Calculate return values
	return tuple(
	self.extremes_transformer.transform(value)
	for value in self.model.get_return_value(
	exceedance_probability=exceedance_probability, alpha=alpha, **kwargs
	)
	)

Regarding your second question, extreme events are assumed to be random IID (Independent and Identically Distributed) so it doesn't matter which year you look at. The question you can answer though is What is the probability of wave height exceeding X at least once in N years?. That would be inverse of probability of not observing it N years in a row.

[nucflash]

Thank you George!

I see, I should divide the wave height by the extreme_rates. I use BM, and I use the same block_size as the return_period_size, I think I shouldn't be affected too much by it, but it is very good to know about it!

With this in mind and for my first question, "What is the probability of wave height exceeding X every year", I write:

wave_height = 7  # in feet
return_period = 0.01  # ~1yr in logscale
extreme_rates = return_period_size / self.extremes_kwargs["block_size"]   # 365/365 = 1
1 - model.model.cdf(wave_height / extreme_rates)

And thank you for the pointer for my second question: "What is the probability of wave height exceeding X at least once in N years?". It feels like I should work along the lines of $1 - \prod_{n} P(X < x)$.