旅館及びホテルにおける日本人・外国人宿泊客の
都道府県別増減から考える旅館の復活
―岡山県湯原温泉の視点からインバウンド需要を旅館に取り込む方策―

import os
import numpy as np
import pandas as pd
import matplotlib
matplotlib.use('Agg')
import matplotlib.pyplot as plt
from scipy import stats

plt.rcParams['font.family'] = 'Hiragino Sans'
plt.rcParams['axes.unicode_minus'] = False
plt.rcParams['figure.dpi'] = 150

FIG_DIR = 'html/figures'
DATA_B  = 'data/raw/SSDSE-B-2026.csv'
os.makedirs(FIG_DIR, exist_ok=True)

df = pd.read_csv(DATA_B, encoding='cp932', header=1)
df = df[df['地域コード'].str.match(r'^R\d{5}', na=False)].copy()
df['年度'] = df['年度'].astype(int)

cols = ['旅館営業施設数（ホテルを含む）', '旅館営業施設客室数（ホテルを含む）',
        '延べ宿泊者数', '外国人延べ宿泊者数']
print("収録年度:", df['年度'].min(), "〜", df['年度'].max())
print("2019年の都道府県数:", (df['年度'] == 2019).sum())
print("使用する列:", cols)

piv_shisetsu = df.pivot(index='都道府県', columns='年度',
                        values='旅館営業施設数（ホテルを含む）')
piv_kyakushitsu = df.pivot(index='都道府県', columns='年度',
                           values='旅館営業施設客室数（ホテルを含む）')

chg_shisetsu    = (piv_shisetsu[2023]    - piv_shisetsu[2012])    / piv_shisetsu[2012]    * 100
chg_kyakushitsu = (piv_kyakushitsu[2023] - piv_kyakushitsu[2012]) / piv_kyakushitsu[2012] * 100

print("=== 2012→2023年 増減率（ホテルを含む合算・%） ===")
print(f"施設数が減少した都道府県: {(chg_shisetsu < 0).sum()} / 47")
print(f"客室数が減少した都道府県: {(chg_kyakushitsu < 0).sum()} / 47")
print(f"施設数増減率  中央値: {chg_shisetsu.median():+.1f}%")
print(f"客室数増減率  中央値: {chg_kyakushitsu.median():+.1f}%")
print(f"客室数増減率  最大: {chg_kyakushitsu.idxmax()} {chg_kyakushitsu.max():+.1f}% / "
      f"最小: {chg_kyakushitsu.idxmin()} {chg_kyakushitsu.min():+.1f}%")
print(f"岡山県: 施設数 {chg_shisetsu['岡山県']:+.1f}% / 客室数 {chg_kyakushitsu['岡山県']:+.1f}%")

fig, ax = plt.subplots(figsize=(8, 6))
bp = ax.boxplot([chg_shisetsu.values, chg_kyakushitsu.values],
                tick_labels=['施設数の増減率', '客室数の増減率'],
                patch_artist=True, widths=0.45,
                medianprops=dict(color='#C62828', linewidth=2))
for patch, c in zip(bp['boxes'], ['#90CAF9', '#1565C0']):
    patch.set_facecolor(c)
    patch.set_alpha(0.85)
ax.axhline(0, color='#555', linewidth=1, linestyle='--')
ax.set_ylabel('2012→2023年 増減率（%）')
ax.set_title('図1（参考再現）: 宿泊施設数・客室数の増減率のばらつき（47都道府県）\n'
             'SSDSE-B-2026「旅館営業施設数／客室数（ホテルを含む）」\n'
             '※原論文（図1）は SSDSE-2019B で旅館とホテルを分離・2007→2016年',
             fontsize=11)
ax.grid(axis='y', alpha=0.3)
plt.tight_layout()
plt.savefig(f'{FIG_DIR}/2019_H5_1_fig1.png', bbox_inches='tight')
plt.close()
print('図1 保存: 2019_H5_1_fig1.png')

