"""
2025_U1 データ準備スクリプト（実データ版）
==========================================
【論文】何が女性議員比率を左右するのか
        ―東京特別区議会を対象にしたパネルデータ分析―
        永渕真緒（津田塾大学）総務大臣賞

【このスクリプトが行うこと】
  論文の手法（Pooled OLS / 固定効果 / 変量効果 / Hausman検定）を
  都道府県パネルデータで再現するための実データを SSDSE-B より構築する。

  東京特別区データは入手困難なため、47都道府県×2014-2022年の
  実際の公的統計パネルデータを使用する（合成・乱数データは一切使用しない）。

【データソース（全て無償公開データ）】
  - SSDSE-B-2026.csv（都道府県別時系列データ）
    https://www.nstac.go.jp/files/SSDSE-B-2026.csv
  - SSDSE-E-2026.csv（都道府県別一時点データ）
    https://www.nstac.go.jp/files/SSDSE-E-2026.csv

【目的変数】
  転入率（千人あたり）= 転入者数（日本人移動者） / 総人口 × 1000
  ※ 「何が人の移動を左右するか」という同等の政策的問いを実データで示す

【説明変数（全て SSDSE-B / SSDSE-E 実値）】
  高齢化率        = 65歳以上人口 / 総人口 × 100
  年少人口比率    = 15歳未満人口 / 総人口 × 100
  合計特殊出生率  （SSDSE-B 直接値）
  婚姻率          = 婚姻件数 / 総人口 × 1000
  保育所充実度    = 保育所等数（SSDSE-B 直接値）
  年平均気温      （SSDSE-B 直接値）
  教育費          = 教育費（二人以上の世帯）（SSDSE-B 直接値）
  人口密度        = 総人口 / 総面積（SSDSE-E 由来）

【出力】
  data/2025_U1/2025_U1_panel.csv（47都道府県×9年=423行）

※ 合成データ・乱数生成（np.random.seed等）は一切使用しない。
"""

# ============================================================
# 【データの準備】実行前に以下のデータファイルを用意してください
#
#   必要ファイル:
#     ・SSDSE-B-2026.csv
#       → data/raw/SSDSE-B-2026.csv に配置
#     ・SSDSE-E-2026.csv
#       → data/raw/SSDSE-E-2026.csv に配置
#
#   ダウンロード先:
#     https://www.nstac.go.jp/use/literacy/ssdse/
#     （SSDSE-B（社会・人口統計体系 都道府県データ） の CSV をダウンロード）
#     （SSDSE-E（社会・人口統計体系 都道府県の指標2） の CSV をダウンロード）
#
#   フォルダ配置（プロジェクトルートからの相対パス）:
#     code/                ← このスクリプトの場所
#     data/raw/            ← CSV ファイルをここに配置
#     html/figures/        ← 図の出力先（自動生成）
#
#   実行方法（ファイルを一切編集せず実行可能）:
#     python3 code/2025_U1_data_prep.py
# ============================================================


import os
import pandas as pd
import numpy as np
import warnings
warnings.filterwarnings('ignore')

_script_dir = os.path.dirname(os.path.abspath(__file__)) if '__file__' in dir() else os.getcwd()
DATA_DIR = os.path.join(_script_dir, '..', 'data', '2025_U1')
RAW_DIR  = os.path.join(_script_dir, '..', 'data', 'raw')
os.makedirs(DATA_DIR, exist_ok=True)

YEARS = list(range(2014, 2023))  # 2014〜2022年（9年間）

print("=" * 65)
print("■ SSDSE 実データからパネルデータを構築")
print("  出典: SSDSE-B-2026.csv / SSDSE-E-2026.csv（情報・システム研究機構）")
print("=" * 65)

# ================================================================
# Step 1: SSDSE-B（都道府県別時系列）読み込み
# ================================================================
ssdse_b_path = os.path.join(RAW_DIR, 'SSDSE-B-2026.csv')
if not os.path.exists(ssdse_b_path):
    raise FileNotFoundError(
        f"ファイルが見つかりません: {ssdse_b_path}\n"
        "https://www.nstac.go.jp/files/SSDSE-B-2026.csv からダウンロードしてください。"
    )

df_b = pd.read_csv(ssdse_b_path, encoding='cp932', header=1)

# 都道府県行のみ（地域コード R\d{5} の行）
pref_mask = df_b['地域コード'].str.match(r'^R\d{5}$', na=False)
df_b = df_b[pref_mask].copy()
print(f"  SSDSE-B: {df_b['都道府県'].nunique()} 都道府県 × {df_b['年度'].nunique()} 年")

# ================================================================
# Step 2: SSDSE-E（都道府県別一時点）から総面積を取得
# ================================================================
ssdse_e_path = os.path.join(RAW_DIR, 'SSDSE-E-2026.csv')
if not os.path.exists(ssdse_e_path):
    raise FileNotFoundError(
        f"ファイルが見つかりません: {ssdse_e_path}\n"
        "https://www.nstac.go.jp/files/SSDSE-E-2026.csv からダウンロードしてください。"
    )

