# -*- coding: utf-8 -*- """L55 特定盛土等規制区域 単独 3 研究例分析 — 広島県の「特定盛土等規制区域」 16 ファイル・約 4,625 km² (県の 55%) を 3 つの独立した研究角度で読み解く カバー宣言: 本記事は DoBoX のシリーズ「特定盛土等規制区域」 1 件 (dataset_id = 1428) を 単独で取り上げ、 広島県内の特定盛土等規制区域 (Shapefile 16 ファイル, polygon 数百, 約 4,625 km²) を 3 つの独立した研究角度 (RQ1 / RQ2 / RQ3) で並列に分析する。 Phase 2 盛土規制系 4 本目 (= L52 許可 / L53 届出 / L54 宅造区域 / L55 特定区域) の 最終ピース。L54 の双子姉妹だが、地理的射程は山間地・農地・林地まで含む県下全域。 「特定盛土等規制区域」 とは: 盛土規制法 (= 宅地造成及び特定盛土等規制法, 2023-05-26 施行) 第 30 条第 1 項に基づき、知事が指定する区域。 L54 (12 条 1 項 = 都市計画区域内宅地) と異なり、 山間地・農地・林地等を含む県下全域を対象に指定可能。 2021 年 7 月 静岡県熱海市 伊豆山土石流災害 (盛土崩落) の上流域盛土を 規制するために新設された区域類型で、L54 が捕捉しきれない山間地大規模盛土を 法的監督下に置くための「最後のピース」。 本データは広島県全 23 市町中 16 市町・union 面積 約 4,625 km² (広島県総面積 8,479 km² の 約 55%) で、 L54 (~24%) の 2 倍以上の面積を占める。 「30 条 1 項」 (本記事独自定義の中核条文): 『都道府県知事は、宅地造成等工事規制区域以外の土地の区域であって、 土地の傾斜度、渓流の位置その他の自然的条件及び周辺地域における 土地利用の状況その他の社会的条件からみて、 当該区域内の土地において特定盛土等又は土石の堆積が行われた場合には、 これに伴う災害により市街地等区域又は市街地等区域となることが 確実であると認められる区域 (中略) の居住者その他の者の生命又は身体に 危害を生ずるおそれが特に大きいと認められるものを、 特定盛土等規制区域として指定することができる』 つまり宅造規制区域(12 条) 以外で、地形条件 + 周辺人口条件を 満たす場所を指定する。山間地に新規盛土ができても下流市街地に影響する場所が対象。 「許可閾値」 (盛土規制法施行令 第 28 条): 本区域内で許可必要となる工事の規模 (12 条より厳しい): - 盛土で高さ 1m 超の崖を生ずる工事 - 盛土で高さ 2m 超 (崖を生じない場合) - 切土で高さ 2m 超の崖を生ずる工事 - 高さ閾値以下でも面積500m² 超の盛土・切土 - 土石の堆積で面積 300m² 超 + 高さ 2m 超または面積 500m² 超 閾値以下でも本区域内なら40 条 1 項届出義務 (= L53 で扱った届出制の 30 条版)。 「30 条系運用」 (本記事独自定義): L52 (許可) / L53 (届出) は条文ごとに 12 条系と 30 条系の 2 系統を持つ: - L52 許可: 12 条 1 項 (宅造区域内) または 30 条 1 項 (特定区域内) - L53 届出: 21 条 1 項 (宅造区域内) または 40 条 1 項 (特定区域内) 本記事では L52 / L53 のうち本特定区域内に立地するものを 30 条系として識別し、 山間地監督の運用実態を定量化する。 「山間地監督の最後のピース」 (本記事独自定義): L54 (宅造区域) は都市計画区域内 = 主に沿岸都市を対象とした。 その結果 L54 で「警戒のみ」 残存エリア = 全警戒区域の 61.7% (98 km²) という 制度間ギャップが判明 (L54 RQ3)。 L55 (特定区域) はそのギャップ = 山間地警戒区域を埋めるために設計された。 本記事では L54 + L55 を合算した警戒区域カバー率を計算し、 「警戒のみ」 残存エリアの最終同定を試みる。 研究の問い (3 RQ): RQ1 (主研究): 広島県内の特定盛土等規制区域の地理範囲と山間地カバレッジ。 16 市町指定 polygon の面積分布、市町別ランキング、地形タイプ (中山間市町 / 都市部辺縁 / 林地・農地中心) の偏在を読む。 広島県 8,479 km² のうち本区域指定は何 % を占めるか? L54 との地理的補完関係 (L54 沿岸都市偏重 vs L55 山間地偏重) を可視化。 RQ2 (副研究 1): 土砂災害警戒区域 (L10/L11既知) との重複。 本区域 vs 警戒区域の重複率 (L54 の 2.9% との対比)、 L54 + L55 を合算した警戒区域カバー率、 「警戒のみ」 残存エリアの最終同定 (= 二制度合算でもなお抜ける場所)、 二重リスク (規制 + 警戒) の地理集中。 RQ3 (副研究 2): 許可・届出盛土 (L52/L53) との関係 — 30 条系の運用実態。 L52 許可盛土 152 件のうち本区域内 (= 30 条許可) の何件、 L53 届出盛土 284 件のうち本区域内、 30 条系 vs 12 条系の運用差 (規模・件数・地理)、 山間地盛土の監督強度 (件数 / 区域 km²) を診断。 仮説 (5): H1 (山間地偏重, RQ1): 16 区域は山間中山間市町に偏在し、 上位 5 市町 (庄原市・三次市・北広島町・安芸高田市・神石高原町等) で 全体の 50% 以上を占める。L54 (沿岸都市偏重) と地理的に補完。 H2 (大カバー率, RQ1): L55 単独で県の 40% 以上を覆う = L54 (24%) を超える広域指定。 これは 30 条 1 項が地形条件 + 下流人口を満たす広い領域を対象とすることに整合。 H3 (警戒重複高, RQ2): 本区域 ∩ 警戒区域の重複率は L54 (2.9%) より高い 10-30%。 なぜなら 30 条 1 項自体が地形条件 (傾斜度・渓流) + 災害リスク (下流市街地への影響)を 考慮した区域指定で、警戒区域指定基準との空間相関が高いはず。 H4 (制度間ギャップ ≪ L54 単独, RQ2): L54 + L55 合算カバーで「警戒のみ」 は L54 単独 98 km² から 30 km² 以下に縮小。L55 が L54 のギャップを埋める設計通り。 H5 (30 条系の運用は控えめ, RQ3): L52 (許可) のうち本区域内立地は 5-15% 程度。 12 条系 (L54 86%) に比べ少数派。30 条系は新法施行 (2023-05-26) 以降の指定で、 許可案件の蓄積が始まったばかり。L53 (届出) も同傾向。 要件 S 準拠 (1分以内完走): - 区域 Shapefile は MakeValid + STRtree で前処理してから dissolve - 警戒区域 (L54 同様) は dissolve で先に統合済 → overlay 高速 - L52 / L53 / L54 の処理済 CSV を参照することで比較研究は追加コスト ~0 - 期待実行時間: <60 秒 要件 T 準拠 (位置情報あり=地図必須): - RQ1: 県全体マップ (L54 と並列重ね) + 市町別 choropleth - RQ1: L54 + L55 合体マップ (制度補完の可視化) - RQ2: 規制 × 警戒 4 区分マップ (3 制度合算) - RQ3: L52 / L53 × 区域内/外 オーバーレイ 要件 Q 準拠: 図 9 / 表 12 (3 RQ × 多角度)。 データ仕様: - dataset 1428: 特定盛土等規制区域 - リソース 51146 (Shapefile, ~6 MB): 県知事指定 統合版 - リソース 51147-51161 (Shapefile, 計 15 ファイル): 市町別版 (15 市町長指定) - リソース 51162 (XLSX, 23.1 KB): 盛土規制法関係資料 (L52/L53/L54 と共通) - 形式: Shapefile (EPSG:6668 → 6671 変換) + XLSX (2 シート) - 取得日: 2026-05-09 - ライセンス: クリエイティブ・コモンズ表示 4.0 - 作成主体: 広島県土木建築局 都市環境整備課 """ from __future__ import annotations import sys, time, re, json, zipfile, shutil from pathlib import Path sys.path.insert(0, str(Path(__file__).parent)) from _common import ROOT, ASSETS, LESSONS, render_lesson, code, figure import numpy as np import pandas as pd import geopandas as gpd import shapely from shapely import make_valid from shapely.ops import unary_union import requests 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("=== L55 特定盛土等規制区域 単独 3 研究例分析 ===", flush=True) # ============================================================================= # 0. 定数・パス # ============================================================================= TARGET_CRS = "EPSG:6671" # JGD2011 平面直角座標系第 III 系 (m 単位) SOURCE_CRS = "EPSG:6668" # JGD2011 緯度経度 (本データの原 CRS) DATA_DIR = ROOT / "data" / "extras" / "L55_specific_filling_zone" DATA_DIR.mkdir(parents=True, exist_ok=True) ZONE_ZIP = DATA_DIR / "340006_specific_filling_zone_pref_shp_2023-04-27.zip" ZONE_SHP_DIR = DATA_DIR / "shp" XLSX_PATH = DATA_DIR / "340006_filling_regulation_related_documents.xlsx" CITY_DIR = DATA_DIR / "all_cities" CITY_DIR.mkdir(parents=True, exist_ok=True) # 市町別 Shapefile (RID 51147-51161 = 15 市町長指定) CITY_SHP_RIDS = [ (51147, "342050", "三原市"), (51148, "342068", "尾道市"), (51149, "342084", "府中市"), (51150, "342092", "三次市"), (51151, "342106", "庄原市"), (51152, "342114", "大竹市"), (51153, "342122", "東広島市"), (51154, "342131", "廿日市市"), (51155, "342149", "安芸高田市"), (51156, "342157", "江田島市"), (51157, "343684", "安芸太田町"), (51158, "343692", "北広島町"), (51159, "344621", "世羅町"), (51160, "345458", "神石高原町"), (51161, "344541", "大崎上島町"), ] # 関連 Shapefile (土砂災害警戒区域 = L10/L11 既存) WARN_KYUKEI = ROOT / "data" / "extras" / "sediment_shp" / "kyukeisha" / \ "340006_sediment_disaster_hazard_area_steep_slope_20260427" / \ "73_031krp_34_20260427.shp" WARN_DOSEKI = ROOT / "data" / "extras" / "sediment_shp" / "doseki" / \ "340006_sediment_disaster_hazard_area_debris_flow_20260427" / \ "73_031drp_34_20260427.shp" # 行政区域 (L15 既キャッシュ) ADMIN_PREF_ZIP = ROOT / "data" / "extras" / "L15_admin_zones" / "admin_922_広島県.zip" # L52 / L53 (許可・届出盛土) 既処理 CSV (RQ3 比較研究用) L52_PROCESSED = ASSETS / "L52_permits_processed.csv" L53_PROCESSED = ASSETS / "L53_notifications_processed.csv" # L54 (宅造区域) 既処理ジオメトリ (RQ2 制度補完研究用) — Shapefile zip を再読込 L54_DATA_DIR = ROOT / "data" / "extras" / "L54_residential_filling_zone" L54_SHP_PATH = L54_DATA_DIR / "shp" / \ "340006_residential_land_construction_regulation_area_shp_20230427.shp" L54_CITY_DIR = L54_DATA_DIR / "all_cities" # DoBoX DOBOX_BASE = "https://hiroshima-dobox.jp" HDR = {"User-Agent": "DoBoX-MDASH-textbook/1.0 (educational use)"} # 主要日付 LAW_ENFORCED = pd.Timestamp("2023-05-26") ZONE_DESIGNATED = pd.Timestamp("2023-04-27") # 本データの指定日 (= 新法施行直前) PREF_AREA_KM2 = 8479.4 # 広島県総面積 (公的値) # 市町コード → 名前 (L52/L53/L54 共通) 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: "坂町", 369: "安芸太田町", 462: "世羅町", 364: "北広島町", 545: "神石高原町", 454: "大崎上島町", } NAME_TO_CITY_CD = {v: k for k, v in CITY_NAME.items()} # ============================================================================= # 1. データ取得 (なければダウンロード) # ============================================================================= print("\n[1] データ取得", flush=True) t1 = time.time() def ensure_resource(rid: int, dst: Path, label: str): if dst.exists() and dst.stat().st_size > 1000: print(f" cache HIT {label}: {dst.name}") return url = f"{DOBOX_BASE}/resource_download/{rid}" print(f" DL <- {url}") r = requests.get(url, headers=HDR, timeout=120) r.raise_for_status() dst.write_bytes(r.content) print(f" saved -> {dst.name} ({len(r.content)} bytes)") ensure_resource(51146, ZONE_ZIP, "特定盛土区域 県知事指定 Shapefile") ensure_resource(51162, XLSX_PATH, "規制関係 XLSX") # Shapefile 展開 ZONE_SHP_DIR.mkdir(parents=True, exist_ok=True) shp_path = ZONE_SHP_DIR / \ "340006_specific_earth_filling_regulation_area_shp_20230427.shp" if not shp_path.exists(): print(f" unzipping -> {ZONE_SHP_DIR}") with zipfile.ZipFile(ZONE_ZIP) as zf: zf.extractall(ZONE_SHP_DIR) # 市町別 Shapefile (15 件) もダウンロード print(f" 15 市町別 Shapefile も取得 (= 各市町長指定の特定区域)...", flush=True) EXPECTED_FNAMES = { 51147: "342050_specific_filling_zone_mihara_shp_2023-09-29.zip", 51148: "342068_specific_filling_zone_onomichi_shp_2023-09-29.zip", 51149: "342084_specific_filling_zone_fuchu_shp_2023-09-29.zip", 51150: "342092_specific_filling_zone_miyoshi_shp_2023-09-29.zip", 51151: "342106_specific_filling_zone_shobara_shp_2023-09-29.zip", 51152: "342114_specific_filling_zone_otake_shp_2023-09-29.zip", 51153: "342122_specific_filling_zone_higashihiroshima_shp_2023-09-29.zip", 51154: "342131_specific_filling_zone_hatsukaichi_shp_2023-09-28.zip", 51155: "342149_specific_filling_zone_akitakata_shp_2023-09-29.zip", 51156: "342157_specific_filling_zone_etajima_shp_2023-09-29.zip", 51157: "343684_specific_filling_zone_akiota_shp_2023-09-29.zip", 51158: "343692_specific_filling_zone_kitahiroshima_shp_2023-09-29.zip", 51159: "344621_specific_filling_zone_sera_shp_2023-09-29.zip", 51160: "345458_specific_filling_zone_jinsekikogen_shp_2023-09-29.zip", 51161: "344541_specific_filling_zone_osakikamijima_shp_2023-09-29.zip", } for rid, code5, name in CITY_SHP_RIDS: fname = EXPECTED_FNAMES.get(rid, f"{code5}.zip") ensure_resource(rid, CITY_DIR / fname, f"市町 {code5} {name}") print(f" ({time.time()-t1:.1f}s)", flush=True) def make_valid_safe(g): if g is None or g.is_empty: return g if g.is_valid: return g try: return make_valid(g) except Exception: return g.buffer(0) # ============================================================================= # 2. 特定区域 Shapefile の読込・前処理 # 重要: 県知事ファイル (51146) は 16 市町すべての polygon を統合済み = 主要源データ # 市町別ファイル (51147-51161) は同じ polygon を市町別に切り出した 派生版 # (確認済: 庄原市の市町別ファイル 1057.94 km² = 県知事ファイル中の対応 polygon と完全一致) # → 本研究では県知事ファイル単独を主源データとする (重複二重カウント回避) # → 市町別ファイルは「市町長確認済」 のメタ情報として活用 # ============================================================================= print("\n[2] 特定区域 Shapefile 読込・前処理 (県知事ファイル 1 件 = 主源)", flush=True) t1 = time.time() # 2-1 県知事指定 (主源) zone_pref_raw = gpd.read_file(shp_path) print(f" 県知事指定: source CRS {zone_pref_raw.crs}, rows {len(zone_pref_raw)}") zone_pref_raw = zone_pref_raw.to_crs(TARGET_CRS) zone_pref_raw["geometry"] = zone_pref_raw.geometry.apply(make_valid_safe) zone_pref_raw["poly_area_km2"] = zone_pref_raw.geometry.area / 1e6 # 主要 polygon 識別: 面積 > 0.001 km² (= 1000 m²) のものを採用 THRESHOLD_KM2 = 0.001 pref_named = zone_pref_raw[zone_pref_raw["poly_area_km2"] > THRESHOLD_KM2].copy() pref_named = pref_named.reset_index(drop=True) print(f" 県知事: 全 {len(zone_pref_raw)} -> 主要 {len(pref_named)} polygon " f"(>{THRESHOLD_KM2*1e6:.0f} m^2)") print(f" 県知事 主要面積合計: {pref_named['poly_area_km2'].sum():.2f} km^2") # N03_004 で直接市町名取得 (国土数値情報の行政区域連携属性) # N03_004 が NaN の polygon は後で sjoin で補完 print(f" N03_004 (市町名) で直接同定可能: " f"{pref_named['N03_004'].notna().sum()} / {len(pref_named)}") # 2-2 市町別ファイル (15 件) はメタ情報のみ確認 — 同じ polygon の派生版なので統合せず city_meta_list = [] for rid, code5, name in CITY_SHP_RIDS: zip_files = list(CITY_DIR.glob(f"{code5}_*.zip")) if not zip_files: print(f" WARN: {code5} {name} no zip") continue zip_path = zip_files[0] extract_dir = CITY_DIR / f"tmp_{code5}" extract_dir.mkdir(exist_ok=True) if not any(extract_dir.glob("*.shp")): with zipfile.ZipFile(zip_path) as zf: zf.extractall(extract_dir) shp_files = list(extract_dir.glob("*.shp")) if not shp_files: continue g = gpd.read_file(shp_files[0]).to_crs(TARGET_CRS) g["geometry"] = g.geometry.apply(make_valid_safe) a_km2 = g.geometry.area.sum() / 1e6 n_poly = len(g[g.geometry.area > THRESHOLD_KM2 * 1e6 / 10]) city_meta_list.append({ "市町コード": code5, "市町名": name, "市町長ファイル_polygon数": n_poly, "市町長ファイル_面積_km2": round(a_km2, 2), }) city_meta_df = pd.DataFrame(city_meta_list) print(f" 市町別ファイル ({len(CITY_SHP_RIDS)}): メタ情報のみ取得 (重複統合は行わない)") print(f" 市町別 polygon 計 (= pref と重複): " f"{city_meta_df['市町長ファイル_polygon数'].sum()} polygon, " f"{city_meta_df['市町長ファイル_面積_km2'].sum():.0f} km^2") # 2-3 zones (主源 = pref) の最終データ整備 pref_named["指定主体"] = "県知事" pref_named["指定主体_code"] = "340006" zones = pref_named.reset_index(drop=True) zones["poly_id"] = np.arange(len(zones)) n_pref_polys = len(zones) n_city_polys = int(city_meta_df["市町長ファイル_polygon数"].sum()) print(f" zones (主源 = 県知事): {len(zones)} polygon, {zones['poly_area_km2'].sum():.2f} km^2") print(f" ({time.time()-t1:.1f}s)", flush=True) # ============================================================================= # 3. 行政区域・警戒区域・L54 区域の読込 # ============================================================================= print("\n[3] 行政区域・警戒区域・L54 区域の読込", flush=True) t1 = time.time() with zipfile.ZipFile(ADMIN_PREF_ZIP) as zf: name = [n for n in zf.namelist() if n.endswith(".geojson")][0] admin = gpd.read_file(f"zip://{ADMIN_PREF_ZIP}!{name}").to_crs(TARGET_CRS) admin_diss = admin.dissolve(by="CITY_CD", as_index=False) print(f" admin: {len(admin_diss)} 市町ブロック") # L54 (宅造区域) を読込 (RQ2 で制度補完合算用) print(f" L54 (宅造規制区域) 読込中...", flush=True) l54_pref_gdf = gpd.read_file(L54_SHP_PATH).to_crs(TARGET_CRS) l54_pref_gdf["geometry"] = l54_pref_gdf.geometry.apply(make_valid_safe) l54_pref_named = l54_pref_gdf[l54_pref_gdf["区域整理番"].notna()].copy() l54_city_geoms = [] for zp in sorted(L54_CITY_DIR.glob("*.zip")): code5 = zp.name[:6] extract = L54_CITY_DIR / f"tmp_{code5}" extract.mkdir(exist_ok=True) if not any(extract.glob("*.shp")): with zipfile.ZipFile(zp) as zf: zf.extractall(extract) shp = list(extract.glob("*.shp")) if shp: g = gpd.read_file(shp[0]).to_crs(TARGET_CRS) g["geometry"] = g.geometry.apply(make_valid_safe) l54_city_geoms.extend(list(g.geometry)) l54_all_geoms_list = list(l54_pref_named.geometry) + l54_city_geoms l54_union = unary_union(l54_all_geoms_list) l54_area_km2 = l54_union.area / 1e6 print(f" L54 union 面積: {l54_area_km2:.2f} km^2 ({l54_area_km2/PREF_AREA_KM2*100:.1f}% of 県)") # 警戒区域 warn_kyukei = gpd.read_file(WARN_KYUKEI).to_crs(TARGET_CRS) warn_doseki = gpd.read_file(WARN_DOSEKI).to_crs(TARGET_CRS) print(f" 警戒区域: 急傾斜 {len(warn_kyukei):,} polygon, 土石流 {len(warn_doseki):,} polygon") print(f" 警戒区域 dissolve... ", flush=True) warn_all_geom = unary_union(list(warn_kyukei.geometry) + list(warn_doseki.geometry)) warn_all = gpd.GeoDataFrame(geometry=[warn_all_geom], crs=TARGET_CRS) warn_total_km2 = warn_all_geom.area / 1e6 print(f" warn_all 統合面積: {warn_total_km2:.2f} km^2") print(f" ({time.time()-t1:.1f}s)", flush=True) # ============================================================================= # 4. polygon -> 市町紐付け (N03_004 主, sjoin で補完) # ============================================================================= print("\n[4] polygon -> 市町紐付け (N03_004 + sjoin 補完)", flush=True) t1 = time.time() # 4-1 N03_004 (国土数値情報 行政区域名) を直接使う zones["市町名"] = zones["N03_004"] # 4-2 N03_004 が NaN の polygon は sjoin (intersects) で補完 unfilled = zones[zones["市町名"].isna()].copy() print(f" N03_004 由来: {zones['市町名'].notna().sum()} / {len(zones)}") print(f" N03_004 NaN (sjoin 補完対象): {len(unfilled)} polygon, " f"{unfilled['poly_area_km2'].sum():.1f} km^2") if len(unfilled) > 0: # representative point + admin sjoin unfilled["rep_pt"] = unfilled.geometry.representative_point() upts = gpd.GeoDataFrame( unfilled[["poly_id", "poly_area_km2"]], geometry=unfilled["rep_pt"], crs=TARGET_CRS) j = gpd.sjoin(upts, admin_diss[["CITY_CD", "geometry"]], how="left", predicate="within").drop_duplicates("poly_id") pid_to_cd = j.set_index("poly_id")["CITY_CD"].to_dict() pid_to_name = {pid: CITY_NAME.get(int(cd), None) if pd.notna(cd) else None for pid, cd in pid_to_cd.items()} # 不足分は intersects 最大重複で割り当て (代表点が境界外に出る大規模 polygon 用) for pid, n in list(pid_to_name.items()): if n is None: poly = zones.loc[zones["poly_id"] == pid, "geometry"].iloc[0] best_cd = None best_area = 0 for _, ar in admin_diss.iterrows(): try: a = poly.intersection(ar.geometry).area except Exception: a = 0 if a > best_area: best_area = a best_cd = int(ar["CITY_CD"]) if best_cd is not None: pid_to_name[pid] = CITY_NAME.get(best_cd, None) # 反映 for pid, n in pid_to_name.items(): if n: zones.loc[zones["poly_id"] == pid, "市町名"] = n zones["市町名"] = zones["市町名"].fillna("不明") zones["CITY_CD"] = zones["市町名"].map(NAME_TO_CITY_CD).fillna(-1).astype(int) n_resolved = (zones["市町名"] != "不明").sum() print(f" 最終解決: {n_resolved}/{len(zones)} polygon") print(zones["市町名"].value_counts().head(20)) # 市町別に dissolve zone_dissolved = zones.dissolve(by="市町名", as_index=False) zone_dissolved["geometry"] = zone_dissolved.geometry.apply(make_valid_safe) zone_dissolved["zone_area_km2"] = zone_dissolved.geometry.area / 1e6 zone_dissolved["zone_id"] = zone_dissolved["市町名"].map(NAME_TO_CITY_CD).fillna(-1).astype(int) print(f" 市町別 dissolve 後: {len(zone_dissolved)} 市町") print(zone_dissolved[["市町名", "zone_area_km2"]].sort_values("zone_area_km2", ascending=False).head(20)) print(f" ({time.time()-t1:.1f}s)", flush=True) # ============================================================================= # 5. RQ1 集計: 地理範囲と山間地カバレッジ # ============================================================================= print("\n[5] RQ1 集計", flush=True) t1 = time.time() # 5-1 市町別指定面積ランキング city_summary = zones.groupby("市町名").agg( 指定polygon数=("poly_id", "count"), 指定面積_km2=("poly_area_km2", "sum"), ).reset_index().sort_values("指定面積_km2", ascending=False) city_summary["順位"] = np.arange(1, len(city_summary) + 1) city_summary["県シェア_%"] = (city_summary["指定面積_km2"] / PREF_AREA_KM2 * 100).round(2) city_summary["指定面積シェア_%"] = (city_summary["指定面積_km2"] / city_summary["指定面積_km2"].sum() * 100).round(2) top5_share = city_summary.head(5)["指定面積_km2"].sum() / city_summary["指定面積_km2"].sum() total_zone_area_sum = city_summary["指定面積_km2"].sum() print(f" 全県指定面積合計: {total_zone_area_sum:.2f} km^2 " f"({total_zone_area_sum/PREF_AREA_KM2*100:.1f}% of 県)") print(f" Top 5 シェア = {top5_share*100:.1f}%") print(f" Top 5: {city_summary.head(5)['市町名'].tolist()}") n_cities_designated = (city_summary["市町名"] != "不明").sum() # 5-2 polygon 個別の規模分布 poly_size_stats = zones["poly_area_km2"].describe() print(f" polygon 規模: 中央値 {poly_size_stats['50%']:.3f} km^2, " f"max {poly_size_stats['max']:.2f} km^2") # 5-3 地形タイプ分類 (本研究独自定義) — L55 は L54 と地理が逆なので分類定義も補完的に再定義 def zone_type_l55(r): cn = r["市町名"] # A_中山間 (大規模山間地市町 = 本制度の主舞台): 庄原・三次・北広島・安芸高田・神石高原・安芸太田・世羅 if cn in ("庄原市", "三次市", "北広島町", "安芸高田市", "神石高原町", "安芸太田町", "世羅町", "府中市"): return "A_中山間" # B_都市辺縁 (都市計画区域外辺縁を含む沿岸都市): 三原・尾道・東広島・廿日市 if cn in ("三原市", "尾道市", "東広島市", "廿日市市"): return "B_都市辺縁" # C_島嶼/小規模沿岸: 大崎上島・江田島・大竹 if cn in ("大崎上島町", "江田島市", "大竹市", "竹原市"): return "C_島嶼小規模" return "D_その他" zone_dissolved["地形タイプ"] = zone_dissolved.apply(zone_type_l55, axis=1) type_summary = zone_dissolved.groupby("地形タイプ").agg( 市町数=("zone_id", "count"), 指定面積_km2=("zone_area_km2", "sum"), ).reset_index().sort_values("指定面積_km2", ascending=False) type_summary["シェア_%"] = (type_summary["指定面積_km2"] / type_summary["指定面積_km2"].sum() * 100).round(1) print(f"\n 地形タイプ集計:") print(type_summary) print(f" ({time.time()-t1:.1f}s)", flush=True) # ============================================================================= # 6. RQ2 集計: 警戒区域との重ね合わせ + L54 制度補完 # ============================================================================= print("\n[6] RQ2 集計: 警戒区域 × L54 + L55 制度補完", flush=True) t1 = time.time() # 6-1 L55 union 計算 zone_pref_geom = unary_union(list(zone_dissolved.geometry)) total_zone_area_union = zone_pref_geom.area / 1e6 pref_share_union = total_zone_area_union / PREF_AREA_KM2 * 100 print(f" L55 union 面積: {total_zone_area_union:.2f} km^2 ({pref_share_union:.2f}% of 県)") # 6-2 L55 ∩ 警戒 warn_pref_geom = warn_all_geom overlap_geom_l55 = zone_pref_geom.intersection(warn_pref_geom) overlap_area_km2 = overlap_geom_l55.area / 1e6 overlap_in_zone_pct = overlap_area_km2 / total_zone_area_union * 100 overlap_in_warn_pct = overlap_area_km2 / (warn_pref_geom.area / 1e6) * 100 print(f" L55 ∩ 警戒 = {overlap_area_km2:.2f} km^2") print(f" H3: 規制区域内に占める警戒比率 = {overlap_in_zone_pct:.1f}% (予想 10-30%)") # 6-3 L54 ∪ L55 合算カバレッジ combined_geom = unary_union([l54_union, zone_pref_geom]) combined_area_km2 = combined_geom.area / 1e6 combined_share = combined_area_km2 / PREF_AREA_KM2 * 100 print(f" L54 ∪ L55 合算 = {combined_area_km2:.2f} km^2 ({combined_share:.1f}% of 県)") # L54 と L55 の重複 l54_l55_overlap = l54_union.intersection(zone_pref_geom) l54_l55_overlap_km2 = l54_l55_overlap.area / 1e6 print(f" L54 ∩ L55 = {l54_l55_overlap_km2:.2f} km^2 (= 制度間共通エリア)") # 6-4 L54 ∪ L55 ∩ 警戒 = 二制度合算で警戒区域をどれだけカバーするか combined_warn_overlap = combined_geom.intersection(warn_pref_geom) combined_warn_km2 = combined_warn_overlap.area / 1e6 combined_warn_pct = combined_warn_km2 / (warn_pref_geom.area / 1e6) * 100 print(f" (L54 ∪ L55) ∩ 警戒 = {combined_warn_km2:.2f} km^2 " f"(= 警戒の {combined_warn_pct:.1f}% を 2 制度でカバー)") # 6-5 警戒のみ (L54 + L55 でも捕捉できない警戒) = 制度間ギャップの最終残存 warn_only_combined = warn_pref_geom.difference(combined_geom) warn_only_combined_km2 = warn_only_combined.area / 1e6 warn_only_combined_pct = warn_only_combined_km2 / (warn_pref_geom.area / 1e6) * 100 print(f" 警戒のみ (L54+L55 でも残る) = {warn_only_combined_km2:.2f} km^2 " f"(全警戒の {warn_only_combined_pct:.1f}%)") # 参考: L54 単独のときの warn_only (= 98 km², 61.7%) warn_only_l54 = warn_pref_geom.difference(l54_union) warn_only_l54_km2 = warn_only_l54.area / 1e6 warn_only_l54_pct = warn_only_l54_km2 / (warn_pref_geom.area / 1e6) * 100 print(f" 参考: L54 単独で 警戒のみ = {warn_only_l54_km2:.2f} km^2 " f"({warn_only_l54_pct:.1f}%)") print(f" H4: L55 追加で警戒のみ縮小 {warn_only_l54_km2:.0f} -> {warn_only_combined_km2:.0f} km^2") # 6-6 市町別 L55 二重リスク double_risk_rows = [] for _, r in zone_dissolved.iterrows(): cn = r["市町名"] z_geom = r.geometry inter = z_geom.intersection(warn_pref_geom) inter_km2 = inter.area / 1e6 if not inter.is_empty else 0.0 double_risk_rows.append({ "市町名": cn, "zone_area_km2": r["zone_area_km2"], "二重リスク_km2": round(inter_km2, 3), "区域内警戒比_%": round(inter_km2 / r["zone_area_km2"] * 100, 1) if r["zone_area_km2"] > 0 else 0, }) double_risk = pd.DataFrame(double_risk_rows).sort_values("二重リスク_km2", ascending=False) print(f"\n 二重リスク 上位 5 区域:") print(double_risk.head(5)) top3_double = double_risk.head(3) top3_double_share = (top3_double["二重リスク_km2"].sum() / double_risk["二重リスク_km2"].sum()) print(f" 上位 3 シェア = {top3_double_share*100:.1f}%") # 規制のみ (L55 のみ, 警戒区域外) zone_only_geom = zone_pref_geom.difference(warn_pref_geom) zone_only_area_km2 = zone_only_geom.area / 1e6 print(f" L55 のみ (警戒区域外) = {zone_only_area_km2:.2f} km^2") # 警戒のみ (L55 単独) warn_only_l55_geom = warn_pref_geom.difference(zone_pref_geom) warn_only_l55_km2 = warn_only_l55_geom.area / 1e6 print(f" L55 単独基準: 警戒のみ = {warn_only_l55_km2:.2f} km^2") print(f" ({time.time()-t1:.1f}s)", flush=True) # ============================================================================= # 7. RQ3 集計: 許可・届出盛土 (L52/L53) との関係 — 30 条系の運用実態 # ============================================================================= print("\n[7] RQ3 集計: L52/L53 30 条系運用", flush=True) t1 = time.time() L52_OK = L52_PROCESSED.exists() L53_OK = L53_PROCESSED.exists() if L52_OK: df_l52 = pd.read_csv(L52_PROCESSED, encoding="utf-8-sig") print(f" L52 (許可盛土) 読込: {len(df_l52)} 件") g52 = gpd.GeoDataFrame( df_l52, geometry=gpd.points_from_xy(df_l52["lon"], df_l52["lat"]), crs="EPSG:4326").to_crs(TARGET_CRS) g52["in_l55"] = g52.geometry.within(zone_pref_geom) g52["in_l54"] = g52.geometry.within(l54_union) n_l52_in_l55 = int(g52["in_l55"].sum()) n_l52_in_l54 = int(g52["in_l54"].sum()) n_l52_in_both = int((g52["in_l54"] & g52["in_l55"]).sum()) n_l52_in_neither = int((~g52["in_l54"] & ~g52["in_l55"]).sum()) n_l52_only_l55 = int((~g52["in_l54"] & g52["in_l55"]).sum()) n_l52_only_l54 = int((g52["in_l54"] & ~g52["in_l55"]).sum()) print(f" L52 152 件:") print(f" L55 内: {n_l52_in_l55} ({n_l52_in_l55/len(g52)*100:.1f}%)") print(f" L54 内: {n_l52_in_l54} ({n_l52_in_l54/len(g52)*100:.1f}%)") print(f" L54 ∩ L55: {n_l52_in_both} 件") print(f" L55 のみ (L54 外): {n_l52_only_l55} 件 = 30 条系単独許可") print(f" L54 のみ (L55 外): {n_l52_only_l54} 件 = 12 条系単独許可") print(f" どちらも外: {n_l52_in_neither} 件") in_l55_share = n_l52_in_l55 / len(g52) print(f" H5: L52 区域内 (L55) = {in_l55_share*100:.1f}% (予想 5-15%)") # L55 内/外の規模比較 l52_in_med = g52[g52["in_l55"]][["height_max_m", "area_m2", "fill_volume_m3"]].median() l52_out_med = g52[~g52["in_l55"]][["height_max_m", "area_m2", "fill_volume_m3"]].median() print(f" L55 内 中央値: 高さ {l52_in_med['height_max_m']:.2f}m / " f"面積 {l52_in_med['area_m2']:.0f}m² / 土量 {l52_in_med['fill_volume_m3']:.0f}m³") print(f" L55 外 中央値: 高さ {l52_out_med['height_max_m']:.2f}m / " f"面積 {l52_out_med['area_m2']:.0f}m² / 土量 {l52_out_med['fill_volume_m3']:.0f}m³") else: df_l52 = None g52 = None n_l52_in_l55 = n_l52_in_l54 = n_l52_in_both = n_l52_in_neither = 0 n_l52_only_l55 = n_l52_only_l54 = 0 in_l55_share = np.nan if L53_OK: df_l53 = pd.read_csv(L53_PROCESSED, encoding="utf-8-sig") print(f"\n L53 (届出盛土) 読込: {len(df_l53)} 件") g53 = gpd.GeoDataFrame( df_l53, geometry=gpd.points_from_xy(df_l53["lon"], df_l53["lat"]), crs="EPSG:4326").to_crs(TARGET_CRS) g53["in_l55"] = g53.geometry.within(zone_pref_geom) g53["in_l54"] = g53.geometry.within(l54_union) n_l53_in_l55 = int(g53["in_l55"].sum()) n_l53_in_l54 = int(g53["in_l54"].sum()) n_l53_only_l55 = int((~g53["in_l54"] & g53["in_l55"]).sum()) n_l53_only_l54 = int((g53["in_l54"] & ~g53["in_l55"]).sum()) n_l53_neither = int((~g53["in_l54"] & ~g53["in_l55"]).sum()) in_l55_share_53 = n_l53_in_l55 / len(g53) print(f" L53 284 件:") print(f" L55 内: {n_l53_in_l55} ({in_l55_share_53*100:.1f}%) = 40 条系届出") print(f" L54 内: {n_l53_in_l54} ({n_l53_in_l54/len(g53)*100:.1f}%) = 21 条系届出") print(f" L55 のみ: {n_l53_only_l55}, L54 のみ: {n_l53_only_l54}, " f"両外: {n_l53_neither}") else: df_l53 = None g53 = None n_l53_in_l55 = n_l53_in_l54 = n_l53_only_l55 = n_l53_only_l54 = n_l53_neither = 0 in_l55_share_53 = np.nan # 区域別の許可・届出件数 if L52_OK and L53_OK: zd_for_join = zone_dissolved[["市町名", "zone_area_km2", "geometry"]].rename( columns={"市町名": "市町名_zone"}) g52_zone = gpd.sjoin(g52, zd_for_join, how="left", predicate="within") l52_count = (g52_zone.dropna(subset=["市町名_zone"]) .groupby("市町名_zone").size().reset_index(name="L52_30条許可")) l52_count = l52_count.rename(columns={"市町名_zone": "市町名"}) g53_zone = gpd.sjoin(g53, zd_for_join, how="left", predicate="within") l53_count = (g53_zone.dropna(subset=["市町名_zone"]) .groupby("市町名_zone").size().reset_index(name="L53_40条届出")) l53_count = l53_count.rename(columns={"市町名_zone": "市町名"}) zone_density = zone_dissolved[["市町名", "zone_area_km2"]].merge( l52_count, on="市町名", how="left").merge( l53_count, on="市町名", how="left").fillna(0) zone_density["L52_30条許可"] = zone_density["L52_30条許可"].astype(int) zone_density["L53_40条届出"] = zone_density["L53_40条届出"].astype(int) zone_density["L52_密度_件km2"] = (zone_density["L52_30条許可"] / zone_density["zone_area_km2"]).round(4) zone_density["L53_密度_件km2"] = (zone_density["L53_40条届出"] / zone_density["zone_area_km2"]).round(4) zone_density = zone_density.sort_values("L53_40条届出", ascending=False) print(f"\n 区域内 L52/L53 件数 上位 5:") print(zone_density.head(5)[["市町名", "zone_area_km2", "L52_30条許可", "L53_40条届出"]]) else: zone_density = pd.DataFrame() print(f" ({time.time()-t1:.1f}s)", flush=True) # ============================================================================= # 8. CSV 出力 (再現用中間データ) # ============================================================================= print("\n[8] 中間 CSV 出力", flush=True) t1 = time.time() # 市町別 dissolve サマリ zone_summary_out = zone_dissolved[["zone_id", "市町名", "zone_area_km2", "地形タイプ"]].copy() zone_summary_out.to_csv(ASSETS / "L55_zones_summary.csv", index=False, encoding="utf-8-sig") # 市町別ランキング city_summary.to_csv(ASSETS / "L55_city_designation_ranking.csv", index=False, encoding="utf-8-sig") # polygon 個別 (主源 = 県知事ファイル分) zones_csv = pd.DataFrame(zones.drop(columns=["geometry"])) rep_pts = zones.geometry.representative_point().to_crs("EPSG:4326") zones_csv["lon_rep"] = rep_pts.x.values zones_csv["lat_rep"] = rep_pts.y.values # 大きすぎる列は削減 keep_cols = ["指定主体", "指定主体_code", "poly_area_km2", "poly_id", "市町名", "CITY_CD", "lon_rep", "lat_rep", "N03_001", "N03_003", "N03_004", "N03_007"] zones_csv = zones_csv[[c for c in keep_cols if c in zones_csv.columns]] zones_csv.to_csv(ASSETS / "L55_zones_polygons.csv", index=False, encoding="utf-8-sig") # 市町メタ (15 市町長ファイル) city_meta_df.to_csv(ASSETS / "L55_city_meta.csv", index=False, encoding="utf-8-sig") # RQ2 二重リスク double_risk.to_csv(ASSETS / "L55_double_risk_by_zone.csv", index=False, encoding="utf-8-sig") # RQ3 区域内 L52/L53 if L52_OK and L53_OK: zone_density.to_csv(ASSETS / "L55_zone_30jou_count.csv", index=False, encoding="utf-8-sig") g52_with = pd.DataFrame(g52.drop(columns=["geometry"])) g52_with.to_csv(ASSETS / "L55_l52_in_l55_l54.csv", index=False, encoding="utf-8-sig") g53_with = pd.DataFrame(g53.drop(columns=["geometry"])) g53_with.to_csv(ASSETS / "L55_l53_in_l55_l54.csv", index=False, encoding="utf-8-sig") # 元 Shapefile zip コピー shutil.copy(ZONE_ZIP, ASSETS / "L55_zones_pref_raw.zip") print(f" CSV 6 ファイル + Shapefile zip 1 ファイル") print(f" ({time.time()-t1:.1f}s)", flush=True) # ============================================================================= # 9. 図の生成 (9 図) # ============================================================================= print("\n[9] 図の生成 (9 図)", flush=True) t1 = time.time() def save_fig(name, dpi=120): p = ASSETS / name plt.savefig(p, dpi=dpi, bbox_inches="tight", facecolor="white") plt.close('all') return p PREF_BOUNDS = admin_diss.total_bounds PXMIN, PYMIN, PXMAX, PYMAX = PREF_BOUNDS L55_COLOR = "#7c3aed" # 特定区域: 紫 (山間部 = L54 の緑と対比) L54_COLOR = "#2da44e" # 宅造区域: 緑 WARN_KYUKEI_COLOR = "#cf222e" WARN_DOSEKI_COLOR = "#bf8700" PERMIT_IN_COLOR = "#1f6feb" PERMIT_OUT_COLOR = "#a045a0" TYPE_COLOR = { "A_中山間": "#2da44e", "B_都市辺縁": "#1f6feb", "C_島嶼小規模": "#fa8c16", "D_その他": "#888", } # ---- 図 1 (RQ1): L55 全区域マップ + 市町別指定面積 choropleth ---- fig, axes = plt.subplots(1, 2, figsize=(16, 7)) # 1-a 全区域マップ ax = axes[0] admin_diss.boundary.plot(ax=ax, color="#888", linewidth=0.5) admin_diss.plot(ax=ax, facecolor="#f8f8f5", edgecolor="#aaa", linewidth=0.4, alpha=0.8) zone_dissolved.plot(ax=ax, color=L55_COLOR, edgecolor="#3b1a64", linewidth=0.4, alpha=0.55, label=f"L55 特定区域 ({n_cities_designated} 市町)") # 上位 5 にラベル for _, r in city_summary.head(5).iterrows(): z = zone_dissolved[zone_dissolved["市町名"] == r["市町名"]] if len(z) == 0: continue g = z.geometry.iloc[0] cx, cy = g.representative_point().x, g.representative_point().y ax.annotate(f"{r['市町名']}\n{r['指定面積_km2']:.0f} km²", (cx, cy), fontsize=9, ha="center", va="center", color="black", weight="bold", bbox=dict(boxstyle="round,pad=0.2", facecolor="white", alpha=0.85)) ax.set_xlim(PXMIN, PXMAX) ax.set_ylim(PYMIN, PYMAX) ax.set_aspect("equal") ax.set_title(f"図 1a: L55 特定盛土等規制区域 ({len(zones)} polygon, " f"{total_zone_area_union:.0f} km², 県の {pref_share_union:.1f}%)", fontsize=11) ax.set_xlabel("X (m, EPSG:6671)", fontsize=9) ax.set_ylabel("Y (m, EPSG:6671)", fontsize=9) ax.legend(loc="upper right", fontsize=10) # 1-b 市町別指定面積 choropleth ax = axes[1] admin_area = admin_diss.copy() city_to_area = dict(zip(city_summary["市町名"], city_summary["指定面積_km2"])) admin_area["指定面積_km2"] = admin_area["CITY_CD"].map( {NAME_TO_CITY_CD.get(n, -1): v for n, v in city_to_area.items()}).fillna(0) admin_area[admin_area["指定面積_km2"] == 0].plot(ax=ax, facecolor="#f0f0f0", edgecolor="#888", linewidth=0.4) admin_area[admin_area["指定面積_km2"] > 0].plot(ax=ax, column="指定面積_km2", cmap="Purples", edgecolor="#444", linewidth=0.4, legend=True, legend_kwds={"label": "指定面積 (km²)", "shrink": 0.7}) for _, r in city_summary.head(5).iterrows(): cd = NAME_TO_CITY_CD.get(r["市町名"], -1) if cd > 0: a = admin_diss[admin_diss["CITY_CD"] == cd] if len(a) == 0: continue cx, cy = a.geometry.iloc[0].centroid.x, a.geometry.iloc[0].centroid.y ax.annotate(f"{r['市町名']}\n{r['指定面積_km2']:.0f} km²", (cx, cy), fontsize=9, ha="center", va="center", color="black", weight="bold", bbox=dict(boxstyle="round,pad=0.2", facecolor="white", alpha=0.85)) ax.set_xlim(PXMIN, PXMAX) ax.set_ylim(PYMIN, PYMAX) ax.set_aspect("equal") ax.set_title(f"図 1b: 市町別指定面積 choropleth — Top 5 シェア {top5_share*100:.0f}%", fontsize=11) ax.set_xlabel("X (m)", fontsize=9) ax.set_ylabel("Y (m)", fontsize=9) fig.suptitle(f"図 1 (RQ1): L55 特定盛土等規制区域の県内地理範囲 — " f"{n_cities_designated} 市町 / {total_zone_area_union:.0f} km² / 県の {pref_share_union:.1f}%", fontsize=13, y=1.02) plt.tight_layout() save_fig("L55_fig1_zone_map.png") # ---- 図 2 (RQ1): L54 vs L55 制度補完マップ + 合算カバレッジ ---- fig, axes = plt.subplots(1, 2, figsize=(16, 7)) # 2-a L54 vs L55 重ね地図 (補完性の可視化) ax = axes[0] admin_diss.boundary.plot(ax=ax, color="#888", linewidth=0.5) admin_diss.plot(ax=ax, facecolor="#f8f8f5", edgecolor="#aaa", linewidth=0.4, alpha=0.6) # L54 (緑) を先に gpd.GeoDataFrame(geometry=[l54_union], crs=TARGET_CRS).plot( ax=ax, color=L54_COLOR, alpha=0.5, edgecolor="#1a4d2e", linewidth=0.4, label=f"L54 宅造区域 ({l54_area_km2:.0f} km²)") # L55 (紫) を後で zone_dissolved.plot(ax=ax, color=L55_COLOR, alpha=0.45, edgecolor="#3b1a64", linewidth=0.4, label=f"L55 特定区域 ({total_zone_area_union:.0f} km²)") ax.set_xlim(PXMIN, PXMAX) ax.set_ylim(PYMIN, PYMAX) ax.set_aspect("equal") ax.set_title(f"図 2a: L54 (緑) × L55 (紫) — 制度補完の地理 (L54 ∩ L55 = {l54_l55_overlap_km2:.1f} km²)", fontsize=11) legend_h = [ Patch(facecolor=L54_COLOR, alpha=0.5, label=f"L54 宅造 ({l54_area_km2:.0f} km², 沿岸都市)"), Patch(facecolor=L55_COLOR, alpha=0.45, label=f"L55 特定 ({total_zone_area_union:.0f} km², 山間地)"), ] ax.legend(handles=legend_h, loc="upper right", fontsize=10, title="制度") # 2-b polygon 規模分布 (log) ax = axes[1] poly_a = zones["poly_area_km2"][zones["poly_area_km2"] > 0] ax.hist(poly_a, bins=np.logspace(np.log10(poly_a.min()), np.log10(poly_a.max()), 30), color=L55_COLOR, edgecolor="#333", alpha=0.85) ax.set_xscale("log") ax.set_xlabel("polygon 個別面積 (km², log)", fontsize=11) ax.set_ylabel("件数", fontsize=11) ax.axvline(poly_a.median(), color="#cf222e", linestyle="--", linewidth=2, label=f"中央値 {poly_a.median():.3f} km²") ax.set_title(f"図 2b: polygon 規模分布 (n={len(zones)}) — 中央値 " f"{poly_a.median():.3f} km², max {poly_a.max():.0f} km²", fontsize=11) ax.legend(fontsize=9) ax.grid(alpha=0.3, which="both") fig.suptitle(f"図 2 (RQ1): 制度補完マップ + polygon 規模分布 — " f"L54 ∪ L55 = {combined_area_km2:.0f} km² ({combined_share:.1f}% of 県)", fontsize=13, y=1.02) plt.tight_layout() save_fig("L55_fig2_l54_l55_complementarity.png") # ---- 図 3 (RQ1): 地形タイプ別構成 + 市町別ランキング ---- fig, axes = plt.subplots(1, 2, figsize=(16, 6)) # 3-a 地形タイプパイ ax = axes[0] type_present = type_summary[type_summary["指定面積_km2"] > 0] labels = [f"{t}\n{int(n)} 市町 / {a:.0f} km² ({s:.0f}%)" for t, n, a, s in zip(type_present["地形タイプ"], type_present["市町数"], type_present["指定面積_km2"], type_present["シェア_%"])] colors = [TYPE_COLOR.get(t, "#888") for t in type_present["地形タイプ"]] ax.pie(type_present["指定面積_km2"], labels=labels, colors=colors, startangle=90, wedgeprops=dict(edgecolor="white", linewidth=2), textprops={"fontsize": 10}) ax.set_title(f"図 3a: 地形タイプ別 指定面積構成 (全 {total_zone_area_union:.0f} km²)", fontsize=11) # 3-b 市町別ランキング棒 ax = axes[1] top15 = city_summary[city_summary["市町名"] != "不明"].head(15).iloc[::-1] y = np.arange(len(top15)) bar_colors = [TYPE_COLOR.get(zone_dissolved[zone_dissolved["市町名"] == cn] ["地形タイプ"].iloc[0], "#888") if (zone_dissolved["市町名"] == cn).any() else "#888" for cn in top15["市町名"]] ax.barh(y, top15["指定面積_km2"], color=bar_colors, edgecolor="#333") ax.set_yticks(y) ax.set_yticklabels(top15["市町名"], fontsize=9) ax.set_xlabel("指定面積 (km²)", fontsize=11) ax.set_title(f"図 3b: 市町別指定面積 (Top 15) — 全県 {total_zone_area_sum:.0f} km²", fontsize=11) ax.grid(axis="x", alpha=0.3) for i, (cn, v, sh) in enumerate(zip(top15["市町名"], top15["指定面積_km2"], top15["指定面積シェア_%"])): ax.text(v + 5, i, f"{v:.0f} km² ({sh:.1f}%)", va="center", fontsize=8) legend_handles = [Patch(facecolor=c, edgecolor="#333", label=t) for t, c in TYPE_COLOR.items() if t != "D_その他"] ax.legend(handles=legend_handles, loc="lower right", fontsize=8, title="地形タイプ") fig.suptitle(f"図 3 (RQ1): 地形タイプ別構成 + 市町別ランキング — " f"H1: Top 5 シェア {top5_share*100:.0f}%", fontsize=13, y=1.02) plt.tight_layout() save_fig("L55_fig3_type_ranking.png") # ---- 図 4 (RQ2): L55 × 警戒区域 重ね合わせ ---- fig, ax = plt.subplots(figsize=(13, 8.5)) admin_diss.boundary.plot(ax=ax, color="#888", linewidth=0.4) admin_diss.plot(ax=ax, facecolor="#fafafa", edgecolor="#bbb", linewidth=0.3, alpha=0.6) warn_kyukei.plot(ax=ax, color=WARN_KYUKEI_COLOR, alpha=0.30, edgecolor="none") warn_doseki.plot(ax=ax, color=WARN_DOSEKI_COLOR, alpha=0.30, edgecolor="none") zone_dissolved.plot(ax=ax, color=L55_COLOR, alpha=0.40, edgecolor="#3b1a64", linewidth=0.4) # 上位 3 二重リスク市町ラベル for _, r in top3_double.iterrows(): z = zone_dissolved[zone_dissolved["市町名"] == r["市町名"]] if len(z) > 0: cx = z.geometry.iloc[0].representative_point().x cy = z.geometry.iloc[0].representative_point().y ax.annotate(f"{r['市町名']}\n二重リスク {r['二重リスク_km2']:.1f} km²", (cx, cy), fontsize=9, ha="center", va="center", color="black", weight="bold", bbox=dict(boxstyle="round,pad=0.2", facecolor="white", alpha=0.85)) ax.set_xlim(PXMIN, PXMAX) ax.set_ylim(PYMIN, PYMAX) ax.set_aspect("equal") ax.set_title(f"図 4 (RQ2): L55 (紫) × 警戒区域 (赤=急傾斜, 橙=土石流) — " f"重複 {overlap_area_km2:.1f} km² (L55 内 {overlap_in_zone_pct:.1f}%, " f"警戒の {overlap_in_warn_pct:.1f}%)", fontsize=12) legend_handles = [ Patch(facecolor=L55_COLOR, alpha=0.5, label="L55 特定盛土区域"), Patch(facecolor=WARN_KYUKEI_COLOR, alpha=0.4, label="急傾斜警戒"), Patch(facecolor=WARN_DOSEKI_COLOR, alpha=0.4, label="土石流警戒"), ] ax.legend(handles=legend_handles, loc="upper right", fontsize=10, title="凡例") ax.set_xlabel("X (m, EPSG:6671)", fontsize=10) ax.set_ylabel("Y (m, EPSG:6671)", fontsize=10) plt.tight_layout() save_fig("L55_fig4_overlay_warn.png") # ---- 図 5 (RQ2): L54 単独 vs L54+L55 合算の警戒のみ縮小 ---- fig, axes = plt.subplots(1, 2, figsize=(16, 7)) # 5-a 4 区分パイ (L54 ∪ L55 合算 vs 警戒) ax = axes[0] both_in_combined = combined_warn_km2 # (L54∪L55) ∩ 警戒 combined_only = combined_area_km2 - both_in_combined # 規制のみ (合算) warn_only_final = warn_only_combined_km2 # 警戒のみ (合算後) neither = PREF_AREA_KM2 - (both_in_combined + combined_only + warn_only_final) parts = [ (f"規制 ∩ 警戒\n(L54+L55 合算)", both_in_combined, "#cf222e"), (f"L54+L55 のみ\n(警戒外)", combined_only, "#7c3aed"), (f"警戒のみ\n(両制度外 = 残るギャップ)", warn_only_final, "#fa8c16"), ("どちらも該当外", max(neither, 0), "#dddddd"), ] labels = [f"{n}\n{a:.0f} km² ({a/PREF_AREA_KM2*100:.1f}%)" for n, a, _ in parts] colors = [c for _, _, c in parts] ax.pie([a for _, a, _ in parts], labels=labels, colors=colors, startangle=90, wedgeprops=dict(edgecolor="white", linewidth=2), textprops={"fontsize": 9}) ax.set_title(f"図 5a: 県全土 {PREF_AREA_KM2:.0f} km² の 4 区分 (L54+L55 合算 × 警戒)", fontsize=11) # 5-b L54 単独 vs L54+L55 合算: 警戒のみ縮小の比較棒 ax = axes[1] data = [ ("L54 単独", warn_only_l54_km2, warn_only_l54_pct, L54_COLOR), ("L54+L55 合算", warn_only_combined_km2, warn_only_combined_pct, "#7c3aed"), ] xs = np.arange(len(data)) ax.bar(xs, [d[1] for d in data], color=[d[3] for d in data], edgecolor="#333") ax.set_xticks(xs) ax.set_xticklabels([d[0] for d in data], fontsize=11) ax.set_ylabel("警戒のみ 残存面積 (km²)", fontsize=11) ax.set_title(f"図 5b: 警戒のみ残存 — L55 追加で {warn_only_l54_km2:.0f} → " f"{warn_only_combined_km2:.0f} km² " f"({(warn_only_l54_km2-warn_only_combined_km2)/warn_only_l54_km2*100:.0f}% 縮小)", fontsize=11) for i, (n, v, p, _) in enumerate(data): ax.text(i, v + 1, f"{v:.0f} km²\n(全警戒の {p:.0f}%)", ha="center", va="bottom", fontsize=10, weight="bold") ax.grid(axis="y", alpha=0.3) # 削減幅を示す矢印 ax.annotate("", xy=(1, warn_only_combined_km2), xytext=(0, warn_only_l54_km2), arrowprops=dict(arrowstyle="->", color="#cf222e", lw=2)) ax.text(0.5, (warn_only_l54_km2 + warn_only_combined_km2)/2, f"削減 {warn_only_l54_km2-warn_only_combined_km2:.0f} km²", color="#cf222e", fontsize=11, weight="bold", ha="center", bbox=dict(boxstyle="round,pad=0.2", facecolor="white", alpha=0.85)) fig.suptitle(f"図 5 (RQ2): 警戒区域の制度カバレッジ — H4 L55 がギャップを埋める", fontsize=13, y=1.02) plt.tight_layout() save_fig("L55_fig5_warn_coverage_l54_l55.png") # ---- 図 6 (RQ2): 市町別二重リスク choropleth ---- fig, axes = plt.subplots(1, 2, figsize=(16, 7)) # 6-a choropleth ax = axes[0] admin_dr = admin_diss.copy() city_to_dr = dict(zip(double_risk["市町名"], double_risk["二重リスク_km2"])) admin_dr["二重リスク_km2"] = admin_dr["CITY_CD"].map( {NAME_TO_CITY_CD.get(n, -1): v for n, v in city_to_dr.items()}).fillna(0) admin_dr[admin_dr["二重リスク_km2"] == 0].plot(ax=ax, facecolor="#f0f0f0", edgecolor="#888", linewidth=0.4) admin_dr[admin_dr["二重リスク_km2"] > 0].plot(ax=ax, column="二重リスク_km2", cmap="OrRd", edgecolor="#444", linewidth=0.4, legend=True, legend_kwds={"label": "二重リスク (km²)", "shrink": 0.7}) for _, r in top3_double.iterrows(): cd = NAME_TO_CITY_CD.get(r["市町名"], -1) if cd > 0: a = admin_diss[admin_diss["CITY_CD"] == cd] if len(a) == 0: continue cx, cy = a.geometry.iloc[0].centroid.x, a.geometry.iloc[0].centroid.y ax.annotate(f"{r['市町名']}\n{r['二重リスク_km2']:.1f} km²", (cx, cy), fontsize=9, ha="center", va="center", color="black", weight="bold", bbox=dict(boxstyle="round,pad=0.2", facecolor="white", alpha=0.85)) ax.set_xlim(PXMIN, PXMAX) ax.set_ylim(PYMIN, PYMAX) ax.set_aspect("equal") ax.set_title(f"図 6a: 市町別 L55 × 警戒 二重リスク面積", fontsize=11) ax.set_xlabel("X (m)", fontsize=9) ax.set_ylabel("Y (m)", fontsize=9) # 6-b ランキング ax = axes[1] dr_top = double_risk[double_risk["二重リスク_km2"] > 0].head(12).iloc[::-1] y = np.arange(len(dr_top)) ax.barh(y, dr_top["zone_area_km2"], color="#cccccc", alpha=0.4, edgecolor="#888", zorder=0, label="区域面積 (km²)") ax.barh(y, dr_top["二重リスク_km2"], color="#cf222e", edgecolor="#333", label="二重リスク面積 (km²)") ax.set_yticks(y) ax.set_yticklabels(dr_top["市町名"], fontsize=9) ax.set_xlabel("面積 (km²)", fontsize=11) ax.set_title(f"図 6b: L55 二重リスク Top 12 — 上位 3 シェア {top3_double_share*100:.0f}%", fontsize=11) ax.grid(axis="x", alpha=0.3) ax.legend(fontsize=9, loc="lower right") for i, (cn, dr_v, z_v, pct) in enumerate(zip(dr_top["市町名"], dr_top["二重リスク_km2"], dr_top["zone_area_km2"], dr_top["区域内警戒比_%"])): ax.text(z_v + 5, i, f"{dr_v:.1f} ({pct:.1f}%)", va="center", fontsize=8) fig.suptitle(f"図 6 (RQ2): L55 二重リスク市町別ランキング — 警戒区域内 0.9% (低重複)", fontsize=13, y=1.02) plt.tight_layout() save_fig("L55_fig6_double_risk.png") # ---- 図 7 (RQ3): L52/L53 × L55 in/out オーバーレイ ---- if L52_OK and L53_OK: fig, axes = plt.subplots(1, 2, figsize=(16, 7)) # 7-a L52 (許可) を L55 内 (青) / L54 内 (緑) / 両外 (灰) で色分け ax = axes[0] admin_diss.boundary.plot(ax=ax, color="#888", linewidth=0.4) admin_diss.plot(ax=ax, facecolor="#fafafa", edgecolor="#bbb", linewidth=0.3, alpha=0.6) # L54, L55 区域も薄く gpd.GeoDataFrame(geometry=[l54_union], crs=TARGET_CRS).plot( ax=ax, color=L54_COLOR, alpha=0.20, edgecolor="none") zone_dissolved.plot(ax=ax, color=L55_COLOR, alpha=0.20, edgecolor="none") # L52 点 g52_l54 = g52[g52["in_l54"] & ~g52["in_l55"]] g52_l55 = g52[g52["in_l55"] & ~g52["in_l54"]] g52_neither = g52[~g52["in_l54"] & ~g52["in_l55"]] ax.scatter(g52_l54.geometry.x, g52_l54.geometry.y, s=30, c=L54_COLOR, edgecolor="black", linewidth=0.4, alpha=0.9, marker="o", label=f"L54 内 (12 条系) {len(g52_l54)}") ax.scatter(g52_l55.geometry.x, g52_l55.geometry.y, s=50, c=L55_COLOR, edgecolor="black", linewidth=0.4, alpha=0.95, marker="^", label=f"L55 内 (30 条系) {len(g52_l55)}") ax.scatter(g52_neither.geometry.x, g52_neither.geometry.y, s=30, c="#888", edgecolor="black", linewidth=0.4, alpha=0.7, marker="x", label=f"両外 {len(g52_neither)}") ax.set_xlim(PXMIN, PXMAX) ax.set_ylim(PYMIN, PYMAX) ax.set_aspect("equal") ax.set_title(f"図 7a: L52 許可盛土 152 件 × L54/L55 in/out — " f"L55 内 {n_l52_in_l55} 件 ({in_l55_share*100:.1f}%)", fontsize=11) ax.legend(loc="upper right", fontsize=9) # 7-b L53 (届出) も同様 ax = axes[1] admin_diss.boundary.plot(ax=ax, color="#888", linewidth=0.4) admin_diss.plot(ax=ax, facecolor="#fafafa", edgecolor="#bbb", linewidth=0.3, alpha=0.6) gpd.GeoDataFrame(geometry=[l54_union], crs=TARGET_CRS).plot( ax=ax, color=L54_COLOR, alpha=0.20, edgecolor="none") zone_dissolved.plot(ax=ax, color=L55_COLOR, alpha=0.20, edgecolor="none") g53_l54 = g53[g53["in_l54"] & ~g53["in_l55"]] g53_l55 = g53[g53["in_l55"] & ~g53["in_l54"]] ax.scatter(g53_l54.geometry.x, g53_l54.geometry.y, s=20, c=L54_COLOR, edgecolor="black", linewidth=0.3, alpha=0.85, marker="o", label=f"L54 内 (21 条系) {len(g53_l54)}") ax.scatter(g53_l55.geometry.x, g53_l55.geometry.y, s=30, c=L55_COLOR, edgecolor="black", linewidth=0.3, alpha=0.9, marker="^", label=f"L55 内 (40 条系) {len(g53_l55)}") ax.set_xlim(PXMIN, PXMAX) ax.set_ylim(PYMIN, PYMAX) ax.set_aspect("equal") ax.set_title(f"図 7b: L53 届出盛土 284 件 × L54/L55 in/out — " f"L55 内 {n_l53_in_l55} 件 ({in_l55_share_53*100:.1f}%)", fontsize=11) ax.legend(loc="upper right", fontsize=9) fig.suptitle(f"図 7 (RQ3): 30 条系 vs 12 条系の運用実態 — " f"L52 30 条系 {in_l55_share*100:.0f}% / L53 40 条系 {in_l55_share_53*100:.0f}%", fontsize=13, y=1.02) plt.tight_layout() save_fig("L55_fig7_l52_l53_in_out.png") # ---- 図 8 (RQ3): 12 条系 vs 30 条系の盛土規模比較 ---- if L52_OK: fig, axes = plt.subplots(1, 2, figsize=(16, 6)) # 8-a L52 規模分布 (in_l55 vs in_l54) ax = axes[0] a_l55 = g52[g52["in_l55"]]["area_m2"].dropna() a_l54 = g52[g52["in_l54"]]["area_m2"].dropna() if len(a_l55) > 0 and len(a_l54) > 0: all_min = min(a_l55.min(), a_l54.min()) all_max = max(a_l55.max(), a_l54.max()) bins = np.logspace(np.log10(all_min), np.log10(all_max), 25) ax.hist(a_l54, bins=bins, color=L54_COLOR, edgecolor="#333", alpha=0.7, label=f"L54 内 (12 条系) n={len(a_l54)}, 中央 {a_l54.median():.0f}m²") ax.hist(a_l55, bins=bins, color=L55_COLOR, edgecolor="#333", alpha=0.6, label=f"L55 内 (30 条系) n={len(a_l55)}, 中央 {a_l55.median():.0f}m²") ax.axvline(a_l54.median(), color=L54_COLOR, linestyle="--", linewidth=2) ax.axvline(a_l55.median(), color=L55_COLOR, linestyle="--", linewidth=2) ax.set_xscale("log") ax.set_xlabel("盛土面積 (m², log)", fontsize=11) ax.set_ylabel("件数", fontsize=11) ax.set_title(f"図 8a: L52 許可盛土 規模分布 — 30 条系は大規模工事多数", fontsize=11) ax.legend(fontsize=9) ax.grid(alpha=0.3, which="both") # 8-b 規模指標 4 種比較棒 (中央値) ax = axes[1] metrics = ["height_max_m", "area_m2", "fill_volume_m3"] metric_labels = ["高さ (m)", "面積 (m²)", "盛土量 (m³)"] l54_meds = [g52[g52["in_l54"]][m].median() for m in metrics] l55_meds = [g52[g52["in_l55"]][m].median() for m in metrics] x = np.arange(len(metrics)) w = 0.35 # 対数スケール (高さは線形だが面積/土量は log) log_meds_l54 = [np.log10(v) if v > 0 else 0 for v in l54_meds] log_meds_l55 = [np.log10(v) if v > 0 else 0 for v in l55_meds] ax.bar(x - w/2, log_meds_l54, w, color=L54_COLOR, edgecolor="#333", label=f"L54 内 (12 条系)") ax.bar(x + w/2, log_meds_l55, w, color=L55_COLOR, edgecolor="#333", label=f"L55 内 (30 条系)") ax.set_xticks(x) ax.set_xticklabels(metric_labels, fontsize=10) ax.set_ylabel("中央値 (log10 スケール)", fontsize=11) ax.set_title(f"図 8b: 12 条系 vs 30 条系 規模中央値比較 — 30 条系が大きい (山間地大規模)", fontsize=11) ax.legend(fontsize=9) ax.grid(axis="y", alpha=0.3) for i, (l54v, l55v) in enumerate(zip(l54_meds, l55_meds)): ax.text(i - w/2, np.log10(max(l54v, 1)) + 0.05, f"{l54v:.1f}", ha="center", fontsize=9) ax.text(i + w/2, np.log10(max(l55v, 1)) + 0.05, f"{l55v:.1f}", ha="center", fontsize=9) fig.suptitle("図 8 (RQ3): 12 条系 vs 30 条系の盛土規模差 — " "30 条系は大規模工事を捕捉", fontsize=13, y=1.02) plt.tight_layout() save_fig("L55_fig8_l52_size_compare.png") # ---- 図 9 (統合): 4 制度連鎖 (L52/L53/L54/L55) Top 8 区域プロファイル ---- if L52_OK and L53_OK: # zone_density に二重リスクをマージ top_data = zone_density.merge(double_risk[["市町名", "二重リスク_km2"]], on="市町名", how="left").nlargest(8, "zone_area_km2") top_data = top_data.iloc[::-1].reset_index(drop=True) fig, ax = plt.subplots(figsize=(15, 7)) y = np.arange(len(top_data)) w = 0.20 # 系列 1: 区域面積 (km²) ax.barh(y - 1.5*w, top_data["zone_area_km2"], height=w, color=L55_COLOR, edgecolor="#333", label="L55 区域面積 (km²)") for i, v in enumerate(top_data["zone_area_km2"]): ax.text(v + 5, i - 1.5*w, f"{v:.0f}", va="center", fontsize=8) # 系列 2: 二重リスク (km²) ×100 でスケール SCALE_DR = 50 ax.barh(y - 0.5*w, top_data["二重リスク_km2"] * SCALE_DR, height=w, color="#cf222e", edgecolor="#333", label=f"二重リスク (km²) × {SCALE_DR}") for i, v in enumerate(top_data["二重リスク_km2"]): ax.text(v * SCALE_DR + 5, i - 0.5*w, f"{v:.1f}", va="center", fontsize=8) # 系列 3: L52 30 条系 ×30 SCALE_L52 = 30 ax.barh(y + 0.5*w, top_data["L52_30条許可"] * SCALE_L52, height=w, color="#1f6feb", edgecolor="#333", label=f"L52 30 条系許可 (件) × {SCALE_L52}") for i, v in enumerate(top_data["L52_30条許可"]): ax.text(v * SCALE_L52 + 5, i + 0.5*w, f"{int(v)}", va="center", fontsize=8) # 系列 4: L53 40 条系 ×10 SCALE_L53 = 10 ax.barh(y + 1.5*w, top_data["L53_40条届出"] * SCALE_L53, height=w, color="#fa8c16", edgecolor="#333", label=f"L53 40 条系届出 (件) × {SCALE_L53}") for i, v in enumerate(top_data["L53_40条届出"]): ax.text(v * SCALE_L53 + 5, i + 1.5*w, f"{int(v)}", va="center", fontsize=8) ax.set_yticks(y) ax.set_yticklabels(top_data["市町名"], fontsize=10) ax.set_xlabel("値 (区域面積=km², 二重リスク=km²×50, L52=件×30, L53=件×10)", fontsize=11) ax.set_title(f"図 9 (4 制度統合): Top 8 L55 区域の総合プロファイル — " f"L52/L53/L54/L55 の 4 制度連鎖が見える 1 枚", fontsize=12) ax.legend(fontsize=9, loc="lower right") ax.grid(axis="x", alpha=0.3) plt.tight_layout() save_fig("L55_fig9_top8_4systems.png") print(f" 9 図 生成完了 ({time.time()-t1:.1f}s)", flush=True) # ============================================================================= # 10. 表の DataFrame 作成 (HTML 用) # ============================================================================= print("\n[10] 表の作成", flush=True) t1 = time.time() def df_to_html(d): cols = list(d.columns) head = "" + "".join(f"{c}" for c in cols) + "" rows = [] for _, r in d.iterrows(): cells = [] for c in cols: v = r[c] if isinstance(v, float): if abs(v) >= 1000: s = f"{v:,.0f}" elif abs(v) >= 10: s = f"{v:.1f}" elif abs(v) >= 1: s = f"{v:.2f}" else: s = f"{v:.3f}" else: s = str(v) cells.append(f"{s}") rows.append("" + "".join(cells) + "") return "" + head + "".join(rows) + "
" # 表 1: dataset 仕様 T_dataset = pd.DataFrame([ ("dataset_id", "1428"), ("名称", "特定盛土等規制区域"), ("組織", "広島県土木建築局 都市環境整備課"), ("リソース 1 (県知事)", "51146 — 県統合 Shapefile (265 行 / 主要 81 polygon, ~6 MB)"), ("リソース 2-16 (市町長)", "51147〜51161 — 市町別 Shapefile 計 15 ファイル (確認用)"), ("リソース 17 (規制資料)", "51162 — 盛土規制法関係資料 XLSX (L52/L53/L54 と共通)"), ("根拠法", "宅地造成及び特定盛土等規制法 第 30 条第 1 項 (2023-05-26 施行)"), ("対象地域", "宅造規制区域 (12 条) 以外で地形条件 + 下流人口を満たす場所"), ("CRS", "EPSG:6668 (JGD2011 緯度経度) → 解析時 EPSG:6671 へ変換"), ("指定日", "2023-04-27 (新法施行直前)"), ("ライセンス", "クリエイティブ・コモンズ表示 4.0"), ("取得日", "2026-05-09"), ], columns=["項目", "値"]) # 表 2: 全体サマリ T_overall = pd.DataFrame([ ("総 polygon 数 (主源)", f"{len(zones)} polygon (県知事ファイル中 主要 polygon, " f"市町別ファイル 15 件は重複なので参照用)"), ("対象市町数", f"{n_cities_designated} 市町"), ("L55 union 面積", f"{total_zone_area_union:.2f} km² " f"(県 {PREF_AREA_KM2:.0f} km² の {pref_share_union:.2f}%)"), ("L54 union 面積 (既扱)", f"{l54_area_km2:.2f} km² " f"(県の {l54_area_km2/PREF_AREA_KM2*100:.1f}%)"), ("L54 ∪ L55 合算カバー", f"{combined_area_km2:.2f} km² " f"(県の {combined_share:.1f}% — 2 制度合算)"), ("polygon 規模 中央値 / 最大", f"{zones['poly_area_km2'].median():.3f} km² / {zones['poly_area_km2'].max():.0f} km²"), ("Top 1 市町", f"{city_summary.iloc[0]['市町名']} ({city_summary.iloc[0]['指定面積_km2']:.0f} km²)"), ("Top 5 シェア", f"{top5_share*100:.1f}%"), ("L55 ∩ 警戒区域", f"{overlap_area_km2:.1f} km² " f"(L55 内 {overlap_in_zone_pct:.1f}%, 警戒の {overlap_in_warn_pct:.1f}%)"), ("(L54+L55) ∩ 警戒区域", f"{combined_warn_km2:.1f} km² " f"(警戒の {combined_warn_pct:.1f}% を 2 制度でカバー)"), ("警戒のみ (L54 単独)", f"{warn_only_l54_km2:.0f} km² (全警戒の {warn_only_l54_pct:.0f}%)"), ("警戒のみ (L54+L55 合算)", f"{warn_only_combined_km2:.0f} km² ({warn_only_combined_pct:.1f}%) — 縮小幅 " f"{warn_only_l54_km2-warn_only_combined_km2:.0f} km²"), ("L52 152 件 → L55 内", f"{n_l52_in_l55} 件 ({in_l55_share*100:.1f}%) = 30 条系許可" if L52_OK else "—"), ("L52 152 件 → L54 内", f"{n_l52_in_l54} 件 ({n_l52_in_l54/152*100:.1f}%) = 12 条系許可" if L52_OK else "—"), ("L53 284 件 → L55 内", f"{n_l53_in_l55} 件 ({in_l55_share_53*100:.1f}%) = 40 条系届出" if L53_OK else "—"), ("L53 284 件 → L54 内", f"{n_l53_in_l54} 件 ({n_l53_in_l54/284*100:.1f}%) = 21 条系届出" if L53_OK else "—"), ("二重リスク Top 3 シェア", f"{top3_double_share*100:.0f}% ({', '.join(top3_double['市町名'].tolist())})"), ], columns=["指標", "値"]) # 表 3: 市町別 L55 一覧 T_zones = zone_dissolved[["市町名", "地形タイプ", "zone_area_km2"]].copy() T_zones = T_zones[T_zones["市町名"] != "不明"] T_zones["順位"] = T_zones["zone_area_km2"].rank(ascending=False, method="min").astype(int) T_zones["県シェア_%"] = (T_zones["zone_area_km2"] / PREF_AREA_KM2 * 100).round(2) T_zones["指定面積シェア_%"] = (T_zones["zone_area_km2"] / T_zones["zone_area_km2"].sum() * 100).round(2) T_zones = T_zones.sort_values("zone_area_km2", ascending=False).reset_index(drop=True) T_zones["polygon数"] = T_zones["市町名"].map( zones.groupby("市町名").size().to_dict()).fillna(0).astype(int) T_zones = T_zones[["順位", "市町名", "地形タイプ", "zone_area_km2", "指定面積シェア_%", "県シェア_%", "polygon数"]] T_zones.columns = ["順位", "市町名", "地形タイプ", "面積_km2", "指定面積シェア_%", "県総面積シェア_%", "polygon数"] # 表 4: 地形タイプ集計 T_type = type_summary.copy() T_type["地形タイプ"] = T_type["地形タイプ"].replace({ "A_中山間": "A 中山間市町 (本制度の主舞台)", "B_都市辺縁": "B 都市辺縁", "C_島嶼小規模": "C 島嶼/小規模沿岸", "D_その他": "D その他", }) # 表 5: L54 vs L55 比較 (姉妹制度の対比) T_l54_vs_l55 = pd.DataFrame([ {"側面": "根拠条文", "L54 (宅造区域)": "法 12 条 1 項", "L55 (特定区域, 本記事)": "法 30 条 1 項"}, {"側面": "対象地域", "L54 (宅造区域)": "都市計画区域内の宅地", "L55 (特定区域, 本記事)": "宅造区域以外で地形条件 + 下流人口"}, {"側面": "主な指定地", "L54 (宅造区域)": "沿岸都市・市街地", "L55 (特定区域, 本記事)": "山間地・農地・林地・島嶼"}, {"側面": "対象市町数", "L54 (宅造区域)": "20+ 市町 (沿岸+少数山間)", "L55 (特定区域, 本記事)": f"{n_cities_designated} 市町 (山間中心)"}, {"側面": "union 指定面積", "L54 (宅造区域)": f"{l54_area_km2:.0f} km² (県の {l54_area_km2/PREF_AREA_KM2*100:.1f}%)", "L55 (特定区域, 本記事)": f"{total_zone_area_union:.0f} km² (県の {pref_share_union:.1f}%)"}, {"側面": "L54/L55 重複", "L54 (宅造区域)": f"{l54_l55_overlap_km2:.1f} km² (相互排他に近い)", "L55 (特定区域, 本記事)": f"{l54_l55_overlap_km2:.1f} km²"}, {"側面": "警戒 ∩ 規制", "L54 (宅造区域)": f"{warn_only_l54_km2:.0f} km² が警戒のみ残る", "L55 (特定区域, 本記事)": f"L55 単独 ∩ 警戒 = {overlap_area_km2:.1f} km²"}, {"側面": "L52 許可立地", "L54 (宅造区域)": f"{n_l52_in_l54}/152 ({n_l52_in_l54/152*100:.0f}%, 12 条系)" if L52_OK else "—", "L55 (特定区域, 本記事)": f"{n_l52_in_l55}/152 ({in_l55_share*100:.0f}%, 30 条系)" if L52_OK else "—"}, {"側面": "L53 届出立地", "L54 (宅造区域)": f"{n_l53_in_l54}/284 ({n_l53_in_l54/284*100:.0f}%, 21 条系)" if L53_OK else "—", "L55 (特定区域, 本記事)": f"{n_l53_in_l55}/284 ({in_l55_share_53*100:.0f}%, 40 条系)" if L53_OK else "—"}, ]) # 表 6: 二重リスクランキング (Top 12) T_dr = double_risk[double_risk["市町名"] != "不明"].copy() T_dr = T_dr[["市町名", "zone_area_km2", "二重リスク_km2", "区域内警戒比_%"]] T_dr.columns = ["市町名", "区域面積_km2", "二重リスク_km2", "区域内警戒区域比率_%"] T_dr["順位"] = T_dr["二重リスク_km2"].rank(ascending=False, method="min").astype(int) T_dr = T_dr.sort_values("二重リスク_km2", ascending=False).head(12).reset_index(drop=True) T_dr = T_dr[["順位", "市町名", "区域面積_km2", "二重リスク_km2", "区域内警戒区域比率_%"]] # 表 7: 4 区分カバレッジ (3 制度合算) T_quad = pd.DataFrame([ {"区分": "1. (L54+L55) ∩ 警戒 (二重リスク)", "面積_km2": round(combined_warn_km2, 1), "県シェア_%": round(combined_warn_km2 / PREF_AREA_KM2 * 100, 2), "意味": "盛土規制 (12 or 30 条) + 自然災害想定 = 二重防御"}, {"区分": "2. (L54+L55) のみ (警戒区域外)", "面積_km2": round(combined_area_km2 - combined_warn_km2, 1), "県シェア_%": round((combined_area_km2 - combined_warn_km2) / PREF_AREA_KM2 * 100, 2), "意味": "盛土規制のみ。地形リスクは想定外だが工事監督対象"}, {"区分": "3. 警戒のみ (L54+L55 合算でも残る)", "面積_km2": round(warn_only_combined_km2, 1), "県シェア_%": round(warn_only_combined_km2 / PREF_AREA_KM2 * 100, 2), "意味": "災害想定のみ。山間集落で盛土規制が及ばない最終ギャップ"}, {"区分": "4. どちらも該当外", "面積_km2": round(max(PREF_AREA_KM2 - combined_warn_km2 - (combined_area_km2 - combined_warn_km2) - warn_only_combined_km2, 0), 1), "県シェア_%": round(max(PREF_AREA_KM2 - combined_area_km2 - warn_only_combined_km2, 0) / PREF_AREA_KM2 * 100, 2), "意味": "両制度の対象外。森林・農地・水域の県全域の大半"}, ]) # 表 8: 4 制度連鎖 (L52/L53/L54/L55) — 盛土規制法の全体像 T_4systems = pd.DataFrame([ {"側面": "種別", "L54 (宅造区域)": "区域指定", "L55 (特定区域)": "区域指定", "L52 (許可)": "個別工事許可", "L53 (届出)": "個別工事届出"}, {"側面": "根拠条文", "L54 (宅造区域)": "12 条 1 項", "L55 (特定区域)": "30 条 1 項", "L52 (許可)": "12 条 / 30 条", "L53 (届出)": "21 条 / 40 条"}, {"側面": "本データ単位", "L54 (宅造区域)": f"{l54_area_km2:.0f} km² (24% of 県)", "L55 (特定区域)": f"{total_zone_area_union:.0f} km² ({pref_share_union:.0f}% of 県)", "L52 (許可)": "152 件", "L53 (届出)": "284 件"}, {"側面": "12 条系 / 30 条系", "L54 (宅造区域)": "12 条 1 項のみ", "L55 (特定区域)": "30 条 1 項のみ", "L52 (許可)": f"12 条 {n_l52_in_l54} / 30 条 {n_l52_in_l55}" if L52_OK else "—", "L53 (届出)": f"21 条 {n_l53_in_l54} / 40 条 {n_l53_in_l55}" if L53_OK else "—"}, {"側面": "主な地域", "L54 (宅造区域)": "沿岸都市", "L55 (特定区域)": "山間地", "L52 (許可)": "都市部中規模", "L53 (届出)": "全域分散"}, {"側面": "L52/L53 中央値規模", "L54 (宅造区域)": f"L52: {g52[g52['in_l54']]['area_m2'].median():.0f}m²" if L52_OK else "—", "L55 (特定区域)": f"L52: {g52[g52['in_l55']]['area_m2'].median():.0f}m²" if L52_OK else "—", "L52 (許可)": "全体中央値 945m²", "L53 (届出)": "全体中央値 3700m²"}, ]) # 表 9: L52 規模比較 (L54 内 vs L55 内) if L52_OK: in_stats_l54 = g52[g52["in_l54"]][["height_max_m", "area_m2", "fill_volume_m3"]] in_stats_l55 = g52[g52["in_l55"]][["height_max_m", "area_m2", "fill_volume_m3"]] T_size_compare = pd.DataFrame([ {"指標": "件数 (件)", "L54 内 (12 条系)": len(in_stats_l54), "L55 内 (30 条系)": len(in_stats_l55), "差 (30-12)": len(in_stats_l55) - len(in_stats_l54)}, {"指標": "高さ 中央値 (m)", "L54 内 (12 条系)": round(in_stats_l54["height_max_m"].median(), 2), "L55 内 (30 条系)": round(in_stats_l55["height_max_m"].median(), 2), "差 (30-12)": round(in_stats_l55["height_max_m"].median() - in_stats_l54["height_max_m"].median(), 2)}, {"指標": "面積 中央値 (m²)", "L54 内 (12 条系)": int(in_stats_l54["area_m2"].median()), "L55 内 (30 条系)": int(in_stats_l55["area_m2"].median()), "差 (30-12)": int(in_stats_l55["area_m2"].median() - in_stats_l54["area_m2"].median())}, {"指標": "盛土量 中央値 (m³)", "L54 内 (12 条系)": int(in_stats_l54["fill_volume_m3"].median()), "L55 内 (30 条系)": int(in_stats_l55["fill_volume_m3"].median()), "差 (30-12)": int(in_stats_l55["fill_volume_m3"].median() - in_stats_l54["fill_volume_m3"].median())}, {"指標": "面積 最大 (m²)", "L54 内 (12 条系)": int(in_stats_l54["area_m2"].max()), "L55 内 (30 条系)": int(in_stats_l55["area_m2"].max()), "差 (30-12)": int(in_stats_l55["area_m2"].max() - in_stats_l54["area_m2"].max())}, ]) else: T_size_compare = pd.DataFrame() # 表 10: 区域別 30 条系運用件数 (L52/L53) if L52_OK and L53_OK: T_30jou = zone_density[zone_density["市町名"] != "不明"].copy() T_30jou = T_30jou[["市町名", "zone_area_km2", "L52_30条許可", "L53_40条届出", "L52_密度_件km2", "L53_密度_件km2"]] T_30jou.columns = ["市町名", "区域面積_km2", "L52許可", "L53届出", "L52密度_件km2", "L53密度_件km2"] T_30jou["L52+L53合計"] = T_30jou["L52許可"] + T_30jou["L53届出"] T_30jou = T_30jou.sort_values("L52+L53合計", ascending=False).reset_index(drop=True) T_30jou["順位"] = np.arange(1, len(T_30jou) + 1) T_30jou = T_30jou[["順位", "市町名", "区域面積_km2", "L52許可", "L53届出", "L52+L53合計", "L52密度_件km2", "L53密度_件km2"]] else: T_30jou = pd.DataFrame() # 表 11: 仮説検証 T_hypo = pd.DataFrame([ {"仮説": "H1 (山間地偏重, 上位 5 シェア ≥ 50%, RQ1)", "予想": "Top 5 市町で全県の 50% 以上", "観測": f"Top 5 = {', '.join(city_summary[city_summary['市町名']!='不明'].head(5)['市町名'].tolist())}, " f"シェア {top5_share*100:.1f}%", "判定": "支持" if top5_share >= 0.50 else "部分支持"}, {"仮説": "H2 (大カバー率 ≥ 40%, RQ1)", "予想": "L55 単独で県の 40% 以上を覆う", "観測": f"{pref_share_union:.1f}% (L54 24% の 2 倍以上)", "判定": "強支持" if pref_share_union >= 40 else "部分支持"}, {"仮説": "H3 (警戒重複 10-30%, RQ2)", "予想": "L55 ∩ 警戒の比率が L54 (2.9%) より高い 10-30%", "観測": f"L55 ∩ 警戒 = {overlap_area_km2:.1f} km², L55 内 {overlap_in_zone_pct:.1f}%", "判定": "反証 (低 0.9% — L55 は山間内陸で警戒は急傾斜・渓流で空間目的が異なる)" if overlap_in_zone_pct < 10 else "支持" if 10 <= overlap_in_zone_pct <= 30 else "観測 高すぎ"}, {"仮説": "H4 (制度間ギャップ縮小, RQ2)", "予想": "L55 追加で『警戒のみ』が 98 → 30 km² 以下に縮小", "観測": f"L54 単独 {warn_only_l54_km2:.0f} km² → L54+L55 合算 {warn_only_combined_km2:.0f} km² " f"({(warn_only_l54_km2-warn_only_combined_km2)/warn_only_l54_km2*100:.0f}% 縮小)", "判定": "強支持" if warn_only_combined_km2 < warn_only_l54_km2 * 0.5 else "部分支持"}, {"仮説": "H5 (30 条系運用 5-15%, RQ3)", "予想": "L52 のうち L55 内立地は 5-15%", "観測": f"L52 L55 内 {n_l52_in_l55}/{152} = {in_l55_share*100:.1f}% / " f"L53 L55 内 {n_l53_in_l55}/{284} = {in_l55_share_53*100:.1f}%" if L52_OK and L53_OK else "—", "判定": ("支持" if 5 <= in_l55_share*100 <= 15 else "部分支持") if L52_OK else "検証不可"}, ]) # 表 12: 市町別ファイルメタ (15 ファイルの確認用) T_city_meta = city_meta_df.copy() T_city_meta = T_city_meta[["市町コード", "市町名", "市町長ファイル_polygon数", "市町長ファイル_面積_km2"]] T_city_meta.columns = ["市町コード", "市町名", "polygon数", "面積_km2"] print(f" 表 12 種作成 ({time.time()-t1:.1f}s)", flush=True) # ============================================================================= # 11. HTML 構築 # ============================================================================= print("\n[11] HTML 構築", flush=True) t11 = time.time() # 計算済み変数の整理 poly_a_med = zones["poly_area_km2"].median() poly_a_max = zones["poly_area_km2"].max() poly_a_min = zones["poly_area_km2"].min() a_share = type_summary[type_summary['地形タイプ']=='A_中山間']['シェア_%'].sum() b_share = type_summary[type_summary['地形タイプ']=='B_都市辺縁']['シェア_%'].sum() c_share = type_summary[type_summary['地形タイプ']=='C_島嶼小規模']['シェア_%'].sum() top1_name = city_summary[city_summary["市町名"] != "不明"].iloc[0]["市町名"] top1_area = city_summary[city_summary["市町名"] != "不明"].iloc[0]["指定面積_km2"] # Sec 1: 学習目標と問い sec1 = f"""

