🐥

pythonで道路ネットワーク分析

2024/12/22に公開

この記事はRICORA Advent Calender 2024の19日目の記事になります。

はじめに

学校のグループ開発で、pythonのOSMnxというパッケージを使用した道路ネットワークの経路探索を行ったため、基本的な使い方について紹介します。

OSMnxについて

OSMnxについて、OSMnxのドキュメントの冒頭部分google翻訳したものを引用します。
https://osmnx.readthedocs.io/en/stable/index.html

OSMnx は、OpenStreetMap から街路網やその他の地理空間フィーチャを簡単にダウンロード、モデル化、分析、視覚化できる Python パッケージです。1 行のコードで、歩行、運転、または自転車のネットワークをダウンロードしてモデル化し、分析して視覚化できます。都市のアメニティ/興味のあるポイント、建物のフットプリント、交通機関の停留所、標高データ、街路の方向、速度/移動時間、ルーティングも簡単に操作できます。

今回は、OpenStreetMapのデータを使用できるpythonのパッケージ、OSMnxを使用した以下のやり方について説明します。

  • 道路ネットワークをグラフとして取得
  • 道路ネットワークの描画
  • 道路ネットワーク上での経路探索

道路ネットワークの取得

OSMnxでは、graphモジュールの関数を使用することで、OpenStreetMapの道路データをグラフとして取得することができます。
グラフは、pythonのパッケージであるNetworkXという、グラフやネットワークを扱えるパッケージの多重有向グラフの型であるMultiDiGraphとして渡されます。

プログラム

import osmnx as ox

center = (35.9181569050815, 139.90789773749773) # 東京理科大学 野田キャンパスの緯度経度
g = ox.graph_from_point(center_point=center, dist=1000)

print(type(g))
# <class 'networkx.classes.multidigraph.MultiDiGraph'>

説明

まず、1行目のimportでOSMnxをインポートし、使えるようにします。
次に、graph_from_pointメソッドを使用してグラフを取得します。
graph_from_pointメソッドでは、引数center_pointにグラフを取得する範囲の中心となる点の緯度経度をタプルで指定し、distに取得するグラフの範囲を表す四角形の中心からの距離[m]で指定します。
center_pointに与える緯度経度は、GoogleMapで右クリックすることで、クリックした位置の緯度経度をコピーすることができます。

グラフの範囲指定は、中心点からの距離を指定する以外にも、

  • 長方形として範囲指定
  • 多角形として範囲指定
  • 中心点を住所で指定

などでできます。
詳しくは、以下のドキュメントを参照してください。
https://osmnx.readthedocs.io/en/stable/user-reference.html#module-osmnx.graph

グラフに含まれる情報の取得

OSMnxで取得できるグラフには、頂点や辺などに様々な情報が入ったものが渡されます。

プログラム

import osmnx as ox

# グラフを取得
center = (35.9181569050815, 139.90789773749773) # 東京理科大学 野田キャンパスの緯度経度
g = ox.graph_from_point(center_point=center, dist=1000)

# 頂点の情報
print(g.nodes)
# [278880597, 278880600, 278880601, 278880604, 278880607, 339577868, …, 12328163480]

# 頂点の詳細情報
print(g.nodes(data=True))
# [(278880597, {'y': 35.9172941, 'x': 139.9062309, 'street_count': 4}), (278880600, {'y': 35.9167047, 'x': 139.9052927, 'highway': 'crossing', 'street_count': 4}), …, (12328163480, {'y': 35.9201241, 'x': 139.9109177, 'street_count': 3})]

# 辺の情報
print(g.edges)
# [(278880597, 9295881005, 0), (278880597, 9295880884, 0), (278880597, 5290910420, 0), …, (12328163480, 5290904441, 0)]

# 辺の詳細情報
print(g.edges(data=True))
# [(278880597, 9295881005, {'osmid': 25584888, 'highway': 'unclassified', 'oneway': False, 'reversed': False, 'length': 9.078}), …
# (12328163480, 5290904441, {'osmid': [822225952, 1332559252], 'highway': 'service', 'oneway': False, 'reversed': False, 'length': 64.492, 'geometry': <LINESTRING (139.911 35.92, 139.911 35.92, 139.911 35.92)>})]

