稀疏矩阵的存储

如果将整个稀疏矩阵都存储进入内存中的话, 那将会占用相当大的空间, 根据稀疏矩阵中非0元素在矩阵中的分布以及个数, 可以使用不同的数据结构来储存矩阵, 这样能大大减小占用的内存数量, 但这种轻量化存储的劣势是如果要访问其中某个特定位置的元素, 往往更复杂, 而且需要额外的structure从而将原始矩阵清晰地还原, 再去访问其中特定的元素。

大体来说, 促成稀疏矩阵的格式共有两大类:

  • 支持对原始稀疏矩阵高效修改的数据结构, 例如DOK (Dictionary of keys), LIL (List of lists), or COO (Coordinate list)等, 用这种数据结构可以很方便地对往原始矩阵中假如新的元素, 或者从原始矩阵中删去元素 (增加元素指的是将原来为0的元素变为非0, 同样删除同理)

  • 支持快速访问到特定位置的元素的数据结构以及支持基于矩阵形式的操作, 即两个稀疏矩阵相加或者相减等, 例如CSR (Compressed Sparse Row) 和 CSC (Compressed Sparse Column).

Dictionary of keys (DOK)

DOK 是一个字典, 其中每个key是(row_index, column_index), 而value则是在该位置上非0元素的值。这种格式的储存方式非常适合增删改, 即如果要将增加某个位置的非0元素直接向字典中加一个键值对就行, 同样删除同理, 但如果要查某一个位置的元素是否是非0元素, 那么就需要遍历整个字典的keys, 较为耗时, 而且字典的key是无序的, 如下为Scipy中实现DOK的代码, 注意DOK代码的意思是我们输入进去一个稀疏矩阵, scipy帮你把这个大的稀疏矩阵压缩成dok_matrix中的形式, 输出的mtx对象就是压缩后的DOK矩阵对象, 访问DOK矩阵中的元素的话需要用(row_index, col_index)的元组进行索引。

from scipy import sparse
import numpy as np
#除此之外dok_matrix还接收数组形式的array作为输入
mtx = sparse.dok_matrix((5, 5), dtype=np.float64)
print(mtx)

for ir in range(5):
    for ic in range(5):
        mtx[ir, ic] = 1.0 * (ir != ic)

print(mtx)
print(mtx.todense())

#结果是(注意第一个print输出的是空值, 因为其就构造了一个(5,5)的零矩阵, 其中没有非0元素自然没有非零元素键值对):


  (0, 1)    1.0
  (0, 2)    1.0
  (0, 3)    1.0
  (0, 4)    1.0
  (1, 0)    1.0
  (1, 2)    1.0
  (1, 3)    1.0
  (1, 4)    1.0
  (2, 0)    1.0
  (2, 1)    1.0
  (2, 3)    1.0
  (2, 4)    1.0
  (3, 0)    1.0
  (3, 1)    1.0
  (3, 2)    1.0
  (3, 4)    1.0
  (4, 0)    1.0
  (4, 1)    1.0
  (4, 2)    1.0
  (4, 3)    1.0

[[0. 1. 1. 1. 1.]
 [1. 0. 1. 1. 1.]
 [1. 1. 0. 1. 1.]
 [1. 1. 1. 0. 1.]
 [1. 1. 1. 1. 0.]]

如果想从mtx中取键值对的话, 可以采用如下方法:

print(mtx[1, 1])
>>> 0.0

print(mtx[1, 1:3])
>>> (0, 1)  1.0

print(mtx[1, 1:3].todense())
>>> [[0. 1.]]

print(mtx[[2,1], 1:3].todense())
#注意此处是先取列表中第一个元素2, 取得[2,1:3]的值, 然后再取第二个元素1, 取得[1,1:3]的值
>>> [[1. 0.]
    [0. 1.]]

List of lists (LIL)

LIL的数据结构是列表, 每一行对应着的矩阵中该行, 而LIL中每一行都是一个列表, 包含在该行出现非0元素的列位置以及非0元素的值, 可能的实现方式是每一行都包含一个列表, 该列表是个升序排序后的列表, 每个元素代表在原矩阵中在该行非0元素出现的列位置

Advantages of the LIL format

  • supports flexible slicing
  • changes to the matrix sparsity structure are efficient

