0
0

Delete article

Deleted articles cannot be recovered.

Draft of this article would be also deleted.

Are you sure you want to delete this article?

Pythonで学ぶ経路生成②-2|A* のヒューリスティック関数はどう設計するか — 距離の測り方で探索は16倍変わる

0
Posted at

速いのに、壊れている

前回のA*(エースター)は、同じ最短コスト26.97を保ったまま、探索を 369 → 113セルに減らしました。ではヒューリスティック(heuristic:ゴールまであとどれくらいか、の見積もりを返す関数)をマンハッタン距離(縦の差+横の差。碁盤の目を縦横だけで進む距離)に差し替えると、どうなるか。

展開ノードは 113 → 23セル。ダイクストラ比で 369 ÷ 23 = 約16倍の削減です。コストは 26.97 のままで、最短経路も出ています。4つ試した中でいちばん速く、しかも答えも合っていました。

ところが、この設定は壊れています。

スタート地点での見積もりは $h$(start) = 34.00。本当のコストは 26.97。見積もりが実コストを 7.03 も上回っています。A*が最短を保証するための条件(許容性)を、このマップの全381セルのうち271セルで破っていました。最短経路が出たのは、たまたまです。

この記事は、その「速いのに間違っている」を入口に、h に何を入れるべきかを実測で詰めます。読み終えると次の3つが手に入ります。

  • 距離の測り方4種(マンハッタン/ユークリッド/オクタイル/重み付き)の使い分けが、移動モデル(4近傍か8近傍か)から機械的に決まること
  • 許容性(admissibility)がなぜ最短保証に必要かを、全セル数え上げで確認する手順
  • 重み付きA*($w>1$)は「保証が消える」だけで「実際に壊れる」とは別であること。実測では $w=2$ でも最短は壊れませんでした

本記事は連載②(A*)の補講です。②の章末に「$h$ を盛って、やりすぎの閾値を探す」という宿題を残したので、その答え合わせにあたります。連載の本線である③はポテンシャル法で、そちらは予定どおり別記事で公開します。

まだ①②を読んでいない方は、先にこちらからどうぞ。本記事は②のコードをそのまま使います。

前提:連載のどこにいるか、何を書き足すか

図:連載のマップ。本記事は本線から枝分かれした②の補講で、次に来るのは③ポテンシャル法。

  • 想定読者:A*を一度実装したものの、heuristic() の中身を何にすればいいのか決めきれない方
  • 環境:マップを持つ NumPy だけ(優先度キューの heapq、時間計測の timeit、math は標準ライブラリ。図を描くときだけ Matplotlib)。実測環境は Apple M1 Pro / macOS 26.3 / Python 3.10.9 / NumPy 1.26.4
  • マップ:連載①②と完全に同じもの(20×20・中央に縦の壁と1マスの隙間・ぬかるみ帯は重み3)。同じ土俵で比べるためです

ファイルは②のものを流用し、本筋で新しく書くのは2つだけです(残り4本は検証用で、折りたたみに入れてあります)。

ファイル ②から 役割
constant.py / grid_map.py / dijkstra.py / astar.py / plot.py 流用(変更なし) 定数・マップ・ダイクストラ・A*本体・可視化の下地
heuristics.py 新規 距離の測り方4種と、重み付きを作る関数
compare_heuristics.py 新規 $h$ を差し替えて横並び比較する(本記事の主役)
count_touched.py / admissible_audit.py / sweeps.py / plot_heuristics.py 新規(折りたたみ) 数え分け・許容性の全数チェック・総当たり走査・作図

⚠️ A*本体(astar.py)は1文字も変えません。 変えるのは「どの関数を h として使うか」だけです。

h の仕事は「残りコストの下限」を答えること

②の復習です。A*は優先度キューから取り出す順番を、$g$(スタートからの実コスト)ではなく $f$ で決めます。

f(v) = g(v) + h(v)

$g$ はもう分かっている確定値なので、設計の自由度があるのは $h$ だけです。つまり A*の性能を決めるつまみは $h$ ひとつしかありません。これは実測で確かめられます。heuristic() が常に 0.0 を返すようにすると、A*の探索はダイクストラと完全に同じになります。展開ノード数が369で一致するだけでなく、返ってきた explored(確定させた順に並んだセルのリスト)を == で比べても True。確定の順序まで1つもずれません。

では $h$ は大きければ大きいほど良いのか。ここに一線があります。

許容性(admissibility):$h$ は、そのセルからゴールまでの本当のコストを、超えてはいけない。

超えると何が起きるか。A*は「あっちは遠そうだ」と判断して、まだ短い可能性が残っている道を早々に捨てます。捨てた先に本当の近道があれば、最短経路を取りこぼす。逆に $h$ が本当のコスト以下である限り、取りこぼしは起きません。

図:$h$ の大きさと A* の性質。右へ行くほど探索は速くなり、真の残りコストを追い越した瞬間に最短保証を失う。ねらいは「真値のすぐ下」。

$h$ の設計とは、この線を越えない範囲で、できるだけ大きな値を返す関数を作ることです。「大きいほど速い、でも越えたら保証が消える」というトレードオフを、以降は数値で見ていきます。

