SC.NumPy 05 | 如何实现本地地球科学数据集与NumPy间交互？

Introduction

前面我们讲述了NumPy数组的概念、创建、计算、统计等基本内容，以及索引、切片与拼接等较为高阶的应用。并对应地提出了一些可能的地球科学中的应用场景。

但是其中就涉及到一个问题，我们所有用到的案例数据都是使用NumPy随机生成，或自定义了一些特殊数组。而我们的数据都是保留在本地的，如何打通本地储存与运行内存间的这一环就很重要。(说人话：就是怎么把本地文件读取到NumPy中，或从NumPy中保存到本地？)

因此，本文将探索NumPy自带的几种文件读取(Input)和输出(Output)，即文件IO。

❝

需要注意，在更广阔的应用场景中，我们的数据更多地是以矢量、栅格的形式存在的，而非本文提到的这几种数据类型(可能仅有.txt文件会是较为常见的存储类型，在一些气象站资料存储中较为常见)。

这需要用到一些更专门的库来实现它们的IO，相关内容将在后续不断补充、扩展。但这几种类型作为临时的数据存储方式，倒也还是不错的。

但其数组最终读取出来后的数据组织形式，以及处理方式终究还是大多基于NumPy和我们接下来将要开启的Pandas进行的(Pandas的更新内容中，我们将对常见的.csv和Excel表格等格式的文件处理详细展开)。

因此，理解NumPy数组的原理和运行机制对后续至关重要。在可预见的未来，我们所有地球科学类型数据的处理都将围绕NumPy展开。

二进制文件IO

NumPy中可以使用save和load函数将数组保存为.npy格式的二进制文件。

import numpy as np

arr = np.array([[1, 2], [3, 4]])

np.save('my_array.npy', arr)
loaded_arr = np.load('my_array.npy')

print(loaded_arr)

[[1 2]
 [3 4]]

文本文件IO

虽然二进制文件具有更快的读取速度，但有时候我们希望检查具体输出的值。而.npy格式文件不被Windows所原生支持，所以这里就要提到我们最常打交道的格式了——文本文件(.txt)。

它也可以使用NumPy简单地实现读写：

arr = np.array([[1, 2, 3], [4, 5, 6], [7, 8, 9]])

np.savetxt('my_file.txt', arr)
loaded_arr = np.loadtxt('my_file.txt')

print(loaded_arr)

[[1. 2. 3.]
 [4. 5. 6.]
 [7. 8. 9.]]

存储文件时，除了指定文件名之外，还有更多可供自定义的参数。这里主要介绍我们最常用到的分隔符和数据类型。

arr = np.random.rand(3, 4)
print(arr, '\n')

np.savetxt('random_array.txt', arr, fmt='%.2f', delimiter='\t')         # 使用tab作为分隔符
loaded_arr = np.loadtxt('random_array.txt', delimiter='\t')

print(loaded_arr)

[[0.94374423 0.0944434  0.85795381 0.75437499]
 [0.80656512 0.60834451 0.79539787 0.61392404]
 [0.92506516 0.0246624  0.55445595 0.19859673]] 

[[0.94 0.09 0.86 0.75]
 [0.81 0.61 0.8  0.61]
 [0.93 0.02 0.55 0.2 ]]

可以看到，我们存储时由于指定了fmt参数(format)，原来的float64类型的数组被储存为保留两位小数(这里用到了%.2f来格式化数据，其含义是保留两位小数的浮点型。所涉及的模块内容我们将在后续系统讲解)。

这对于我们储存空间有限，而压缩小数位数并不会过多地影响计算精度时极为有效，能极大地节省存储空间。而指定分隔符更多地是便于阅读，或用于其他指定了数据组织形式的应用中。

如果我们使用了与指定分隔符不一致的方式打开文件，数据读取将失败：

arr = np.random.rand(3, 4)

np.savetxt('random_array.txt', arr, fmt='%.2f', delimiter=',')         # 使用tab作为分隔符
loaded_arr = np.loadtxt('random_array.txt', delimiter='\t')

