{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 논문 1 · 챕터 5 — 시나리오 시뮬레이션\n",
    "\n",
    "[![Open In Colab](https://colab.research.google.com/assets/colab-badge.svg)](https://colab.research.google.com/github/sdkparkforbi/gva-scenario-2025-2035/blob/main/ch5/code/ch5_scenarios.ipynb)\n",
    "\n",
    "> **Google Colab에서 바로 실행**됩니다. 배지를 누른 뒤 `런타임 → 모두 실행`.\n",
    "> 데이터는 아래 셀에서 자동으로 내려받고, 그림은 셀 실행 시 바로 표시됩니다.\n",
    "\n",
    "각 셀에는 **셀 번호 · 셀 제목**이 달려 있고, 코드에는 줄별 **주석**을 넣었습니다."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 0. 준비 — Colab 환경과 데이터"
   ]
  },
  {
   "cell_type": "code",
   "metadata": {},
   "execution_count": null,
   "outputs": [],
   "source": [
    "# ===== 셀 1 · 한글 폰트 (koreanize-matplotlib) =====\n",
    "# 그림의 한글이 깨지지 않도록 koreanize-matplotlib 를 설치하고 불러옵니다.\n",
    "# import 하는 순간 matplotlib 한글 폰트가 자동 설정됩니다.\n",
    "_font_ok = False\n",
    "try:\n",
    "    import koreanize_matplotlib                       # 이미 설치돼 있으면 바로 사용\n",
    "    _font_ok = True\n",
    "except Exception:\n",
    "    try:\n",
    "        import subprocess, sys\n",
    "        subprocess.run([sys.executable,'-m','pip','install','-q','koreanize-matplotlib'], check=False)\n",
    "        import koreanize_matplotlib                   # 설치 후 다시 불러오기\n",
    "        _font_ok = True\n",
    "    except Exception:\n",
    "        _font_ok = False\n",
    "if not _font_ok:                                      # 폴백: 나눔폰트 직접 설치·등록\n",
    "    import subprocess\n",
    "    subprocess.run(['apt-get','-qq','-y','install','fonts-nanum'], check=False)\n",
    "    import matplotlib.font_manager as fm, matplotlib.pyplot as plt\n",
    "    _p = '/usr/share/fonts/truetype/nanum/NanumGothic.ttf'\n",
    "    try:\n",
    "        fm.fontManager.addfont(_p)\n",
    "        plt.rcParams['font.family'] = fm.FontProperties(fname=_p).get_name()  # 실제 등록 이름 사용\n",
    "    except Exception: pass\n",
    "import matplotlib.pyplot as plt\n",
    "plt.rcParams['axes.unicode_minus'] = False           # 마이너스 기호 깨짐 방지\n",
    "print('한글 폰트 준비 완료')"
   ]
  },
  {
   "cell_type": "code",
   "metadata": {},
   "execution_count": null,
   "outputs": [],
   "source": [
    "# ===== 셀 2 · 데이터 다운로드 =====\n",
    "# 이 챕터의 입력 데이터를 서버에서 로컬 data/ 폴더로 내려받습니다.\n",
    "import os, urllib.request, urllib.parse\n",
    "BASE = \"https://aiforalab.com/gva-scenario-2025-2035/ch5/data/\"          # 서버의 이 챕터 데이터 폴더\n",
    "FILES = [\"dlog_rates_sa.csv\", \"dsge_exog_paths.csv\", \"exog_quarterly_filled.csv\", \"levels_sa.csv\", \"scenario_industry_2035.csv\"]\n",
    "os.makedirs(\"data\", exist_ok=True)\n",
    "for f in FILES:                                      # 한글 파일명은 URL 인코딩 후 요청\n",
    "    urllib.request.urlretrieve(BASE + urllib.parse.quote(f), os.path.join('data', f))\n",
    "print('다운로드 완료:', len(FILES), '개 ->', sorted(os.listdir('data')))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 1. 방법 구현과 계산"
   ]
  },
  {
   "cell_type": "code",
   "metadata": {},
   "execution_count": null,
   "outputs": [],
   "source": [
    "# ===== 셀 3 · 라이브러리 임포트·설정 =====\n",
    "import os, json\n",
    "import numpy as np, pandas as pd\n",
    "import matplotlib\n",
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "cell_type": "code",
   "metadata": {},
   "execution_count": null,
   "outputs": [],
   "source": [
    "# ===== 셀 4 · Colab 경로: 데이터=data/, 그림=figures/ =====\n",
    "# ----- Colab 경로: 데이터=data/, 그림=figures/ -----\n",
    "import os\n",
    "DATA = \"data\"; PREP = \"data\"; FIG = \"figures\"; OUT = \"figures\"\n",
    "os.makedirs(\"data\", exist_ok=True); os.makedirs(\"figures\", exist_ok=True)"
   ]
  },
  {
   "cell_type": "code",
   "metadata": {},
   "execution_count": null,
   "outputs": [],
   "source": [
    "# ===== 셀 5 · 계산 =====\n",
    "plt.rcParams[\"savefig.dpi\"]=130\n",
    "P,Q,LAM,LAMX,WIN=4,2,0.1,1.0,30.0\n",
    "EXO=[\"oil_g\",\"trd_g\",\"d_rrate\",\"rfx_g\",\"lab_g\"]\n",
    "SCN=[\"S1\",\"S2\",\"S3\",\"S4\",\"S5\"]; SCN_KO={\"S1\":\"기준\",\"S2\":\"공급제약+충격\",\"S3\":\"제약완화\",\"S4\":\"기술도약\",\"S5\":\"지정학\"}\n",
    "CLR={\"S1\":\"#333\",\"S2\":\"#c0392b\",\"S3\":\"#1f4e79\",\"S4\":\"#2e7d32\",\"S5\":\"#b8892f\"}"
   ]
  },
  {
   "cell_type": "code",
   "metadata": {},
   "execution_count": null,
   "outputs": [],
   "source": [
    "# ===== 셀 6 · 데이터·공통표본 (챕터3과 동일) =====\n",
    "# 데이터·공통표본 (챕터3과 동일)\n",
    "rate=pd.read_csv(os.path.join(DATA,\"dlog_rates_sa.csv\"),encoding=\"utf-8-sig\",index_col=0); rate.index=pd.PeriodIndex(rate.index,freq=\"Q\")\n",
    "g=rate.drop(columns=[\"총량GDP\"]).dropna().clip(-WIN,WIN)\n",
    "ex=pd.read_csv(os.path.join(DATA,\"exog_quarterly_filled.csv\"),encoding=\"utf-8-sig\"); ex.index=pd.PeriodIndex(ex[\"분기\"],freq=\"Q\"); ex=ex[EXO].astype(float)\n",
    "common=g.index.intersection(ex.index).sort_values(); g=g.loc[common]; ex=ex.loc[common]\n",
    "names=list(g.columns); n=len(names); ne=len(EXO); gv=g.values; xv=ex.values"
   ]
  },
  {
   "cell_type": "code",
   "metadata": {},
   "execution_count": null,
   "outputs": [],
   "source": [
    "# ===== 셀 7 · 함수 정의: ar_sigma, design =====\n",
    "def ar_sigma(m,p):\n",
    "    s=np.zeros(m.shape[1])\n",
    "    for i in range(m.shape[1]):\n",
    "        y=m[:,i]; Y=y[p:]; Z=np.column_stack([y[p-l:len(y)-l] for l in range(1,p+1)]+[np.ones(len(Y))])\n",
    "        b,*_=np.linalg.lstsq(Z,Y,rcond=None); s[i]=(Y-Z@b).std(ddof=Z.shape[1])\n",
    "    return s\n",
    "def design(gvv,xvv,p,q):\n",
    "    T0=max(p,q); TT=len(gvv)-T0; blocks=[]\n",
    "    for l in range(1,p+1): blocks.append(gvv[[t-l for t in range(T0,len(gvv))]])\n",
    "    for s in range(0,q+1): blocks.append(xvv[[t-s for t in range(T0,len(gvv))]])\n",
    "    return gvv[T0:], np.column_stack(blocks+[np.ones(TT)])"
   ]
  },
  {
   "cell_type": "code",
   "metadata": {},
   "execution_count": null,
   "outputs": [],
   "source": [
    "# ===== 셀 8 · BVARX 추정 (전표본) =====\n",
    "# BVARX 추정 (전표본)\n",
    "sy=ar_sigma(gv,P); sx=xv.std(0)\n",
    "d=[]\n",
    "for l in range(1,P+1):\n",
    "    for j in range(n): d.append((l**2)*(sy[j]**2)/LAM**2)\n",
    "for s in range(0,Q+1):\n",
    "    for m in range(ne): d.append(((s+1)**2)*(sx[m]**2)/LAMX**2)\n",
    "d.append(1e-6); d=np.array(d)\n",
    "Y,X=design(gv,xv,P,Q); B=np.linalg.solve(X.T@X+np.diag(d),X.T@Y)\n",
    "\n",
    "# DSGE 시나리오 외생 경로\n",
    "dsge=pd.read_csv(os.path.join(DATA,\"dsge_exog_paths.csv\"),encoding=\"utf-8-sig\")\n",
    "Hf=dsge[dsge[\"시나리오\"]==\"S1\"].shape[0]      # 40분기"
   ]
  },
  {
   "cell_type": "code",
   "metadata": {},
   "execution_count": null,
   "outputs": [],
   "source": [
    "# ===== 셀 9 · 함수 정의: scen_forecast =====\n",
    "def scen_forecast(scn):\n",
    "    xs=dsge[dsge[\"시나리오\"]==scn][EXO].values          # (Hf, ne)\n",
    "    exog_full=np.vstack([xv, xs]); start=len(xv)\n",
    "    buf=[r.copy() for r in gv[-P:]]; out=[]\n",
    "    for h in range(Hf):\n",
    "        ti=start+h; row=[]\n",
    "        for l in range(1,P+1): row.append(buf[-l])\n",
    "        for s in range(0,Q+1): row.append(exog_full[ti-s])\n",
    "        y=np.concatenate(row+[[1.0]])@B; buf.append(y); out.append(y)\n",
    "    return np.array(out)                                   # (Hf, n) 성장률\n",
    "\n",
    "# 가중치(2025Q4 실질수준)\n",
    "lev=pd.read_csv(os.path.join(DATA,\"levels_sa.csv\"),encoding=\"utf-8-sig\",index_col=0)\n",
    "w=lev[names].iloc[-1].values; w=w/w.sum()"
   ]
  },
  {
   "cell_type": "code",
   "metadata": {},
   "execution_count": null,
   "outputs": [],
   "source": [
    "# ===== 셀 10 · 계산 =====\n",
    "gdp_idx={}; ind2035={}\n",
    "for s in SCN:\n",
    "    gr=scen_forecast(s)                                    # (Hf,n) %\n",
    "    agg=(gr*w[None,:]).sum(1)                              # 가중 총량 성장 %\n",
    "    gdp_idx[s]=100*np.exp(np.cumsum(agg/100))\n",
    "    ind2035[s]=100*np.exp(np.cumsum(gr/100,0))[-1]         # 산업 2035 지수(2025=100)\n",
    "ind_df=pd.DataFrame(ind2035,index=names)[SCN]\n",
    "ind_df.round(1).to_csv(os.path.join(DATA,\"scenario_industry_2035.csv\"),encoding=\"utf-8-sig\")"
   ]
  },
  {
   "cell_type": "code",
   "metadata": {},
   "execution_count": null,
   "outputs": [],
   "source": [
    "# ===== 셀 11 · 함수 정의: short =====\n",
    "def short(s): return (s.replace(\" 제조업\",\"\").replace(\" 및 \",\"·\").replace(\"서비스업\",\"서비스\"))[:9]"
   ]
  },
  {
   "cell_type": "code",
   "metadata": {},
   "execution_count": null,
   "outputs": [],
   "source": [
    "# ===== 셀 12 · 결과 출력 =====\n",
    "print(\"챕터5 — 시나리오 시뮬레이션 (BVARX×DSGE 조건화)\")\n",
    "print(f\"{'시나리오':<16}{'2035 GDP지수(2025=100)':>22}\")\n",
    "for s in SCN: print(f\"{s+' '+SCN_KO[s]:<16}{gdp_idx[s][-1]:>22.1f}\")\n",
    "labor=(ind_df[\"S3\"]-ind_df[\"S2\"]).sort_values()          # 제약완화-공급제약 = 노동효과\n",
    "print(\"\\n[노동효과(S3−S2) 큰 산업]\")\n",
    "for nm,v in labor.tail(4)[::-1].items(): print(f\"  +{v:.1f}  {short(nm)}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 2. 결과와 그림"
   ]
  },
  {
   "cell_type": "code",
   "metadata": {},
   "execution_count": null,
   "outputs": [],
   "source": [
    "# ===== 셀 13 · 그림 그리기 =====\n",
    "t=[pd.Period(q,\"Q\").to_timestamp(how=\"end\") for q in dsge[dsge['시나리오']=='S1']['분기']]\n",
    "# fig1: GDP 지수 시나리오 팬\n",
    "fig,ax=plt.subplots(figsize=(9,5))\n",
    "for s in SCN: ax.plot(t,gdp_idx[s],color=CLR[s],lw=2,label=f\"{s} {SCN_KO[s]} ({gdp_idx[s][-1]:.0f})\")\n",
    "ax.set_title(\"그림 1. 시나리오별 총량 실질 GDP 경로 (2025Q4=100)\",fontweight=\"bold\")\n",
    "ax.set_ylabel(\"지수 (2025=100)\"); ax.legend(fontsize=9); ax.grid(alpha=.3); ax.margins(x=.02)\n",
    "fig.tight_layout(); fig.savefig(os.path.join(FIG,\"fig1_gdp_fan.png\"),bbox_inches=\"tight\"); plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "metadata": {},
   "execution_count": null,
   "outputs": [],
   "source": [
    "# ===== 셀 14 · fig2: 산업별 노동효과(S3-S2) =====\n",
    "# fig2: 산업별 노동효과(S3-S2)\n",
    "fig,ax=plt.subplots(figsize=(8,9)); ls=labor\n",
    "cols=[\"#c0392b\" if x<0 else \"#1f4e79\" for x in ls.values]\n",
    "ax.barh(range(n),ls.values,color=cols); ax.axvline(0,color=\"#333\",lw=.8)\n",
    "ax.set_yticks(range(n)); ax.set_yticklabels([short(x) for x in ls.index],fontsize=7)\n",
    "ax.set_title(\"그림 2. 노동공급 효과 = 제약완화(S3) - 공급제약(S2)\\n2035 산업 지수 차이 (클수록 인구에 민감)\")\n",
    "ax.set_xlabel(\"2035 지수 차이 (%p)\")\n",
    "fig.tight_layout(); fig.savefig(os.path.join(FIG,\"fig2_industry.png\"),bbox_inches=\"tight\"); plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "metadata": {},
   "execution_count": null,
   "outputs": [],
   "source": [
    "# ===== 셀 15 · 결과 출력 =====\n",
    "json.dump({\"gdp2035\":{s:round(float(gdp_idx[s][-1]),1) for s in SCN},\n",
    "           \"labor_top\":[[short(k),round(float(v),1)] for k,v in labor.tail(4)[::-1].items()]},\n",
    "          open(os.path.join(DATA,\"ch5_results.json\"),\"w\",encoding=\"utf-8\"),ensure_ascii=False,indent=2)\n",
    "print(\"\\n저장: ch5/figures/(fig1,fig2), ch5/data/(scenario_industry_2035.csv, ch5_results.json)\")"
   ]
  }
 ],
 "metadata": {
  "colab": {
   "name": "ch5_scenarios.ipynb",
   "provenance": []
  },
  "kernelspec": {
   "name": "python3",
   "display_name": "Python 3"
  },
  "language_info": {
   "name": "python"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}