re_estimate function

Re-estimate the posterior distributions with different burn-in

After running the MCMC chain for the given number of steps, the trace plots may indicate that too small value of burn-in was used in the first place. This function enables re-estimating the posterior distributions with a different value of burn-in, without the need to run the MCMC chain again.

re_estimate(Output, BurnIn = 0)

Arguments

• Output: list, output of the main function estintp.
• BurnIn: new value of burn-in.

Returns

List containing the parameter estimates along with the 2.5% and 97.5% quantiles of the posterior distributions, along with auxiliary objects needed for printing and plotting the outputs.

Details

The output of the main function binspp contains all the intermediate states of the chain (sampled with the required frequency) no matter what the original value of burn-in was. This enables simple and quick re-estimation of the posterior distributions with either higher or lower value of burn-in than the one used originally. The output of this function has the same structure as the output of the main function estintp().

Examples

library(spatstat)
# Prepare the dataset:
X = trees_N4
x_left = x_left_N4
x_right = x_right_N4
y_bottom = y_bottom_N4
y_top = y_top_N4

z_beta = list(refor = cov_refor, slope = cov_slope)
z_alpha = list(tmi = cov_tmi, tdensity = cov_tdensity)
z_omega = list(slope = cov_slope, reserv = cov_reserv)

# Determine the union of rectangles:
W = owin(c(x_left[1], x_right[1]), c(y_bottom[1], y_top[1]))
if (length(x_left) >= 2) {
  for (i in 2:length(x_left)) {
    W2 = owin(c(x_left[i], x_right[i]), c(y_bottom[i], y_top[i]))
    W = union.owin(W, W2)
  }
}

# Dilated observation window:
W_dil = dilation.owin(W, 100)

# Default parameters for prior distributions:
control = list(NStep = 100, BurnIn = 50, SamplingFreq = 5)

# MCMC estimation:
Output = estintp(X, control, x_left, x_right, y_bottom, y_top, W_dil,
                 z_beta, z_alpha, z_omega, verbose = FALSE)

# Text output + series of figures:
print_outputs(Output)
plot_outputs(Output)

# Recompute the outputs when another value of burn-in is desired,
# without running the chain again:
Out2 <- re_estimate(Output, BurnIn = 80)
print_outputs(Out2)
plot_outputs(Out2)