# R/plot_cdf.R In spesenti/SWIM: Scenario Weights for Importance Measurement

#### Documented in plot_cdf

```#' Plotting the Empirical Distribution Functions of a Stressed Model
#'
#' Plots the empirical distribution function of a stressed model
#'     component (random variable) under the scenario weights.
#'
#' @inheritParams  sensitivity
#' @inheritParams  plot_sensitivity
#' @inheritParams  summary.SWIM
#' @param xCol     Numeric or character, (name of) the column of the underlying data
#'                 of the \code{object} (\code{default = 1}).
#' @param n        Integer, the number of points used to plot
#'                 \code{stat_ecdf} in \code{ggplot} (\code{default
#'                 = 500}).
#' @param x_limits Vector, the limits of the x-axis of the plot, the
#'                 value for \code{xlim} in the \code{coord_cartesian}
#'                 function in \code{ggplot}.
#' @param y_limits Vector, the limits of the y-axis of the plot, the
#'                 value for \code{ylim} in the \code{coord_cartesian}
#'                 function in \code{ggplot}.
#
#' @return If \code{displ = TRUE}, a plot displaying the empirical
#'     distribution function of the stochastic model under the
#'     scenario weights.
#'
#'     If \code{displ = FALSE}, a data.frame for customised plotting with
#'     \code{ggplot}. The data.frame contains the columns: the column,
#'     \code{xCol}, of the data of the stressed model,
#'     \code{stress} (the stresses) and \code{value} (the values). \cr
#'     Denote by \code{res} the return of the function call, then
#'     \code{ggplot} can be called via:
#'     \deqn{ggplot(res, aes(x = res[ ,1], w = value))}
#'     \deqn{ + stat_ecdf(aes(color = factor(stress)), n = n).}
#'     Note that the ggplot2 default of \code{stat_ecdf} does not
#'     take \code{weight} as an aesthetic. We use the workaround
#'     by Nicolas Woloszko, see Note below.
#'
#' @note This function is based on the ggplot \code{stat_ecdf}
#'     function. However, the \code{stat_ecdf} does not allow for
#'     specifying \code{weights}, thus the function is based on
#'     the workaround by Nicolas Woloszko, see
#'     \url{https://github.com/NicolasWoloszko/stat_ecdf_weighted}.
#'
#' @examples
#' ## example with a stress on VaR
#' set.seed(0)
#' x <- as.data.frame(cbind(
#'   "normal" = rnorm(10^5),
#'   "gamma" = rgamma(10^5, shape = 2)))
#' res1 <- stress(type = "VaR", x = x,
#'   alpha = c(0.75, 0.95), q_ratio = 1.15)
#' plot_cdf(res1, xCol = 1, wCol = 1:2, base = TRUE)
#' plot_cdf(res1, xCol = 1, wCol = 1:2, base = TRUE,
#'   x_limits = c(0, 5), y_limits = c(0.5, 1))
#'
#' @seealso See \code{\link{cdf}} for the empirical distribution function
#'     of a stressed model and \code{\link{quantile_stressed}} for
#'     sample quantiles of a stressed model.
#'
#' @export

plot_cdf <- function(object, xCol = 1, wCol = "all", base = FALSE, n = 500,                             x_limits, y_limits, displ = TRUE){
if (!is.SWIM(object)) stop("Object not of class SWIM")
if (anyNA(object\$x)) warning("x contains NA")
x_data <- get_data(object)[, xCol]
if(is.character(xCol)) x_name <- xCol
if(is.null(colnames(get_data(object)))) x_name <- paste("X", xCol, sep = "") else if(!is.character(xCol)) x_name <- colnames(get_data(object))[xCol]
if (is.character(wCol) && wCol == "all") wCol <- 1:ncol(get_weights(object))
plot_data <- data.frame(x_data, get_weights(object)[ , wCol])
names(plot_data) <- c(x_name, paste("stress", wCol, sep = " "))
if (base == TRUE){
plot_data <- cbind(plot_data, base = rep(1, length(x_data)))
}
plot_data <- reshape2::melt(plot_data, id.var = x_name, variable.name = "stress", value.name = "value")

if (displ == TRUE){
if (missing(x_limits)) x_limits <- c(min(x_data)-0.1, max(x_data))
if (missing(y_limits)) y_limits <- c(0,1)
ggplot2::ggplot(plot_data, ggplot2::aes_(x = plot_data[,1], weight = ~value)) +
stat_ecdf(ggplot2::aes(color = factor(stress)), n = n) +
ggplot2::labs(x = x_name, y = "ecdf") +
ggplot2::coord_cartesian(xlim = x_limits, ylim = y_limits) +
ggplot2::theme_minimal() +
ggplot2::theme(legend.title = ggplot2::element_blank(), legend.key = ggplot2::element_blank(), legend.text = ggplot2::element_text(size = 10))
} else {
return(plot_data)
}
}

```
spesenti/SWIM documentation built on Jan. 11, 2020, 3:33 a.m.