MonteCarlo {mistral} | R Documentation |
Estimate a failure probability using a crude Monte Carlo method.
MonteCarlo(dimension, lsf, N_max = 5e+05, N_batch = 1000, q = 0, lower.tail = TRUE, precision = 0.05, plot = FALSE, output_dir = NULL, verbose = 0)
dimension |
the dimension of the input space. |
lsf |
the function defining safety/failure domain. |
N_max |
maximum number of calls to the |
N_batch |
number of points onte evalutae the |
q |
the quantile |
lower.tail |
as for pxxxx functions, TRUE for estimating P(lsf(X) < q), FALSE for P(lsf(X) > q) |
precision |
a targeted maximum value for the coefficient of variation. |
plot |
to plot the contour of the |
output_dir |
to save a copy of the plot in a pdf. This name will be pasted with "_Monte_Carlo_brut.pdf" |
verbose |
to control the level of outputs in the console; either 0 or 1 or 2 for almost no outputs to a high level output. |
This implementation of the crude Monte Carlo method works with evaluating batchs of points sequentialy until a given precision is reached on the final estimator
An object of class list
containing the failure probability and some
more outputs as described below:
p |
the estimated probabilty. |
ecdf_MC |
the empiracal cdf got with the generated samples. |
cov |
the coefficient of variation of the Monte Carlo estimator. |
Ncall |
the total numnber of calls to the |
X |
the generated samples. |
Y |
the value |
Problem is supposed to be defined in the standard space. If not, use UtoX
to do so. Furthermore, each time a set of vector is defined as a matrix, ‘nrow’
= dimension
and ‘ncol’ = number of vector to be consistent with
as.matrix
transformation of a vector.
Algorithm calls lsf(X) (where X is a matrix as defined previously) and expects a vector in return. This allows the user to optimise the computation of a batch of points, either by vectorial computation, or by the use of external codes (optimised C or C++ codes for example) and/or parallel computation.
Clement WALTER clement.walter@cea.fr
R. Rubinstein and D. Kroese:
Simulation and the Monte Carlo method
Wiley (2008)
#First some considerations on the usage of the lsf. #Limit state function defined by Kiureghian & Dakessian : # Remember you have to consider the fact that the input will be a matrix ncol >= 1 lsf_wrong = function(x, b=5, kappa=0.5, e=0.1) { b - x[2] - kappa*(x[1]-e)^2 # work only with a vector of lenght 2 } lsf_correct = function(x){ apply(x, 2, lsf_wrong) } lsf = function(x, b=5, kappa=0.5, e=0.1) { x = as.matrix(x) b - x[2,] - kappa*(x[1,]-e)^2 # vectorial computation, run fast } y = lsf(X <- matrix(rnorm(20), 2, 10)) #Compare running time ## Not run: require(microbenchmark) X = matrix(rnorm(2e5), 2) microbenchmark(lsf(X), lsf_correct(X)) ## End(Not run) #Example of parallel computation require(doParallel) lsf_par = function(x){ foreach(x=iter(X, by='col'), .combine = 'c') %dopar% lsf(x) } #Try Naive Monte Carlo on a given function with different failure level ## Not run: res = list() res[[1]] = MonteCarlo(2,lsf,q = 0,plot=TRUE) res[[2]] = MonteCarlo(2,lsf,q = 1,plot=TRUE) res[[3]] = MonteCarlo(2,lsf,q = -1,plot=TRUE) ## End(Not run) #Try Naive Monte Carlo on a given function and change number of points. ## Not run: res = list() res[[1]] = MonteCarlo(2,lsf,N_max = 10000) res[[2]] = MonteCarlo(2,lsf,N_max = 100000) res[[3]] = MonteCarlo(2,lsf,N_max = 500000) ## End(Not run)