SciPy稀疏矩阵实战：从格式选型到性能优化，告别内存爆炸

做数据密集型任务的人迟早会碰到同一个问题矩阵太大内存装不下计算还慢得离谱。我最早遇到这个状况是在构建用户行为特征矩阵的时候——几百万行、几千万列如果老老实实用 NumPy 的稠密二维数组存光一个矩阵就能吃掉几十个 GB 内存还没开始训练模型机器就卡死了。当时一位老工程师随口提醒我一句你这不是稠密数据换成 SciPy 稀疏矩阵吧。就是那次之后我才真正理解数据结构选型比优化代码重要得多这句话的分量。这篇内容面向所有需要处理大规模表格化数据的人不管你做的是推荐系统、文本挖掘、图结构分析还是有限元仿真这类典型的数值计算只要你发现手里的矩阵大部分元素是 0scipy.sparse就是绕不开的核心工具。它能把内存占用砍掉几个数量级运算速度也能明显提升关键是使用起来并没有想象中复杂。接下来我会把格式选型、几种常见构建方式、运算技巧和实际项目里踩过的坑一次讲透希望能帮你少走点弯路。1. 先看痛点为什么矩阵内存会“爆炸”1.1 稠密矩阵的内存账本很多刚接触稀疏矩阵的同学第一个困惑是同样是放着这么多数据凭什么换个存储格式就能省下几十 GB 内存要理解这一点得先算一笔账。一个很常见的场景是商品推荐。假设我们有 50 万个用户、200 万件商品要构建一个“用户-商品交互矩阵”行代表用户列代表商品格子里的数字表示用户对商品的打分或者点击次数。如果用 NumPy 的稠密二维数组来存元素个数是50 万 × 200 万 1 万亿个元素如果每个元素是 8 字节的float64内存需求大约是 8 TB。这个数字对一个单体机器来说基本是灾难级别的存在现实里根本没法直接操作。但仔细想想一个普通用户这辈子能接触的商品数量才多少个就算买了 500 件商品这 500 万个格子中间也只有一个格子有值其余 499 万多个格子全是 0。整张矩阵的稀疏度可能高达 99.99% 以上。稀疏矩阵想做的事情说白了就一句话别再傻傻地给 0 元素留位置只记录那些真正有数据的格子。依然拿上面这个例子来看平均每个用户有 500 个交互记录那么实际需要存储的数据量是50 万用户 × 500 个记录 ≈ 2.5 亿个元素就算每个元素加上行索引、列索引的额外开销内存也只是几百 MB 的级别。一句话总结稀疏矩阵的核心思想用索引信息换存储空间。1.2 什么时候才值得转成稀疏矩阵不是所有矩阵都适合用稀疏存储。如果矩阵比较小或者非零元素比例很高强行用稀疏矩阵反而可能会更慢、更占内存。这里给一个我常用的粗略判断标准矩阵规模超过 1000 × 1000非零元素占比低于 10%逻辑上是天然稀疏的数据例如用户行为矩阵、文本词频矩阵、图邻接矩阵、有限元刚度矩阵。真实项目里我倾向直接看非零元素占比。如果矩阵中 0 的比例超过 90%那基本就不用犹豫了。但如果是一个 100 × 100 的小矩阵哪怕里面全是 0转不转稀疏矩阵意义都不大直接用 NumPy 处理反而方便。有一点我想强调稀疏矩阵在科学计算里几乎是必备品。比如有限元模拟里刚度矩阵和自由度编号相关绝大多数节点只跟相邻节点有交互矩阵稀疏度非常高再比如 PageRank 这类图算法图中每个节点连出去的边都很有限。这些场景中如果没有稀疏矩阵计算规模会直接收敛到无法进行。2. SciPy 稀疏矩阵格式大盘点选对格式才高效scipy.sparse这个模块里有六七种常见的稀疏矩阵格式刚上手的人很容易被这些缩写搞晕。我不打算把文档翻译一遍只想讲清楚每种格式的存储逻辑、优缺点和适合的场合这在实际项目中才是取舍的关键。2.1 COO 与 CSR最常用的两个主力格式COOCoordinate Format坐标格式的思路最直白同时记录三个等长数组——行索引row、列索引col、数据值data。每一个“非零元素”对应三元组(行号, 列号, 数值)。这种格式非常接近我们的直觉尤其在从文件读数据的时候天然就长这个样子。import numpy as np from scipy import sparse # 三个等长数组分别存储行号、列号、数值 row np.array([0, 1, 2, 2, 1]) col np.array([1, 0, 3, 1, 2]) data np.array([4.0, 5.0, 6.0, 7.0, 8.0]) coo sparse.coo_matrix((data, (row, col)), shape(3, 4)) print(coo.toarray())CSRCompressed Sparse Row压缩行存储是工程实践里最常见的格式。它把 COO 中大量重复的行号信息进行压缩。存储上分为三块data所有非零元素的数值indices每个非零元素对应的列索引indptr表示每一行第一个非零元素在data中的起始位置。很多人看到indptr就发怵其实它就是一个分段点数组。比如indptr [0, 2, 4, 5]表示第 0 行的元素在data[0:2]第 1 行的元素在data[2:4]第 2 行的元素在data[4:5]。这个设计让 CSR 做“按行切分”和“矩阵乘法”特别顺手。为什么强调 CSR 好用因为大多数机器学习算法和数值求解器在读取矩阵时底层循环都是按照“一行一行”的思路去访问的。CSR 把行连续地存放在内存里CPU 的缓存命中率高访问速度自然就快。所以我的经验是多数场景直接无脑选 CSR等遇到列访问频繁的场景再考虑 CSC。2.2 CSC、LIL、DOK、BSR、DIA各自解决什么问题CSCCompressed Sparse Column是 CSR 的“转置版本”压缩的是列信息。如果你经常需要按列切片、按列统计比如计算文本矩阵中每个词的逆文档频率CSC 会更快。LILList of Lists的底层用 Python 列表来存储每行的非零元素。它的优势是支持灵活的切片赋值非常方便逐步构建矩阵。工作室里调试小规模稀疏矩阵时我经常用lil_matrix因为可以直接给某个位置赋值lil sparse.lil_matrix((3, 4)) lil[0, 2] 10DOKDictionary of Keys基于 Python 字典实现键是(行, 列)值是数据。它的优势同样是赋值方便特别适合小规模动态构建。缺点也很明显Python 循环和字典在这些场景下比较慢不适合大规模计算。BSRBlock Sparse Row适合非零元素以小块矩阵形式出现的场景比如有限元计算中的单元刚度矩阵直接按块存储能大幅减少索引开销。DIADiagonal Format专门用来存对角带状矩阵比如偏微分方程里常见的五点差分格式矩阵。如果矩阵的非零元素都集中在对角线附近DIA 的效率极高但这在通用场景里比较少见所以我用得不多。下面这个表是项目里经常拿来做选型参考的建议收藏格式存储结构优点适用场景COO三元组构建简单适合从数据源直接构造数据处理流水线的中间格式CSR压缩行行访问快矩阵乘法高效通用计算、机器学习、求解器CSC压缩列列访问快按列统计、列切片LIL嵌套列表支持切片赋值动态构建阶段DOK字典单点赋值灵活小规模动态构建BSR分块压缩分块运算内存友好有限元类分块矩阵DIA对角存储对角线密集时高效偏微分方程离散矩阵一个关键认知格式之间不是彼此替代而是配合使用。我个人的固定操作路线一直是先用 COO 或者 LIL/DOK 把数据构建好然后统一转换成 CSR交给后续计算链路。这样既发挥了 COO 构建快的优势又在后续矩阵运算时享受到 CSR 的高性能。3. 构建与运算实操别把稀疏矩阵用成稠密矩阵3.1 三种构建方式与性能对比方式一从稠密数组直接转换新手最容易写出这种代码dense np.array([[1, 0, 0], [0, 2, 0], [0, 0, 3]]) sp sparse.csr_matrix(dense)这种方式有个大问题你的稠密矩阵本身就占满了内存转过去之后内存峰值依然很高。我只有在验证算法正确性的时候会拿小规模数据这么做生产环境几乎不碰。方式二先坐标三元组再转 CSR这是我在真实项目里最喜欢的方式。不管数据来自数据库还是日志文件最终都能整理成(行号, 列号, 值)的三列形式。直接把这三列喂给coo_matrix构建很快然后再转成 CSRrows np.array([0, 0, 1, 2, 2]) cols np.array([1, 3, 0, 2, 3]) vals np.array([1.0, 2.0, 3.0, 4.0, 5.0]) coo sparse.coo_matrix((vals, (rows, cols)), shape(1000, 1000)) csr coo.tocsr()这里有个我踩过的坑如果同一个(行, 列)位置出现多次coo_matrix不会自动做搞笑相加它会先保留重复项直到你调用tocsr()的时候才把重复项合并。2018 年刚上手那会儿我完全不知道这个细节导致排查了整整一下午最后发现是数据源里同一对用户和商品产生了多次行为记录。所以现在只要是从明细数据构建矩阵我一定会对三元组先做一次去重聚合避免后面求解结果莫名其妙地出错。方式三动态构建用 LIL/DOK如果是在循环里逐步填充矩阵比如在迭代算法中实时更新某个权重矩阵用lil_matrix和dok_matrix会更顺手mat sparse.lil_matrix((5, 5)) for i in range(5): mat[i, (i 1) % 5] i result mat.tocsr()这里要注意动态构建只是构建阶段用 LIL/DOK一旦构建完成请立刻转为 CSR/CSC。因为 LIL 和 DOK 做算术运算性能很差会把 Python 层的开销拉满。3.2 高效运算加法、乘法、求和、切片稀疏矩阵最爽的地方在于常规的线性代数操作基本都有现成实现而且会主动避开 0 元素的无效计算。矩阵乘法是让我经验最深刻的一环。比如计算用户行为矩阵 A 和它的转置相乘用来得到用户之间的相似度矩阵# A用户-物品矩阵CSR 格式 similarity A A.T这种操作如果先转成稠密矩阵再算内存立刻爆炸但用稀疏矩阵的话由于矩阵本身很“空”中间结果依然会保持稀疏结构。我自己在一个千万级别的用户评分矩阵上做过尝试用A A.T得到的相似度矩阵依然是稀疏的内存开销比稠密版本低了几个量级。你有没有想过为什么 Sparsity 能带来这种收益因为稀疏矩阵乘法在底层只对非零项做操作。比如计算一行乘以一列时它只在两个非零项都存在的条件下才做乘积其他情况直接跳过。这本质上是在用稀疏性换取线性时间的缩减效果特别明显。切片操作也有一些门道。CSR 格式做行切片非常直接因为它本来就是按行压缩的row_slice csr[10:20, :] # 取 10 到 19 行非常快但如果用 CSR 做列切片性能就差很多了因为它需要遍历列索引并重构矩阵。这种时候可以先将矩阵转成 CSC再做列切片csc_mat csr.tocsc() col_slice csc_mat[:, 5] # 取某一列我现在的习惯是预先把矩阵按照使用方式选好格式而不是在每步操作里反复转换。后面的排查环节我会再展开讲。求和与统计记住一个 API——sum()返回的是另一个矩阵不是 NumPy 数组。比如对行求和row_sums csr.sum(axis1) # 得到的是 np.matrix 类型 row_sums_array np.asarray(row_sums).ravel()不少同事第一次写类似代码直接拿着这个返回值去和数组做拼接结果报维度错误非常浪费时间。这里的细节值得警惕。3.3 什么时候该转回稠密矩阵经常有人问我“既然稀疏矩阵这么好是不是所有环节都要用它”答案是否定的。稀疏矩阵的索引访问、随机访问、逐元素复杂计算性能都很差。比如要对矩阵的每一个非零元素做log(1x)变换因为没法用向量化的 NumPy 操作直接对整个数据块处理性能就出不来虽然 SciPy 在部分元素级操作上做了优化但总归不如连续内存的稠密数组。我这里有个比较务实的判断逻辑矩阵很大、很稀疏而且只需要做线性代数运算或整体统计时保持稀疏数据已经降到几千行、几千列且后续逻辑需要频繁随机索引访问、需要和机器学习模型的特征矩阵对齐时输出前转回稠密矩阵只要还在大数据链路中转稠密就是最后一步不该提前发生。4. 实战案例线性求解、特征分解和机器学习特征矩阵4.1 稀疏线性方程组求解在我参与的很多科学计算项目中最核心的问题最后都落在求解大型稀疏线性方程组上形状一般是A x b其中 A 是稀疏矩阵b 是已知向量。如果 A 是几十万行、几十万列的稠密矩阵scipy.linalg.solve会直接把内存耗尽但使用scipy.sparse.linalg.spsolve它可以利用 A 的稀疏结构快速完成 LU 分解并求解。from scipy.sparse import csr_matrix from scipy.sparse.linalg import spsolve import numpy as np # 假设这是一个 5 阶稀疏矩阵 A_dense np.array([ [4.0, 1.0, 0, 0, 0], [1.0, 4.0, 1.0, 0, 0], [0, 1.0, 4.0, 1.0, 0], [0, 0, 1.0, 4.0, 1.0], [0, 0, 0, 1.0, 4.0], ]) A csr_matrix(A_dense) b np.array([1.0, 2.0, 3.0, 4.0, 5.0]) x spsolve(A, b) print(x)如果是很大的矩阵默认的直接法可能不够快这时候可以考虑迭代法例如bicgstab、gmres等from scipy.sparse.linalg import bicgstab x, info bicgstab(A, b, tol1e-8, maxiter1000)迭代法的优势是内存占用更低且可以方便地设置收敛精度。但迭代法对矩阵的条件数敏感条件数大的时候可能收敛很慢。这时候可以先做预处理例如用scipy.sparse.linalg.spilu做不完全 LU 分解作为预处理矩阵。这块的调参经验不少往往需要根据具体问题多试几次。4.2 特征值计算在图和网络分析里的应用大型网络的谱聚类、PageRank 等算法核心步骤经常是求解稀疏矩阵的部分特征值和特征向量。实际数据里网络邻接矩阵的规模非常庞大稠密方式求特征值会全面失落。scipy.sparse.linalg.eigs支持只求最大的几个特征值且内部使用迭代法不会生成完整稠密矩阵。from scipy.sparse.linalg import eigs # 假设 A 是邻接矩阵转换的稀疏矩阵 eigenvalues, eigenvectors eigs(A.astype(float), k10, whichLM) print(eigenvalues)这里有一个经常出错的点eigs默认处理的是复数特征值问题如果确认矩阵是实数对称的记得使用eigsh否则计算和输出都可能不符合预期。eigsh对内存和时间的消耗都更友好别选错了入口。我在做网络聚类时流程基本是这样构建图的邻接矩阵或 Laplacian 矩阵使用 COO 构建后转 CSR使用eigsh求前 k 个特征向量把特征向量作为节点的特征用简单聚类完成切割。4.3 机器学习中的稀疏特征矩阵如果你做文本分类一定会碰到 TF-IDF 特征矩阵。假设我们有两百万篇文档、词汇表大小是一百万那么特征矩阵在数学上就是 2 × 10 的六次方行、10 的六次方列稠密存储毫无希望但用稀疏矩阵存储 TF-IDF每篇文档只有几百上千个非零项整体占内存极小。在 sklearn 里TfidfVectorizer输出的本身就是scipy.sparse的 CSR 格式from sklearn.feature_extraction.text import TfidfVectorizer texts [这是第一句话, 这是第二句话包含一些词, 另一个完全不同的小文档] vec TfidfVectorizer() X vec.fit_transform(texts) print(type(X)) # class scipy.sparse._csr.csr_matrix这时候千万不要X.toarray()去转稠密除非你已经确认数据量很小。常规做法就是把X直接传给模型比如 LR、SVM、XGBoost 都能接入稀疏矩阵。训练时模型内部也会针对稀疏结构做专门的优化速度比稠密输入更快。遇到需要做特征预处理的情况比如归一化直接使用sklearn.preprocessing里的稀疏友好类即可尽量别在中间步骤转回稠密。5. 性能陷阱排查这些坑我替你踩过了5.1 toarray() 是把内存炸弹很多新手有个直觉操作就是为了看得到矩阵长什么样直接调用toarray()。我特别理解这种需求但在大矩阵上这是绝对禁忌。一个 50 万 × 200 万的真实稀疏矩阵转换过程相当于一口气申请几 TB 内存结果必然是 OOM进程被杀掉。真要是想看某个局部数据只需要切片出小区域再转稠密small csr[:10, :10].toarray()原则就是只稠密化你真正需要的那一小块数据。5.2 隐式转换和格式错配稀疏矩阵之间做运算时SciPy 会自动检查格式兼容性。如果格式不一致它会在底层进行转换这种隐式转换单次看可能只要几十毫秒但在循环里被反复触发的时候性能损耗会呈指数级放大。我写过一段入门的笨拙代码在循环里反复对同一个 CSR 矩阵取列再和另一个 CSC 矩阵做计算结果每轮都在转格式几分钟的任务跑了 40 分钟。后来我在循环外面把两个矩阵统一成相同格式时间直接降到 3 分钟。所以我的建议是在进入循环前就做一次格式抉择。关于indptr越界等问题也值得注意。如果你手动修改稀疏矩阵的内部数组破坏了indptr、indices、data的一致性行为会非常诡异尽量不要去直接操作内部结构。5.3 稀疏矩阵的广播机制是硬伤SciPy 稀疏矩阵不支持和广播语义匹配的普通逐元素运算。举个例子如果你试图把一个(3, 1)的数组减去一个(3, 4)的稀疏矩阵结果大概率会出现错误或者返回一个你根本不需要的稠密结果。广播是 NumPy 里非常自然的习惯但稀疏矩阵并不遵循。我的经验是先做一次彻底的形状分析之后用矩阵乘法和显式构造来替代广播场景。比如对稀疏矩阵按列标准化我不会去对每列做广播除法而是先计算列和向量再左乘一个对角稀疏矩阵或者使用更底层的操作。5.4 稀疏矩阵求逆几乎总是错误选择我知道有同学新学线性代数的时候习惯性遇到Ax b就写x A^-1 b。但在大规模稀疏矩阵场景下千万别直接A.inv()。稀疏矩阵的逆往往可能是稠密矩阵这意味着你可能求出了一个几百 GB 的稠密结果而且求逆的时间复杂度也远高于直接求解。正确做法是上一节提到的spsolve或迭代法。只有在极小矩阵、对运算速度要求不高的时候求逆才是一个可用选项。下面把这个环节里的常见问题整理成速查表方便你对照排查问题可能原因解决方案内存 OOM某一步触发了toarray()或隐式稠密化全链路保持稀疏格式切片看数据矩阵运算结果错误构建时存在重复位置没有做聚合对三元组先去重汇总再构建 COO 转 CSR切片特别慢CSR 做列切片先tocsc()再做列切片循环运算慢每轮都在做格式转换循环前统一格式全局只转一次eigs返回复数矩阵是实数对称但用了eigs换成eigsh稀疏矩阵和数组广播出错稀疏矩阵不支持 NumPy 广播语义用矩阵乘法或显式构造替代6. 项目实践中的性能调优细节6.1 数据加载阶段尽量批量构建我见过很多同学在构建稀疏矩阵时喜欢在循环里逐个mat[i, j] val。如果对象是 LIL 或 DOK勉强能跑但数据量一大依然很慢。更推荐的方式是先把源数据解析成三个 NumPy 数组然后一次性调用coo_matrix。瓶颈在 IO 和解析上不在矩阵构建本身。如果把构建写成循环会额外引入 Python 层的解释开销数据量越大越明显。实际项目中我一般用生成器或者分块解析的方式处理原始日志每读完一个 batch 就把三元组追加到内存列表最后将整个列表拼成 NumPy 数组再统一构建 COO转 CS。如果源文件过大也可以先落盘到临时文件但“批量构建”的内核始终是统一的。6.2 留意数据类型和精度稀疏矩阵的dtype默认从数据推导。如果我们只插入整数矩阵的数据类型会是整数一旦后续做某些浮点运算容易触发类型转换影响性能。我建议在构建时显式指定dtypenp.float64尤其当后续要传给求解器和模型时。此外如果你确实只需要 0、1 的布尔型矩阵可以考虑dtypebool。这在某些图算法里能节省很多内存但要注意做算术运算时结果可能不再是布尔。6.3 对象引用与内存释放大型项目里中间对象特别多稍不注意就会同时持有多个大稀疏矩阵。sparse矩阵本身占用内存不小除了用del删除外最好配合gc.collect()。不过更推荐的做法是尽量把流程设计成单向管道同一个内存位置可以被后续对象复用。在 Jupyter Notebook 环境里调试大矩阵我把这些经验归结为四句话没有产生大对象的中间变量不要轻易保存副本输入的原始数据如果已经转换成稀疏矩阵原始数组可以尽早释放不同阶段的矩阵使用函数作用域隔离避免全局变量引用链拖拽生命周期长时间运行的服务里注意缓存别把一个早该释放的矩阵偷偷留在内存里。6.4 优先使用稀疏感知的算法库除了 SciPy 本身很多第三方库对稀疏矩阵有专门优化。特征工程领域一些库在内部会自动识别 scipy.sparse 输入图算法领域也有大量库直接基于 scipy.sparse 存储邻接矩阵。如果你在某个场景发现自己的稀疏计算很慢先不要怀疑 SciPy而是检查自己是不是在“稀疏的外壳”里做了“稠密的操作”比如不必要的toarray()、不必要的格式化运算。这些才是性能下降的主要来源。我自己写了几年稀疏矩阵相关代码最大的体会是稀疏矩阵不是银弹但它一定是你处理高维稀疏数据的第一选择。只要在构建阶段选对格式在全链路中保持稀疏严格控制稠密化的时机内存和性能的问题基本能控制得很好。最后再补充一个小技巧。当你第一次拿到一个大数据矩阵时不管它是不是稀疏先打印一下shape、nnz和dtype心里有个底。这批信息能帮你判断后面该走哪条路、会不会出现内存紧张。我在项目里习惯用一行代码快速检查print(csr.shape, csr.nnz, csr.dtype)nnz是非零元素个数这个数字趁早掌握很多设计决策其实都取决于它。