numpy：有效地读取一个大数组答案

【问题标题】：numpy: efficiently reading a large arraynumpy：有效地读取一个大数组
【发布时间】：2011-05-20 22:33:44
【问题描述】：

我有一个二进制文件，其中包含一个密集的n*m 32 位浮点矩阵。将它读入 Fortran 排序的 numpy 数组的最有效方法是什么？

文件大小为数 GB。我可以控制格式，但它必须紧凑（即长度约为 4*n*m 字节）并且必须易于从非 Python 代码生成。

编辑：该方法必须直接生成一个Fortran有序矩阵（由于数据的大小，我无法创建一个C有序矩阵然后对其进行转换到一个单独的 Fortran 有序副本中。）

【问题讨论】：

scipy.org/Cookbook/InputOutput 会回答您的问题吗？（参见“二进制文件”部分）
@nimrodm 谢谢。事实上，我已经尝试过那里描述的一些方法。我问这个问题是希望有人能挺身而出，他要么有做我想做的事情的第一手经验，要么熟悉numpy 内部结构并可以从这个角度提出建议。
通常我在将非常大的数组读入 numpy 时发现，我需要提前知道大小，以便预先分配适当的数组来保存数据。事先知道尺寸吗？如果没有，请尝试使用两遍方法：首先扫描以发现数据的大小/维度，然后分配数组，然后读取/解析该数组。
@Peter 好点，谢谢。我事先知道大小（我控制数据格式，所以我可以将大小写成文件头的一部分。）

标签： python performance numpy scipy large-files

【解决方案1】：

NumPy 提供fromfile() 来读取二进制数据。

a = numpy.fromfile("filename", dtype=numpy.float32)

将创建一个包含您的数据的一维数组。要将其作为二维 Fortran 有序 n x m 矩阵访问，您可以对其进行整形：

a = a.reshape((n, m), order="FORTRAN")

[编辑：reshape() 在这种情况下实际上复制了数据（请参阅 cmets）。要在不复制的情况下执行此操作，请使用

a = a.reshape((m, n)).T

感谢 Joe Kingtion 指出这一点。]

但老实说，如果您的矩阵有几 GB，我会选择像 h5py 或 PyTables 这样的 HDF5 工具。这两个工具都有常见问题条目，将工具与另一个工具进行比较。我一般更喜欢 h5py，虽然 PyTables 似乎更常用（而且两个项目的范围略有不同）。

HDF5 文件可以用数据分析中使用的大多数编程语言编写。链接的维基百科文章中的接口列表不完整，例如还有一个R interface。但我其实不知道你想用哪种语言来写数据……

【讨论】：

@Sven 谢谢。与简单的n*m 数组相比，使用 HDF5 存储一个大的密集浮点矩阵有什么优势？
@aix：仅举几例：在给定时间的内存中透明地仅保留部分矩阵的能力，透明压缩和处理矩阵维度的元数据的能力，消除错误的重要来源。而且由于它真的很容易使用，我没有看到真正的缺点。
@Sven 你能否在你的例子中澄清一些事情。如果我在 4GB 文件上执行 a = numpy.fromfile("filename", dtype=numpy.float32) 后跟 a = a.reshape((n, m), order="FORTRAN")，这是否可能会在内存中创建一个 4GB“C”矩阵，然后立即制作另一个 4GB 内存副本以将其转换为 Fortran 格式？跨度>
@aix - Rehshaping 不会复制数组，它只是返回它的新视图。在这种情况下没有内存重复。
@aix - 我假设您正在从磁盘读取 Fortran 有序数组。你得到一个平面数组，它实际上是一个 ixj Fortran 有序数组。这与 jxi C 有序数组相同，只是它是转置的。因此，我们将其重塑为 jxi，然后将其转置为 ixj。 Numpy 认为是 C 排序的 jxi 被视为 ixj，但这相当于直接将 Fortran 排序的 ixj 存储在内存中。

【解决方案2】：

Numpy 基本上将数组存储为平面向量。多个维度只是 Numpy 迭代器使用的不同视图和跨步创建的错觉。

有关 Numpy 内部工作原理的全面但易于理解的说明，请参阅优秀的 chapter 19 on The Beatiful Code book。

至少 Numpy array() 和 reshape() 有 C ('C')、Fortran ('F') 或保留顺序 ('A') 的参数。另见问题How to force numpy array order to fortran style?

使用默认 C 索引的示例 (row-major order)：

>>> a = np.arange(12).reshape(3,4) # <- C order by default
>>> a
array([[ 0,  1,  2,  3],
       [ 4,  5,  6,  7],
       [ 8,  9, 10, 11]])
>>> a[1]
array([4, 5, 6, 7])

>>> a.strides
(32, 8)

使用 Fortran 顺序进行索引 (column-major order)：

>>> a = np.arange(12).reshape(3,4, order='F')
>>> a
array([[ 0,  3,  6,  9],
       [ 1,  4,  7, 10],
       [ 2,  5,  8, 11]])
>>> a[1]
array([ 1,  4,  7, 10])

>>> a.strides
(8, 24)

另一种观点

此外，您始终可以使用数组的参数 T 获得另一种视图：

>>> a = np.arange(12).reshape(3,4, order='C')
>>> a.T
array([[ 0,  4,  8],
       [ 1,  5,  9],
       [ 2,  6, 10],
       [ 3,  7, 11]])

>>> a = np.arange(12).reshape(3,4, order='F')
>>> a.T
array([[ 0,  1,  2],
       [ 3,  4,  5],
       [ 6,  7,  8],
       [ 9, 10, 11]])

您也可以手动设置步幅：

>>> a = np.arange(12).reshape(3,4, order='C')
>>> a
array([[ 0,  1,  2,  3],
       [ 4,  5,  6,  7],
       [ 8,  9, 10, 11]])
>>> a.strides
(32, 8)
>>> a.strides = (8, 24)
>>> a
array([[ 0,  3,  6,  9],
       [ 1,  4,  7, 10],
       [ 2,  5,  8, 11]])

【讨论】：