Disadvantages of the LIL format

  • arithmetic operations LIL + LIL are slow (consider CSR or CSC)
  • slow column slicing (consider CSC)
  • slow matrix vector products (consider CSR or CSC)

Intended Usage

  • LIL is a convenient format for constructing sparse matrices
  • once a matrix has been constructed, convert to CSR or CSC format for fast arithmetic and matrix vector operations
  • consider using the COO format when constructing large matrices

Data Structure

  • An array (self.rows) of rows, each of which is a sorted list of column indices of non-zero elements.
  • The corresponding nonzero values are stored in similar fashion in self.data.

Scipy代码如下:

from scipy import sparse
from numpy.random import rand
import numpy as np
mtx = sparse.lil_matrix((4, 5))
data = np.round(rand(2, 3))

print(data)
>>> [[0. 1. 0.]
    [0. 1. 0.]]
mtx[:2, [1, 2, 3]] = data

print(mtx)
>>> (0, 2)  1.0
    (1, 2)  1.0

print(mtx.todense())
>>> [[0. 0. 1. 0. 0.]
     [0. 0. 1. 0. 0.]
     [0. 0. 0. 0. 0.]
     [0. 0. 0. 0. 0.]]

print(mtx.toarray())
>>> [[0. 0. 1. 0. 0.]
     [0. 0. 1. 0. 0.]
     [0. 0. 0. 0. 0.]
     [0. 0. 0. 0. 0.]]
mtx = sparse.lil_matrix([[0, 1, 2, 0], [3, 0, 1, 0], [1, 0, 0, 1]])

print(mtx.todense()) 
>>> [[0 1 2 0]
     [3 0 1 0]
     [1 0 0 1]]
     
print(mtx)
>>>   (0, 1)    1
      (0, 2)    2
      (1, 0)    3
      (1, 2)    1
      (2, 0)    1
      (2, 3)    1
      
print(mtx[:2, :])
>>> (0, 1)  1
    (0, 2)  2
    (1, 0)  3
    (1, 2)  1

print(mtx[:2, :].todense())
>>> [[0 1 2 0]
     [3 0 1 0]]

print(mtx[1:2, [0,2]].todense())
>>> [[3 1]]

print(mtx.todense())
>>> [[0 1 2 0]
     [3 0 1 0]
     [1 0 0 1]]

Coordinate list (COO)

COO是一个列表, 其元素为元组, 每个元组的格式是(row_index,column_index,value), 整个列表首先按row_index排序, 之后在此基础上再按column_index排序, 排序的意义在于可以更快的查找到特定元素, 同样这种数据结构很适合在构建matrix的时候往matrix中一步一步加元素, 即incremental matrix construction。

Advantages of the COO format

  • facilitates fast conversion among sparse formats
  • permits duplicate entries (see example)
  • very fast conversion to and from CSR/CSC formats

Disadvantages of the COO format

  • does not directly support:

    • arithmetic operations
    • slicing

Intended Usage

  • COO is a fast format for constructing sparse matrices
  • Once a matrix has been constructed, convert to CSR or CSC format for fast arithmetic and matrix vector operations
  • By default when converting to CSR or CSC format, duplicate (i,j) entries will be summed together. This facilitates efficient construction of finite element matrices and the like. (see example)

代码如下:

from scipy import sparse
import numpy as np
mtx = sparse.coo_matrix((3, 4), dtype=np.int8)
print(mtx.todense())
>>> [[0 0 0 0]
     [0 0 0 0]
     [0 0 0 0]]

row = np.array([0, 3, 1, 0])
col = np.array([0, 3, 1, 2])
data = np.array([4, 5, 7, 9])
mtx = sparse.coo_matrix((data, (row, col)), shape=(4, 4))
print(mtx)
>>>   (0, 0)    4
      (3, 3)    5
      (1, 1)    7
      (0, 2)    9

print(mtx.todense())
>>> [[4 0 9 0]
     [0 7 0 0]
     [0 0 0 0]
     [0 0 0 5]]

