segunda derivada numpy de uma matriz ndimensional

Eu tenho um conjunto de dados de simulação em que gostaria de encontrar a menor inclinação em n dimensões. O espaçamento dos dados é constante ao longo de cada dimensão, mas não é o mesmo (eu poderia mudar isso por uma questão de simplicidade).

Eu posso viver com alguma imprecisão numérica, especialmente nas bordas. Eu preferiria não gerar um spline e usar esse derivado; apenas os valores brutos seriam suficientes.

É possível calcular a primeira derivada comnumpy usando onumpy.gradient() função.

import numpy as np

data = np.random.rand(30,50,40,20)
first_derivative = np.gradient(data)
# second_derivative = ??? <--- there be kudos (:

Este é um comentário sobre laplace versus a matriz hessiana; isso não é mais uma pergunta, mas serve para ajudar a entender os futuros leitores.

Uso como caso de teste uma função 2D para determinar a área 'mais plana' abaixo de um limite. As figuras a seguir mostram a diferença nos resultados entre usar o mínimo desecond_derivative_abs = np.abs(laplace(data)) e o mínimo de:

second_derivative_abs = np.zeros(data.shape)
hess = hessian(data)
# based on the function description; would [-1] be more appropriate? 
for i in hess[0]: # calculate a norm
    for j in i[0]:
        second_derivative_abs += j*j

A escala de cores representa os valores das funções, as setas representam a primeira derivada (gradiente), o ponto vermelho o ponto mais próximo de zero e a linha vermelha o limite.

A função geradora dos dados foi( 1-np.exp(-10*xi**2 - yi**2) )/100.0 com xi, yi sendo gerado comnp.meshgrid.

Laplace:

Hessian:

questionAnswers(3)

yourAnswerToTheQuestion