| Title: | Extreme Value Analysis for Circular Data |
|---|---|
| Description: | General functions for performing extreme value analysis on a circular domain as part of the statistical methodology in the paper by Konzen, E., Neves, C., and Jonathan, P. (2021). Modeling nonstationary extremes of storm severity: Comparing parametric and semiparametric inference. Environmetrics, 32(4), e2667. |
| Authors: | Evandro Konzen |
| Maintainer: | Evandro Konzen <[email protected]> |
| License: | GPL-3 |
| Version: | 0.1.1 |
| Built: | 2024-11-28 06:54:54 UTC |
| Source: | CRAN |
Calculate T-year levels for spline ML model
CalcRLsplineML( Data, drc, h, xiBoot, sigBoot, TTs = c(100, 10000), thetaGrid = 1:360, timeRange, thr )CalcRLsplineML( Data, drc, h, xiBoot, sigBoot, TTs = c(100, 10000), thetaGrid = 1:360, timeRange, thr )
Data |
Response variable |
drc |
Directional covariate |
h |
Bandwidth value |
xiBoot |
Bootstrap estimates for EVI |
sigBoot |
Bootstrap estimates for shape |
TTs |
T-year levels. For example, TTs = c(100, 10000). |
thetaGrid |
Grid values at which the estimation is performed |
timeRange |
Time range of the sample |
thr |
Threshold values along thetaGrid |
List including bootstrap estimates of T-year levels.
SplineML for examples.
## See also examples in vignettes: # vignette("localMethods", package = "circularEV") # vignette("splineML", package = "circularEV")## See also examples in vignettes: # vignette("localMethods", package = "circularEV") # vignette("splineML", package = "circularEV")
A vector of length 1521.
drcdrc
A vector of length 1521.
Directional covariate of HsSP data used by Reistad, M., Breivik, Ø., Haakenstad, H., Aarnes, O. J., Furevik, B. R., and Bidlot, J.-R. (2011), A high-resolution hindcast of wind and waves for the North Sea, the Norwegian Sea, and the Barents Sea, J. Geophys. Res., 116:1-18.
A vector of length 1521.
HsSPHsSP
A vector of length 1521.
HsSP data used by Reistad, M., Breivik, Ø., Haakenstad, H., Aarnes, O. J., Furevik, B. R., and Bidlot, J.-R. (2011), A high-resolution hindcast of wind and waves for the North Sea, the Norwegian Sea, and the Barents Sea, J. Geophys. Res., 116:1-18.
Local bootstrap estimation of EVI, scale and T-year levels
LocalEstim( Data, drc, thr = NULL, thetaGrid, nBoot = 100, EVIestimator = "Mom", h = 30, useKernel = TRUE, concent = 10, movThr = TRUE, TTs = NULL, timeRange = NULL )LocalEstim( Data, drc, thr = NULL, thetaGrid, nBoot = 100, EVIestimator = "Mom", h = 30, useKernel = TRUE, concent = 10, movThr = TRUE, TTs = NULL, timeRange = NULL )
Data |
Response variable |
drc |
Directional covariate |
thr |
Threshold values along thetaGrid |
thetaGrid |
Grid values at which the estimation is performed |
nBoot |
Number of bootstrap resamples. Default to 100. |
EVIestimator |
It can be either "ML" or "Mom" |
h |
Bandwidth value |
useKernel |
Logical. If TRUE (default), use kernel to assign weights depending on the directional distance. |
concent |
Concentration parameter value for von Mises kernel |
movThr |
Logical. If TRUE (default), moving threshold within the window used. |
TTs |
T-year levels. For example, TTs = c(100, 10000). |
timeRange |
Time range of the sample |
See Konzen, E., Neves, C., and Jonathan, P. (2021). Modeling nonstationary extremes of storm severity: Comparing parametric and semiparametric inference. Environmetrics, 32(4), e2667.
List including bootstrap estimates of EVI, scale and T-year levels.
data(HsSP) data(drc) timeRange <- 54.5 idx <- order(drc) drc <- drc[idx] Data <- HsSP[idx] set.seed(1234) Data <- Data + runif(length(Data), -1e-4, 1e-4) thetaVec <- 1:360 data(thresholdExampleMom) # loads threshold example thrResultMom <- thresholdExampleMom h <- 60 useKernel <- TRUE concent <- 10 movThr <- TRUE nBoot <- 30 set.seed(1234) output <- LocalEstim(Data=Data, drc=drc, thr=thrResultMom, thetaGrid=thetaVec, nBoot=nBoot, EVIestimator="Mom", h=h, useKernel=useKernel, concent=concent, movThr=movThr, TTs=c(100, 10000), timeRange=timeRange) RLBoot <- output$RLBoot PlotParamEstim(bootEstimates=output$xiBoot, thetaGrid=thetaVec, ylab=bquote(hat(xi)), alpha=0.05, ylim=NULL, cex.axis=15, cex.lab=2, thrWidth=2) PlotParamEstim(bootEstimates=output$sigBoot, thetaGrid=thetaVec, ylab=bquote(hat(sigma)), alpha=0.05, ylim=NULL, cex.axis=15, cex.lab=2, thrWidth=2) # 100-year level PlotRL(RLBootList=RLBoot, thetaGrid=thetaVec, Data=Data, drc=drc, TTs=c(100, 10000), whichPlot=1, alpha=0.05, ylim=NULL, pointSize=1, cex.axis=15, cex.lab=2, thrWidth=2) PolarPlotRL(RLBootList=RLBoot, thetaGrid=thetaVec, Data=Data, drc=drc, TTs=c(100, 10000), whichPlot=1, alpha=0.05, ylim=NULL, pointSize=4, fontSize=12, lineWidth=2) ## See examples in vignette: # vignette("localMethods", package = "circularEV")data(HsSP) data(drc) timeRange <- 54.5 idx <- order(drc) drc <- drc[idx] Data <- HsSP[idx] set.seed(1234) Data <- Data + runif(length(Data), -1e-4, 1e-4) thetaVec <- 1:360 data(thresholdExampleMom) # loads threshold example thrResultMom <- thresholdExampleMom h <- 60 useKernel <- TRUE concent <- 10 movThr <- TRUE nBoot <- 30 set.seed(1234) output <- LocalEstim(Data=Data, drc=drc, thr=thrResultMom, thetaGrid=thetaVec, nBoot=nBoot, EVIestimator="Mom", h=h, useKernel=useKernel, concent=concent, movThr=movThr, TTs=c(100, 10000), timeRange=timeRange) RLBoot <- output$RLBoot PlotParamEstim(bootEstimates=output$xiBoot, thetaGrid=thetaVec, ylab=bquote(hat(xi)), alpha=0.05, ylim=NULL, cex.axis=15, cex.lab=2, thrWidth=2) PlotParamEstim(bootEstimates=output$sigBoot, thetaGrid=thetaVec, ylab=bquote(hat(sigma)), alpha=0.05, ylim=NULL, cex.axis=15, cex.lab=2, thrWidth=2) # 100-year level PlotRL(RLBootList=RLBoot, thetaGrid=thetaVec, Data=Data, drc=drc, TTs=c(100, 10000), whichPlot=1, alpha=0.05, ylim=NULL, pointSize=1, cex.axis=15, cex.lab=2, thrWidth=2) PolarPlotRL(RLBootList=RLBoot, thetaGrid=thetaVec, Data=Data, drc=drc, TTs=c(100, 10000), whichPlot=1, alpha=0.05, ylim=NULL, pointSize=4, fontSize=12, lineWidth=2) ## See examples in vignette: # vignette("localMethods", package = "circularEV")
Plot of circular data
PlotData( Data, drc, thr = NULL, ylim = NULL, pointSize = 4, cex.axis = 15, cex.lab = 2, thrWidth = 2, thrColor = "#D45E1A", thrLineType = 1, ylab = NULL )PlotData( Data, drc, thr = NULL, ylim = NULL, pointSize = 4, cex.axis = 15, cex.lab = 2, thrWidth = 2, thrColor = "#D45E1A", thrLineType = 1, ylab = NULL )
Data |
Response variable |
drc |
Directional covariate |
thr |
Threshold values along thetaGrid |
ylim |
Range of values |
pointSize |
Size of points (observations) |
cex.axis |
Graphical parameter |
cex.lab |
Graphical parameter |
thrWidth |
Threshold width |
thrColor |
Threshold colour |
thrLineType |
Graphical parameter |
ylab |
y-axis label |
Plot of circular data, possibly including a threshold.
data(HsSP) data(drc) PlotData(Data=HsSP, drc=drc, thr=NULL, pointSize=1, cex.axis=15, cex.lab=2, thrWidth=2) data(thresholdExampleML) # loads threshold example PlotData(Data=HsSP, drc=drc, thr=thresholdExampleML, pointSize=1, cex.axis=15, cex.lab=2, thrWidth=2)data(HsSP) data(drc) PlotData(Data=HsSP, drc=drc, thr=NULL, pointSize=1, cex.axis=15, cex.lab=2, thrWidth=2) data(thresholdExampleML) # loads threshold example PlotData(Data=HsSP, drc=drc, thr=thresholdExampleML, pointSize=1, cex.axis=15, cex.lab=2, thrWidth=2)
Plot of parameter estimates with bootstrap confidence intervals
PlotParamEstim( bootEstimates, thetaGrid = 1:360, alpha = 0.05, ylim = NULL, cex.axis = 15, cex.lab = 2, thrWidth = 2, ylab = NULL, thrColor = "#D45E1A" )PlotParamEstim( bootEstimates, thetaGrid = 1:360, alpha = 0.05, ylim = NULL, cex.axis = 15, cex.lab = 2, thrWidth = 2, ylab = NULL, thrColor = "#D45E1A" )
bootEstimates |
Bootstrap estimates (for example, shape or scale) |
thetaGrid |
Grid values at which the estimation is performed |
alpha |
Significance level for the confidence intervals. Default to 0.05. |
ylim |
Range for the y-axis |
cex.axis |
Graphical parameter |
cex.lab |
Graphical parameter |
thrWidth |
Threshold width |
ylab |
y-axis label |
thrColor |
Threshold colour |
Plot of parameter estimates.
SplineML and LocalEstim for examples.
## See examples in vignettes: # vignette("localMethods", package = "circularEV") # vignette("splineML", package = "circularEV")## See examples in vignettes: # vignette("localMethods", package = "circularEV") # vignette("splineML", package = "circularEV")
Plot of T-year levels
PlotRL( RLBootList, Data, drc, thetaGrid = 1:360, TTs, whichPlot, alpha = 0.05, ylim = NULL, pointSize = 1, cex.axis = 15, cex.lab = 2, thrWidth = 2, thrColor = "#D45E1A", ylab = NULL )PlotRL( RLBootList, Data, drc, thetaGrid = 1:360, TTs, whichPlot, alpha = 0.05, ylim = NULL, pointSize = 1, cex.axis = 15, cex.lab = 2, thrWidth = 2, thrColor = "#D45E1A", ylab = NULL )
RLBootList |
List containing bootstrap estimates of T-year levels |
Data |
Response variable |
drc |
Directional covariate |
thetaGrid |
Grid values at which the estimation is performed |
TTs |
T-year levels. For example, TTs = c(100, 10000). |
whichPlot |
Index identifying which T-year level should be plotted from TTs. If TTs = c(100, 10000), then whichPlot=2 produces a plot for the 10000-year level |
alpha |
Significance level for the confidence intervals. Default to 0.05. |
ylim |
Range for the y-axis |
pointSize |
Size of points (observations) |
cex.axis |
Graphical parameter |
cex.lab |
Graphical parameter |
thrWidth |
Threshold width |
thrColor |
Threshold colour |
ylab |
y-axis label |
Plot of T-year levels.
SplineML and LocalEstim for examples.
## See also examples in vignettes: # vignette("localMethods", package = "circularEV") # vignette("splineML", package = "circularEV")## See also examples in vignettes: # vignette("localMethods", package = "circularEV") # vignette("splineML", package = "circularEV")
Polar plot of circular data
PolarPlotData( Data, drc, thr = NULL, ylim = NULL, pointSize = 1, fontSize = 12, thrWidth = 4, thrColor = "#D45E1A" )PolarPlotData( Data, drc, thr = NULL, ylim = NULL, pointSize = 1, fontSize = 12, thrWidth = 4, thrColor = "#D45E1A" )
Data |
Response variable |
drc |
Directional covariate |
thr |
Threshold values along thetaGrid |
ylim |
Range of values |
pointSize |
Size of points (observations) |
fontSize |
Font size |
thrWidth |
Threshold width |
thrColor |
Threshold colour |
Polar plot of circular data, possibly including a threshold
data(HsSP) data(drc) PolarPlotData(Data=HsSP, drc=drc, thr=NULL, pointSize=4, fontSize=14, thrWidth=4, ylim=c(0,max(HsSP))) data(thresholdExampleML) # loads threshold example PolarPlotData(Data=HsSP, drc=drc, thr=thresholdExampleML, pointSize=4, fontSize=12, thrWidth=4, ylim=c(0,max(HsSP)))data(HsSP) data(drc) PolarPlotData(Data=HsSP, drc=drc, thr=NULL, pointSize=4, fontSize=14, thrWidth=4, ylim=c(0,max(HsSP))) data(thresholdExampleML) # loads threshold example PolarPlotData(Data=HsSP, drc=drc, thr=thresholdExampleML, pointSize=4, fontSize=12, thrWidth=4, ylim=c(0,max(HsSP)))
Polar plot of T-year levels
PolarPlotRL( RLBootList, Data, drc, thetaGrid = 1:360, TTs, whichPlot, alpha = 0.05, ylim = NULL, pointSize = 4, fontSize = 12, lineWidth = 4 )PolarPlotRL( RLBootList, Data, drc, thetaGrid = 1:360, TTs, whichPlot, alpha = 0.05, ylim = NULL, pointSize = 4, fontSize = 12, lineWidth = 4 )
RLBootList |
List containing bootstrap estimates of T-year levels |
Data |
Response variable |
drc |
Directional covariate |
thetaGrid |
Grid values at which the estimation is performed |
TTs |
T-year levels. For example, TTs = c(100, 10000). |
whichPlot |
Index identifying which T-year level should be plotted from TTs. If TTs = c(100, 10000), then whichPlot=2 produces a plot for the 10000-year level |
alpha |
Significance level for the confidence intervals. Default to 0.05. |
ylim |
Range for the y-axis |
pointSize |
Size of points (observations) |
fontSize |
Font size |
lineWidth |
Threshold width |
Polar plot of T-year levels.
SplineML and LocalEstim for examples.
## See also examples in vignettes: # vignette("localMethods", package = "circularEV") # vignette("splineML", package = "circularEV")## See also examples in vignettes: # vignette("localMethods", package = "circularEV") # vignette("splineML", package = "circularEV")
Spline ML fitting
SplineML( excesses, drc, thetaVec = 0:360, nBoot = 100, numIntKnots = 10, knotsType = "eqSpaced", lambda = seq(0, 2, by = 0.5), kappa = seq(0, 2, by = 0.5), nCandidatesInit = 1000, numCores = 2 )SplineML( excesses, drc, thetaVec = 0:360, nBoot = 100, numIntKnots = 10, knotsType = "eqSpaced", lambda = seq(0, 2, by = 0.5), kappa = seq(0, 2, by = 0.5), nCandidatesInit = 1000, numCores = 2 )
excesses |
Excesses data |
drc |
Directional covariate |
thetaVec |
Grid values at which the threshold will be evaluated |
nBoot |
Number of bootstrap resamples |
numIntKnots |
Number of internal knots |
knotsType |
Position of knots. Default to "eqSpaced". Otherwise, the knots will be placed at the quantiles of observed directions. |
lambda |
Penalty parameter values for lambda |
kappa |
Penalty parameter values for kappa |
nCandidatesInit |
Number of initial parameter vectors. Optimisation will start with the best. |
numCores |
Number of CPU cores to be used |
See Konzen, E., Neves, C., and Jonathan, P. (2021). Modeling nonstationary extremes of storm severity: Comparing parametric and semiparametric inference. Environmetrics, 32(4), e2667.
List of bootstrap estimates of shape and scale, and optimal values of lambda and kappa.
data(HsSP) data(drc) timeRange <- 54.5 idx <- order(drc) drc <- drc[idx] Data <- HsSP[idx] set.seed(1234) Data <- Data + runif(length(Data), -1e-4, 1e-4) thetaVec <- 1:360 data(thresholdExampleML) # loads threshold example thrResultML <- thresholdExampleML lambda <- 100 kappa <- 40 thrPerObs <- thrResultML[drc] excess <- Data - thrPerObs drcExcess <- drc[excess>0] excess <- excess[excess>0] splineFit <- SplineML(excesses = excess, drc = drcExcess, nBoot = 30, numIntKnots = 16, lambda=lambda, kappa=kappa, numCores=2) xiBoot <- splineFit$xi sigBoot <- splineFit$sig PlotParamEstim(bootEstimates=xiBoot, thetaGrid=0:360, ylab=bquote(hat(xi)), alpha=0.05, ylim=NULL, cex.axis=15, cex.lab=2, thrWidth=2) PlotParamEstim(bootEstimates=sigBoot, thetaGrid=0:360, ylab=bquote(hat(sigma)), alpha=0.05, ylim=NULL, cex.axis=15, cex.lab=2, thrWidth=2) h <- 60 # needed for calculating local probability of exceedances RLBoot <- CalcRLsplineML(Data=Data, drc=drc, xiBoot=xiBoot, sigBoot=sigBoot, h=h, TTs=c(100, 10000), thetaGrid=thetaVec, timeRange=timeRange, thr=thrResultML) # 100-year level PlotRL(RLBootList=RLBoot, thetaGrid=thetaVec, Data=Data, drc=drc, TTs=c(100, 10000), whichPlot=1, alpha=0.05, ylim=NULL, pointSize=1, cex.axis=15, cex.lab=2, thrWidth=2) PolarPlotRL(RLBootList=RLBoot, thetaGrid=thetaVec, Data=Data, drc=drc, TTs=c(100, 10000), whichPlot=1, alpha=0.05, ylim=c(0, 25), pointSize=4, fontSize=12, lineWidth=2) ## See also examples in vignette: # vignette("splineML", package = "circularEV")data(HsSP) data(drc) timeRange <- 54.5 idx <- order(drc) drc <- drc[idx] Data <- HsSP[idx] set.seed(1234) Data <- Data + runif(length(Data), -1e-4, 1e-4) thetaVec <- 1:360 data(thresholdExampleML) # loads threshold example thrResultML <- thresholdExampleML lambda <- 100 kappa <- 40 thrPerObs <- thrResultML[drc] excess <- Data - thrPerObs drcExcess <- drc[excess>0] excess <- excess[excess>0] splineFit <- SplineML(excesses = excess, drc = drcExcess, nBoot = 30, numIntKnots = 16, lambda=lambda, kappa=kappa, numCores=2) xiBoot <- splineFit$xi sigBoot <- splineFit$sig PlotParamEstim(bootEstimates=xiBoot, thetaGrid=0:360, ylab=bquote(hat(xi)), alpha=0.05, ylim=NULL, cex.axis=15, cex.lab=2, thrWidth=2) PlotParamEstim(bootEstimates=sigBoot, thetaGrid=0:360, ylab=bquote(hat(sigma)), alpha=0.05, ylim=NULL, cex.axis=15, cex.lab=2, thrWidth=2) h <- 60 # needed for calculating local probability of exceedances RLBoot <- CalcRLsplineML(Data=Data, drc=drc, xiBoot=xiBoot, sigBoot=sigBoot, h=h, TTs=c(100, 10000), thetaGrid=thetaVec, timeRange=timeRange, thr=thrResultML) # 100-year level PlotRL(RLBootList=RLBoot, thetaGrid=thetaVec, Data=Data, drc=drc, TTs=c(100, 10000), whichPlot=1, alpha=0.05, ylim=NULL, pointSize=1, cex.axis=15, cex.lab=2, thrWidth=2) PolarPlotRL(RLBootList=RLBoot, thetaGrid=thetaVec, Data=Data, drc=drc, TTs=c(100, 10000), whichPlot=1, alpha=0.05, ylim=c(0, 25), pointSize=4, fontSize=12, lineWidth=2) ## See also examples in vignette: # vignette("splineML", package = "circularEV")
A vector of threshold values at directions 1,...,360. It is used for spline ML and local ML examples.
thresholdExampleMLthresholdExampleML
A vector of 360 values.
It has been generated as follows:
data(HsSP)
data(drc)
timeRange <- 54.5
idx <- order(drc)
drc <- drc[idx]
Data <- HsSP[idx]
set.seed(1234)
Data <- Data + runif(length(Data), -1e-4, 1e-4)
thetaVec <- 1:360
thresholdExampleML <- ThrSelection(Data=Data, drc=drc, h=60, b=0.35,
thetaGrid=thetaVec, EVIestimator="ML",
useKernel=T, concent=10, bw=30, numCores=2)$thr
A vector of threshold values at directions 1,...,360. It is used for local Moment examples.
thresholdExampleMomthresholdExampleMom
A vector of 360 values.
It has been generated as follows:
data(HsSP)
data(drc)
timeRange <- 54.5
idx <- order(drc)
drc <- drc[idx]
Data <- HsSP[idx]
set.seed(1234)
Data <- Data + runif(length(Data), -1e-4, 1e-4)
thetaVec <- 1:360
thresholdExampleMom <- ThrSelection(Data=Data, drc=drc, h=60, b=0.35,
thetaGrid=thetaVec, EVIestimator="Mom",
useKernel=T, concent=10, bw=30, numCores=2)$thr
This function selects a moving threshold for circular data using an automatic procedure for selecting the local number of exceedances
ThrSelection( Data, drc, h = 30, b = 0.35, thetaGrid, EVIestimator = "ML", useKernel = TRUE, concent = 10, bw = 30, numCores = 2 )ThrSelection( Data, drc, h = 30, b = 0.35, thetaGrid, EVIestimator = "ML", useKernel = TRUE, concent = 10, bw = 30, numCores = 2 )
Data |
Response variable |
drc |
Directional covariate |
h |
Bandwidth value |
b |
Parameter used in the automatic procedure for selection of local number of exceedances |
thetaGrid |
Grid values at which the estimation is performed |
EVIestimator |
It can be either "ML" or "Mom" |
useKernel |
Logical. If TRUE (default), use kernel to assign weights depending on the directional distance. |
concent |
Concentration parameter value for von Mises kernel |
bw |
Bandwidth parameter value for smoothing the sample path of the selected threshold |
numCores |
Number of CPU cores to be used |
See Konzen, E., Neves, C., and Jonathan, P. (2021). Modeling nonstationary extremes of storm severity: Comparing parametric and semiparametric inference. Environmetrics, 32(4), e2667.
List containing the selected threshold and selected number of local exceedances at each direction in the grid.
PlotData and PolarPlotData to see how
the threshold can be visualised.
data(HsSP) data(drc) timeRange <- 54.5 idx <- order(drc) drc <- drc[idx] Data <- HsSP[idx] set.seed(1234) Data <- Data + runif(length(Data), -1e-4, 1e-4) thetaVec <- 1:360 thrResultMom <- ThrSelection(Data=Data, drc=drc, h=60, b=0.35, thetaGrid=thetaVec, EVIestimator="Mom", useKernel=T, concent=10, bw=30, numCores=2)$thr thrResultML <- ThrSelection(Data=Data, drc=drc, h=60, b=0.35, thetaGrid=thetaVec, EVIestimator="ML", useKernel=T, concent=10, bw=30, numCores=2)$thr ## See also examples in vignettes: # vignette("localMethods", package = "circularEV") # vignette("splineML", package = "circularEV")data(HsSP) data(drc) timeRange <- 54.5 idx <- order(drc) drc <- drc[idx] Data <- HsSP[idx] set.seed(1234) Data <- Data + runif(length(Data), -1e-4, 1e-4) thetaVec <- 1:360 thrResultMom <- ThrSelection(Data=Data, drc=drc, h=60, b=0.35, thetaGrid=thetaVec, EVIestimator="Mom", useKernel=T, concent=10, bw=30, numCores=2)$thr thrResultML <- ThrSelection(Data=Data, drc=drc, h=60, b=0.35, thetaGrid=thetaVec, EVIestimator="ML", useKernel=T, concent=10, bw=30, numCores=2)$thr ## See also examples in vignettes: # vignette("localMethods", package = "circularEV") # vignette("splineML", package = "circularEV")