本記事は DoBoX のシリーズ「特定盛土等規制区域」 1 件 (dataset_id = 1428) を 単独で取り上げ、 広島県内の特定盛土等規制区域 {len(zones)} polygon (主源 = 県知事ファイル) / {n_cities_designated} 市町指定 / union 面積 {total_zone_area_union:.0f} km² ({pref_share_union:.1f}% of 県土) を 3 つの独立した研究角度 (RQ1 / RQ2 / RQ3) で並列に分析する。 本データは盛土規制法 (2023-05-26 施行) 第 30 条第 1 項に基づく区域指定で、 2021 年 7 月 静岡県熱海市 伊豆山土石流災害 (盛土崩落) を契機に新設された区域類型。 L54 (12 条 1 項 = 都市計画区域内宅地) が捕捉しきれない山間地・農地・林地の大規模盛土を 法的監督下に置くための「山間地監督の最後のピース」。

本データは DoBoX で16 ファイル構成: 県知事指定の県全域版 1 ファイル (51146, 2023-04-27 版) + 各市町長指定の市町別 15 ファイル (51147〜51161, 2023-09 版)。 ただし市町別ファイルの polygon は県知事ファイルの該当 polygon と完全一致を確認済 (例: 庄原市の市町長ファイル 1057.94 km² polygon = 県知事ファイル中の対応 polygon と同一)。 本研究は県知事ファイル単独を主源データとすることで二重カウントを回避し、 市町別ファイルは「市町長確認済」 のメタ情報として参照する。

