加入收藏 | 设为首页 | 会员中心 | 我要投稿 李大同 (https://www.lidatong.com.cn/)- 科技、建站、经验、云计算、5G、大数据,站长网!
当前位置: 首页 > 编程开发 > Python > 正文

python – numpy多维索引和对角线对称

发布时间:2020-12-20 13:42:28 所属栏目:Python 来源:网络整理
导读:我有一个非常大的numpy阵列…… power = ...print power.shape (3,10,19,75,10) 这是对称的w.r.t. 10×10部分,即以下2-d矩阵是对称的 power[i,:,j,k,l,m]power[i,m,:] 对于i,m的所有值 我可以利用这个因子-4的增益吗?例如.将矩阵保存到文件时(使用savez_com
我有一个非常大的numpy阵列……

power = ...
print power.shape
>>> (3,10,19,75,10)

这是对称的w.r.t. 10×10部分,即以下2-d矩阵是对称的

power[i,:,j,k,l,m]
power[i,m,:]

对于i,m的所有值

我可以利用这个因子-4的增益吗?例如.将矩阵保存到文件时(使用savez_compressed 50 mb)

我的尝试:

size = 10
row_idx,col_idx = np.tril_indices(size)
zip_idx = zip(row_idx,col_idx)
print len(zip_idx),zip_idx[:5]
>>> 55 [(0,0),(1,1),(2,1)]
all_idx = [(r0,c0,r1,c1) for (r0,c0) in zip_idx for (r1,c1) in zip_idx]
print len(all_idx),all_idx[:5]
>>> 3025 [(0,(0,1,2,1)]
a,b,c,d = zip(*all_idx)
tril_part = np.transpose(s.power,3,4,5,6))[:,a,d]
print tril_part.shape
>>> (3,3025)

这看起来很难看,但“有效”……一旦我也可以从tril_part恢复权力……

我想这会产生两个问题:

>从权力转向tril_part的更好方式?
>如何从tril_part转向权力?

编辑:“大小”评论显然是有效的,但请忽略它:-)恕我直言,问题的索引部分是独立的.我一直发现自己想要为较小的矩阵做类似的索引.

解决方法

你是在正确的道路上.使用np.tril_indices,您确实可以巧妙地索引这些较低的三角形.还有待改进的是数据的实际索引/切片.

请试试这个(复制和粘贴):

import numpy as np
shape = (3,10)
p = np.arange(np.prod(shape)).reshape(shape)  # this is not symmetric,but not important

ix,iy = np.tril_indices(10)
# In order to index properly,we need to add axes. This can be done by hand or with this
ix1,ix2 = np.ix_(ix,ix)
iy1,iy2 = np.ix_(iy,iy)

p_ltriag = p[:,ix1,iy1,ix2,iy2]
print p_ltriag.shape  # yields (55,55,75),axis order can be changed if needed

q = np.zeros_like(p)
q[:,iy2] = p_ltriag  # fills the lower triangles on both sides
q[:,iy2,ix2] = p_ltriag  # fills the lower on left,upper on right
q[:,iy2] = p_ltriag  # fills the upper on left,lower on right
q[:,ix2] = p_ltriag  # fills the upper triangles on both sides

数组q现在包含p的对称化版本(其中上三角形被下三角形的内容替换).注意,最后一行包含反转顺序的iy和ix索引,实质上是创建下三角矩阵的转置.

下三角形的比较
为了比较,我们将所有上三角形设置为0

ux,uy = np.triu_indices(10)
p[:,ux,uy] = 0
q[:,uy] = 0
p[:,uy] = 0

print ((p - q) ** 2).sum()  # euclidean distance is 0,so p and q are equal

print ((p ** 2).sum(),(q ** 2).sum())  # prove that not all entries are 0 ;) - This has a negative result due to an overflow

(编辑:李大同)

【声明】本站内容均来自网络,其相关言论仅代表作者个人观点,不代表本站立场。若无意侵犯到您的权利,请及时与联系站长删除相关内容!

    推荐文章
      热点阅读