对称函数的梯度在所有维度上都应具有相同的导数。
numpy.gradient提供了不同的组件。
这是MWE。
import numpy as np
x = (-1,0,1)
y = (-1,0,1)
X,Y = np.meshgrid(x,y)
f = 1/(X*X + Y*Y +1.0)
print(f)
>> [[0.33333333 0.5 0.33333333]
[0.5 1. 0.5 ]
[0.33333333 0.5 0.33333333]]
这两个维度的值相同。
但是
np.gradient(f)
给[array([[ 0.16666667, 0.5 , 0.16666667],
[ 0. , 0. , 0. ],
[-0.16666667, -0.5 , -0.16666667]]),
array([[ 0.16666667, 0. , -0.16666667],
[ 0.5 , 0. , -0.5 ],
[ 0.16666667, 0. , -0.16666667]])]
梯度的两个分量都不同。
为什么这样?
我在解释输出时缺少什么?
最佳答案
让我们逐步完成此步骤。所以首先,正如喵喵狗正确提到的
numpy计算一个方向上的导数。
Numpy计算梯度的方式
重要的是要注意np.gradient
使用中心差的含义(为简单起见,我们只看一个方向):
grad_f[i] = (f[i+1] - f[i])/2 + (f[i] - f[i-1])/2 = (f[i+1] - f[i-1])/2
在边界处numpy计算(以
min
为例)grad_f[min] = f[min+1] - f[min]
grad_f[max] = f[max] - f[max-1]
在您的情况下,边界为
0
和2
。2D盒
如果使用多个维,则需要考虑导数的方向。
np.gradient
计算所有可能方向上的导数。让我们重现您的结果:让我们沿着列移动,以便使用行向量进行计算
f[1,:] - f[0,:]
输出量
array([0.16666667, 0.5 , 0.16666667])
恰好是渐变的第一个元素的第一行。
该行是使用居中导数计算的,因此:
(f[2,:]-f[1,:])/2 + (f[1,:]-f[0,:])/2
输出量
array([0., 0., 0.])
第三行:
f[2,:] - f[1,:]
输出量
array([-0.16666667, -0.5 , -0.16666667])
对于另一个方向,只需交换
:
和数字,并记住您现在正在计算列向量。在对称函数的情况下(如您的情况),这直接导致转置导数。3D外壳
x_ = (-1,0,4)
y_ = (-3,0,1)
z_ = (-1,0,12)
x, y, z = np.meshgrid(x_, y_, z_, indexing='ij')
f = 1/(x**2 + y**2 + z**2 + 1)
np.gradient(f)[1]
输出量
array([[[ *2.50000000e-01, 4.09090909e-01, 3.97702165e-04*],
[ 8.33333333e-02, 1.21212121e-01, 1.75554093e-04],
[-8.33333333e-02, -1.66666667e-01, -4.65939801e-05]],
[[ **4.09090909e-01, 9.00000000e-01, 4.03045231e-04**],
[ 1.21212121e-01, 2.00000000e-01, 1.77904287e-04],
[-1.66666667e-01, -5.00000000e-01, -4.72366556e-05]],
[[ ***1.85185185e-02, 2.03619910e-02, 3.28827183e-04***],
[ 7.79727096e-03, 8.54700855e-03, 1.45243282e-04],
[-2.92397661e-03, -3.26797386e-03, -3.83406181e-05]]])
此处给出的梯度是沿着行计算的(
0
沿着矩阵,1
沿着行,2
沿着列)。这可以通过
(f[:,1,:] - f[:,0,:])
输出量
array([[*2.50000000e-01, 4.09090909e-01, 3.97702165e-04*],
[**4.09090909e-01, 9.00000000e-01, 4.03045231e-04**],
[***1.85185185e-02, 2.03619910e-02, 3.28827183e-04***]])
我添加了星号,以便可以清楚地在哪里找到相应的行向量。由于我们计算了方向
1
上的梯度,因此我们必须寻找行向量。如果要重现整个渐变,可以通过以下方法完成
np.stack(((f[:,1,:] - f[:,0,:]), (f[:,2,:] - f[:,0,:])/2, (f[:,2,:] - f[:,1,:])), axis=1)
n暗箱
我们可以概括我们在这里学到的东西,以计算沿方向的任意函数的梯度。
def grad_along_axis(f, ax):
f_grad_ind = []
for i in range(f.shape[ax]):
if i == 0:
f_grad_ind.append(np.take(f, i+1, ax) - np.take(f, i, ax))
elif i == f.shape[ax] -1:
f_grad_ind.append(np.take(f, i, ax) - np.take(f, i-1, ax))
else:
f_grad_ind.append((np.take(f, i+1, ax) - np.take(f, i-1, ax))/2)
f_grad = np.stack(f_grad_ind, axis=ax)
return f_grad
哪里
np.take(f, i, ax) = f[:,...,i,...,:]
i
位于索引ax
。关于python - 对称函数的numpy.gradient的组件不同,我们在Stack Overflow上找到一个类似的问题:https://stackoverflow.com/questions/54891699/