# -*- coding: utf-8 -*- """L79 河川トンネル基本情報 単独 3 研究例分析 — 広島県管理 3 河川トンネル放水路 / 3 水系 / 2 規模クラス を 3 角度で深掘り カバー宣言: 本記事は DoBoX のシリーズ「河川トンネル基本情報」 1 件 (dataset_id = 1270) を 単独で取り上げ、 広島県が管理する 河川トンネル (= 治水トンネル放水路) 3 本を 3 つの独立した研究角度 (RQ1 / RQ2 / RQ3) で並列に分析する。 「県の隠れた治水主役 — 都市河川を地中で逃がす装置」 を、 構造・治水機能・道路トンネルとの用途差の 3 軸から立体的に描く。 「河川トンネル」 とは: 河川 (の一部または分派) を地中の管路で導く土木構造物。 本データの 3 本はすべて放水路 (= バイパス放流路) であり、 本川の洪水時に支川や本川下流部を地中ショートカットさせ、 下流の人家・市街地への洪水到達を抑える治水トンネルとして機能する。 河川法 (1964) 第 6 条 第 2 項 が規定する河川管理施設であり、 道路トンネル (道路法管轄) とは制度・目的・構造がまったく異なる。 本データは広島県が管理する 3 本のみ収録 (= 県管理。 国管理の太田川 放水路本線等は含まれない)。 本記事は L67 (道路トンネル 157 件単独) と厳密に区別: L67 = 道路法管轄、 道路を山岳貫通で結ぶ。 平均延長 314m。 L79 = 河川法管轄、 河川を地中バイパスで逃がす。 全 3 本平均 1,173m。 両者は「トンネル」 という共通形態を持つが、 用途 (交通 vs 治水)・規模・整備の希少性が異なる。 本記事 RQ3 で「道路 157 vs 河川 3」 = 1 : 52 の件数比 + 規模比較を実証。 本記事は L4-L11 (河川浸水想定区域シリーズ) と厳密に区別: L4-L11 = 河川浸水想定区域の結果分析 (= 浸水深・浸水時間・複合災害)。 L79 = 浸水想定区域を引き起こす本川を地中で逃がす対策装置側の分析。 RQ2 で「3 本の各トンネル下流端から半径 2km 圏内の浸水想定区域面積」 を計算し、 「治水トンネルがどれだけの恩恵範囲を持つか」 を量化する。 研究の問い (3 RQ): RQ1 (主研究): 広島県の河川トンネルの構造 — 規模・水系・地理分布はどう描けるか? 3 本 × 11 列 (名称 / 事務所 / 上下流所在地・緯度経度 / 海岸河川名 / 延長 / 流下能力) を単位流下能力 (m³/s/m)などの本記事独自指標も 含めて 9 軸で集計し、 「県内河川トンネルの物理形状」 を初めて系統的に 記述する。 H1 = 全 3 本が放水路仮説、 H2 = 延長は400-2000m 帯に集中し道路トンネル中央値 (206m) より長い仮説。 RQ2 (副研究 1): トンネルの治水機能 — 浸水想定区域への恩恵範囲 はどこまで届くか? 各トンネルの上流端 (取水) と下流端 (放流) の2 点 GeoDataFrameを 生成し、 下流端から半径 2 km 圏内のL72 既キャッシュ浸水想定区域 (浸水深最大値ポリゴン, 613 件 / 829.7 km²) との重なり面積を計算。 「下流恩恵」 = 下流端 2 km 圏内の浸水想定区域面積 (km²) と本記事独自 定義し、 3 本をランキング。 H3 = 八幡川 (広島市直下流) が 下流恩恵 県内最大仮説、 H4 = 流下能力 (m³/s) と下流恩恵範囲が 正の相関仮説。 RQ3 (副研究 2): L67 道路トンネル 157 件との用途差 — 道路 vs 河川 はどう現れるか? L67 既知 道路トンネル 157 件 vs L79 河川トンネル 3 件の件数比 52:1、 平均延長、 整備の希少性、 制度的位置付け (道路法 vs 河川法) を比較。 「道路 = 山岳貫通のための交通インフラ、 河川 = 都市内の隠れた治水主役」 の用途差の二極構造を実証する。 H5 = 河川トンネル平均延長は 道路トンネル平均の3 倍以上仮説。 仮説 (5): H1 (RQ1, 全件 放水路): 河川トンネル 3 本は全件 が放水路 (= 本川バイパス治水トンネル) であり、 「導水トンネル」 「砂防トンネル」 等の 非バイパス用途は 0 件。 これは「広島県の河川トンネル整備史 = 都市河川の 洪水対策史」 という制度的選択の物理的証拠仮説。 H2 (RQ1, 延長 400-2000m 帯): 延長は400-2000m 帯に 3/3 = 100% が集中。 道路トンネルでは延長中央値 206m に対し、 河川トンネルは中央値が 4 倍以上。 これは「河川トンネルは大規模な工事を経て例外的に整備される」 という希少性の量的証拠仮説。 H3 (RQ2, 八幡川 下流恩恵 最大): 八幡川トンネル放水路は広島市直下流で 下流端 2 km 圏内の浸水想定区域面積 (= 下流恩恵) が県内 3 本中 最大。 これは「都市直下流の治水トンネル = 最大の治水投資効果」 という 立地戦略の物理的証拠仮説。 H4 (RQ2, 流下能力 ↔ 下流恩恵 正相関): 流下能力 (m³/s) と下流恩恵範囲 (km²) の Pearson r >= 0.5。 大規模トンネルほど大きな治水効果を持つという 直観の量的検証仮説。 H5 (RQ3, 河川 vs 道路 規模比 ≥ 3): 河川トンネル平均延長 / 道路トンネル 平均延長 (= L67 既知 314m) ≥ 3。 これは「河川トンネル = 大規模 例外整備、 道路トンネル = 中規模日常整備」 の規模二極構造の量的検証仮説。 要件 S 準拠 (1 分以内完走): - 全データ 3 件 / 11 列 → 超軽量 - POINT geometry × 2 系列 (上流端/下流端) で 6 点のみ → sjoin 高速 - L72 浸水想定区域キャッシュ (613 件) を再利用 → 大規模 IO 不要 - 本スクリプト 1 本で完結 (~10 秒目標) 要件 T 準拠 (位置情報あり = 地図必須): - RQ1: 県全域 配置マップ (上流端〜下流端を線で結び、 流下能力でバブル) - RQ2: 浸水想定区域 + 河川トンネル + 下流恩恵 2 km バッファの重ね合わせ - RQ3: L67 道路トンネル 157 件 (背景) + L79 河川トンネル 3 本 (前景) の二極構造マップ 要件 Q 準拠: 図 7 / 表 11+ (3 RQ × 多角度: 構造 / 治水 / 道路対比) データ仕様: - dataset 1270: 河川トンネル基本情報 (CSV) - 形式: CSV (UTF-8 BOM), 3 行 × 11 列 - 列: 河川トンネル名称, 事務所名, 上流端 所在地, 下流端 所在地, 上流端 緯度(10進), 上流端 経度(10進), 下流端 緯度(10進), 下流端 経度(10進), 海岸・河川名, 延長(m), 流下能力(m3/s) - 全件 河川法上の河川管理施設 (放水路)区分 - 全件 管理者 = 広島県 - ライセンス: クリエイティブ・コモンズ表示 (CC-BY) メモリ対策: Figure ごとに plt.close('all') で確実に解放。 """ from __future__ import annotations import sys, time from pathlib import Path sys.path.insert(0, str(Path(__file__).parent)) from _common import (ROOT, ASSETS, LESSONS, render_lesson, code, figure, ensure_dataset) import numpy as np import pandas as pd import geopandas as gpd from shapely.geometry import Point, LineString import matplotlib matplotlib.use("Agg") import matplotlib.pyplot as plt from matplotlib.patches import Patch from matplotlib.lines import Line2D plt.rcParams["font.family"] = "Yu Gothic" plt.rcParams["axes.unicode_minus"] = False t_all = time.time() print("=== L79 河川トンネル基本情報 単独 3 研究例分析 ===", flush=True) # ============================================================================= # 0. 定数・パス # ============================================================================= TARGET_CRS = "EPSG:6671" # JGD2011 平面直角第 III 系 (m 単位) DATASET_ID = 1270 RESOURCE_ID = 39771 DATA_DIR = ROOT / "data" / "extras" / "L79_river_tunnels" DATA_DIR.mkdir(parents=True, exist_ok=True) SRC_CSV = DATA_DIR / "river_tunnel_basic.csv" # 既存キャッシュの再利用 ADMIN_GPKG = ROOT / "data" / "extras" / "L44_storm_surge" / "_cache" / "admin_diss.gpkg" FLOOD_GPKG = ROOT / "data" / "extras" / "L72_emergency_road" / "_cache_flood_max.gpkg" ROAD_TUNNEL_CSV = ROOT / "data" / "extras" / "tunnel_basic.csv" # 下流恩恵バッファ (本記事独自) DOWNSTREAM_BUFFER_M = 2000 # 下流端から 2 km 圏内 = 治水恩恵範囲 # 規模クラス境界 (本記事独自, 河川トンネル用) def length_class(L): if pd.isna(L): return "(不明)" if L >= 1500: return "大 (>=1500m)" if L >= 800: return "中 (800-1500m)" return "小 (<800m)" # 流下能力クラス (本記事独自) def discharge_class(Q): if pd.isna(Q): return "(不明)" if Q >= 100: return "大 (>=100 m³/s)" if Q >= 30: return "中 (30-100 m³/s)" return "小 (<30 m³/s)" # 水系の表示色 (S52 / L78 と統一思想) WS_COLORS = { "太田川": "#0969da", "沼田川": "#1f883d", "江の川": "#cf222e", } # 事務所別 色 OFFICE_COLORS = { "西部建設事務所": "#0969da", "西部建設事務所東広島支所": "#1f883d", "北部建設事務所": "#cf222e", } # 道路トンネル比較定数 (L67 既知) ROAD_TUNNEL_MEAN_LEN = 314.5 ROAD_TUNNEL_MEDIAN_LEN = 206.0 ROAD_TUNNEL_COUNT = 157 # ============================================================================= # 1. データ取得 (なければ DoBoX から自動取得) # ============================================================================= print("\n[1] データ取得", flush=True) t1 = time.time() ensure_dataset(SRC_CSV, dataset_id=DATASET_ID, resource_id=RESOURCE_ID, label="河川トンネル基本情報") df = pd.read_csv(SRC_CSV, encoding="utf-8-sig") df = df.dropna(subset=["河川トンネル名称"]).reset_index(drop=True) df = df[df["河川トンネル名称"].astype(str).str.strip() != ""].reset_index(drop=True) print(f" 河川トンネル件数: {len(df)}", flush=True) assert len(df) == 3, f"想定 3 件、 実際 {len(df)} 件" # 数値列の型変換 for c in ["上流端　緯度(10進数)", "上流端　経度(10進数)", "下流端　緯度(10進数)", "下流端　経度(10進数)", "延長(m)", "流下能力(m3/s)"]: df[c] = pd.to_numeric(df[c], errors="coerce") # 列名を短縮 (扱いやすく) df = df.rename(columns={ "上流端　所在地": "上流端_所在地", "下流端　所在地": "下流端_所在地", "上流端　緯度(10進数)": "上流緯度", "上流端　経度(10進数)": "上流経度", "下流端　緯度(10進数)": "下流緯度", "下流端　経度(10進数)": "下流経度", "海岸・河川名": "本川_放水路", "延長(m)": "延長_m", "流下能力(m3/s)": "流下能力_m3s", }) # 派生列 df["水系"] = df["本川_放水路"].astype(str).str.split("　").str[0] df["放水路名"] = df["本川_放水路"].astype(str).str.split("　").str[1].fillna("") df["用途"] = df["河川トンネル名称"].apply( lambda s: "放水路" if "放水路" in str(s) else "その他") df["延長クラス"] = df["延長_m"].map(length_class) df["流下能力クラス"] = df["流下能力_m3s"].map(discharge_class) # 本記事独自指標 1: 単位延長あたり流下能力 (m³/s/m) = トンネル断面の # 「水を流す力 / 1m あたり」 効率 df["単位流下効率_m3s_per_m"] = (df["流下能力_m3s"] / df["延長_m"]).round(4) # 本記事独自指標 2: 「治水ポテンシャル指数」 = 流下能力 × log(延長) # 「短くて細い vs 長くて太い」 を 1 軸で比較する合成指標 df["治水ポテンシャル指数"] = (df["流下能力_m3s"] * np.log10(df["延長_m"])).round(2) n_tunnels = len(df) print(f" ({time.time()-t1:.2f}s)", flush=True) # ============================================================================= # 2. GeoDataFrame 化 (上流端 + 下流端 + LineString) # ============================================================================= print("\n[2] GeoDataFrame 化", flush=True) t2 = time.time() # 上流端 POINT gdf_up = gpd.GeoDataFrame( df.assign(端点="上流端").copy(), geometry=[Point(lon, lat) for lon, lat in zip(df["上流経度"], df["上流緯度"])], crs="EPSG:4326").to_crs(TARGET_CRS) # 下流端 POINT gdf_dn = gpd.GeoDataFrame( df.assign(端点="下流端").copy(), geometry=[Point(lon, lat) for lon, lat in zip(df["下流経度"], df["下流緯度"])], crs="EPSG:4326").to_crs(TARGET_CRS) # トンネル経路 LineString (上流→下流の直線近似) gdf_line = gpd.GeoDataFrame( df.copy(), geometry=[LineString([Point(ux, uy), Point(dx, dy)]) for ux, uy, dx, dy in zip(df["上流経度"], df["上流緯度"], df["下流経度"], df["下流緯度"])], crs="EPSG:4326").to_crs(TARGET_CRS) # 直線近似距離 (m) を追加 df["直線距離_m"] = gdf_line.geometry.length.round(1).tolist() df["迂回率"] = (df["延長_m"] / df["直線距離_m"]).round(3) gdf_line["直線距離_m"] = df["直線距離_m"].values gdf_line["迂回率"] = df["迂回率"].values print(f" ({time.time()-t2:.2f}s)", flush=True) # ============================================================================= # 3. 行政界 sjoin (上流/下流 各端の市町判定) # ============================================================================= print("\n[3] 行政界 sjoin", flush=True) t3 = time.time() admin = gpd.read_file(ADMIN_GPKG).to_crs(TARGET_CRS) CITY_NAME = { 101: "広島市中区", 102: "広島市東区", 103: "広島市南区", 104: "広島市西区", 105: "広島市安佐南区", 106: "広島市安佐北区", 107: "広島市安芸区", 108: "広島市佐伯区", 202: "呉市", 203: "竹原市", 204: "三原市", 205: "尾道市", 207: "福山市", 208: "府中市", 209: "三次市", 210: "庄原市", 211: "大竹市", 212: "東広島市", 213: "廿日市市", 214: "安芸高田市", 215: "江田島市", 302: "府中町", 304: "海田町", 307: "熊野町", 309: "坂町", 368: "安芸太田町", 369: "安芸太田町", 412: "北広島町", 462: "世羅町", 545: "神石高原町", 543: "神石高原町", } admin["市町名"] = admin["CITY_CD"].map(CITY_NAME).fillna( admin["CITY_CD"].astype(str)) def _sjoin_city(gdf_pt, idcol): gdf_pt[idcol] = [f"L79_{i:02d}" for i in range(len(gdf_pt))] j = gpd.sjoin(gdf_pt[[idcol, "geometry"]], admin[["市町名", "geometry"]], how="left", predicate="within").drop_duplicates(idcol) return gdf_pt[idcol].map(j.set_index(idcol)["市町名"]).fillna("不明") df["上流端_市町"] = _sjoin_city(gdf_up.copy(), "seg_id_up").values df["下流端_市町"] = _sjoin_city(gdf_dn.copy(), "seg_id_dn").values print(f" 上流端 市町: {df['上流端_市町'].tolist()}", flush=True) print(f" 下流端 市町: {df['下流端_市町'].tolist()}", flush=True) print(f" ({time.time()-t3:.2f}s)", flush=True) # ============================================================================= # 4. RQ1 集計: 構造 (規模・水系・地理分布) # ============================================================================= print("\n[4] RQ1 集計 — 構造", flush=True) t4 = time.time() # (1) 用途別 T_use = (df.groupby("用途") .agg(本数=("河川トンネル名称", "count"), 延長平均=("延長_m", "mean"), 延長合計=("延長_m", "sum"), 流下能力平均=("流下能力_m3s", "mean")) .round(1).reset_index()) T_use["シェア_%"] = (T_use["本数"] / n_tunnels * 100).round(1) T_use = T_use.sort_values("本数", ascending=False).reset_index(drop=True) n_bypass = int(T_use.loc[T_use["用途"] == "放水路", "本数"].iloc[0]) \ if (T_use["用途"] == "放水路").any() else 0 share_bypass = round(100 * n_bypass / n_tunnels, 1) h1_ok = share_bypass == 100.0 # (2) 延長クラス別 T_len = (df.groupby("延長クラス") .agg(本数=("河川トンネル名称", "count"), 延長平均=("延長_m", "mean"), 流下能力平均=("流下能力_m3s", "mean"), 代表トンネル=("河川トンネル名称", lambda s: " / ".join(s.tolist()))) .round(1).reset_index()) T_len["シェア_%"] = (T_len["本数"] / n_tunnels * 100).round(1) len_order = ["大 (>=1500m)", "中 (800-1500m)", "小 (<800m)"] T_len["順"] = T_len["延長クラス"].apply( lambda x: len_order.index(x) if x in len_order else 99) T_len = T_len.sort_values("順").drop(columns="順").reset_index(drop=True) n_400_2000 = int(((df["延長_m"] >= 400) & (df["延長_m"] <= 2000)).sum()) h2_ok = n_400_2000 == n_tunnels # (3) 流下能力クラス別 T_q = (df.groupby("流下能力クラス") .agg(本数=("河川トンネル名称", "count"), 延長平均=("延長_m", "mean"), 流下能力平均=("流下能力_m3s", "mean")) .round(1).reset_index()) T_q["シェア_%"] = (T_q["本数"] / n_tunnels * 100).round(1) q_order = ["大 (>=100 m³/s)", "中 (30-100 m³/s)", "小 (<30 m³/s)"] T_q["順"] = T_q["流下能力クラス"].apply( lambda x: q_order.index(x) if x in q_order else 99) T_q = T_q.sort_values("順").drop(columns="順").reset_index(drop=True) # (4) 水系別 T_water = (df.groupby("水系") .agg(本数=("河川トンネル名称", "count"), 延長合計=("延長_m", "sum"), 流下能力合計=("流下能力_m3s", "sum"), 代表トンネル=("河川トンネル名称", lambda s: " / ".join(s.tolist()))) .round(1).reset_index()) T_water["シェア_%"] = (T_water["本数"] / n_tunnels * 100).round(1) T_water = T_water.sort_values("本数", ascending=False).reset_index(drop=True) # (5) 事務所別 T_office = (df.groupby("事務所名") .agg(本数=("河川トンネル名称", "count"), 延長合計=("延長_m", "sum"), 流下能力合計=("流下能力_m3s", "sum"), 代表トンネル=("河川トンネル名称", lambda s: " / ".join(s.tolist()))) .round(1).reset_index()) T_office = T_office.sort_values("本数", ascending=False).reset_index(drop=True) # (6) 個体台帳 (RQ1 主役) T_overview = df[["河川トンネル名称", "事務所名", "水系", "放水路名", "上流端_所在地", "下流端_所在地", "上流端_市町", "下流端_市町", "延長_m", "直線距離_m", "迂回率", "流下能力_m3s", "単位流下効率_m3s_per_m", "治水ポテンシャル指数", "延長クラス", "流下能力クラス"]].copy() T_overview = T_overview.sort_values("延長_m", ascending=False).reset_index(drop=True) print(f" 全 {n_tunnels} 本 / 用途 = 放水路 {n_bypass} ({share_bypass}%)", flush=True) print(f" 延長 400-2000m: {n_400_2000}/{n_tunnels}", flush=True) print(f" H1 (全件 放水路): {h1_ok}", flush=True) print(f" H2 (延長 400-2000m 帯集中): {h2_ok}", flush=True) print(f" ({time.time()-t4:.2f}s)", flush=True) # ============================================================================= # 5. RQ2 集計: 治水機能 (浸水想定区域との関係) # ============================================================================= print("\n[5] RQ2 集計 — 治水機能", flush=True) t5 = time.time() # (1) L72 浸水想定区域キャッシュを読込み # 元データは 4 億超頂点の超詳細ポリゴン → 本記事用に # pre-clip + simplify したサブセットを 1 度だけキャッシュ作成 L79_FLOOD_CACHE = DATA_DIR / "_flood_clipped_for_L79.gpkg" # まずトンネル下流端 4 km bbox (= 2 km バッファに余裕付け) を作って # それで先に flood を pre-clip してから処理する gdf_dn_bbox_geom = gpd.GeoDataFrame( geometry=[Point(lon, lat) for lon, lat in zip(df["下流経度"], df["下流緯度"])], crs="EPSG:4326").to_crs(TARGET_CRS).geometry.buffer( DOWNSTREAM_BUFFER_M + 500).union_all().envelope if L79_FLOOD_CACHE.exists(): print(f" flood cache ヒット: {L79_FLOOD_CACHE.name}", flush=True) flood = gpd.read_file(L79_FLOOD_CACHE).to_crs(TARGET_CRS) else: print(" flood cache 未作成 → 初回ビルド (~1-2 分)", flush=True) raw = gpd.read_file(FLOOD_GPKG).to_crs(TARGET_CRS) # 全 bbox で先に絞り込み (sindex 利用) → 大半を捨てる cand_idx = list(raw.sindex.query(gdf_dn_bbox_geom, predicate="intersects")) raw = raw.iloc[cand_idx].copy() print(f" pre-clip: {len(raw)} 件", flush=True) # 2D 化 + make_valid + simplify (10m 許容) で頂点数を 1/100 程度に from shapely import make_valid from shapely.ops import transform as _shp_transform raw["geometry"] = raw.geometry.apply( lambda g: _shp_transform(lambda x, y, z=None: (x, y), g)) raw["geometry"] = raw.geometry.apply(make_valid) raw["geometry"] = raw.geometry.simplify(10.0, preserve_topology=True) raw["geometry"] = raw.geometry.buffer(0) # bbox で 2 度クリップ (シンプル化で範囲微変更したため) raw = raw[raw.geometry.intersects(gdf_dn_bbox_geom)].copy() raw.to_file(L79_FLOOD_CACHE, driver="GPKG") flood = raw print(f" cache 書込: {len(flood)} 件", flush=True) flood_total_km2 = round(flood.geometry.area.sum() / 1e6, 2) print(f" 浸水想定区域 (clip後): {len(flood)} 件 / " f"{flood_total_km2} km² (2km バッファ近傍)", flush=True) # (2) 各トンネルの下流端から半径 2 km バッファを生成 gdf_dn["seg_id"] = [f"L79_{i:02d}" for i in range(len(gdf_dn))] gdf_dn["buffer_2km"] = gdf_dn.geometry.buffer(DOWNSTREAM_BUFFER_M) gdf_buf = gpd.GeoDataFrame( gdf_dn[["seg_id", "河川トンネル名称", "水系", "下流端_所在地", "流下能力_m3s", "延長_m"]].copy(), geometry=gdf_dn["buffer_2km"], crs=TARGET_CRS) # (3) 各バッファ内の浸水想定区域面積を計算 # 速度対策: 全 union は重く topology error も出やすい。 候補のみ抽出 → # 候補 1 件ずつ intersection し、 次に candidates 同士の重なりを # 「intersection 結果集合の dissolve」 で取り除く 2 段階方式。 sindex = flood.sindex benefit_rows = [] for _, r in gdf_buf.iterrows(): cand_idx = list(sindex.query(r.geometry, predicate="intersects")) benefit_area = 0.0 if cand_idx: # 各候補と バッファ の intersection を個別計算 → list に cand = flood.iloc[cand_idx].copy() inters = cand.geometry.intersection(r.geometry) # 候補間の重なりを除くため dissolve (= union) # ここで対象は バッファ範囲内に切り抜かれた小さな小片 (= 軽量) # なので union_all が現実的時間で完了する valid = inters[~inters.is_empty] if len(valid): try: u = valid.union_all() benefit_area = float(u.area) except Exception: # union が失敗した場合は単純合計 (重なりは過大評価) benefit_area = float(valid.area.sum()) benefit_km2 = benefit_area / 1e6 buf_area = r.geometry.area / 1e6 benefit_rows.append({ "河川トンネル名称": r["河川トンネル名称"], "水系": r["水系"], "下流端_所在地": r["下流端_所在地"], "流下能力_m3s": r["流下能力_m3s"], "延長_m": r["延長_m"], "下流バッファ面積_km2": round(buf_area, 2), "下流恩恵_km2": round(benefit_km2, 3), "下流恩恵率_%": round(100 * benefit_km2 / buf_area, 2), }) T_benefit = pd.DataFrame(benefit_rows).sort_values( "下流恩恵_km2", ascending=False).reset_index(drop=True) T_benefit["順位"] = np.arange(1, len(T_benefit) + 1) # (4) 八幡川の下流恩恵ランク yahata_rank = int(T_benefit.loc[ T_benefit["河川トンネル名称"].str.contains("八幡川"), "順位"].iloc[0]) yahata_benefit = float(T_benefit.loc[ T_benefit["河川トンネル名称"].str.contains("八幡川"), "下流恩恵_km2"].iloc[0]) h3_ok = yahata_rank == 1 # (5) 流下能力 vs 下流恩恵 の Pearson 相関 q_arr = T_benefit["流下能力_m3s"].astype(float).values b_arr = T_benefit["下流恩恵_km2"].astype(float).values if len(q_arr) >= 2 and np.std(q_arr) > 0 and np.std(b_arr) > 0: r_qb = float(np.corrcoef(q_arr, b_arr)[0, 1]) else: r_qb = float("nan") h4_ok = (not np.isnan(r_qb)) and r_qb >= 0.5 # (6) RQ2 の集計サマリ T_rq2_stats = pd.DataFrame([ {"統計": "下流恩恵 平均 (km²)", "値": round(T_benefit["下流恩恵_km2"].mean(), 3)}, {"統計": "下流恩恵 中央 (km²)", "値": round(T_benefit["下流恩恵_km2"].median(), 3)}, {"統計": "下流恩恵 最大 (km²)", "値": round(T_benefit["下流恩恵_km2"].max(), 3)}, {"統計": "下流恩恵 最小 (km²)", "値": round(T_benefit["下流恩恵_km2"].min(), 3)}, {"統計": "下流恩恵率 平均 (%)", "値": round(T_benefit["下流恩恵率_%"].mean(), 2)}, {"統計": "流下能力 vs 下流恩恵 Pearson r", "値": round(r_qb, 3)}, {"統計": "下流バッファ半径 (km)", "値": DOWNSTREAM_BUFFER_M / 1000.0}, {"統計": "解析対象 浸水想定区域 (km²)", "値": flood_total_km2}, ]) print(f" 八幡川 下流恩恵ランク: {yahata_rank} 位 ({yahata_benefit:.3f} km²)", flush=True) print(f" 流下能力 vs 下流恩恵 r = {r_qb:.3f}", flush=True) print(f" H3 (八幡川 下流恩恵 1 位): {h3_ok}", flush=True) print(f" H4 (流下能力 ↔ 下流恩恵 r >= 0.5): {h4_ok}", flush=True) print(f" ({time.time()-t5:.2f}s)", flush=True) # ============================================================================= # 6. RQ3 集計: L67 道路トンネルとの用途差対比 # ============================================================================= print("\n[6] RQ3 集計 — 道路 vs 河川対比", flush=True) t6 = time.time() # (1) L67 道路トンネルデータ読込み (地図プロット用 + 統計用) road = pd.read_csv(ROAD_TUNNEL_CSV, encoding="utf-8-sig") road = road.dropna(subset=["施設名", "緯度（10進数）", "経度（10進数）"]).copy() # 緯度経度の異常値 (緯度 > 50) は経度と入れ替え (L67 と同じ補正) swap_mask = road["緯度（10進数）"] > 50 road.loc[swap_mask, ["緯度（10進数）", "経度（10進数）"]] = \ road.loc[swap_mask, ["経度（10進数）", "緯度（10進数）"]].values road["延長(m)"] = pd.to_numeric(road["延長(m)"], errors="coerce") road["建設年度"] = pd.to_numeric(road["建設年度"], errors="coerce") road = road.dropna(subset=["延長(m)"]).reset_index(drop=True) road_total = len(road) road_mean_len = float(road["延長(m)"].mean()) road_median_len = float(road["延長(m)"].median()) road_year_median = float(road["建設年度"].median()) river_mean_len = float(df["延長_m"].mean()) river_median_len = float(df["延長_m"].median()) river_total = n_tunnels # (2) 件数比 + 規模比較 ratio_count = round(road_total / river_total, 1) ratio_mean_len = round(river_mean_len / road_mean_len, 2) ratio_median_len = round(river_median_len / road_median_len, 2) h5_ok = ratio_mean_len >= 3.0 T_compare = pd.DataFrame([ {"指標": "件数 (本)", "道路トンネル (L67)": road_total, "河川トンネル (L79)": river_total, "比率 (河川/道路 倍)": round(river_total / road_total, 4), "備考": f"件数比 {ratio_count}:1 (= 道路の希少性)"}, {"指標": "平均延長 (m)", "道路トンネル (L67)": round(road_mean_len, 1), "河川トンネル (L79)": round(river_mean_len, 1), "比率 (河川/道路 倍)": ratio_mean_len, "備考": "河川トンネルは道路の" + f"{ratio_mean_len:.1f} 倍長い"}, {"指標": "中央値延長 (m)", "道路トンネル (L67)": round(road_median_len, 1), "河川トンネル (L79)": round(river_median_len, 1), "比率 (河川/道路 倍)": ratio_median_len, "備考": "中央値でも河川が大きく上回る"}, {"指標": "最長 (m)", "道路トンネル (L67)": round(road["延長(m)"].max(), 1), "河川トンネル (L79)": round(df["延長_m"].max(), 1), "比率 (河川/道路 倍)": round(df["延長_m"].max() / road["延長(m)"].max(), 3), "備考": "道路の最長トンネルがやや上回る"}, {"指標": "最短 (m)", "道路トンネル (L67)": round(road["延長(m)"].min(), 1), "河川トンネル (L79)": round(df["延長_m"].min(), 1), "比率 (河川/道路 倍)": round(df["延長_m"].min() / road["延長(m)"].min(), 2), "備考": "河川の最短でも道路 100m 級と同等"}, ]) # (3) 用途・制度的位置付け比較 T_purpose = pd.DataFrame([ {"観点": "管轄法", "道路トンネル (L67)": "道路法 (1952)", "河川トンネル (L79)": "河川法 (1964)"}, {"観点": "主目的", "道路トンネル (L67)": "交通 (山岳貫通の時間短縮)", "河川トンネル (L79)": "治水 (洪水バイパス放流)"}, {"観点": "管理者", "道路トンネル (L67)": "県 (国道 + 県道) 道路河川管理課", "河川トンネル (L79)": "県 河川課 (建設事務所)"}, {"観点": "点検制度", "道路トンネル (L67)": "5 年に 1 回 (道路法施行規則 改正 2014)", "河川トンネル (L79)": "河川管理施設点検要領 (河川法施行令)"}, {"観点": "通常時の状態", "道路トンネル (L67)": "常時供用 (車両通行)", "河川トンネル (L79)": "通常時 = 待機 / 洪水時のみ通水"}, {"観点": "整備の希少性", "道路トンネル (L67)": "中規模・多数 (157 本)", "河川トンネル (L79)": "大規模・例外 (3 本のみ)"}, {"観点": "幾何学的形状", "道路トンネル (L67)": "馬蹄形断面 + 換気・照明設備", "河川トンネル (L79)": "円形管路 + 流速制御構造"}, ]) # (4) 道路トンネル長距離 Top 10 と河川トンネル全件の混合比較表 top10_road = road.nlargest(10, "延長(m)")[["施設名", "延長(m)", "建設年度", "道路種別", "住所(市町)"]].copy() top10_road = top10_road.rename(columns={ "施設名": "トンネル名", "延長(m)": "延長_m", "建設年度": "整備年", "住所(市町)": "所在市町"}) top10_road["分類"] = "道路 (L67)" river_short = df[["河川トンネル名称", "延長_m", "下流端_市町"]].copy() river_short = river_short.rename(columns={ "河川トンネル名称": "トンネル名", "下流端_市町": "所在市町"}) river_short["道路種別"] = "—" river_short["整備年"] = "—" river_short["分類"] = "河川 (L79)" T_long_top = pd.concat( [top10_road[["トンネル名", "延長_m", "整備年", "道路種別", "所在市町", "分類"]], river_short[["トンネル名", "延長_m", "整備年", "道路種別", "所在市町", "分類"]]], ignore_index=True) T_long_top = T_long_top.sort_values("延長_m", ascending=False).reset_index(drop=True) T_long_top["全国順位 (本記事 13 件中)"] = np.arange(1, len(T_long_top) + 1) print(f" L67 道路 {road_total} 件 / 平均延長 {road_mean_len:.1f} m", flush=True) print(f" L79 河川 {river_total} 件 / 平均延長 {river_mean_len:.1f} m", flush=True) print(f" 件数比 = {ratio_count}:1, 規模比 (平均) = {ratio_mean_len}", flush=True) print(f" H5 (規模比 ≥ 3): {h5_ok}", flush=True) print(f" ({time.time()-t6:.2f}s)", flush=True) # ============================================================================= # 7. 仮説検証総合 # ============================================================================= print("\n[7] 仮説検証", flush=True) t7 = time.time() T_hyp = pd.DataFrame([ {"仮説": "H1 全件 放水路 (RQ1)", "観測値": f"観測 = {share_bypass}% ({n_bypass}/{n_tunnels})", "判定": "強支持" if h1_ok else "反証", "解釈": (f"H1 {'強支持' if h1_ok else '反証'}: 全件 (3/3) が" f"放水路 = 都市河川の洪水バイパス。" f" 「導水トンネル」 「砂防トンネル」 等の非バイパス用途は 0 件。" f" 広島県の河川トンネル整備史 = 都市河川の洪水対策史" f"であることが量的に確認された。")}, {"仮説": "H2 延長 400-2000m 帯集中 (RQ1)", "観測値": f"観測 = {n_400_2000}/{n_tunnels} ({round(100*n_400_2000/n_tunnels,1)}%)", "判定": "強支持" if h2_ok else "反証", "解釈": (f"H2 {'強支持' if h2_ok else '反証'}: 延長" f" {df['延長_m'].min():.0f}〜{df['延長_m'].max():.0f}m。" f" 中央値 {river_median_len:.0f}m は道路トンネル中央値" f" ({ROAD_TUNNEL_MEDIAN_LEN:.0f}m) の {round(river_median_len/ROAD_TUNNEL_MEDIAN_LEN, 1)} 倍。" f" 「河川トンネル = 大規模工事を経て例外的に整備される」" f"希少性の量的証拠。")}, {"仮説": "H3 八幡川 下流恩恵 1 位 (RQ2)", "観測値": f"観測 = {yahata_rank} 位 ({yahata_benefit:.3f} km²)", "判定": "強支持" if h3_ok else "反証", "解釈": (f"H3 {'強支持' if h3_ok else '反証'}: 八幡川トンネル放水路" f"は下流端 2 km 圏内の浸水想定区域面積 (= 下流恩恵)" f" {yahata_benefit:.3f} km² で県内 {yahata_rank} 位。" f" {'広島市直下流の都市治水トンネル = 最大の治水投資効果という立地戦略の物理的証拠。' if h3_ok else '都市直下流でも下流恩恵が必ずしも最大とは限らず、 浸水想定区域の地理的偏在が支配的。'}")}, {"仮説": "H4 流下能力 ↔ 下流恩恵 正相関 (RQ2)", "観測値": f"観測 Pearson r = {r_qb:.3f}", "判定": "強支持" if h4_ok else ("反証" if not np.isnan(r_qb) else "判定不能"), "解釈": (f"H4 {'強支持' if h4_ok else ('反証' if not np.isnan(r_qb) else '判定不能')}: " f"流下能力 (m³/s) と下流恩恵 (km²) の Pearson r = {r_qb:.3f}。" f" {'大規模トンネルほど大きな治水効果という直観の量的検証。' if h4_ok else 'n=3 と小さく、 都市直下流の有無 = 立地条件が流下能力以上に支配的。'}" f" n=3 のため統計的検出力は限定。")}, {"仮説": "H5 河川 vs 道路 規模比 ≥ 3 (RQ3)", "観測値": f"観測 = 河川 {river_mean_len:.0f}m / 道路 {road_mean_len:.0f}m / 比 {ratio_mean_len}", "判定": "強支持" if h5_ok else "反証", "解釈": (f"H5 {'強支持' if h5_ok else '反証'}: 河川トンネル平均延長は" f"道路トンネルの {ratio_mean_len} 倍。" f" 件数比は逆転して {ratio_count}:1 (= 道路が圧倒)、" f"「道路 = 中規模・日常整備、 河川 = 大規模・例外整備」" f"の規模二極構造が定量化された。")}, ]) print(T_hyp[["仮説", "観測値", "判定"]].to_string(index=False), flush=True) print(f" ({time.time()-t7:.2f}s)", flush=True) # ============================================================================= # 8. CSV 出力 # ============================================================================= print("\n[8] CSV 出力", flush=True) t8 = time.time() T_overview.to_csv(ASSETS / "L79_tunnel_overview.csv", index=False, encoding="utf-8-sig") T_use.to_csv(ASSETS / "L79_use_summary.csv", index=False, encoding="utf-8-sig") T_len.to_csv(ASSETS / "L79_length_summary.csv", index=False, encoding="utf-8-sig") T_q.to_csv(ASSETS / "L79_discharge_summary.csv", index=False, encoding="utf-8-sig") T_water.to_csv(ASSETS / "L79_watershed_summary.csv", index=False, encoding="utf-8-sig") T_office.to_csv(ASSETS / "L79_office_summary.csv", index=False, encoding="utf-8-sig") T_benefit.to_csv(ASSETS / "L79_downstream_benefit.csv", index=False, encoding="utf-8-sig") T_rq2_stats.to_csv(ASSETS / "L79_rq2_stats.csv", index=False, encoding="utf-8-sig") T_compare.to_csv(ASSETS / "L79_road_vs_river_compare.csv", index=False, encoding="utf-8-sig") T_purpose.to_csv(ASSETS / "L79_purpose_compare.csv", index=False, encoding="utf-8-sig") T_long_top.to_csv(ASSETS / "L79_top_long_combined.csv", index=False, encoding="utf-8-sig") T_hyp.to_csv(ASSETS / "L79_hypothesis_check.csv", index=False, encoding="utf-8-sig") print(f" ({time.time()-t8:.2f}s)", flush=True) # ============================================================================= # 9. 図の生成 # ============================================================================= print("\n[9] 図の生成", flush=True) t9 = time.time() # 県全域 表示 bbox XMIN, YMIN = -15000, -220000 XMAX, YMAX = 125000, -90000 def save_fig(name, dpi=120): p = ASSETS / name plt.savefig(p, dpi=dpi, bbox_inches="tight", facecolor="white") plt.close('all') return p # ---- 図 1 (RQ1): 県全域 配置マップ — 上下流端 + LineString + バブル ---- print(" fig1: 配置マップ", flush=True) fig, ax = plt.subplots(figsize=(13, 8)) admin.plot(ax=ax, color="#fff8e8", edgecolor="#888", linewidth=0.4, alpha=0.55) # トンネル経路 LineString for _, r in gdf_line.iterrows(): ws = r["水系"] color = WS_COLORS.get(ws, "#666") ax.plot([pt[0] for pt in r.geometry.coords], [pt[1] for pt in r.geometry.coords], color=color, linewidth=4.0, alpha=0.7, solid_capstyle="round", zorder=3) # 上流端 (空丸) + 下流端 (塗丸) for _, r in gdf_up.iterrows(): ws = r["水系"] ax.scatter(r.geometry.x, r.geometry.y, s=80, facecolor="white", edgecolor=WS_COLORS.get(ws, "#666"), linewidth=2.0, zorder=5) for _, r in gdf_dn.iterrows(): ws = r["水系"] q = float(r["流下能力_m3s"]) ax.scatter(r.geometry.x, r.geometry.y, s=q * 5 + 80, c=WS_COLORS.get(ws, "#666"), edgecolor="#222", linewidth=1.0, alpha=0.85, zorder=6) # トンネル名ラベル for _, r in gdf_dn.iterrows(): name_short = r["河川トンネル名称"].replace("放水路", "") label = (f"{name_short}\n" f"延長 {r['延長_m']:.0f}m / 流下 {r['流下能力_m3s']:.0f}m³/s") ax.annotate(label, (r.geometry.x, r.geometry.y), xytext=(10, 8), textcoords="offset points", fontsize=9, color="#333", bbox=dict(boxstyle="round,pad=0.3", facecolor="#fff8c5", edgecolor="#d4a72c", alpha=0.8)) ax.set_xlim(XMIN, XMAX) ax.set_ylim(YMIN, YMAX) ax.set_aspect("equal") ax.set_title(f"図 1 (RQ1): 県管理 河川トンネル {n_tunnels} 本の地理配置\n" f"線 = トンネル経路 (色 = 水系), " f"白丸 = 上流端 (取水), 塗丸 = 下流端 (放流, サイズ = 流下能力)", fontsize=11) ax.set_xlabel("X (m, EPSG:6671)") ax.set_ylabel("Y (m, EPSG:6671)") handles = [Patch(facecolor=c, edgecolor="#222", label=f"{ws} 水系") for ws, c in WS_COLORS.items()] handles.extend([ Line2D([0], [0], marker='o', color='w', markerfacecolor='white', markeredgecolor='#666', markersize=10, markeredgewidth=2.0, label='上流端 (取水点)'), Line2D([0], [0], marker='o', color='w', markerfacecolor='#666', markeredgecolor='#222', markersize=12, label='下流端 (放流点)') ]) ax.legend(handles=handles, loc="lower left", fontsize=9, framealpha=0.95) ax.grid(True, linestyle="--", alpha=0.3) plt.tight_layout() save_fig("L79_fig1_tunnel_map.png") # ---- 図 2 (RQ1): 構造プロファイル 4 panels (用途/延長/流下能力/単位効率) ---- print(" fig2: 構造プロファイル", flush=True) fig, axes = plt.subplots(2, 2, figsize=(14, 10)) # (a) 用途別 円グラフ ax = axes[0, 0] use_counts = T_use.set_index("用途")["本数"] ax.pie(use_counts.values, labels=[f"{u}\n({c} 本)" for u, c in use_counts.items()], colors=["#0969da"], autopct="%d%%", startangle=90, textprops={"fontsize": 11}, wedgeprops={"edgecolor": "#222", "linewidth": 1.5}) ax.set_title("(a) 用途別 — 100% 放水路", fontsize=11) # (b) 延長クラス別 ax = axes[0, 1] xs = np.arange(len(len_order)) counts = [int(T_len.loc[T_len["延長クラス"] == s, "本数"].iloc[0]) if (T_len["延長クラス"] == s).any() else 0 for s in len_order] cols = ["#cf222e", "#cf6f00", "#1a7f37"] ax.bar(xs, counts, color=cols, edgecolor="#222", linewidth=0.6) for x, v in zip(xs, counts): if v > 0: ax.text(x, v + 0.05, f"{v} 本", ha="center", fontsize=11, fontweight="bold") ax.set_xticks(xs) ax.set_xticklabels([s.split(" ")[0] for s in len_order], fontsize=10) ax.set_ylabel("本数") ax.set_ylim(0, max(counts) * 1.3 if max(counts) > 0 else 1) ax.set_title("(b) 延長クラス別 (3 区分)", fontsize=11) ax.grid(True, axis="y", alpha=0.3) # (c) 個別 延長 + 流下能力 散布 ax = axes[1, 0] for _, r in df.iterrows(): ws = r["水系"] ax.scatter(r["延長_m"], r["流下能力_m3s"], s=400, c=WS_COLORS.get(ws, "#666"), edgecolor="#222", linewidth=1.2, alpha=0.85, label=f"{ws} 水系") ax.annotate(r["河川トンネル名称"].replace("放水路", ""), (r["延長_m"], r["流下能力_m3s"]), xytext=(8, 5), textcoords="offset points", fontsize=10, color="#333", fontweight="bold") ax.set_xlabel("延長 (m)") ax.set_ylabel("流下能力 (m³/s)") ax.set_title("(c) 延長 × 流下能力 — 個体散布", fontsize=11) ax.grid(True, linestyle="--", alpha=0.3) # 凡例 (重複削除) handles, labels = ax.get_legend_handles_labels() seen = set() ulab = [(h, l) for h, l in zip(handles, labels) if l not in seen and not seen.add(l)] ax.legend([h for h, l in ulab], [l for h, l in ulab], loc="lower right", fontsize=9) # (d) 単位流下効率 (m³/s/m) ax = axes[1, 1] df_sorted = df.sort_values("単位流下効率_m3s_per_m", ascending=True).reset_index(drop=True) ys = np.arange(len(df_sorted)) labels_y = [n.replace("放水路", "") for n in df_sorted["河川トンネル名称"]] vals = df_sorted["単位流下効率_m3s_per_m"].astype(float).tolist() cols = [WS_COLORS.get(w, "#888") for w in df_sorted["水系"].tolist()] ax.barh(ys, vals, color=cols, edgecolor="#222", linewidth=0.5) for y, v in zip(ys, vals): ax.text(v + max(vals) * 0.02, y, f"{v:.4f}", va="center", fontsize=10, fontweight="bold") ax.set_yticks(ys) ax.set_yticklabels(labels_y, fontsize=10) ax.set_xlabel("単位延長あたり流下能力 (m³/s/m)") ax.set_title("(d) 単位流下効率 — 「1m あたり何 m³/s 流せるか」", fontsize=11) ax.grid(True, axis="x", alpha=0.3) fig.suptitle(f"図 2 (RQ1): 構造プロファイル — 用途 × 延長 × " f"流下能力 × 単位効率", fontsize=12.5, y=1.00) plt.tight_layout() save_fig("L79_fig2_structure_profile.png") # ---- 図 3 (RQ2): 浸水想定区域 + 河川トンネル + 下流恩恵バッファ重ね合わせ ---- print(" fig3: 浸水想定 + バッファ重ね合わせ", flush=True) fig, ax = plt.subplots(figsize=(13, 8)) admin.plot(ax=ax, color="#fff8e8", edgecolor="#888", linewidth=0.4, alpha=0.55) # 浸水想定区域 (薄青で背景) flood.plot(ax=ax, color="#7cb9e8", edgecolor="none", alpha=0.45, zorder=2) # 下流恩恵バッファ (赤枠点線) for _, r in gdf_buf.iterrows(): gpd.GeoSeries([r.geometry], crs=TARGET_CRS).plot( ax=ax, facecolor="#ffe0b2", edgecolor="#cf222e", linewidth=1.5, linestyle="--", alpha=0.5, zorder=3) # トンネル経路 for _, r in gdf_line.iterrows(): ws = r["水系"] color = WS_COLORS.get(ws, "#666") ax.plot([pt[0] for pt in r.geometry.coords], [pt[1] for pt in r.geometry.coords], color=color, linewidth=4.0, alpha=0.95, solid_capstyle="round", zorder=4) # 下流端 for _, r in gdf_dn.iterrows(): ws = r["水系"] ax.scatter(r.geometry.x, r.geometry.y, s=350, c=WS_COLORS.get(ws, "#666"), edgecolor="#222", linewidth=1.5, alpha=0.95, zorder=5) # ラベル for _, r in gdf_dn.iterrows(): name = r["河川トンネル名称"].replace("放水路", "") benefit = T_benefit.loc[ T_benefit["河川トンネル名称"] == r["河川トンネル名称"], "下流恩恵_km2"].iloc[0] rank = T_benefit.loc[ T_benefit["河川トンネル名称"] == r["河川トンネル名称"], "順位"].iloc[0] label = f"{name}\n下流恩恵 {benefit:.3f} km²\n(ランク {rank})" ax.annotate(label, (r.geometry.x, r.geometry.y), xytext=(15, 12), textcoords="offset points", fontsize=9, color="#333", bbox=dict(boxstyle="round,pad=0.3", facecolor="#fff8c5", edgecolor="#d4a72c", alpha=0.85)) ax.set_xlim(XMIN, XMAX) ax.set_ylim(YMIN, YMAX) ax.set_aspect("equal") ax.set_title(f"図 3 (RQ2): 浸水想定区域 (青) × 河川トンネル × " f"下流恩恵 {DOWNSTREAM_BUFFER_M//1000} km バッファ\n" f"L72 既キャッシュ浸水ポリゴン {len(flood)} 件 / " f"{flood_total_km2} km² と重ね合わせ", fontsize=11) ax.set_xlabel("X (m, EPSG:6671)") ax.set_ylabel("Y (m, EPSG:6671)") handles = [ Patch(facecolor="#7cb9e8", edgecolor="none", label="浸水想定区域 (L72 cache)"), Patch(facecolor="#ffe0b2", edgecolor="#cf222e", linestyle="--", label=f"下流恩恵 {DOWNSTREAM_BUFFER_M//1000} km バッファ"), ] handles.extend([Patch(facecolor=c, edgecolor="#222", label=f"{ws}") for ws, c in WS_COLORS.items()]) ax.legend(handles=handles, loc="lower left", fontsize=9, framealpha=0.95) ax.grid(True, linestyle="--", alpha=0.3) plt.tight_layout() save_fig("L79_fig3_flood_overlay.png") # ---- 図 4 (RQ2): 下流恩恵ランキング + 流下能力散布 ---- print(" fig4: 下流恩恵ランキング", flush=True) fig, axes = plt.subplots(1, 2, figsize=(15, 6.5)) # 左: 下流恩恵ランキング ax = axes[0] T_b_sorted = T_benefit.sort_values("下流恩恵_km2", ascending=True).reset_index(drop=True) ys = np.arange(len(T_b_sorted)) labels_y = [n.replace("放水路", "") for n in T_b_sorted["河川トンネル名称"]] vals = T_b_sorted["下流恩恵_km2"].astype(float).tolist() cols = [WS_COLORS.get(w, "#888") for w in T_b_sorted["水系"].tolist()] bars = ax.barh(ys, vals, color=cols, edgecolor="#222", linewidth=0.5) for y, v in zip(ys, vals): ax.text(v + max(vals) * 0.02, y, f"{v:.3f} km²", va="center", fontsize=10, fontweight="bold") # 1 位 (= 下流恩恵最大) を赤枠で強調 top_idx = T_b_sorted["下流恩恵_km2"].idxmax() ax.barh(top_idx, vals[top_idx], facecolor="none", edgecolor="#cf222e", linewidth=2.5) ax.set_yticks(ys) ax.set_yticklabels(labels_y, fontsize=10.5) ax.set_xlabel("下流恩恵 (km²) — 下流端 2 km バッファ × 浸水想定区域") ax.set_title(f"下流恩恵ランキング — " f"1 位 = {T_b_sorted.iloc[top_idx]['河川トンネル名称'].replace('放水路','')} " f"({vals[top_idx]:.3f} km²)", fontsize=10.5) ax.grid(True, axis="x", alpha=0.3) # 右: 流下能力 vs 下流恩恵 散布 ax = axes[1] for _, r in T_benefit.iterrows(): ws = r["水系"] ax.scatter(r["流下能力_m3s"], r["下流恩恵_km2"], s=r["延長_m"] * 0.4 + 200, c=WS_COLORS.get(ws, "#666"), edgecolor="#222", linewidth=1.0, alpha=0.85, label=f"{ws} 水系") ax.annotate(r["河川トンネル名称"].replace("放水路", ""), (r["流下能力_m3s"], r["下流恩恵_km2"]), xytext=(8, 5), textcoords="offset points", fontsize=9.5, color="#333", fontweight="bold") ax.set_xlabel("流下能力 (m³/s)") ax.set_ylabel("下流恩恵 (km²)") ax.set_title(f"流下能力 × 下流恩恵 — Pearson r = {r_qb:.3f}\n" f"(円サイズ = 延長, 色 = 水系)", fontsize=10.5) handles, labels = ax.get_legend_handles_labels() seen = set() ulab = [(h, l) for h, l in zip(handles, labels) if l not in seen and not seen.add(l)] ax.legend([h for h, l in ulab], [l for h, l in ulab], loc="upper left", fontsize=9) ax.grid(True, linestyle="--", alpha=0.3) fig.suptitle("図 4 (RQ2): 下流恩恵ランキング × 流下能力相関 " "— 「治水投資効果」 の量化", fontsize=12.5, y=1.02) plt.tight_layout() save_fig("L79_fig4_downstream_benefit.png") # ---- 図 5 (RQ3): L67 道路トンネル (背景) + L79 河川トンネル (前景) ---- print(" fig5: 道路 + 河川 二極構造マップ", flush=True) # 道路トンネル GeoDataFrame 化 road_geom = [Point(lon, lat) for lon, lat in zip(road["経度（10進数）"], road["緯度（10進数）"])] gdf_road = gpd.GeoDataFrame( road, geometry=road_geom, crs="EPSG:4326").to_crs(TARGET_CRS) fig, ax = plt.subplots(figsize=(13, 8)) admin.plot(ax=ax, color="#fff8e8", edgecolor="#888", linewidth=0.4, alpha=0.55) # L67 道路トンネル (背景, 灰色小さい点) gdf_road.plot(ax=ax, color="#888", markersize=8, edgecolor="none", alpha=0.55, label=f"L67 道路トンネル ({road_total} 件)", zorder=2) # L79 河川トンネル LineString for _, r in gdf_line.iterrows(): ws = r["水系"] color = WS_COLORS.get(ws, "#666") ax.plot([pt[0] for pt in r.geometry.coords], [pt[1] for pt in r.geometry.coords], color=color, linewidth=5.0, alpha=0.95, solid_capstyle="round", zorder=4) # 下流端 (赤大きい) for _, r in gdf_dn.iterrows(): ws = r["水系"] ax.scatter(r.geometry.x, r.geometry.y, s=400, c=WS_COLORS.get(ws, "#cf222e"), edgecolor="#222", linewidth=2.0, alpha=0.95, zorder=5) name = r["河川トンネル名称"].replace("放水路", "") ax.annotate(name, (r.geometry.x, r.geometry.y), xytext=(12, 10), textcoords="offset points", fontsize=10, color="#cf222e", fontweight="bold", bbox=dict(boxstyle="round,pad=0.3", facecolor="#fff", edgecolor="#cf222e", alpha=0.9)) ax.set_xlim(XMIN, XMAX) ax.set_ylim(YMIN, YMAX) ax.set_aspect("equal") ax.set_title(f"図 5 (RQ3): L67 道路 {road_total} 件 (灰背景) + " f"L79 河川 {n_tunnels} 本 (前景) — " f"件数比 {ratio_count}:1 の二極構造", fontsize=11) ax.set_xlabel("X (m, EPSG:6671)") ax.set_ylabel("Y (m, EPSG:6671)") handles = [ Line2D([0], [0], marker='o', color='w', markerfacecolor='#888', markersize=6, label=f'L67 道路トンネル ({road_total} 件)'), ] handles.extend([Line2D([0], [0], color=c, linewidth=4, label=f'L79 {ws} 系河川トンネル') for ws, c in WS_COLORS.items()]) ax.legend(handles=handles, loc="lower left", fontsize=9, framealpha=0.95) ax.grid(True, linestyle="--", alpha=0.3) plt.tight_layout() save_fig("L79_fig5_road_vs_river_map.png") # ---- 図 6 (RQ3): 規模分布対比 (ヒスト + 箱ひげ) ---- print(" fig6: 規模分布対比", flush=True) fig, axes = plt.subplots(1, 2, figsize=(15, 6)) # 左: ヒストグラム (両者対数軸) ax = axes[0] bins_log = np.logspace(np.log10(15), np.log10(2500), 20) ax.hist(road["延長(m)"].astype(float), bins=bins_log, color="#888", edgecolor="#222", alpha=0.65, label=f"L67 道路 ({road_total} 件)") # 河川トンネル位置 (3 本) を赤縦線で for _, r in df.iterrows(): ax.axvline(r["延長_m"], color="#cf222e", linewidth=2.0, alpha=0.85, linestyle="--", zorder=5) ax.text(r["延長_m"], ax.get_ylim()[1] * 0.85, r["河川トンネル名称"].replace("放水路", "")[:4], color="#cf222e", fontsize=9, ha="center", rotation=90, fontweight="bold") ax.axvline(road_median_len, color="#0969da", linewidth=2.0, alpha=0.7, label=f"道路 中央値 {road_median_len:.0f}m") ax.axvline(river_median_len, color="#cf222e", linewidth=2.0, alpha=0.95, label=f"河川 中央値 {river_median_len:.0f}m") ax.set_xscale("log") ax.set_xlabel("延長 (m, 対数軸)") ax.set_ylabel("道路トンネル 本数") ax.set_title(f"延長分布対比 — 河川 (赤縦線) は道路の長尺帯に位置", fontsize=10.5) ax.legend(loc="upper right", fontsize=9) ax.grid(True, which="both", linestyle="--", alpha=0.3) # 右: 箱ひげ + 散布 ax = axes[1] data_box = [road["延長(m)"].astype(float).values, df["延長_m"].astype(float).values] bp = ax.boxplot(data_box, tick_labels=["L67 道路\n(n=157)", "L79 河川\n(n=3)"], patch_artist=True, widths=0.6, boxprops=dict(linewidth=1.0), medianprops=dict(color="#cf222e", linewidth=2.0)) for patch, color in zip(bp["boxes"], ["#888", "#cf222e"]): patch.set_facecolor(color) patch.set_alpha(0.5) # 散布点重ね xs_road = np.random.normal(1, 0.05, len(road)) ax.scatter(xs_road, road["延長(m)"].astype(float), s=8, c="#888", alpha=0.4, zorder=2) xs_riv = np.array([2, 2, 2]) ax.scatter(xs_riv, df["延長_m"], s=200, c="#cf222e", edgecolor="#222", linewidth=1.0, alpha=0.85, zorder=4) for x, r in zip(xs_riv, df.itertuples()): name = r.河川トンネル名称.replace("放水路", "") ax.annotate(name, (x, r.延長_m), xytext=(10, 0), textcoords="offset points", fontsize=9, color="#cf222e", fontweight="bold") ax.set_yscale("log") ax.set_ylabel("延長 (m, 対数軸)") ax.set_title(f"箱ひげ — 平均比 {ratio_mean_len:.2f} / 中央値比 {ratio_median_len:.2f}", fontsize=10.5) ax.grid(True, which="both", linestyle="--", alpha=0.3, axis="y") fig.suptitle("図 6 (RQ3): 道路 vs 河川 — 規模分布の二極構造", fontsize=12.5, y=1.02) plt.tight_layout() save_fig("L79_fig6_size_compare.png") # ---- 図 7 (RQ3): 用途・制度比較 + 件数 vs 規模 二軸 ---- print(" fig7: 用途・制度比較", flush=True) fig, axes = plt.subplots(1, 2, figsize=(15, 7)) # 左: 件数 vs 規模 (棒 + 線) ax_b = axes[0] labels_x = ["L67 道路\n(道路法)", "L79 河川\n(河川法)"] counts = [road_total, river_total] mean_lens = [road_mean_len, river_mean_len] xs = np.arange(2) b1 = ax_b.bar(xs - 0.2, counts, width=0.4, color="#888", edgecolor="#222", linewidth=0.6, label="件数") ax_b.set_ylabel("件数", color="#888") ax_b.tick_params(axis='y', labelcolor='#888') for x, v in zip(xs - 0.2, counts): ax_b.text(x, v + max(counts) * 0.02, f"{v}", ha="center", fontsize=11, fontweight="bold", color="#444") ax_b.set_yscale("log") ax_b.set_xticks(xs) ax_b.set_xticklabels(labels_x, fontsize=11) ax_b.set_title(f"件数 vs 規模 — 件数比 {ratio_count}:1 / 規模比 {ratio_mean_len:.2f}", fontsize=10.5) ax_r = ax_b.twinx() b2 = ax_r.bar(xs + 0.2, mean_lens, width=0.4, color="#cf222e", edgecolor="#222", linewidth=0.6, label="平均延長 (m)") ax_r.set_ylabel("平均延長 (m)", color="#cf222e") ax_r.tick_params(axis='y', labelcolor='#cf222e') for x, v in zip(xs + 0.2, mean_lens): ax_r.text(x, v + max(mean_lens) * 0.02, f"{v:.0f}m", ha="center", fontsize=11, fontweight="bold", color="#cf222e") ax_b.grid(True, axis="y", alpha=0.3, linestyle="--") # 右: 用途・制度比較 表 (テキスト描画) ax = axes[1] ax.axis("off") ax.set_title("制度・用途の二極構造", fontsize=11, pad=12) table_data = [["観点", "L67 道路トンネル", "L79 河川トンネル"]] for _, r in T_purpose.iterrows(): table_data.append([r["観点"], r["道路トンネル (L67)"][:22], r["河川トンネル (L79)"][:22]]) tbl = ax.table(cellText=table_data, cellLoc="left", loc="center", colWidths=[0.18, 0.40, 0.40]) tbl.auto_set_font_size(False) tbl.set_fontsize(9) tbl.scale(1, 1.5) # ヘッダ装飾 for i in range(3): cell = tbl[0, i] cell.set_facecolor("#0969da") cell.set_text_props(color="white", fontweight="bold") # 道路 vs 河川 列を色分け for r_i in range(1, len(table_data)): tbl[r_i, 1].set_facecolor("#f0f3f6") tbl[r_i, 2].set_facecolor("#fff8c5") fig.suptitle("図 7 (RQ3): 件数 × 規模 + 用途・制度比較", fontsize=12.5, y=0.99) plt.tight_layout() save_fig("L79_fig7_purpose_compare.png") print(f" ({time.time()-t9:.2f}s)", flush=True) # ============================================================================= # 10. HTML 出力 # ============================================================================= print("\n[10] HTML 出力", flush=True) t10 = time.time() def df_to_html(d): return d.to_html(index=False, classes="t", border=1, escape=False, na_rep="—") # ----- セクション 1: 学習目標と問い ----- sec1 = f"""