⚠️ 厳密には、許容性だけでは足りない場合があります。 本記事のA*は、いちど確定させたセル(visited)を二度と見直しません。この実装で最短を保証するには、許容性に加えて一貫性(consistency:隣のセルへ1歩動いたときの $h$ の減り方が、その1歩のコストを超えないこと。$h(u) - h(v) \le c(u,v)$)も必要です。本記事で許容的だった3つ($h=0$・ユークリッド・オクタイル)は、全セル・全辺で一貫性も満たしていました(違反辺0)。破っていたのはマンハッタン(301辺)と重み付き($w=2$ で624辺)で、こちらはそもそも許容的でもありません。

グリッドで使える4つの距離

$h$ の中身は「ゴールまでまっすぐ進めたら何コストか」の見積もりです。どう進めるか(移動モデル)が変われば、正しい見積もりも変わります。

距離 式($dr$=縦の差、$dc$=横の差) 整合する移動モデル 8近傍で許容的か
マンハッタン $dr + dc$ 4近傍(上下左右のみ) ❌(斜めを2と数えるので過大)
ユークリッド $\sqrt{dr^2 + dc^2}$ 連続空間・任意の角度に進める ✅(ただし甘い)
オクタイル $\sqrt{2}\min(dr,dc) + \bigl(\max(dr,dc)-\min(dr,dc)\bigr)$ 8近傍(斜めコスト$\sqrt{2}$) ✅(いちばんタイト)
重み付き $w \times$(上のどれか) 最短でなくてよい場面 $w>1$ なら ❌

この4つを1ファイルにまとめます。すべて同じシグネチャ(セルとゴールを受け取って float を返す)にしておくのが要点で、こうしておけば後から差し替えるだけで比較できます。

# ヒューリスティック関数(heuristic:ゴールまであとどれくらいか、の見積もりを返す関数)の詰め合わせ

# ライブラリの読み込み
import math                        # 平方根・斜辺の計算
from typing import Callable        # 「関数を受け取る/返す」ことを型で表す

# 自作モジュールの読み込み
from constant import COST_STRAIGHT, COST_DIAGONAL    # 1マス移動の基本コスト(縦横1・斜め√2)

# ヒューリスティック関数の型(セルとゴールを受け取って推定コストを返す関数)
Heuristic = Callable[[tuple, tuple], float]


def h_zero(cell: tuple, goal: tuple) -> float:
    """
    常に 0 を返す(=ゴールまでを一切見積もらない)

    このとき f = g + 0 = g となり、A* はダイクストラ法そのものに戻る。

    パラメータ
        cell: 現在セル (row, col)
        goal: ゴールセル (row, col)

    戻り値
        h: つねに 0.0
    """
    return 0.0


def h_manhattan(cell: tuple, goal: tuple) -> float:
    """
    マンハッタン距離(縦の差+横の差。碁盤の目を縦横だけで進む距離)を返す

    4近傍(上下左右のみ)の移動モデルに整合したヒューリスティック。
    斜めに進める8近傍で使うと、斜め1歩(√2≒1.41)を「縦1+横1=2」と数えるため
    実コストを過大評価する(=許容的でなくなる)。

    パラメータ
        cell: 現在セル (row, col)
        goal: ゴールセル (row, col)

    戻り値
        h: ゴールまでの推定コスト
    """
    dr = abs(cell[0] - goal[0])
    dc = abs(cell[1] - goal[1])
    return COST_STRAIGHT * (dr + dc)


def h_euclidean(cell: tuple, goal: tuple) -> float:
    """
    ユークリッド距離(定規で測るまっすぐな直線距離)を返す

    グリッドの上を斜めにしか進めない状況でも必ず直線距離以上かかるので、
    8近傍でも許容的(過大評価しない)。ただし見積もりが甘いぶん探索は広がる。

    パラメータ
        cell: 現在セル (row, col)
        goal: ゴールセル (row, col)

    戻り値
        h: ゴールまでの推定コスト
    """
    dr = abs(cell[0] - goal[0])
    dc = abs(cell[1] - goal[1])
    return math.hypot(dr, dc)    # hypot: 直角三角形の斜辺の長さ √(dr² + dc²)


def h_octile(cell: tuple, goal: tuple) -> float:
    """
    オクタイル距離(斜めに進める前提での格子上の最短距離)を返す

    「斜めで詰められるぶんは斜め(√2)、残りはまっすぐ(1)」と見積もる。
    8近傍・斜めコスト√2の移動モデルにぴったり整合する。連載②で採用したもの。

    パラメータ
        cell: 現在セル (row, col)
        goal: ゴールセル (row, col)

    戻り値
        h: ゴールまでの推定コスト
    """
    dr = abs(cell[0] - goal[0])
    dc = abs(cell[1] - goal[1])
    diagonal = min(dr, dc)                # 斜めで詰められるマス数
    straight = max(dr, dc) - diagonal     # 残りをまっすぐ進むマス数
    return COST_DIAGONAL * diagonal + COST_STRAIGHT * straight


