抱歉标题含糊不清。我有两个相关的问题。
首先，假设我有一个函数“hessian”，给定两个参数 (x, y) 返回一个矩阵。我现在想计算在二维空间上运行的 (x,y) 矩阵。我想做类似的事情:
x = linspace(1, 4, 100).reshape(-1,1)
y = linspace(1, 4, 100).reshape(1,-1)
H = vectorize(hessian)(x, y)
生成的 H 形状为 (100,100,2,2)。以上不起作用(ValueError:使用序列设置数组元素)。我唯一想出的是
H = array([ hessian(xx, yy) for xx in x.flat for yy in y.flat ]).reshape(100,100,2,2)
有没有更好、更直接的方法？
其次，现在 H 的形状为 (100,100,2,2) 并且 dominant_eigenvector(X) 完全符合您的想法。
U, V = hsplit(array(map(dominant_eigenvector, H.reshape(10000,2,2))), 2)
我再次需要使用列表理解来进行迭代并将结果重新打包到一个手动指定形状的数组中。有没有更直接的方法可以达到同样的效果？
谢谢!
编辑:按照 Paul 和 JoshAdel 的建议，我实现了一个适用于数组的 hessian 版本，就在这里
def hessian(w1, w2):
w1 = atleast_1d(w1)[...,newaxis,newaxis]
w2 = atleast_1d(w2)[...,newaxis,newaxis]
o1, o2 = ix_(*map(xrange, Z.shape))
W = Z * pow(w1, o1) * pow(w2, o2)
A = (W).sum()
B = (W * o1).sum()
C = (W * o2).sum()
D = (W * o1 * o1).sum()
E = (W * o1 * o2).sum()
F = (W * o2 * o2).sum()
return array([[ D/A - B/A*B/A, E/A - B/A*C/A ],
[ E/A - B/A*C/A, F/A - C/A*C/A ]])
Z 可以被认为是一个大约 250x150 的全局数组。 o1 和 o2 索引 Z 的两个维度来计算 像 $\sum_{i,j} Z_{ij} * i * j$ 这样的东西。
这个版本的问题是中间数组 太大了。如果 w1 和 w2 是像 w1_k w2_l 这样的数组 W 变为 W_{k,l,i,j}，numpy 对此给出 ValueError: too big。
请您参考如下方法:
您可以尝试使用 meshgrid，也许您必须将 xn、yn 展平:
x = linspace(1, 4, 100)
y = linspace(1, 4, 100)
xn,yn=meshgrid(x,y)
H = vectorize(hessian)(xn, yn)