# 以下の相関係数はすべて【原論文の報告値】であり、本スクリプトで再計算したものではない。
# 原データ: SSDSE-2019B（衛生行政報告例）・宿泊旅行統計調査（2008年・2016年）
paper_t1 = {   # 原論文 表1: 旅館営業施設客室数増減（2007→2016）× 宿泊タイプ別稼働率
    '2008年 旅館稼働率':         (0.320, False),
    '2008年 リゾートホテル稼働率': (0.530, True),
    '2008年 ビジネスホテル稼働率': (0.477, True),
    '2008年 シティホテル稼働率':   (0.492, True),
    '2016年 旅館稼働率':         (0.225, False),
    '2016年 リゾートホテル稼働率': (0.662, True),
    '2016年 ビジネスホテル稼働率': (0.438, True),
    '2016年 シティホテル稼働率':   (0.492, True),
}
paper_t2 = {   # 原論文 表2: 外国人延べ宿泊者数増減（2008→2016）× 宿泊タイプ別の増減
    '旅館':          (0.150, False),
    'リゾートホテル': (0.325, False),
    'ビジネスホテル': (0.620, True),
    'シティホテル':   (0.545, True),
}

fig, axes = plt.subplots(1, 2, figsize=(12.5, 5.5))
for ax, d, ttl in zip(
        axes, [paper_t1, paper_t2],
        ['表1: 旅館営業施設客室数増減と\n宿泊タイプ別稼働率の相関',
         '表2: 外国人延べ宿泊者数増減と\n宿泊タイプ別増減の相関']):
    labels = list(d.keys())[::-1]
    vals   = [d[k][0] for k in labels]
    sig    = [d[k][1] for k in labels]
    colors = ['#1565C0' if s else '#B0BEC5' for s in sig]
    bars = ax.barh(labels, vals, color=colors)
    for b, v, s in zip(bars, vals, sig):
        ax.text(v + 0.012, b.get_y() + b.get_height()/2,
                f'{v:.3f}' + ('**' if s else ''), va='center', fontsize=10)
    ax.set_xlim(0, 0.8)
    ax.set_xlabel('相関係数（原論文の報告値）')
    ax.set_title(ttl, fontsize=11)
    ax.grid(axis='x', alpha=0.3)
fig.suptitle('図2: 原論文（表1・表2）が報告した相関係数の可視化 '
             '（**＝有意確率1%未満、灰色＝有意でない）\n'
             '※本図は原論文の報告値をそのまま図示したものであり、再計算ではない',
             fontsize=12)
plt.tight_layout(rect=[0, 0, 1, 0.90])
plt.savefig(f'{FIG_DIR}/2019_H5_1_fig2.png', bbox_inches='tight')
plt.close()
print('図2 保存: 2019_H5_1_fig2.png（原論文の報告値の可視化）')

nat = df.groupby('年度')[['延べ宿泊者数', '外国人延べ宿泊者数']].sum()
nat['日本人延べ宿泊者数'] = nat['延べ宿泊者数'] - nat['外国人延べ宿泊者数']
nat_man = nat / 10000  # 万人泊

print("=== 全国の延べ宿泊者数（万人泊） ===")
for y in [2012, 2019, 2020, 2023]:
    print(f"{y}年  日本人 {nat_man.loc[y, '日本人延べ宿泊者数']:,.0f}  "
          f"外国人 {nat_man.loc[y, '外国人延べ宿泊者数']:,.0f}")
r1219 = nat.loc[2019, '外国人延べ宿泊者数'] / nat.loc[2012, '外国人延べ宿泊者数']
print(f"外国人延べ宿泊者数 2012→2019年: {r1219:.2f}倍（+{(r1219-1)*100:.0f}%）")

fig, ax = plt.subplots(figsize=(9.5, 5.5))
ax.plot(nat_man.index, nat_man['日本人延べ宿泊者数'], marker='o',
        color='#1565C0', label='日本人（延べ−外国人）')
ax.plot(nat_man.index, nat_man['外国人延べ宿泊者数'], marker='s',
        color='#E65100', label='外国人')
ax.axvspan(2020, 2021, color='#FFCDD2', alpha=0.4)
ax.text(2020.5, nat_man['日本人延べ宿泊者数'].max()*0.55, 'コロナ禍',
        ha='center', color='#C62828', fontsize=11)
ax.set_xlabel('年')
ax.set_ylabel('延べ宿泊者数（万人泊）')
ax.set_title('図3（実データ再現）: 日本人・外国人の延べ宿泊者数の推移（全国計）\n'
             'SSDSE-B-2026・2012〜2023年 ※原論文（図5）は宿泊旅行統計調査 2007〜2018年',
             fontsize=11)
