Multilevel/Hierarchical models and clustered-robust standard errors

[John Dohn]

Multilevel/hierarchical model with clustered-robust standard errors. I have reviewed various posts on this topic, including this post pointing at a cross-nested hierarchical specification, this discussion on hierarchical probit models mentioning that the panel variable must be nested within the cluster variable, and this post showing how mixing hierarchical models and clustered-robust standard errors is usually a bad idea. However, I believe I am in a situation where a multilevel model with clustered-robust standard errors may be justified.

The features of my data structure are the following:
- Individual-level data, in particular, a repeated cross-section dataset (i.e. not a panel) across many countries (variable COUNTRY)
- Dichotomous dependent variable (DV)
- Multiple scores for each individual for one independent variable (IV1 below), to various possible choices (P, ranging from 4 to 15, depending on country and year)
- ther variables defined at the individual level (e.g. see IV2)

The following table provides an idea of the data structure at hand:

ID	COUNTRY	DV	IV1	IV2	P
1	1	0	2	1	1
1	1	0	3	1	2
1	1	1	10	1	3
2	1	0	2	3	1
2	1	1	7	3	2

In normal circumstances, I would estimate a three-level hierarchical logit model with IDxP-level observations nested in ID, nested in a COUTRYxYEAR indicator.
However, I only have a few observations within individuals (always less than 16) and this would not satisfy the asymptotic assumptions needed to estimate random effects in a frequentist setting.

Thus, I would like to estimate a hierarchical logit model with observations nested at the COUNTRY-YEAR level (adding a random intercept to account for this), but then I also want to account for the fact that observations are repeated within individuals introducing clustered-robust standard errors at the ID level. This way, I want to avoid unrealistic standard errors for the IV1 coefficient.

My code:

melogit DV IV1 IV2 || COUNTRY:, vce(cluster ID)

gets the following error:

highest-level groups are not nested within ID
r(459);

Questions:
1) Is there some flaw in the model that I cannot see?
2) How can I estimate such a model in STATA?

Thank you all.
John

[Clyde Schechter]

You cannot use vce(cluster ID) when you have a COUNTRY: level in your model because COUNTRY: is not nested within ID. It is the other way around.

I also don't grasp your concern that you only have a modest number of observations within individual (by which I assume you mean ID). That's not a problem if you have a large number of individuals. You will still have plenty of degrees of freedom for estimating the residual variance.

You repeatedly refer to COUNTRY-YEAR, but your data has no YEAR variable, so I don't know what you are talking about there. Similarly you say that variable P ranges from 4 to 15, but in your example it takes on only the values 1, 2, or 3.

I'm a little hesitant to go out on a limb here because, as you can see, I don't think I fully understand your data set. But off hand it looks like I would do this as:

melogit DV IV1 IV2 || COUNTRY: || ID:, vce(cluster COUNTRY)

By using vce(cluster COUNTRY) you are accounting for the non-independence of observations both within COUNTRY and within ID. And since COUNTRY: is also the top-level in the model, Stata will not object to this choice. I think it is important to include the ID: level in the model, and I would regard it as mis-specified without it. What matters here is not whether the number of observations per ID is small, but whether the number of distinct values of ID itself is small. If you have only a handful of people, then take out the ID: level and put i.ID into the fixed effects instead. But one way or another you need to account for the fact that there are likely to be individual-level effects.