Notice: On March 31, it was **announced** that Statalist is moving from an email list to a **forum**. The old list will shut down on April 23, and its replacement, **statalist.org** is already up and running.

[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

From |
"Mak, Timothy" <timothy.mak07@imperial.ac.uk> |

To |
"statalist@hsphsun2.harvard.edu" <statalist@hsphsun2.harvard.edu> |

Subject |
st: What's wrong with my differentiation? |

Date |
Tue, 26 Oct 2010 13:08:36 +0100 |

Dear list, Sorry if this is not much of a Stata question, but it has been bugging me for hours. Let's say we input the following: input yesno freq 1 4 0 6 end logit yesno [fw=freq] matrix list e(V) We get a variance of 0.416667, suggesting the hessian = -1/0.416667 = -2.4 However, for some reason, when I tried to get this hessian analytically, I fail. I know it must be to do with my calculations, but I haven't been able to spot the mistake for 5 hours, and I checked it with Wolframalpha, and it's still wrong. So here's my derivation: y=4 z=6 loglikelihood = lik = y*log(p) + z*log(1-p) + constant p = 1/(1+exp(-x)) = exp(x)/(1+exp(x)) dlik/dx = dlik/dp * dp/dx (chain rule) dlik/dp = y/p - z/(1-p) dp/dx = exp(x)/(1+exp(x))^2 = exp(-x)/(1+exp(-x))^2 Therefore, dlik/dx = (y/p - z/(1-p)) * exp(x)/(1+exp(x))^2 Denote: A=exp(x)/(1+exp(x))^2 = dp/dx B=(y/p - z/(1-p)) ******************** dlik/dx = A*B ******************** d2lik/dx2 = d/dx (dlik/dx) = A dB/dx + B dA/dx (product rule) dB/dx = dB/dp * dp/dx = dB/dp * A (chain rule) dB/dp = -yp^(-2) + z(1-p)^(-2) = z/(1-p)^2 - y/p^2 ******************************** dB/dx = A* (z/(1-p)^2 - y/p^2 ) ******************************** dA/dx = exp(x) d/dx (1+exp(x))^(-2) + (1+exp(x))^(-2) * d/dx (exp(x)) (product rule) d/dx(1+exp(x))^(-2) = -2exp(x)(1+exp(x))^(-3) (chain rule) d/dx(exp(x)) = exp(x) Therefore, ***************************************************************************** dA/dx = -2 * exp(x) * exp(x) * (1+exp(x))^(-3) + (1+exp(x))^(-2) * exp(x) = A * (-2exp(x)/(1+exp(x)) + 1) = A * (1 - 2/(1+exp(-x))) ***************************************************************************** Therefore, ******************************************************************** d2lik/dx2 = A*A*( z/(1-p)^2 - y/p^2 ) + B * A * (1 - 2/(1+exp(-x))) ******************************************************************** Now, for the above toy example, the ML of p = 0.4, which means on the logit scale, x=logit(0.4)=-0.405 Plugging x=-0.405 into the above equations: A=0.24 B=0 d2lik/dx2 = 0.24^2 * (6/0.36 - 4/0.16) = -0.48 Hessian = -0.48 != -2.4, which is what the answer should be! Thank you very much for your help if you've read so far! As I mentioned before I tried putting the above into Wolframalpha, but it ends up with such a complicated formula that after 1 or 2 hours, I still haven't been able to work out what's gone wrong. Following Wolframalpha's formula, however, I did manage to get -2.4, so obviously my formula above contains a mistake. (The reason I need the analytic formula is so that I can do write an estimation procedure for something, for which the above represents a special case. What I'm writing is a small extension to -logit-, although I can't exactly use -logit- in my engine, hence my wish to write down the analytic formula.) Thanks again, Tim * * For searches and help try: * http://www.stata.com/help.cgi?search * http://www.stata.com/support/statalist/faq * http://www.ats.ucla.edu/stat/stata/

**References**:**st: problem using catplot to graph vertical bars***From:*"Leung, Adelaine K.W." <Adelaine_Leung@hms.harvard.edu>

**Re: st: problem using catplot to graph vertical bars***From:*Ulrich Kohler <kohler@wzb.eu>

- Prev by Date:
**st: RE: RE: RE: nl @ within ado** - Next by Date:
**st: RE: problem using catplot to graph vertical bars** - Previous by thread:
**Re: st: problem using catplot to graph vertical bars** - Next by thread:
**st: RE: problem using catplot to graph vertical bars** - Index(es):