Learn R Programming

HMMHSMM (version 0.1.0)

exceedanceplotHMMgev: Plot Exceedance Probabilities from GEV-HMM

Description

Computes and plots exceedance probabilities of a given threshold from a fitted Hidden Markov Model (HMM) with generalized extreme value (GEV) state-dependent distributions. Bootstrap confidence intervals are included. Probabilities are aggregated over user-defined time periods.

Usage

exceedanceplotHMMgev(x, HMM, threshold, varcov = NULL, B = 10000, level = 0.95,
                     time_structure = NULL,
                     plot_title = "Exceedance Probabilities Over Time",
                     save_plot = FALSE, filename = NULL,
                     width = 12, height = 8, dpi = 300,
                     verbose = TRUE)

Value

Invisibly returns a list containing:

exceedance_probs

Matrix of exceedance probabilities over time (3 rows: estimate, lower CI, upper CI).

time_info

List containing time axis information and labels.

threshold

The threshold value used.

A time-series plot is also produced, showing aggregated exceedance probabilities with confidence intervals.

Arguments

x

Numeric vector. The observed data sequence.

HMM

A fitted HMM object (output from findmleHMM), with state-dependent GEV parameters loc, scale, shape, and estimated transition probabilities.

threshold

Numeric. The exceedance threshold. Probabilities of exceeding this value are computed.

varcov

Optional variance–covariance matrix of parameter estimates. If NULL, computed internally via HMMVarianceMatrix with obsdist = "gev".

B

Integer. Number of bootstrap replicates used for confidence intervals. Default is 10000.

level

Numeric between 0 and 1. Confidence level for intervals. Default is 0.95.

time_structure

Optional list describing the time scale. Includes unit, observations_per_unit, and optionally start_point, end_point. Supported units: "year", "day", "hour", "week", "month", or "custom" with conversion_factor and unit_name.

plot_title

Character string. Title for the plot. Default is "Exceedance Probabilities Over Time".

save_plot

Logical. If TRUE, the plot is saved to file. Default is FALSE.

filename

Character string or NULL. Filename for saved plot (if save_plot = TRUE). Must be specified when save_plot = TRUE. Default is NULL.

width

Numeric. Width of saved plot in inches. Default is 12.

height

Numeric. Height of saved plot in inches. Default is 8.

dpi

Integer. Resolution of saved plot in dots per inch. Default is 300.

verbose

Logical. If TRUE, progress messages are printed to the console. Default is TRUE.

Author

Aimee Cody

Details

For each state, the exceedance probability of threshold is computed from the fitted GEV parameters. Uncertainty is estimated by parametric bootstrap using draws from a multivariate normal approximation of parameter estimates. Observation-level probabilities are aggregated over time periods defined by time_structure (e.g., years, months, weeks). When verbose = TRUE, progress messages are displayed during computation.

See Also

findmleHMM, for fitting HMMs. generateHMM, for simulating HMM data. globaldecodeHMM, for decoding the most likely state sequence. HMMVarianceMatrix, for variance–covariance estimation.

Examples

Run this code
set.seed(99)
J <- 2
Pi <- matrix(c(0.9, 0.1,
               0.2, 0.8), nrow = J, byrow = TRUE)
delta <- c(0.5, 0.5)
obspar <- list(loc = c(0, 5), scale = c(1, 2), shape = c(0.1, -0.1))

sim <- generateHMM(n = 200, J = J, obsdist = "gev",
                   obspar = obspar, Pi = Pi, delta = delta)

HMM_fit <- findmleHMM(x = sim$x, J = J, obsdist = "gev",
                      obspar = obspar, Pi = Pi, delta = delta)

time_struct <- list(unit = "week", observations_per_unit = 10, start_point = 1)

# \donttest{
result <- exceedanceplotHMMgev(x = sim$x, HMM = HMM_fit,
                               threshold = 20, B = 10000,
                               time_structure = time_struct,
                               plot_title = "GEV-HMM Exceedance Probabilities")

result <- exceedanceplotHMMgev(x = sim$x, HMM = HMM_fit,
                               threshold = 20, B = 10000,
                               time_structure = time_struct,
                               save_plot = TRUE,
                               filename = tempfile(fileext = ".png"),
                               verbose = FALSE)
# }

Run the code above in your browser using DataLab