Stationary

Load the data

Loading the pseudo-observations into an object of type Pseudoensemble:

filename = "../../../test/data/A2020_Analyse_Historique_QMA_SLSO00003.nc"
pensemble = ErrorsInVariablesExtremes.load_discharge_distribution(filename)

6-element Vector{Pseudodata}:
 Pseudodata:
  name: LN24HA
  year: Vector{Int64}[60]
  value:Vector{Distributions.LogNormal{Float64}}[60]

 Pseudodata:
  name: MG24HA
  year: Vector{Int64}[60]
  value:Vector{Distributions.LogNormal{Float64}}[60]

 Pseudodata:
  name: MG24HI
  year: Vector{Int64}[60]
  value:Vector{Distributions.LogNormal{Float64}}[60]

 Pseudodata:
  name: MG24HK
  year: Vector{Int64}[60]
  value:Vector{Distributions.LogNormal{Float64}}[60]

 Pseudodata:
  name: MG24HQ
  year: Vector{Int64}[60]
  value:Vector{Distributions.LogNormal{Float64}}[60]

 Pseudodata:
  name: MG24HS
  year: Vector{Int64}[60]
  value:Vector{Distributions.LogNormal{Float64}}[60]

The vector of years can be extracted:

years = pensemble.value[1].year;

Loading the covariable:

co2data = CSV.read("../../../test/data/RCPdata.csv", DataFrame);
filter!(row->row.Year in pensemble.value[1].year, co2data)
locationcov= Extremes.buildVariables(co2data, [:RCP85])

MCMC

res = gevfitbayes(pensemble, locationcov = locationcov)

Display of annual maximum estimates

# Compute the estimates
ŷ = vec(mean( res.maxima[:, :, 1].value, dims=1))

# Compute the 95% credible interval
ymin = [quantile(vec(res.maxima[:, j, 1].value), .025) for j in 1:length(ŷ)]
ymax= [quantile(vec(res.maxima[:, j, 1].value), .975) for j in 1:length(ŷ)]

df = DataFrame(Year = years, Discharge=ŷ, ymin = ymin, ymax=ymax)

# plot the ponctual estimates and the intervals
set_default_plot_size(12cm, 8cm)
plot(df, x=:Year, y=:Discharge, Geom.line, Geom.point,
    ymin=ymin, ymax=ymax, Geom.ribbon
)

Fit of the GEV law adjusted to the estimates of the annual maxima

set_default_plot_size(21cm ,16cm)
diagnosticplots(res.parameters)