I have a log likelihood function of a model and I want to find $\mu$ and $\sigma^2$ which maximize the log likelihood. Since the log lik function is quite complex, I decided to use Nelder-Mead algorithm from scipy module. What I would like to do now is to estimate the errors on $\mu$ and $\sigma^2$ retrieved by Nelder-Mead algorithm.

How to do it without calculating the Hessian of the function?

EDIT: I found a way to approximate the Hessian (using numdifftools.core.Hessian) but I'm not sure how much can I rely on it.

EDIT 2: In some cases I manage to use the approximated Hessian (from numdifftools) but in other cases, it doesn't work. And it seems like there is a number precision problem. So I'm still looking for a solution with which I can efficiently asses the approximation error without using Hessian or an implementation of a Hessian approximation which works good in high number precision.

  • $\begingroup$ Does the data follow some specific distribution? $\endgroup$ – HelloGoodbye Aug 22 '18 at 10:04
  • $\begingroup$ Yes, it follows a Gaussian distribution. $\endgroup$ – Đorđe Relić Aug 22 '18 at 10:32
  • $\begingroup$ Have you considered bootstrapping? $\endgroup$ – Denziloe Aug 22 '18 at 13:16
  • $\begingroup$ @Denziloe I haven't. Mostly because of the reason it would not give me one error value for each parameter and also it would depend on the number of samples I choose, right? The thing is that I repeat the process of inferring mean and variance thousands of times (for many subsets of my given data set) and I would like to have a more convenient way of getting the error rather than sampling. $\endgroup$ – Đorđe Relić Aug 22 '18 at 14:04
  • $\begingroup$ Sure, if optimisation is important then it makes sense to rule out bootstrapping. Choosing the number of samples isn't an issue though, you just need to choose "enough". The more the better but after a moderate number the answer should start to converge. And you'd get errors for each parameter. $\endgroup$ – Denziloe Aug 22 '18 at 14:07

You could assume that the real mean and variance are μ* and σ2* respectively, and use the formulas you use to estimate the estimators μ and σ2 (which are themselves stochastic variables, since they depend on your sampled values which are stochastic values) to calculate the variances of μ and σ2 as functions of the variance in the sampled values, i.e. σ2*.

Since you don't know σ2* , you could assume that σ2 ≈ σ2* and plug in the value for σ2 instead of σ2* when evaluating the expressions for the variances of μ and σ2. The square root of those variances, i.e. the standard deviations of the estimators, give you an indication of the error of the estimators.


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.