説明

頂点や辺の情報は、nodesやedgesで取得することができ、さらに詳細の情報は関数として呼び出して引数dataにTrueを指定することで取得できます。

頂点の情報としては、頂点に割り振られているIDが取得でき、詳細には、頂点の緯度経度や隣接している頂点の個数などが含まれます。
辺の情報としては、辺の始点と終点の頂点のIDが取得でき、詳細には、辺の長さやその辺が属する道路区分などが含まれます。

辺のhighwayタグに割り当てられている値については、以下のOpenStreetMapのwikiを確認してください。
https://wiki.openstreetmap.org/wiki/JA:Key:highway#値
また、主要な区分については、motorway, primary, secondalyなどがあり、これらは国ごとに割り当てが異なっています。日本では、motorwayは高速道路、primaryは国道、secondalyは二桁以下の県道というように割り当てられています。日本でのその他の割り当てやほかの国での割り当てについては、以下を参照してください。
https://wiki.openstreetmap.org/wiki/International_highway_classification_equivalence

グラフに含まれるデータをDataFrameに変換

convertモジュールを使用することで、取得した多重有向グラフを多重無向グラフにしたり、単純有向グラフにしたり、地理的データを扱うGeoPandasパッケージに含まれるDataFrameのような型GeoDataFrameに変換したりを行うことができます。
今回は、GeoDataFrameに変換して、分析を行います。

プログラム

import osmnx as ox

# グラフを取得
center = (35.9181569050815, 139.90789773749773) # 東京理科大学 野田キャンパスの緯度経度
g = ox.graph_from_point(center_point=center, dist=1000)

# 頂点のDataFrame
print(ox.graph_to_gdfs(g, nodes=True, edges=False))
"""
                     y           x  ...          highway                    geometry
osmid                               ...                                             
278880597    35.917294  139.906231  ...              NaN  POINT (139.90623 35.91729)
278880600    35.916705  139.905293  ...         crossing  POINT (139.90529 35.91670)
278880601    35.916687  139.905241  ...  traffic_signals  POINT (139.90524 35.91669)
278880604    35.914339  139.905492  ...              NaN  POINT (139.90549 35.91434)
278880607    35.911852  139.905066  ...              NaN  POINT (139.90507 35.91185)
...                ...         ...  ...              ...                         ...
12208885445  35.913591  139.899317  ...              NaN  POINT (139.89932 35.91359)
12270428266  35.915181  139.914740  ...              NaN  POINT (139.91474 35.91518)
12270428267  35.915220  139.914833  ...              NaN  POINT (139.91483 35.91522)
12328163479  35.920099  139.911092  ...              NaN  POINT (139.91109 35.92010)
12328163480  35.920124  139.910918  ...              NaN  POINT (139.91092 35.92012)

[2060 rows x 5 columns]
"""

# 辺のDataFrame
print(ox.graph_to_gdfs(g, nodes=False, edges=True))
"""
                                               osmid       highway  ...  width  ref
u           v           key                                         ...            
278880597   9295881005  0                   25584888  unclassified  ...    NaN  NaN
            9295880884  0                 1007605084       service  ...    NaN  NaN
            5290910420  0                 1337112287      tertiary  ...    NaN  NaN
            9310691628  0                 1337112287      tertiary  ...    NaN  NaN
278880600   278880601   0                   25584851  unclassified  ...    NaN  NaN
...                                              ...           ...  ...    ...  ...
12328163479 12328163479 0                 1332559253       service  ...    NaN  NaN
                        1                 1332559253       service  ...    NaN  NaN
12328163480 9298384785  0                  658351494       footway  ...    NaN  NaN
            12328163479 0                 1332559254       service  ...    NaN  NaN
            5290904441  0    [822225952, 1332559252]       service  ...    NaN  NaN

[5515 rows x 14 columns]
"""

説明

