Announcement

Collapse
No announcement yet.
X
  • Filter
  • Time
  • Show
Clear All
new posts

  • Creating coefplots using lincom output

    Dear Stata expert,
    I am estimating the below model where the coefficients of interest are "POSTGRAD" + "POSTGRAD#SUBJECT" which aims to capture the wage advantage to a postgraduate of studying a particular subject. The omitted subject dummy is subject 11(=Education). To capture the overall effect, I am using the "lincom" command to add the relevant coefficients together. My goal is to display these combined point estimates and confidence intervals in a coefplot. This a categorical by categorical interaction and I am using the most up to date version of Stata.

    Code:
    reg LOGHOURPAY i.POSTGRAD##ib11.SUBJECT AGE AGESQ MALE i.REGION FULLTIME i.TENURE i.BIRTHCOHORT i.MARRIED ETHNICITY BADHEALTH i.year, robust
    lincom 1.POSTGRAD + 1.POSTGRAD#1.SUBJECT
    lincom 1.POSTGRAD + 1.POSTGRAD#2.SUBJECT
    lincom 1.POSTGRAD + 1.POSTGRAD#2.SUBJECT
    lincom 1.POSTGRAD + 1.POSTGRAD#3.SUBJECT
    lincom 1.POSTGRAD + 1.POSTGRAD#4.SUBJECT
    lincom 1.POSTGRAD + 1.POSTGRAD#5.SUBJECT
    lincom 1.POSTGRAD + 1.POSTGRAD#6.SUBJECT
    lincom 1.POSTGRAD + 1.POSTGRAD#7.SUBJECT
    lincom 1.POSTGRAD + 1.POSTGRAD#8.SUBJECT
    lincom 1.POSTGRAD + 1.POSTGRAD#9.SUBJECT
    lincom 1.POSTGRAD + 1.POSTGRAD#10.SUBJECT
    lincom 1.POSTGRAD + 0

    Any help would be greatly appreicated

  • #2
    Use margins which allows you to post these estimates to e(b) and e(V)

    Code:
    qui reg LOGHOURPAY i.POSTGRAD##ib11.SUBJECT AGE AGESQ MALE i.REGION FULLTIME i.TENURE i.BIRTHCOHORT i.MARRIED ETHNICITY BADHEALTH i.year, robust
    estimates store full
    forval i=1/10{
        estimates restore full
        margins, expression(_b[1.POSTGRAD]+ _b[1.POSTGRAD#`i'.SUBJECT]) post
        estimates store sub`i'
    }
    Last edited by Andrew Musau; 20 Apr 2021, 08:42.

    Comment


    • #3
      Hi Andrew,
      Thank you for your quick reply. What code would I need to convert these saved estimates into a coefplot with confidence intervals?

      Thank you

      Comment


      • #4
        You can install coefplot from SSC and graph directly.

        Code:
        ssc install coefplot, replace
        local coefs
        forval i=1/9{
            local coefs "`coefs'  sub`i' \"
        }
        local coefs `coefs' sub10
        di "`coefs'"    
        coefplot (`coefs'),  aseq swapnames
        So in the code in #2, we named these estimates sub1, sub2, ,,,, sub10. I use a loop to generate these, but in the command, these will be entered as

        Code:
        coefplot (sub1 \  sub2 \  sub3 \  sub4 \  sub5 \  sub6 \  sub7 \  sub8 \  sub9 \ sub10),  aseq swapnames
        You can relabel the coefficients using the option -ylab()-, e.g.,

        Code:
        coefplot (`coefs'),  aseq swapnames ylab(1 "my lab 1" 2 "my lab 2" 3 "my lab 3")
        Last edited by Andrew Musau; 20 Apr 2021, 10:10.

        Comment


        • #5
          Dear Andrew,
          Thank you so much for your help - the code has allowed me to create the graph I was aiming for. I just had one more question: Would I be able to alter the code to include a term for subject 11 so that it is represented in the graph? The coefficient for this subject is just b[1.POSTGRAD] seeing as it is the omitted dummy variable for subject.

          Thank you

          Comment


          • #6
            So you should have:

            Code:
            estimates restore full
            margins, expression(_b[1.POSTGRAD]+ _b[1.POSTGRAD#11.SUBJECT]) post
            est sto sub11
            Then

            Code:
            coefplot (`coefs' \ sub11), aseq swapnames omitted

            Comment


            • #7
              Dear Andrew,
              Thank you for providing this code, I have been able to create the graph I was aiming for. I have since thought of another question, and I was wondering if you could assit me. I am conducting the above regressions separately for men and women but have since realised it may be useful to display the results for each subject (for both gender) on one graph. Is there a way I could combine the results of two regressions onto one graph?
              Regards,
              Ella Garrett

              Comment


              • #8
                For coefplot, it does not matter where the estimates come from. So estimates 1-10 can be from one regression and 11-20 from another. If you want to differentiate the coefficients, store the coefficients separately and label them the same across regressions.

                Code:
                sysuse auto, clear
                eststo m1: regress mpg weight turn trunk displacement if foreign
                eststo m2: regress mpg weight turn trunk displacement if !foreign
                coefplot m1 m2,  drop(_cons) leg(order(2 "Foreign" 4 "Domestic")) scheme(s1color)
                Click image for larger version

Name:	Graph.png
Views:	1
Size:	14.4 KB
ID:	1605956

                Comment


                • #9
                  Dear Andrew,

                  The code I am currently using is as follows:

                  eststo m1: qui reg LOGHOURPAY i.MAS##ib4.BROAD AGE AGESQ i.REGION FULLTIME i.TENURE i.BIRTHCOHORT i.MARRIED ETHNICITY BADHEALTH i.year if SEX==1, robust
                  estimates store full
                  forval i=1/3{
                  estimates restore full
                  margins, expression(_b[1.MAS]+ _b[1.MAS#`i'.BROAD]) post
                  estimates store sub`i'
                  }
                  estimates restore full
                  margins, expression(_b[1.MAS]+ _b[1.MAS#4.BROAD]) post
                  est sto sub4
                  ssc install coefplot, replace
                  local coefs
                  forval i=1/2{
                  local coefs "`coefs' sub`i' "
                  }
                  local coefs `coefs' sub3
                  di "`coefs'"
                  eststo m2: qui reg LOGHOURPAY i.PHD##ib4.BROAD AGE AGESQ i.REGION FULLTIME i.TENURE i.BIRTHCOHORT i.MARRIED ETHNICITY BADHEALTH i.year if SEX==1, robust
                  estimates store full
                  forval i=1/3{
                  estimates restore full
                  margins, expression(_b[1.PHD]+ _b[1.PHD#`i'.BROAD]) post
                  estimates store sub`i'
                  }
                  estimates restore full
                  margins, expression(_b[1.PHD]+ _b[1.PHD#4.BROAD]) post
                  est sto sub4
                  ssc install coefplot, replace
                  local coefs
                  forval i=1/2{
                  local coefs "`coefs' sub`i' "
                  }
                  local coefs `coefs' sub3
                  di "`coefs'"
                  coefplot m1(`coefs' \ sub4) m2(`coefs'\ sub4), vertical recast(bar) ciopts(recast(rcap))

                  But Stata produces an error message saying "invalid something: unmatched open parenthesis or bracket". Here BROAD refers to a broader subject grouping and MAS and PHD refer to different levels of postgraduate study. I was aiming to produce a coefplot displaying the BROAD coefficients for both MAS and PHD. Would you be able to point out where I have made a mistake in my code?
                  Regards,
                  Ella

                  Comment


                  • #10
                    You need to write a program to rename the margins coefficients. The following example uses erepost from SSC.

                    Code:
                    ssc install erepost, replace

                    Name estimates uniquely, e.g., coef1, coef2, ..., but name the coefficients the same across models. By default, margins names all coefficients "_cons", but you can replace this with your own names.

                    Code:
                    cap program drop repostb
                    program repostb, eclass
                    erepost b = b, rename
                    end
                    
                    sysuse auto
                    *FIRST REGRESSION
                    eststo m1: regress mpg price weight i.rep78
                    local j 1
                    forval i=2/4{
                    est restore m1
                    margins, expression(_b[price]+_b[`i'.rep78]) post
                    mat b=e(b)
                    mat colname b = "coef`i'"
                    repostb
                    est sto coef`j'
                    local ++j
                    }
                    local j `=`j'-1'
                    *SECOND REGRESSION
                    eststo m2: regress turn price weight i.rep78
                    forval i=2/4{
                    est restore m2
                    margins, expression(_b[price]+_b[`i'.rep78]) post
                    mat b=e(b)
                    mat colname b = "coef`i'"
                    local ++j
                    repostb
                    est sto coef`j'
                    }
                    coefplot (coef1\ coef2\coef3) (coef4\ coef5\coef6), aseq swapnames ylab("") legend(order(2 "mpg" 4 "turn"))
                    Res.:
                    Click image for larger version

Name:	Graph.png
Views:	1
Size:	14.3 KB
ID:	1606164

                    Last edited by Andrew Musau; 27 Apr 2021, 05:18.

                    Comment


                    • #11
                      Dear Andrew,

                      Applying the code provided above to my dataset gives me:
                      ssc install erepost, replace

                      cap program drop repostb
                      program repostb, eclass
                      erepost b = b, rename
                      end

                      eststo m1: qui reg LOGHOURPAY i.MAS##ib4.BROAD AGE AGESQ i.REGION FULLTIME i.TENURE i.BIRTHCOHORT i.MARRIED ETHNICITY BADHEALTH i.year if SEX==1, robust
                      local j 1
                      forval i=1/3{
                      est restore m1
                      margins, expression(_b[1.MAS]+_b[1.MAS#`i'.BROAD]) post
                      mat b=e(b)
                      mat colname b = "coef`i'"
                      repostb
                      est sto coef`j'
                      local ++j
                      }
                      local j `=`j'-1'

                      eststo m2: qui reg LOGHOURPAY i.PHD##ib4.BROAD AGE AGESQ i.REGION FULLTIME i.TENURE i.BIRTHCOHORT i.MARRIED ETHNICITY BADHEALTH i.year if SEX==1, robust
                      forval i=1/3{
                      est restore m2
                      margins, expression(_b[1.PHD]+_b[1.PHD#`i'.BROAD]) post
                      mat b=e(b)
                      mat colname b = "coef`i'"
                      local ++j
                      repostb
                      est sto coef`j'
                      }

                      coefplot (coef1\coef2\coef3) (coef4\coef5\coef6), aseq swapnames omitted xlab(1 "Medicine & Related" 2 "STEM" 3 "Business & Finance") vertical recast(bar) legend(order(2 "Masters" 4 "PhD"))

                      which produces the below graph:

                      Would you be able to point out where I have made a mistake in my code so that the labels on the x axis match up? I also would like the graph to display the omitted broad subject category:
                      (_b[1.MAS]+ _b[1.MAS#4.BROAD]) and
                      (_b[1.PHD]+ _b[1.PHD#11.BROAD]) for both masters and PhD level. Thank you for your continued support

                      Comment


                      • #12
                        I cannot see any graph in your post #11. Upload it in .PNG format.

                        Comment


                        • #13
                          Click image for larger version

Name:	Screenshot 2021-04-28 at 10.57.13.png
Views:	1
Size:	208.6 KB
ID:	1606428
                          Dear Stata expert, I am estimating the below model where the coefficients of interest are "POSTGRAD" + "POSTGRAD#SUBJECT" which aims to

                          Comment


                          • #14
                            Would you be able to point out where I have made a mistake in my code so that the labels on the x axis match up?
                            The axis is continuous, so it looks to me like you need:

                            Code:
                            xlab(1.5 "Medicine & Related" 4.5 "STEM" 7.5 "Business & Finance")
                            You may adjust these values so that they align with the middle of the bars. However, I think what you want is:


                            Code:
                            xlab("") eqlabels("Medicine & Related"  "STEM" "Business & Finance")


                            I also would like the graph to display the omitted broad subject category:
                            See #6 on saving the coefficients and using the options -omitted- and possibly -baselevels-. Note that the omitted coefficients are equal to 0 and have no confidence intervals. So with your design, nothing will be displayed apart from the category label. Here is an example that shows this.


                            Code:
                            sysuse auto, clear
                            regress mpg i.foreign
                            coefplot, omitted drop(_cons) baselevels recast(bar) bcolor(gray%40) vertical scheme(s1mono)

                            Res.:
                            Click image for larger version

Name:	Graph.png
Views:	1
Size:	24.3 KB
ID:	1606435

                            Last edited by Andrew Musau; 28 Apr 2021, 06:04.

                            Comment

                            Working...
                            X