def make_weighted(base: Heuristic, weight: float) -> Heuristic:
    """
    既存のヒューリスティックを weight 倍にした「重み付きヒューリスティック」を作る

    weight > 1 にすると h が実コストを超えうるので、最適性の保証は失われる
    (=保証が消える。実際に最短が壊れるかどうかはマップ次第)。

    パラメータ
        base: もとにするヒューリスティック関数
        weight: h に掛ける倍率(1.0 でもとの関数と同じ)

    戻り値
        weighted: weight 倍した値を返すヒューリスティック関数
    """
    def weighted(cell: tuple, goal: tuple) -> float:
        return weight * base(cell, goal)

    return weighted

差し替えは1行です。astar.py はモジュール内の heuristic を呼んでいるので、外から代入すれば中身が入れ替わります。

# ヒューリスティックを1行だけ差し替えて、A* の探索がどう変わるか見る

# 自作モジュールの読み込み
import astar as A                        # モジュールごと読む(中の heuristic を差し替えるため)
from constant import CONNECTIVITY
from grid_map import make_sample_map
from heuristics import h_manhattan

grid_map = make_sample_map()

A.heuristic = h_manhattan    # ← これが「差し替え1行」。astar.py 自体はいじらない

path, explored, cost = A.astar(grid_map, CONNECTIVITY.EIGHT)
print(f"展開ノード={len(explored)}  コスト={cost:.2f}  通過点={len(path)}")
print(f"h(start)={h_manhattan(grid_map.start, grid_map.goal):.2f}")
展開ノード=23  コスト=26.97  通過点=23
h(start)=34.00

⚠️ この差し替え方(モジュールの関数を外から代入する=モンキーパッチ)は比較実験向けの近道です。 実務のコードなら astar(grid_map, connectivity, heuristic=h_octile) のように引数で渡すほうが素直です。ここでは②のコードを1文字も変えずに比較したいので、あえてこの手を使っています。

実測比較:距離の測り方で展開ノードは16倍変わる

比較スクリプトです。8種類の設定を同じマップで順に解き、展開ノード数・コスト・実行時間・h(start)・許容性・最適性を並べます。

# ヒューリスティックを差し替えて A* を横並び比較する(同一マップ・8近傍)

# ライブラリの読み込み
import timeit    # 実行時間の計測(同じ処理を何回も走らせて最速値を取る)

# 自作モジュールの読み込み
import astar as A                        # A* 本体(モジュールごと読む=あとで heuristic を差し替えるため)
from constant import CONNECTIVITY        # 近傍の種類(4近傍 / 8近傍)
from dijkstra import dijkstra            # 比較の基準(最適コストの正解を出す)
from grid_map import make_sample_map     # 連載①②と同じサンプルマップ
from heuristics import h_zero, h_manhattan, h_euclidean, h_octile, make_weighted

REPEAT = 30    # 実行時間は30回まわして最速値を採用する(他プロセスの影響を減らす)


def measure(grid_map, heuristic) -> tuple:
    """
    ヒューリスティックを差し替えて A* を1回解き、結果と最速実行時間を返す

    パラメータ
        grid_map: 占有格子地図(スタート・ゴール設定済み)
        heuristic: 使うヒューリスティック関数

    戻り値
        explored_n: 展開ノード数(確定させたセルの数)
        cost: 得られた経路のコスト
        best_ms: REPEAT 回まわしたうちの最速実行時間 [ミリ秒]
        h_start: スタート地点での h の値
    """
    A.heuristic = heuristic    # astar.py が参照している heuristic を差し替える(1行)

    path, explored, cost = A.astar(grid_map, CONNECTIVITY.EIGHT)
    times = timeit.repeat(lambda: A.astar(grid_map, CONNECTIVITY.EIGHT),
                          number=1, repeat=REPEAT)
    h_start = heuristic(grid_map.start, grid_map.goal)

    return len(explored), cost, min(times) * 1000, h_start


def main() -> None:
    grid_map = make_sample_map()

    # 基準:ダイクストラ法(最適コストの正解と、展開ノード数の上限)
    d_path, d_explored, d_cost = dijkstra(grid_map, CONNECTIVITY.EIGHT)
    print(f"[基準] ダイクストラ  展開={len(d_explored)}  コスト={d_cost:.2f}")
    print(f"       (以降の「最適?」は、このコスト {d_cost:.2f} と一致するかで判定する)\n")

    # 比較するヒューリスティック(名前, 関数)
    settings = [
        ("h_zero",       h_zero),
        ("manhattan",    h_manhattan),
        ("euclidean",    h_euclidean),
        ("octile",       h_octile),
        ("octile w=1.5", make_weighted(h_octile, 1.5)),
        ("octile w=2",   make_weighted(h_octile, 2.0)),
        ("octile w=3",   make_weighted(h_octile, 3.0)),
        ("octile w=5",   make_weighted(h_octile, 5.0)),
    ]

    # 列見出しは半角のみ(全角文字を混ぜると桁がずれるため)
    header = (f"{'h':<14}{'expand':>7}{'ratio':>8}{'cost':>8}"
              f"{'best_ms':>9}{'h(start)':>10}{'admissible':>12}{'optimal':>10}")
    print(header)
    print("-" * len(header))

    for name, heuristic in settings:
        explored_n, cost, best_ms, h_start = measure(grid_map, heuristic)
        ratio      = len(d_explored) / explored_n                 # ダイクストラ比の展開ノード削減倍率
        admissible = "OK" if h_start <= d_cost + 1e-9 else "NG"   # h(start) が実コストを超えていないか
        optimal    = "OK" if abs(cost - d_cost) < 1e-9 else f"NG(+{cost - d_cost:.2f})"
        print(f"{name:<14}{explored_n:>7}{ratio:>7.1f}x{cost:>8.2f}"
              f"{best_ms:>9.3f}{h_start:>10.2f}{admissible:>12}{optimal:>10}")


