KDTreeで最寄り点を探すときは、distance_upper_boundで許容距離を指定し、未対応の返り値を除いてから添字として使う。見つからなかった点には距離infと点数と同じ添字が返るため、そのまま配列へ渡すと範囲外参照になる。遠すぎる点を無理に対応付けない設計が重要である。
最小例で確かめる
以下は説明用に作った小さなデータである。実測データや実行速度の測定結果ではない。コード全体をexample.pyとして保存すれば、入力ファイルを別途用意せずに実行できる。assertは、この例で成り立つべき形や値を確認するために入れてある。
import numpy as np
import matplotlib.pyplot as plt
from scipy.spatial import KDTree
points = np.array([[0.,0.], [1.,0.], [0.,1.]]) # m
values = np.array([10.,20.,30.])
queries = np.array([[0.1,0.1], [0.9,0.1], [2.,2.]])
limit = 0.25
assert np.isfinite(points).all() and np.isfinite(queries).all()
tree = KDTree(points, copy_data=True)
distance, index = tree.query(queries, k=1, eps=0., distance_upper_bound=limit)
matched = np.isfinite(distance) & (index < tree.n)
result = np.full(len(queries), np.nan)
result[matched] = values[index[matched]]
print("distances:", np.round(distance, 6).tolist())
print("indices:", index.tolist())
print("matched:", matched.tolist())
print("values:", result.tolist())
assert index.tolist() == [0,1,tree.n]
np.testing.assert_allclose(distance[:2], np.sqrt(0.02))
assert np.isinf(distance[2]) and not matched[2]
np.testing.assert_allclose(result, [10.,20.,np.nan], equal_nan=True)
# 全点との直接距離でも、小さい例の対応を検査する。
full_distance = np.linalg.norm(queries[:, None, :] - points[None, :, :], axis=2)
np.testing.assert_array_equal(index[matched], full_distance[matched].argmin(axis=1))
fig, ax = plt.subplots(figsize=(6, 5), layout="constrained")
ax.scatter(points[:,0], points[:,1], s=80, color="#26778c", label="reference points")
ax.scatter(queries[matched,0], queries[matched,1], marker="x", s=90,
color="#a83262", label="matched queries")
ax.scatter(queries[~matched,0], queries[~matched,1], marker="x", s=90,
color="#e38619", label="unmatched query")
for point in points:
ax.add_patch(plt.Circle(point, limit, fill=False, color="#26778c", alpha=0.4))
for q, i in zip(queries[matched], index[matched]):
ax.plot([q[0], points[i,0]], [q[1], points[i,1]], color="#a83262")
ax.set(xlabel="x [m]", ylabel="y [m]", title="Nearest point within 0.25 m",
xlim=(-0.4,2.4), ylim=(-0.4,2.4))
ax.set_aspect("equal", adjustable="box")
ax.legend(loc="upper left", fontsize=9)
fig.savefig("PY090-neighbors.png", dpi=160)
plt.close(fig)
実行結果
distances: [0.141421, 0.141421, inf]
indices: [0, 1, 3]
matched: [True, True, False]
values: [10.0, 20.0, nan]
未対応には特別な値が返る
tree.nは登録した点の数で、この例では3である。有効な添字は0、1、2なので、未対応のindex=3はvalues[3]として参照できない。距離がinfであることと添字が範囲内であることをマスクで確認し、対応した行だけ値を取り出す。残りの結果は初期値のNaNとして残している。
np.where(matched, values[index], np.nan)のように書くと、引数を評価する段階でvalues[index]が実行され、未対応の添字でも参照しようとする。先にindex[matched]へ絞ってから取り出すのが安全である。未対応の意味を表す値と、配列へ使える添字を混同しないようにしたい。
座標の単位と距離を合わせる
この例は同じ平面のx、yをメートルで表し、既定のユークリッド距離を使っている。図の円は半径0.25mの範囲である。xがメートルでyが秒など異なる尺度を持つ特徴量なら、そのまま距離を計算する意味を先に決める必要がある。標準化も距離の定義を変更する操作になる。
緯度経度をそのまま平面座標として扱い、distance_upper_boundをメートルのつもりで渡してはいけない。地理的な距離が必要なら、適切な座標変換や球面上の距離を使う方法を選ぶ。利用する距離の種類、座標系、単位をセットで記録しておこう。
入力と出力のshapeを確認する
pointsは点数×次元数、queriesは問い合わせ数×次元数である。k=1では距離と添字の最後の軸が省略され、この例はどちらも長さ3の1次元になる。k=[1]や複数近傍を指定するとshapeが変わるので、同じマスク代入をそのまま流用する前に返り値を確認する。
eps=0は近似探索の許容幅を追加しない指定である。近似探索を使う場合は、結果の距離保証と採用条件の関係を確認したい。距離が上限とほぼ等しい点には浮動小数点や境界の扱いも関係するため、許容距離の境界を含めるという業務規則が重要なら、候補を取得した後の距離判定を明示的に設計する。
最近傍は一対一対応ではない
複数のqueryが同じ基準点を最寄りとして返すことは正常である。全点を一対一に割り当てたい問題なら、KDTree.queryだけではその制約を満たさない。同距離の候補がある場合も、特定の添字が選ばれることを前提に意味付けしない。距離だけで十分か、IDや属性条件が必要かを考える。
copy_data=Trueは木の作成に使った座標の独立したコピーを保持させる指定である。元の配列を後から書き換えて探索構造と内容が食い違うのを避けるため、ここでは明示した。大規模データではメモリ負担との兼ね合いがあるので、元配列を変更しない契約にする方法も含めて選ぶ。
例では小さい配列に限り全組合せの直接距離を計算し、対応した点の添字を照合した。これはKDTreeの高速性を測った結果ではない。未対応件数、最大採用距離、代表点の位置を確認すれば、入力座標の単位違いや遠すぎる対応を早めに見つけやすい。
動作確認環境と参考資料
Linux・CPython 3.12.14、NumPy 2.3.5、pandas 2.2.3、SciPy 1.17.0、Matplotlib 3.10.8の環境で掲載コードを実行した。使用するライブラリはコード冒頭のimportを参照してほしい。公式資料の最新版と、この実行確認版は区別している。数値の末尾や表の表示幅は環境によって変わることがある。