[Numpy-discussion] Calculation of a hessian

Jose Gomez-Dans jgomezdans at gmail.com
Fri Aug 8 05:51:41 EDT 2014


Your function looks fairly simple to differentiate by hand, but if you have
access to the gradient (or you estimate it numerically using scipy...),
this function might do the job:

def hessian ( x, the_func, epsilon=1e-8):
    """Numerical approximation to the Hessian
    Parameters
    ------------
    x: array-like
        The evaluation point
    the_func: function
        The function. We assume that the function returns the function
value and
        the associated gradient as the second return element
    epsilon: float
        The size of the step
    """

    N = x.size
    h = np.zeros((N,N))
    df_0 = the_func ( x )[1]
    for i in xrange(N):
        xx0 = 1.*x[i]
        x[i] = xx0 + epsilon
        df_1 = the_func ( x )[1]
        h[i,:] = (df_1 - df_0)/epsilon
        x[i] = xx0
    return h

Jose


On 8 August 2014 08:31, Kiko <kikocorreoso at gmail.com> wrote:

> Hi all,
>
> I am trying to calculate a Hessian. I am using numdifftools for this (
> https://pypi.python.org/pypi/Numdifftools).
>
> My question is, is it possible to make it using pure numpy?.
>
> The actual code is like this:
>
>
> *import numdifftools as nd*
> *import numpy as np*
>
> *def log_likelihood(params):*
> *    sum1 = 0; sum2 = 0*
> *    mu = params[0]; sigma = params[1]; xi = params[2]*
> *    for z in data:*
> *        x = 1 + xi * ((z-mu)/sigma)*
> *        sum1 += np.log(x)*
> *        sum2 += x**(-1.0/xi)*
> *    return -((-len(data) * np.log(sigma)) - (1 + 1/xi)*sum1 - sum2) #
> negated so we can use 'minimum'*
>
> *kk = nd.Hessian(log_likelihood)*
>
> Thanks in advance.
>
> _______________________________________________
> NumPy-Discussion mailing list
> NumPy-Discussion at scipy.org
> http://mail.scipy.org/mailman/listinfo/numpy-discussion
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mail.python.org/pipermail/numpy-discussion/attachments/20140808/0261d55a/attachment.html>


More information about the NumPy-Discussion mailing list