how to plot cefficients of quantiel regression for panel data

[bruce fan]

Hi guys,

I am using the following code to run a quantile regression for panel data:

forvalues i=0.1(0.1)0.9{
qregpd gini carrier_mmc2 lnhhi, quantile(`i') id(id) fix(quarter) optimize(mcmc) ///
draws(10000) burn(1000) arate(0.5)
}

two questions: 1) I want to include id fixed effects and time(quarter) fixed effects, can I do this with the above codes?

2) I want to plot coefficients for regressor carrier_mmc2 by quantiles with 95% confidence interval, I am wondering if there are any codes can do this?


Thank you! Any suggest will be appreciated!

[Clyde Schechter]

Caution: I do not know -qregpd-; it is not a native Stata command and it is the etiquette of this forum that you identify its source and credit its author when referring to it. Anyway, I can't comment on your first question. But as for graphing the coefficients as a function of the quantile, assuming that like official Stata estimation commands -qregpd- returns a _b[] matrix with its coefficients, you can do it like this:

capture postutil clear
tempfile holding
postfile handle quantile coefficient using `holding'

forvalues i=0.1(0.1)0.9{
    qregpd gini carrier_mmc2 lnhhi, quantile(`i') id(id) fix(quarter) optimize(mcmc) ///
        draws(10000) burn(1000) arate(0.5)
    post handle (`i') (_b[carrier_mmc2])
}
postclose handle

use `holding', clear
graph twoway connected coefficient quantile

As an aside, -forvalues i = 0.1(0.1)0.9- is not a good programming practice. In this particular instance it will do you no harm, but as a general principle you should avoid using floating point numbers as increments (unless the non-integer part is an exact negative power of 2). The reason is that they are not representable exactly in finite precision binary, just as 1/3 has no exact representation in finite-precision decimals. So Stata is not actually using 0.1 as the increment; it is using the closest approximation to 0.1 it can find that has a finite binary representation. Consequently, if you had a longer series of iterations, the error caused by that approximation would add up, and you would eventually be substantially off from your intended series of numbers, and perhaps the stopping point (0.9) would be missed entirely. The safer way to code this is:

forvalues i = 1(1)/9 {
    local ii = `i'/10
    // command(s) involving reference to `ii' here
}

As I say, in this particular case, this problem does not bite and what you wrote will work just fine. But you should get in the habit of doing it the safer way, because someday you might trip over this and be utterly baffled by what went wrong.

[bruce fan]

Hi, Clyde (@Clyde Schechter)

Thank you very much for your suggestions. I changed my codes as exactly you said, including the floating point numbers. One more question: I cannot get the confidence interval on the graph. Do you know any way to do it? Thanks!

[Clyde Schechter]

Well, this is a little harder for me to answer as I do not know -qregpd-. Official Stata estimation commands leave behind a matrix r(table), and the confidence limits are found in the fifth and sixth rows of that matrix. So, if -qregpd- does that, you can modify the code as follows:

capture postutil clear
tempfile holding
postfile handle quantile coefficient lb ub using `holding'

forvalues i=1/9 {
    local ii = `i'/10
    qregpd gini carrier_mmc2 lnhhi, quantile(`ii') id(id) fix(quarter) optimize(mcmc) ///
        draws(10000) burn(1000) arate(0.5)
    matrix M = r(table)
    post handle (`ii') (_b[carrier_mmc2]) (M[1, 5]) (M[1,6])
}
postclose handle

use `holding', clear
graph twoway connected coefficient lb ub quantile

(Changes related to grabbing confidence bounds in bold face.)

Now, if -qregpd- doesn't leave r(table) behind, then the alternative is to grab the standard error out of _se[], analogous to the way you grabbed the coefficient out of _b[] and then calculate the confidence bounds yourself from the coefficient and its standard error (I'm assuming here that the confidence bounds come from a normal or t-distribution in -qregpd-.)

Note also that you might prefer the appearance of a graph created with -ciplot- or -rcap-. (-ciplot-, by Nick Cox, is available from SSC.)

[bruce fan]

Hey Clyde,

Thank you so much! It is a normal distribution, and I figured out how to plot the shaded graph with CI using graph twoway. Thanks!

[Matthew J. Baker]

Clyde Schechter and bruce fan --

I am one of the authors of the qregpd command, along with Travis A. Smith and David Powell. First, thanks Clyde Schechter for providing this code as a lot of people ask about it and now I can direct them here!

Pursuant to that, bruce fan - you mentioned that you were able to modify the code to include confidence intervals. Could you provide an example just so it is up for posterity's sake?

Best,

Matthew J. Baker

[bruce fan]

Hi Mat, @Matthew J. Baker Thanks for the qregpd command. Here are my codes to plot the graph with CI.

**********************
set seed 8675309
capture postutil clear
tempfile holding
postfile handle quantile coefficient se using `holding'

forvalues i=1(1)9{
local ii=`i'/10
qregpd gini carrier_mmc2 lnhhi, quantile(`ii') id(id) fix(quarter) optimize(mcmc) ///
draws(10000) burn(1000) arate(0.5)
post handle (`ii') (_b[carrier_mmc2]) (_se[carrier_mmc2])
}
postclose handle

use `holding', clear
gen lb=coefficient-se*1.65 //90% Confidence interval
gen ub=coefficient+se*1.65 //90% Confidence interval
graph twoway (rarea lb ub quantile) || (line coefficient quantile)
***********************

I also have one question for you, as I asked at the very beginning of this post, when I used your qregpd command as shown in the loop, did I include both id fixed effects and time fixed effects? Thank you!

[River Huang]

Dear Clyde, I followed your code

webuse grunfeld, clear

// Clyde Schechter
set seed 1234567
capture postutil clear
tempfile holding
postfile handle quantile coefficient lb ub using `holding'

forvalues i=1/9 {
    local ii = `i'/10
    qregpd invest mvalue kstock, quantile(`ii') id(company) fix(year) optimize(mcmc) ///
        draws(2000) burn(1000) arate(0.5)
    matrix M = r(table)
    post handle (`ii') (_b[mvalue]) (M[1, 5]) (M[1,6])
}
postclose handle

use `holding', clear
graph twoway connected coefficient lb ub quantile

but encountered the following error message

Adaptive MCMC optimization


Quantile Regression for Panel Data (QRPD)
     Number of obs:               200
     Number of groups:             10
     Min obs per group:            20
     Max obs per group:            20
------------------------------------------------------------------------------
      invest |      Coef.   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
      mvalue |   .0308365   .0016018    19.25   0.000      .027697    .0339761
      kstock |   .0857145   .0070731    12.12   0.000     .0718514    .0995775
------------------------------------------------------------------------------
No excluded instruments - standard QRPD estimation.


MCMC diagonstics:
     Mean acceptance rate:      0.381
     Total draws:                2000
     Burn-in draws:              1000
     Draws retained:             1000
     Value of objective function:   
             Mean:              -5.7679
             Min:               -8.1767
             Max:               -5.0766
MCMC notes:
     *Point estimates correspond to mean of draws.
     *Standard errors are derived from variance of draws.
Adaptive MCMC optimization


Quantile Regression for Panel Data (QRPD)
     Number of obs:               200
     Number of groups:             10
     Min obs per group:            20
     Max obs per group:            20
------------------------------------------------------------------------------
      invest |      Coef.   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
      mvalue |   .0465982   .0033938    13.73   0.000     .0399465      .05325
      kstock |   .0875357   .0055371    15.81   0.000     .0766832    .0983882
------------------------------------------------------------------------------
No excluded instruments - standard QRPD estimation.


MCMC diagonstics:
     Mean acceptance rate:      0.228
     Total draws:                2000
     Burn-in draws:              1000
     Draws retained:             1000
     Value of objective function:   
             Mean:              -1.4330
             Min:               -4.6860
             Max:               -0.1909
MCMC notes:
     *Point estimates correspond to mean of draws.
     *Standard errors are derived from variance of draws.
Adaptive MCMC optimization


Quantile Regression for Panel Data (QRPD)
     Number of obs:               200
     Number of groups:             10
     Min obs per group:            20
     Max obs per group:            20
------------------------------------------------------------------------------
      invest |      Coef.   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
      mvalue |   .0721535   .0029434    24.51   0.000     .0663846    .0779225
      kstock |   .0272715   .0113692     2.40   0.016     .0049882    .0495548
------------------------------------------------------------------------------
No excluded instruments - standard QRPD estimation.


MCMC diagonstics:
     Mean acceptance rate:      0.268
     Total draws:                2000
     Burn-in draws:              1000
     Draws retained:             1000
     Value of objective function:   
             Mean:              -4.2031
             Min:               -5.2876
             Max:               -2.3942
MCMC notes:
     *Point estimates correspond to mean of draws.
     *Standard errors are derived from variance of draws.
Adaptive MCMC optimization


Quantile Regression for Panel Data (QRPD)
     Number of obs:               200
     Number of groups:             10
     Min obs per group:            20
     Max obs per group:            20
------------------------------------------------------------------------------
      invest |      Coef.   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
      mvalue |   .0846515   .0019593    43.20   0.000     .0808114    .0884917
      kstock |   .1420837   .0132979    10.68   0.000     .1160202    .1681472
------------------------------------------------------------------------------
No excluded instruments - standard QRPD estimation.


MCMC diagonstics:
     Mean acceptance rate:      0.196
     Total draws:                2000
     Burn-in draws:              1000
     Draws retained:             1000
     Value of objective function:   
             Mean:              -1.5938
             Min:              -10.0185
             Max:               -0.2784
MCMC notes:
     *Point estimates correspond to mean of draws.
     *Standard errors are derived from variance of draws.
Adaptive MCMC optimization


Quantile Regression for Panel Data (QRPD)
     Number of obs:               200
     Number of groups:             10
     Min obs per group:            20
     Max obs per group:            20
------------------------------------------------------------------------------
      invest |      Coef.   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
      mvalue |   .0619978   .0278343     2.23   0.026     .0074436    .1165521
      kstock |   .3867696   .0972512     3.98   0.000     .1961608    .5773784
------------------------------------------------------------------------------
No excluded instruments - standard QRPD estimation.


MCMC diagonstics:
     Mean acceptance rate:      0.310
     Total draws:                2000
     Burn-in draws:              1000
     Draws retained:             1000
     Value of objective function:   
             Mean:              -0.9249
             Min:               -6.0812
             Max:               -0.1187
MCMC notes:
     *Point estimates correspond to mean of draws.
     *Standard errors are derived from variance of draws.
Adaptive MCMC optimization


Quantile Regression for Panel Data (QRPD)
     Number of obs:               200
     Number of groups:             10
     Min obs per group:            20
     Max obs per group:            20
------------------------------------------------------------------------------
      invest |      Coef.   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
      mvalue |   -2849956   3.46e+07    -0.08   0.934    -7.07e+07    6.50e+07
      kstock |  -1.15e+07   1.39e+08    -0.08   0.934    -2.85e+08    2.62e+08
------------------------------------------------------------------------------
No excluded instruments - standard QRPD estimation.


MCMC diagonstics:
     Mean acceptance rate:      0.956
     Total draws:                2000
     Burn-in draws:              1000
     Draws retained:             1000
     Value of objective function:   
             Mean:              -0.5726
             Min:               -1.4513
             Max:               -0.0000
MCMC notes:
     *Point estimates correspond to mean of draws.
     *Standard errors are derived from variance of draws.
Adaptive MCMC optimization


Quantile Regression for Panel Data (QRPD)
     Number of obs:               200
     Number of groups:             10
     Min obs per group:            20
     Max obs per group:            20
------------------------------------------------------------------------------
      invest |      Coef.   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
      mvalue |   5.40e+22   6.45e+23     0.08   0.933    -1.21e+24    1.32e+24
      kstock |   2.28e+23   2.72e+24     0.08   0.933    -5.10e+24    5.56e+24
------------------------------------------------------------------------------
No excluded instruments - standard QRPD estimation.


MCMC diagonstics:
     Mean acceptance rate:      0.994
     Total draws:                2000
     Burn-in draws:              1000
     Draws retained:             1000
     Value of objective function:   
             Mean:              -0.0722
             Min:               -0.1247
             Max:               -0.0000
MCMC notes:
     *Point estimates correspond to mean of draws.
     *Standard errors are derived from variance of draws.
Adaptive MCMC optimization


Quantile Regression for Panel Data (QRPD)
     Number of obs:               200
     Number of groups:             10
     Min obs per group:            20
     Max obs per group:            20
------------------------------------------------------------------------------
      invest |      Coef.   Std. Err.      z    P>|z|     [95% Conf. Interval]
-------------+----------------------------------------------------------------
      mvalue |   1.88e+12   3.22e+12     0.58   0.559    -4.43e+12    8.20e+12
      kstock |  -2.16e+12   3.71e+12    -0.58   0.560    -9.44e+12    5.11e+12
------------------------------------------------------------------------------
No excluded instruments - standard QRPD estimation.


MCMC diagonstics:
     Mean acceptance rate:      0.233
     Total draws:                2000
     Burn-in draws:              1000
     Draws retained:             1000
     Value of objective function:   
             Mean:              -0.4507
             Min:               -3.4836
             Max:               -0.0000
MCMC notes:
     *Point estimates correspond to mean of draws.
     *Standard errors are derived from variance of draws.
Adaptive MCMC optimization
estimates post: matrix has missing values
r(504);

Do you have any idea why this is happening? Thanks.

[Benjamin Sas]

Just had to do the same, and used different code which makes it easier to extract more coefficients of variables, and easier to understand

** VARIABLES WHERE TO PUT RESULTS
local ill_list malaria diarrhea dos toux peau ent oeil dentaire blessure diabetes mauxventre autre 

foreach var in `ill_list'{
    gen `var' =.
    gen se`var' = .
    gen ll`var' = .
    gen ul`var' =.
}



local quantiles "5 10 15 20 25 30 35 40 45 50 55 60 65 75 80 85 90 95"
local n: word count `quantiles'

gen q =.
forvalues i = 1/`n'{
    local qq: word `i' of `quantiles'
    replace q = `qq' in `i'
}


local i = 1

foreach q in `quantiles'{
    
    * Panel Unconditional Panel Reg
    qregpd quarter_sante     $control_vars ///
                            , q(`q') id(hhid) fix(passage)
                            
                mat aux = r(table) // Output of reg

    * Fill Estiamtes from uqreg
    local col = 9 // column where variables of interest start
        
    foreach var in `ill_list'{
        
        replace `var'     = aux[1,`col'] in `i'
        replace se`var' = aux[2,`col'] in `i'
        replace ll`var' = aux[5,`col'] in `i'    
        replace ul`var' = aux[6,`col'] in `i'
        
        local col = `col' + 1
        
    }

    
    local i = `i' + 1
}                        



* Graph
graph twoway (rarea lldiabetes uldiabetes q) || (line diabetes q) ///
    , xlabel(5(5)95)

[Dwi Utami]

Well, this is a little harder for me to answer as I do not know -qregpd-. Official Stata estimation commands leave behind a matrix r(table), and the confidence limits are found in the fifth and sixth rows of that matrix. So, if -qregpd- does that, you can modify the code as follows:

capture postutil clear
tempfile holding
postfile handle quantile coefficient lb ub using `holding'

