Xarrayで出力するnetCDF-4/HDF5のチャンク・フィルタ(圧縮)・データパッキング

## はじめに

本記事では、Xarray を用いて出力した netCDF‑4 ファイルのチャンクや圧縮に関する情報をまとめる。Xarrayの使い方だけでなく、その基礎となるnetCDF‑4やHDF5の仕様についても解説する。

:::note warn
この記事では主に用語解説を ChatGPT で生成し、手動で修正・確認を行っている。コード断片はすべて実行し、結果を確認済みである。
:::

:::note
今後、Bit-groomingなどの量子化、bloscの詳細やベンチマーク結果を追記する予定
:::

## 背景知識

この節では、netCDFの圧縮に関する基礎となるnetCDFやHDF5の仕様について説明する。Xarrayの使い方を知りたい方は読み飛ばして、必要に応じて拾い読みするのでも問題ない。

### HDF5のチャンクとフィルタ

netCDF-4を理解するにあたり、内部形式として利用しているHDF5の機能のチャンクとフィルタについて理解する必要がある。

**チャンク(chunk)** は、変数データを等サイズの小さな多次元ブロックに分割してディスクへ保存する仕組みである。必要なチャンクだけを読み書きできるため部分アクセスが高速化される。

**フィルタ(filter)** とはデータセットをチャンク単位で読み書きするときに適用される入出力パイプラインで、圧縮・変換・誤り検出などを自動化する。書き込み時にフィルタを順に適用し、読み込み時には逆順で解除するため、ユーザーは通常の API を介して透過的に利用できる。フィルタには、HDF5 本体にあらかじめ組み込まれた **定義済みフィルタ（predefined filters）** に加えて、**サードパーティフィルター（third-party filters)** をプラグインとして組み込むことができる。

定義済みフィルタは以下の6種類である。

|フィルタ|ID|説明|
|---|---:|---|
|**Deflate (zlib)**|1|zlib を用いたdeflate圧縮 |
|**Shuffle**|2|バイト単位のシャッフル前処理|
|**Fletcher32**|3|32 bit のチェックサムで誤り検出|
|**Szip**|4|NASA SZIP 圧縮|
|N-Bit|5|任意のデータ型に対し、無効領域を除去してビット幅を最適化|
|Scale-Offset|6|浮動小数点／整数データのビット幅を縮小する非可逆圧縮|

このうちnetCDF-4で用いられるのはShuffle, Fletcher32, Deflate, Szipの4つである。圧縮フィルタのうち、Szipは特許制限があることからオプションであり、Deflateが主に使われる。Scale-OffsetおよびN-Bitフィルタは利用されず、非可逆圧縮を行う場合はnetCDF-4の量子化やCF規約のパッキングを用いる。

