Common random numbers in multivariate simulations

We describe major improvements to the testing capabilities of SimOpt, a library of simulation-optimization problems and solvers. Foremost among these improvements is a transition to GitHub that makes SimOpt easier to use and maintain. We also design two new wrapper functions that facilitate empirical comparisons of solvers. The wrapper functions make extensive use of common random numbers (CRN) both within and across solvers for various purposes; e.g., identifying random initial solutions and running simulation replications. We examine some of the intricacies of using CRN to compare simulation-optimization solvers.

We present Q-SASS, a quasi-Newton method for unconstrained stochastic optimization that does not rely on common random numbers. Most existing quasi-Newton approaches leverage common random numbers to construct second-order updates. However, motivated by challenges in variational quantum algorithms—where such coordination is not possible—we consider the setting in which function values and gradients are accessible only through noisy probabilistic zeroth- and first-order oracles, and no common random numbers can be exploited. We derive high-probability tail bounds on the iteration complexity of our algorithm for nonconvex, convex, and strongly convex (more generally, those satisfying the PL condition) objective functions. Finally, we demonstrate the empirical benefits of our quasi-Newton updating scheme on both synthetic and quantum chemistry problems.

We consider noisy optimization and some traditional variance reduction techniques aimed at improving the convergence rate, namely (i) common random numbers (CRN), which is relevant for population-based noisy optimization and (ii) stratified sampling, which is relevant for most noisy optimization problems. We present artificial models of noise for which common random numbers are very efficient, and artificial models of noise for which common random numbers are detrimental. We then experiment on a desperately expensive unit commitment problem. As expected, stratified sampling is never detrimental. Nonetheless, in practice, common random numbers nonetheless provided, by far, most of the improvement.

In many practical problems, simulation models are used to support complex decision-making processes comparing hundreds or thousands of solutions. These problems typically have a key objective but the final decision may be dependent on other factors, which cannot be incorporated into the simulation model. In such cases, decision-makers may request a short list of ‘good’ solutions, which work well for the main objective and satisfy one or more chance constraints. While fully sequential ranking and selection procedures can be effective at solving these problems, surveys of experimentation practice suggest that they are under-utilized, potentially due to difficulties automating commercial software. We develop an approach with just two stages of replications. The approach, which has been designed to cope with the use of common random numbers, draws on ideas from indifference zones and makes use of bootstrapping to find a subset of high quality solutions. A Python implementation is freely available.

We describe a practical two-stage algorithm, BootComp, for multi-objective optimization via simulation. Our algorithm finds a subset of good designs that a decision-maker can compare to identify the one that works best when considering all aspects of the system, including those that cannot be modeled. BootComp is designed to be straightforward to implement by a practitioner with basic statistical knowledge in a simulation package that does not support sequential ranking and selection. These requirements restrict us to a two-stage procedure that works with any distributions of the outputs and allows for the use of common random numbers. Comparisons with sequential ranking and selection methods suggest that it performs well, and we also demonstrate its use analyzing a real simulation aiming to determine the optimal ward configuration for a UK hospital.

