多维数组库Xarray 在数据处理中的应用

曹永正 王立鹏 高 山 温连杰 王 翠

（1、山东省海洋生态环境与防灾减灾重点实验室，山东 青岛266033 2、国家海洋局北海预报中心，山东 青岛 266033）

1 概述

地球科学领域一直与“大数据”密切相关，高分辨率的观测或模拟资料能够更准确地刻画不同尺度的物理过程，但也给数据分析处理带来了困难。近年来随着地球观测手段和物理模型的进步，人类能够获取的气象、海洋和陆地要素越来越多，其中再分析和模式预报结果多以正交网格形式呈现，而船舶、浮标和测站观测记录则是散点形式。为了更方便处理这些不同要素、不同形式的数据，Xarray 提供了一套较为统一的格式和操作工具。

Xarray 是一个开源项目，面向需要N 维数组，尤其是需要数据分析的科研工作，旨在让多维数据有关的工作更加简单、高效。Xarray 在类似NumPy的原始数组之上引入了维度、坐标和属性形式的标签，气象、海洋数据的地理空间和时间、单位、变量名称、注释等信息在处理过程中都得以保留。Xarray 通过数据引擎的方式，对netCDF、GRIB 和ZARR 等研究中常用格式都提供了相应支持。

Xarray 最早开发与2014 年，期间经过数十次更新迭代，功能已经日趋成熟和强大，目前最新版本为v0.19.0，本文内容基于该版本。

2 数据的查看与获取

针对不同需求，Xarray 提供了两种数据类型：

DataArray：带有标签或名称的多维数组。DataArray的基础是一个多维数组矩阵，在此基础上添加了用以描述该数组的元数据信息，既维度名称、坐标和属性等。

Dataset：可以看做是多个DataArray的集合，这些DataArray具有相同的维度，是彼此“对齐”的，他们有各自的元数据信息，但也可以有公用的部分。大多数可以在单个DataArray 上执行的操作都可以在Dataset 上执行。

对于只包含单个要素的数据文件，构建DataArray 类型对象即可，但实际工作中为了方便分发，一个数据文件中往往包含相同时空范围的多个要素信息，因此Dataset 类型更常使用。

不管是DataArray 还是Dataset，其中的数组都以Variable 对象的形式呈现，各种运算可以通过维度名称实现数组广播。

我们以包含2020 年12 月一个月内10m 风速u、v 分量和海平面气压3 个要素的ERA5 再分析数据（文件名202012_atm.nc）为例，首先在Python 环境下导入Xarray 库和Python 环境下的高效矩阵计算库numpy、数据分析库pandas：

可以看到，Xarray 以202012_atm.nc 数据文件为基础，创建了一个Dataset 对象。对象有3 个维度（Dimensions），分别为纬度（latitude）、经度（longitude）和时间（time），各维度方向的格点数为721、1440 和744。

所有数值数组均为Variable 对象格式，其中各维度数值归类到“坐标（Coordinates）”中，而变量10 米经向风矢量（u10）、10米纬向风矢量（v10）和平均海平面气压（msl）则归属于“数据变量（Data variables）”。其中纬度、经度和时间为1 维变量，即沿自身方向延伸的1 维数组，而u10、v10 和msl 为3 维变量，第0 维为时间，第1 维为纬度，第2 维为经度。除时间外，各项数值以32bits 浮点数格式存储，而时间则转换为专有的datetime64 格式以方便计算。

本数据集中含有2 项共有属性（Attributes），Conventions 和history，描述了数据来源和生成操作历史。要想进一步查看每个变量的属性信息，可以在环境中输入“file1.msl”，显示如下：

属性中即显示变量msl的单位为Pa，全名为“Mean sea level pressure”，标准名称为“air_pressure_at_mean_sea_level”。

为了进行复杂的分析计算，我们可以把数据从Dataset 或DataArray 中提取为numpy 支持的无属性数组格式：

msl_data=file1.msl.data

提取出的数组不再含有坐标和属性信息，为单纯的3 维格点气压值。

3 子集选取与合并

知道数据集基本结构后，我们就可以根据需要提取相应子集。对于DataArray 对象，Xarray 提供的loc()方法可以让用户像对待numpy 中Array 对象一样进行数据切片操作。例如从全球范围分辨率0.25°的原始气压数据中选取中国渤、黄海范围1°分辨率气压场，通过以下命令就能得到：

target_field=file1.msl.loc[:,40:30:4,117:127:4]

但该方法无法用于Dataset 对象。对于Dataset，需要借助sel()或isel()方法。其中sel()方法的输入参数为坐标值或坐标值列表，isel()方法的输入参数为坐标位置索引。两种方法用法相似，但sel()对大多数使用者来说更直观。例如选取u10、v10 和msl全部3 个变量的渤、黄海范围1°分辨率气压场，可以：

除了提取子集，在实际工作中也经常会遇到需要将两个或多个较小数据集合并为一个大数据集，Xarray 中可以用concat()和merge()方法实现。这里增加2021 年1 月数据集：

file2=xr.open_dataset('202101_atm.nc')

concat()方法是维度上的合并，需要指定操作维度，需要合并的数据集以列表方式传入：

target_file=xr.concat([file1,file2],dim="time")

而merge()方法是变量上的合并，需要合并的数据集同样以列表方式传入：

target_file=xr.merge([file1,file2])

merge()方法不需要指定操作维度，倘若file1与file2 坐标上存在差异，其维度坐标会自动补全，但变量场的缺失值会自动设置为缺省值Nan。

4 数据插值与计算

为了方便数据分析，Xarray 还提供了一些较基础的计算工具，包括时间、空间插值和简单的统计计算。

