[Rd] proposed simulate.glm method

Martin Maechler maechler at stat.math.ethz.ch
Sat Feb 14 20:10:26 CET 2009


>>>>> "NicLK" == Nicholas Lewin-Koh <nikko at hailmail.net>
>>>>>     on Sat, 14 Feb 2009 08:34:45 -0800 writes:

    NicLK> Hi, For extended glms such as gams, gnm or other
    NicLK> distributions such as negative binomial, would there
    NicLK> need to be a separate simulate method?  

Not necessarily,  as I said, the "glm"s are now also dealt with
in simulate.lm() and Heather more or less confirmed that this
gives correct results for "gnm" objects.

For gam(), I'd strongly expect the same to apply, but there
maybe sophisticated gam() models where the result is currently
not correct.  That's, BTW, also true for  
    simulate(lm(...., weights), ...)

    NicLK>  Or, could the current framework, rather than
    NicLK> stopping with an error look for the appropriate model
    NicLK> matrix, coefficients, distribution function and
    NicLK> family object to simulate from?

What do you mean?
A situation where there's no supported 'family'
or a situation where  predict(<obj>) does not work as it's
supposed in the current framework,
or ????

If there are such cases, we'd have to consider them together
with the corresponding package author.  It may often make sense
fthen that the author changes his methods {predict(), ..} such
that the (now) extended simulate.lm() will work automatically;
Alternatively, the author can provide  simulate.<myclass>().

But I'm not sure I'm answering the question you've asked..
Martin

    NicLK> Nicholas


    >> Message: 9 Date: Fri, 13 Feb 2009 21:27:57 +0100 From:
    >> Martin Maechler <maechler at stat.math.ethz.ch> Subject: Re:
    >> [Rd] proposed simulate.glm method To: Heather Turner
    >> <Heather.Turner at warwick.ac.uk> Cc: r-devel at r-project.org,
    >> Martin Maechler <maechler at stat.math.ethz.ch> Message-ID:
    >> <18837.55245.15158.29378 at cmath-5.math.ethz.ch>
    >> Content-Type: text/plain; charset=us-ascii
    >> 
    >> Thank you, Heather and Ben,
    >> 
    >> >>>>> "HT" == Heather Turner
    >> <Heather.Turner at warwick.ac.uk> >>>>> on Fri, 13 Feb 2009
    >> 15:52:37 +0000 writes:
    >> 
    HT> Yes, thanks to Ben for getting the ball rolling. His
    HT> code was more streamlined than mine, pointing to further
    HT> simplifications which I've included in the extended
    HT> version below.
    >> 
    HT> The code for the additional families uses functions from
    HT> MASS and SuppDists - I wasn't sure about the best way to
    HT> do this, so have just used :: for now.
    >> 
    HT> It appears to be working happily for both glm and gnm
    HT> objects (no gnm-specific code used).
    >> 
    HT> Best wishes,
    >> 
    HT> Heather
    >> 
    >> [....]
    >> 
    >> I have now followed Brian Ripley's suggetion to just
    >> extend simulate.lm() to also deal with "glm" objects, but
    >> using Heather's suggestions for the different families;
    >> I've just commited src/library/stats/R/lm.R with the new
    >> code.  (get it from svn.r-project.org/R/trunk/ or this
    >> night's R-devel tarball).
    >> 
    >> One difference to your propsal: Instead of just
    >> object$fitted , the code is using fitted(object)
    >> ... something which should properly to the na.action
    >> used.
    >> 
    >> For the (MASS and) SuppDists package requirement, I'm
    >> using the following
    >> 
    >> if(is.null(tryCatch(loadNamespace("SuppDists"), error =
    >> function(e) NULL))) stop("Need CRAN package 'SuppDists'
    >> for 'inverse.gaussian' family")
    >> 
    >> 
    >> I've not yet updated the help page for simulate(), and
    >> have only tested relatively few cases for binomial,
    >> poisson and Gamma.  I've wanted to expose this to you, so
    >> you can provide more feedback and possibly even a patch
    >> to
    >> svn.r-project.org/R/trunk/src/library/stats/man/simulate.Rd
    >> 
    >> Martin
    >> 
    >> 
    >> 

______________________________________________
    NicLK> R-devel at r-project.org mailing list
    NicLK> https://stat.ethz.ch/mailman/listinfo/r-devel



More information about the R-devel mailing list