sumsim {dvmisc}

Summarize Simulation Results

Description

Creates table summarizing results of statistical simulations, providing common metrics of performance like mean bias, standard deviation, mean standard error, mean squared error, and confidence interval coverage.

Usage

sumsim(estimates, ses = NULL, truth = NULL, theta_0 = 0,
  statistics = c("mean_bias", "sd", "mean_se", "mse", "coverage"),
  alpha = 0.05, digits = 3, listwise_deletion = TRUE)

Arguments

estimates

Numeric matrix where each column gives the point estimates for a particular method across multiple trials.

ses

Numeric matrix where each column gives the standard errors for a particular method across multiple trials.

truth

Numeric value specifying the true value of the parameter being estimated.

theta_0

Numeric value specifying null value for hypothesis test H_0: theta = theta_0. Only used for calculating empirical power.

statistics

Numeric vector specifying which performance metrics should be calculated. Possible values are "n" for number of trials, "mean", "median", "mean_bias", "median_bias", "sd", "iqr", "mean_se" (for mean standard error), "mse" (for mean squared error), "coverage" (for confidence interval coverage), "ci_width" for median confidence interval width, and "power" for empirical power.

alpha

Numeric value specifying alpha for confidence interval. Set to 0.05 for the usual 95% CI, 0.1 for a 90% CI, and so forth.

digits

Numeric value or vector specifying the number of decimal places to include.

listwise_deletion

Logical value for whether to remove trials in which any of the estimators have missing values.

Value

Numeric matrix.

Examples

# For X ~ N(mu, sigma^2), the MLE for sigma^2 is the sample variance with n 
# in the denominator, but the unbiased version with (n - 1) is typically used 
# for its unbiasedness. Compare these estimators in 1,000 trials with n = 25. 
MLE <- c()
Unbiased <- c()
for (ii in 1: 1000) {
   x <- rnorm(n = 25)
   MLE[ii] <- sum((x - mean(x))^2) / 25
   Unbiased[ii] <- sum((x - mean(x))^2) / 24
 }
sumsim(estimates = cbind(MLE, Unbiased), truth = 1)