row = np.array([0, 0, 1, 3, 1, 0, 0])
col = np.array([0, 2, 1, 3, 1, 0, 0])
data = np.array([1, 1, 1, 1, 1, 1, 1])
mtx = sparse.coo_matrix((data, (row, col)), shape=(4, 4))
print(mtx.todense())
>>> [[3 0 1 0]
     [0 2 0 0]
     [0 0 0 0]
     [0 0 0 1]]

#不支持索引
print(mtx[2, 3])
>>> Traceback (most recent call last):
  File "/Users/shenyi/Documents/Untitled.py", line 21, in <module>
    print(mtx[2, 3])
TypeError: 'coo_matrix' object is not subscriptable

可以用scipy中的sparse中的csr_matrix和csc_matrix, 其中CSR代表Compressed Sparse Row matrix而CSC代表Compressed Sparse Column matrix

Compressed sparse row和Compressed sparse column

CSR或者 compressed row storage (CRS)用三个列表来表达稀疏矩阵, 其跟COO较为相似, 但比COO更进一步的是其还将每一行的空间压缩了, 这种类型的数据结构支持快速row access和矩阵的相加, 同理 CSC支持快速column access和矩阵的相加:

我们用CSR举例, 首先第一种最简单的即生成全是0的矩阵:

import numpy as np
from scipy.sparse import csr_matrix
csr_matrix((3, 4), dtype=np.int8).toarray()

#结果是:
>>> array([[0, 0, 0, 0],
           [0, 0, 0, 0],
           [0, 0, 0, 0]], dtype=int8)

如果想往这里边加非0元素, 可以传入这些非0元素的数值以及其所以在行和列的下标位置:

row = np.array([0, 0, 1, 2, 2, 2])
col = np.array([0, 2, 2, 0, 1, 2])
data = np.array([1, 2, 3, 4, 5, 6])
#即第(0,0)位置上的元素是1, 第(0,2)位置上的元素是2, 以此类推
#注意此处调用csr_matrix的时候如果是传入行和列的列表的话要以元组形式传入, 然后再将该元组与data组成一个新的元祖
csr_matrix((data, (row, col)), shape=(3, 3)).toarray()

#结果是:
>>> array([[1, 0, 2],
           [0, 0, 3],
           [4, 5, 6]])

另外一种调用方式是更复杂一些的, 首先传入数据列表, 接着传入indices和indptr, 如下:

indptr = np.array([0, 2, 3, 6])
indices = np.array([0, 2, 1, 0, 1, 2])
data = np.array([1, 2, 3, 4, 5, 6])
#注意此处的csr_matrix调用时传参方法与上述不同, 上面用col和row的方法中是将col和row打包成一个元组, 此处没有如此操作
csr_matrix((data, indices, indptr), shape=(3, 3)).toarray()

#结果是:
>>> array([[1, 0, 2],
           [0, 3, 0],
           [4, 5, 6]])

其中我们首先看indices, indices代表每一行的非0元素的位置, indptr列表中最后一个元素代表整个矩阵中所有非0元素的数量, 除了最后一个元素外, 其他每个元素从第一个元素开始依次指代每一行的数据在indices中开始的位置。

我们来看indptr中第一个元素, 第一个元素代表最后生成的矩阵的第一行, 而该元素等于0代表矩阵中第一行是从indices中第0个元素开始; indptr中第二个元素代表矩阵中的第二列, 而该元素等于2代表在indices这个列表中第二个元素是矩阵中第二列的起始元素

也就是indptr中第二个元素代表其对应着矩阵中第二列, 该元素等于2, 代表indices[2]是第二列的第一个非零元素的位置, 发现indices[2]=1, 说明在矩阵中该行第一个非0元素出现在该行的第一个位置(以0为起始下标), 因为indices和data的下标是一一对应的, 所以此时data[2]=3, 所以矩阵的(1,1)的值为3

附上一张图:

CSR

其中用csr矩阵的形式的优势和劣势为:

Advantages of the CSR format

  • efficient arithmetic operations CSR + CSR, CSR * CSR, etc.
  • efficient row slicing
  • fast matrix vector products

Disadvantages of the CSR format

  • slow column slicing operations (consider CSC)
  • changes to the sparsity structure are expensive (consider LIL or DOK)