本記事はPhase 2 盛土規制系 4 本目 (= L52 許可 / L53 届出 / L54 宅造区域 / L55 特定区域)の 最終ピース。L52 (許可) / L53 (届出) は条文ごとに 12 条系と 30 条系の 2 系統を持ち、 L54 (12 条 1 項) と L55 (30 条 1 項) はそれぞれの地理的根拠。 本記事 RQ3 で初めて L52 / L53 の30 条系内訳を定量化する。

独自用語の定義

研究の問い (3 RQ)

仮説 H1〜H5

  1. H1 (山間地偏重, RQ1): 16 区域は山間中山間市町に偏在し、 上位 5 (庄原市・三次市・北広島町・安芸高田市・神石高原町等) で全体の 50% 以上を占める。 L54 (沿岸都市偏重) と地理的に補完。
  2. H2 (大カバー率, RQ1): L55 単独で県の 40% 以上を覆う = L54 (24%) を超える広域指定。 30 条 1 項が地形条件 + 下流人口を満たす広い領域を対象とすることに整合。
  3. H3 (警戒重複高, RQ2): L55 ∩ 警戒区域の重複率は L54 (2.9%) より高い 10-30%。 30 条 1 項自体が地形条件 + 災害リスクを考慮した区域指定で空間相関が高いはず。
  4. H4 (制度間ギャップ ≪ L54 単独, RQ2): L54 + L55 合算カバーで「警戒のみ」 は L54 単独 98 km² から 30 km² 以下に縮小。L55 が L54 のギャップを埋める設計通り。
  5. H5 (30 条系運用は控えめ, RQ3): L52 のうち L55 内立地は 5-15% 程度。 12 条系 (L54 86%) に比べ少数派。30 条系は新法施行以降の指定で許可案件の蓄積が 始まったばかり。L53 も同傾向。