由于Xarray 将时间、空间维度都看做不同名称命名的坐标，因此各维度的插值均使用interp()方法实现。使用时，仅需将要调整的维度坐标值传入，未传入的坐标将保持原样。例如要插值得到50-55°N、120-125°E 范围内分辨率0.1°，2020 年12 月15 日08-14 时期间间隔0.5 小时的气压场，可以通过以下命令：latitude=lat_new,time=time_new)

图1 气压变化折线图（（a）为默认样式，（b）为添加颜色和标记样式参数）

Xarray 内置了几种基础的统计方法，有最大值（max()）、最小值（min()）、平均值（mean()）、中值（median()）、总和（sum()）和标准差（std()）等。但来自外部的数学运算多数也可以直接进行，例如四则运算、三角函数计算等。这些运算会直接作用在“Data variables”下的每个Variable 对象中，相当于对多维矩阵进行的计算。而最终结果会还原输入形式，即DataArray 或Dataset。

5 可视化

Matplotlib 是Python 环境下功能最全面、应用最广泛的可视化工具库，借助Matplotlib 强大的绘图功能，Xarray 可以对DataArray 和Dataset 数据进行多种形式的可视化。由于集成了Matplotlib 应用接口，使用Xarray 时不必先将数据导出到Matplotlib，而是直接调用plot()函数。这样不仅避免了繁琐的数据转移流程，还可以方便地将DataArray 或Dataset 中的属性和维度等信息应用在绘图中。

目前Xarray 支持的绘图样式涵盖了一维、二维和三维形式下多种常见的图表。plot()函数默认情况下为折线图，例如，查看2020 年12 月期间（120°E，36°N）处气压变化情况，只需直接使用plot()命令即可得到相应图像，通过添加绘图参数可进一步控制图表样式：

得到的图像中，Xarray 自动为横轴采用了时间格式标签，为纵轴标注了变量名称和单位，并将经纬度坐标作为图像标题。

图2 Xarray 支持的常用图表样式

除最基本的折线图外，Xarray 还可以绘制直方图、阶梯图、二维等值线图、矢量图、流线图等，示例如下：

直方图：file1.sel(longitude=120,latitude=36).msl.plot.hist()

阶梯图：file1.sel(longitude=120,latitude=36).msl[:10].plot.step()

填色等值线图：file1.sel(longitude=np.linspace(120,130,11),latitude=np.linspace(30,40,11),time='2020-12-01T12').msl.plot.contour()

网格填色图：file1.sel (longitude=np.linspace (120,130,11),latitude=np.linspace(30,40,11),time='2020-12-01T12').msl.plot.pcolormesh()

矢量图：file1.sel (longitude=np.linspace (120,130,11),latitude=np.linspace (30,40,11),time='2020-12-01T12').plot.quiver(x='longitude',y='latitude',u='u10',v='v10')

流线图：file1.sel (longitude=np.linspace (120,130,11),latitude=np.linspace (30,40,11),time='2020-12-01T12').plot.streamplot(x='longitude',y='latitude',u='u10',v='v10')

需要注意的是，DataArray 和Dataset 分别对应着不同的图表类型，即有的图表类型以DataArray 为输入，而另一些则以Dataset 为输入，这是由图表性质决定的。上述例子中，折线图、直方图、阶梯图、填色等值线图和网格填色图的输入为DataArray，矢量图和流线图的输入为Dataset。

6 数据的修改、创建与保存

对于读取为DataArray 或Dataset 格式的数据，可以使用Xarray 对其数值和维度、属性等信息进行修改。由于这些信息在DataArray 或Dataset 中以对象的形式呈现，因此修改实质上是对目标对象的操作。Xarray 为不同对象提供了不同操作方式，其中变量数值可以直接用赋值的方式修改，维度需要通过专门函数assign_coords()，属性则为“键-值对”的字典形式。例如，将file1 中的时间维由2020 年12 月改为1 月，变量msl的数值用随机数替代，并将属性中的“history”内容改为其它文本，可以通过以下几个步骤实现。首先借助pandas 生成新时间轴序列：

使用numpy 生成与变量msl 具有相同维度范围的随机数矩阵：

用生成的新数据为目标对象赋值：

最后修改属性中history 键对应值：

file1.attrs['history']='Edited by XXX at 2021-9-8'

在已有其它来源数据的情况下，还可以使用Xarray 从头创建一个DataArray 或Dataset，进而保存为netCDF 等文件格式。这里我们使用刚才生成的随机数矩阵,可以通过如下方式构建一个新DataArray 并以netCDF 格式保存：

其中da_new.nc 为此次保存的文件名。修改过的file1 也可以用同样的方法保存在本地。

7 结论

本文探索和介绍了Xarray 在地球物理数据处理中的常用功能和使用方法，对于科研中常见的资料格式，例如netCDF、GRIB 和ZARR 等，Xarray 提供了一套强大的工具和命令，除了本文中介绍的数据读取、查看、选取与合并、插值计算、可视化、修改保存等功能外，还可以借助Lamda 表达式实现自定义函数，以及与Numpy、Pandas、Dask、matplotlib 等工具结合使用，实现更复杂的数据分析、可视化，甚至并行计算。在数据量日益增大的今天，Python 环境下的数据科学生态系统成为了不可或缺的生产力来源，而Xarray 就为地球科学领域研究者提供了进入这套生态系统最方便的门户和工具。

Xarray 具有使用简单、集成度高、可拓展性强的优点，但也存在着如对复杂格式文件支持不完善、部分函数的计算效率低、与其它工具同时使用偶现Bug 等，目前仍在开发中。当前版本Xarray 功能已经覆盖了大部分数据操作需求，能够为数据研究分析工作带来很大方便。