Average Marginal Effects via runmlwin

Welcome to the forum for runmlwin users. Feel free to post your question about runmlwin here. The Centre for Multilevel Modelling take no responsibility for the accuracy of these posts, we are unable to monitor them closely. Do go ahead and post your question and thank you in advance if you find the time to post any answers!

Go to runmlwin: Running MLwiN from within Stata >> http://www.bristol.ac.uk/cmm/software/runmlwin/
Forum rules
Post Reply
leap
Posts: 12
Joined: Mon Jul 11, 2011 11:29 am

Average Marginal Effects via runmlwin

Post by leap »

Dear everyone,

My question is the following:

How can I get average marginal effects and theirs standard errors from MLwiN via runmlwin?

A previous post on how to get predictions from MLwiN via runmlwin has already been answered. The post is extremely useful to obtain marginal effects at average and/or specific values of the covariates (http://www.cmm.bristol.ac.uk/forum/view ... a2112#p452).

Using a similar method than the one explained in the post above, I can easily obtained the average marginal effects but not the standard errors. Without being able to use of the command margins, my understanding is that I would need to calculate the Jacobian matrix. However, this is not so straightforward in a software such as STATA, which is not meant for mathematics but applied statistics.

I have been struggling with this issue and still cannot find a solution, so any help would be gladly appreciated.

Thank you,

Léa Pessin
Top
GeorgeLeckie
Site Admin
Posts: 432
Joined: Fri Apr 01, 2011 2:14 pm

Re: Average Marginal Effects via runmlwin

Post by GeorgeLeckie »

Hi Léa,

Calculating the AME is not so hard (see below). But I don't think there is any easy way to calculate the standard error of this average.

Sorry not to be of more help

Best wishes

George

Code: Select all

use "http://www.bristol.ac.uk/cmm/media/runmlwin/tutorial.dta", clear

generate pass = (normexam>0)

melogit pass standlrt i.girl

margins, dydx(girl)

runmlwin pass cons standlrt girl, ///
	level1(student:) ///
	discrete(distribution(binomial) link(logit) denom(cons) pql2) ///
	nopause

predictnl ame = ///
	invlogit(_b[cons] + _b[standlrt]*standlrt + _b[girl]*1) ///
	- invlogit(_b[cons] + _b[standlrt]*standlrt + _b[girl]*0)

tabstat ame, stat(mean)
with output...

Code: Select all

. use "http://www.bristol.ac.uk/cmm/media/runmlwin/tutorial.dta", clear

. 
. generate pass = (normexam>0)

. 
. melogit pass standlrt i.girl

Iteration 0:   log likelihood = -2274.4478  
Iteration 1:   log likelihood = -2271.5217  
Iteration 2:   log likelihood = -2271.5149  
Iteration 3:   log likelihood = -2271.5149  

Logistic regression                             Number of obs      =      4059

                                                Wald chi2(2)       =    734.81
Log likelihood = -2271.5149                     Prob > chi2        =    0.0000
------------------------------------------------------------------------------
        pass |      Coef.   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
    standlrt |    1.28272   .0475777    26.96   0.000      1.18947    1.375971
      1.girl |   .2207443   .0737735     2.99   0.003     .0761509    .3653376
       _cons |  -.0851353   .0575463    -1.48   0.139    -.1979239    .0276533
------------------------------------------------------------------------------

. 
. margins, dydx(girl)

Average marginal effects                          Number of obs   =       4059
Model VCE    : OIM

Expression   : Predicted mean, predict()
dy/dx w.r.t. : 1.girl

------------------------------------------------------------------------------
             |            Delta-method
             |      dy/dx   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
      1.girl |    .042072   .0140612     2.99   0.003     .0145126    .0696314
------------------------------------------------------------------------------
Note: dy/dx for factor levels is the discrete change from the base level.

. 
. runmlwin pass cons standlrt girl, ///
>         level1(student:) ///
>         discrete(distribution(binomial) link(logit) denom(cons) pql2) ///
>         nopause
 
MLwiN 2.29 multilevel model                     Number of obs      =      4059
Binomial logit response model
Estimation algorithm: IGLS, PQL2

Run time (seconds)   =       1.39
Number of iterations =          5
------------------------------------------------------------------------------
        pass |      Coef.   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
        cons |  -.0851353   .0575431    -1.48   0.139    -.1979177    .0276472
    standlrt |    1.28272   .0475705    26.96   0.000     1.189484    1.375956
        girl |   .2207442   .0737696     2.99   0.003     .0761584    .3653299
------------------------------------------------------------------------------


. 
. predictnl ame = ///
>         invlogit(_b[cons] + _b[standlrt]*standlrt + _b[girl]*1) ///
>         - invlogit(_b[cons] + _b[standlrt]*standlrt + _b[girl]*0)

. 
. tabstat ame, stat(mean)

    variable |      mean
-------------+----------
         ame |   .042072
------------------------
Best wishes

George
Top
GeorgeLeckie
Site Admin
Posts: 432
Joined: Fri Apr 01, 2011 2:14 pm

Re: Average Marginal Effects via runmlwin

Post by GeorgeLeckie »

Hi Léa,

Just had another thought.

If you can specify the same model using one of Stata's own multilevel modelling commands, but are avoiding Stata because it is computationally slow with your particular model and data, then you could use runmlwin to obtain the point estimates then plug them into the Stata multilevel modelling command as starting values. Then fit the model. The Stata multilevel modelling command should take far fewer iterations than normal and so converge far more quickly. So really this is just away of giving sensible starting values to Stata. You can then of course use the margins command however you want.

Example given below

George

Commands:

Code: Select all

* Load the data
use "http://www.bristol.ac.uk/cmm/media/runmlwin/tutorial.dta", clear

* Generate a binary response
generate pass = (normexam>0)

* Fit the model in runmlwin
runmlwin pass cons standlrt girl, ///
	level1(student:) ///
	discrete(distribution(binomial) link(logit) denom(cons) pql2) ///
	nopause

* Save the vector of runmlwin point estimates
matrix b1 = e(b)

* Specify model using melogit, but do not estiamte iit
melogit pass standlrt i.girl, noestimate

* Save the vector of starting values
matrix b = e(b)

* List the runmlwin point estimates
matrix list b1

* List the melogit starting values
matrix list b

* Replace the melogit starting values with the runmlwin point estimates
matrix b[1,1] = b1[1,2]
matrix b[1,2] = 0
matrix b[1,3] = b1[1,3]
matrix b[1,4] = b1[1,1]

* List the updated melogit starting values
matrix list b

* Note model fits in one iteration, because we provided starting values which were effectively the estimates
melogit pass standlrt i.girl, from(b)

* Note that had we fitted the model with no starting values, it would have taken three iterations
melogit pass standlrt i.girl

* Compute the average marginal effect
margins, dydx(girl)
Output:

Code: Select all

. * Load the data
. use "http://www.bristol.ac.uk/cmm/media/runmlwin/tutorial.dta", clear

. 
. * Generate a binary response
. generate pass = (normexam>0)

. 
. * Fit the model in runmlwin
. runmlwin pass cons standlrt girl, ///
>         level1(student:) ///
>         discrete(distribution(binomial) link(logit) denom(cons) pql2) ///
>         nopause
 
MLwiN 2.29 multilevel model                     Number of obs      =      4059
Binomial logit response model
Estimation algorithm: IGLS, PQL2

Run time (seconds)   =       1.48
Number of iterations =          5
------------------------------------------------------------------------------
        pass |      Coef.   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
        cons |  -.0851353   .0575431    -1.48   0.139    -.1979177    .0276472
    standlrt |    1.28272   .0475705    26.96   0.000     1.189484    1.375956
        girl |   .2207442   .0737696     2.99   0.003     .0761584    .3653299
------------------------------------------------------------------------------


. 
. * Save the vector of runmlwin point estimates
. matrix b1 = e(b)

. 
. * Specify model using melogit, but do not estiamte iit
. melogit pass standlrt i.girl, noestimate

Posting starting values:

Logistic regression                             Number of obs      =      4059

 ( 1)  [pass]0b.girl = 0
------------------------------------------------------------------------------
        pass |      Coef.  Legend
-------------+----------------------------------------------------------------
    standlrt |   1.170814  _b[standlrt]
      1.girl |   .1997406  _b[1.girl]
       _cons |  -.0626776  _b[_cons]
------------------------------------------------------------------------------
Note: The above coefficient values are starting values and not the result of a fully fitted model.

. 
. * Save the vector of starting values
. matrix b = e(b)

. 
. * List the runmlwin point estimates
. matrix list b1

b1[1,4]
           FP1:        FP1:        FP1:         OD:
          cons    standlrt        girl     bcons_1
y1  -.08513525   1.2827201   .22074416           1

. 
. * List the melogit starting values
. matrix list b

b[1,4]
          pass:       pass:       pass:       pass:
                        0b.          1.            
      standlrt        girl        girl       _cons
y1   1.1708138           0   .19974063  -.06267763

. 
. * Replace the melogit starting values with the runmlwin point estimates
. matrix b[1,1] = b1[1,2]

. matrix b[1,2] = 0

. matrix b[1,3] = b1[1,3]

. matrix b[1,4] = b1[1,1]

. 
. * List the updated melogit starting values
. matrix list b

b[1,4]
          pass:       pass:       pass:       pass:
                        0b.          1.            
      standlrt        girl        girl       _cons
y1   1.2827201           0   .22074416  -.08513525

. 
. * Note model fits in one iteration, because we provided starting values which were effectively the estimates
. melogit pass standlrt i.girl, from(b)

Iteration 0:   log likelihood = -2271.5149  
Iteration 1:   log likelihood = -2271.5149  

Logistic regression                             Number of obs      =      4059

                                                Wald chi2(2)       =    734.81
Log likelihood = -2271.5149                     Prob > chi2        =    0.0000
------------------------------------------------------------------------------
        pass |      Coef.   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
    standlrt |    1.28272   .0475777    26.96   0.000      1.18947    1.375971
      1.girl |   .2207443   .0737735     2.99   0.003     .0761509    .3653376
       _cons |  -.0851353   .0575463    -1.48   0.139    -.1979239    .0276533
------------------------------------------------------------------------------

. 
. * Note that had we fitted the model with no starting values, it would have taken three iterations
. melogit pass standlrt i.girl

Iteration 0:   log likelihood = -2274.4478  
Iteration 1:   log likelihood = -2271.5217  
Iteration 2:   log likelihood = -2271.5149  
Iteration 3:   log likelihood = -2271.5149  

Logistic regression                             Number of obs      =      4059

                                                Wald chi2(2)       =    734.81
Log likelihood = -2271.5149                     Prob > chi2        =    0.0000
------------------------------------------------------------------------------
        pass |      Coef.   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
    standlrt |    1.28272   .0475777    26.96   0.000      1.18947    1.375971
      1.girl |   .2207443   .0737735     2.99   0.003     .0761509    .3653376
       _cons |  -.0851353   .0575463    -1.48   0.139    -.1979239    .0276533
------------------------------------------------------------------------------

. 
. * Compute the average marginal effect
. margins, dydx(girl)

Average marginal effects                          Number of obs   =       4059
Model VCE    : OIM

Expression   : Predicted mean, predict()
dy/dx w.r.t. : 1.girl

------------------------------------------------------------------------------
             |            Delta-method
             |      dy/dx   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
      1.girl |    .042072   .0140612     2.99   0.003     .0145126    .0696314
------------------------------------------------------------------------------
Note: dy/dx for factor levels is the discrete change from the base level.
Top
leap
Posts: 12
Joined: Mon Jul 11, 2011 11:29 am

Re: Average Marginal Effects via runmlwin

Post by leap »

Hello George,

Thank you for your prompt and detailed answer.

Unfortunately (or fortunately thanks to the existence of MLwiN), I am fully using the functionality of MLwiN (MCMC estimation, etc.) so I cannot use STATA as an alternative.

Is there any hope that runmlwin will one day be compatible with margins? It is really a shame because margins and marginsplot really offer some nice ways to illustrate regression output.

Thank you again,

Léa
Top
GeorgeLeckie
Site Admin
Posts: 432
Joined: Fri Apr 01, 2011 2:14 pm

Re: Average Marginal Effects via runmlwin

Post by GeorgeLeckie »

Hi Léa,

I'm afraid we are not planning to make runmlwin compatible with margins. Not because it would not be useful, but because it would require a lot of programming.

If you are using MCMC estimation then things might actually be easier...

(1) Run the model for 1000 iterations and the save the chains as data. Or run it for 10000 iterations and set thinning to 10.

(2) Go back to the original data and expand the estimation sample by 1000 and then merge in the MCMC chains

(3) Calculate the AME in the usual way, but do it separately at each iteration of the MCMC chain

You will get 1000 values for the AME. The mean is your point estimate for the AME, the sd is your standard error.

Best wishes

George



variate normal.
Top
Post Reply

Return to “runmlwin user forum”