I know this kind of question has been asked before, but I can't find anything that clearly elucidates the issue. What is the 'right' way to model the follow situation:
Let's say I pair two people up randomly (on the basis of their sex) such that there are four possible groups: MM, FF, MF, and FM. I have them play some behavioral economics game, and I am interested in their before-and-after differences on some variable, let's say their happiness.
In R, I can think of several potential ways to model this relationship. Assume treatment is a categorical variable containing the above four categories:
lm(happy_score ~ treatment + education + income)
However, this doesn't account for the fact that, within every 2-person team, there will be correlations unaccounted for in the above specification.
So, we can add a fixed effect for team (a unique identifier for each 2-person team)
lm(happy_score ~ treatment + education + income + team)
This seems to not be getting at what I want, as it will compute a separate slope for each team, when I really want to just 'control' for the correlation that is likely to be shared between two individuals on the same team.
I've also been told to address this issue we can cluster standard errors at the team level, so:
happy_score ~ treatment + education + income,
data = data,
clusters = team,
se = "stata"
But I'm not sure what this is doing that is different from adding a fixed effect. Should I have both fixed effects and clustered standard errors?
I've also been pushed in the direction of mixed models, but I don't understand what they are substantively doing that is different from some combination of fixed effects and interactions:
lmer(happy_score ~ treatment + education + income
(1|team), data = data)
Are there intuitive ways of thinking about each of these options?