[R-sig-ME] Cluster-robust SEs & random effects -- seeking some clarification
J.D. Haltigan
jh@|t|g@ @end|ng |rom gm@||@com
Thu Aug 4 05:46:30 CEST 2022
Thought I would bump my last post to see if anyone might weigh in on my
more general statements about re: unions (clusters). Would be most
appreciated as I continue to wrap my head around some things.
Thanks in advance for any thoughts. I appreciate your time.
On Sat, Jul 30, 2022 at 8:43 PM J.D. Haltigan <jhaltiga using gmail.com> wrote:
> And to take one more step: the inference is that the single
> population-level treatment effect is drawing from a 'random sample' of
> unions (clusters). That is, one can not assume that union's are the same.
> They are not interchangeable. They are drawn from a random population. This
> is the point of the exercise for me. If we remove the random effect for
> union in the model I shared, we end up with a model in which only the fixed
> effects of pairID (treatment-control pairs for each pair of treatment
> control villages) are estimated (albeit using cluster-robust SEs). So, if
> by adding that random effect for union the treatment intervention is no
> longer significant (as opposed to a model in which there is no random
> effect of union modeled), what is that telling us? That some of the between
> cluster (union) variance in intercepts is contributing to variation in the
> response variable, yes?
>
> I realize you have not read the paper, nor are necessarily interested in
> this discourse, but any remarks are greatly appreciated.
>
>
> On Sat, Jul 30, 2022 at 8:36 PM Ben Bolker <bbolker using gmail.com> wrote:
>
>> Yes.
>>
>> On 2022-07-30 8:12 p.m., J.D. Haltigan wrote:
>> > Thanks, Ben. So in the model you remarked on, would that be a
>> > 'random-intercepts only' model?
>> >
>> >
>> > On Sat, Jul 30, 2022 at 7:53 PM Ben Bolker <bbolker using gmail.com
>> > <mailto:bbolker using gmail.com>> wrote:
>> >
>> > I haven't been following the whole thread that carefully, but I
>> want to
>> > emphasize that
>> >
>> > posXsymp~treatment + pairID + (1 | union)
>> >
>> > is *not*, by any definition I'm familiar with, a "random-slopes
>> model";
>> > that is, it only estimates a single population-level treatment
>> > effect/doesn't allow the effect of treatment to vary across groups
>> > defined by 'union'. You would need a random-effect term of the form
>> > (treatment | union).
>> >
>> > Reasons why you might *not* want to do this:
>> >
>> > * if treatment only varies across and not within levels of union
>> > ("union is nested within treatment" according to some terminology),
>> > then
>> > this variation is unidentifiable
>> > * maybe you have decided that you don't have enough data/want a
>> more
>> > parsimonious model.
>> >
>> > Schielzeth and Forstmeier, among many others (this is the
>> example I
>> > know of), have cautioned about the consequences of leaving out
>> > random-slopes terms.
>> >
>> > Schielzeth, Holger, and Wolfgang Forstmeier. “Conclusions beyond
>> > Support: Overconfident Estimates in Mixed Models.” Behavioral
>> Ecology
>> > 20, no. 2 (March 1, 2009): 416–20.
>> > https://doi.org/10.1093/beheco/arn145
>> > <https://doi.org/10.1093/beheco/arn145>.
>> >
>> >
>> > On 2022-07-30 7:43 p.m., J.D. Haltigan wrote:
>> > > Addendum:
>> > >
>> > > It just occurred to me on my walk that I think I am getting a bit
>> > lost in
>> > > some of the differences in nomenclature across scientific silos.
>> > In the
>> > > original model that they specified, which treated the 'pairID'
>> > variable as
>> > > a control variable for which they controlled for 'fixed effects'
>> of
>> > > control/treatment villages (in their own language in the paper)
>> using
>> > > cluster-robust SEs, I think this is indeed a 'random-intercepts
>> > only' model
>> > > in the language of Hamaker et al. They implement the 'absorb'
>> > command in
>> > > STATA which I believe aggregates across the pairIDs to generate
>> an
>> > > 'omnibus' F-test of sorts for the pairID variable (in the ANOVA
>> > > nomenclature). I say this as when I specify the pairID variable
>> > in the lmer
>> > > model I shared (or in a fixest model I conducted to replicate the
>> > original
>> > > Abalauck results in R), I get the estimates for all the pairs
>> > (i.e., there
>> > > is no way to aggregate across them--though I think formally the
>> > models are
>> > > the same if we are unconcerned about any one pairID
>> > [treatment/control
>> > > village pair].
>> > >
>> > > So, in the lmer model I shared where I specify a specific random
>> > effects
>> > > term for the 'cluster' variable, I think this indeed is allowing
>> > for random
>> > > slopes across the clusters which implies the treatment effect may
>> > vary
>> > > across the clusters (and we might anticipate it will for various
>> > reasons I
>> > > can elaborate on). More generally: we are generalizing to *any*
>> > universe of
>> > > villages (say in the entire world) where the treatment
>> > intervention (masks)
>> > > may vary across villages. This is the crux of invoking the random
>> > effects
>> > > model (i.e., random slopes model).
>> > >
>> > > I realize this is a mouthful, but I think the way these terms
>> (e.g.,
>> > > random/fixed effects models etc.) are used across disciplines
>> > makes things
>> > > a bit confusing.
>> > >
>> > > On Sat, Jul 30, 2022 at 5:25 PM J.D. Haltigan <
>> jhaltiga using gmail.com
>> > <mailto:jhaltiga using gmail.com>> wrote:
>> > >
>> > >> This is a very helpful walkthrough, James. My responses are
>> > italicized
>> > >> under yours to maintain thread readability. The key is
>> > Generalizability
>> > >> here and (as I also note in my last reply) the idea is to
>> > Generalize to a
>> > >> universe of "any villages or clusters." That is, the target
>> > population we
>> > >> are generalizing to is *any* random population.
>> > >>
>> > >> On Sat, Jul 30, 2022 at 3:01 PM James Pustejovsky
>> > <jepusto using gmail.com <mailto:jepusto using gmail.com>>
>> > >> wrote:
>> > >>
>> > >>> Hi J.D.,
>> > >>> A few comments/reactions inline below.
>> > >>> James
>> > >>>
>> > >>> On Wed, Jul 27, 2022 at 5:37 PM J.D. Haltigan
>> > <jhaltiga using gmail.com <mailto:jhaltiga using gmail.com>> wrote:
>> > >>>
>> > >>>> ...
>> > >>>>
>> > >>> In the original investigation, the authors did not invoke a
>> random
>> > >>>> effects model (but did use the pairIDs to control for fixed
>> > effects as
>> > >>>> noted and with robust SEs). Thus, in the original
>> > investigation there was
>> > >>>> *no* specification of a random effects model for the 'cluster'
>> > variable. We
>> > >>>> know from some other work there were some biases in village
>> > mapping and
>> > >>>> other possible sources of between-cluster variation that
>> might be
>> > >>>> anticipated to have influence--at the random intercepts
>> > level--so we are
>> > >>>> looking into how specifying 'cluster' as a random effect might
>> > change the
>> > >>>> fixed effects estimates for the treatment intervention effect.
>> > In the
>> > >>>> Hamaker et al. language, it is indeed a 'random intercepts'
>> > only model.
>> > >>>>
>> > >>>
>> > >>> I don't follow how using a random intercepts model improves the
>> > >>> generalizability warrant here. The random intercepts model is
>> > essentially
>> > >>> just a re-weighted average of the pair-specific effects in the
>> > original
>> > >>> analysis, where the weights are optimally efficient if the
>> model is
>> > >>> correctly specified. That last clause carries a lot of weight
>> > here--correct
>> > >>> specification means 1) treatment assignment is unrelated to the
>> > random
>> > >>> effects, 2) the treatment effect is constant across clusters,
>> 3)
>> > >>> distributional assumptions are valid (i.e., homoskedasticity at
>> > each level
>> > >>> of the model).
>> > >>>
>> > >>> If the effects are heterogeneous, then I would think that
>> including
>> > >>> random slopes on the treatment indicator would provide a better
>> > basis for
>> > >>> generalization. But even then, the warrant is still pretty
>> > vague---what is
>> > >>> the hypothetical population of villages from which the observed
>> > villages
>> > >>> are sampled?
>> > >>>
>> > >>
>> > >> *In the most basic model (without baseline controls) the model
>> > takes the
>> > >> form: myModel = lmer(posXsymp~treatment + pairID + (1 | union),
>> > data =
>> > >> myData). I believe--correct me if I am wrong--that this
>> reflects a
>> > >> random-intercepts only model, but I may be mistaken. If I am,
>> > and this is
>> > >> allowing for random slopes on the treatment indicator, then I
>> > will need to
>> > >> rethink my statements. *
>> > >>
>> > >>>
>> > >>>
>> > >>>> Given this, however, does it also make sense to include the
>> > cluster
>> > >>>> robust SEs for the fixed effects which would account for
>> possible
>> > >>>> heterogeneity of treatment effects (i.e., slopes) across
>> > clusters?s
>> > >>>>
>> > >>>> If you're committed to the random intercepts model, then yes I
>> > think so
>> > >>> because using cluster robust SEs at least acknowledges the
>> > possibility of
>> > >>> heterogeneous treatment effects.
>> > >>>
>> > >>
>> > >> *If the above model does allow for both random intercepts and
>> > slopes, then
>> > >> perhaps the use of cluster robust SEs is redundant in some sense
>> > since the
>> > >> random slopes would be modeling the heterogeneity in treatment
>> > effects?*
>> > >>
>> > >>>
>> > >>>
>> > >>>
>> > >>>> Bottom line: in their original analyses, clusters are seen as
>> > >>>> interchangeable from a conceptual perspective (rather than
>> > drawn from a
>> > >>>> random universe of observations). When one scales up evidence
>> > to a universe
>> > >>>> of observations that are random (as they would be in the
>> > intended universe
>> > >>>> of inference in the real-world), then we are better
>> > positioned, I think, to
>> > >>>> adjudicate whether the mask intervention effect is
>> 'practically
>> > >>>> significant' (in addition to whether the focal effect remains
>> > marginally
>> > >>>> significant from a frequentist perspective).
>> > >>>>
>> > >>> As noted above, this argument is a bit vague to me. If there's
>> > concern
>> > >>> about generalizability, then my first question would be: what
>> > is the target
>> > >>> population to which you are trying to generalize?
>> > >>>
>> > >>
>> > >> *Essentially, the target population we are trying to generalize
>> > to is a
>> > >> random selection of villages. Any random selection of villages.
>> > In other
>> > >> words, villages should not be seen as interchangeable. We are
>> > interested in
>> > >> whether the effects generalize to any randomly selected
>> village. *
>> > >>
>> > >>>
>> > >>>
>> > >>
>> > >
>> > > [[alternative HTML version deleted]]
>> > >
>> > > _______________________________________________
>> > > R-sig-mixed-models using r-project.org
>> > <mailto:R-sig-mixed-models using r-project.org> mailing list
>> > > https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>> > <https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models>
>> >
>> > --
>> > Dr. Benjamin Bolker
>> > Professor, Mathematics & Statistics and Biology, McMaster University
>> > Director, School of Computational Science and Engineering
>> > (Acting) Graduate chair, Mathematics & Statistics
>> >
>> > _______________________________________________
>> > R-sig-mixed-models using r-project.org
>> > <mailto:R-sig-mixed-models using r-project.org> mailing list
>> > https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>> > <https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models>
>> >
>>
>> --
>> Dr. Benjamin Bolker
>> Professor, Mathematics & Statistics and Biology, McMaster University
>> Director, School of Computational Science and Engineering
>> (Acting) Graduate chair, Mathematics & Statistics
>>
>
[[alternative HTML version deleted]]
More information about the R-sig-mixed-models
mailing list