日本の標準地域メッシュコードと PostGIS で都市環境を定量化する

CREATE TABLE mesh_metrics (
    mesh_code VARCHAR(10) PRIMARY KEY,  -- 250mメッシュ
    geom GEOMETRY(Polygon, 4326),       -- メッシュ境界ポリゴン
    population INTEGER,                  -- 人口（国勢調査）
    elevation_avg REAL,                  -- 平均標高
    hazard_flags JSONB,                  -- 災害リスクフラグ
    env_scores JSONB,                    -- 環境スコア
    updated_at TIMESTAMPTZ DEFAULT now()
);
 
-- 空間インデックス
CREATE INDEX idx_mesh_geom ON mesh_metrics USING GIST(geom);
-- コード前方一致検索用
CREATE INDEX idx_mesh_code ON mesh_metrics (mesh_code text_pattern_ops);

CREATE OR REPLACE FUNCTION lat_lng_to_mesh_250m(
    lat DOUBLE PRECISION,
    lng DOUBLE PRECISION
) RETURNS VARCHAR(10) AS $$
DECLARE
    lat_idx INTEGER; lng_idx INTEGER;
    m1_lat INTEGER; m1_lng INTEGER;
    m2_lat INTEGER; m2_lng INTEGER;
    m3_lat INTEGER; m3_lng INTEGER;
    m4 INTEGER; m5 INTEGER;
    rem_lat DOUBLE PRECISION; rem_lng DOUBLE PRECISION;
BEGIN
    -- 1次メッシュ
    m1_lat := floor(lat * 1.5)::INTEGER;
    m1_lng := (floor(lng) - 100)::INTEGER;
    rem_lat := lat * 1.5 - m1_lat;
    rem_lng := lng - floor(lng);
 
    -- 2次メッシュ (8分割)
    m2_lat := floor(rem_lat * 8)::INTEGER;
    m2_lng := floor(rem_lng * 8)::INTEGER;
    rem_lat := rem_lat * 8 - m2_lat;
    rem_lng := rem_lng * 8 - m2_lng;
 
    -- 3次メッシュ (10分割)
    m3_lat := floor(rem_lat * 10)::INTEGER;
    m3_lng := floor(rem_lng * 10)::INTEGER;
    rem_lat := rem_lat * 10 - m3_lat;
    rem_lng := rem_lng * 10 - m3_lng;
 
    -- 4次 (2分割)
    m4 := floor(rem_lat * 2)::INTEGER * 2 + floor(rem_lng * 2)::INTEGER + 1;
    rem_lat := rem_lat * 2 - floor(rem_lat * 2);
    rem_lng := rem_lng * 2 - floor(rem_lng * 2);
 
    -- 5次 (2分割)
    m5 := floor(rem_lat * 2)::INTEGER * 2 + floor(rem_lng * 2)::INTEGER + 1;
 
    RETURN lpad(m1_lat::TEXT, 2, '0')
        || lpad(m1_lng::TEXT, 2, '0')
        || m2_lat::TEXT || m2_lng::TEXT
        || m3_lat::TEXT || m3_lng::TEXT
        || m4::TEXT || m5::TEXT;
END;
$$ LANGUAGE plpgsql IMMUTABLE STRICT;

SELECT env_scores, hazard_flags
FROM mesh_metrics
WHERE mesh_code = lat_lng_to_mesh_250m(35.6812, 139.7671);

-- 1km メッシュ単位（3次メッシュ）で人口を集計
SELECT
    left(mesh_code, 8) AS mesh_1km,
    SUM(population) AS total_pop,
    AVG((env_scores->>'convenience')::REAL) AS avg_convenience
FROM mesh_metrics
WHERE mesh_code LIKE '5339%'  -- 東京周辺の1次メッシュ
GROUP BY left(mesh_code, 8);

-- 各メッシュの指標を百分位ランクに変換
WITH ranked AS (
    SELECT
        mesh_code,
        (env_scores->>'convenience')::REAL AS convenience,
        PERCENT_RANK() OVER (
            ORDER BY (env_scores->>'convenience')::REAL
        ) AS pct_rank
    FROM mesh_metrics
    WHERE mesh_code LIKE '5339%'
)
SELECT mesh_code, convenience, round(pct_rank * 100) AS percentile
FROM ranked
WHERE mesh_code = '5339451234';

SELECT
    mesh_code,
    (
        (env_scores->>'convenience')::REAL * 0.3 +
        (env_scores->>'safety')::REAL * 0.25 +
        (env_scores->>'childcare')::REAL * 0.2 +
        (env_scores->>'transit')::REAL * 0.25
    ) AS livability_score
FROM mesh_metrics
ORDER BY livability_score DESC
LIMIT 100;

-- 行政区画 × メッシュの空間結合で区別集計
SELECT
    a.city_name,
    AVG((m.env_scores->>'convenience')::REAL) AS avg_convenience,
    AVG((m.env_scores->>'safety')::REAL) AS avg_safety,
    COUNT(*) AS mesh_count
FROM mesh_metrics m
JOIN admin_boundary a ON ST_Contains(a.geom, m.geom)
WHERE a.prefecture = '東京都'
GROUP BY a.city_name
ORDER BY avg_convenience DESC;