if __name__ == "__main__":
    main()

実行します。

python compare_heuristics.py
[基準] ダイクストラ  展開=369  コスト=26.97
       (以降の「最適?」は、このコスト 26.97 と一致するかで判定する)

h              expand   ratio    cost  best_ms  h(start)  admissible   optimal
------------------------------------------------------------------------------
h_zero            369    1.0x   26.97    1.821      0.00          OK        OK
manhattan          23   16.0x   26.97    0.136     34.00          NG        OK
euclidean         138    2.7x   26.97    0.774     24.04          OK        OK
octile            113    3.3x   26.97    0.663     24.04          OK        OK
octile w=1.5       25   14.8x   26.97    0.166     36.06          NG        OK
octile w=2         24   15.4x   26.97    0.159     48.08          NG        OK
octile w=3         23   16.0x   26.97    0.154     72.12          NG        OK
octile w=5         19   19.4x   35.94    0.133    120.21          NG NG(+8.97)

読みやすく並べ替えたのが下の表です。「速いか」と「最短が保証されるか」は別の列である、という一点に注目してください。

$h$(距離の測り方) 展開ノード ダイクストラ比 経路コスト $h$(start) 許容的か このマップの結果
$h=0$(=ダイクストラ) 369 1.0倍 26.97 0.00 ✅ ✅ 最短(保証あり)
ユークリッド 138 2.7倍 26.97 24.04 ✅ ✅ 最短(保証あり)
オクタイル(②で採用) 113 3.3倍 26.97 24.04 ✅ ✅ 最短(保証あり)
マンハッタン 23 16.0倍 26.97 34.00 ❌ ⚠️ 最短だが偶然(保証なし)
オクタイル×1.5 25 14.8倍 26.97 36.06 ❌ ⚠️ 最短だが保証なし
オクタイル×2 24 15.4倍 26.97 48.08 ❌ ⚠️ 最短だが保証なし
オクタイル×3 23 16.0倍 26.97 72.12 ❌ ⚠️ 最短だが保証なし
オクタイル×5 19 19.4倍 35.94 120.21 ❌ ❌ 壊れた(+8.97)

表:同一マップ・8近傍での実測。「許容的か」は $h$(start) が実コスト26.97以下かで判定している。

探索範囲を絵にすると、この差がそのまま見えます。

ヒューリスティック別の探索範囲の比較
図:水色が展開ノード(確定させたセル)、橙が得られた経路。$h$ を強くするほど水色が細り、マンハッタンでは経路の上だけになる。4枚とも経路コストは26.97で同じ。

実測から言えることを3つに絞ります。

1. 展開ノードの削減はユークリッド2.7倍、オクタイル3.3倍、マンハッタン16.0倍。 ユークリッドとオクタイルは $h$(start) がどちらも24.04で同じですが、探索量は138と113で違います。理由はスタート以外のセルにあります。オクタイルは全381セルで必ずユークリッド以上の値を返し(うち323セルで厳密に大きい)、平均は 14.10 対 13.45。真の残りコストの平均が16.48なので、オクタイルのほうが天井に近い=タイトです。許容性を守りながら値を上げると、それだけ探索が細る。この関係が数値で見えます。

2. 実行時間は展開ノード数にほぼ比例しますが、環境依存です。 私の手元では 1.821ms → 0.136ms(約13倍)。best-of-30(30回まわして最速値)で測っていますが、機種やPythonのバージョンで変わるので、倍率のほうを見てください。

3. 「展開ノード16倍」は、いちばん都合のよい数え方です。 展開ノードは「優先度キューから取り出して確定させたセル」で、キューに入れただけのセルは含みません。両方数えるとこうなります。

$h$ 展開(確定) 触ったセル(コストを計算した) キューへの push 回数
$h=0$ 369 379 421
マンハッタン 23 87 90
ユークリッド 138 193 289
オクタイル 113 172 251

表:数え方を変えると倍率も変わる。触ったセルで見ると 379 → 87 で約4.4倍。「16倍」は展開ノードという指標での話。

この数え分けをするスクリプト(count_touched.py・クリックで展開)
# 「展開ノード(確定させたセル)」と「触ったセル(キューに入れたセル)」を数え分ける

# ライブラリの読み込み
import heapq    # 優先度キュー

# 自作モジュールの読み込み
from constant import CONNECTIVITY
from grid_map import GridMap, make_sample_map
from heuristics import Heuristic, h_zero, h_manhattan, h_euclidean, h_octile, make_weighted