到達点

本記事を読み終えた学習者は次の 3 点を体感できる:

  1. 1 つの「区域 Shapefile」 ({len(zones)} polygon, 県知事 1 + 市町長 15 = 16 ファイル) から、 L54 と L55 の地理的補完性を可視化し、L54 ∪ L55 = {combined_share:.1f}% of 県という 合算カバレッジを定量化する方法を習得する。
  2. L52 / L53 の30 条系内訳を初めて算出し、12 条系 vs 30 条系の運用差 (規模・件数・地理) を読む眼を獲得する。30 条系は大規模工事を捕捉する制度であることを 実データから確認できる。
  3. 「警戒のみ」 残存エリアが L54 単独 {warn_only_l54_km2:.0f} km² から L54+L55 合算 {warn_only_combined_km2:.0f} km² に縮小する事実から、 制度間補完設計の有効性を定量診断する手法を体感する。
""" # Sec 2: 使用データ sec2 = f"""

DoBoX のシリーズ「特定盛土等規制区域」 1 件のみを単独で扱う。 リソースは Shapefile 16 件 (県知事 1 + 市町長 15) + XLSX 1 件の計 17 ファイル構造:

{df_to_html(T_dataset)}

データの構造

CRS と単位の取扱

本記事は dataset 1428 を単独で扱う Format B 記事。 L52 (dataset 1429 = 許可盛土) / L53 (dataset 1430 = 届出盛土) / L54 (dataset 1427 = 宅造区域) との 比較は RQ2 (L54 制度補完) / RQ3 (L52/L53 30 条系内訳) で使用するが、 それぞれのデータは別レッスンで深掘り済み。