{
 "name": "ValueError",
 "message": "could not convert string '0.25,0.85,0.28,0.77' to float64 at row 0, column 1.",
 "stack": "---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
Cell In[22], line 4
      1 arr = np.random.rand(3, 4)
      3 np.savetxt('random_array.txt', arr, fmt='%.2f', delimiter=',')         # 使用tab作为分隔符
----> 4 loaded_arr = np.loadtxt('random_array.txt', delimiter='\\t')

File c:\\Users\\RitasCake\\.conda\\envs\\eb-py310\\lib\\site-packages\
umpy\\lib\
pyio.py:1318, in loadtxt(fname, dtype, comments, delimiter, converters, skiprows, usecols, unpack, ndmin, encoding, max_rows, quotechar, like)
   1315 if isinstance(delimiter, bytes):
   1316     delimiter = delimiter.decode('latin1')
-> 1318 arr = _read(fname, dtype=dtype, comment=comment, delimiter=delimiter,
   1319             converters=converters, skiplines=skiprows, usecols=usecols,
   1320             unpack=unpack, ndmin=ndmin, encoding=encoding,
   1321             max_rows=max_rows, quote=quotechar)
   1323 return arr

File c:\\Users\\RitasCake\\.conda\\envs\\eb-py310\\lib\\site-packages\
umpy\\lib\
pyio.py:979, in _read(fname, delimiter, comment, quote, imaginary_unit, usecols, skiplines, max_rows, converters, ndmin, unpack, dtype, encoding)
    976     data = _preprocess_comments(data, comments, encoding)
    978 if read_dtype_via_object_chunks is None:
--> 979     arr = _load_from_filelike(
    980         data, delimiter=delimiter, comment=comment, quote=quote,
    981         imaginary_unit=imaginary_unit,
    982         usecols=usecols, skiplines=skiplines, max_rows=max_rows,
    983         converters=converters, dtype=dtype,
    984         encoding=encoding, filelike=filelike,
    985         byte_converters=byte_converters)
    987 else:
    988     # This branch reads the file into chunks of object arrays and then
    989     # casts them to the desired actual dtype.  This ensures correct
    990     # string-length and datetime-unit discovery (like `arr.astype()`).
    991     # Due to chunking, certain error reports are less clear, currently.
    992     if filelike:

ValueError: could not convert string '0.25,0.85,0.28,0.77' to float64 at row 0, column 1."
}

因此，使用NumPy读取第三方数据时，需要首先明确该数据源的存储形式。

上面我们讲到了存储时能存储为特定的数据类型，可节省存储空间。同样，当我们面对一个已有的大型数据集时，它本身已具有特定的数据类型。而我们受制于自己的硬件条件，可能需要压缩读取。我们也可以在读取时直接指定将数据读取成什么类型：

arr = np.random.rand(3, 4) * 100
print(arr, '\n')

np.savetxt('random_array.txt', arr, delimiter='\t')

loaded_arr = np.loadtxt('random_array.txt', delimiter='\t')
print(loaded_arr, '\n')

loaded_arr = np.loadtxt('random_array.txt', delimiter='\t', dtype=np.float16)
print(loaded_arr)

[[72.52279016 27.60526614 73.44863686 21.84277912]
 [32.77975841 96.19431062  7.34798216  5.90261337]
 [44.2292748  34.33757862 81.68238205 43.03481254]] 

[[72.52279016 27.60526614 73.44863686 21.84277912]
 [32.77975841 96.19431062  7.34798216  5.90261337]
 [44.2292748  34.33757862 81.68238205 43.03481254]] 

[[72.5   27.61  73.44  21.84 ]
 [32.78  96.2    7.348  5.902]
 [44.22  34.34  81.7   43.03 ]]

可以看到，我们原封不动地储存了一个数组。而通过指定数据类型dtype决定了最终读取出的数据精度。

而另一个更常用的操作则是，有些时候前几行内容是无效的，或者我们只需要某几列的数据。此时，我们可以进一步对需要特殊处理的位置使用参数限定：

arr = np.random.rand(5, 5) * 100
print(arr, '\n')

np.savetxt('random_array.txt', arr, delimiter='\t')         # 使用tab作为分隔符

loaded_arr = np.loadtxt('random_array.txt', delimiter='\t', skiprows=2, usecols=[0, 2, 4])
print(loaded_arr, '\n')

[[87.38476251 55.83061031 89.12212303 12.96685324 82.89666845]
 [79.13914198 91.36254096 47.18566036 68.06809675 58.00357039]
 [70.04846718 42.88226646 44.15204615 92.01741989 89.78729344]
 [17.31092689 70.17482622 25.05784428 34.69986909 31.42033019]
 [ 2.0699918  27.87860504 60.35116244 99.816091   70.98181825]] 

[[70.04846718 44.15204615 89.78729344]
 [17.31092689 25.05784428 31.42033019]
 [ 2.0699918  60.35116244 70.98181825]] 

这里我们使用skiprows参数跳过了前两行，并使用usecols参数指定了只读取[0, 2, 4]三列的内容。

这里你应该发现了，列数也是从0开始数的，尽管它还不是数组；而skiprows则是正常的自然计数方式。因此，Python中很多设定我们都需要注意，可能一些细节纰漏就会导致很严重的错误。

更自由地数组IO——npz

上面无论是看不见的二进制文件也好，看得见的文本文件也好，每次我们终究只能存储一个数组。

这也就意味着，假设我们有一个区域的长期气温，我们可以把它转换为(t, x, y)的三维数组进行储存。但是如果我们有多个区域，而且可能它们空间范围不一样，因而数组维度长度也不一样，把它们压缩到一个数组中进行储存就会很麻烦。

