exceedance_prob.py

def calcualte_exceedance_probability(exceed_values, posterior):
    n = posterior.shape[1]
    ix = []
    probs_m = []
    probs_05 = []
    probs_95 = []
    for i in exceed_values:
        p = ((posterior>i).sum(1)/n)
        p05 = np.quantile(p, 0.01)
        p95 = np.quantile(p, 0.99)
        probs_m.append(p.mean())
        probs_05.append(p05)
        probs_95.append(p95)
        ix.append(i)
    return pd.DataFrame(dict(exceed_value=ix, ep_mean=probs_m, ep_lb=probs_05, ep_ub=probs_95))

exceedance = calcualte_exceedance_probability(list(range(1, 60)), ev_post)
plt.figure(figsize=(20, 12))
plt.plot(exceedance.exceed_value, exceedance.ep_mean, c="r", label="ev")
plt.fill_between(exceedance.exceed_value, exceedance.ep_lb, exceedance.ep_ub, color="b", alpha=.3)
plt.legend()