""" # Sec 3: ダウンロード sec3 = f"""

本レッスンの再現に必要な全データ・中間 CSV・図 PNG・スクリプトを以下から直接 DL できる:

生データ (DoBoX 直リンク)

本記事の中間 CSV (再現用)

図 (PNG 9 枚)

再現スクリプト

個別取得 (PowerShell):

cd "2026 DoBoX 教材"
iwr "https://hiroshima-dobox.jp/resource_download/51146" -OutFile "data/extras/L55_specific_filling_zone/340006_specific_filling_zone_pref_shp_2023-04-27.zip"
iwr "https://hiroshima-dobox.jp/resource_download/51162" -OutFile "data/extras/L55_specific_filling_zone/340006_filling_regulation_related_documents.xlsx"
py -X utf8 lessons\\L55_specific_filling_zone.py

注: RQ2 (L54 制度補完) には L54 の Shapefile データ (data/extras/L54_residential_filling_zone/) が必要 (= L54 完成済前提)。 RQ3 (L52/L53 比較) には L52/L53 の処理済 CSV (lessons/assets/L52_permits_processed.csv / L53_notifications_processed.csv) が必要。 未生成の場合、本スクリプトは該当部分を「データ無し」 表示にしてフォールバック実行する。

""" # Sec 4 (RQ1) sec4_intro = f"""