def astar_counted(grid_map: GridMap, heuristic: Heuristic,
                  connectivity: CONNECTIVITY = CONNECTIVITY.EIGHT) -> tuple:
    """
    A* を解きながら「確定させたセル数」と「キューに入れたセル数」を数える

    astar.py の astar() に数え上げだけを足したもの(探索の中身は同じ)。

    パラメータ
        grid_map: 占有格子地図(スタート・ゴール設定済み)
        heuristic: 使うヒューリスティック関数
        connectivity: 近傍の種類(4近傍 / 8近傍)

    戻り値
        explored_n: 展開ノード数(優先度キューから取り出して確定させたセル数)
        touched_n: 一度でもコストを計算してキューに入れたセル数
        pushes: キューへ入れた回数(同じセルを何度も入れ直すのでセル数より多い)
    """
    start = grid_map.start
    goal  = grid_map.goal

    dist     = {start: 0.0}
    visited  = set()
    explored = []
    pushes   = 0

    queue = [(heuristic(start, goal), 0.0, start)]

    while queue:
        f, cost, current = heapq.heappop(queue)
        if current in visited:
            continue
        visited.add(current)
        explored.append(current)

        if current == goal:
            break

        for n_cell, move_cost in grid_map.neighbors(current, connectivity):
            new_cost = cost + move_cost
            if n_cell not in dist or new_cost < dist[n_cell]:
                dist[n_cell] = new_cost
                heapq.heappush(queue, (new_cost + heuristic(n_cell, goal), new_cost, n_cell))
                pushes += 1

    return len(explored), len(dist), pushes


def main() -> None:
    grid_map = make_sample_map()

    settings = [
        ("h_zero",     h_zero),
        ("manhattan",  h_manhattan),
        ("euclidean",  h_euclidean),
        ("octile",     h_octile),
        ("octile w=2", make_weighted(h_octile, 2.0)),
    ]

    header = f"{'h':<12}{'explored':>10}{'touched':>9}{'pushes':>8}"
    print(header)
    print("-" * len(header))
    for name, heuristic in settings:
        explored_n, touched_n, pushes = astar_counted(grid_map, heuristic)
        print(f"{name:<12}{explored_n:>10}{touched_n:>9}{pushes:>8}")


if __name__ == "__main__":
    main()

マンハッタンはなぜ8近傍で反則になるのか

原因は1マスの数え方です。8近傍では斜めに1歩進めて、そのコストは $\sqrt{2} \fallingdotseq 1.41$。ところがマンハッタン距離はその1歩を「縦1+横1=2」と数えます。

斜め1歩の見積もりの差(実コスト√2 と マンハッタン 2)
図:斜め1歩の見積もりの差。斜めに進むぶんだけマンハッタンは上に外れ、真斜めに向かうときは実コストの $\sqrt{2}$ 倍まで盛る。

実測すると、マンハッタンの値はオクタイルの 1.000〜1.4142倍(真横・真下なら等しく、真斜めで最大の $\sqrt{2}$ 倍)に収まりました。真の残りコストとの比も最大1.414倍。つまり 8近傍でマンハッタンを使うのは、気づかないうちに「$w$ が最大 $\sqrt{2}$ の重み付きA*」をやっているのと同じです。表でマンハッタン(23)とオクタイル×1.5(25)が近い数字になったのは偶然ではありません。

破れの範囲を全数で確認しましょう。ゴールから逆向きにダイクストラを1回まわせば、全セルの「本当の残りコスト」が一度に手に入ります。それと $h$ を突き合わせるだけです。

許容性の全数チェック(admissible_audit.py・クリックで展開)
# h が「本当の残りコスト」を超えているセルを数える(許容性の全数チェック)

# ライブラリの読み込み
import heapq    # 優先度キュー

# 自作モジュールの読み込み
from constant import CONNECTIVITY
from grid_map import GridMap, make_sample_map
from heuristics import h_manhattan, h_euclidean, h_octile, make_weighted


def true_cost_to_goal(grid_map: GridMap,
                      connectivity: CONNECTIVITY = CONNECTIVITY.EIGHT) -> dict:
    """
    全セルについて「そのセルからゴールまでの本当の最小コスト」を求める

    ゴールから逆向きにダイクストラを1回まわすだけで全セル分が一度に得られる。
    移動コストは「入る先のセルの地形」で決まる(=向きで値が変わる)ので、
    辺を明示的に反転させてから探索する。

    パラメータ
        grid_map: 占有格子地図(ゴール設定済み)
        connectivity: 近傍の種類(4近傍 / 8近傍)

    戻り値
        dist: セル -> ゴールまでの最小コスト の辞書(到達できないセルは含まない)
    """
    # 反転した辺を作る: v から u へ、重みは「u から v へ入るコスト」
    reverse = {}
    for row in range(grid_map.rows):
        for col in range(grid_map.cols):
            u = (row, col)
            if not grid_map.is_free(u):
                continue
            for v, move_cost in grid_map.neighbors(u, connectivity):
                reverse.setdefault(v, []).append((u, move_cost))

    dist  = {grid_map.goal: 0.0}
    queue = [(0.0, grid_map.goal)]
    done  = set()

    while queue:
        cost, current = heapq.heappop(queue)
        if current in done:
            continue
        done.add(current)
        for prev, move_cost in reverse.get(current, []):
            new_cost = cost + move_cost
            if prev not in dist or new_cost < dist[prev]:
                dist[prev] = new_cost
                heapq.heappush(queue, (new_cost, prev))

    return dist


