1. Objective:

Generate a PIT histogram for the ETAS model applied to the Hector Mine sequence using the pixel partition. The histogram, with 90% simulated confidence bands, is modeled after figure 4 in Thorarinsdottir (2013).

## deldir 0.0-22

2. Components

The Hector Mine sequence has been loaded in as catalog, the simulated empirical distribution of \(N(A)\) under the ETAS model fit to this sequence has been loaded in as N_vec.

3. Randomized PIT Histogram

We construct the histogram of the randomized PIT values associated with the observed counts in each of 400 pixels.

n_pixels <- 400
polygon_list <- pix_to_polygon(n_pixels, rw = c(0, 1, 0, 1))
counts <- counts_in_pix(catalog[, 1], catalog[,1], polygon_list)

# Calculate randomized PIT values
x_zero <- counts[counts == 0]
x_nzero <- counts[counts != 0]
v <- runif(n_pixels)

# Load in Empirical ETAS cdf
N_vec <- unlist(read.csv("../simulations/N_vec.csv", header = T))
Fx <- pETAS(x_nzero, N_vec)
Fxm1 <- pETAS(x_nzero - 1, N_vec)
u_nzero <- Fxm1 + v[counts != 0] * (Fx - Fxm1)
u_zero <- v[counts == 0] * pETAS(0, N_vec)
u <- c(u_nzero, u_zero)

90% Bands

Before plotting the histogram, we load in the simulated distribution of the counts in the five histogram bins under the null hypothesis, qmat. The 5% and the 95% quantiles of these empirical distributions will serve as the confidence bands.

umat <- as.matrix(read.csv("../simulations/umat.csv"))
qmat <- apply(umat, 1, quantile, probs = c(.05, .95))
ymax <- 1.1 * max(table(cut(u, breaks = seq(0, 1, .2))), max(qmat[2, ]))
hist(u, breaks = 4, col = "gray", border = "gray", xlab = "PIT value", ylab = "Frequency",
     ylim = c(0, ymax), main = "ETAS model")
break_x <- seq(0, 1, .2)
for(i in 1:dim(qmat)[2]) {
  lines(x = c(break_x[i], break_x[i + 1]), y = rep(qmat[1, i], 2), lty = 2)
  lines(x = c(break_x[i], break_x[i + 1]), y = rep(qmat[2, i], 2), lty = 2)
}