How to extract parameters from a non-stationary generalized Pareto (GP) model in R?

Question

I am now using the extRemes package to build a non-stationary GP model, and I find it difficult to extract the parameters.

Non-stationary scale parameter

library(extRemes)
data(Fort)
fit1 <- fevd(Prec, Fort, threshold=0.395,
        scale.fun=~sin(2 * pi * (year - 1900)/365.25) + cos(2 * pi * (year - 1900)/365.25),
        type="GP", use.phi=TRUE, verbose=TRUE)

According to the fevd help page, log(scale(y)) = phi(y) = phi0 + phi1 * g1(y) + phi2 * g2(y) + ...

Now, we have phi0, phi1, and phi2 from the results, but what are g1(y) and g2(y) in the above function?

Also, how can we understand the scale.fun in fit1? What does scale.fun=~sin(2 * pi * (year - 1900)/365.25) + cos(2 * pi * (year - 1900)/365.25) stand for? For example, if we use scale.fun=~Fort$year, we assume that year has a linear influence on the scale parameter.

Non-stationary threshold

fit2 <- fevd(Prec, Fort, threshold=0.475, threshold.fun=~I(-0.15 * cos(2 * pi * month / 12)),
        type="GP", verbose=TRUE)

From fit2, how can we compute the changing threshold values based on threshold.fun = ~I(-0.15 * cos(2 * pi * month/12))? Thanks for any help.

Answer 1

The g1(y) etc. are the functions you gave to fevd via the scale.fun argument. So, phi0 is like an intercept term, ph1 is the coefficient for sin(2 * pi * (year - 1900)/365.25), etc. Because you used use.phi = TRUE, your model for the scale parameter is estimated to be: log( scale ) = -0.84 - 0.23 * sin(2 * pi * (year - 1900)/365.25) - 0.25 * cos(2 * pi * (year - 1900)/365.25), and your estimated shape parameter is about 0.21 (so heavy-tailed).