graph_to_gdfsメソッドを使用することで、グラフに含まれる情報をDataFrameのような形に変換することができます。引数nodes、edgesに引数を指定しない場合、頂点と辺が両方出力されて、間に改行がなくその部分が見ずらいので、プログラムでは分けて出力しています。

これを利用して、DataFrameとして操作することで、「グラフの中に含まれるhighwayタグの種類」などを知ることができます。

プログラム

import osmnx as ox
from osmnx import plot

# グラフを取得
center = (35.9181569050815, 139.90789773749773) # 東京理科大学 野田キャンパスの緯度経度
g = ox.graph_from_point(center_point=center, dist=1000)

# 辺のDataFrame
gdf = ox.graph_to_gdfs(g, nodes=False, edges=True)

# highwayの中のユニークな要素を出力
list_highway = gdf["highway"].tolist()
list_highway = list(map(str, list_highway))
print(set(list_highway))
# {'cycleway', 'pedestrian', 'tertiary', "['footway', 'service']", "['footway', 'steps']", 'unclassified', "['residential', 'footway', 'steps']", "['path', 'service']", "['residential', 'footway', 'service']", "['residential', 'footway']", "['residential', 'unclassified']", "['unclassified', 'cycleway']", "['pedestrian', 'footway']", "['unclassified', 'tertiary']", 'trunk', "['path', 'footway', 'steps']", "['residential', 'track']", "['track', 'path']", 'path', "['path', 'steps']", "['unclassified', 'service']", 'track', "['residential', 'path']", 'footway', 'steps', 'service', 'residential', "['residential', 'service']", "['track', 'path', 'service']"}

# highwayを1つ持っている中でユニークな要素を出力
list_highway_mono = gdf["highway"].tolist()
list_highway_mono = [e for e in list_highway_mono if type(e) == str]
print(set(list_highway_mono))
# {'track', 'cycleway', 'footway', 'steps', 'service', 'residential', 'path', 'pedestrian', 'tertiary', 'unclassified', 'trunk'}

説明

DataFrameのhighwayの列を取り出し、tolistを使用してリストに変換します。辺の中には複数のタグを与えられている辺があり、それらはリストの形で保管されています。setを使用して、重複した要素を削除してユニークな要素のみにしますが、リストを含んでいるとsetに変換できないため、先にリストの全ての要素を文字列型に変換しておきます。

ただ、リストの形が含まれていると見づらいため、下では文字列型ではない要素を事前に取り除いてsetに変換し、出力しています。

グラフの描画

取得した道路ネットワークの描画は、plotモジュールを使用することでできます。
今回はMatplotlibを使用して、描画を行います。

プログラム

import osmnx as ox

# グラフを取得
center = (35.9181569050815, 139.90789773749773) # 東京理科大学 野田キャンパスの緯度経度
g = ox.graph_from_point(center_point=center, dist=1000)

# グラフの描画
ox.plot_graph(g, save=True, filepath="road.png")

出力

説明

グラフの取得までは、前のプログラムと同じものです。
グラフの描画は、plot_graphメソッドを使用して行っています。描画自体は、描画するグラフを渡すだけでできます。サンプルプログラムでは、画像として出力するために、引数saveとfilepathを指定しています。
その他にも、

  • 頂点や辺の色の指定
  • 頂点の大きさ
  • 辺の太さ

などを指定することができます。「頂点や辺の色の指定」については、後の項目で説明します。

グラフ上の経路の描画

グラフにおける経路の取得は、routingモジュールを使用してできます。
今回は、ランダムに選んだ二つの頂点間の最短距離を求めて描画します。

プログラム

import random
import osmnx as ox

# グラフを取得
center = (35.9181569050815, 139.90789773749773) # 東京理科大学 野田キャンパスの緯度経度
g = ox.graph_from_point(center_point=center, dist=1000)

# ランダムに頂点を二つ選ぶ
nodes = list(g.nodes)
node_start = nodes[random.randint(0, len(nodes) - 1)]
node_goal = nodes[random.randint(0, len(nodes) - 1)]

# 最短距離の経路
route = ox.shortest_path(g, node_start, node_goal)
print(type(route))
# <class 'list'>
print(route)
# [9927298470, 9927298472, 9927298469, 1346855894, 2892639894,…, 6165206661]

