geostatistics-interpolation

地球統計学・補間：知らないことを補間する。プロ向け製品と見た目のよいラスターの違いは、不確実性の面と誠実なクロスバリデーションにあります。ここではどちらも必須です。

手法の選択

用途	手法
滑らかな現象、迅速	IDW（べき乗パラメータを調整、1～3）
空間構造をもつ物理現象、不確実性が必要	通常クリギング（プロの既定の選択）
明確な勾配、共変量あり	普遍クリギングまたは回帰クリギング
共変量が豊富（DEM、土地被覆）	回帰クリギング（ランダムフォレスト残差クリギング）
カテゴリカル	指示クリギング
ハニカム状のテッセレーション、不確実性が必要	最近傍 IDW は妥当ですが、制約がありません。IDW と同程度の見た目になる製品には注意してください。

探索的（補間の前に）

• 点を値とともにプロットする。重複を確認する（重複はクリギング行列を特異にします）。
• ヒストグラム：強く歪んだ変数（降雨、濃度）は通常対数変換が必要です。予測値は適切に逆変換してください（バイアス補正付き対数正規クリギング）。
• トレンドを確認：x と y に対する値のプロット。明確なトレンド → 普遍/回帰クリギング。
• サンプリングが偏っている場合は再クラスタリング（密集した値が重みを持ちすぎないように）。

バリオグラム

import gstools as gs

bin_center, gamma = gs.vario_estimate((x, y), values, max_dist=dmax)  # dmax ≈ half extent
model = gs.Exponential(dim=2)
model.fit_variogram(bin_center, gamma, nugget=True)
print(model)   # report: nugget, sill, range — plus the fitted plot

• ラグ間隔 ≈ ドメインの 1/3、1ラグあたり30～50ペア。ビンあたりのペア数が少なすぎるとノイズになります。
• 方向性バリオグラム（0/45/90/135°）で異方性を確認します。地質学的・環境場は異方性であることが多く、無視するのではなく異方性モデルを当てはめてください。
• パラメータを解釈して報告します：ナゲット（測定誤差 + 微小スケール分散）、シル、レンジ（相関距離）。ナゲットがシルに近い場合、データはほとんど情報を支えていません。正直にそう報告してください。
• バリオグラムのレンジを超えて意味のある補間をしないこと。超える場合はそのセルをマスクするか、フラグを立ててください。

クリギングの実行

from pykrige.ok import OrdinaryKriging

ok = OrdinaryKriging(x, y, values, variogram_model="exponential",
                     variogram_parameters={"sill": s, "range": r, "nugget": n},
                     coordinates_type="euclidean")
z, ss = ok.execute("grid", gridx, gridy)   # ss = kriging VARIANCE — keep it!

• 投影 CRS で作業します（明示的に地理対応のモデルを使わない限り、度単位のクリギング距離は誤りです）。
• 2つのラスターを提供します：予測値 AND クリギング標準偏差（np.sqrt(ss)）。不確実性マップは、次にどこでサンプリングすべきか、どこでマップを信用すべきでないかを示します。
• 探索近傍：最近傍16～32点が典型です。文書化してください。

バリデーション — 一個抜き法とその先

• 小さい n には LOOCV：ME（バイアス ≈ 0？）、RMSE、標準化 RMSE（クリギング分散が誠実なら ≈ 1）を報告します。
• クラスタ化したサンプルでは LOOCV は良く見えすぎます。空間ブロック CV（標準レシピ：ml-experiment-standards → references/spatial-cv-protocol.md）も実行し、両方を報告してください。
• 実測値 vs 予測値の散布図を 1:1 線とともに。CV 残差をマップします — 空間的にクラスタ化した残差は、欠落したトレンド/共変量を示します。
• 素朴なベースライン（大域平均、IDW）と比較します：クリギングはその複雑さに見合う成果を上げなければなりません。

レポートテンプレート

## Surface: <variable>
- n = <>, extent, CRS, transform applied: <log/none>
- Variogram: <model>, nugget/sill/range = <...>, anisotropy: <...>
- Method: <OK/UK/RK + covariates>
- LOOCV: ME <>, RMSE <>, RMSE_std <>; block-CV RMSE <>
- Uncertainty: kriging SD raster delivered; area beyond reliable range masked

落とし穴チェックリスト

• プロットを見ずに自動フィットしたバリオグラムをそのまま受け入れる。
• 広い範囲で EPSG:4326 のままクリギングする。
• 分散/SD を伴わない予測ラスターだけを納品する。
• 対数クリギングを単純な exp() で逆変換する（バイアス！）。
• データの凸包をはるかに超えた外挿を同等の信頼度で提示する。
• 観測点座標の重複がフィットをクラッシュさせる、または黙って歪める。

実行契約

• ワークフロー: サンプリングと変換のニーズを調査し、空間依存性をモデル化し、候補サーフェスを適合させ、単純なベースラインに対して検証し、不確実性を定量化し、サポートされない外挿をマスクします。
• 判断ルール: 透明なベースラインには IDW またはトレンド法を選択し、防御可能なバリオグラムが存在する場合はクリギングを、共変量がブロック検証を改善する場合は回帰クリギングを選択します。
• 検証プロトコル: バリオグラム診断、LOOCV および空間ブロック CV 残差、ベースラインの改善、不確実性のキャリブレーション、サンプルハルを超えた挙動を検証します。
• 失敗モード: サンプリングが疎すぎる場合、依存性が存在しない場合、重複が支配的な場合、距離単位が無効な場合、または不確実性がキャリブレーションされていない場合は、自信のあるサーフェスを保留します。
• 成果物: 予測サーフェス、不確実性サーフェス、バリオグラムとパラメータ、検証テーブル、サンプリングサポートマスク、および手法の限界。
• ソースの鮮度: パッケージの API やデフォルトに依存する前に、権威あるソースレジストリを参照し、確認日を記録してください。