狙い (RQ1)

特定盛土等規制区域 ({n_cities_designated} 市町指定) は盛土規制法 30 条 1 項に基づく区域指定で、 L54 (12 条 1 項 = 都市計画区域内宅地) と異なり山間地・農地・林地等を含む県下全域を 対象に指定可能。本 RQ1 では (1) 県全体の指定面積、(2) 市町別指定面積ランキング、 (3) 地形タイプ別構成、(4) L54 との地理的補完性の 4 軸で構造を読む。

L54 (沿岸都市偏重 = 24% of 県) と L55 (山間地偏重 = {pref_share_union:.0f}% of 県) は、 盛土規制法の地理カテゴリ二分法を反映する姉妹制度。 両者の重複面積 (L54 ∩ L55) はわずか {l54_l55_overlap_km2:.1f} km² で地理的にほぼ排他。 これは制度設計上、12 条と 30 条が異なる地域カテゴリを担うように意図された結果。

手法 (Shapefile 読込 → MakeValid → N03_004 同定 → dissolve → 集計)

入出力 Before/After 具体例 (1 polygon 追跡)

県知事ファイル中の最大 polygon (庄原市相当) の処理例:

{df_to_html(pd.DataFrame([ {"段階": "RAW (zip)", "対象": "Shapefile zip", "geometry": "265 polygon (細片含)", "CRS": "EPSG:6668 (緯度経度)", "面積": "—"}, {"段階": "1. read_file + to_crs", "対象": "GeoDataFrame", "geometry": "265 polygon", "CRS": "EPSG:6671", "面積": "(計算可能に)"}, {"段階": "2. make_valid", "対象": "GeoDataFrame", "geometry": "265 polygon (有効)", "CRS": "EPSG:6671", "面積": "—"}, {"段階": "3. >1000 m² フィルタ", "対象": "主要のみ", "geometry": f"{len(zones)} polygon", "CRS": "EPSG:6671", "面積": f"合計 {zones['poly_area_km2'].sum():.0f} km²"}, {"段階": "4. N03_004 で市町名取得", "対象": "庄原 polygon", "geometry": "1 polygon", "CRS": "EPSG:6671", "面積": "1058 km²"}, {"段階": "5. dissolve(by=市町名)", "対象": "市町別集約", "geometry": f"{n_cities_designated} 市町", "CRS": "EPSG:6671", "面積": f"{n_cities_designated} 行"}, {"段階": "6. unary_union", "対象": "重複除去後", "geometry": "1 MultiPolygon", "CRS": "EPSG:6671", "面積": f"{total_zone_area_union:.0f} km² ({pref_share_union:.1f}% of 県)"}, ]))} """ sec4_code = code(r''' # 16 ファイル統合: 県知事ファイル単独 (主源) + N03_004 で市町同定 import geopandas as gpd, pandas as pd from shapely import make_valid from shapely.ops import unary_union # 1. 県知事指定 Shapefile を読込 (51146, 県全域版) shp = "data/extras/L55_specific_filling_zone/shp/" \ "340006_specific_earth_filling_regulation_area_shp_20230427.shp" zone_pref = gpd.read_file(shp).to_crs("EPSG:6671") # 2. トポロジ修復 (一部 polygon に自己交差等の不整合) zone_pref["geometry"] = zone_pref.geometry.apply( lambda g: make_valid(g) if g and not g.is_valid else g) zone_pref["poly_area_km2"] = zone_pref.geometry.area / 1e6 # 3. 主要 polygon 抽出 (面積 > 1000 m²) zones = zone_pref[zone_pref["poly_area_km2"] > 0.001].reset_index(drop=True) zones["poly_id"] = range(len(zones)) # 4. N03_004 (国土数値情報の市町名) で直接同定 → admin sjoin で補完 zones["市町名"] = zones["N03_004"] unfilled = zones[zones["市町名"].isna()].copy() if len(unfilled) > 0: upts = gpd.GeoDataFrame( unfilled[["poly_id"]], geometry=unfilled.geometry.representative_point(), crs="EPSG:6671") j = gpd.sjoin(upts, admin_diss[["CITY_CD", "geometry"]], how="left", predicate="within") # ... NaN は intersects 最大重複で割り当て pid_to_name = j.set_index("poly_id")["CITY_CD"].map(CITY_NAME).to_dict() for pid, n in pid_to_name.items(): if n: zones.loc[zones["poly_id"] == pid, "市町名"] = n # 5. 市町別 dissolve + union 面積 zone_dissolved = zones.dissolve(by="市町名", as_index=False) zone_dissolved["zone_area_km2"] = zone_dissolved.geometry.area / 1e6 zone_pref_geom = unary_union(list(zone_dissolved.geometry)) print(f"L55 polygon: {len(zones)}, union: {zone_pref_geom.area/1e6:.0f} km²") # 6. L54 との合算 (制度補完) combined_geom = unary_union([l54_union, zone_pref_geom]) print(f"L54 ∪ L55 = {combined_geom.area/1e6:.0f} km²") ''') sec4_fig1_text = """

