我np.einsum
用来乘以概率表,例如:
np.einsum('ijk,jklm->ijklm', A, B)
问题是我正在处理超过 26 个随机变量(轴),所以如果我为每个随机变量分配一个字母,我的字母就会用完。有没有另一种方法可以指定上述操作来避免这个问题,而无需使用混乱的np.sum
操作np.dot
?
我np.einsum
用来乘以概率表,例如:
np.einsum('ijk,jklm->ijklm', A, B)
问题是我正在处理超过 26 个随机变量(轴),所以如果我为每个随机变量分配一个字母,我的字母就会用完。有没有另一种方法可以指定上述操作来避免这个问题,而无需使用混乱的np.sum
操作np.dot
?
简短的回答是,您可以使用 52 个字母(大写和小写)中的任何一个。这就是英语中的所有字母。任何更高级的轴名称都必须映射到这 52 个或一组等效的数字上。实际上,您将希望在任何一次einsum
通话中使用这 52 个中的一小部分。
@kennytm
建议使用替代输入语法。一些示例运行表明这不是解决方案。26 仍然是实际限制(尽管有可疑的错误消息)。
In [258]: np.einsum(np.ones((2,3)),[0,20],np.ones((3,4)),[20,2],[0,2])
Out[258]:
array([[ 3., 3., 3., 3.],
[ 3., 3., 3., 3.]])
In [259]: np.einsum(np.ones((2,3)),[0,27],np.ones((3,4)),[27,2],[0,2])
---------------------------------------------------------------------------
ValueError Traceback (most recent call last)
<ipython-input-259-ea61c9e50d6a> in <module>()
----> 1 np.einsum(np.ones((2,3)),[0,27],np.ones((3,4)),[27,2],[0,2])
ValueError: invalid subscript '|' in einstein sum subscripts string, subscripts must be letters
In [260]: np.einsum(np.ones((2,3)),[0,100],np.ones((3,4)),[100,2],[0,2])
---------------------------------------------------------------------------
ValueError Traceback (most recent call last)
<ipython-input-260-ebd9b4889388> in <module>()
----> 1 np.einsum(np.ones((2,3)),[0,100],np.ones((3,4)),[100,2],[0,2])
ValueError: subscript is not within the valid range [0, 52]
我不完全确定为什么您需要超过 52 个字母(大写和小写),但我确定您需要进行某种映射。您不想一次einsum
使用超过 52 个轴编写字符串。生成的迭代器将太大(对于内存或时间)。
我正在描绘某种可以用作的映射函数:
astr = foo(A.names, B.names)
# foo(['i','j','k'],['j','k','l','m'])
# foo(['a1','a2','a3'],['a2','a3','b4','b5'])
np.einsum(astr, A, B)
https://github.com/hpaulj/numpy-einsum/blob/master/einsum_py.py
是 Python 版本的einsum
. 粗略地说einsum
,解析下标字符串,创建一个op_axes
列表,可用于np.nditer
设置所需的积和计算。使用此代码,我可以查看翻译是如何完成的:
从__name__
块中的一个例子:
label_str, op_axes = parse_subscripts('ik,kj->ij', Labels([A.ndim,B.ndim]))
print op_axes
# [[0, -1, 1], [-1, 1, 0], [0, 1, -1]] fine
# map (4,newaxis,3)(newaxis,3,2)->(4,2,newaxis)
print sum_of_prod([A,B],op_axes)
您的示例,具有完整的诊断输出是
In [275]: einsum_py.parse_subscripts('ijk,jklm->ijklm',einsum_py.Labels([3,4]))
jklm
{'counts': {105: 1, 106: 2, 107: 2, 108: 1, 109: 1},
'strides': [],
'num_labels': 5,
'min_label': 105,
'nop': 2,
'ndims': [3, 4],
'ndim_broadcast': 0,
'shapes': [],
'max_label': 109}
[('ijk', [105, 106, 107], 'NONE'),
('jklm', [106, 107, 108, 109], 'NONE')]
('ijklm', [105, 106, 107, 108, 109], 'NONE')
iter labels: [105, 106, 107, 108, 109],'ijklm'
op_axes [[0, 1, 2, -1, -1], [-1, 0, 1, 2, 3], [0, 1, 2, 3, 4]]
Out[275]:
(<einsum_py.Labels at 0xb4f80cac>,
[[0, 1, 2, -1, -1], [-1, 0, 1, 2, 3], [0, 1, 2, 3, 4]])
使用'ajk,jkzZ->ajkzZ'
更改标签,但结果相同op_axes
。
这是翻译功能的初稿。它应该适用于任何列表列表(可散列项):
def translate(ll):
mset=set()
for i in ll:
mset.update(i)
dd={k:v for v,k in enumerate(mset)}
x=[''.join([chr(dd[i]+97) for i in l]) for l in ll]
# ['cdb', 'dbea', 'cdbea']
y=','.join(x[:-1])+'->'+x[-1]
# 'cdb,dbea->cdbea'
In [377]: A=np.ones((3,1,2),int)
In [378]: B=np.ones((1,2,4,3),int)
In [380]: ll=[list(i) for i in ['ijk','jklm','ijklm']]
In [381]: y=translate(ll)
In [382]: y
Out[382]: 'cdb,dbea->cdbea'
In [383]: np.einsum(y,A,B).shape
Out[383]: (3, 1, 2, 4, 3)
使用set
映射索引对象意味着最终的索引字符是无序的。只要您指定不应该成为问题的 RHS。我也忽略了ellipsis
。
==================
输入的列表版本einsum
被转换为einsum_list_to_subscripts()
(in numpy/core/src/multiarray/multiarraymodule.c
) 中的下标字符串版本。它替换ELLIPSIS
为“...”。如果( s < 0 || s > 2*26)
其中s
一个子列表中的数字在哪里,则会引发 [0,52] 错误消息。并转换s
为字符串
if (s < 26) {
subscripts[subindex++] = 'A' + s;
}
else {
subscripts[subindex++] = 'a' + s;
但看起来第二种情况不起作用;我收到类似 26 的错误。
ValueError: invalid subscript '{' in einstein sum subscripts string, subscripts must be letters
这'a'+s
是错误的,如果s>26
:
In [424]: ''.join([chr(ord('A')+i) for i in range(0,26)])
Out[424]: 'ABCDEFGHIJKLMNOPQRSTUVWXYZ'
In [425]: ''.join([chr(ord('a')+i) for i in range(0,26)])
Out[425]: 'abcdefghijklmnopqrstuvwxyz'
In [435]: ''.join([chr(ord('a')+i) for i in range(26,52)])
Out[435]: '{|}~\x7f\x80\x81\x82\x83\x84\x85\x86\x87\x88\x89\x8a\x8b\x8c\x8d\x8e\x8f\x90\x91\x92\x93\x94'
那'a'+s
是错误的;应该是:
In [436]: ''.join([chr(ord('a')+i-26) for i in range(26,52)])
Out[436]: 'abcdefghijklmnopqrstuvwxyz'
我提交了https://github.com/numpy/numpy/issues/7741
这么长时间后这个bug的存在表明子列表格式并不常见,并且在该列表中使用大数字的频率更低。
您可以使用einsum(op0, sublist0, op1, sublist1, ..., [sublistout])
表单代替i,j,ik->ijk
,API 不限于 52 个轴*。文档中显示了这种详细形式如何对应于 ijk 形式。
OP的
np.einsum('ijk,jklm->ijklm', A, B)
将被写为
np.einsum(A, [0,1,2], B, [1,2,3,4], [0,1,2,3,4])
(* 注意:实现仍然限于 26 轴。请参阅@hpaulj 的答案和他的错误报告以获取解释)
来自 numpy 示例的等价物:
>>> np.einsum('ii', a)
>>> np.einsum(a, [0,0])
>>> np.einsum('ii->i', a)
>>> np.einsum(a, [0,0], [0])
>>> np.einsum('ij,j', a, b)
>>> np.einsum(a, [0,1], b, [1])
>>> np.einsum('ji', c)
>>> np.einsum(c, [1,0])
>>> np.einsum('..., ...', 3, c)
>>> np.einsum(3, [...], c, [...])
>>> np.einsum('i,i', b, b)
>>> np.einsum(b, [0], b, [0])
>>> np.einsum('i,j', np.arange(2)+1, b)
>>> np.einsum(np.arange(2)+1, [0], b, [1])
>>> np.einsum('i...->...', a)
>>> np.einsum(a, [0, ...], [...])
>>> np.einsum('ijk,jil->kl', a, b)
>>> np.einsum(a, [0,1,2], b, [1,0,3], [2,3])