Profile Likelihood Confidence Intervals for the Quantiles of the Maximum on a given Period

Source: R/PLODE.R

Compute the lower and upper limits of a confidence interval on the quantiles \(q_{M^\star}(p)\) of the maximum \(M^\star\) on a given period of time. This computation is made for a number of probabilities \(p\).

quantMaxPLODE(
  object,
  date,
  probIni = 0.7,
  level = 0.95,
  out = c("data.frame", "array"),
  trace = 1
)

Arguments

object

An object with class "TVGEV".

date

A vector that can be coerced to the class "Date" giving the beginnings of the years in the period.

probIni

The probability that will be used for the initialisation of the tangential ODE.

level

The confidence level \(1 - \alpha\).

out

The type of output wanted. The choice "array" leads to an array with its dimensions corresponding to the probability \(p\), the type of result (lower/upper bound, quantile) and the confidence level. With the choice "data.frame" we get an object inheriting from this class, actually a quantMax.TVGEV object for which some methods like autoplot are available.

trace

Integer level of verbosity.

Value

An object with class "quantMax.TVGEV" inheriting from "data.frame"

Details

For a given probability \(p\), the profile likelihood confidence interval for the quantile \(\eta(\boldsymbol{\psi},\,p):=q_{M^\star}(p;\,\boldsymbol{\psi})\) can be obtained by minimising and maximising \(\eta(\boldsymbol{\psi},\,p)\) w.r.t. \(\boldsymbol{\psi}\) under the constraint \(\ell(\boldsymbol{\psi}) \geq \ell_{\max} - \delta\) where \(\ell(\boldsymbol{\psi})\) is the log-likelihood function, and \(\delta := q_{\chi^2(1)}(1 -\alpha) /2\). By differentiating the Lagrangian for this problem w.r.t. the probability \(p\), we get an ODE for the vector \(\boldsymbol{\psi}(p)\). This ODE may be called a tangential ODE inasmuch as it describes the evolution of the vector \(\boldsymbol{\psi}(p)\) on the likelihood surface defined in the parameter space by the equation \(\ell(\boldsymbol{\psi}) = \ell_{\max} - \delta\). The initial ODE can be found either by solving a constrained optimisation problem or by solving another ODE, the radial ODE that describes the evolution of a point moving from an origin at the ML estimate \(\hat{\boldsymbol{\psi}}\) to a destination located on the surface likelihood.

Caution

This function is experimental.

References

Yves Deville (2024) "Profile Likelihood via Optimisation and Differential Equations" arXiv doi:10.48550/arXiv.2404.02774 .

Examples

if (FALSE) { # \dontrun{
   library(deSolve)
   df <- within(TXMax_Dijon, Date <- as.Date(sprintf("%4d-01-01", Year)))
   ## fit a TVGEV model. Only the location parameter is TV.
   object <- TVGEV(data = df, response = "TXMax", date = "Date",
                   design = breaksX(date = Date, breaks = "1970-01-01", degree = 1),
                   estim = "nloptr",
                   loc = ~ t1 + t1_1970)
   newDate <- as.Date(paste0(2025:2034, "-01-01"))
   st1 <- system.time(
      plODE <- quantMaxPLODE(object, date = newDate, level = c(0.70, 0.95))
   )
   st2 <- system.time(
       pl <- quantMax(object, date = newDate, conf = "proflik", trace = 1,
                  level = c(0.70, 0.95))
   )
   g <- autoplot(pl)
   g <- g + geom_line(aes(x = 1 - Prob, y = L, group = Level), data = plODE) +
      geom_line(aes(x = 1 - Prob, y = U, group = Level), data = plODE)
   g
} # }