# After running a probit, how can I generate the margins for the whole distribution?

I'm using Stata.

I ran a probit of the form $$\text{outcome}_i = \beta \ f(\text{income}_i) + \gamma\text{ Controls}_i$$

Where $$f(\text{income}_i)$$ is a fractional polynomial. I'm interested on the marginal probability $$\frac{\delta \ P \ (\text{outcome}_i | \text{mean (Controls}_i))}{\delta \ \text{income}_i}$$ for all values of $$i$$. When I run  margins, dydx(varname) I only get the average marginal effect, not for the whole distribution. Moreover, I have to calcualte bootstrap standard errors, so I need some way to get the marginal probability that doesn't rely on specifying the functional form beforehand (i.e. _b[var1] + 2*_b[var2]*var2  for a 2nd order polynomial), as $$f(\cdot)$$ will change for each of the 500 bootstrap estimations.

A bit of context in case this is an instance of an xy problem.

I'm implementing an estimator of the form:

$$\beta(\mathbf{x}) = \triangledown \Psi(\mathbf{x}) + \frac{\Psi(\mathbf{x})}{\mathrm{G_M}(\mathbf{x})}\triangledown \mathrm{G_M}(\mathbf{x})$$

Where $$\mathbf{x}$$ is a vector of size $$N$$ (number of observations), $$\Psi(\mathbf{x})$$ is the prediction of the conditional mean, $$\mathrm{G_M}(\mathbf{x})$$ is the prediction of the probability that the outcome is greater than zero, and $$\triangledown \Psi(\mathbf{x})$$, $$\triangledown \mathrm{G_M}(\mathbf{x})$$ refer to their respective marginal effects.

The paper I'm following specifies:  I hope I have been clear enough. Thanks in advance.

• margins has a badly documented gen(me_income) option. Does that do the trick for storing the marginal effects (rather than just their average)? – Dimitriy V. Masterov Oct 18 at 1:12
• That does seem to do the trick regarding storing the marginal effects! I read the documentation for margins all over and couldn't find that! Thanks! – Pablo Derbez Oct 18 at 1:28
• Try help undocumented to see the other goodies. Can I ask how you specify the derivative with respect to the fp-ized variable? – Dimitriy V. Masterov Oct 18 at 1:59
• Thanks! I'm not sure I understood your question but I think I've found a solution to the derivatives issue: After running margins, xb(income_1 income_2) gen(me) it generates me1 and me2. These are the marginal effects wrt the variable transformed according to fp (i.e. raised to some power) so by the chain rule the marginal effect with respect to the non transformed variable should be me1*dydx where dydx is the derivative of the transformed variable wrt to the original. I wrote some code to do that, if you're interested I can post it. – Pablo Derbez Oct 18 at 3:19
• Also you can post your answer and I will mark it as solved, that turned out to be the main issue. – Pablo Derbez Oct 18 at 3:22

First of all, thanks to Dimitriy V. Masterov for pointing me towards the fact that I can save the margins for the whole distribution using margins, dydx(var1 var2), gen(newvar), which will generate newvar1 and newvar2 for var1 and var2 respectively.

Concerning the automation of marginal probability from the fp command, I wrote some code that can handle it, although only for two fp terms.

A fractional polynomial $$x^{f(p)}+x^{f(q)}$$ of $$x$$ is simply $$x^p + x^q$$, with the excepton that $$x^{f(0)}=\ln(x)$$. Also, if any term is repeated, say $$x^{f(p)}+x^{f(p)}$$ then the second term is $$ln(x^{f(p)})$$.

Our outcome is $$y= x^{f(p)}+x^{f(q)}$$. Then, $$\frac{\delta y}{\delta x} = \frac{\delta y}{\delta x^{f(p)}}(\frac{\delta x^{f(p)}}{\delta x})+ \frac{\delta y}{\delta x^{f(q)}}(\frac{\delta x^{f(q)}}{\delta x})$$

Where $$\frac{\delta y}{\delta x^{f(p)}}$$ and $$\frac{\delta y}{\delta x^{f(q)}}$$ are the results from margins, dydx(var1 var2) on the transformed variables,

and $$\frac{\delta x^{f(p)}}{\delta x}$$, $$\frac{\delta x^{f(q)}}{\delta x}$$ are the derivatives of the transformed variables with respect to $$x$$. Thus the code is as follows:

fp <ingpre_noL>, scale(1 1) replace: probit anytransfer <ingpre_noL> \$controls i.id_entidad [pw=factor], robust cluster(upm) iterate(30)

margins, dydx(ingpre_noL_1 ingpre_noL_2) gen(me)

tempvar beta1p beta2p

scalar r = e(fp_fp)[1,1]

if r==0 {
gen beta1p' = 1/(ingpre_noL)
}
else if r!=0 {
gen beta1p' = r*(ingpre_noL)^(r-1)
}

scalar t = e(fp_fp)[1,2]

if r!=t & t==0{
gen beta2p' = 1/(ingpre_noL+1)
}
else if r!=t & t!=0 {
gen beta2p' = t*(ingpre_noL+1)^(t-1)
}
else if r==t & t==0 {
gen beta2p' = 1/(ingpre_noL+1)*ln(ingpre_noL)
}
else if r==t & t!=0{
gen beta2p' = t/(ingpre_noL+1)
}

gen beta_p = me1*beta1p' + me2*beta2p'