在numpy中，从另一个多维数组中的二维指数数组创建多维数组。

Question

我有一个问题，我正在努力解决。

我有一个形状为(3,4,5)的数组，其中包含4个指标和5个区域的3个时间步长的数值。在这个例子中，让我们说这个数组是。

values = np.arange(60).reshape(3,4,5)

我还有一个 "地图"，它显示了5个区域在5x5网格上的位置，例如：

map = np.array([[1,1,1,2,2],
               [1,1,2,2,2],
               [3,3,3,2,2],
               [3,3,4,4,4],
               [3,3,5,5,4]])

我想最终得到一个数组(result)，在5x5网格上显示数值（类似于 map)，每一个时间步长(3)和指标(4)，所以

>>>result.shape 
(3,4,5,5)

我希望这有意义。

Answer 1

如果我的理解正确，5个区域中的 values.shape 不涉及 (5, 5) 在 map.shape中的数值范围，但只限于 map. 首先要做的是使值在 map 遵循传统的 python 索引惯例：从零开始计数，而不是从一开始。

下面是如何在map与值的形状不同的情况下做到这一点。

import numpy as np

m, n, p = 2, 3, 4
values = np.arange(m*n*p).reshape(m, n, p)

map = np.array([
     [0, 0, 0, 3, 3, 3],
     [0, 0, 1, 3, 3, 2],
     [1, 1, 1, 3, 2, 2],
     [2, 1, 2, 2, 2, 1],
     [2, 2, 2, 2, 1, 1]
     ])

q, r = map.shape

map_r = map.reshape(q*r)

# result_r: shape (q*r, m, n)
# If you want to keep the 1-based numbering, index as [:, :, i-1].
result_r = np.array([values[:, :, i] for i in map_r])
result = result_r.transpose(1, 2, 0).reshape(m, n, q, r)

# test
jq, jr = 4, 5
mapval = map[jq, jr]
print(f'map[{jq}, {jr}] = {mapval}')
print(f'result[:, :, {jq}, {jr}] =\n{result[:, :, jq, jr]}')
print(f'values[:, :, {mapval}] =\n{values[:, :, mapval]}')

输出

map[4, 5] = 1
result[:, :, 4, 5] =
[[ 1  5  9]
 [13 17 21]]
values[:, :, 1] =
[[ 1  5  9]
 [13 17 21]]

旁白： 多维数组如果定义最左边的索引来代表一个更高层次的概念，那么它的工作通常会更简洁。如果你可以将其表述为 "M个小部件的N个列表"，那么自然的数组形状将是 (N, M)，不 (M,N). 所以，你会有 values.shape = (p, m, n) 和 result.shape (q,r,m,n) (也许还可以交换 m 和 n 以及，这取决于您以后如何处理数据）。) 您将避免大量的 [:, :] 指数和 transpose 操作。

Answer 2

下面是一个使用np.pad的实现。

values = np.arange(60).reshape(3,4,5)

# I renamed map to mapping, because map is a basic function in python

mapping = np.array([[1,1,1,2,2],
                  [1,1,2,2,2],
                  [3,3,3,2,2],
                  [3,3,4,4,4],
                  [3,3,5,5,4]])

v0,v1,v2 = values.shape
m0,m1 = mapping.shape

new = np.zeros((v0,v1,m0,m1))

for a in range(v2):
    q = values[:,:,a].reshape(v0,v1,1)
    r,c = np.where((mapping-1)==a)
    npad = ((0,0), (0,0), (0,len(c)-1))
    new[:,:,r,c] = np.pad(q, pad_width=npad, mode='symmetric')