ボランティア・社会参加活動と
地域コミュニティ指標の分析

pop = df['総人口']

# (1) 保育所数／千人 … 子育てコミュニティの代理変数
df['保育所数_千人'] = df['保育所等数'] / pop * 1000

# (2) 病院・診療所数／千人 … 医療・地域福祉の代理変数
df['医療施設数_千人'] = (df['一般病院数'] + df['一般診療所数']) / pop * 1000

# (3) 幼稚園数／千人 … 就学前教育コミュニティ
df['幼稚園数_千人'] = df['幼稚園数'] / pop * 1000

# (4) 教養娯楽費（社会参加コストの代理）
df['教養娯楽費'] = df['教養娯楽費（二人以上の世帯）']

# (5) ごみのリサイクル率 … 環境意識・市民参加の代理変数
df['リサイクル率'] = df['ごみのリサイクル率']

# (6) 転出率 … 低いほど地域定着（地域凝集性の逆指標）
df['転出率_千人'] = df['転出者数（日本人移動者）'] / pop * 1000

# (7) 高齢化率 … コミュニティ構造の指標
df['高齢化率'] = df['65歳以上人口'] / pop * 100

# (8) 合計特殊出生率 … 地域活力・若年層定着の指標
df['出生率'] = df['合計特殊出生率']

# 目的変数: 教養娯楽費（社会参加・余暇活動への支出 = 社会参加指標の代理）
Y_NAME = '教養娯楽費'

# 説明変数
X_NAMES = ['保育所数_千人', '医療施設数_千人', '幼稚園数_千人',
           'リサイクル率', '転出率_千人', '高齢化率', '出生率']

# 地域区分マップ
region_map = {
    '北海道': '北海道・東北', '青森県': '北海道・東北', '岩手県': '北海道・東北',
    '宮城県': '北海道・東北', '秋田県': '北海道・東北', '山形県': '北海道・東北',
    '福島県': '北海道・東北', '茨城県': '関東', '栃木県': '関東', '群馬県': '関東',
    '埼玉県': '関東', '千葉県': '関東', '東京都': '関東', '神奈川県': '関東',
    '新潟県': '中部', '富山県': '中部', '石川県': '中部', '福井県': '中部',
    '山梨県': '中部', '長野県': '中部', '岐阜県': '中部', '静岡県': '中部', '愛知県': '中部',
    '三重県': '近畿', '滋賀県': '近畿', '京都府': '近畿', '大阪府': '近畿',
    '兵庫県': '近畿', '奈良県': '近畿', '和歌山県': '近畿',
    '鳥取県': '中国・四国', '島根県': '中国・四国', '岡山県': '中国・四国',
    '広島県': '中国・四国', '山口県': '中国・四国', '徳島県': '中国・四国',
    '香川県': '中国・四国', '愛媛県': '中国・四国', '高知県': '中国・四国',
    '福岡県': '九州・沖縄', '佐賀県': '九州・沖縄', '長崎県': '九州・沖縄',
    '熊本県': '九州・沖縄', '大分県': '九州・沖縄', '宮崎県': '九州・沖縄',
    '鹿児島県': '九州・沖縄', '沖縄県': '九州・沖縄'
}
region_colors = {
    '北海道・東北': '#4e9af1', '関東': '#e05c5c', '中部': '#f0a500',
    '近畿': '#5cb85c', '中国・四国': '#9b59b6', '九州・沖縄': '#f39c12'
}
df['地域'] = df['都道府県'].map(region_map)
df['地域色'] = df['地域'].map(region_colors)

# 分析用サブデータ（欠損除去）
cols_use = [Y_NAME] + X_NAMES
dfa = df[['都道府県', '地域', '地域色'] + cols_use].dropna().copy()
print(f"分析対象: {len(dfa)}都道府県")

fig1, ax1 = plt.subplots(figsize=(9, 7))
n = len(corr_mat)
im = ax1.imshow(corr_mat.values, cmap='RdBu_r', vmin=-1, vmax=1, aspect='auto')
plt.colorbar(im, ax=ax1, fraction=0.04, pad=0.04)
ax1.set_xticks(range(n))
ax1.set_yticks(range(n))
ax1.set_xticklabels(corr_mat.columns, fontsize=9, rotation=45, ha='right')
ax1.set_yticklabels(corr_mat.columns, fontsize=9)
for i in range(n):
    for j in range(n):
        val = corr_mat.values[i, j]
        color = 'white' if abs(val) > 0.55 else 'black'
        ax1.text(j, i, f'{val:.2f}', ha='center', va='center',
                 fontsize=8, color=color, fontweight='bold')