サードパーティーの圧縮プラグインは[HDF5 Compression Plugins](https://github.com/HDFGroup/hdf5_plugins)にある。これらのプラグインは、活発にメンテナンスされるMaintainedとそうでないCommunityの2つのカテゴリに分けられる。netCDF-4では以下に示すMaintainedのうちBzip2, Blosc, Zstandardフィルタが利用される。

|フィルタ|ID|説明|
|---|---:|---|
|**Bzip2**|307|Burrows–Wheeler + Huffman ベースの高圧縮フィルタ|
|LZF|32000| DEFLATE 風の高速圧縮フィルタ
|**Blosc**|32001|マルチスレッド／SIMD 対応メタコンプレッサ|
|LZ4|32004| 非常に高速なロスレス圧縮。|
|JPEG-XR|32007|JPEG-XR 画像圧縮フィルタ|
|Bitshuffle (BSHUF)|32008|ビットレベルのシャッフルで圧縮率向上|
|**Zstandard (Zstd)**|32015|Facebook製圧縮アルゴリズム|
|ZFP|32013|浮動小数点配列向け高スループット可逆／非可逆圧縮|
|Blosc2|32026|64-bit super-chunk 対応の次世代 Blosc。Blosc1 を拡張|

### netCDF-4の上位仕様: NASA ESDS-RFC-022

netCDF-4/HDF5の仕様は2011年に[NASA ESDS-RFC-022](https://www.earthdata.nasa.gov/about/esdis/esco/standards-practices/netcdf-4hdf5)として定められた。これはnetCDF I/O ライブラリに依存しない独立した標準仕様であり、netCDF-4 が HDF5 のサブセットであることを明文化している。

netCDF-4で利用されるフィルタに関して、当時存在した標準フィルタ（主に deflate と szip）に言及しているが、それ以外のフィルタについては言及はない。

### netCDF-4の下位仕様: UnidataのNetCDF Documentation

Unidata は NSFからの助成金で運営される UCAR のコミュニティプログラムで、netCDF ライブラリの標準を策定・維持している。

Unidataが提供するnetCDF-Cライブラリのドキュメントの[Appendix B. File Format Specifications](https://docs.unidata.ucar.edu/netcdf-c/current/file_format_specifications.html)はRFC-022の下位仕様にあたり、具体的なファイルバイナリの文法と実装上の注意点を網羅している。netCDF-C ライブラリのバージョンアップに合わせて随時改訂される。

netCDF-4の標準フィルタについては、[Appendix D. NetCDF-4 Filter Support](https://docs.unidata.ucar.edu/netcdf-c/current/filters.html)に書かれており、HDF5の定義済みフィルタのうちShuffle, Fletcher32, Deflate (zlib) , szipを継承している。それに加えて、サードパーティフィルタのうちZstandardおよびBzip2も標準フィルタに含まれる。Bloscについては、netCDF-4/HDF5の標準フィルタではないが、NCZarr向けのオプションとしてサポートされる。

validなnetCDFは標準フィルタのみを利用する必要がある。逆にそれ以外のHDF5のフィルタを使用した場合、netCDFとしてはinvalidであり、フィルタが正しく処理されずポータビリティや互換性の問題が発生する可能性がある。

### netCDF-Cライブラリが提供する圧縮フィルタAPI

netCDF-C ライブラリは Unidata が公式に提供するnetCDF データモデルとファイルフォーマットの参照実装である。

HDF5のフィルタを適用する汎用の`nc_def_var_filter`に加えて、各標準フィルタに対応した専用の(dedicated) APIを提供している。

[deflate](https://docs.unidata.ucar.edu/netcdf-c/current/group__variables.html)であれば
```C
int nc_def_var_deflate(int ncid, int varid,
                       int shuffle, int deflate, int deflate_level);
```
のように圧縮レベルに加えてshuffleフィルタの適用の有無を指定できる。shuffleフィルタを利用する場合、出力ファイルはHDFのshuffleフィルタ(ID 2)とdeflateフィルタ(ID 1)が適用されることになる。

一方、[Zstandard](https://docs.unidata.ucar.edu/netcdf-c/current/dfilter_8c.html)の場合
```C
int nc_def_var_zstandard(int ncid, int varid,
                         int level);
```
圧縮レベルのみでshuffleフィルタの適用オプションは存在しない。shuffleフィルタを使用したい場合は、汎用の`nc_def_var_filter`で指定する必要がある。これはBzip2 (`nc_def_var_bzip2`) についても同様である。

[Blosc](https://docs.unidata.ucar.edu/netcdf-c/current/dfilter_8c.html)は標準フィルタではないが、netCDF-4/HDF5の出力時にも利用可能である。
```C
int nc_def_var_blosc(int ncid, int varid,
                     unsigned subcompressor,
                     unsigned level,
                     unsigned blocksize,
                     unsigned addshuffle);
```
`addshuffle`で(byte)shuffleおよびbitshuffleを指定できるほか、`subcompressor`（Blosc コーデック）として以下を指定可能である。
| 定数名             | 値  | 説明                         |
|------------------|----:|----------------------------|
| BLOSC_LZ     |   0 | Blosc 独自の LZ ベース（FastLZ）コーデック  |
| BLOSC_LZ4    |   1 | LZ4（非常に高速）    |
| BLOSC_LZ4HC  |   2 | LZ4HC（LZ4 の高圧縮モード） |
| BLOSC_SNAPPY |   3 | Snappy（Google 製の高速圧縮) |
| BLOSC_ZLIB   |   4 | zlib（Deflate）   |
| BLOSC_ZSTD   |   5 | Zstandard（高圧縮) |

注意点として、HDF5のフィルタパイプラインではBlosc (ID 32001)単独が出現する。
addshuffleやsubcompressorでBitshuffleやLZ4を指定しても、HDF5のBitshuffle (ID 32008)やLZ4 (ID 32004)などは出現しない（逆にこれらのフィルタを使用したHDF5はnetCDF-4としてはinvalidである）。

### CFメタデータ規約

Climate and Forecast（CF）メタデータ規約は、気候や大気海洋のnetCDFデータセットを対象として、変数名・次元・属性などのメタデータ構造を定めるものである。

[CF 規約の8.1 章](https://cfconventions.org/cf-conventions/cf-conventions.html#packed-data)では、浮動小数点データを `int8` や `int16` に **整数量子化**して保存するパッキングについて定められている。具体的には読み込み時に

```python
unpacked = packed * scale_factor + add_offset
```

で元の値へ復元する。`scale_factor` と `add_offset` は対象変数の **属性 (attribute)** として記録される（ 参照）。


## Xarrayを用いたnetCDF出力

### 検証環境について

本記事の検証環境は、MacBook Air M3, メモリ24GB, macOS Sonoma 14.7上のcondaにより構築したPython環境である。
`$ conda list` から主要ライブラリのみ抜粋して示す。パッケージのバージョンアップに伴い、将来的に動作が変わる可能性がある。

```
python                    3.12.4          h30c5eda_0_cpython    conda-forge
libnetcdf                 4.9.2           nompi_he469be0_114    conda-forge
hdf5                      1.14.3          nompi_hec07895_105    conda-forge
zlib                      1.3.1                hfb2fe0b_1    conda-forge
zstd                      1.5.6                hb46c0d2_0    conda-forge
zstandard                 0.22.0          py312h721a963_1    conda-forge
netcdf4                   1.7.1           nompi_py312hec02768_101    conda-forge
h5py                      3.11.0          nompi_py312h903599c_102    conda-forge
h5netcdf                  1.6.1              pyhd8ed1ab_0    conda-forge
hdf5plugin                5.1.0           py312h5518634_0    conda-forge
xarray                    2024.6.0           pyhd8ed1ab_1    conda-forge
numpy                     1.26.4          py312h8442bc7_0    conda-forge
dask                      2024.6.2           pyhd8ed1ab_0    conda-forge
```

上記パッケージは、conda であれば `conda install <パッケージ名>` でインストールできる。

以下のインポート文は、以降省略する。
```Python
import xarray as xr
import h5py
```

### バックエンドエンジン
Xarray の `Dataset.to_netcdf()` および `xr.open_dataset()` では、netCDF-4/HDF5の出力に次の 2 種類のバックエンドエンジンを指定できる。

|エンジン名|概要|
|---|---|
|**`netcdf4`**|Unidata 製の C ライブラリ（netCDF‑C/HDF5）を Python から利用する _netCDF4‑Python_ バインディング|
|**`h5netcdf`**|純粋 Python 実装の h5py 経由で HDF5 ファイルを読み書きする軽量インターフェース|

エンジンは以下のように指定する（省略時は `netcdf4` が選択される）。


```python
# ds は xarray.Dataset
ds.to_netcdf("file.nc", engine="netcdf4")
ds.to_netcdf("file.nc", engine="h5netcdf")
```

どちらのエンジンでもフォーマットのデフォルトは `format="NETCDF4"` なので、本記事では明示しない。

### チャンクサイズの指定

チャンクサイズは読み込み速度だけでなく、書き込み性能や圧縮効率も左右しますが、とりわけ「部分読み込み」の高速化に大きく寄与する。圧縮フィルタはチャンク単位で適用されるため、欲しいデータ領域が一つ（あるいは少数）のチャンクに収まるようにレイアウトを設計すれば、そのチャンクだけを展開（解凍）するだけで済み、入出力量と処理時間を大幅に削減できる。

チャンクサイズはエンジンに関係なく `encoding` で指定できる。

```python
# ds は xarray.Dataset
# 2 次元変数 "var" を例とする
ds.to_netcdf(
    "file.nc",
    encoding={"var": {"chunksizes": (72, 144)}},
)
```

生成された netCDF‑4/HDF5 ファイルのチャンクは、コマンドラインなら

```shell
h5ls -v file.nc
```

Python では

```python
with h5py.File("file.nc", "r") as f:
    dset = f["var"]
    print(dset.chunks)
```

で確認できる。


#### デフォルトのチャンクサイズ

チャンクサイズを明示しない場合、挙動はエンジンによって異なる（配列形状 `(721, 1440)` の例）：

- **`netcdf4` エンジン**：`(721, 1440)`（全体が 1 チャンク）
    
- **`h5netcdf` エンジン**：`(46, 180)`（自動設定）
    
#### 混乱しやすいポイント：DaskとnetCDF-4/HDF5

Xarray でファイルを読み込むと、既定では `numpy.ndarray`（NumPy 配列）になる。すべての演算は即時評価で実行され、全データがメモリに読み込まれる。これに対し、Dask を使うと遅延評価・並列計算が可能な `dask.array.Array`（Dask 配列）をバックエンドにでき、Dask 側でも「チャンク」という分割単位を持つ。

```python
da_chunked = ds["var"].chunk({"latitude": 72, "longitude": 144})
print(da_chunked.chunks)
```

**Dask 配列のチャンクと netCDF‑4/HDF5 のチャンクは無関係** である。Dask 配列を `to_netcdf()` しても、ファイル側のチャンクは `encoding`（未指定ならデフォルト）に従う。逆に、既存ファイルのチャンクがどうであれ `xr.open_dataset()` で `chunks=` オプションを与えない限り、NumPy 配列として読み込まれる。

（ちなみに Zarr 形式では Dask のチャンクがそのまま書き出しチャンクに反映されるため、netCDF とは挙動が異なる。）

### フィルタ(圧縮)


定義済みフィルタのうち、Xarray の netCDF 出力で使える主なものはDeflate, Shuffle, Fletcher32の 3 種である。サードパーティーフィルタでは、Deflate より高い圧縮率や速度を得やすい **Zstandard (zstd)** が代表的で、こちらも Xarray から利用できる。

#### Deflate (zlib)

netCDF‑4/HDF5 で最も一般的なロスレス圧縮が Deflate (zlib) である。以下では _Deflate_・_zlib_・_gzip_ を同一視してかまわない[^1]。

両エンジンとも次の指定で Deflate を有効化できる。

```python
ds.to_netcdf(
    "file.nc",
    encoding={"var": {"zlib": True, "complevel": 4}},
)
```

`h5netcdf` エンジンでは h5py 互換の指定も可能。

```python
ds.to_netcdf(
    "file.nc",
    engine="h5netcdf",
    encoding={"var": {"compression": "gzip", "compression_opts": 4}},
)
```


圧縮レベル `complevel` は 1–9（高いほど高圧縮・低速）。未指定時のデフォルトは両エンジンとも **4**[^2]。参考までに `gzip` コマンドや Python 標準 `zlib` のデフォルトは 6[^3]である。

Deflate フィルタの有無は `h5ls -v` または
```Python
with h5py.File("file.nc", "r") as f:
    dset = f["var"]
    print(dset.compression,　dset.compression_opts)
```
で確認できる。

#### Shuffle

Shuffle フィルタは Deflate などと組み合わせ、圧縮率を向上させる前処理用フィルタである。

```Python
ds.to_netcdf("file.nc", encoding={
    "var": {"zlib": True, "complevel": 4, "shuffle": True}
})
```
`shuffle` のデフォルトはエンジンで異なり、`netcdf4` は **True**、`h5netcdf` は **False** なので注意。

Shuffleフィルタが使用されているかは、`h5ls -v`または
```Python
with h5py.File("file.nc", "r") as f:
    dset = f["var"]
    print(dset.shuffle)
```
で確認できる。

#### Zstandard (zstd)

Zstandard の指定方法はエンジンで異なる。

`netcdf4` エンジンの場合、

```Python
ds.to_netcdf("file.nc", engine="netcdf4", encoding={
    "var": {"compression": "zstd", "complevel": 4}
})
```

`complevel` は 1–22（デフォルト 4）。記事執筆時点では **`shuffle=True` を指定しても shuffle が有効にならない**。これは元のnetCDF-Cライブラリの`nc_def_var_zstandard`にshuffle引数がないことを引き継いでいるためだと考えられる。

`h5netcdf` エンジンの場合、

```Python
import hdf5plugin
ds.to_netcdf("file.nc", engine="h5netcdf", encoding={
    "var": {"compression": hdf5plugin.Zstd(3), "shuffle": True}
})
```
`hdf5plugin.Zstd()` の引数を省略するとレベル 3（Zstandard のデフォルト）になる。

Zstandard の有無は `h5ls -v` で確認するか、h5py の低レベル API でフィルタを列挙する。

```Python
with h5py.File("file.nc", "r") as f:
    dset = f["var"]

    # 低レベルで作成プロパティリストを取得
    plist = dset.id.get_create_plist()

    # 登録されているフィルター数を取得
    nfilters = plist.get_nfilters()
    print(f"Number of filters: {nfilters}")

    # 各フィルターの詳細を出力
    for idx in range(nfilters):
        filter_id, flags, cd_values, name = plist.get_filter(idx)
        print(f"Filter-{idx}:")
        print(f"  id        : {filter_id}")
        print(f"  name      : {name}")
        print(f"  flags     : {flags}")
        print(f"  parameters: {cd_values}")
```

### データパッキング

Xarray では `xr.open_dataset()` がデフォルトでアンパックを行うが、生の整数値を取得したい場合は `mask_and_scale=False` を指定する。

この方法だけで **`float32` と比べてサイズを int16 なら 1/2，int8 なら 1/4 に削減**できる。  
圧縮フィルタ（zlib／zstd など）を併用しなくても一定の節約効果が得られる点が利点である。

```python
ds.to_netcdf(
    "file.nc",
    encoding={
        "var": {
            "dtype": "int16",
            "scale_factor": 0.01,
            "add_offset": 273.15,
            "_FillValue": -32767,
        }
    },
)
```

NetCDF-C の既定では **`int16` の `_FillValue` は −32767**[^5]。特に理由がなければこの値を採用するとよい。

パッキングの有無は

```bash
ncdump -h file.nc
```

で `scale_factor` と `add_offset` 属性の有無を確認すれば判別できる。

## 圧縮設定の性能比較

本節では、配列形状 **(721 × 1440)** の単一変数（海面更正気圧）を題材に、次の 3 指標を計測・比較した。変数種別や配列サイズが変われば結果も変動し得るため、**目安として参照**してほしい。

1. **ファイルサイズ** — HDF5 の `get_storage_size()` を MB 換算  
2. **書き込み時間** — `to_netcdf()` を 10 回連続実行した合計  
3. **読み込み時間** — `xr.open_dataset().load()` を 10 回連続実行した合計  

すべての実験は同一の入力ファイル `input.nc` を基に、出力ファイル `output.nc` を上書きして行った。

---

### チャンク 1 つ／複数分割 × zlib 圧縮の効果

| 比較ケース | チャンクサイズ | 圧縮 |  
|------------|---------------|------|  
| A | (721, 1440) | なし |  
| B | (46, 180)  | zlib (level 4) |  
| C | (721, 1440) | zlib (level 4) |

共通設定：`netcdf4` エンジン／shuffle = True／パッキングなし  

![figure1.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/829116/d78e88e0-2e6b-4aee-9134-e7844963f7f4.png)

* **zlib 圧縮**でファイルサイズは約 **80 % 減**。  
* 本データでは **チャンクを分割しないほうが 10 % 程度小さく**なった。東西方向に値がなだらかな変数のためで、降水量のように局所性が大きい変数では逆転する可能性がある。  
* 圧縮に伴い書き込み／読み込み時間はおおむね **2 倍**になった。

---

### Deflate (zlib) vs. Zstandard のレベル比較

| 圧縮方式 | 評価レベル |
|----------|-----------|
| zlib | 1, 4, 9 |
| zstd | 1, 3, 10, 16, 22 |

共通設定：`h5netcdf` エンジン／チャンク (721, 1440)／shuffle = True／パッキングなし  

![figure2.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/829116/edc1afc6-1ac5-4552-93f7-8294888ba8da.png)

* **同等の圧縮速度**では *zstd* がわずかに高圧縮。  
* zlib は **level 4** が速度と圧縮率のバランス良好。  
* zstd はデフォルトの **level 3** ではやや圧縮不足。本データでは **level 10** がバランス良い。

![figure2-2.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/829116/8b2b6895-e29e-4780-917c-88f33a119d09.png)

- 読み込み速度は圧縮レベルにほぼ依存せず、**zstdはzlibの約半分の時間**

---

### Shuffle フィルタの有無

共通設定：`netcdf4` エンジン／チャンク (721, 1440)／zlib level 4／パッキングなし  

![figure3.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/829116/015df7ad-18ed-41ee-9685-adf8b9bcfaa7.png)

* shuffle = **True** でファイルサイズが **約 50 % 減**（zstd でも同傾向）。  
* 書き込み時間は **約 10 % 減**、読み込み時間はほぼ不変。

---

### 圧縮 × パッキングの組合せ効果

| ケース | パッキング | 圧縮 | shuffle |
|--------|-----------|------|---------|
| A | なし（float32） | なし | — |
| B | int16 (`scale_factor=1`, `add_offset=100000`) | なし | — |
| C | なし | zlib 4 | True |
| D | int16 | zlib 4 | True |

共通設定：`netcdf4` エンジン／チャンク (721, 1440)  

![figure4.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/829116/9dbae738-8054-4c91-8293-50bcfcf68d9e.png)

* **パッキングのみ**（A→B）でサイズは **2 分の 1**。  
* **圧縮あり** ではパッキング効果は小さい
    * ただし降水量のような変数で0.01mmなどの小さな値がパッキングにより0に丸められるのであれば、**圧縮＋パッキング** の相乗効果が期待できる。
---

### バックエンドエンジンの速度比較

共通設定：チャンク (721, 1440)／zlib level 4／shuffle = True／パッキングなし  

![figure5.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/829116/9e49e7a4-6c69-491a-98ba-64174faae54a.png)

* **書き込み** — `netcdf4`（C バインディング）が約 **20 %** 高速。  
* **読み込み** — 両エンジンで差はほぼなし。

---

### 結論と推奨設定

| 目的 | 推奨設定 | 備考 |
|------|----------|------|
| **バランス重視** | `netcdf4` エンジン / チャンク 1 つ / zlib(level 4) / shuffle = True | デフォルト値をほぼそのまま利用 |
| **圧縮率最優先** | `h5netcdf` エンジン / zstd(level 10〜16) / shuffle = True | 書き込みは遅くなるがサイズ最小 |
| **さらに削減したい** | 上記＋ **int16 パッキング** | 許容誤差が明確な変数に限る |

これらを基準に、**データ特性と許容時間**に応じて設定を調整すると良いだろう。

[^1]: **Deflate** は LZ77 とハフマン符号化を組み合わせたロスレス圧縮アルゴリズム（RFC 1951）。**zlib** は Deflate 流を RFC 1950 形式で包むストリーム／ライブラリ。**gzip** は Deflate 圧縮データに独自ヘッダーと CRC トレーラーを付加したファイル形式（RFC 1952）。Xarray の `compression="gzip"` も内部的に Deflate フィルタを呼び出すだけで、`.gz` コンテナを直接ラップするわけではない。
[^2]: `engine="netcdf4"`（netCDF4-Python ライブラリ）を通じて変数を作成する際、内部では `Dataset.createVariable(..., zlib=True)` が呼ばれ、圧縮レベル `complevel` のデフォルト値が4に設定される。https://unidata.github.io/netcdf4-python/ 
`engine="h5netcdf"` を使う場合も、`Group.createVariable(..., zlib=True)` の `complevel` デフォルトが4になる。https://h5netcdf.org/generated/h5netcdf.legacyapi.Group.html
[^3]: Linux の `gzip` コマンドでは、圧縮レベルの指定なしで実行した場合のデフォルトが6である。https://refspecs.linuxbase.org/LSB_3.0.0/LSB-PDA/LSB-PDA/gzip.html
Python の `zlib.compress()` や `zlib.compressobj()` における `level` のデフォルト（`Z_DEFAULT_COMPRESSION`）も現在6である。 https://docs.python.org/3/library/zlib.html
[^5]: https://docs.unidata.ucar.edu/netcdf/documentation/4.8.0/netcdf_8h.html