Fully sequential selection procedures have been developed in the field of stochastic simulation to find the simulated system with the best expected performance when the number of alternatives is finite. Kim and Nelson proposed the procedure to allow for unknown and unequal variances and the use of common random numbers. approximates the raw sum of differences between observations from two systems as a Brownian motion process with drift and uses a triangular continuation region to decide the stopping time of the selection process. In this paper new fully sequential selection procedures are derived that employ a more effective sum of differences, which is called a controlled sum. Two provably valid procedures and an approximate procedure are described. Empirical results and a realistic illustration are provided to compare the efficiency of these procedures with other procedures that solve the same problem. [Supplemental materials are available for this article. Go to the publisher's online edition of IIE Transactions for the following supplemental resources: Proofs and guidelines to choose appropriate parameters.]

The design and analysis of Monte Carlo experiments, with special reference to structural equation modelling, is discussed in this article. These topics merit consideration, since the validity of the conclusions drawn from a Monte Carlo study clearly hinges on these features. It is argued that comprehensive Monte Carlo experiments can be implemented on a PC if the experiments are adequately designed. This is especially important when investigating modern computer intensive methodologies like resampling and Markov Chain Monte Carlo methods. We are faced with three fundamental challenges in Monte Carlo experimentation. The first problem is statistical precision, which concerns the reliability of the obtained results. External validity, on the other hand, depends on the number of experimental conditions, and is crucial for the prospects of generalising the results beyond the specific experiment. Finally, we face the constraint on available computer resources. The conventional wisdom in designing and analysing Monte Carlo experiments embodies no explicit specification of meta-model for analysing the output of the experiment, the use of case studies or full factorial designs as experimental plans, no use of variance reduction techniques, a large number of replications, and "eyeballing" of the results. A critical examination of the conventional wisdom is presented in this article. We suggest that the following alternative procedures should be considered. First of all, we argue that it is profitable to specify explicit meta-models, relating the chosen performance statistics and experimental conditions. Regarding the experimental plan, we recommend the use of incomplete designs, which will often result in considerable savings. We also consider the use of common random numbers in the simulation phase, since this may enhance the precision in estimating meta-models. The use of fewer replications per trial, enabling us to investigate an increased number of experimental conditions, should also be considered in order to improve the external validity at the cost of the conventionally excessive precision.

To analyze simulation runs which use the same random numbers, the blocking concept of experimental design is not needed. Instead, this paper applies a linear regression model with a nondiagonal covariance matrix. This covariance matrix does not need to have a specific pattern such as constant covariances. A simple example yields surprising results. The paper proposes a new framework for the error analysis. This framework consists of three factors (namely, common random numbers, replication, model validity), each with three levels.

Parametric sensitivity of biochemical networks is an indispensable tool for studying system robustness properties, estimating network parameters, and identifying targets for drug therapy. For discrete stochastic representations of biochemical networks where Monte Carlo methods are commonly used, sensitivity analysis can be particularly challenging, as accurate finite difference computations of sensitivity require a large number of simulations for both nominal and perturbed values of the parameters. In this paper we introduce the common random number (CRN) method in conjunction with Gillespie's stochastic simulation algorithm, which exploits positive correlations obtained by using CRNs for nominal and perturbed parameters. We also propose a new method called the common reaction path (CRP) method, which uses CRNs together with the random time change representation of discrete state Markov processes due to Kurtz to estimate the sensitivity via a finite difference approximation applied to coupled reaction paths that emerge naturally in this representation. While both methods reduce the variance of the estimator significantly compared to independent random number finite difference implementations, numerical evidence suggests that the CRP method achieves a greater variance reduction. We also provide some theoretical basis for the superior performance of CRP. The improved accuracy of these methods allows for much more efficient sensitivity estimation. In two example systems reported in this work, speedup factors greater than 300 and 10,000 are demonstrated.

We consider the problem of stochastic simulation optimization with common random numbers over a numerical search domain. We propose the Knowledge Gradient for Common Random Numbers (KG-CRN) sequential sampling algorithm, a simple elegant modification to the Knowledge Gradient that incorporates the use of correlated noise in simulation outputs with Gaussian Process meta-models. We compare this method against the standard Knowledge Gradient and a more recently proposed variation that allows for pairwise sampling. Our method significantly outperforms both baselines under identical laboratory conditions while greatly reducing computational cost compared to pairwise sampling.

The article describes some recently-proposed procedures that identify the best simulated system when common random numbers are used. The procedures are based on a Bayesian average-case analysis, rather than a worst-case indifference zone formulation. The procedures allow decision makers to focus on reducing either the expected opportunity cost loss associated with potentially selecting an inferior system, or the probability of incorrect selection. Numerical experiments indicate that the new procedures outperform two existing procedures with respect to several criteria for a well-known selection problem.

Chi-square tests based on maximum likelihood (ML) estimation of covariance structures often incorrectly over-reject the null hypothesis: Σ = Σ θ when the sample size is small. Reweighted least squares (RLS) avoids this problem. In some models, the vector of parameter must contain means, variances, and covariances, yet whether RLS also works in mean and covariance structures remains unexamined. This research extends RLS to mean and covariance structures, evaluating a generalized least squares function with ML parameter estimates. A Monte Carlo simulation study was carried out to examine the statistical performance of ML vs RLS with multivariate normal data. Based on empirical rejection frequencies and empirical averages of test statistics, this study shows that RLS performs much better than ML in mean and covariance structure models when sample sizes are small, whereas it does not perform better than ML to reject misspecified models.

ABSTRACTThe central composite design (CCD) is often used to estimate the second-order linear model. This paper uses a correlation induction strategy of common random numbers (CRN) in simulation experiment and utilizes the induced correlations to obtain better estimates for the second-order linear model. This strategy assigns the CRN to all design points in the CCD. An appropriate selection of the axial points in CCD makes the weighted least squares (WLS) estimator be equivalent to ordinary least squares (OLS) estimator in estimating the linear model parameters of CCD. We analytically investigate the efficiency of this strategy in estimation of model parameters. Under certain conditions, this correlation induction strategy yields better results than independent random number strategy in estimating model parameters except intercept. The simulation experiment on a selected model supports such results. We expect a suggested random number assignment is useful in application of CCD in simulation experiments.

Traditionally, exactness, numerical stability and speed are the three main criteria for evaluating algorithms for random variate generation. However, it is sometimes required that the algorithms provide correlation between generated variates for the purpose of inducing dependence among the output of simulation runs. The inverse transformation, which produces optimal correlation induction, often performs poorly in terms of the first three criteria. Algorithms based on composition, rejection, and special properties which often excel in terms of the first three criteria, tend to scramble the use of random numbers, causing many attempts at common random numbers, antithetic variates and external control variates to fail. The concept of obtaining correlation via algorithms other than the inverse transformation is examined here. To demonstrate feasibility, previously developed algorithms for Poisson and binomial random variate generation are modified to obtain both positive and negative correlation between runs....

We consider optimization via simulation over a finite set of alternatives. We employ a Bayesian value-ofinformation approach in which we allow both correlated prior beliefs on the sampling means and correlated sampling. Correlation in the prior belief allow us to learn about an alternative’s value from samples of similar alternatives. Correlation in sampling, achieved through common random numbers, allows us to reduce the variance in comparing one alternative to another. We allow for a more general combination of both types of correlation than has been offered previously in the Bayesian ranking and selection literature. We do so by giving an exact expression for the value of information for sampling the difference between a pair of alternatives, and derive new knowledge-gradient methods based on this valuation. 1