Efficient computation of multivariate empirical distribution functions at the observed values

Abstract

Consider the evaluation of model-based functions of cumulative distribution functions that are integrals. When the cumulative distribution function does not have a tractable form but simulation of the multivariate distribution is easily feasible, we can evaluate the integral via a Monte Carlo sample, replacing the model-based distribution function by the empirical distribution function. Given a simulation sample of size N, the naive method uses O(N^2) comparisons to compute the empirical distribution function at all N sample vectors. To obtain faster computational speed when N needs to be large to achieve a desired accuracy, we propose methods modified from the popular merge sort and quicksort algorithms that preserve their average O(N*logN) complexity in the bivariate case. The modified merge sort algorithm can be extended to the computation of a d-dimensional empirical distribution function at the observed values with O(N*(logN)^(d-1)) complexity. Simulation studies suggest that the proposed algorithms provide substantial time savings when N is large.