# グラフと経路の描画
ox.plot_graph_route(g, route, save=True, filepath="road.png")

出力

説明

まず、nodesをリストに変換することで頂点のIDからなるリストを作り、その中からランダムに二つ取り出すことで頂点を二つ選びます。
最短経路は、shortest_pathメソッドに、グラフ、始点の頂点ID、終点の頂点IDを渡すことで取得することができます。経路は、頂点IDの列を表すint型のリストの形で渡されます。
グラフと経路の描画は、plot_graph_routeメソッドに、引数としてグラフと経路を渡すことでできます。また、引数node_colorやedge_colorに文字列型で色を指定することにより、頂点や辺の色を変えることができます。

頂点と辺に色を付けたグラフの描画

plotモジュールのメソッドには、引数を指定することで、頂点や辺ごとに色を付けて描画することができます。

頂点に色を付ける

プログラム

import random
import osmnx as ox
from osmnx import plot

# グラフを取得
center = (35.9181569050815, 139.90789773749773) # 東京理科大学 野田キャンパスの緯度経度
g = ox.graph_from_point(center_point=center, dist=1000)

# ランダムに頂点を二つ選ぶ
nodes = list(g.nodes)
node_start = nodes[random.randint(0, len(nodes) - 1)]
node_goal = nodes[random.randint(0, len(nodes) - 1)]

# 最短距離の経路
route = ox.shortest_path(g, node_start, node_goal)

# 経路から頂点を色付け
get_colornum = {n: 0 for n in g.nodes}
for i in range(len(route)):
    get_colornum[route[i]] = i + 1 # 経路に含まれる頂点は始点から1,2,…,それ以外は0
color_scale = ["white"] + plot.get_colors(len(route), cmap="hsv") # 色を配列で取得
colors = [color_scale[get_colornum[n]] for n in g.nodes] # get_colornumの番号で頂点に色を割り当て

# グラフの描画
ox.plot_graph(g, node_color=colors, save=True, filepath="road.png")

出力

説明

このプログラムでは、頂点に色を付けることによって、経路を表現しています。

頂点に割り当てる色はcolorsでリストの形で定義していて、「nodesのi番目のIDを持つ頂点にはcolorsのi番目の色で塗られる」という形になっています。そのcolorsをplot_graphの引数node_colorに渡すことで、色を反映させることができます。

色については、「red」「white」といった単語や「#FFFF00」といった色コードで指定することができます。今回は、経路をグラデーションで表示するため、get_colorsメソッドを使用して、カラーマップから等間隔な色を取得しています。
get_colorsメソッドでは、第1引数で等間隔でいくつの色を取得するかを指定していて、cmapではどのようなカラーマップを使用するかを指定しています。

使用できるカラーマップは、matplotlibのもので以下から確認できます。

https://matplotlib.org/stable/users/explain/colors/colormaps.html

経路を表示は、先のようにplot_graph_routeメソッドを使用することでもできますが、「経路が密集している」「進んでいる方向を知りたい」といった場合はこちらの方がわかりやすく表示できます。
また、plot_graph_routeメソッドでは、経路がちゃんとしている(隣接した頂点を移動し続けている)ものでないと描画ができないため、自分のプログラムで導いた経路のデバッグとしても使用できます。

辺に色を付ける

プログラム

import osmnx as ox

# グラフを取得
center = (35.9181569050815, 139.90789773749773) # 東京理科大学 野田キャンパスの緯度経度
g = ox.graph_from_point(center_point=center, dist=1000)

# highwayタグによって色付け
ec = ["red" if type(g.edges[e]["highway"]) == str and g.edges[e]["highway"] == "unclassified" else "white" for e in g.edges]

# グラフと経路の描画
ox.plot_graph(g, edge_color=ec, save=True, filepath="road.png")

出力

説明

頂点のときと同様に、edgesで繰り返すことによって辺の色を持つリストを作成します。今回は、辺の持つhighwayタグが「unclassified」なら赤、そうでないなら白で描画されるように作成しました。ただし、複数のタグが割り当てられている辺には文字列のリストが入っていて、文字列と判定ができないため、最初に型の判定を行っています。
作成したリストをplot_graphメソッドの引数edge_colorに渡すことによって反映させることができます。

