Skip to content

marginal effect lmer

4 messages · Andigo, John Fox

Andigo
# 
Hi everybody,

I try to calculate and display the marginal effect(s) in a hierarchical
model using lmer. Here is my model:

m1<- lmer(vote2011~ Catholic + Attendance+ logthreshold + West +
Catholicproportion+
  (Catholic * Catholicproportion) + (Attendance*Catholicproportion) +
Catholicproportion?+ (Catholic *Catholicproportion?)+
  (Attendance* Catholicproportion?) + (1 + Attendance+ Catholic | canton),
data=dat2011, verbose = T)

I want to display the me of the individual level variable Catholic depending
on the contextual variable Catholicproportion (showing also the 95% ci). So
far I tried a bit with the "effects" package, but without success.

Does anybody know how to display that?

Thanks in advance for your help!

Andigo



--
View this message in context: http://r.789695.n4.nabble.com/marginal-effect-lmer-tp4637421.html
Sent from the R help mailing list archive at Nabble.com.
John Fox
# 
Dear Andigo,

You don't say what problems you encountered, and I don't know how to interpret the superscript 2s in your lmer() command, but if the intention is to fit a quadratic in Catholicproportion, then you can use poly(Catholicproportion, 2) in formulating the model. That way, effect() will be able to figure out the structure of the model. BTW, the "effects" computed by effect() are not what are commonly called "marginal effects," but rather fitted values under the model for particular combinations of predictors.

I hope this helps,
 John

------------------------------------------------
John Fox
Sen. William McMaster Prof. of Social Statistics
Department of Sociology
McMaster University
Hamilton, Ontario, Canada
http://socserv.mcmaster.ca/jfox/
	
On Mon, 23 Jul 2012 02:46:32 -0700 (PDT)
Andigo <Andreas.Goldberg at unige.ch> wrote:
Andigo
# 
Dear John,

yes, the superscript shall say that the quadratic term of Catholicproportion
is included as well (plus the interaction terms with it). In my dataset the
quadratic Catholicproportion is already included as a variable of its own,
so it should be fine or do I have to use the syntax with "poly..." you
mentioned? If so how do I have to rewrite the command for the whole model?

Out of your description of the package I just tried the example you
mentioned for a lmer model, so I tried:

plot(effect("Catholic:Catholicproportion", m1), grid=TRUE)

The result are ten (dk why 10 and not just 1) little plots, which look
rather linear and not the curvilinear pattern it should be. 
So far I always calculated the marginal effects in Stata, which works fine
for rather simple hierarchical models. But when the models become more
complicated (more random effects, more interaction terms,..) Stata often
fails to calculate the models, so that I cannot use it for displaying the
marginal effects. Thats why I try to find a solution in R to calculate the
marginal effects exactly the same way as I do in Stata. In Stata I followed
the syntax by Brambor, Clark and Golder
(https://files.nyu.edu/mrg217/public/interaction.html). Now I just wonder if
there is a way to calculate the marginal effects in R exactly the same way
as they do in Stata?

Best,
Andigo



--
View this message in context: http://r.789695.n4.nabble.com/marginal-effect-lmer-tp4637421p4637452.html
Sent from the R help mailing list archive at Nabble.com.
John Fox
# 
Dear Andigo,

On Mon, 23 Jul 2012 07:42:49 -0700 (PDT)
Andigo <Andreas.Goldberg at unige.ch> wrote:
It's hard to know where to begin, particularly since you deleted the earlier discussion. Retrieving that from my reply to you, the model you fit was
I don't really know what these variables are but I'll assume that all are numeric variables and not, e.g., 0/1 dummy regressors that you put in the model instead of factors in the same way that you generated the quadratic regressor for Catholicproportion -- both of these are bad ideas in R because they disguise the structure of the model. 

Anyway, I suppose that what you want in the fixed-effects part of the model is

logthreshold + West + (Catholic + Attendance)*poly(Catholicproportion, 2)

This will give you an orthogonal quadradic in Catholicproportion that interacts both with Catholic and Attendance and will allow effect() to figure out the structure of the model; if you want a raw quadratic (which would provide the same fit to the data but be somewhat less numerically stable), you could substitute poly(Catholicproportion, 2, raw=TRUE).
I'm not sure what you read, but what you got is what you asked for; with the respecified model, you could ask for the "Catholic:poly(Catholicproportion, 2)" effect, or more simply plot(allEffects(m1)) to plot all high-order fixed-effects terms in the model.
Catholic and Catholicproportion are presumably numeric variables, and as documented in ?effect, each will by default be set to 10 levels spanning the range of the variable. There are then 10*10 == 100 combinations of values of these variables at which the effect will be estimated. That you expected one plot suggests to me that you don't quite understand what effect() computes (though all the lines could be put on one plot by setting the argument multiline=TRUE -- ?plot.eff). The lines are straight because effect() couldn't read your mind to know that Catholicproportion2 and Catholicproportion are related, which is why using poly() instead is a good idea.
I already explained that effect() doesn't compute marginal effects. There may well be an R package that does, but I'm not aware of one. It should, however, be a simple matter to compute marginal effects from the lmer() results if you find marginal effects interesting.

Best,
 John