当然也可以把它们储存为多个不同的文件，但同类型的数据，我们只需要管理一个文件很多情况下还是比要管理一堆文件要来的方便。

那么有没有一种，不管数组形状，我们只需要记住它的名字，一股脑打包压缩就行了。然后，读取时将特定名字的数据提取出来就行的方式？

那就不得不提NumPy自带的.npz格式了，它其实和MATLAB的.mat与R语言的.rds等有点像。但是因为NumPy只专注于数组本身，内容和形式没有那么丰富(也许可以考虑融合Pandas的结构)。

arr0 = np.array([1, 2, 3, 4, 5])
arr1 = np.array([[1, 2, 3], [4, 5, 6], [7, 8, 9]])
arr2 = np.random.rand(3, 3, 3)

np.savez('my_array.npz', a=arr0, b=arr1, c=arr2)

data = np.load('my_array.npz')
data0 = data['a']
data1 = data['b']

data2 = np.load('my_array.npz')['c']                # 当然也可以直接一步到位

print(data0, '\n')
print(data1, '\n')
print(data2)

[1 2 3 4 5] 

[[1 2 3]
 [4 5 6]
 [7 8 9]] 

[[[0.90735062 0.0610163  0.57842624]
  [0.15831279 0.3820223  0.28797242]
  [0.1336372  0.18660548 0.67522842]]

 [[0.10486749 0.7858729  0.41940512]
  [0.86897401 0.60569692 0.99331381]
  [0.62838982 0.9566743  0.90177265]]

 [[0.3436089  0.2955513  0.99697909]
  [0.02338764 0.09861164 0.58682854]
  [0.26079274 0.03435842 0.33366696]]]

以上就是.npz的简单用法，我们把三个不同类型的数组存入了同一.npz文件中，并给各自分配了一个别名。在读取时，我们只需要根据对应名字读取出需要的内容即可。

Supplementary

之前的内容中，感觉还有一个比较关键的概念没有提及，这里进行补充。其他可能存在的一些遗漏，也会在后续更新中以这种形式进行增补。

这里主要想提到的概念就是浅拷贝和深拷贝的问题。

❝

浅拷贝与深拷贝

• 浅拷贝

浅拷贝是指在拷贝对象时，只拷贝对象的引用，而不拷贝对象本身。也就是说，拷贝的对象和原对象共享同一块内存空间，修改其中一个对象，另一个对象也会受到影响。

• 深拷贝

深拷贝会创建一个新的数组，并复制原始数组的所有数据到新的数组中。这意味着原始数组和拷贝数组拥有不同的内存空间。对其中一个数组进行修改，不会影响另一个数组。

所以，如果我们需要对一个数据进行不同的两种计算，而我们两个数组是使用浅拷贝连接的话。两个数据的计算结果将是一致的，可能导致非常巨大的错误，需要特别注意。

我们首先给出一个浅拷贝的例子：

a = np.array([1, 2, 3])
b = a

print(a)
print(b)

a[0] = 10

print(a)
print(b)

[1 2 3]
[1 2 3]
[10  2  3]
[10  2  3]

可以看到，我们使用赋值符号=直接建立起了数组a和b间的联系，但是这是一种共享内存地址的浅层连接，不会重复占用内存。意味着，其中一个数值发生变化，另一个也将改变。

a = np.array([1, 2, 3])
b = a.copy()

print(a)
print(b)

a[0] = 10
b[1] = 20

print(a)
print(b)

[1 2 3]
[1 2 3]
[10  2  3]
[ 1 20  3]

而如果我们使用copy将数组完全复制一遍到另一块内存中，它们的存储相互独立，后续的操作也就各行其是，互不影响了。

后记

关于NumPy的主体内容在这里就告一段落了，后面应该会有涉及到它的诸多案例和应用。NumPy所支撑的数组确实是科学计算最核心的东西，需要深入地理解才便于我们后续更高阶地建模和分析。

不得不说，从0开始建立起一个知识体系是很艰难的，真正动起笔来远比想象中困难和复杂很多倍。但是写了几篇之后也逐渐开始轻车熟路，没有之前那么大的压力了。

也很感谢大家的支持，一个月不到的时间就突破了一百粉，比我预期的快了很多倍。那么，向着五百粉、千粉的目标，还需要持续地输出更多内容。

后续我们将开始更新Pandas，这个工具我觉得不仅是科学计算本身，就算是用于简化我们平时的Excel拖表也能省去相当多的时间和麻烦，是涉足Python不得不一试的强大工具。

但我觉得在提出了Pandas基本的概念后(数据类型DataFrame和Series)，可能需要花上两、三期对Python基础语法快速补完(包括控制流、函数等核心内容)，以便于我们后续更复杂多样的代码实现。

以上就是这期的全部内容和对于接下来内容的设想了。那么，我们下期见。

Manuscript: RitasCake

Proof: Philero; RitasCake

获取更多资讯，欢迎订阅微信公众号：Westerlies