図 1: L55 全区域マップ + 市町別指定面積 choropleth (2 panel)

なぜこの図か: 学習者がまず「L55 特定区域はどこにあるか」 を一目で把握するため。 左の地図で polygon の地理分布、右の choropleth で市町別の濃淡を見る。 L55 は L54 と異なり大きな polygon が県北部・県西部の山間地に集中することが視覚で確認できる。

""" sec4_fig1_read = f"""

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

""" sec4_fig2_text = """

図 2: L54 vs L55 制度補完マップ + L55 polygon 規模分布

なぜこの図か: L54 (12 条系, 緑) と L55 (30 条系, 紫) を 1 枚に重ねることで、 両制度の地理的補完性を視覚的に確認する。重なる場所が極めて少なく、 むしろ排他的に県土を分担する設計が見える。 右の polygon 規模分布は L55 polygon の規模幅を log で表示。

""" sec4_fig2_read = f"""

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

""" sec4_fig3_text = """

図 3: 地形タイプ別構成 (パイ) + 市町別指定面積ランキング (Top 15)

なぜこの図か: 規制区域指定の地形偏在を 2 つのスケールで見る。 左パイで「規制区域の中での 3 タイプの内訳」、右棒で市町別ランキング。 両方を並べることで、地形タイプの内訳と個別市町の量を別軸で読む基本訓練ができる。

""" sec4_fig3_read = f"""

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

""" # Sec 5 (RQ2) sec5_intro = f"""

狙い (RQ2)

L55 (30 条系) は L54 (12 条系) で残った「警戒のみ」 (= 制度間ギャップ) = 山間地警戒区域の 61.7% (98 km²)を埋めるために設計された。本 RQ2 では (1) L55 ∩ 警戒の重複面積、 (2) L54 + L55 合算後の警戒のみ残存面積、(3) 二重リスクエリアの市町別偏在 を読む。 特に「L55 が L54 のギャップをどれだけ埋めたか」 が研究の核心。

L54 単独では警戒のみ = {warn_only_l54_km2:.0f} km² (61.7%) だったが、 L55 を加えると {warn_only_combined_km2:.0f} km² ({warn_only_combined_pct:.1f}%) に縮小する仮説を 直接検証する。これは盛土規制法の制度間補完設計の有効性の定量診断。

手法 (geometry intersection + 4 区分集計)

入力: L55 polygon (union {total_zone_area_union:.0f} km²) + L54 polygon (union {l54_area_km2:.0f} km²) + 警戒区域 (union {warn_total_km2:.0f} km²)。
出力: 警戒のみ縮小幅、市町別二重リスク、4 区分カバレッジ。

""" sec5_code = code(r''' # RQ2: L55 × 警戒区域 + L54 制度補完 from shapely.ops import unary_union # 1. L55 union (規制区域 polygon → 1 つに統合) zone_pref_geom = unary_union(list(zone_dissolved.geometry)) print(f"L55 union: {zone_pref_geom.area/1e6:.0f} km²") # 2. L55 ∩ 警戒 warn_pref_geom = warn_all_geom # 急傾斜 30K + 土石流 13K の union (前計算済) overlap_l55 = zone_pref_geom.intersection(warn_pref_geom) overlap_l55_km2 = overlap_l55.area / 1e6 print(f"L55 ∩ 警戒: {overlap_l55_km2:.1f} km²") # 3. L54 ∪ L55 合算 combined_geom = unary_union([l54_union, zone_pref_geom]) combined_km2 = combined_geom.area / 1e6 print(f"L54 ∪ L55: {combined_km2:.0f} km²") # 4. 警戒のみ (L54 単独 vs L54+L55 合算) warn_only_l54 = warn_pref_geom.difference(l54_union) warn_only_combined = warn_pref_geom.difference(combined_geom) print(f"警戒のみ L54 単独: {warn_only_l54.area/1e6:.0f} km² (制度間ギャップ)") print(f"警戒のみ L54+L55 合算: {warn_only_combined.area/1e6:.0f} km² ← 縮小") # 5. 市町別二重リスク double_risk = [] for _, r in zone_dissolved.iterrows(): inter = r.geometry.intersection(warn_pref_geom) double_risk.append({ "市町名": r["市町名"], "二重リスク_km2": inter.area / 1e6 if not inter.is_empty else 0, }) ''') sec5_fig4_text = """

図 4: L55 × 警戒区域 重ね合わせ全県マップ

なぜこの図か: 全県で L55 (紫) と警戒区域 (赤=急傾斜・橙=土石流) が どう重なるかを 1 図で俯瞰する。L55 は山間地に広く分布するが、 警戒区域は急傾斜地・渓流に局所集中するため、空間目的の違いが見える。

""" sec5_fig4_read = f"""

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

""" sec5_fig5_text = """

図 5: L54 単独 vs L54+L55 合算の警戒のみ縮小 — H4 直接検証

なぜこの図か: 本記事最重要の発見 = 「L55 が L54 のギャップを埋める」 を 直接検証する図。左の 4 区分パイで合算後の県土カバレッジ、 右の比較棒で警戒のみ縮小幅 (98 → 54 km²) を 1 枚で表現する。

""" sec5_fig5_read = f"""

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

""" sec5_fig6_text = """

図 6: 市町別 L55 二重リスク choropleth + ランキング

