# Comparing estimation procedures via Monte Carlo simulations

I am looking for references of studies that show how to use Monte Carlo simulations to compare different estimators for any given parameter of any probability distribution (for example, comparing the MLE vs. the method-of-moments estimator, for the mean in an exponential distribution). I would like to know which properties are usually compared (like bias, MSE, etc.), and in general what things can be studied to compare their performance.

A Monte Carlo experiment towards comparing estimators is nothing special as it relies on simulation to approximate expectations, as any other Monte Carlo experiments. Hence, if you have a collection of procedures or estimators $\delta_1,\ldots,\delta_p$ to be compared under a loss or penalty function $\text{L}(\delta,\theta)$, the expected loss of procedure $\delta_i$ $$\mathbb{E}_\theta[\text{L}(\delta,\theta)]=\int_\mathfrak{X} \text{L}(\delta(x),\theta)f(x;\theta)\,\text{d}x$$where $f(x;\theta)$ denotes the density of the observation(s), is replaced with a Monte Carlo approximation$$\mathbb{E}_\theta[\text{L}(\delta_i,\theta)]\approx\frac{1}{T}\sum_{t=1}^T \omega_t \delta_i(x_t)$$where the $x_t$'s are generated according to a certain stochastic process and the $\omega_t$'s are appropriate weights to make the Monte Carlo approximation either (approximately) unbiased or with a finite variance. Both the process and the weights are usually dependent on $\theta$. Once the Monte Carlo experiment is designed, graphs of the different expected losses can be plotted against $\theta$ to assess the properties of those different estimators. It is usually recommended to correlate (or recycle) the simulated variates to reduce the variability in the graphs. Here is an example taken from our book, Introducing Monte Carlo Methods with R (2010, Fig. 4.8), comparing Monte Carlo approximations to the expected quadratic loss (or risk) $$\mathbb{E}_\theta[||\delta(X)-\theta||^2]$$ of several James-Stein estimators as the norm $||\theta||$ of a Normal mean $\theta\in\mathbb{R}^5$ varies from $0$ to $5$:
The above curves of the expected quadratic losses are smooth because they are based on the same $\mathrm{N}_5(0,I_5)$ simulations for all values of $\theta$.