ax.legend()
ax.grid(alpha=0.3)
plt.tight_layout()
plt.savefig(f'{FIG_DIR}/2019_H5_1_fig3.png', bbox_inches='tight')
plt.close()
print('図3 保存: 2019_H5_1_fig3.png')

piv_gaikoku = df.pivot(index='都道府県', columns='年度', values='外国人延べ宿泊者数')
chg_gaikoku = ((piv_gaikoku[2019] - piv_gaikoku[2012]) / piv_gaikoku[2012] * 100).sort_values()

print("=== 外国人延べ宿泊者数 増減率 2012→2019年（コロナ前・%） ===")
print("上位5県:")
for p, v in chg_gaikoku.tail(5)[::-1].items():
    print(f"  {p}: +{v:,.0f}%")
print("下位3県:")
for p, v in chg_gaikoku.head(3).items():
    print(f"  {p}: {v:+,.0f}%")
print(f"岡山県: +{chg_gaikoku['岡山県']:,.0f}%（全国 {int((chg_gaikoku.rank(ascending=False))['岡山県'])} 位）")

# 原論文 表2 の趣旨の実データ版: 外国人宿泊の増減と客室数の増減の相関（同期間）
# 「増減数（絶対数）」と「増減率」の両方で計算し、結果の違いを確認する
diff_gaikoku = piv_gaikoku[2019] - piv_gaikoku[2012]
diff_kyaku   = piv_kyakushitsu[2019] - piv_kyakushitsu[2012]
r_n, p_n = stats.pearsonr(diff_gaikoku, diff_kyaku.reindex(diff_gaikoku.index))
chg_kyaku1219 = (piv_kyakushitsu[2019] - piv_kyakushitsu[2012]) / piv_kyakushitsu[2012] * 100
r_r, p_r = stats.pearsonr(chg_gaikoku.reindex(chg_kyaku1219.index), chg_kyaku1219)
print("外国人宿泊の増減 × 客室数の増減（2012→2019, ホテルを含む合算）:")
print(f"  増減数どうし: r = {r_n:.3f}, p = {p_n:.5f}  ← 県の規模に引っ張られて強い正")
print(f"  増減率どうし: r = {r_r:.3f}, p = {p_r:.4f}  ← 規模の影響を除くとほぼ無相関")

highlight = {'香川県': '#C62828', '鳥取県': '#C62828',
             '岡山県': '#E65100', '福島県': '#6A1B9A'}
colors = [highlight.get(p, '#90A4AE') for p in chg_gaikoku.index]
fig, ax = plt.subplots(figsize=(9, 12))
ax.barh(chg_gaikoku.index, chg_gaikoku.values, color=colors)
for p in highlight:
    v = chg_gaikoku[p]
    ax.text(v + 15, list(chg_gaikoku.index).index(p), f'{p} +{v:,.0f}%',
            va='center', fontsize=10, color=highlight[p], fontweight='bold')
ax.set_xlabel('外国人延べ宿泊者数 増減率 2012→2019年（%）')
ax.set_title('図4（実データ再現）: 都道府県別 外国人延べ宿泊者数の増減率ランキング\n'
             'SSDSE-B-2026・2012→2019年（コロナ前） '
             '※原論文（図6）は宿泊旅行統計調査 2008→2016年の階級区分図',
             fontsize=11)
ax.grid(axis='x', alpha=0.3)
ax.tick_params(axis='y', labelsize=9)
plt.tight_layout()
plt.savefig(f'{FIG_DIR}/2019_H5_1_fig4.png', bbox_inches='tight')
plt.close()
print('図4 保存: 2019_H5_1_fig4.png')

# 湯原温泉の宿泊者数は「湯原町旅館協同組合調べ」の非公開データのため再現不可。
# 原論文本文が明記した2時点の値のみを示す（推移グラフは原論文の図7を参照）。
yubara_1989 = 220832   # 原論文の報告値（バブル期のピーク）
yubara_2018 = 128274   # 原論文の報告値
chg = (yubara_2018 - yubara_1989) / yubara_1989 * 100
print("=== 湯原温泉の宿泊者数（原論文の報告値・湯原町旅館協同組合調べ） ===")
print(f"1989年: {yubara_1989:,} 人")
print(f"2018年: {yubara_2018:,} 人")
print(f"増減率: {chg:.1f}%（約30年でほぼ半減。推移グラフは原論文の図7を参照）")
print("\n完了: 図は html/figures/ に保存されました")