これを使用することで、どの道路にどのタグが割り当てられているのかを見ることができます。

自分で経路探索をする

OSMnxを使用して取得したグラフ上で、ダイクストラ法を用いて最短経路を求めます。

プログラム

import heapq
import random
import osmnx as ox

def get_shortest_route(g, node_start, node_goal):
    # 頂点が探索済みであるか
    node_used = {node: False for node in list(g.nodes)}
    # 開始頂点から各頂点までの最短距離
    node_dist = {node: float("inf") for node in list(g.nodes)}
    node_dist[node_start] = 0
    # その頂点にたどり着く直前の頂点
    node_prev = {node: -1 for node in list(g.nodes)}
    # 頂点と距離のタプルを含むキュー(優先度付き)
    queue = [(0, node_start)]

    # 各頂点への最短距離を求める
    while len(queue) > 0:
        totaldist, node = heapq.heappop(queue)

        # 目的地に到達した
        if node == node_goal:
            break
        
        # すでに確定済みの頂点
        if node_used[node]:
            continue
        
        # 探索済みに変える
        node_used[node] = True
        # 隣接する頂点すべてに対して繰り返す
        for node_next in g.neighbors(node):
            dist = float(g.edges[node, node_next, 0]["length"])
            # node_nextまでの距離が今までの最小値よりも小さい場合はnode_distを更新しキューに追加する
            if totaldist + dist < node_dist[node_next]:
                node_dist[node_next] = totaldist + dist
                node_prev[node_next] = node
                heapq.heappush(queue, (node_dist[node_next], node_next))

    # スタートとゴールがつながっていない
    if(node_dist[node_goal] == float("inf")):
        return []
    
    # 現在の頂点(最初はゴール)
    node_cur = node_goal
    # 復元された経路
    route = [node_goal]
    # 最短経路を復元する
    while node_cur != node_start:
        # ひとつ前に戻って経路に追加
        node_cur = node_prev[node_cur]
        route.append(node_cur)
 
    # 目的地->出発地の順になっているため逆にする
    route.reverse()
    return route

# グラフを取得
center = (35.9181569050815, 139.90789773749773) # 東京理科大学 野田キャンパスの緯度経度
g = ox.graph_from_point(center_point=center, dist=1000)

# ランダムに頂点を二つ選ぶ
nodes = list(g.nodes)
node_start = nodes[random.randint(0, len(nodes) - 1)]
node_goal = nodes[random.randint(0, len(nodes) - 1)]

# 最短距離の経路
route = get_shortest_route(g, node_start, node_goal)

# # グラフと経路の描画
ox.plot_graph_route(g, route, save=True, filepath="road.png")

出力

説明

shortest_pathメソッドを使用して最短経路を取得していた部分を、作成したget_shortest_routeに変更しました。
get_shortest_routeでは、ダイクストラ法を使用した最短経路の探索とその復元を行っています。
経路探索に必要となる「隣接頂点の取得」は、neighborsメソッドで行うことができます。また、辺の長さについては、辺のlengthタグに入っているため、辺を指定してそのlengthを取り出しています。
多重有向グラフでは始点と終点を決定してもその頂点間に複数の辺が存在することが考えられるため、辺の指定は「始点の頂点」「終点の頂点」「そのような頂点を結ぶ辺のうち何番目か」の3つで行います。今回は一番最初の辺として3つ目の値は0としています。

まとめ

OSMnxを使用することでどんなことができるのかや使い方について、私が実際に使ったものを中心に紹介しました。

リンク

OSMnxのドキュメント
https://osmnx.readthedocs.io/en/stable/index.html

OSMnxで得られるNetworkXの多重有向グラフ(MultiDiGraph)のドキュメント
https://networkx.org/documentation/stable/reference/classes/multidigraph.html

OpenStreetMapのWikiの日本語ページ
https://wiki.openstreetmap.org/wiki/JA:Main_Page

Discussion