Python中的稀疏3d矩阵/数组？

Question

在scipy中,我们可以使用scipy.sparse.lil_matrix()等构造一个稀疏矩阵.但矩阵在2d.

我想知道在Python中是否存在稀疏3d矩阵/数组(张量)的现有数据结构？

ps我在3d中有很多稀疏数据,需要一个张量来存储/执行乘法.如果没有现有的数据结构,是否有任何建议来实现这样的张量？

Answer 1

很高兴建议一个(可能是显而易见的)实现,如果你有新的依赖项的时间和空间,可以在纯Python或C/Cython中进行,并且需要它更快.

N维中的稀疏矩阵可以假设大多数元素都是空的,因此我们使用键入元组的字典:

class NDSparseMatrix:
  def __init__(self):
    self.elements = {}

  def addValue(self, tuple, value):
    self.elements[tuple] = value

  def readValue(self, tuple):
    try:
      value = self.elements[tuple]
    except KeyError:
      # could also be 0.0 if using floats...
      value = 0
    return value

你会像这样使用它:

sparse = NDSparseMatrix()
sparse.addValue((1,2,3), 15.7)
should_be_zero = sparse.readValue((1,5,13))

您可以通过验证输入实际上是一个元组,并且它只包含整数来使这个实现更加健壮,但这只会减慢速度,所以除非您以后将代码发布到全世界,否则我不会担心.

编辑 - 矩阵乘法问题的Cython实现,假设其他张量是N维NumPy数组(numpy.ndarray)可能如下所示:

#cython: boundscheck=False
#cython: wraparound=False

cimport numpy as np

def sparse_mult(object sparse, np.ndarray[double, ndim=3] u):
  cdef unsigned int i, j, k

  out = np.ndarray(shape=(u.shape[0],u.shape[1],u.shape[2]), dtype=double)

  for i in xrange(1,u.shape[0]-1):
    for j in xrange(1, u.shape[1]-1):
      for k in xrange(1, u.shape[2]-1):
        # note, here you must define your own rank-3 multiplication rule, which
        # is, in general, nontrivial, especially if LxMxN tensor...

        # loop over a dummy variable (or two) and perform some summation:
        out[i,j,k] = u[i,j,k] * sparse((i,j,k))

  return out

虽然你总是需要手动解决这个问题,因为(如代码注释中所述)你需要定义你正在总结的索引,并且要小心数组长度或事情不起作用!

编辑2 - 如果另一个矩阵也稀疏,那么你不需要进行三向循环:

def sparse_mult(sparse, other_sparse):

  out = NDSparseMatrix()

  for key, value in sparse.elements.items():
    i, j, k = key
    # note, here you must define your own rank-3 multiplication rule, which
    # is, in general, nontrivial, especially if LxMxN tensor...

    # loop over a dummy variable (or two) and perform some summation 
    # (example indices shown):
    out.addValue(key) = out.readValue(key) + 
      other_sparse.readValue((i,j,k+1)) * sparse((i-3,j,k))

  return out

我对C实现的建议是使用一个简单的结构来保存索引和值:

typedef struct {
  int index[3];
  float value;
} entry_t;

然后,您需要一些函数来分配和维护这些结构的动态数组,并根据需要快速搜索它们; 但是你应该在担心这些东西之前测试Python实现的性能.

Answer 2

截至 2017 年，另一种答案是sparse包。根据包本身，它在 NumPy 之上并scipy.sparse通过概括scipy.sparse.coo_matrix布局实现了稀疏多维数组。

这是从文档中获取的示例：

import numpy as np
n = 1000
ndims = 4
nnz = 1000000
coords = np.random.randint(0, n - 1, size=(ndims, nnz))
data = np.random.random(nnz)

import sparse
x = sparse.COO(coords, data, shape=((n,) * ndims))
x
# <COO: shape=(1000, 1000, 1000, 1000), dtype=float64, nnz=1000000>

x.nbytes
# 16000000

y = sparse.tensordot(x, x, axes=((3, 0), (1, 2)))

y
# <COO: shape=(1000, 1000, 1000, 1000), dtype=float64, nnz=1001588>

Answer 3

看看Python中的sparray - 稀疏n维数组(作者Jan Erik Solem).也可以在github上找到.