ax1.set_title('社会参加関連変数のPearson相関行列（2022年・47都道府県）', fontsize=13, pad=15)
fig1.tight_layout()
fig1.savefig(os.path.join(FIG_DIR, '2020_U5_1_fig1.png'), dpi=150, bbox_inches='tight')
plt.close(fig1)
print("fig1 保存完了")

fig2, ax2 = plt.subplots(figsize=(9, 6))
for region, grp in dfa.groupby('地域'):
    ax2.scatter(grp['リサイクル率'], grp[Y_NAME],
                color=region_colors[region], s=60, label=region, zorder=3, alpha=0.85)
    for _, row in grp.iterrows():
        ax2.annotate(row['都道府県'], (row['リサイクル率'], row[Y_NAME]),
                     fontsize=6.5, ha='center', va='bottom',
                     xytext=(0, 4), textcoords='offset points', color='#333')

# 回帰線
x_plot = dfa['リサイクル率'].values
y_plot = dfa[Y_NAME].values
slope, intercept, r_val, p_val, _ = stats.linregress(x_plot, y_plot)
x_line = np.linspace(x_plot.min(), x_plot.max(), 100)
ax2.plot(x_line, slope * x_line + intercept, color='#333', linewidth=1.5,
         linestyle='--', zorder=2)
ax2.set_xlabel('ごみのリサイクル率（%）', fontsize=11)
ax2.set_ylabel('教養娯楽費（円，二人以上世帯）', fontsize=11)
ax2.set_title(f'リサイクル率と社会参加指標の関係\n（r = {r_val:.3f}, p = {p_val:.4f}, 47都道府県, 2022年）',
              fontsize=12)
ax2.legend(loc='upper left', fontsize=9, framealpha=0.85)
ax2.grid(True, alpha=0.3, linestyle='--')
fig2.tight_layout()
fig2.savefig(os.path.join(FIG_DIR, '2020_U5_1_fig2.png'), dpi=150, bbox_inches='tight')
plt.close(fig2)
print("fig2 保存完了")
print(f"  リサイクル率 vs 教養娯楽費: r={r_val:.3f}, p={p_val:.4f}")

x_labels = [
    '保育所数／千人', '医療施設数／千人', '幼稚園数／千人',
    'リサイクル率', '転出率／千人', '高齢化率', '合計特殊出生率'
]
sorted_idx = np.argsort(coefs)
sorted_coefs = coefs[sorted_idx]
sorted_labels = [x_labels[i] for i in sorted_idx]
sorted_pvals  = pvals[sorted_idx]
sorted_ci_lo  = conf.values[sorted_idx, 0]
sorted_ci_hi  = conf.values[sorted_idx, 1]

colors_bar = ['#e05c5c' if p < 0.05 else '#aab4c8' for p in sorted_pvals]

fig3, ax3 = plt.subplots(figsize=(8, 5))
bars = ax3.barh(range(len(sorted_coefs)), sorted_coefs,
                color=colors_bar, height=0.55, zorder=3)
for i, (lo, hi) in enumerate(zip(sorted_ci_lo, sorted_ci_hi)):
    ax3.plot([lo, hi], [i, i], color='#555', linewidth=1.5, zorder=4)
    ax3.plot([lo, lo], [i - 0.15, i + 0.15], color='#555', linewidth=1.5, zorder=4)
    ax3.plot([hi, hi], [i - 0.15, i + 0.15], color='#555', linewidth=1.5, zorder=4)
ax3.axvline(0, color='#333', linewidth=1, linestyle='-')
ax3.set_yticks(range(len(sorted_labels)))
ax3.set_yticklabels(sorted_labels, fontsize=10)
ax3.set_xlabel('標準化偏回帰係数 (β)', fontsize=11)
ax3.set_title(f'重回帰分析：標準化偏回帰係数\n（目的変数：教養娯楽費, R²={ols_model.rsquared:.3f}）', fontsize=12)
ax3.grid(True, axis='x', alpha=0.3, linestyle='--')
# 有意・非有意の凡例
from matplotlib.patches import Patch
legend_elements = [Patch(facecolor='#e05c5c', label='p < 0.05（有意）'),
                   Patch(facecolor='#aab4c8', label='p ≥ 0.05（非有意）')]
ax3.legend(handles=legend_elements, loc='lower right', fontsize=9)
fig3.tight_layout()
fig3.savefig(os.path.join(FIG_DIR, '2020_U5_1_fig3.png'), dpi=150, bbox_inches='tight')
plt.close(fig3)
print("fig3 保存完了")