本記事は DoBoX のシリーズ「河川トンネル基本情報」 1 件 (dataset_id = {DATASET_ID}) を単独で取り上げ、 広島県管理の 河川トンネル {n_tunnels} 本を 3 つの独立した研究角度 (RQ1 / RQ2 / RQ3) で並列に分析する。

本記事の独自定義 (本セクションで初出):

研究の問い (3 RQ)

  1. RQ1 (主研究) — 広島県の河川トンネルの構造 — 規模・水系・ 地理分布はどう描けるか? {n_tunnels} 本 × 11 列を、 独自指標 (単位流下効率 / 治水ポテンシャル) を含む 9 軸で集計し、 「県内河川トンネルの物理形状」 を初めて系統的に記述する。
  2. RQ2 (副研究 1) — トンネルの治水機能 — 浸水想定区域への 恩恵範囲はどこまで届くか? 各トンネル下流端から {DOWNSTREAM_BUFFER_M//1000} km 圏内の浸水想定区域面積 (= 下流恩恵) を計算し、 治水投資効果をランキングする。
  3. RQ3 (副研究 2) — L67 道路トンネル {road_total} 件との 用途差はどう現れるか? 件数比 {ratio_count}:1、 平均延長比 {ratio_mean_len}、 制度差 (道路法 vs 河川法) を比較し、 「同じ『トンネル』 でも目的・規模・希少性が全く異なる」 ことを 実証する。

仮説 (5 個)

到達点

本記事を読み終えると、 学習者は以下を達成する:

""" # ----- セクション 2: 使用データ ----- sec2 = f"""

本記事が使うのは DoBoX dataset {DATASET_ID}「河川トンネル基本情報」 のみ。 ただし RQ2 では L72 既キャッシュの浸水想定区域ポリゴン ({len(flood)} 件 / {flood_total_km2} km²) を、 RQ3 では L67 既取得の道路トンネル CSV ({road_total} 件) を従属参照する。

主データ (本記事の研究対象)

項目内容備考
データセット 河川トンネル基本情報 (DoBoX #{DATASET_ID}) 広島県 河川課 提供
形式CSV (UTF-8 BOM) ~0.8 KB の超軽量データセット
件数{n_tunnels} 件 × 11 列 「個体名で全件を語れる」 粒度
主キー河川トンネル名称 八幡川/宮領川/片丘川 トンネル放水路
位置情報上流端 / 下流端 各緯度経度 (10 進) 2 端点の POINT を生成可能 ({n_tunnels*2} 点)
主要属性 事務所名 / 所在地 / 海岸河川名 / 延長(m) / 流下能力(m³/s) 「規模 + 機能」 が 2 列で表現される
延長範囲{df['延長_m'].min():.0f} 〜 {df['延長_m'].max():.0f} m 道路トンネル中央値 {ROAD_TUNNEL_MEDIAN_LEN:.0f}m の 2〜10 倍
流下能力範囲{df['流下能力_m3s'].min():.0f} 〜 {df['流下能力_m3s'].max():.0f} m³/s 1 桁台〜 100m³/s 超の 20 倍幅
ライセンスクリエイティブ・コモンズ表示 (CC-BY) DoBoX オープンデータ

形式特性の注意点

""" # ----- セクション 3: ダウンロード ----- sec3 = f"""

本記事の再現に必要なすべてを直リンクで提供する。 HTML だけ読めば学習者が完全再現できることが目標 (要件 A)。

生データ (DoBoX 1 dataset, 1 リソース)

このスクリプト本体

従属参照 (本記事は読込みのみ、 取得不要)

中間 CSV (本記事生成、 再利用可)

図 (PNG, 直 DL 可)

""" # ----- セクション 4: RQ1 ----- sec4_code = ''' import pandas as pd import numpy as np import geopandas as gpd from shapely.geometry import Point, LineString # (1) CSV 読込み (UTF-8 BOM, 3 行 × 11 列) df = pd.read_csv("data/extras/L79_river_tunnels/river_tunnel_basic.csv", encoding="utf-8-sig") # (2) 数値列の型変換 for c in ["上流端　緯度(10進数)", "上流端　経度(10進数)", "下流端　緯度(10進数)", "下流端　経度(10進数)", "延長(m)", "流下能力(m3/s)"]: df[c] = pd.to_numeric(df[c], errors="coerce") # (3) 列名短縮 + 派生列 df = df.rename(columns={ "上流端　緯度(10進数)": "上流緯度", "上流端　経度(10進数)": "上流経度", "下流端　緯度(10進数)": "下流緯度", "下流端　経度(10進数)": "下流経度", "海岸・河川名": "本川_放水路", "延長(m)": "延長_m", "流下能力(m3/s)": "流下能力_m3s", }) df["水系"] = df["本川_放水路"].str.split("　").str[0] df["放水路名"] = df["本川_放水路"].str.split("　").str[1] df["用途"] = df["河川トンネル名称"].apply( lambda s: "放水路" if "放水路" in s else "その他") # (4) 本記事独自指標 df["単位流下効率_m3s_per_m"] = df["流下能力_m3s"] / df["延長_m"] df["治水ポテンシャル指数"] = df["流下能力_m3s"] * np.log10(df["延長_m"]) # (5) GeoDataFrame 化 — 上流端/下流端 POINT + LineString 経路 gdf_up = gpd.GeoDataFrame(df, geometry=[Point(lon, lat) for lon, lat in zip(df["上流経度"], df["上流緯度"])], crs="EPSG:4326").to_crs("EPSG:6671") gdf_dn = gpd.GeoDataFrame(df, geometry=[Point(lon, lat) for lon, lat in zip(df["下流経度"], df["下流緯度"])], crs="EPSG:4326").to_crs("EPSG:6671") gdf_line = gpd.GeoDataFrame(df, geometry=[LineString([Point(ux, uy), Point(dx, dy)]) for ux, uy, dx, dy in zip(df["上流経度"], df["上流緯度"], df["下流経度"], df["下流緯度"])], crs="EPSG:4326").to_crs("EPSG:6671") # (6) 用途・延長クラス別 集計 print(df.groupby("用途").size()) print(df.groupby("水系")["延長_m"].agg(["count", "mean", "sum"])) ''' sec4 = f"""

狙い (RQ1)

RQ1 では「県の隠れた治水主役、 河川トンネル {n_tunnels} 本の物理形状」を 初めて系統的に記述する。 具体的には {n_tunnels} 本を用途 × 延長 × 流下能力 × 単位流下効率 × 治水ポテンシャル指数 × 水系 × 事務所の 9 軸で集計し、 「規模はどう分布するか / 水系はどこに集中するか / 1 m あたりの 流下効率は何 m³/s か」 を 1 枚で俯瞰できるようにする。 H1 (全件 放水路) は 「県の河川トンネル整備史 = 都市河川の洪水対策史」 仮説、 H2 (延長 400-2000m 帯集中) は「大規模工事の例外整備」 という規模クラスタの希少性仮説。

手法 — 6 ステップ

  1. STEP 1: CSV 読込み + 型変換
    CSV (UTF-8 BOM, {n_tunnels} 行 × 11 列) を read_csv() で読込み、 数値列 6 個 (上下流端 緯度経度 / 延長 / 流下能力) を pd.to_numeric(errors="coerce") で数値化。 トンネル名空欄行を除外。
  2. STEP 2: 列名短縮 + 派生列追加
    長い列名 (例: 上流端　緯度(10進数)) を短く改名。 「海岸・河川名」 列を str.split("\\u3000") で 水系列 + 放水路名列に分解 (例: 「太田川　八幡川放水路」 → 水系=「太田川」、 放水路名=「八幡川放水路」)。
  3. STEP 3: 独自指標 3 個の生成
    (a) 単位流下効率 = 流下能力 / 延長 (m³/s/m)、 (b) 治水ポテンシャル指数 = 流下能力 × log₁₀(延長)、 (c) 迂回率 = 延長 / 直線距離。 これらは公式データには無い、 本記事独自の比較指標。
  4. STEP 4: GeoDataFrame 化 (3 系列)
    上流端 POINT + 下流端 POINT + 経路 LineString の 3 つを生成、 to_crs("EPSG:6671") で平面直角第 III 系に投影。
  5. STEP 5: 行政界 sjoin で上下流端の市町判定
    L44 既キャッシュ admin_diss.gpkg を再利用。 上流端と下流端それぞれ別々に sjoin し、 「ある県内で上流→下流が どの市町を貫通するか」 を可視化。
  6. STEP 6: 9 軸集計
    用途別・延長クラス別・流下能力クラス別・水系別・事務所別の 5 グループ集計、 および個体台帳 (15 列 × {n_tunnels} 行) を生成。

実装

狙いと方法を踏まえた実装コード。 列名短縮 + 派生列 + GeoDataFrame 3 系列の 3 段構成。

{code(sec4_code)}

結果と読み取り

(a) {n_tunnels} 本 配置マップ (図 1)

なぜこの図か: {n_tunnels} 本という小さな母集団は個体名で全件を語れる 唯一の DoBoX シリーズ。 棒グラフだけで終わらせず、 全 {n_tunnels} 本を 水系色 + 経路 LineString + 上流端白丸 + 下流端塗丸 + 流下能力バブルサイズ で 1 枚に描くことで、 「県の隠れた治水トンネルがどこに分布しているか / どの方向に水を逃がすか」 を地理的に直感する (要件 T)。

{figure("assets/L79_fig1_tunnel_map.png", f"図 1: 県管理河川トンネル {n_tunnels} 本の配置 (上下流端 + 経路)")}

{df_to_html(T_water)}

図 1 / 表 (水系別) から読み取れること:

(b) 構造プロファイル 4 panels (図 2)

なぜこの図か: H1 (全件 放水路) と H2 (延長 400-2000m 帯) を 1 枚で示す ための 4 panel 構成。 円グラフ (用途) + 棒グラフ (延長クラス) + 散布 (延長 × 流下能力) + 棒 (単位流下効率) で「県管理河川トンネルの構造プロファイル」 を 4 角度から記述する。

{figure("assets/L79_fig2_structure_profile.png", f"図 2: 用途 × 延長 × 流下能力 × 単位効率 — 4 panels")}

{df_to_html(T_use)}

{df_to_html(T_len)}

{df_to_html(T_q)}

図 2 / 表から読み取れること:

(c) 個体台帳 (3 本 全項目)

なぜこの表か: {n_tunnels} 本という小規模データセットなら、 集計表だけでなく 「全件・全列」 を 1 表で読み下せる のが学習者の理解を最大化する。 延長で降順ソートし、 派生指標 (単位流下効率・治水ポテンシャル指数・迂回率) も含めた 15 列の個体台帳を提示。

{df_to_html(T_overview)}

表から読み取れること:

(d) 事務所別管轄

なぜこの表か: 維持管理の制度的単位は「事務所」 であり、 県の組織構造と治水投資配分を反映する。

{df_to_html(T_office)}

表から読み取れること:

""" # ----- セクション 5: RQ2 ----- sec5_code = ''' import geopandas as gpd # (1) L72 既キャッシュ浸水想定区域を読込み (613 件 / 829.7 km²) flood = gpd.read_file( "data/extras/L72_emergency_road/_cache_flood_max.gpkg" ).to_crs("EPSG:6671") flood["geometry"] = flood.geometry.buffer(0) # 自交差修正 # (2) 各トンネル下流端から半径 2 km バッファを生成 gdf_dn["buffer_2km"] = gdf_dn.geometry.buffer(2000) # (3) 各バッファ × 浸水想定区域 の重なり面積を計算 (= 下流恩恵) flood_union = flood.unary_union benefits = [] for _, r in gdf_dn.iterrows(): inter = r["buffer_2km"].intersection(flood_union) benefits.append({ "トンネル名": r["河川トンネル名称"], "下流恩恵_km2": inter.area / 1e6, "バッファ面積_km2": r["buffer_2km"].area / 1e6, }) T_benefit = pd.DataFrame(benefits).sort_values("下流恩恵_km2", ascending=False) print(T_benefit) # (4) 流下能力 vs 下流恩恵 の Pearson 相関 import numpy as np q = T_benefit["流下能力_m3s"].values b = T_benefit["下流恩恵_km2"].values r = np.corrcoef(q, b)[0, 1] print(f"Pearson r = {r:.3f}") ''' sec5 = f"""

狙い (RQ2)

RQ2 では「河川トンネルがどれだけの治水恩恵を下流に届けているか」を 量的に測る。 公開データには「下流人家数」 「下流被災想定」 等の直接指標が 無いため、 本記事独自に「下流恩恵 = 下流端から {DOWNSTREAM_BUFFER_M//1000} km 圏内の浸水想定区域面積 (km²)」 を定義し、 L72 既キャッシュ浸水ポリゴン ({len(flood)} 件 / {flood_total_km2} km²) と intersection で重ね合わせる。 H3 (八幡川 1 位) は「都市直下流の最大投資効果」 仮説、 H4 (流下能力 ↔ 下流恩恵 r ≥ 0.5) は「規模 ↔ 効果」 直観の量的検証仮説。

限界 (要件 J): 「下流恩恵」 は下流端 2 km 圏内の浸水想定区域面積 の代理指標。 実際の治水効果は (a) トンネル流下能力 (m³/s)、 (b) 下流河道の ボトルネック解消、 (c) 上流流域面積、 (d) 想定降雨確率 など多変数の関数。 本記事の値は「下流のどれだけのエリアが浸水想定下にあるか = 治水ニーズの大きさ」 を反映するもので、 「実際に何 % の浸水を抑制したか」 ではない。 n=3 のため統計検定力も限定的。

手法 — 4 ステップ

  1. STEP 1: L72 浸水想定区域キャッシュ読込み
    L72 (緊急輸送道路) 構築時に作成した _cache_flood_max.gpkg ({len(flood)} polygon / {flood_total_km2} km²) を再利用。 これは河川浸水想定区域の 最大想定 (= 想定最大規模) のディゾルブ済みポリゴン集合。
  2. STEP 2: 下流端 2 km バッファ生成
    geometry.buffer({DOWNSTREAM_BUFFER_M}) で半径 {DOWNSTREAM_BUFFER_M//1000} km の円形バッファを {n_tunnels} 個生成。 投影座標系が EPSG:6671 (m 単位) なので、 buffer の引数は m 単位 の数値で渡せる。
  3. STEP 3: バッファ × 浸水想定区域の intersection
    浸水ポリゴンを unary_union で 1 ジオメトリに統合してから、 各バッファと intersection。 結果ジオメトリの .area / 1e6 で km² 換算。
  4. STEP 4: 流下能力との相関 + ランキング
    下流恩恵で降順ソートして 3 本ランク付け。 流下能力との Pearson r を計算 (n=3)。

実装

狙いと方法を踏まえた実装コード。 浸水ポリゴン読込み + バッファ + intersection の 3 段構成。

{code(sec5_code)}

結果と読み取り

(a) 浸水想定 + 河川トンネル + 下流恩恵バッファ重ね合わせ (図 3)

なぜこの図か: H3 (下流恩恵) を直感的に示すには、 「浸水想定の海に トンネルが矢を放つ」 ような重ね合わせ図が最も雄弁。 浸水想定 (薄青) + 下流バッファ (赤点線) + トンネル経路 (色線) + 下流端 (大円) を 4 layer で 重ねることで、 各トンネルがどのエリアを「逃がし先」 にしているかが見える (要件 T)。

{figure("assets/L79_fig3_flood_overlay.png", f"図 3: 浸水想定 (青) × トンネル × 下流 {DOWNSTREAM_BUFFER_M//1000} km バッファ")}

{df_to_html(T_benefit)}

{df_to_html(T_rq2_stats)}

図 3 / 表から読み取れること:

(b) 下流恩恵ランキング × 流下能力相関 (図 4)

なぜこの図か: H4 (流下能力 ↔ 下流恩恵 正相関) を可視化するには、 左ランキング (= 順位) と右散布 (= 相関) の 2 panel 構成が最適。 ランキングで 「どこが大きいか」、 散布で「規模と効果はどう連動するか」 を見せる。

{figure("assets/L79_fig4_downstream_benefit.png", f"図 4: 下流恩恵ランキング × 流下能力相関")}

図 4 / 統計から読み取れること:

""" # ----- セクション 6: RQ3 ----- sec6_code = ''' # RQ3: L67 道路トンネル 157 件 vs L79 河川トンネル 3 本 の用途差比較 # (1) L67 道路トンネル CSV 読込み (16 列 × 157 件) road = pd.read_csv("data/extras/tunnel_basic.csv", encoding="utf-8-sig") road = road.dropna(subset=["施設名", "緯度（10進数）", "経度（10進数）"]) # 緯度経度の異常値補正 (緯度 > 50 の行は経度と入れ替え, L67 既知) swap_mask = road["緯度（10進数）"] > 50 road.loc[swap_mask, ["緯度（10進数）", "経度（10進数）"]] = \\ road.loc[swap_mask, ["経度（10進数）", "緯度（10進数）"]].values road["延長(m)"] = pd.to_numeric(road["延長(m)"], errors="coerce") # (2) 規模・件数比較 print(f"L67 道路 件数: {len(road)} / 平均延長 {road['延長(m)'].mean():.1f}m") print(f"L79 河川 件数: {len(df)} / 平均延長 {df['延長_m'].mean():.1f}m") print(f"件数比 (道路/河川): {len(road)/len(df):.1f}:1") print(f"規模比 (河川/道路 平均): {df['延長_m'].mean()/road['延長(m)'].mean():.2f}") # (3) 用途・制度比較は表で記述 (T_purpose) # 道路法 (1952) vs 河川法 (1964)、 5 年点検 vs 河川管理施設点検要領、 # 常時供用 vs 待機 (洪水時のみ通水), 馬蹄形断面 vs 円形管路 等 ''' sec6 = f"""

狙い (RQ3)

RQ3 では「同じ『トンネル』 でも道路と河川では用途・制度・規模が全く異なる」 ことを量的に実証する。 L67 既知 道路トンネル {road_total} 件 vs L79 河川トンネル {n_tunnels} 件の件数比 {ratio_count}:1、 規模比 (平均延長) {ratio_mean_len}、 さらに管轄法 (道路法 vs 河川法)・主目的・整備の希少性 の 7 観点で比較する。 H5 (規模比 ≥ 3) は「道路 = 中規模日常整備、 河川 = 大規模例外整備」 の規模二極構造仮説。

手法 — 4 ステップ

  1. STEP 1: L67 既取得 道路トンネル CSV 読込み
    data/extras/tunnel_basic.csv ({road_total} 件 × 16 列) を読込み。 緯度経度の異常値 (緯度 > 50 = 緯度経度入れ替え誤入力) を L67 既知の補正処理で復旧。
  2. STEP 2: 件数 + 規模統計の対比
    L67 と L79 の件数 / 平均延長 / 中央値延長 / 最長 / 最短を クロス計算。 比率を 5 行 × 4 列の表に整理。
  3. STEP 3: 用途・制度比較表の作成
    管轄法・主目的・管理者・点検制度・通常時の状態・整備の希少性・ 幾何学的形状の 7 観点でテキスト比較表を作成。 これは公式データには無い、 制度学の知識を入れた付加情報。
  4. STEP 4: 重ね合わせマップ + 規模分布対比 + 件数 × 規模 棒図
    L67 道路 {road_total} 件 (灰背景) + L79 河川 {n_tunnels} 本 (赤前景) の 重ね合わせマップ、 ヒストグラム + 箱ひげの規模対比、 件数 × 規模の 二軸対比の 3 図を生成。

実装

L67 道路トンネルデータの再利用 + 件数規模比較 + 制度比較表の 3 段構成。

{code(sec6_code)}

結果と読み取り

(a) 道路 + 河川 二極構造マップ (図 5)

なぜこの図か: 「件数比 {ratio_count}:1」 を地図 1 枚で直感する には、 道路トンネルを灰色小点で背景レイヤに、 河川トンネルを赤色 LineString + 大円で前景レイヤに重ねるのが最適。 道路の数の多さと、 河川の例外的な希少性が視覚的に対比される (要件 T)。

{figure("assets/L79_fig5_road_vs_river_map.png", f"図 5: L67 道路 {road_total} 件 (背景) + L79 河川 {n_tunnels} 件 (前景)")}

図 5 から読み取れること:

(b) 規模分布対比 — ヒストグラム + 箱ひげ (図 6)

なぜこの図か: H5 (規模比 ≥ 3) を可視化するには、 道路トンネルのヒストグラムに河川トンネル位置を縦線で重ねる のが 分かりやすい。 「道路の長尺帯 (右裾) に河川が位置する」 ことが一目で分かる。

{figure("assets/L79_fig6_size_compare.png", f"図 6: 道路 vs 河川 — ヒストグラム + 箱ひげ")}

{df_to_html(T_compare)}

図 6 / 表から読み取れること:

(c) 件数 × 規模 + 用途・制度比較 (図 7)

なぜこの図か: 「件数の二極構造」 と「用途・制度の二極構造」 を 1 枚で 示すため、 左を件数 vs 規模 二軸棒、 右を用途・制度比較表 の 2 panel 構成。 量的差 (左) と質的差 (右) を並置することで、 河川トンネルと 道路トンネルが「同名異質」であることを総合的に伝える。

{figure("assets/L79_fig7_purpose_compare.png", f"図 7: 件数 × 規模 + 用途・制度比較")}

{df_to_html(T_purpose)}

図 7 / 表から読み取れること:

""" # ----- セクション 7: 仮説検証総合 ----- sec7 = f"""

仮説検証総合表

{df_to_html(T_hyp)}

結果の総合解釈

3 RQ × 5 仮説の検証結果から、 広島県の河川トンネル {n_tunnels} 本について以下の 3 つの実証的知見が得られた:

  1. (RQ1 — 構造) 「全件 放水路」 の都市河川治水戦略
    {n_tunnels} 本すべてが放水路 (= 洪水バイパス治水トンネル) であり、 延長は{df['延長_m'].min():.0f}〜{df['延長_m'].max():.0f}mの範囲、 平均 {river_mean_len:.0f}m は道路トンネル中央値 ({ROAD_TUNNEL_MEDIAN_LEN:.0f}m) の {round(river_mean_len/ROAD_TUNNEL_MEDIAN_LEN, 1)} 倍。 H1, H2 ともに{'支持' if h1_ok else '不支持'}/{'支持' if h2_ok else '不支持'}。 これは広島県の河川トンネル整備史 = 都市河川の洪水対策史 (= 大規模例外整備) であることの量的証拠。
  2. (RQ2 — 治水機能) 立地が下流恩恵を支配する非線形関係
    下流端 2 km 圏内の浸水想定区域面積 (= 下流恩恵) で 3 本をランキング: 最大 {T_benefit['下流恩恵_km2'].max():.3f} km² 、 最小 {T_benefit['下流恩恵_km2'].min():.3f} km²。 流下能力との Pearson r = {r_qb:.3f}。 H3, H4 = {'支持' if h3_ok else '不支持'}/{'支持' if h4_ok else ('反証' if not np.isnan(r_qb) else '判定不能')}。 規模 (流下能力) は単独では下流恩恵を予測しきれず、 立地 (都市直下流かどうか) が支配的要因として浮上した。 これは「治水投資の評価軸は流下能力だけではない」 ことを示す重要な発見。
  3. (RQ3 — 用途差) 件数比 {ratio_count}:1 + 規模比 {ratio_mean_len} の二極構造
    L67 道路 {road_total} 件 vs L79 河川 {n_tunnels} 件 = 件数比 {ratio_count}:1、 平均延長比 {ratio_mean_len}。 道路 = 中規模 ({road_mean_len:.0f}m) × 多数 ({road_total} 件) = 日常整備、 河川 = 大規模 ({river_mean_len:.0f}m) × 希少 ({n_tunnels} 件) = 例外整備。 管轄法 (道路法 vs 河川法)・主目的 (交通 vs 治水)・通常時の状態 (常時供用 vs 待機) の 3 軸で完全分離。 H5 = {'支持' if h5_ok else '不支持'}。

3 RQ を統合した「県管理河川トンネル」 の見立て

RQ1 〜 RQ3 を統合すると、 広島県の河川トンネル {n_tunnels} 本は「全件 放水路 × 立地が下流恩恵を支配 × 道路との二極構造」 の 3 重特性として描ける: (1) 用途は全件 放水路の単一構造、 (2) 流下能力ではなく立地が下流恩恵を決め、 (3) 件数比 {ratio_count}:1 + 規模比 {ratio_mean_len} で道路トンネルと完全分離。 これは単なる「インフラ台帳」 ではなく、 県の都市河川治水戦略の物理的形態であり、 3 本それぞれが県の主要都市帯 (広島市西区 / 東広島市 / 三次市) を地下から守る 隠れた治水主役として価値を持つ。

3 つの研究角度 (構造 / 治水機能 / 道路対比) は完全に独立で、 単一の RQ では 見えないものを 3 RQ 並列で初めて立体的に見せる。 {n_tunnels} 本という 極端に小さな母集団でも、 個体名で全件を語れる粒度と独自指標 (下流恩恵 / 単位流下効率 / 治水ポテンシャル指数) 導入で 研究水準の深掘りが可能であることを実証した。

""" # ----- セクション 8: 発展課題 ----- sec8 = f"""

結果から導かれる新たな問い

発展課題 1: 国管理河川トンネルとの統合分析

結果 X: 本記事は県管理 {n_tunnels} 本のみを扱った。 一方、 同じ広島県内には 国土交通省管理の太田川放水路 (国直轄、 高瀬堰〜広島湾 を本川分派で バイパス) や、 中国地方整備局管理の他放水路が存在する可能性が高い。
新仮説 Y: 県管理河川トンネル (本記事 {n_tunnels} 本) と国管理放水路を 統合した「広島県内全河川バイパス施設母集団」 で同じ 3 RQ を再走したら、 H1 (全件 放水路) は支持継続、 H2 (延長 400-2000m) は太田川放水路 (全長 1.8 km, 県内最大級) が含まれることで強支持になる仮説。
課題 Z: (1) 国土交通省 中国地方整備局のオープンデータ (河川局 河川トンネル諸量) を取得。 (2) 県管理 {n_tunnels} 本と統合し、 管理者軸を加えた 4 軸クロス (用途 × 管理者 × 規模 × 水系) を生成。 (3) 「県管理 vs 国管理」 で構造プロファイルがどう違うかを量的比較し、 役割分担の制度的境界を可視化する。

発展課題 2: 実放流データを使った下流恩恵の検証

結果 X: RQ2 で「下流恩恵」 を浸水想定区域面積で代用した。 これは 静的な治水ニーズの指標であって、 「実際にどれだけの洪水流量を逃がしたか」 の動的指標ではない。
新仮説 Y: 過去 10 年の実洪水イベント (2018 西日本豪雨等) における 各トンネルの実放流量 × 時間積分と、 本記事の「下流恩恵」 は r ≥ 0.7 の強い正相関を示す仮説。 つまり静的指標が動的効果の良い代理指標である仮説。
課題 Z: (1) 県の河川管理データから各トンネルの洪水時放流流量 時系列を取得 (公開されていれば)。 (2) イベント別積算放流量を計算し、 本記事の「下流恩恵」 と相関分析。 (3) 静的指標と動的指標の整合性を 量的に評価し、 公開データから治水効果を推定する手法の信頼性を検証。

発展課題 3: 河川トンネル下流部の建物被害想定との重ね合わせ

結果 X: 「下流恩恵」 は浸水想定区域面積のみを使った。 実際の社会的価値は「面積 × 人口 × 建物資産価値」 の関数であるはず。
新仮説 Y: 八幡川トンネル放水路 (広島市直下流) は下流端 2 km 圏内の 推定資産価値が県内最高仮説。 これは「都市直下流の治水トンネル = 最大の経済効果」 という投資効率の物理的証拠。
課題 Z: (1) 国勢調査小地域 (町丁字) の人口データと 固定資産税課税台帳の家屋評価額を取得 (匿名化済の集計値)。 (2) 各トンネルのバッファ内人口・推定資産を計算。 (3) 「面積ベース下流恩恵」 と「資産ベース下流恩恵」 の差を比較し、 都市直下流ダム = 最高 ROI 投資の量的証拠を提示する。

発展課題 4: 河川トンネル整備史の年代分析

結果 X: 本データには建設年情報が無く、 RQ1 で年代分析を 割愛した。 一方、 道路トンネル (L67) は 1927-2018 の建設年を全件保持し、 1980-90s 集中という整備史を量的に示せる。
新仮説 Y: 河川トンネル {n_tunnels} 本は1980-2000s に整備された 都市河川対策第 2 期の産物仮説。 県の都市化進展で河川氾濫被害が増えた 時期に対応している。
課題 Z: (1) 県のダム史・河川史資料 (年史・記念誌等) から各トンネル の建設年・契約金額を調査。 (2) 県の都市化指標 (DID 拡大史・市街地面積年次) と整備史の相関を分析。 (3) 「都市河川対策の 3 期区分」 (S30s ダム期 / S50s 〜 H10s 河川改修期 / H10s〜トンネル放水路期) を量的に検証する。

発展課題 5: 都道府県別 河川トンネル整備密度比較

結果 X: H5 で県内の件数比 {ratio_count}:1 (道路 vs 河川) を 量的確認した。 一方、 大都市圏 (東京・大阪) は地下河川を多数整備しており、 {ratio_count}:1 とは大きく異なる比率を持つ可能性が高い。
新仮説 Y: 河川トンネル整備密度 (= 件数 / 県面積) は都市化率と正相関 仮説。 東京都・大阪府は地下河川群で密度が広島県より1 桁以上高い。 一方、 中山間県 (島根・鳥取) は密度ほぼゼロ。
課題 Z: (1) 各都道府県のオープンデータポータルから河川トンネル 整備状況を集約。 (2) 都市化率 (DID 人口比) ・1 級水系数・面積で 正規化した整備密度を計算。 (3) 全国 47 都道府県の整備密度ランキング を可視化し、 「広島県は中位」 「都市県は高位」 という相対位置を実証。 県の制度的選択を全国比較で相対化する。

""" # ----- セクション組み立て ----- sections = [ ("学習目標と問い", sec1), ("使用データ", sec2), ("ダウンロード", sec3), ("【RQ1】河川トンネル構造研究 — 規模 × 水系 × 地理分布", sec4), ("【RQ2】治水機能研究 — 下流恩恵 × 流下能力", sec5), ("【RQ3】L67 道路トンネルとの対比研究 — 用途差・規模差", sec6), ("仮説検証総合", sec7), ("発展課題", sec8), ] html = render_lesson( num=79, title="L79 河川トンネル基本情報 単独 3 研究例分析", tags=["河川トンネル", "放水路", "治水", "GIS", "オープンデータ", "下流恩恵", "道路対比"], time="30〜45 分", level="リテラシ〜中級", data_label=(f'河川トンネル基本情報 (1 dataset / 1 リソース) — ' f'{n_tunnels} 本 / {df["水系"].nunique()} 水系 / ' f'放水路 100%'), sections=sections, script_filename="L79_river_tunnels.py", ) OUT_HTML = LESSONS / "L79_river_tunnels.html" OUT_HTML.write_text(html, encoding="utf-8") print(f" HTML 出力: {OUT_HTML}", flush=True) print(f" ({time.time()-t10:.2f}s)", flush=True) # ============================================================================= # 11. 完了 # ============================================================================= print(f"\n=== L79 完了: {time.time()-t_all:.2f}s ===", flush=True) print(f" CSV (中間): {len(list(ASSETS.glob('L79_*.csv')))} ファイル", flush=True) print(f" PNG (図): {len(list(ASSETS.glob('L79_*.png')))} ファイル", flush=True)