I'm looking for the code for calculating marginal effects for the logit and count part separately after fitting a hurdle poisson model.
-
Login or Register
- Log in with
- You are not logged in. You can browse but not post. Login or Register by clicking 'Login or Register' at the top-right of this page. For more information on Statalist, see the FAQ.
-
-
Your post neither shows what commands you used nor specifies what marginal effects you want. Assuming you are referring to the command hplogit from the Stata Journal, you can specify the equation number in margins using "\(\text{eq}(\#n)\)" where \(n \in \mathbb{Z}^+\). Check out what marginal effects can be computed after Poisson and their formulas from the documentation.
A simple example from the documentation:Code:help poisson_postestimation
Options for predict
+------+
----+ Main +--------------------------------------------------------------------------------------------------------------------------------------------------------
n, the default, calculates the predicted number of events, which is exp(xb) if neither offset() nor exposure() was specified when the model was fit;
exp(xb + offset) if offset() was specified; or
exp(xb)*exposure if exposure() was specified.
Res.:Code:webuse ships, clear hplogit accident co_70_74 co_75_79 op_75_79, nolog *DEFAULT IS XB (LINEAR PREDICTION). NO OFFSET NOR EXPOSURE SPECIFIED. margins, expression(exp(predict(eq(#2)))) dydx(*)
For the logit equation, refer to https://www.statalist.org/forums/for...cient-estimate.Code:. hplogit accident co_70_74 co_75_79 op_75_79, nolog Poisson-Logit Hurdle Regression Number of obs = 34 Wald chi2(3) = 0.53 Log likelihood = -252.85156 Prob > chi2 = 0.9113 ------------------------------------------------------------------------------ | Coefficient Std. err. z P>|z| [95% conf. interval] -------------+---------------------------------------------------------------- logit | co_70_74 | -15.34799 785.9749 -0.02 0.984 -1555.83 1525.135 co_75_79 | -14.96603 1124.472 -0.01 0.989 -2218.89 2188.958 op_75_79 | -.6931546 .9486858 -0.73 0.465 -2.552544 1.166235 _cons | -.0000216 .6324555 -0.00 1.000 -1.239612 1.239568 -------------+---------------------------------------------------------------- poisson | co_70_74 | -.4402382 .1158772 -3.80 0.000 -.6673534 -.213123 co_75_79 | -1.037688 .1894416 -5.48 0.000 -1.408987 -.6663894 op_75_79 | .1547408 .1127558 1.37 0.170 -.0662566 .3757382 _cons | 2.827939 .0962007 29.40 0.000 2.639389 3.016489 ------------------------------------------------------------------------------ AIC Statistic = 15.109 . . *DEFAULT IS XB (LINEAR PREDICTION). NO OFFSET NOR EXPOSURE SPECIFIED. . . margins, expression(exp(predict(eq(#2)))) dydx(*) Average marginal effects Number of obs = 34 Model VCE: OIM Expression: exp(predict(eq(#2))) dy/dx wrt: co_70_74 co_75_79 op_75_79 ------------------------------------------------------------------------------ | Delta-method | dy/dx std. err. z P>|z| [95% conf. interval] -------------+---------------------------------------------------------------- co_70_74 | -6.470658 1.817137 -3.56 0.000 -10.03218 -2.909135 co_75_79 | -15.25203 2.994656 -5.09 0.000 -21.12144 -9.38261 op_75_79 | 2.274393 1.649322 1.38 0.168 -.9582189 5.507005 ------------------------------------------------------------------------------Last edited by Andrew Musau; 05 Jul 2024, 18:27. -
Thank you for your assistance Andrew. Much appreciated. Just to find out if the default reference variable for the logit part is 1? I have negative coefficients on my logit part of the model wher I know for sure that they should be positive.Last edited by Teresa Mungazi; 09 Jul 2024, 17:36.Comment
-
That appears to be the case, i.e., it estimates \(Pr(y=0)\). From
you haveCode:viewsource jhpoi_logit_ll.ado
You can also see that comparing it to logit.* LOGIT-POISSON HURDLE LOG-LIKELIHOOD: Joseph Hilbe: 30Sep2005
program jhpoi_logit_ll
version 9.1
args lnf beta1 beta2
tempvar pi mu
qui gen double `pi' = exp(`beta1')
qui gen double `mu' = exp(`beta2')
qui replace `lnf' = cond($ML_y1==0, ln(`pi'/(1+`pi')), ///
ln(1/(1+`pi')) -`mu' + $ML_y1 * `beta2' ///
- lngamma($ML_y1 + 1) - ln(1 - exp(-`mu')) )
end
Res.:Code:webuse ships, clear hplogit accident co_65_69 op_75_79, nolog logit accident co_65_69 op_75_79, nolog
Code:. hplogit accident co_65_69 op_75_79, nolog Poisson-Logit Hurdle Regression Number of obs = 34 Wald chi2(2) = 1.59 Log likelihood = -269.21136 Prob > chi2 = 0.4507 ------------------------------------------------------------------------------ | Coefficient Std. err. z P>|z| [95% conf. interval] -------------+---------------------------------------------------------------- logit | co_65_69 | .4266586 .8733348 0.49 0.625 -1.285046 2.138363 op_75_79 | -.9567469 .8379033 -1.14 0.254 -2.599007 .6855134 _cons | -.8417876 .6349363 -1.33 0.185 -2.08624 .4026647 -------------+---------------------------------------------------------------- poisson | co_65_69 | .4818266 .1097559 4.39 0.000 .2667091 .6969442 op_75_79 | .0018699 .1088363 0.02 0.986 -.2114452 .2151851 _cons | 2.461543 .0960013 25.64 0.000 2.273384 2.649702 ------------------------------------------------------------------------------ AIC Statistic = 16.012 . . logit accident co_65_69 op_75_79, nolog Logistic regression Number of obs = 34 LR chi2(2) = 1.67 Prob > chi2 = 0.4349 Log likelihood = -17.717553 Pseudo R2 = 0.0449 ------------------------------------------------------------------------------ accident | Coefficient Std. err. z P>|z| [95% conf. interval] -------------+---------------------------------------------------------------- co_65_69 | -.4266584 .8733347 -0.49 0.625 -2.138363 1.285046 op_75_79 | .9567466 .8379032 1.14 0.254 -.6855135 2.599007 _cons | .8417872 .6349362 1.33 0.185 -.4026649 2.086239 ------------------------------------------------------------------------------ .Last edited by Andrew Musau; 09 Jul 2024, 18:40.Comment
-
Thank you Andrew. Any way of changing it? I'd like 0 to be the reference category for the logit part, just for a clean intepretation of my coefficients. Can I manipulate the code so that I can reverse the default?Comment
-
I think it is easier to recode the outcome rather than mess with the program files.
Code:gen newoutcome=!outcome hplogit newoutcome ...
Last edited by Andrew Musau; 11 Jul 2024, 16:54.Comment
-
Sorry for #6, I was only thinking about the logit part of the model and forgetting about the Poisson part. Recoding the outcome won't work. Estimate with the default settings and then use margins to get the results that you want.
Res.:Code:webuse ships clear hplogit accident co_70_74 co_75_79 op_75_79, nolog margins, expression(predict(xb)*-1) dydx(*)
Code:. hplogit accident co_70_74 co_75_79 op_75_79, nolog Poisson-Logit Hurdle Regression Number of obs = 34 Wald chi2(3) = 0.53 Log likelihood = -252.85156 Prob > chi2 = 0.9113 ------------------------------------------------------------------------------ | Coefficient Std. err. z P>|z| [95% conf. interval] -------------+---------------------------------------------------------------- logit | co_70_74 | -15.34799 785.9749 -0.02 0.984 -1555.83 1525.135 co_75_79 | -14.96603 1124.472 -0.01 0.989 -2218.89 2188.958 op_75_79 | -.6931546 .9486858 -0.73 0.465 -2.552544 1.166235 _cons | -.0000216 .6324555 -0.00 1.000 -1.239612 1.239568 -------------+---------------------------------------------------------------- poisson | co_70_74 | -.4402382 .1158772 -3.80 0.000 -.6673534 -.213123 co_75_79 | -1.037688 .1894416 -5.48 0.000 -1.408987 -.6663894 op_75_79 | .1547408 .1127558 1.37 0.170 -.0662566 .3757382 _cons | 2.827939 .0962007 29.40 0.000 2.639389 3.016489 ------------------------------------------------------------------------------ AIC Statistic = 15.109 . . margins, expression(predict(xb)*-1) dydx(*) Average marginal effects Number of obs = 34 Model VCE: OIM Expression: predict(xb)*-1 dy/dx wrt: co_70_74 co_75_79 op_75_79 ------------------------------------------------------------------------------ | Delta-method | dy/dx std. err. z P>|z| [95% conf. interval] -------------+---------------------------------------------------------------- co_70_74 | 15.34799 785.9749 0.02 0.984 -1525.135 1555.83 co_75_79 | 14.96603 1124.472 0.01 0.989 -2188.958 2218.89 op_75_79 | .6931546 .9486858 0.73 0.465 -1.166235 2.552545 ------------------------------------------------------------------------------Comment
-
Thank you so much Andrew. Much appreciated!Comment
-
Oh, I'm only seeing the second part of your response now. So, is your expression for marginals for joint marginals or for the poisson part? How about the logit part? Why would recording not work?Comment
-
You cannot recode the outcome as it is 0/positive count and you need the actual counts for poisson. The coefficients from the poisson part are fine. Consider:
Res.:Code:webuse ships, clear hplogit accident co_65_69 op_75_79, nolog poisson accident co_65_69 op_75_79 if accident, nolog
If you want the base to be zero for the logit part, what I recommend from #7 is that you estimate with default settings then get the wanted estimates from margins. Here again is how you can get these estimates from margins. Alternatively, install cmp from SSC and estimate this model with that command. Below, I use estout from SSC to output the estimates,Code:. hplogit accident co_65_69 op_75_79, nolog Poisson-Logit Hurdle Regression Number of obs = 34 Wald chi2(2) = 1.59 Log likelihood = -269.21136 Prob > chi2 = 0.4507 ------------------------------------------------------------------------------ | Coefficient Std. err. z P>|z| [95% conf. interval] -------------+---------------------------------------------------------------- logit | co_65_69 | .4266586 .8733348 0.49 0.625 -1.285046 2.138363 op_75_79 | -.9567469 .8379033 -1.14 0.254 -2.599007 .6855134 _cons | -.8417876 .6349363 -1.33 0.185 -2.08624 .4026647 -------------+---------------------------------------------------------------- poisson | co_65_69 | .4818266 .1097559 4.39 0.000 .2667091 .6969442 op_75_79 | .0018699 .1088363 0.02 0.986 -.2114452 .2151851 _cons | 2.461543 .0960013 25.64 0.000 2.273384 2.649702 ------------------------------------------------------------------------------ AIC Statistic = 16.012 . . poisson accident co_65_69 op_75_79 if accident, nolog Poisson regression Number of obs = 26 LR chi2(2) = 18.41 Prob > chi2 = 0.0001 Log likelihood = -251.49396 Pseudo R2 = 0.0353 ------------------------------------------------------------------------------ accident | Coefficient Std. err. z P>|z| [95% conf. interval] -------------+---------------------------------------------------------------- co_65_69 | .4818187 .1097539 4.39 0.000 .266705 .6969325 op_75_79 | .0018698 .1088331 0.02 0.986 -.2114391 .2151788 _cons | 2.461551 .0959977 25.64 0.000 2.2734 2.649703 ------------------------------------------------------------------------------
Res.:Code:webuse ships clear hplogit accident co_70_74 co_75_79 op_75_79, nolog est sto hplogit *LOGIT margins, expression(-predict(eq(#1))) dydx(*) post est sto logit *POISSON est restore hplogit margins, expression(predict(eq(#2))) dydx(*) post est sto poisson *ORIGINAL ESTIMATES esttab hplogit, unstack mlab(none) *TRANSFORMED: LOGIT WITH BASE 0 esttab logit poisson, mlab(logit poisson)
Code:. esttab hplogit, unstack mlab(none) -------------------------------------------- (1) logit poisson -------------------------------------------- co_70_74 -15.35 -0.440*** (-0.02) (-3.80) co_75_79 -14.97 -1.038*** (-0.01) (-5.48) op_75_79 -0.693 0.155 (-0.73) (1.37) _cons -0.0000216 2.828*** (-0.00) (29.40) -------------------------------------------- N 34 -------------------------------------------- t statistics in parentheses * p<0.05, ** p<0.01, *** p<0.001 . . *TRANSFORMED: LOGIT WITH BASE 0 . . esttab logit poisson, mlab(logit poisson) -------------------------------------------- (1) (2) logit poisson -------------------------------------------- co_70_74 15.35 -0.440*** (0.02) (-3.80) co_75_79 14.97 -1.038*** (0.01) (-5.48) op_75_79 0.693 0.155 (0.73) (1.37) -------------------------------------------- N 34 34 -------------------------------------------- t statistics in parentheses * p<0.05, ** p<0.01, *** p<0.001 .Comment
-
Thank you. I'll give this a try.Comment
Comment