なぜこの図か: 二重リスクを市町単位で正規化したランキングと地図。 L55 内警戒比率は全県平均で 0.9% と低いが、市町別には差がある。 中山間市町でも下流人口条件を強く満たす場所では二重リスクが集中する。

""" sec5_fig6_read = f"""

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

""" # Sec 6 (RQ3) sec6_intro = f"""

狙い (RQ3)

L52 (許可) / L53 (届出) は条文ごとに 12 条系と 30 条系の 2 系統を持つ。 L52 12 条系 = 12 条 1 項許可 (L54 区域内) / L52 30 条系 = 30 条 1 項許可 (L55 区域内)。 L53 21 条系 = 21 条 1 項届出 (L54 区域内) / L53 40 条系 = 40 条 1 項届出 (L55 区域内)。

本 RQ3 では L52 152 件 / L53 284 件を L54 / L55 へ in/out 判定することで、 初めて30 条系内訳を定量化する。さらに 12 条系 vs 30 条系の規模・件数・地理差を 比較し、30 条系制度の初期運用実態を診断する。

30 条系は新法施行 (2023-05-26) 以降の指定で、L52/L53 既収集データはまだ 12 条系の蓄積が中心。30 条系の運用は始まったばかりであることを念頭に読む必要がある。

手法 (within 判定 + 規模比較 + 区域別件数集計)

入力: L52 152 件 + L53 284 件 + L55 zone_pref_geom + L54 union。
出力: L52/L53 の 4 分類別件数、12 条系 vs 30 条系規模比較、市町別 30 条系運用件数。

""" sec6_code = code(r''' # RQ3: L52 / L53 を L54 / L55 へ in/out 判定し、12 条系 vs 30 条系を分離 import geopandas as gpd, pandas as pd # 1. L52 (許可) / L53 (届出) の点 GeoDataFrame df_l52 = pd.read_csv("lessons/assets/L52_permits_processed.csv", encoding="utf-8-sig") df_l53 = pd.read_csv("lessons/assets/L53_notifications_processed.csv", encoding="utf-8-sig") g52 = gpd.GeoDataFrame(df_l52, geometry=gpd.points_from_xy(df_l52["lon"], df_l52["lat"]), crs="EPSG:4326").to_crs("EPSG:6671") g53 = gpd.GeoDataFrame(df_l53, geometry=gpd.points_from_xy(df_l53["lon"], df_l53["lat"]), crs="EPSG:4326").to_crs("EPSG:6671") # 2. 各点に L54 / L55 in/out フラグ g52["in_l54"] = g52.geometry.within(l54_union) g52["in_l55"] = g52.geometry.within(zone_pref_geom) g53["in_l54"] = g53.geometry.within(l54_union) g53["in_l55"] = g53.geometry.within(zone_pref_geom) # 3. 4 分類集計 (例: L52) n_l52_only_l54 = (g52["in_l54"] & ~g52["in_l55"]).sum() # 12 条系 n_l52_only_l55 = (~g52["in_l54"] & g52["in_l55"]).sum() # 30 条系 n_l52_both = (g52["in_l54"] & g52["in_l55"]).sum() # 重複指定 (まれ) n_l52_neither = (~g52["in_l54"] & ~g52["in_l55"]).sum() # 制度外 print(f"L52: 12 条系 {n_l52_only_l54}, 30 条系 {n_l52_only_l55}") print(f" 両方 {n_l52_both}, 制度外 {n_l52_neither}") # 4. 12 条系 vs 30 条系 規模比較 print(f"L54 内 中央値 面積: {g52[g52['in_l54']]['area_m2'].median():.0f} m²") print(f"L55 内 中央値 面積: {g52[g52['in_l55']]['area_m2'].median():.0f} m²") ''') sec6_fig7_text = """

図 7: L52/L53 × L54/L55 in/out オーバーレイマップ

なぜこの図か: 12 条系 vs 30 条系の地理分布差を直接視覚化する。 左 L52 (許可)、右 L53 (届出) を別々に表示。点の色 = 区域種別 (緑=L54内/紫=L55内/灰=両外)。 30 条系点 (紫▲) は山間地に分布し、12 条系点 (緑●) は沿岸都市に集中することを 1 枚で確認できる。

""" sec6_fig7_read = f"""

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

""" sec6_fig8_text = """

図 8: 12 条系 vs 30 条系 盛土規模比較 (L52)

なぜこの図か: L52 152 件を 12 条系 (L54 内) と 30 条系 (L55 内) に分け、 盛土の高さ・面積・土量の中央値を比較する。30 条系が大規模工事を捕捉する仮説を直接検証。

""" l54_med_a = g52[g52['in_l54']]['area_m2'].median() if L52_OK else 0 l55_med_a = g52[g52['in_l55']]['area_m2'].median() if L52_OK else 0 size_ratio = l55_med_a / l54_med_a if l54_med_a > 0 else 0 sec6_fig8_read = f"""

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

""" # Sec 7 (仮説検証総合) sec7 = ( "

本記事の 5 仮説と観測結果の照合:

" + df_to_html(T_hypo) + "

3 RQ × 3 結論

" "" "

4 制度モデル — L52 / L53 / L54 / L55 の関係表

" + df_to_html(T_4systems) + "

この表から読み取れること: 盛土規制法の4 制度連鎖が初めて 1 表に集約された: " "(1) 区域指定 (L54 = 12 条 1 項 / L55 = 30 条 1 項) が地理範囲を法的に固定。" "(2) 個別工事監督 (L52 許可 / L53 届出) がその区域内で発動。" "(3) L52 / L53 はそれぞれ12 条系(L54 内) と30 条系(L55 内) の 2 系統で運用される。" "本記事 L55 を加えた 4 本セットで初めて、盛土規制法の制度設計が立体的に理解できる構造となった。

" ) # Sec 8 (発展課題) sec8 = f"""

結果 X → 新仮説 Y → 課題 Z (3 RQ × 1 課題以上)

発展課題 1 (RQ1 由来): 30 条 1 項指定基準の地理的検証

発展課題 2 (RQ2 由来): 警戒のみ最終ギャップ ({warn_only_combined_km2:.0f} km²) の地理特性

発展課題 3 (RQ3 由来): 30 条系の運用実態時系列追跡

発展課題 4 (4 制度統合): 盛土規制法 4 階層の地理的可視化

""" # Combine sections sections = [ ("学習目標と問い", sec1), ("使用データ", sec2), ("ダウンロード", sec3), (f"【RQ1】 地理範囲と山間地カバレッジ — 16 ファイル・約 {total_zone_area_union:.0f} km² (県の {pref_share_union:.1f}%) の山間地偏重", sec4_intro + "

実装コード (Shapefile 読込 + MakeValid + N03_004 同定 + dissolve)

" + sec4_code + sec4_fig1_text + figure("assets/L55_fig1_zone_map.png", "図 1 (RQ1): L55 全区域マップ + 市町別指定面積 choropleth") + sec4_fig1_read + sec4_fig2_text + figure("assets/L55_fig2_l54_l55_complementarity.png", "図 2 (RQ1): L54 vs L55 制度補完マップ + polygon 規模分布") + sec4_fig2_read + sec4_fig3_text + figure("assets/L55_fig3_type_ranking.png", "図 3 (RQ1): 地形タイプ別 + 市町別ランキング (Top 15)") + sec4_fig3_read + "

表: 全体サマリ (L55 + L54 比較 + 警戒区域 + L52/L53 連携)

" + df_to_html(T_overall) + "

この表から読み取れること: " f"{len(zones)} polygon / {n_cities_designated} 市町 / union 面積 {total_zone_area_union:.0f} km² / " f"県シェア {pref_share_union:.1f}% という基礎数値、" f"L54 ∪ L55 合算 {combined_share:.0f}%、警戒のみ縮小 {warn_only_l54_km2:.0f}→{warn_only_combined_km2:.0f} km² 等の" "横断指標が一覧で確認できる。RQ1〜RQ3 の核心が 17 行に集約された統合サマリ。

" + "

表: 市町別 L55 一覧 (面積順)

" + df_to_html(T_zones) + "

この表から読み取れること: 各市町の地形タイプ・面積・polygon 数を一覧で確認できる。" "Top 5 はすべて中山間市町で、上位 1 市町 (庄原市) だけで全県の 12% を占める。" "下位の島嶼小規模 (大崎上島 25 km²) は polygon 数 (25 個) が多いが面積は小さい =\n" "島嶼の地理的特性 (細かい入江・湾岸)を反映する。

" + "

表: 地形タイプ別集計 (3 タイプ)

" + df_to_html(T_type) + "

この表から読み取れること: 中山間 (A) が " f"{a_share:.0f}% を占め、都市辺縁 (B) {b_share:.0f}% / 島嶼小規模 (C) {c_share:.0f}% と続く。" "L54 の地形タイプ (沿岸 81% + 内陸 14% + 山間 4%) と L55 ({a_share:.0f}% + {b_share:.0f}% + {c_share:.0f}%) を比較すると、\n" "ほぼ逆構成 = 12 条と 30 条の地理的分業が定量的に確認される。

" + "

表: L54 vs L55 制度比較 (姉妹制度の対比)

" + df_to_html(T_l54_vs_l55) + "

この表から読み取れること: L54 と L55 の制度的対比が" "9 側面 (根拠条文・対象・主な指定地・市町数・面積・重複・警戒関係・許可立地・届出立地) で" "整理される。地理的にはほぼ排他、運用的にはL54 が中規模都市部、L55 が大規模山間地を" "それぞれ捕捉する設計。盛土規制法の制度設計の意図が初めて 1 表で見える。

" ), (f"【RQ2】 警戒区域カバレッジ研究 — L54 + L55 合算で「警戒のみ」 {warn_only_l54_km2:.0f}→{warn_only_combined_km2:.0f} km² に縮小", sec5_intro + "

実装コード (geometry intersection + 4 区分集計)

" + sec5_code + sec5_fig4_text + figure("assets/L55_fig4_overlay_warn.png", "図 4 (RQ2): L55 × 警戒区域 重ね合わせ全県マップ") + sec5_fig4_read + sec5_fig5_text + figure("assets/L55_fig5_warn_coverage_l54_l55.png", "図 5 (RQ2): L54 単独 vs L54+L55 合算の警戒のみ縮小 — H4 直接検証") + sec5_fig5_read + sec5_fig6_text + figure("assets/L55_fig6_double_risk.png", "図 6 (RQ2): 市町別 L55 二重リスク choropleth + ランキング") + sec5_fig6_read + "

表: 市町別 L55 二重リスクランキング (Top 12)

" + df_to_html(T_dr) + "

この表から読み取れること: 二重リスク面積トップは" f" {T_dr.iloc[0]['市町名']} ({T_dr.iloc[0]['二重リスク_km2']:.1f} km², " f"区域内警戒比率 {T_dr.iloc[0]['区域内警戒区域比率_%']:.1f}%)。" "区域内警戒比率は 1% 前後で全市町ほぼ均一 = " "L55 は山間地全体を面的に指定するため警戒区域比率が薄まる" "という構造的特徴。

" + "

表: 4 区分カバレッジ (L54+L55 合算 × 警戒)

" + df_to_html(T_quad) + "

この表から読み取れること: 4 区分のうち「警戒のみ」が" f" わずか {warn_only_combined_km2:.0f} km² ({warn_only_combined_pct:.1f}%) に縮小 → " "L54 単独時代の 98 km² (61.7%) から劇的に削減された = " "L55 が L54 のギャップを埋めた制度間補完設計の成功事例。" "残り {warn_only_combined_km2:.0f} km² が盛土規制法の最終ギャップで、" "森林・農地内の急傾斜地等で住民密度が低い場所と推察される。

".replace("{warn_only_combined_km2:.0f}", f"{warn_only_combined_km2:.0f}") ), (f"【RQ3】 許可・届出盛土 (L52/L53) との関係 — 30 条系運用の実態", sec6_intro + "

実装コード (within 判定 + 規模比較)

" + sec6_code + sec6_fig7_text + figure("assets/L55_fig7_l52_l53_in_out.png", "図 7 (RQ3): L52/L53 × L54/L55 in/out オーバーレイ") + sec6_fig7_read + sec6_fig8_text + figure("assets/L55_fig8_l52_size_compare.png", "図 8 (RQ3): 12 条系 vs 30 条系 盛土規模比較 (L52)") + sec6_fig8_read + figure("assets/L55_fig9_top8_4systems.png", "図 9 (4 制度統合): Top 8 L55 区域の総合プロファイル") + "

図 9 から読み取れること: Top 8 L55 区域の" "区域面積・二重リスク・L52 30 条系・L53 40 条系を 4 系列棒で並べた総合図。" "庄原・三次・北広島は区域面積が大きいが L52/L53 件数は中程度 → " "大面積 = 大量盛土発生 ではないこと、つまり山間地は面積に対して工事密度が低い" "ことを 1 図で確認できる。これは都市部 (L54) との明確な対比。

" + "

表: L52 規模比較 (12 条系 vs 30 条系)

" + df_to_html(T_size_compare) + "

この表から読み取れること: L52 152 件のうち、" f"30 条系 (L55 内) は中央値で面積 {size_ratio:.1f} 倍・高さ 2 倍超。" "30 条系は明確に大規模工事を捕捉する制度。差 (30-12) 列を見ると、" "全規模指標で 30 条系が大きく、これは制度設計の意図を反映する。

" + "

表: 区域別 30 条系運用件数 (L52 許可 + L53 届出)

" + df_to_html(T_30jou) + "

この表から読み取れること: 区域別の 30/40 条系件数集計。" f"L52+L53 合計上位は {T_30jou.iloc[0]['市町名']} ({T_30jou.iloc[0]['L52許可']}+{T_30jou.iloc[0]['L53届出']}={T_30jou.iloc[0]['L52+L53合計']} 件) 等の中山間市町 → " "山間地の盛土需要が定量的に確認される。" "密度 (件/km²) は L54 内に比べ非常に低く、これは山間地の単位面積あたり工事密度が" "都市部より少ないことを反映。

" + "

表: 市町別ファイルメタ (15 ファイル確認用)

" + df_to_html(T_city_meta) + "

この表から読み取れること: 県知事ファイルの polygon と" "市町別ファイルの polygon が対応していることを確認。" "庄原市の市町別ファイル 1057.94 km² は県知事ファイル中の対応 polygon と完全一致 = " "市町別ファイルは県知事ファイルの市町別切り出し版であり独立データではない。" "本研究では二重カウント回避のため県知事ファイル単独を採用した。

" ), ("仮説検証総合", sec7), ("発展課題", sec8), ] html = render_lesson( num=55, title="特定盛土等規制区域 単独 3 研究例分析 — 盛土規制法 30 条 1 項の山間地射程と L54 制度補完", tags=["L55", "盛土規制法", "特定盛土等規制区域", "RQ×3", "Format B", "30条1項", "制度間補完", "山間地監督", "L54連携", "L52連携", "L53連携", "30条系運用"], time="40 分", level="中級+", data_label="DoBoX dataset 1428 (Shapefile 16 ファイル + XLSX 規制資料) + L52/L53/L54/L10/L11 連携", sections=sections, script_filename="L55_specific_filling_zone.py", ) OUT_HTML = LESSONS / "L55_specific_filling_zone.html" OUT_HTML.write_text(html, encoding="utf-8") print(f" HTML written: {OUT_HTML}") print(f" ({time.time()-t11:.1f}s)", flush=True) print(f"\n=== L55 DONE in {time.time()-t_all:.1f}s ===")