simulate                package:stats                R Documentation

_S_i_m_u_l_a_t_e _R_e_s_p_o_n_s_e_s

_D_e_s_c_r_i_p_t_i_o_n:

     Simulate one or more responses from the distribution corresponding
     to a fitted model object.

_U_s_a_g_e:

     simulate(object, nsim, seed, ...)

_A_r_g_u_m_e_n_t_s:

  object: an object representing a fitted model.

    nsim: number of response vectors to simulate.  Defaults to '1'.

    seed: an object specifying if and how the random number generator
          should be initialized ('seeded').
           For the "lm" method, either 'NULL' or an integer that will
          be used in a call to 'set.seed' before simulating the
          response vectors.  If set, the value is saved as the '"seed"'
          attribute of the returned value.  The default, 'NULL' will
          not change the random generator state, and return
          '.Random.seed' as the '"seed"' attribute, see 'Value'. 

     ...: additional optional arguments.

_D_e_t_a_i_l_s:

     This is a generic function.  Consult the individual modeling
     functions for details on how to use this function.

     Package 'stats' has a method for '"lm"' objects which is used for
     'lm' and 'glm' fits.  There is a method for fits from 'glm.nb' in
     package 'MASS', and hence the case of negative binomial families
     is not covered by the '"lm"' method.

     The methods for linear models fitted by 'lm' or 'glm(family =
     "gaussian")' assume that any weights which have been supplied are
     inversely proportional to the error variance.  For other GLMs the
     (optional) 'simulate' component of the 'family' object is
     used-there is no appropriate simulation method for 'quasi' models
     as they are specified only up to two moments.

     For binomial and Poisson GLMs the dispersion is fixed at one. 
     Integer prior weights w_i can be interpreted as meaning that
     observation i is an average of w_i observations, which is natural
     for binomials specified as proportions but less so for a Poisson,
     for which prior weights are ignored with a warning.

     For a gamma GLM the shape parameter is estimated by maximum
     likelihood (using function 'gamma.shape' in package 'MASS').  The
     interpretation of weights is as multipliers to a basic shape
     parameter, since dispersion is inversely proportional to shape.

     For an inverse gauasian GLM the model assumed is IG(mu_i, lambda
     w_i) (see <URL:
     http://en.wikipedia.org/wiki/Inverse_Gaussian_distribution>) where
     lambda is estimated by the inverse of the dispersion estimate for
     the fit.  The variance is mu_i^3/(lambda w_i) and hence inversely
     proprortional to the prior weights.  The simulation is done by
     function 'rinvGaussian' from the 'SuppDists' package, which must
     be installed.

_V_a_l_u_e:

     Typically, a list of length 'nsim' of simulated responses.  Where
     appropriate the result can be a data frame (which is a special
     type of list).

     For the '"lm"' method, the result is a data frame with an
     attribute '"seed"' containing the 'seed' argument if not 'NULL'
     with '"kind"' attributs the vaue of 'as.list(RNGkind())',
     otherwise (the default) the value of '.Random.seed' before the
     simulation was started.

_S_e_e _A_l_s_o:

     'fitted.values' and 'residuals' for related methods; 'glm', 'lm'
     for model fitting.

     There are further examples in the 'simulate.R' tests file in the
     sources for package 'stats'.

_E_x_a_m_p_l_e_s:

     x <- 1:5
     mod1 <- lm(c(1:3,7,6) ~ x)
     S1 <- simulate(mod1, nsim = 4)
     ## repeat the simulation:
     .Random.seed <- attr(S1, "seed")
     identical(S1, simulate(mod1, nsim = 4))

     S2 <- simulate(mod1, nsim = 200, seed = 101)
     rowMeans(S2) # should be about
     fitted(mod1)

     ## repeat identically:
     (sseed <- attr(S2, "seed")) # seed; RNGkind as attribute
     stopifnot(identical(S2, simulate(mod1, nsim = 200, seed = sseed)))

     ## To be sure about the proper RNGkind, e.g., after
     RNGversion("2.7.0")
     ## first set the RNG kind, then simulate
     do.call(RNGkind, attr(sseed, "kind"))
     identical(S2, simulate(mod1, nsim = 200, seed = sseed))

     ## Binomial GLM examples
     yb1 <- matrix(c(4,4,5,7,8,6,6,5,3,2), ncol = 2)
     modb1 <- glm(yb1 ~ x, family = binomial)
     S3 <- simulate(modb1, nsim = 4)
     # each column of S3 is a two-column matrix.

     x2 <- sort(runif(100))
     yb2 <- rbinom(100, prob = plogis(2*(x2-1)), size = 1)
     yb2 <- factor(1 + yb2, labels = c("failure", "success"))
     modb2 <- glm(yb2 ~ x2, family = binomial)
     S4 <- simulate(modb2, nsim = 4)
     # each column of S4 is a factor