[R] Robust SE & Heteroskedasticity-consistent estimation

[Achim Zeileis]

On Tue, 11 May 2010, RATIARISON Eric wrote:

> Hi,
> I think my "meat"function is wrong.

Yes.

> But the "estfun(resMaxlik)" doesn't work;

Because you haven't supplied one.

It's as I wrote in my previous mails: You need to write a bread() and an 
estfun() method as specified in my previous mail. Details are in the 
sandwich vignette I already pointed you to. The meat() function does not 
have to be modified.
Z

> Error message :
> "Erreur dans UseMethod("estfun") :
>  pas de méthode pour 'estfun' applicable pour un objet de classe "c('maxLik', 'maxim', 'list')""
>
>
>
> -----Message d'origine-----
> De : RATIARISON Eric
> Envoyé : mardi 11 mai 2010 11:47
> À : 'Arne Henningsen'; Achim Zeileis; r-help at r-project.org; Ott-Siim Toomet
> Objet : RE: [R] Robust SE & Heteroskedasticity-consistent estimation
>
> Thank you all of you,;
> 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);
>
> n<-NROW(v)
> k<-NROW(resMaxlik$estimate)
>
> meat <-function(obj,nrow,npar,adjust=FALSE,...)
> { psi<-obj$gradient
>  k<-npar
>  n<-nrow
>  rval<-crossprod(as.matrix(psi))/n
>  if (adjust) rval <- n/(n-k)*rval
>  rval
> }
>
> M<-as.vector(meat(resMaxlik,n,k))
> B <-as.matrix(ginv(resMaxlik$hessian))
>
> vcovhc=B*M*B
>
> #new standard errors
> coeftest(resMaxlik,vcov.=vcovhc)
>
>
> do you think it's sound good ?
>
>
>
>
> -----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
>