def main() -> None:
    grid_map = make_sample_map()
    dist = true_cost_to_goal(grid_map)

    settings = [
        ("manhattan",    h_manhattan),
        ("euclidean",    h_euclidean),
        ("octile",       h_octile),
        ("octile w=1.5", make_weighted(h_octile, 1.5)),
        ("octile w=2",   make_weighted(h_octile, 2.0)),
    ]

    print(f"到達可能なセル数: {len(dist)}")
    header = f"{'h':<14}{'violations':>12}{'worst_over':>12}{'worst_cell':>14}"
    print(header)
    print("-" * len(header))

    for name, heuristic in settings:
        violations = 0
        worst_over = 0.0
        worst_cell = None
        for cell, true_cost in dist.items():
            over = heuristic(cell, grid_map.goal) - true_cost    # 見積もりが実コストを超えた分
            if over > 1e-9:
                violations += 1
                if over > worst_over:
                    worst_over = over
                    worst_cell = cell
        print(f"{name:<14}{violations:>12}{worst_over:>12.2f}{str(worst_cell):>14}")


if __name__ == "__main__":
    main()
到達可能なセル数: 381
h               violations  worst_over    worst_cell
----------------------------------------------------
manhattan              271        7.62       (18, 0)
euclidean                0        0.00          None
octile                   0        0.00          None
octile w=1.5           348       10.45      (19, 11)
octile w=2             375       22.70       (18, 0)

マンハッタンは 381セル中271セルで許容性を破っていました。最悪は左下の (18, 0) で 7.62の盛りすぎ。境界ぎりぎりの話ではなく、マップの7割で反則をしている状態です。

それでも最短が出た理由

では、なぜコスト26.97(最適)が返ってきたのか。展開ノードを調べると答えが出ます。マンハッタンの展開ノード23セルは、最終経路の23セルと完全に同一の集合でした(set(explored) == set(path) が True)。

つまりこのA*は、寄り道をゼロで、スタートからゴールまで一直線に降りていったのです。他の候補を検討していない。にもかかわらず答えが合ったのは、その一直線がこのマップでは最適経路とたまたま一致していたからです。壁の隙間がスタートとゴールを結ぶ斜め線の上にあり、ぬかるみを避ける迂回も斜めの延長線上に乗っていました。

ここで気になるのは「たまたま」がどれくらいの確率なのか、です。ゴールを全マスに順番に置いて380通り試し、マンハッタンA*のコストをダイクストラの最適コストと突き合わせました(走査コードは後の節の sweeps.py)。結果は次のとおりです。

  • 最適を外したゴール:380通り中わずか2通り((3, 18) で +0.586、(3, 19) で +1.172)
  • 残りの378通りでは、非許容のまま最短を返していた

悪化量にも見覚えがあります。+0.586 は $2 - \sqrt{2}$、つまり上の図で見た「斜め1歩ぶんの盛りすぎ」ちょうど1回分。+1.172 はその2回分です。盛った量が、そのまま遠回りの量として出てきています。

7割のセルで反則をしていて、それでも実害が出るのは0.5%。 これが「非許容な $h$ はテストで見つからない」という話の実態です。動かして確かめても、ほぼ通ります。壊れたときだけ、静かに数%遠回りした経路が返ってくる。

⚠️ 「速い」と「正しい」を、同じ実行結果から読み取ってはいけません。 1回の実行で最短が出たことは、最短保証の証拠になりません。保証があるかどうかは $h$ の性質(許容性)で決まり、上のように全セルで突き合わせれば機械的に検査できます。実行結果ではなく $h$ を検査する、が正しい向きです。

34.00 という数字は間違いではない

おもしろいのは、マンハッタンが返した h(start) = 34.00 がでたらめな値ではないことです。同じマップを4近傍(上下左右のみ)で解くと、最適コストは ちょうど34.00。マンハッタン距離は4近傍の世界では真の残りコストとぴったり一致していました(違反 0/381セル)。

4近傍での比較 展開ノード 経路コスト 許容性の違反
ダイクストラ 352 34.00 —
マンハッタン 185 34.00 0 / 381
オクタイル 218 34.00 0 / 381

表:4近傍にすると立場が入れ替わる。マンハッタンは許容的なまま探索を1.9倍絞り、オクタイル(斜めを$\sqrt{2}$で数える=4近傍では控えめすぎる)は許容的だが甘いので218セルまで広がる。

同じ関数が、移動モデルを変えるだけで「反則」から「ぴったり正解」に変わる。 ヒューリスティックの正しさは関数そのものではなく、関数と移動モデルの組み合わせで決まります。これが「近傍と $h$ はセットで整合させる」の中身です。

