[R] Robust SE & Heteroskedasticity-consistent estimation

[Arne Henningsen]

On 13 May 2010 00:41, RATIARISON Eric <ERATIARISON at monceauassurances.com> wrote:
> Hi, here my new version:
> I submit you my case:
> ( Pseudo likehood for exponential family with offset )
>
> loglik <- function(param)  {
>   b<-param
>   m=as.vector(of+z%*%b)
>   ll <- sum(v*m-exp(m))   }
>
> gradlik <- function(param) {
>     b<-param
>     m=as.vector(of+z%*%b)
>     gg<-(v-exp(m))*z      }
>
> hesslik <- function(param) {
>     b<-param
>     m=as.vector(of+z%*%b)
>     hh <- -t(exp(m)*z)%*%z }
>
> resMaxlik <- maxLik(loglik,grad=gradlik,hess=hesslik,start=sc,method="nr",tol=1e-4);
>
> resMaxlik$offset<-of
> resMaxlik$x<-z
> resMaxlik$y<-v
>
>
>
> estfun <- function(obj,...)
> {
> m=as.vector(obj$offset+obj$x%*%obj$estimate)
> (obj$y-exp(m))*obj$x
> }
>
> M <- function(obj, adjust = FALSE, ...)
> {
> psi <- estfun(obj)
> k <- NCOL(psi)
> n <- NROW(psi)
> rval <- crossprod(as.matrix(psi))/n
> if(adjust) rval <- n/(n - k) * rval
> rval
> }
>
>
> B <- function(obj, ...)
> { as.matrix(ginv(obj$hessian))
> }
>
>
> So S <-B(resMaxlik)*M(resMaxlik)*B(resMaxlik) directly. It's ok.
>
> But I call sandwich function like this : sandwich(resMaxlik,meat.=M,bread.=B)
> It returns a error message : "Erreur dans UseMethod("estfun") :
>  pas de méthode pour 'estfun' applicable pour un objet de classe "c('maxLik', 'maxim', 'list')"

I have added an estfun() method and a bread() method for objects of
class "maxLik" so that you can use "sandwich" now.

Please note:

a) The methods can be applied only if maxLik() was called with
argument "grad" equal to a gradient function or (if no gradient
function is specified) argument "logLik" equal to a log-likelihood
function that return the gradients or log-likelihood values,
respectively, for *each observation* (i.e. in the BHHH form).

b) The package with the new features is not yet available on CRAN but
on R-Forge: http://r-forge.r-project.org/scm/?group_id=254

Please let me know if you experience any problems.

/Arne


> -----Message d'origine-----
> De : Arne Henningsen [mailto:arne.henningsen at googlemail.com]
> Envoyé : mardi 11 mai 2010 08:25
> À : Achim Zeileis; RATIARISON Eric; r-help at r-project.org; Ott-Siim Toomet
> Objet : Re: [R] Robust SE & Heteroskedasticity-consistent estimation
>
> On 11 May 2010 00:52, Achim Zeileis <Achim.Zeileis at uibk.ac.at> wrote:
>> On Mon, 10 May 2010, RATIARISON Eric wrote:
>>> I'm using maxlik with functions specified (L, his gradient & hessian).
>>>
>>> Now I would like determine some robust standard errors of my estimators.
>>>
>>> So I 'm try to use vcovHC, or hccm or robcov for example
>>>
>>> but in use one of them with my result of maxlik, I've a the following
>>> error message :
>>>
>>> Erreur dans terms.default(object) : no terms component
>>>
>>> Is there some attributes to give to maxlik objet for "fitting" the call
>>> of vcovHC?
>>
>> This is discussed in
>>  vignette("sandwich-OOP", package = "sandwich")
>> one of the vignettes accompanying the "sandwich" package that provides the
>> vcovHC() function. At the very least, you need an estfun() method which
>> extracts the gradient contributions per observation. Then you need a bread()
>> function, typically based on the observed Hessian. Then you can compute the
>> basic sandwich() estimators.
>
> Is it possible to implement the estfun() method and the bread()
> function in the maxLik package so that vcovHC() can be easily used by
> all users of the maxLik package? If yes: should we (Eric, Achim, Ott,
> Arne) implement this feature together?
>
> /Arne
>
> --
> Arne Henningsen
> http://www.arne-henningsen.name
>



-- 
Arne Henningsen
http://www.arne-henningsen.name