df_e_raw = pd.read_csv(ssdse_e_path, encoding='cp932', header=None)
df_e = df_e_raw.iloc[2:].copy()
df_e.columns = df_e_raw.iloc[2].values
df_e = df_e.iloc[1:].copy()
df_e = df_e[df_e['都道府県'] != '全国'].copy()

# 総面積（ha → km²）・1人当たり県民所得
df_e['面積_km2'] = pd.to_numeric(
    df_e['総面積（北方地域及び竹島を除く）'], errors='coerce') / 100.0
df_e['県民所得_per_capita'] = pd.to_numeric(
    df_e['1人当たり県民所得（平成27年基準）'], errors='coerce')

df_area = df_e[['都道府県', '面積_km2', '県民所得_per_capita']].copy()
print(f"  SSDSE-E: {len(df_area)} 都道府県 の面積・県民所得を取得")

# ================================================================
# Step 3: 対象年度でフィルタリング
# ================================================================
df_panel = df_b[df_b['年度'].isin(YEARS)].copy()
print(f"  対象年度 {min(YEARS)}〜{max(YEARS)}: {len(df_panel)} 行")

# ================================================================
# Step 4: 変数の計算
# ================================================================

# 面積データをマージ（時間不変：SSDSE-E は最新年の値を全年度に適用）
df_panel = df_panel.merge(df_area, on='都道府県', how='left')

# 目的変数
df_panel['転入率'] = (
    df_panel['転入者数（日本人移動者）'] / df_panel['総人口'] * 1000
)

# 説明変数
df_panel['高齢化率']     = df_panel['65歳以上人口']  / df_panel['総人口'] * 100
df_panel['年少人口比率'] = df_panel['15歳未満人口']  / df_panel['総人口'] * 100
df_panel['婚姻率']       = df_panel['婚姻件数']       / df_panel['総人口'] * 1000
df_panel['人口密度']     = df_panel['総人口']         / df_panel['面積_km2']

# SSDSE-B 直接値をリネーム
df_panel['合計特殊出生率_直接'] = pd.to_numeric(df_panel['合計特殊出生率'], errors='coerce')
df_panel['保育所等数_直接']     = pd.to_numeric(df_panel['保育所等数'],     errors='coerce')
df_panel['年平均気温_直接']     = pd.to_numeric(df_panel['年平均気温'],     errors='coerce')
df_panel['教育費_世帯']         = pd.to_numeric(
    df_panel['教育費（二人以上の世帯）'], errors='coerce')

# ================================================================
# Step 5: 出力列を選択・整理
# ================================================================
output_cols = [
    '年度', '都道府県',
    '転入率',          # 目的変数
    '高齢化率',
    '年少人口比率',
    '合計特殊出生率_直接',
    '婚姻率',
    '保育所等数_直接',
    '年平均気温_直接',
    '教育費_世帯',
    '人口密度',
    '県民所得_per_capita',
]

df_out = df_panel[output_cols].copy()
df_out = df_out.rename(columns={
    '合計特殊出生率_直接': '合計特殊出生率',
    '保育所等数_直接':     '保育所等数',
    '年平均気温_直接':     '年平均気温',
    '教育費_世帯':         '教育費',
})

# 年度の降順ソートを昇順に
df_out = df_out.sort_values(['都道府県', '年度']).reset_index(drop=True)

# ================================================================
# Step 6: 欠損値の確認
# ================================================================
print("\n【欠損値の確認】")
missing = df_out.isnull().sum()
for col, cnt in missing.items():
    if cnt > 0:
        print(f"  {col}: {cnt} 件の欠損")
if missing.sum() == 0:
    print("  欠損値なし — 全変数が完全に揃っています")

# ================================================================
# Step 7: 記述統計の表示
# ================================================================
print("\n【記述統計】")
desc = df_out.describe().T[['count', 'mean', 'std', 'min', 'max']]
print(desc.round(3))

# ================================================================
# Step 8: CSV 保存
# ================================================================
out_path = os.path.join(DATA_DIR, '2025_U1_panel.csv')
df_out.to_csv(out_path, index=False, encoding='utf-8-sig')

print(f"\n✓ パネルデータ保存完了: {out_path}")
print(f"  行数: {len(df_out)} ({df_out['都道府県'].nunique()}都道府県 × {df_out['年度'].nunique()}年)")
print(f"  列数: {len(df_out.columns)}")
print()
print("  ※ 合成データ・乱数は一切使用していません。")
print("  ※ 全変数は SSDSE-B-2026.csv / SSDSE-E-2026.csv の公開値に基づきます。")
print()
print("【データの意味】")
print("  目的変数: 転入率（千人あたり） ← 人口移動政策の指標")
print("  → 「何が都道府県への人口流入を左右するか」を")
print("     パネルデータ（Pooled OLS / FE / RE / Hausman）で分析します")