重み付きA*:保証が消えることと、実際に壊れること

最後に、$h$ を意図的に盛る重み付きA*(Weighted A*:$h$ を $w$ 倍して探索を急がせる手法)です。$w>1$ なら許容性は壊れます。オクタイル×1.5 の時点で 348/381セルが違反、×2 では 375/381セルが違反。それでも実測はこうでした。

$w$ 展開ノード 経路コスト 許容性の違反 このマップの結果
1.0 113 26.97 0 / 381 最短(保証あり)
1.5 25 26.97 348 / 381 最短だが保証なし
2.0 24 26.97 375 / 381 最短だが保証なし
3.0 23 26.97 380 / 381 最短だが保証なし
5.0 19 35.94 380 / 381 +8.97 悪化

表:重みを上げていったときの実測。保証が消える境目($w>1$)と、実際に壊れる境目(このマップでは $w=3.25$)はまったく別の場所にある。

重み付きA*:wを上げたときの探索範囲と経路
図:左から $w=1.0$(113セル・26.97)、$w=2.0$(24セル・26.97)、$w=5.0$(19セル・35.94)。$w=5$ の橙の経路だけが、右上のぬかるみ帯(薄茶・重み3)を突っ切っている。

壊れ方の中身は、図の右端に出ています。 $w=5$ の経路は通過点19セルで、$w=1$ の23セルより短い。それでいてコストは35.94まで上がっています。ぬかるみ(重み3)を4セル突っ切って距離を稼ぎ、コストで損をしたのです。「見積もりを盛る」とは、A*に「多少高くつく道でもいいから前へ進め」と言うことなので、こういう負け方をします。

ここは書き方に気をつけるところなので、はっきり分けて書きます。

⚠️ 「$w>1$ にすると最適ではなくなる」は言いすぎです。 正しくは「最適性の保証が消える」。保証が消えることと、実際に最短が壊れることは別で、このマップでは $w=2$(375セルで許容性を破っている)でもコストは26.97のままでした。$w$ を1.0から5.0まで0.05刻みで走査すると、悪化が始まったのは $w=3.25$($w=3.20$ までは26.97)。逆に「$w=3.2$ までなら安全」も言えません。それは私のこのマップ・このゴールでの観測にすぎません。

壊れる境目を決めているのは、$h$ が無視している情報の量です。 それを確かめるために、ぬかるみの重みを 3.0 から 1.0 に落として(=地形の凹凸をなくして)同じ走査をしました。最適コストは26.97から24.63に下がり、$w=10$ まで走査しても最短は一度も壊れませんでした。オクタイル距離は地形コストを見ていないので、ぬかるみがあるマップでは $h$ が真値から大きく下振れし、その隙間に $w$ を掛けた瞬間に「ぬかるみを突っ切る」判断が生まれる。地形が平坦なら $h$ はほぼ真値なので、いくら盛っても行き先が変わりません。閾値はマップの性質そのものです。

理論の側から言えることも、実はあまり強くありません。$h$ が許容的なとき $w,h$ を使ったA*が返す経路のコストは、最適の $w$ 倍以内に収まると保証されます。$w=5$ なら $5 \times 26.97 = 134.85$ 以内。実測は 35.94(最適の1.33倍) でした。上限が緩すぎて、$w$ を決める役には立ちません。

「壊れる条件」を総当たりで探すスクリプト(sweeps.py・クリックで展開)

ゴールを全マスに動かす走査(マンハッタンの節の380通り)と、$w$ を刻んで上げる走査を1本にまとめたものです。

# 「壊れる条件」を総当たりで探す(ゴールを全マスに動かす/重み w を刻んで上げる)

# 自作モジュールの読み込み
import astar as A
import constant                              # ぬかるみの重みを書き換えるためモジュールごと読む
from constant import CONNECTIVITY, CELL_OBSTACLE, CELL_MUD
from dijkstra import dijkstra
from grid_map import GridMap, make_sample_map
from heuristics import h_manhattan, h_octile, make_weighted


def sweep_goals(heuristic) -> list:
    """
    ゴールを全マスに順番に置いて、そのヒューリスティックが最適を外すゴールを列挙する

    パラメータ
        heuristic: 検査したいヒューリスティック関数

    戻り値
        broken: (ゴール, 最適からの悪化量) のリスト
    """
    base   = make_sample_map()
    broken = []
    total  = 0

    for row in range(base.rows):
        for col in range(base.cols):
            goal = (row, col)
            if base.grid[goal] == CELL_OBSTACLE or goal == base.start:
                continue
            grid_map = GridMap(base.grid, base.start, goal)

            # 正解(最適コスト)はダイクストラで求める
            _, _, opt_cost = dijkstra(grid_map, CONNECTIVITY.EIGHT)
            if opt_cost == float("inf"):
                continue    # そのゴールへは到達できない
            total += 1

            A.heuristic = heuristic
            _, _, cost = A.astar(grid_map, CONNECTIVITY.EIGHT)
            if cost > opt_cost + 1e-9:
                broken.append((goal, round(cost - opt_cost, 3)))

    print(f"  到達可能なゴール {total} 通り中、最適を外した = {len(broken)} 通り  {broken}")
    return broken


