find.threshold: Find thresholds

In PierreMasselot/hhws: Functions for heat-health warning systems

Description

Evaluate different combinations of threshold and indicator weightings according to their sensitivity and specificity.

Usage

find.threshold(indicators, episodes, u.grid, fixed.alphas = NULL,
  alpha.step = 0.1, decreasing.alphas = TRUE, same.alphas = TRUE,
  order.result = NULL, order.decreasing = TRUE, r = NULL,
  trim = NULL, thinning = c("none", "days", "episodes"),
  progressBar = FALSE)

Arguments

indicators	List of matrices containing the lagged covariates to construct weighted indicators. Can also be a matrix in the case of a single indicator.
episodes	A matrix or data.frame containing the indices of extreme episodes. Must contain the time indices in the first column and the corresponding episode index in the second column as returned by episodes.
u.grid	A list of vectors containing the grid of thresholds to be considered. The list must have the same number of elements than indicators. If missing, the default is to use all percentiles from 0.8 to 1.
fixed.alphas	A list of optional prefixed weightings for a subset of indicators. When provided, must have the same length as indicators and contain NULL for each non-fixed indicator.
alpha.step	Numeric value between 0 and .5. The step of the sequence of weightings tested for the indicator.
decreasing.alphas	Logical. If TRUE (the default), the alpha weightings are constrained to decrease with the lag.
same.alphas	Logical. If TRUE, the weightings are constrained to be identical for each indicator. Note that trying different weightings (same.alpha = FALSE) can lead to significantly longer computation time.
order.result	Character or numeric value indicating a column used to order the returned table. Can also be a vector.
order.decreasing	Logical. If TRUE (the default), ordering of the table is made by decreasing order to the column specified by order.result.
r	Positive integer. Number of consecutive values below threshold following an excess to end the episode. By default, take the attribute r in episodes or 1 if absent.
trim	Positive integer. If not NULL (the default), the first trim lines of the result are returned by the function. Depends on the parameter order.result.
thinning	Character string indicating if the results should be thinned before returning. When thinning = "none" (the default), all tested combinations are returned. When thinning = "days", combinations for which there is at least another combinations such that specificity and sensivity are both higher, are removed from the returned results. When thinning = "episodes" the thinning is performed on the basis of episodes instead of days.
progressBar	Logical indicating if a progressBar is displayed during execution of the function. If TRUE, may greatly increase execution time.

Details

We consider a warning system as a couple indicator/threshold used to launch alerts when forecasts of the indicator exceed the threshold. In the present function, the indicators considered are linear combinations of all matrix columns in the parameter indicators, with the constraint that, for each indicator, the weights sum to 1.

The indicator and threshold are determined by evaluating a large range of different weightings and threshold (given in u.grid and alpha.step). For each combination of indicators/thresholds, the function computes the indices corresponding to alerts in the data and compare them to the actual values given in episodes. If thinning != "none") the function then selects a subset of best candidates for which the number of detected days are maximum while the number of false alarms is minimum. Note that if thinning = "episodes", the subset is selected on the basis of detected and false episodes instead of days. It is left to the user to choose the best combination by a trade-off between specificity and sensitivity.

Value

A data.frame containing a subset (unless thinning = "none") of weightings and thresholds. Weightings correspond to the columns with a name containing "alpha" and threshold to names beginning with "threshold". In addition, several scores are given in each line:

Detected	The number of indices in episodes detected by the combination.
Missed	The number of indices in episodes missed by the combination.
Sensitivity	The proportion of indices in episodes detected by the combination, i.e. Detected / n.
False_alarms	The number of false alarms, i.e. of indices found by the combination which are not in episodes.
Specificity	The proportion of false alarms, i.e. False_alarms / n.
Episodes_found	The number of episodes found. An episode is found when at least one of its days is found.
Episodes_sensitivity	The proportion of episodes found.
False_episodes	The number of false episodes found, i.e. absent from the provides episodes.

References

Chebana F., Martel B., Gosselin P., Giroux J.X., Ouarda T.B.M.J., 2013. A general and flexible methodology to define thresholds for heat health watch and warning systems, applied to the province of Quebec (Canada). International journal of biometeorology 57, 631-644.

Pascal M., Laaidi K., Ledrans M., Baffert E., Caserio-Schonemann C., Le Tertre A., Manach J., Medina S., Rudant J., Empereur-Bissonnet P., 2006. France's heat health watch warning system. International journal of biometeorology 50, 144-153.

See Also

episodes for extracting episodes of extreme values and predict_alarms for alarms prediction.

Examples

   library(dlnm)
   data(chicagoNMMAPS)
   x <- chicagoNMMAPS$death
   dates <- as.POSIXlt(chicagoNMMAPS$date)
   n <- nrow(chicagoNMMAPS)

   # Compute over-mortality
   om <- excess(x, dates = dates, order = 15)
  
   # Extract all days for which om is above 40%
   epis <- episodes(om, u = 40)

   # Prepare indicator based on temperature until lag 2
   indic <- matrix(NA, nrow = n, ncol = 3)
   indic[,1] <- chicagoNMMAPS$temp  # lag 0
   indic[,2] <- c(NA, chicagoNMMAPS$temp[-n]) # Lag 1
   indic[,3] <- c(NA, NA, chicagoNMMAPS$temp[1:(n-2)]) # lag 2

   # Evaluate different threshold/indicators based on these episodes
   find.threshold(indic, epis, u.grid = 20:35, thinning = "episodes",
     order.result = "Episodes_found")

PierreMasselot/hhws documentation built on May 25, 2019, 4:01 a.m.