>>> import numpy as np
>>> from scipy.sparse import csc_matrix
>>> csc_matrix((3, 4), dtype=np.int8).toarray()
array([[0, 0, 0, 0],
       [0, 0, 0, 0],
       [0, 0, 0, 0]], dtype=int8)

>>> row = np.array([0, 2, 2, 0, 1, 2])
>>> col = np.array([0, 0, 1, 2, 2, 2])
>>> data = np.array([1, 2, 3, 4, 5, 6])
>>> csc_matrix((data, (row, col)), shape=(3, 3)).toarray()
array([[1, 0, 4],
       [0, 0, 5],
       [2, 3, 6]])

>>> indptr = np.array([0, 2, 3, 6])
>>> indices = np.array([0, 2, 2, 0, 1, 2])
>>> data = np.array([1, 2, 3, 4, 5, 6])
>>> csc_matrix((data, indices, indptr), shape=(3, 3)).toarray()
array([[1, 0, 4],
       [0, 0, 5],
       [2, 3, 6]])

csc矩阵与csr矩阵的构造方式雷同, 只不过之前是一行一行开始扫描indptr和indices的, 现在变为一列一列了, 如下, 我们看indptr[1]=2, 代表indices中第2个元素开始是矩阵第二列的范畴, 而indices[2]=2, 所以说矩阵第二列中从上往下数第二个元素是data[2]=3, 所以矩阵(2,1)元素为3

>>> import numpy as np
>>> from scipy.sparse import csc_matrix
>>> csc_matrix((3, 4), dtype=np.int8).toarray()
array([[0, 0, 0, 0],
       [0, 0, 0, 0],
       [0, 0, 0, 0]], dtype=int8)

>>> row = np.array([0, 2, 2, 0, 1, 2])
>>> col = np.array([0, 0, 1, 2, 2, 2])
>>> data = np.array([1, 2, 3, 4, 5, 6])
>>> csc_matrix((data, (row, col)), shape=(3, 3)).toarray()
array([[1, 0, 4],
       [0, 0, 5],
       [2, 3, 6]])

>>> indptr = np.array([0, 2, 3, 6])
>>> indices = np.array([0, 2, 2, 0, 1, 2])
>>> data = np.array([1, 2, 3, 4, 5, 6])
>>> csc_matrix((data, indices, indptr), shape=(3, 3)).toarray()
array([[1, 0, 4],
       [0, 0, 5],
       [2, 3, 6]])

而其唯一区别即一个是slow column slicing operation一个是slow row slicing operation

如果想从matrix形式变回原来row和col的形式, 可以用如下方法:

Here's your array mtr:

In [11]: mtr Out[11]: \<3x3 sparse matrix of type '\<class 'numpy.int64'\>' with 6 stored elements in Compressed Sparse Row format\> In [12]: mtr.A Out[12]: array([[1, 0, 2], [0, 0, 3], [4, 5, 6]], dtype=int64)

Convert to COO format, and access the data, row and col attributes.

In [13]: c = mtr.tocoo() In [14]: c.data Out[14]: array([1, 2, 3, 4, 5, 6], dtype=int64) In [15]: c.row Out[15]: array([0, 0, 1, 2, 2, 2], dtype=int32) In [16]: c.col Out[16]: array([0, 2, 2, 0, 1, 2], dtype=int32)

或者用如下方法:

import numpy as np
from scipy.sparse import csr_matrix
row = np.array([0, 0, 1, 2, 2, 2])
col = np.array([0, 2, 2, 0, 1, 2])
data = np.array([1, 2, 3, 4, 5, 6])
mtr = csr_matrix((data, (row, col)))

print(mtr.todense())

rows, cols = mtr.nonzero()
data = mtr[rows, cols]

print(rows, cols, data)

#结果是:
[[1 0 2]
 [0 0 3]
 [4 5 6]]
