Previous: Compatibility with sim()

Up: Compatibility with sim()

Next: Calculating Quantities of Interest

Simulating Parameters

Whether you choose to use the default method, or write a model-specific method for simulating parameters, these functions require the same three inputs:

• object: the estimated model or zelig() output.
• num: the number of simulations.
• bootstrap: either TRUE or FALSE.

The output from param() should be either

• If bootstrap = FALSE (default), an matrix with rows corresponding to simulations and columns corresponding to model parameters. Any ancillary parameters should be included in the output matrix.
• If bootstrap = TRUE, a vector containing all model parameters, including ancillary parameters.

There are two ways to simulate parameters:

1. Use the param.default() function to extract parameters from the model and, if bootstrapping is not selected, simulate coefficients using asymptotic normal approximation. The param.default() function relies on two R functions:
   1. coef(): extracts the coefficients. Continuing the Normal regression example from above, the appropriate coef.normal() function is simply:

      coef.normal <- function(object)
        object$coefficients
      

   2. vcov(): extracts the variance-covariance matrix. Again continuing the Poisson example from above:

      vcov.normal <- function(object)
        object$variance
      

2. Alternatively, you can write your own param.contrib() function. This is appropriate when:
   1. Your model has auxiliary parameters, such as 140#140 in the case of the Normal distribution.
   2. Your model performs some sort of correction to the coefficients or the variance-covariance matrix, which cannot be performed in either the coef.contrib() or the vcov.contrib() functions.
   3. Your model does not rely on asymptotic approximation to the log-likelihood. For Bayesian Markov-chain monte carlo models, for example, the param.contrib() function ( param.MCMCzelig() in this case) simply extracts the model parameters simulated in the model-fitting function.
   Continuing the Normal example,

   param.normal <- function(object, num = NULL, bootstrap = FALSE, 
                      terms = NULL) {
     if (!bootstrap) {
       par <- mvrnorm(num, mu = coef(object), Sigma = vcov(object))
       Beta <- parse.par(par, terms = terms, eqn = "mu")
       sigma2 <- exp(parse.par(par, terms = terms, eqn = "sigma2"))
       res <- cbind(Beta, sigma2)
     }
     else {
       par <- coef(object)
       Beta <- parse.par(par, terms = terms,  eqn = "mu")
       sigma2 <- exp(parse.par(par, terms = terms, eqn = "sigma2"))
       res <- c(coef, sigma2)
     }
     res
   }