def sweep_weight(w_max: float = 5.0, step: float = 0.05) -> float:
    """
    重み w を刻んで上げ、最適コストが崩れ始める w を返す

    パラメータ
        w_max: 走査する w の上限
        step: w の刻み幅

    戻り値
        threshold: 崩れ始めた w(w_max まで崩れなければ 0.0)
    """
    grid_map = make_sample_map()
    _, _, opt_cost = dijkstra(grid_map, CONNECTIVITY.EIGHT)

    w = 1.0
    while w <= w_max + 1e-9:
        A.heuristic = make_weighted(h_octile, w)
        _, explored, cost = A.astar(grid_map, CONNECTIVITY.EIGHT)
        if cost > opt_cost + 1e-9:
            print(f"  最適={opt_cost:.2f} → w={w:.2f} で cost={cost:.2f}(展開={len(explored)})に悪化")
            return w
        w = round(w + step, 2)

    print(f"  最適={opt_cost:.2f} → w={w_max:.2f} まで崩れなかった")
    return 0.0


def main() -> None:
    print("[1] マンハッタン距離(8近傍・非許容)が最適を外すゴールを全数探索")
    sweep_goals(h_manhattan)

    print("[2] オクタイル×w が崩れ始める w(ぬかるみの重み3.0=標準)")
    sweep_weight()

    print("[3] 同じ走査を、ぬかるみの重みを1.0(=地形の凹凸なし)にして再実行")
    constant.TERRAIN_COST[CELL_MUD] = 1.0    # 定数の辞書を書き換える(他モジュールも同じ辞書を見ている)
    sweep_weight(w_max=10.0)
    constant.TERRAIN_COST[CELL_MUD] = 3.0    # あとの実行に影響しないよう戻す


if __name__ == "__main__":
    main()
[1] マンハッタン距離(8近傍・非許容)が最適を外すゴールを全数探索
  到達可能なゴール 380 通り中、最適を外した = 2 通り  [((3, 18), 0.586), ((3, 19), 1.172)]
[2] オクタイル×w が崩れ始める w(ぬかるみの重み3.0=標準)
  最適=26.97 → w=3.25 で cost=35.94(展開=19)に悪化
[3] 同じ走査を、ぬかるみの重みを1.0(=地形の凹凸なし)にして再実行
  最適=24.63 → w=10.00 まで崩れなかった

だから実務での判断はこうなります。最短が絶対に必要か(例:燃料や電力が本当にぎりぎり)なら $w=1$ 以外に選択肢はありません。 数%の遠回りが許されるなら $w$ を上げ、壊れ始める値は自分のマップ群で挟み込んで測る。理論の上限ではなく実測で決める領域です。

使い分けの結論

移動モデルが決まれば $h$ は機械的に決まります。ここが本記事のいちばん持ち帰ってほしい表です。

移動モデル 使う $h$ 理由 本記事の実測(展開ノード)
4近傍(上下左右のみ) マンハッタン 1歩のコストが最低1で、斜めがないので過大評価にならない(違反0/381) 352 → 185
8近傍(斜めコスト$\sqrt{2}$) オクタイル 斜めを$\sqrt{2}$で正しく数える。許容的なまま最もタイト 369 → 113
連続空間・任意角度 ユークリッド どんな進み方でも直線距離を下回れない 369 → 138
速度優先(最短でなくてよい) オクタイル×$w$ 探索は激減する。壊れ始める $w$ は自分のマップで測る $w=2$ で 369 → 24

表:移動モデル別の選択指針。「$h$ が何を数えているか」と「実際に何ができるか」を一致させるのが唯一の規則。

要点を3つだけ。

  • $h$ の条件は許容性(真の残りコストを超えない)。守る範囲で大きいほど速い
  • 許容性は全セル数え上げで機械的に検査できる(ゴールから逆向きにダイクストラを1回)。1回の実行結果を見て「最短が出たからOK」と判断しない
  • $w>1$ は「保証が消える」。実際に壊れるかは別問題で、このマップでは $w=3.20$ まで壊れず、$w=3.25$ で壊れた

終わりに

②で残した宿題(「$h$ に何を入れるべきか」)の答え合わせでした。移動モデルに整合した最大の見積もりを選ぶ。速さは $h$ の大きさで、正しさは許容性で決まる。 この2つが分かっていれば大成功です。

経路生成シリーズ全体は、こちらのまとめから辿れます。

ロボット関連の記事はこちらのまとめにあります。

連載の本線は ③ポテンシャル法。格子を離れ、連続空間を「ゴールへの引力」と「障害物からの斥力」で進みます。今回のように「速いのに壊れている」が、あちらではもっと派手な形で出てきます。

役に立ったら いいね・ストック で応援いただけると、次回の励みになります。

0
0
0

Register as a new user and use Qiita more conveniently

  1. You get articles that match your needs
  2. You can efficiently read back useful information
  3. You can use dark theme
What you can do with signing up
0
0

Delete article

Deleted articles cannot be recovered.

Draft of this article would be also deleted.

Are you sure you want to delete this article?