scikit-gstat

SciKit-GStat - 地球統計学

クイックリファレンス

import skgstat as skg
import numpy as np

# Create variogram
V = skg.Variogram(coordinates=coords, values=values, n_lags=15)

# Fit model
V.model = 'spherical'
print(f"Range: {V.parameters[0]:.2f}, Sill: {V.parameters[1]:.2f}")

# Kriging interpolation
ok = skg.OrdinaryKriging(V)
predictions = ok.transform(grid_coords)

主なクラス

クラス	用途
Variogram	経験バリオグラムと理論バリオグラム
OrdinaryKriging	空間相関を用いた補間
DirectionalVariogram	異方性バリオグラム
SpaceTimeVariogram	時空間解析

基本的な操作

バリオグラムの作成と適合

import skgstat as skg

V = skg.Variogram(
    coordinates=coords,      # (n, 2) array of x, y
    values=values,           # (n,) array of measurements
    n_lags=15,
    maxlag='median'          # or specific distance
)

# Fit model: 'spherical', 'exponential', 'gaussian', 'matern', 'stable'
V.model = 'spherical'

# Get parameters
print(f"Range: {V.parameters[0]:.2f}")
print(f"Sill: {V.parameters[1]:.2f}")
print(f"Nugget: {V.parameters[2]:.2f}")
print(f"RMSE: {V.rmse:.4f}")

通常クリギング

import skgstat as skg
import numpy as np

V = skg.Variogram(coords, values, model='spherical')
ok = skg.OrdinaryKriging(V)

# Create prediction grid
x = np.linspace(0, 100, 50)
y = np.linspace(0, 100, 50)
xx, yy = np.meshgrid(x, y)
grid_coords = np.column_stack([xx.ravel(), yy.ravel()])

# Predict
predictions = ok.transform(grid_coords)
Z = predictions.reshape(xx.shape)

# Get variance
ok.return_variance = True
predictions, variance = ok.transform(grid_coords)

方向性バリオグラム

import skgstat as skg

DV = skg.DirectionalVariogram(
    coordinates=coords,
    values=values,
    azimuth=45,          # Direction in degrees
    tolerance=22.5,      # Angular tolerance
    bandwidth='q33'      # Perpendicular bandwidth
)

# Check anisotropy
for az in [0, 45, 90, 135]:
    DV.azimuth = az
    print(f"Azimuth {az}: Range = {DV.parameters[0]:.2f}")

交差検証

import skgstat as skg
from sklearn.model_selection import cross_val_score

V = skg.Variogram(coords, values, model='spherical')
ok = skg.OrdinaryKriging(V)

scores = cross_val_score(ok, coords, values, cv=5, scoring='neg_mean_squared_error')
print(f"CV RMSE: {np.sqrt(-scores.mean()):.4f}")

ロバスト推定量

import skgstat as skg

# Use robust estimator for noisy data
V = skg.Variogram(
    coords, values,
    estimator='cressie'  # 'matheron', 'cressie', 'dowd', 'genton'
)

モデルのクイックリファレンス

モデル	特性
spherical	最も一般的で、原点付近では線形です
exponential	シルに到達せず、漸近的に近づきます
gaussian	原点付近では放物線状で、滑らかです
matern	滑らかさを柔軟に制御できます

代替手段との使い分け

ユースケース	ツール	理由
バリオグラム解析とクリギング	scikit-gstat	モダンな API、sklearn 互換
GSLIB 形式のシミュレーション（SGSIM）	GeostatsPy	GSLIB シミュレーションエンジンを完備
トレンド／ドリフトを伴うクリギング	pykrige	普遍クリギング、回帰クリギング
確率場の生成	gstools	柔軟な共分散、SRF 生成
時空間バリオグラム	scikit-gstat	SpaceTimeVariogram を組み込み
本番環境の地質モデリング	SGeMS / Petrel	GUI、大規模な 3D モデル
ロバストなバリオグラム推定	scikit-gstat	Cressie、Dowd、Genton 推定量
ML パイプラインとの統合	scikit-gstat	sklearn の fit／transform インターフェース

scikit-gstat を選ぶ場合：バリオグラムの適合とクリギングに、Pythonic で scikit-learn 互換の API が必要な場合に適しています。交差検証や ML パイプラインとの統合を伴う、探索的な地球統計解析に最適です。

GeostatsPy を選ぶ場合：GSLIB 互換のシミュレーションワークフロー（SGSIM、SISIM）が必要な場合や、従来の地球統計学の慣例に沿って作業する場合に適しています。

pykrige を選ぶ場合：外部ドリフト変数を用いた普遍クリギングや、地球統計学と機械学習の予測を組み合わせた回帰クリギングが必要な場合に適しています。

一般的なワークフロー

バリオグラムの適合と通常クリギング

• 空間データを numpy 配列（座標と値）として読み込む
• 適切な n_lags と maxlag を指定して Variogram オブジェクトを作成する
• 推定量を比較する：matheron（デフォルト）と、ノイズの多いデータ向けの cressie（ロバスト）
• 複数のモデル（spherical、exponential、gaussian）を適合させ、RMSE を比較する
• 0、45、90、135 度の DirectionalVariogram で異方性を確認する
• RMSE と視覚的な適合度に基づいて最適なモデルを選択する
• 適合済みのバリオグラムから OrdinaryKriging オブジェクトを作成する
• 予測グリッドを定義し、ok.transform(grid_coords) を実行する
• クリギング分散を取得するには ok.return_variance = True を設定する
• cross_val_score() で交差検証し、予測品質を評価する
• 予測値とクリギング分散をマッピングする

よくある問題

問題	解決策
バリオグラムが平坦または不安定	n_lags と maxlag を調整する（maxlag='median' を試す）
モデルの適合度が低い（RMSE が高い）	別のモデルタイプやネスト構造を試す
クリギングが遅すぎる	条件付けに使う点の数またはグリッド解像度を減らす
ナゲットが大きすぎる	測定誤差を示している可能性があります。ロバスト推定量を試してください
異方性が不明瞭	DirectionalVariogram の角度許容差を小さくする

ヒント

1. Maxlag は調査領域の対角線の約 50% に設定します
2. ノイズの多いデータには ロバスト推定量（cressie、dowd）を使用します
3. 複数のモデルを試し、RMSE を比較します
4. クリギングの前に、方向性バリオグラムで 異方性を確認します

参考資料

• バリオグラムモデル - モデルの式とパラメーター
• クリギング手法 - クリギングの種類と設定

スクリプト

• scripts/variogram_analysis.py - バリオグラム解析の完全なワークフロー