[0 0 1 2 2 2] [0 2 2 0 1 2] [[1 2 3 4 5 6]]
最后编辑于 ：
©著作权归作者所有,转载或内容合作请联系作者
  • 人面猴
    序言：七十年代末，一起剥皮案震惊了整个滨河市，随后出现的几起案子，更是在滨河造成了极大的恐慌，老刑警刘岩，带你破解...
    沈念sama阅读 159,117评论 4赞 362
  • 死咒
    序言：滨河连续发生了三起死亡事件，死亡现场离奇诡异，居然都是意外死亡，警方通过查阅死者的电脑和手机，发现死者居然都...
    沈念sama阅读 67,328评论 1赞 293
  • 救了他两次的神仙让他今天三更去死
    文/潘晓璐 我一进店门，熙熙楼的掌柜王于贵愁眉苦脸地迎上来，“玉大人，你说我怎么就摊上这事。” “怎么了？”我有些...
    开封第一讲书人阅读 108,839评论 0赞 243
  • 道士缉凶录：失踪的卖姜人
    文/不坏的土叔 我叫张陵，是天一观的道长。 经常有香客问我，道长，这世上最难降的妖魔是什么？ 我笑而不...
    开封第一讲书人阅读 44,007评论 0赞 206
  • 港岛之恋（遗憾婚礼）
    正文 为了忘掉前任，我火速办了婚礼，结果婚礼上，老公的妹妹穿的比我还像新娘。我一直安慰自己，他们只是感情好，可当我...
    茶点故事阅读 52,384评论 3赞 287
  • 恶毒庶女顶嫁案：这布局不是一般人想出来的
    文/花漫 我一把揭开白布。 她就那样静静地躺着，像睡着了一般。 火红的嫁衣衬着肌肤如雪。 梳的纹丝不乱的头发上，一...
    开封第一讲书人阅读 40,629评论 1赞 219
  • 城市分裂传说
    那天，我揣着相机与录音，去河边找鬼。 笑死，一个胖子当着我的面吹牛，可吹牛的内容都是我干的。 我是一名探鬼主播，决...
    沈念sama阅读 31,880评论 2赞 313
  • 双鸳鸯连环套：你想象不到人心有多黑
    文/苍兰香墨 我猛地睁开眼，长吁一口气：“原来是场噩梦啊……” “哼！你这毒妇竟也来了？” 一声冷哼从身侧响起，我...
    开封第一讲书人阅读 30,593评论 0赞 198
  • 万荣杀人案实录
    序言：老挝万荣一对情侣失踪，失踪者是张志新（化名）和其女友刘颖，没想到半个月后，有当地人在树林里发现了一具尸体，经...
    沈念sama阅读 34,313评论 1赞 243
  • 护林员之死
    正文 独居荒郊野岭守林人离奇死亡，尸身上长有42处带血的脓包…… 初始之章·张勋 以下内容为张勋视角 年9月15日...
    茶点故事阅读 30,575评论 2赞 246
  • 白月光启示录
    正文 我和宋清朗相恋三年，在试婚纱的时候发现自己被绿了。 大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
    茶点故事阅读 32,066评论 1赞 260
  • 活死人
    序言：一个原本活蹦乱跳的男人离奇死亡，死状恐怖，灵堂内的尸体忽然破棺而出，到底是诈尸还是另有隐情，我是刑警宁泽，带...
    沈念sama阅读 28,392评论 2赞 253
  • 日本核电站爆炸内幕
    正文 年R本政府宣布，位于F岛的核电站，受9级特大地震影响，放射性物质发生泄漏。R本人自食恶果不足惜，却给世界环境...
    茶点故事阅读 33,052评论 3赞 236
  • 男人毒药：我在死后第九天来索命
    文/蒙蒙 一、第九天 我趴在偏房一处隐蔽的房顶上张望。 院中可真热闹，春花似锦、人声如沸。这庄子的主人今日做“春日...
    开封第一讲书人阅读 26,082评论 0赞 8
  • 一桩弑父案，背后竟有这般阴谋
    文/苍兰香墨 我抬头看了看天上的太阳。三九已至，却和暖如春，着一层夹袄步出监牢的瞬间，已是汗流浃背。 一阵脚步声响...
    开封第一讲书人阅读 26,844评论 0赞 195
  • 情欲美人皮
    我被黑心中介骗来泰国打工， 没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留，地道东北人。 一个月前我还...
    沈念sama阅读 35,662评论 2赞 274
  • 代替公主和亲
    正文 我出身青楼，却偏偏与公主长得像，于是被迫代替她去往敌国和亲。 传闻我的和亲对象是个残疾皇子，可洞房花烛夜当晚...
    茶点故事阅读 35,575评论 2赞 270