forvalues i=1/9 {
local ii = `i'/10
qregpd gini carrier_mmc2 lnhhi, quantile(`ii') id(id) fix(quarter) optimize(mcmc) ///
draws(10000) burn(1000) arate(0.5)
matrix M = r(table)
post handle (`ii') (_b[carrier_mmc2]) (M[1, 5]) (M[1,6])
}
postclose handle

use `holding', clear
graph twoway connected coefficient lb ub quantile

(Changes related to grabbing confidence bounds in bold face.)

Now, if -qregpd- doesn't leave r(table) behind, then the alternative is to grab the standard error out of _se[], analogous to the way you grabbed the coefficient out of _b[] and then calculate the confidence bounds yourself from the coefficient and its standard error (I'm assuming here that the confidence bounds come from a normal or t-distribution in -qregpd-.)

Note also that you might prefer the appearance of a graph created with -ciplot- or -rcap-. (-ciplot-, by Nick Cox, is available from SSC.)

Dear @Clyde Schechter,

thank you for this code.

I would like to ask whether we can combine the coefficient from quantile regression from qregpd and the coefficient from OLS.

the visualization is like this:

I would really appreciate your help to provide the code.

Thank you very much.

[Clyde Schechter]

Try this:

capture postutil clear
tempfile holding
postfile handle quantile coefficient lb ub using `holding'

forvalues i=1/9 {
    local ii = `i'/10
    qregpd gini carrier_mmc2 lnhhi, quantile(`ii') id(id) fix(quarter) optimize(mcmc) ///
    draws(10000) burn(1000) arate(0.5)
    matrix M = r(table)
    post handle (`ii') (_b[carrier_mmc2]) (M[1, 5]) (M[1,6])
}
postclose handle

regress gini carrier_mmc2 lnhhi
gen ols_coeff = _b[carrier_mmc2]

label var lb "lower limits"
label var ub "upper limits"
label var coefficient "beta"
label var ols_coef "ols"

use `holding', clear
graph twoway connected coefficient lb ub ols_coeff quantile

[Dwi Utami]

Hi @Clyde Schechter,

thank you very much.

[Henry Nguyen]

Dear

Clyde Schechter ,
I used your code "tempfile holding
postfile handle quantile coefficient lb ub using `holding'

forvalues i=1/9 {
local ii = `i'/10
xtqreg blev efwamb1 lagmb lagA lagprof lagsize_sale, quantile(`ii') id(year) fix(id_firm) optimize(mcmc) ///
draws(10000) burn(1000) arate(0.5)
matrix M = r(table)
post handle (`ii') (_b[efwamb1]) (_c[lagmb]) (_d[lagprof]) (_e[lagsize_sales])(M[1, 5]) (M[1,6])
}
postclose handle

use `holding', clear
graph twoway connected coefficient lb ub quantile"
but data be deleted, and I get the results in the picture attach, and comment "variable lagmb not found
r(111)"
Can you tell me what's wrong, please?
Thank you very much,

[Clyde Schechter]

So, two things come to mind.

First, you are getting an error message saying variable lagmb not found. So this variable apparently does not exist in your data set. So investigate whether you have misnamed the variable in your data set, or if you have mistyped the variable's name in your code.

Second, I think you are running this code in pieces rather than all at once. I say that because you could not get any plot at all if you were running it all at once. Once the -xtqreg- command quit with an error because of not finding the variable lagmb, there would be no further execution. So it appears that you ran that block of code separately, and then tried to just ignore the error message about variable lagmb. You then tried to run the rest of the code. But that will just leave you with an empty data set because the loop hasn't actually run. Also, since tempfile holding is defined before the error about lagmb, when that error is reached, tempfile holding disappears. So you must be going through considerable contortions to even get the empty graph you are showing.

Fix up the problem about variable lagmb not being found. And run the code all at once: do not try to continue execution after encountering an error message--you will not get usable results that way. Only when the code runs without error from top to bottom can you expect to have a meaningful graph.

That said, your -postfile- and -post handle- commands are not compatible with each other, and you will get an error when you reach the -post handle- command. In -postfile- you specify four variables, but in -post handle- you specify 7 expressions. There must be an exact one-to-one correspondence between the variables in -postfile- and the expressions in -post handle-.

Finally, apart from the apparent non-existence of variable lagmb, I am unable to understand your -xtqreg- command. The syntax you are using is nothing like the syntax specified in the help file for the version of -xtqreg- on my system. Perhaps there is a different -xtqreg.ado- that you are using. If so, I can't say whether you have errors in that command or not.

[Henry Nguyen]

Thanks so much Clyde Schechter. I will try again.