Simulate the increasing random walk associated with a real-valued continuous random variable.
IRW(dimension, lsf, N = 10, q = Inf, Nevent = Inf, particles,
LSF_particles = lsf(particles), K, burnin = 20, sigma = 0.3,
last.return = TRUE, use.potential = TRUE, plot = FALSE,
print_plot = FALSE, output_dir = NULL)dimension of the input space.
limit state function.
number of particules.
level until which the randow walk is to be generated.
the number of desired events.
to start with some given particles.
value of the lsf on these particles.
kernel transition for conditional generations.
burnin parameter.
radius parameter for K.
if the last event should be returned.
tu use a ‘potential’ matrix to select starting point not directly related to the sample to be moved with the MH algorithm.
if TRUE, the algorithm plots the evolution of the particles. This
requieres to evaluate the lsf on a grid and is only for visual purpose.
if TRUE, print the updated plot after each iteration. This might
be slow; use with a small N. Otherwise it only prints the final plot.
if plots are to be saved in pdf in a given directory. This will
be pasted with ‘_IRW.pdf’. Together with print_plot==TRUE this will
produce a pdf with a plot at each iteration, enabling ‘video’ reconstitution
of the algorithm.
An object of class list containing the following data:
the events of the random walk.
the total number of iterations.
the total number of calls to the lsf.
a matrix containing the final particles.
the value of lsf on the particles.
the threshold considered when generating the random walk.
the target number of events when generating the random walk.
the number of rejected transitions, ie when the proposed point was not stricly greater/lower than the current state.
a vector containing the acceptance rate for each use of the MH algorithm.
This function lets generate the increasing random walk associated with a continous
real-valued random variable of the form Y = lsf(X) where X is
vectorial random variable.
This random walk can be associated with a Poisson process with parameter
N and hence the number of iterations before a given threshold q
is directly related to P[ lsf(X) > q]. It is the core tool of algorithms
such as nested sampling, Last Particle Algorithm or Tootsie Pop Algorithm.
Bascially for N = 1, it generates a sample \(Y = lsf(X)\) and iteratively
regenerates greater than the sought value: \(Y_{n+1} \sim \mu^Y( \cdot \mid Y > Y_n\). This
regeneration step is done with a Metropolis-Hastings algorithm and that is why it is usefull
to consider generating several chains all together (N > 1).
The algorithm stops when it has simulated the required number of events Nevent or when
it has reached the sought threshold q.
C. Walter: Moving Particles: a parallel optimal Multilevel Splitting method with application in quantiles estimation and meta-model based algorithms Structural Safety, 55, 10-25.
C. Walter: Point Process-based Monte Carlo estimation arXiv preprint arXiv:1412.6368.
J. Skilling: Nested sampling for general Bayesian computation Bayesian Analysis, 1(4), 833-859.
M. Huber \& S. Schott: Using TPA for Bayesian inference Bayesian Statistics 9, 9, 257.
A. Guyader, N. Hengartner and E. Matzner-Lober: Simulation and estimation of extreme quantiles and extreme probabilities Applied Mathematics \& Optimization, 64(2), 171-196.
# NOT RUN {
# Get faililng samples for the kiureghian limit state function
# Failure is defined as lsf(X) < 0 so we have to invert the lsf
lsf <- function(x) -1*kiureghian(x)
# }
# NOT RUN {
fail.samp <- IRW(2, lsf, q = 0, N = 10, plot = TRUE)
# }
Run the code above in your browser using DataLab