Robust and clustered standard errrors

I'm kinda hijacking this issue because I think I have discovered a problem.

So,

Is there any easy way to get the sandwich standard errors, or cluster sandwich standard errors after estimation, or do I have to calculate them myself using the gradient and hessian in the results?

Yes, you can use the sandwich and lmtest like so:

m <- maxlik(...)
coeftest(m, vcov = vcovCL, cluster = data$cluster_variable)

You can also use a formula like ~ cluster_variable. vcovCL is for clustered standard errors. Look at the documentation of sandwich for all the options. The reason it works is that maxlik very nicely reports the bread and the meat that sandwich uses.

That being said... I tried to use bootstrapped standard errors (vcovBS). It doesn't work, and I think the reason is that a maxlik model does not report nobs, i.e. if I type nobs(m) I always get 0.

@otoomet would it be possible to make maxlik report nobs too?

This is important because, when looking at the clusters and the data, vcovBS will do

vcovBS.default <- function(x, cluster = NULL, R = 250, start = FALSE, ..., fix = FALSE, use = "pairwise.complete.obs", applyfun = NULL, cores = NULL)
{
  ## set up return value with correct dimension and names
  cf <- coef(x)
  k <- length(cf)
  n <- nobs0(x)

where nobs0 is defined as

nobs0   <- function(x, ...) {
  nobs1 <- if("stats4" %in% loadedNamespaces()) stats4::nobs else stats::nobs
  nobs2 <- function(x, ...) NROW(residuals(x, ...))
  rval <- try(nobs1(x, ...), silent = TRUE)
  if(inherits(rval, "try-error") | is.null(rval)) rval <- nobs2(x, ...)
  return(rval)
}

For a maxlik model this function returns 0, as nobs.

The problem is that we get to

if(NROW(cluster) != n) stop("number of observations in 'cluster' and 'nobs()' do not match")

Since n = 0, it breaks.

For some reason maxlik reports the number of observations with nObs, and nObs(m) reports the same number of observations as NROW(cluster). So if maxlik would also report to nobs, vcovBS would work fine.

arne-henningsen / maxLik

Robust and clustered standard errrors #38