# -*- coding: utf-8 -*- """BOKDPM 거시엔진과 시나리오(보고서1 3장) — 계산 함수 모음. 노트북 01(입수)·02(분석)·03(결과)이 이 모듈의 함수만 호출한다. 경로·상수는 config.py에 둔다. 계산 구조 1. 관측자료 : BOKDPM 관측변수 20개(2000Q1–2025Q4, 104분기) — 추정에 실제로 쓴 Dynare 자료파일(유가 = 두바이유) 2. 재추정 : 베이지언 추정(Metropolis–Hastings 2블록 × 30,000회, 앞 30% 제외) — 행태계수 20개 + 충격 표준편차 34개 3. 기준선 : 과거는 칼만 평활값(SmoothedVariables), 미래는 사후최빈값 모수의 40분기 전망(forecast.Mean) 장기 수렴값(정상상태 13개)은 2026-10 정비에서 현재 여건(공식 추정·목표치, 최근 10년 평균)으로 바꾸었다(SS_REVISION). 4. 시나리오 : 설계 레버 5개(인구 가정 + 감쇠 충격 yus·oil·tcost·dshk) → 캘리브레이션 거시 블록(IS·필립스·테일러·UIP·오쿤 5식, Fair–Taylor 반복해법) → 거시경로 5. 산업 외생 : 기준선(BOKDPM) + 시나리오 편차(거시 블록) → 관측 가능한 조합 6종(수준) → 정상형 변환 2026-10 정비(이전 판 2026-10-08 패키지와 달라진 점) - 정상상태 13개 교체(잠재성장률 4.29→1.8%, 물가목표 3.0→2.0%, 유가 상승률 9.88→2.3% 등) → Dynare 재추정·전망 - 실질금리 투입형을 수준에서 1차 차분으로(수준은 결정적 하락 추세를 지녀 성장률 식과 적분차수가 맞지 않는다) - 경제활동인구 기준 경로를 최근 20분기 추세 연장에서 한국고용정보원 중장기 인력수급 전망(2024–2034)의 증가율로 - 인구 레버를 경제활동인구에 반영(저위·고위 추계의 중위 대비 증가율 차이를 누적) - 교역비용 레버를 IS 식의 실효 해외수요(yus − 0.75·tcost)로 전달(0.75 = 설계 교역량 식의 계수비 0.6/0.8) - S3 세계수요에 설계 레버 yus(+0.6, 12분기)를 더함(모든 시나리오가 레버표의 yus를 그대로 쓴다) 원본 프로그램과의 대응 - 1·2·3: bokdpm_pipeline/01_데이터입수.ipynb, 02_모형추정.ipynb, 03_기준선전망_시나리오.ipynb, mod/*.mod - 4 : repo/ch4/code/ch4_dsge.py - 5 : build_exog_scenarios_v2.py(수준), build_varx_exog.py(정상형) — 통합 번들 industry/inputs/ """ from __future__ import annotations import glob import hashlib import os import re import shutil import subprocess import numpy as np import pandas as pd # ── 공통 상수 ───────────────────────────────────────────────────────── VAROBS = ["UNR_kor", "RS_kor", "LCPI_kor", "LGDP_kor", "Y_kor", "BLT_kor", "GE_kor", "CA_kor", "LS_kor", "LDEF_kor", "RS_us", "LCPI_us", "LGDP_us", "Y_us", "BLT_us", "CD_kor", "CB_kor", "LHOUSE_kor", "LOIL_us", "TB_kor"] SCENARIOS = ["S1", "S2", "S3", "S4", "S5"] SCN_KO = {"S1": "기준", "S2": "공급제약+충격", "S3": "제약완화", "S4": "기술도약", "S5": "지정학"} EXOG_COLS = ["유가", "세계수요", "실질금리", "실질환율", "고용_취업자", "우리경기"] EXOG_LABEL = {"유가": "유가", "세계수요": "세계수요", "실질금리": "실질금리", "실질환율": "실질환율", "고용_취업자": "취업자 수", "우리경기": "국내 경기"} EXOG_DEF = { # 산식(관측 가능한 조합)과 산업모형 투입 형태 "유가": ("LOIL_us = 100·ln(두바이유, 달러/배럴)", "1차 차분"), "세계수요": ("Y_us = 100·ln(미국 실질GDP/잠재GDP)", "수준"), "실질금리": ("CB_kor − (LCPI_kor − LCPI_kor₋₄) = 회사채(3년, AA−) 수익률 − 전년동기 CPI 상승률", "1차 차분"), "실질환율": ("LS_kor + LCPI_us − LCPI_kor = 로그 실질 원/달러 환율", "1차 차분"), "고용_취업자": ("경제활동인구 × (1 − UNR_kor/100), 천 명", "로그 차분(×100)"), "우리경기": ("Y_kor = 한국 산출갭(HP 순환성분)", "수준"), } EXOG_TF = {"유가": "diff", "세계수요": "level", "실질금리": "diff", "실질환율": "diff", "고용_취업자": "dlog", "우리경기": "level"} EXOG_TF_PREV = {**EXOG_TF, "실질금리": "level"} # 이전 판(2026-10-08 이전) 투입형 — 비교용 # ── 2026-10 정비: 정상상태(장기 수렴값) 교체 ───────────────────────────── # (파라미터, 이전 값, 새 값, 의미, 근거). 이전 값은 원 모형(2009–2010) 캘리브레이션으로 2000년대 표본 평균이다. # 근거의 '최근 10년 평균'은 관측자료 2016Q1–2025Q4(물가안정목표 2% 도입 이후 40분기) 평균이다(02 3절에서 다시 계산). SS_REVISION = [ ("growth_kor_ss", 4.2925, 1.80, "한국 잠재성장률(%)", "한국은행(2024.12) 「우리 경제의 잠재성장률과 향후 전망」: 2024–2026년 2% 수준, 2025–2029년 연평균 1.8%"), ("pietar_kor_ss", 3.0000, 2.00, "한국 물가목표(CPI, %)", "한국은행 물가안정목표 2%(2016년 이후)"), ("def_kor_ss", 2.1943, 2.00, "GDP 디플레이터 상승률(%)", "물가목표와 같게(최근 10년 평균 2.08%)"), ("rr_kor_bar_ss", 1.5000, 0.55, "한국 중립 실질 정책금리(%)", "한국은행(2024.11) BOK 경제연구 「한국의 중립금리 추정」: 명목 1.8–3.3%, 중간값 2.55% − 물가목표 2%"), ("tprem_kor_ss", 0.3900, 0.18, "CD(91일) − 콜금리(%p)", "최근 10년 평균"), ("tprem_tb_kor_ss", 0.7645, 0.30, "국고채(3년) − 콜금리(%p)", "최근 10년 평균"), ("tprem_cb_kor_ss", 0.7419, 0.67, "회사채(3년, AA−) − 국고채(3년)(%p)", "최근 10년 평균"), ("house_kor_ss", 4.2000, 2.33, "주택가격 상승률(%)", "최근 10년 평균"), ("ge_kor_ss", 7.5885, 6.57, "통합재정 총지출 증가율(%)", "최근 10년 평균"), ("growth_us_ss", 2.3000, 1.80, "미국 잠재성장률(%)", "미 의회예산처(CBO, 2025.1) 잠재GDP 증가율 2027–2035년 1.8%(세계수요 Y_us의 잠재GDP와 같은 출처)"), ("pietar_us_ss", 2.4982, 2.30, "미국 물가목표(CPI, %)", "CBO 소비자물가(CPI-U) 장기 전망 2.3%(PCE 목표 2%에 CPI 괴리를 더한 값)"), ("rr_us_bar_ss", 2.0000, 0.90, "미국 중립 실질 정책금리(%)", "FOMC 경제전망(2026.9) 장기 연방기금금리 중간값 3.2% − 미국 CPI 목표 2.3%"), ("oil_g_us_ss", 9.8763, 2.30, "유가 상승률(달러 기준, %)", "장기 실질유가 불변(= 미국 CPI 목표 2.3%)"), ] SS_UNCHANGED_NOTE = ("UNR_kor_ss(3.6083)는 바꾸지 않았다. UNR_kor_BAR가 확률보행이라 이 값은 전망·우도에 영향이 없다. " "blt_kor_ss·blt_us_ss·s_kor_ss는 0(갭 변수의 정상상태)이다.") def recent_averages(obs: pd.DataFrame, start: str = "2016Q1") -> pd.Series: """정비표의 '최근 10년 평균' 근거값: 관측자료 start–2025Q4의 전년동기 상승률·금리차·지출 증가율 평균.""" m = obs.index >= start g4 = lambda v: (obs[v] - obs[v].shift(4))[m].mean() return pd.Series({"def_kor_ss": g4("LDEF_kor"), "tprem_kor_ss": (obs["CD_kor"] - obs["RS_kor"])[m].mean(), "tprem_tb_kor_ss": (obs["TB_kor"] - obs["RS_kor"])[m].mean(), "tprem_cb_kor_ss": (obs["CB_kor"] - obs["TB_kor"])[m].mean(), "house_kor_ss": g4("LHOUSE_kor"), "ge_kor_ss": obs["GE_kor"][m].mean(), "참고: CPI 상승률": g4("LCPI_kor"), "참고: 실질GDP 성장률": g4("LGDP_kor"), "참고: 두바이유 상승률": g4("LOIL_us")}, name=f"{start}–2025Q4 평균") def unit_root_battery(series: dict[str, pd.Series]) -> pd.DataFrame: """단위근 검정 묶음: ADF·PP·DF-GLS(귀무 = 단위근), KPSS(귀무 = 정상), 상수(c)·상수+추세(ct) 두 설정의 p값.""" import warnings from arch.unitroot import ADF, DFGLS, KPSS, PhillipsPerron rows = [] with warnings.catch_warnings(): warnings.simplefilter("ignore") for name, s in series.items(): x = s.dropna().astype(float).to_numpy() for tr in ("c", "ct"): rows.append({"계열": name, "결정항": tr, "관측수": len(x), "ADF p": ADF(x, trend=tr, method="aic").pvalue, "PP p": PhillipsPerron(x, trend=tr).pvalue, "DF-GLS p": DFGLS(x, trend=tr, method="aic").pvalue, "KPSS p": KPSS(x, trend=tr).pvalue}) return pd.DataFrame(rows).set_index(["계열", "결정항"]) def steady_state_table(mod_path: str, prev_mod_path: str | None = None) -> pd.DataFrame: """정상상태 교체표 — 이전 값·새 값을 .mod 파일에서 직접 읽어 SS_REVISION과 대조한다.""" names = [r[0] for r in SS_REVISION] new = calibrated_values(mod_path, names) prev = calibrated_values(prev_mod_path, names) if prev_mod_path else {} rows = [] for n, old, val, lab, basis in SS_REVISION: assert abs(new[n] - val) < 1e-9, f"{n}: .mod 값 {new[n]} ≠ 정비표 {val}" if prev: assert abs(prev[n] - old) < 1e-9, f"{n}: 이전 .mod 값 {prev[n]} ≠ 정비표 {old}" rows.append({"파라미터": n, "의미": lab, "이전 값": old, "새 값": val, "근거": basis}) return pd.DataFrame(rows).set_index("파라미터") BK = ["LOIL_us", "Y_us", "Y_kor", "CB_kor", "LCPI_kor", "LCPI_us", "UNR_kor", "LS_kor"] # 외생 구성에 쓰는 BOKDPM 변수 def sha256(path: str) -> str: h = hashlib.sha256() with open(path, "rb") as f: for chunk in iter(lambda: f.read(1 << 20), b""): h.update(chunk) return h.hexdigest() def quarters(start: str, n: int) -> list[str]: """'2000Q1'부터 n개 분기 표기.""" y, q = int(start[:4]), int(start[-1]) - 1 return [f"{y + (q + i) // 4}Q{(q + i) % 4 + 1}" for i in range(n)] def write_exact(df: pd.DataFrame, path: str, **kw) -> None: """부동소수점을 손실 없이 기록한다(17자리).""" df.to_csv(path, float_format="%.17g", encoding="utf-8-sig", **kw) def read_exact(path: str, **kw) -> pd.DataFrame: """write_exact로 기록한 파일을 비트 단위로 같게 읽는다.""" return pd.read_csv(path, float_precision="round_trip", encoding="utf-8-sig", **kw) # ── 1. 입수: 관측 20변수 ─────────────────────────────────────────────── def read_dynare_data(path: str) -> pd.DataFrame: """Dynare 자료파일(.m, 변수명 = [값; …];) → 행=분기, 열=관측변수.""" txt = open(path, encoding="utf-8").read() cols = {} for m in re.finditer(r"^(\w+)\s*=\s*\[(.*?)\];", txt, flags=re.M | re.S): vals = [float(x) for x in m.group(2).replace(";", " ").split()] cols[m.group(1)] = vals df = pd.DataFrame(cols) df.index = quarters("2000Q1", len(df)) df.index.name = "분기" return df def parse_mod(path: str) -> dict: """.mod 파일에서 변수·충격·방정식 수, 관측변수, 추정 대상과 사전분포, 추정 옵션을 읽는다.""" raw = open(path, encoding="utf-8", errors="ignore").read() txt = re.sub(r"/\*.*?\*/", " ", raw, flags=re.S) txt = re.sub(r"//[^\n]*", " ", txt) def block(kw): m = re.search(rf"^\s*{kw}\b(.*?);", txt, flags=re.M | re.S) return [t for t in re.split(r"[\s,]+", m.group(1)) if t] if m else [] model = re.search(r"^\s*model\s*;(.*?)^\s*end\s*;", txt, flags=re.M | re.S).group(1) est = re.search(r"^\s*estimated_params\s*;(.*?)^\s*end\s*;", txt, flags=re.M | re.S).group(1) priors = [] for line in est.split(";"): parts = [p.strip() for p in line.split(",")] if len(parts) < 4 or not parts[0]: continue name = parts[0].replace("stderr", "").strip() priors.append({"파라미터": name, "구분": "shocks_std" if parts[0].startswith("stderr") else "parameters", "사전분포": parts[1].replace("_pdf", ""), "사전평균": float(parts[2]), "사전표준편차": float(parts[3])}) opt = re.search(r"estimation\(([^;]*)\);", txt) return {"내생변수": len(block("var")), "외생충격": len(block("varexo")), "방정식": model.count(";"), "관측변수": block("varobs"), "사전분포": pd.DataFrame(priors), "추정옵션": opt.group(1) if opt else ""} # ── 1-1. 입수(ONLINE=1): ECOS·FRED 재수집 ──────────────────────────────── ECOS_BASE = "https://ecos.bok.or.kr/api" FRED_BASE = "https://api.stlouisfed.org/fred" def fetch_ecos(key: str, stat: str, item: str, freq: str, start: str, end: str | None = None, n: int = 15000) -> pd.Series: """한국은행 ECOS 통계 1계열 → pandas Series(인덱스=원 시점 표기).""" import requests from datetime import datetime if end is None: t = datetime.today() end = {"A": str(t.year), "Q": f"{t.year}Q{(t.month - 1) // 3 + 1}", "M": t.strftime("%Y%m"), "D": t.strftime("%Y%m%d")}[freq] url = f"{ECOS_BASE}/StatisticSearch/{key}/json/kr/1/{n}/{stat}/{freq}/{start}/{end}/{item}" j = requests.get(url, timeout=60).json() if "StatisticSearch" not in j: raise RuntimeError(f"ECOS 응답 오류({stat}/{item}): {j.get('RESULT', j)}") rows = j["StatisticSearch"]["row"] return pd.Series({r["TIME"]: float(r["DATA_VALUE"]) for r in rows if r["DATA_VALUE"] not in (None, "", "-")}).sort_index() def fetch_fred(key: str, series_id: str, start: str = "1999-01-01") -> pd.Series: """FRED 1계열 → pandas Series(인덱스=날짜).""" import requests url = (f"{FRED_BASE}/series/observations?series_id={series_id}&api_key={key}" f"&file_type=json&observation_start={start}") j = requests.get(url, timeout=60).json() if "observations" not in j: raise RuntimeError(f"FRED 응답 오류({series_id}): {j.get('error_message', '')}") s = pd.Series({o["date"]: float(o["value"]) for o in j["observations"] if o["value"] != "."}) s.index = pd.to_datetime(s.index) return s.sort_index() def _q_ts(s: pd.Series) -> pd.Series: s = s.copy() s.index = pd.PeriodIndex(s.index, freq="Q").to_timestamp() return s def _m_to_q(s: pd.Series) -> pd.Series: s = s.copy() fmt = "%Y%m" if len(str(s.index[0])) == 6 else "%Y%m%d" s.index = pd.to_datetime(s.index.astype(str), format=fmt) return s.resample("QS").mean() # 관측 20변수의 원천(추정 자료파일 기준). LOIL_us는 두바이유, BLT_us는 원천 미확정(스냅샷 값 사용). SOURCES = { "LGDP_kor": "ECOS 200Y108/10101 실질GDP(계절조정), 100·ln", "Y_kor": "ECOS 200Y108/10101 실질GDP의 HP 순환성분(λ=1600)", "LDEF_kor": "ECOS 200Y107/10101 명목GDP ÷ 실질GDP, 100·ln", "CA_kor": "ECOS 301Y013/000000 경상수지 ÷ 명목GDP × 100", "BLT_kor": "ECOS 514Y001/AA 대출태도(2002Q1~), 2000–2001은 원본 BOKDPM 자료", "UNR_kor": "ECOS 901Y027/I61BC 실업률", "LCPI_kor": "ECOS 901Y009/0 소비자물가지수(월→분기 평균), 100·ln", "LS_kor": "ECOS 731Y004/0000001 원/달러 환율(월→분기 평균), 100·ln", "LHOUSE_kor": "ECOS 901Y062/P63A 주택매매가격지수(월→분기 평균), 100·ln", "RS_kor": "ECOS 722Y001/0101000 콜금리(월→분기 평균)", "CD_kor": "ECOS 817Y002/010502000 CD(91일)(일→분기 평균)", "TB_kor": "ECOS 817Y002/010200000 국고채(3년)(일→분기 평균)", "CB_kor": "ECOS 817Y002/010300000 회사채(3년, AA−)(일→분기 평균)", "GE_kor": "ECOS 901Y013/B01 통합재정 총지출(누계→분기 흐름), 4분기 로그차분, 2000년은 원본 BOKDPM 자료", "LGDP_us": "FRED GDPC1, 100·ln", "Y_us": "FRED GDPC1·GDPPOT, 100·ln(실질/잠재)", "LCPI_us": "FRED CPIAUCSL(월→분기 평균), 100·ln", "RS_us": "FRED FEDFUNDS(월→분기 평균)", "LOIL_us": "FRED POILDUBUSDM 두바이유(월→분기 평균), 100·ln", "BLT_us": "원천 미확정 합성 대출태도 지수(2026-08-01 자료파일 값 고정; 참고 계열 FRED DRTSCILM)", } def refetch_observables(key_bok: str, key_fred: str, snapshot: pd.DataFrame) -> tuple[pd.DataFrame, pd.Series]: """bokdpm_pipeline/01_데이터입수.ipynb의 수집·변환 절차로 관측 20변수를 다시 만든다. 추정 자료파일과 같은 원천을 쓰도록 유가는 두바이유(POILDUBUSDM)로 받는다. ECOS가 제공하지 않는 구간(KR_BLT 2000–2001, GE_kor 2000년)과 원천이 확정되지 않은 BLT_us는 스냅샷 값을 그대로 쓴다. 반환: (관측 20변수 2000Q1–2025Q4, 참고 계열 FRED DRTSCILM 분기 평균) """ from statsmodels.tsa.filters.hp_filter import hpfilter S = "2000Q1" e = lambda st, it, fr, start, **k: fetch_ecos(key_bok, st, it, fr, start, **k) kr = { "KR_GDP_REAL": _q_ts(e("200Y108", "10101", "Q", S)), "KR_GDP_NOM": _q_ts(e("200Y107", "10101", "Q", S)), "KR_CA": _q_ts(e("301Y013", "000000", "Q", S)), "KR_BLT": _q_ts(e("514Y001", "AA", "Q", "2002Q1")), "KR_UNR": _q_ts(e("901Y027", "I61BC", "Q", S)), "KR_CPI": _m_to_q(e("901Y009", "0", "M", "200001")), "KR_FX": _m_to_q(e("731Y004", "0000001", "M", "200001")), "KR_HOUSE": _m_to_q(e("901Y062", "P63A", "M", "200001")), "KR_CALL": _m_to_q(e("722Y001", "0101000", "M", "200001")), "KR_CD": _m_to_q(e("817Y002", "010502000", "D", "20000101")), "KR_TB": _m_to_q(e("817Y002", "010200000", "D", "20000101")), "KR_CB": _m_to_q(e("817Y002", "010300000", "D", "20000101")), } f = lambda sid: fetch_fred(key_fred, sid, "2000-01-01") us_real, us_pot = _q_ts(f("GDPC1")), _q_ts(f("GDPPOT")) us_pot = us_pot[us_pot.index <= us_real.index.max()] us = {"US_GDP_REAL": us_real, "US_GDP_POT": us_pot, "US_CPI": f("CPIAUCSL").resample("QS").mean(), "US_RATE": f("FEDFUNDS").resample("QS").mean(), "US_OIL": f("POILDUBUSDM").resample("QS").mean()} drtscilm = f("DRTSCILM").resample("QS").mean() raw = pd.concat({**kr, **us}, axis=1).sort_index() lg = 100 * np.log(raw["KR_GDP_REAL"].dropna()) cyc, _ = hpfilter(lg.values, lamb=1600) # 표본 끝까지 쓴 HP 필터(원본과 같은 방식) obs = pd.DataFrame(index=raw.index) obs["LGDP_kor"] = 100 * np.log(raw["KR_GDP_REAL"]) obs["LCPI_kor"] = 100 * np.log(raw["KR_CPI"]) obs["LDEF_kor"] = 100 * np.log(raw["KR_GDP_NOM"] / raw["KR_GDP_REAL"] * 100) obs["LS_kor"] = 100 * np.log(raw["KR_FX"]) obs["LHOUSE_kor"] = 100 * np.log(raw["KR_HOUSE"]) obs["LGDP_us"] = 100 * np.log(raw["US_GDP_REAL"]) obs["LCPI_us"] = 100 * np.log(raw["US_CPI"]) obs["LOIL_us"] = 100 * np.log(raw["US_OIL"]) for c, r in [("RS_kor", "KR_CALL"), ("CD_kor", "KR_CD"), ("TB_kor", "KR_TB"), ("CB_kor", "KR_CB"), ("UNR_kor", "KR_UNR"), ("BLT_kor", "KR_BLT"), ("RS_us", "US_RATE")]: obs[c] = raw[r] gm = e("901Y013", "B01", "M", "199901") # 통합재정 총지출(월, 연초 누계) gm.index = pd.to_datetime(gm.index, format="%Y%m") qe = gm[gm.index.month.isin([3, 6, 9, 12])] flow = {} for ts, v in qe.items(): pv = qe[(qe.index.year == ts.year) & (qe.index.month == ts.month - 3)] flow[pd.Timestamp(ts.year, ts.month - 2, 1)] = v - (pv.iloc[0] if len(pv) else 0) obs["GE_kor"] = (100 * np.log(pd.Series(flow).sort_index()).diff(4)).reindex(obs.index) obs["CA_kor"] = raw["KR_CA"] / raw["KR_GDP_NOM"] * 100 obs["Y_kor"] = pd.Series(cyc, index=lg.index) obs["Y_us"] = 100 * np.log(raw["US_GDP_REAL"] / raw["US_GDP_POT"]) obs = obs.loc["2000-01-01":"2025-10-01"].copy() obs.index = quarters("2000Q1", len(obs)) obs.index.name = "분기" obs["BLT_us"] = snapshot["BLT_us"] # 원천 미확정 — 스냅샷 값 obs.loc["2000Q1":"2001Q4", "BLT_kor"] = snapshot.loc["2000Q1":"2001Q4", "BLT_kor"] obs.loc["2000Q1":"2000Q4", "GE_kor"] = snapshot.loc["2000Q1":"2000Q4", "GE_kor"] drtscilm = drtscilm.loc["2000-01-01":"2025-10-01"] drtscilm.index = quarters("2000Q1", len(drtscilm)) return obs[VAROBS], drtscilm FIXED_FROM_SNAPSHOT = {"BLT_us": "전 구간 스냅샷 값(원천 미확정)", "BLT_kor": "2000–2001 스냅샷 값(원본 BOKDPM)", "GE_kor": "2000년 스냅샷 값(원본 BOKDPM)"} def compare_observables(new: pd.DataFrame, snap: pd.DataFrame) -> pd.DataFrame: """재수집 자료와 스냅샷을 변수별로 비교한다. 증분 상관은 수준 차이와 형태 차이를 구분하려고 함께 낸다.""" rows = [] for v in VAROBS: a, b = snap[v].astype(float), new[v].reindex(snap.index).astype(float) d = (b - a) m = d.notna() dc = np.corrcoef(np.diff(a[m]), np.diff(b[m]))[0, 1] if m.sum() > 8 else np.nan mx = float(d.abs().max()) verdict = "동일" if mx < 1e-4 else ("개정(증분 상관 ≥ 0.99)" if dc >= 0.99 else ("수준·정의 차이" if dc >= 0.9 else "원천 차이")) rows.append({"변수": v, "비교 분기": int(m.sum()), "최대절대차": mx, "평균절대차": float(d.abs().mean()), "증분 상관": dc, "판정": verdict, "스냅샷 고정 구간": FIXED_FROM_SNAPSHOT.get(v, "")}) return pd.DataFrame(rows).set_index("변수") # ── 2. 분석: 사후분포 ───────────────────────────────────────────────── def load_kept_draws(path: str) -> tuple[list[str], np.ndarray, np.ndarray]: """연속 반복을 묶어 저장한 사후 표본(번인 제외분)을 풀어 (이름, 표본행렬, 블록번호)로 반환한다.""" d = read_exact(path) names = [c for c in d.columns if c not in ("블록", "반복수")] rep = d["반복수"].to_numpy() X = np.repeat(d[names].to_numpy(), rep, axis=0) blk = np.repeat(d["블록"].to_numpy(), rep) return names, X, blk def posterior_moments(x: np.ndarray, conf: float = 0.9) -> tuple[float, float, float, float]: """Dynare 7.1 posterior_moments.m과 같은 계산: 평균, 표준편차(n−1), 최단 HPD 구간.""" x = np.sort(np.asarray(x, float)) n = len(x) h = int(round((1 - conf) * n)) kk = x[n - h - 1: n - 1] - x[:h] # MATLAB: kk(ii) = xx(n−h+ii−1) − xx(ii), ii = 1…h i = int(np.argmin(kk)) return float(x.mean()), float(x.std(ddof=1)), float(x[i]), float(x[i] + kk[i]) def posterior_summary(names: list[str], X: np.ndarray, priors: pd.DataFrame, conf: float = 0.9) -> pd.DataFrame: """추정 대상 54개의 사전분포와 사후 평균·표준편차·90% HPD 구간(estimation_results.csv와 같은 배열).""" rows = [] for j, n in enumerate(names): m, s, lo, hi = posterior_moments(X[:, j], conf) rows.append({"파라미터": n, "사후평균": m, "표준편차": s, "HPD_low": lo, "HPD_high": hi}) out = priors.merge(pd.DataFrame(rows), on="파라미터", how="right", validate="one_to_one") out["순서"] = out["구분"].map({"parameters": 0, "shocks_std": 1}) return out.sort_values("순서", kind="stable").drop(columns="순서").reset_index(drop=True) def mcmc_diagnostics(names: list[str], X: np.ndarray, blk: np.ndarray) -> pd.DataFrame: """블록별 이동 비율(유지 표본의 채택률 근사)과 Gelman–Rubin 수렴 지표(PSRF).""" chains = [X[blk == b] for b in np.unique(blk)] n = min(len(c) for c in chains) C = np.stack([c[:n] for c in chains]) # (블록, n, 모수) W = C.var(axis=1, ddof=1).mean(axis=0) B = n * C.mean(axis=1).var(axis=0, ddof=1) psrf = np.sqrt(((n - 1) / n * W + B / n) / W) move = {f"이동비율_블록{b}": float(np.any(np.diff(c, axis=0) != 0, axis=1).mean()) for b, c in zip(np.unique(blk), chains)} out = pd.DataFrame({"파라미터": names, "PSRF": psrf}) out.attrs.update(move) return out # ── 2-1. 분석(FULL=1): Octave + Dynare 실행과 결과 추출 ──────────────────── def run_octave(script: str, cwd: str, octave: str, dynare: str, log_path: str, done_marker: str) -> str: """Octave에서 Dynare 러너를 실행하고 로그를 남긴다. 완료 표지가 없으면 멈춘다.""" res = subprocess.run([octave, "--no-gui", script], cwd=cwd, capture_output=True, text=True, encoding="utf-8", errors="replace", env={**os.environ, "DYNARE_MATLAB": dynare}) with open(log_path, "w", encoding="utf-8") as f: f.write(res.stdout) f.write("\n[stderr]\n" + (res.stderr or "")) if done_marker not in res.stdout: raise RuntimeError(f"Octave/Dynare 실행이 끝나지 않았다({script}). 로그: {log_path}\n{(res.stderr or '')[-1500:]}") return res.stdout def prepare_dynare_workdir(workdir: str, dynare_src: str, datafile: str) -> None: """모형·러너·자료파일을 작업 폴더로 복사한다(이전 Dynare 산출물은 지운다).""" os.makedirs(workdir, exist_ok=True) for name in ("BOKDPM_2025", "BOKDPM_forecast", "+BOKDPM_2025", "+BOKDPM_forecast"): shutil.rmtree(os.path.join(workdir, name), ignore_errors=True) for f in ("BOKDPM_2025.mod", "BOKDPM_forecast.mod", "run_est_local.m", "run_forecast.m"): shutil.copy(os.path.join(dynare_src, f), os.path.join(workdir, f)) shutil.copy(datafile, os.path.join(workdir, "BOKDPM_data_2025.m")) def _loadmat(path): import warnings import scipy.io as sio with warnings.catch_warnings(): warnings.simplefilter("ignore") return sio.loadmat(path, struct_as_record=False, squeeze_me=True) def extract_kept_draws(dyn_dir: str, fname: str = "BOKDPM_2025") -> pd.DataFrame: """metropolis/*.mat에서 번인 이후 표본만 골라 연속 반복을 묶은 표로 만든다(Dynare KeepedDraws 기준).""" mdir = os.path.join(dyn_dir, fname, "metropolis") hist = sorted(glob.glob(os.path.join(mdir, f"{fname}_mh_history_*.mat")))[-1] rec = _loadmat(hist)["record"] first_file, first_line = int(rec.KeepedDraws.FirstMhFile), int(rec.KeepedDraws.FirstLine) names = [str(n).replace("stderr ", "") for n in _loadmat(os.path.join(dyn_dir, fname, "Output", f"{fname}_results.mat"))["bayestopt_"].name] nblck = int(rec.Nblck) parts = [] for b in range(1, nblck + 1): files = sorted(glob.glob(os.path.join(mdir, f"{fname}_mh*_blck{b}.mat")), key=lambda p: int(re.search(r"_mh(\d+)_blck", p).group(1))) xs = [] for p in files: k = int(re.search(r"_mh(\d+)_blck", p).group(1)) if k < first_file: continue x = np.atleast_2d(_loadmat(p)["x2"]) xs.append(x[first_line - 1:] if k == first_file else x) x = np.vstack(xs) new = np.r_[True, np.any(np.diff(x, axis=0) != 0, axis=1)] idx = np.flatnonzero(new) df = pd.DataFrame(x[idx], columns=names) df.insert(0, "반복수", np.diff(np.r_[idx, len(x)])) df.insert(0, "블록", b) parts.append(df) return pd.concat(parts, ignore_index=True) def _endo_names(m) -> list[str]: return [str(n) for n in m["M_"].endo_names][: int(m["M_"].orig_endo_nbr)] def extract_smoothed(results_mat: str, start: str = "2000Q1") -> pd.DataFrame: """추정 결과(.mat)의 칼만 평활값 — 원 모형 내생변수 70개.""" m = _loadmat(results_mat) sv = m["oo_"].SmoothedVariables names = _endo_names(m) df = pd.DataFrame({n: np.asarray(getattr(sv, n), float).ravel() for n in names}) df.index = quarters(start, len(df)) df.index.name = "분기" return df def extract_forecast(results_mat: str, start: str = "2026Q1") -> pd.DataFrame: """전망 결과(.mat)의 forecast.Mean — 원 모형 내생변수 70개, 40분기.""" m = _loadmat(results_mat) fm = m["oo_"].forecast.Mean names = _endo_names(m) df = pd.DataFrame({n: np.asarray(getattr(fm, n), float).ravel() for n in names}) df.index = quarters(start, len(df)) df.index.name = "분기" return df def extract_posterior_mode(results_mat: str) -> pd.Series: """사후최빈값(구조계수 + 충격 표준편차).""" pm = _loadmat(results_mat)["oo_"].posterior_mode out = {} for cat in ("parameters", "shocks_std"): st = getattr(pm, cat) for n in st._fieldnames: out[n] = float(getattr(st, n)) return pd.Series(out, name="사후최빈값") # ── 3. 기준선 요약 ───────────────────────────────────────────────────── def annual_baseline(smoothed: pd.DataFrame, forecast: pd.DataFrame, years=range(2024, 2036)) -> pd.DataFrame: """평활값(과거)과 전망(미래)을 이어 연간 지표를 만든다: 실질GDP 성장률·CPI 상승률(연평균 기준), 콜금리·실업률·산출갭(연평균).""" keep = ["LGDP_kor", "LCPI_kor", "RS_kor", "UNR_kor", "Y_kor", "LOIL_us", "Y_us"] q = pd.concat([smoothed[keep], forecast[keep]]) q["연도"] = [int(s[:4]) for s in q.index] g = q.groupby("연도") lvl = pd.DataFrame({"GDP": g["LGDP_kor"].apply(lambda v: np.exp(v / 100).sum()), "CPI": g["LCPI_kor"].apply(lambda v: np.exp(v / 100).mean()), "OIL": g["LOIL_us"].apply(lambda v: np.exp(v / 100).mean())}) out = pd.DataFrame({ "실질GDP 성장률(%)": 100 * lvl["GDP"].pct_change(), "CPI 상승률(%)": 100 * lvl["CPI"].pct_change(), "콜금리(%)": g["RS_kor"].mean(), "실업률(%)": g["UNR_kor"].mean(), "산출갭(%)": g["Y_kor"].mean(), "두바이유(달러/배럴)": lvl["OIL"], "미국 산출갭(%)": g["Y_us"].mean(), }) out["구분"] = ["실적(평활)" if y <= 2025 else "전망" for y in out.index] return out.loc[list(years)] # ── 4. 시나리오: 설계 레버와 캘리브레이션 거시 블록(ch4_dsge.py 이식) ───────── CALIB = dict(b1=0.50, b2=0.10, b3=0.10, b4=0.30, b6=0.10, lam1=0.55, lam2=0.15, lam3=0.65, g1=0.50, g2=1.50, g3=0.50, a1=0.80, a2=0.30, om4=0.25, om7=0.50, om8=0.20, rrbar=1.5, pietar=2.0) H, BUF = 40, 24 N_SOLVE = H + BUF LAB_ASSUMP = {"S1": "중위", "S2": "저위", "S3": "고위", "S4": "중위", "S5": "중위"} LEVERS = { # (최고치, 감쇠 분기 k) — 2026Q1 최고치에서 k분기에 걸쳐 직선으로 0이 된다 "S1": {}, "S2": {"yus": (-1.2, 12), "oil": (-8.0, 8), "tcost": (+1.5, 12), "dshk": (-1.0, 12)}, "S3": {"yus": (+0.6, 12), "oil": (-3.0, 8), "tcost": (-1.0, 12), "dshk": (+0.5, 12)}, "S4": {"yus": (+0.8, 20), "dshk": (+0.6, 28)}, "S5": {"oil": (+15.0, 8), "tcost": (+2.0, 10), "yus": (-0.6, 10), "dshk": (-0.8, 10)}, } def dec(k: float, n: int = N_SOLVE) -> np.ndarray: """감쇠함수 dec(k)_t = clip(1 − t/k, 0, 1), t = 0(2026Q1), 1, …""" return np.clip(1 - np.arange(n) / k, 0, 1) def scenario_levers(s: str, n: int = N_SOLVE) -> dict[str, np.ndarray]: """시나리오 s의 감쇠 충격 4종(yus 미국수요·oil 유가·tcost 교역비용·dshk 내수) 경로.""" out = {k: np.zeros(n) for k in ("yus", "oil", "tcost", "dshk")} for k, (peak, kk) in LEVERS[s].items(): out[k] = peak * dec(kk, n) if s == "S4": # 원본 식: dshk = +0.6·1·dec(28) out["dshk"] = +0.6 * np.ones(n) * dec(28, n) return out def lever_table() -> pd.DataFrame: rows = [] for s in SCENARIOS: r = {"시나리오": f"{s} {SCN_KO[s]}", "인구(생산연령인구 추계)": LAB_ASSUMP[s]} for k, lab in [("yus", "미국수요 yus"), ("oil", "유가 oil"), ("tcost", "교역비용 tcost"), ("dshk", "내수 dshk")]: r[lab] = f"{LEVERS[s][k][0]:+.1f} ({LEVERS[s][k][1]}분기)" if k in LEVERS[s] else "—" rows.append(r) return pd.DataFrame(rows).set_index("시나리오") def labor_force_growth(pop: pd.DataFrame) -> dict[str, np.ndarray]: """KOSIS 생산연령인구 추계(저·중·고위)의 연간 로그증가율(%)을 분기 증가율로 나눈다.""" pop = pop.set_index("연도") col = {"중위": "중위_천명", "고위": "고위_천명", "저위": "저위_천명"} out = {} for a, c in col.items(): ann = {y: np.log(pop.loc[y, c] / pop.loc[y - 1, c]) * 100 for y in range(2026, 2036)} out[a] = np.array([ann[2026 + h // 4] / 4 for h in range(H)]) return out def population_lever(pop: pd.DataFrame, s: str) -> np.ndarray: """인구 레버(40분기, 100·log): 시나리오 인구 가정과 중위 추계의 분기 증가율 차이를 2026Q1부터 누적한 값. 경제활동인구에 exp(값/100)을 곱해 시나리오별 노동공급 수준 차이를 만든다. 2025Q4에는 모든 시나리오가 같다. """ lfg = labor_force_growth(pop) return np.cumsum(lfg[LAB_ASSUMP[s]] - lfg["중위"]) def keis_growth(labor: pd.DataFrame, keis: pd.DataFrame) -> dict[int, float]: """한국고용정보원 중장기 인력수급 전망(2024–2034)의 경제활동인구 증감 → 연도별 연간 로그증가율(%). 기준 수준은 관측 경제활동인구(ECOS, 계절조정)의 2024년 평균이다. 2024→2029년 증감과 2029→2034년 증감을 각각 5년 균등 로그증가율로 바꾸고, 2035년은 후기(2029–2034) 증가율을 잇는다. """ k = keis.set_index("구간")["경제활동인구_증감_천명"] lf24 = labor["경활"].loc[[f"2024Q{i}" for i in range(1, 5)]].mean() lf29 = lf24 + k["2024-2029"] lf34 = lf29 + k["2029-2034"] g1, g2 = 100 * np.log(lf29 / lf24) / 5, 100 * np.log(lf34 / lf29) / 5 return {y: (g1 if y <= 2029 else g2) for y in range(2025, 2036)} def project_labor_force_keis(labor: pd.DataFrame, keis: pd.DataFrame, q: list[str]) -> np.ndarray: """경제활동인구 기준 경로: 관측 구간(2026Q2까지)은 그대로, 이후는 인력수급 전망의 연간 증가율(분기 1/4씩)로 잇는다.""" s = labor["경활"].reindex(q) y = np.log(s.values) last = np.where(np.isfinite(y))[0][-1] g = keis_growth(labor, keis) for i in range(last + 1, len(q)): y[i] = y[i - 1] + g[int(q[i][:4])] / 400 return np.exp(y) KAPPA_TCOST = 0.6 / 0.8 # 교역비용 1단위 = 미국수요 0.75단위의 교역량 감소(설계 교역량 식 trd_g = 1 + 0.8·yus − 0.6·tcost) def solve_korea_block(yus, oil, dshk, tcost=None, iters: int = 400): """IS·필립스·테일러·UIP 4식을 Fair–Taylor(가우스–자이델 시간 반복)로 푼다. 반환: Y, PIE, RS, RR, S(앞 40분기). tcost가 주어지면 IS 식의 해외수요 항을 실효 해외수요 yus − KAPPA_TCOST·tcost로 바꾼다(2026-10 정비). tcost=None이면 원본(ch4_dsge.py)과 같다. """ c = CALIB b1, b2, b3, b4, b6 = c["b1"], c["b2"], c["b3"], c["b4"], c["b6"] lam1, lam2, lam3, g1, g2, g3 = c["lam1"], c["lam2"], c["lam3"], c["g1"], c["g2"], c["g3"] om4, om7, om8, rrbar, pietar = c["om4"], c["om7"], c["om8"], c["rrbar"], c["pietar"] N = N_SOLVE ext = yus if tcost is None else yus - KAPPA_TCOST * tcost # 실효 해외수요 Y = np.zeros(N); PIE = np.full(N, pietar); RS = np.full(N, rrbar + pietar); S = np.zeros(N) for _ in range(iters): Y0, PIE0, RS0, S0 = Y.copy(), PIE.copy(), RS.copy(), S.copy() for t in range(N): yl = Y[t - 1] if t > 0 else 0.0; yf = Y[t + 1] if t + 1 < N else 0.0 rsl = RS[t - 1] if t > 0 else rrbar + pietar rgap_l = (RS[t - 1] if t > 0 else rrbar + pietar) - pietar - rrbar Y[t] = (b1 * yl + b2 * yf - b3 * rgap_l + b4 * ext[t] - b6 * oil[t] / 10 + dshk[t]) pif = PIE[t + 1] if t + 1 < N else pietar; pil = PIE[t - 1] if t > 0 else pietar PIE[t] = lam3 * pietar + (1 - lam3) * (lam1 * pif + (1 - lam1) * pil + lam2 * 4 * yl) pie_e = PIE[t + 1] if t + 1 < N else pietar RS[t] = g1 * rsl + (1 - g1) * (rrbar + pietar + g2 * (pie_e - pietar) + g3 * 4 * Y[t]) sf = S[t + 1] if t + 1 < N else 0.0; sl = S[t - 1] if t > 0 else 0.0 rgap = RS[t] - pietar - rrbar S[t] = (1 - 0.1) * (om7 * sf + (1 - om7) * sl - om4 * rgap) + om8 * (-Y[t]) if max(np.abs(Y - Y0).max(), np.abs(RS - RS0).max(), np.abs(S - S0).max()) < 1e-8: break RR = RS - PIE return Y[:H], PIE[:H], RS[:H], RR[:H], S[:H] def dsge_paths(pop: pd.DataFrame, lf_base: np.ndarray, drop: str | None = None) -> tuple[pd.DataFrame, pd.DataFrame]: """설계 레버 → 거시 블록 → (시나리오별 외생 성장 경로, 거시경로). ch4_dsge.py의 식·반올림(소수 셋째 자리)에 2026-10 정비를 더했다. lf_base: 경제활동인구 기준 경로(2025Q4–2035Q4, 41개, 천 명) — project_labor_force_keis의 끝 41분기. 정비 내용: ① IS 식 해외수요 = yus − 0.75·tcost, ② 고용 증가율 lab_g = 100·Δln(경제활동인구 × 인구 레버) − Δ실업갭 (원본은 생산연령인구 증가율 + 참가율 추세 − Δ실업갭으로, 외생 '취업자 수'와 다른 노동공급을 썼다). drop: 레버 전달 점검용. '인구'이면 모든 시나리오에 중위 인구를, 'yus'·'oil'·'tcost'·'dshk'이면 그 충격을 0으로 둔다. """ c = CALIB assert len(lf_base) == H + 1 qs = quarters("2026Q1", H) rows, macro = [], [] for s in SCENARIOS: lv = scenario_levers(s) if drop in lv: lv[drop] = np.zeros(N_SOLVE) yus, oil_sh, tcost, dshk = lv["yus"], lv["oil"], lv["tcost"], lv["dshk"] Y, PIE, RS, RR, S = solve_korea_block(yus, oil_sh, dshk, tcost) UNRgap = np.zeros(H) for t in range(H): UNRgap[t] = c["a1"] * (UNRgap[t - 1] if t > 0 else 0) - c["a2"] * Y[t] dUNR = np.diff(np.concatenate([[0], UNRgap])) oil_g = 0.5 + oil_sh[:H] / 8 * 4 trd_g = 1.0 + 0.8 * yus[:H] - 0.6 * tcost[:H] d_rrate = np.diff(np.concatenate([[RR[0]], RR])) rfx_g = np.diff(np.concatenate([[0], S])) plev = np.zeros(H) if drop == "인구" else population_lever(pop, s) lab_g = np.diff(100 * np.log(lf_base)) + np.diff(np.concatenate([[0], plev])) - dUNR for h in range(H): rows.append([s, qs[h], round(oil_g[h], 3), round(trd_g[h], 3), round(d_rrate[h], 3), round(rfx_g[h], 3), round(lab_g[h], 3)]) macro.append([s, qs[h], round(Y[h], 3), round(PIE[h], 3), round(RS[h], 3), round(RR[h], 3), round(S[h], 3), round(UNRgap[h], 3)]) ex = pd.DataFrame(rows, columns=["시나리오", "분기", "oil_g", "trd_g", "d_rrate", "rfx_g", "lab_g"]) mac = pd.DataFrame(macro, columns=["시나리오", "분기", "산출갭Y", "인플레PIE", "정책금리RS", "실질금리RR", "실질환율갭S", "실업갭UNR"]) return ex, mac # ── 5. 산업 외생 6종 ─────────────────────────────────────────────────── def world_demand_dev(s: str, n: int = H) -> np.ndarray: """세계수요 외생에 더하는 편차 = 설계 레버 yus(2026-10 정비: S3 +0.6·dec(12) 포함, 모든 시나리오가 레버표를 따른다).""" return scenario_levers(s, n)["yus"] def world_demand_dev_prev(s: str, n: int = H) -> np.ndarray: """이전 판(build_exog_scenarios_v2.py, 2026-08-28)의 세계수요 편차 — S3를 빠뜨렸다. 비교용.""" z = np.zeros(n) d = lambda k: np.clip(1 - np.arange(n) / k, 0, 1) if s == "S2": return -1.2 * d(12) if s == "S4": return +0.8 * d(20) if s == "S5": return -0.6 * d(10) return z def design_consistency() -> pd.DataFrame: """설계 레버 yus와 세계수요 외생에 더한 편차를 시나리오별로 대조한다(이전 판 편차 함께 표시).""" rows = [] for s in SCENARIOS: lever = scenario_levers(s, H)["yus"] used, prev = world_demand_dev(s), world_demand_dev_prev(s) rows.append({"시나리오": f"{s} {SCN_KO[s]}", "레버 yus 2026Q1": lever[0], "세계수요 편차 2026Q1": used[0], "이전 판 편차 2026Q1": prev[0], "40분기 최대 차이": float(np.abs(lever - used).max()), "일치": bool(np.allclose(lever, used))}) return pd.DataFrame(rows).set_index("시나리오") def project_labor_force(labor: pd.DataFrame, q: list[str], k: int = 20) -> np.ndarray: """(이전 판) 경제활동인구: 관측 구간은 그대로, 이후는 최근 k분기 로그 선형추세로 연장한다. 비교용.""" s = labor["경활"].reindex(q) y = np.log(s.values) obs = np.where(np.isfinite(y))[0] last = obs[-1] kk = min(k, len(obs)) b1, b0 = np.polyfit(np.arange(kk), y[obs[-kk:]], 1) for i in range(last + 1, len(q)): y[i] = b0 + b1 * (kk - 1 + (i - last)) return np.exp(y) def inflation_2000(cpi: pd.DataFrame) -> list[float]: """2000년 전년동기 CPI 상승률(%) — 1999년 지수가 필요하므로 원본 BOKDPM 워크북 지수로 계산(소수 둘째 자리).""" c = cpi.set_index("분기")["CPI지수"] return [round(float(100 * (np.log(c[f"2000Q{k}"]) - np.log(c[f"1999Q{k}"]))), 2) for k in range(1, 5)] def labor_force_base(labor: pd.DataFrame, keis: pd.DataFrame, n: int = 144) -> np.ndarray: """경제활동인구 기준 경로(2000Q1부터 n분기, 천 명) — 관측치 + 인력수급 전망 증가율.""" return project_labor_force_keis(labor, keis, quarters("2000Q1", n)) def build_exog_levels(smoothed: pd.DataFrame, forecast: pd.DataFrame, labor: pd.DataFrame, keis: pd.DataFrame, pop: pd.DataFrame, mac: pd.DataFrame, exg: pd.DataFrame, cpi: pd.DataFrame, use_yus: bool = True, use_pop: bool = True) -> pd.DataFrame: """시나리오별 외생 6종(수준), 2000Q1–2035Q4. build_exog_scenarios_v2.py + 2000년 실질금리 보정(2026-09-02) + 2026-10 정비. 정비: 경제활동인구 = 관측치 + 인력수급 전망 증가율(기준) × exp(인구 레버/100)(시나리오), 세계수요 편차 = 레버 yus(S3 포함). use_yus·use_pop=False는 레버 전달 점검에서 세계수요 편차·인구 레버를 끄는 용도이다. """ bk = {v: np.concatenate([smoothed[v].to_numpy(float), forecast[v].to_numpy(float)]) for v in BK} nH = len(smoothed); N = nH + len(forecast) q = quarters("2000Q1", N) lf = labor_force_base(labor, keis, N) yoy = lambda x: np.array([x[i] - x[i - 4] if i >= 4 else np.nan for i in range(len(x))]) base = {"유가": bk["LOIL_us"].copy(), "세계수요": bk["Y_us"].copy(), "실질금리": bk["CB_kor"] - yoy(bk["LCPI_kor"]), "실질환율": bk["LS_kor"] + bk["LCPI_us"] - bk["LCPI_kor"], "UNR_base": bk["UNR_kor"].copy(), "우리경기": bk["Y_kor"].copy()} col = lambda df, s, c: df[df["시나리오"] == s][c].to_numpy()[:H] RR1, S1, U1, Y1, oil1 = col(mac, "S1", "실질금리RR"), col(mac, "S1", "실질환율갭S"), col(mac, "S1", "실업갭UNR"), col(mac, "S1", "산출갭Y"), col(exg, "S1", "oil_g") f = slice(nH, N) infl00 = inflation_2000(cpi) rr2000 = [round(bk["CB_kor"][k] - infl00[k], 4) for k in range(4)] out = [] for s in SCENARIOS: v = {c: base[c].copy() for c in ("유가", "세계수요", "실질금리", "실질환율", "우리경기")} unr = base["UNR_base"].copy() v["유가"][f] += np.cumsum(col(exg, s, "oil_g") - oil1) # 유가 증가율 편차의 누적 → 수준 편차 if use_yus: v["세계수요"][f] += world_demand_dev(s) v["실질금리"][f] += col(mac, s, "실질금리RR") - RR1 v["실질환율"][f] += col(mac, s, "실질환율갭S") - S1 unr[f] += col(mac, s, "실업갭UNR") - U1 lf_s = lf.copy() if use_pop: lf_s[f] *= np.exp(population_lever(pop, s) / 100.0) v["고용_취업자"] = lf_s * (1 - unr / 100.0) v["우리경기"][f] += col(mac, s, "산출갭Y") - Y1 for i in range(N): out.append([s, q[i]] + [round(float(v[c][i]), 4) for c in EXOG_COLS]) df = pd.DataFrame(out, columns=["시나리오", "분기"] + EXOG_COLS) for k in range(4): # 2000년 실질금리: 1999년 CPI가 필요 df.loc[df["분기"] == f"2000Q{k + 1}", "실질금리"] = rr2000[k] return df def lever_reach(pop: pd.DataFrame, smoothed: pd.DataFrame, forecast: pd.DataFrame, labor: pd.DataFrame, keis: pd.DataFrame, cpi: pd.DataFrame, base: pd.DataFrame) -> pd.DataFrame: """레버 하나씩을 끄고(인구=전부 중위, 충격=0) 외생 6종을 다시 만들어, 2026–2035년 최대 변화폭을 잰다. 0이면 그 레버는 외생 6종에 전달되지 않는다. """ fut = lambda d: d[d["분기"] >= "2026Q1"].set_index(["시나리오", "분기"])[EXOG_COLS] b = fut(base) lf41 = labor_force_base(labor, keis)[-(H + 1):] rows = {} for lever in ("인구", "yus", "oil", "tcost", "dshk"): ex, mac = dsge_paths(pop, lf41, drop=lever) alt = build_exog_levels(smoothed, forecast, labor, keis, pop, mac, ex, cpi, use_yus=(lever != "yus"), use_pop=(lever != "인구")) rows[lever] = (fut(alt) - b).abs().max() out = pd.DataFrame(rows).T out.index.name = "끈 레버" return out def population_transmission(pop: pd.DataFrame, exg: pd.DataFrame, levels: pd.DataFrame, levels_prev: pd.DataFrame | None = None) -> pd.DataFrame: """인구 가정의 크기(설계)와 외생 '취업자 수'의 시나리오 차이(결과)를 2035년 기준으로 대조한다.""" p = pop.set_index("연도") cum = exg.groupby("시나리오")["lab_g"].sum() emp = levels[levels["분기"] == "2035Q4"].set_index("시나리오")["고용_취업자"] rows = [] for s in SCENARIOS: a = LAB_ASSUMP[s] r = {"시나리오": f"{s} {SCN_KO[s]}", "인구 가정": a, "생산연령인구 2035, 중위 대비(%)": 100 * (p.loc[2035, f"{a}_천명"] / p.loc[2035, "중위_천명"] - 1), "인구 레버 2035Q4(100·ln)": float(population_lever(p.reset_index(), s)[-1]), "거시 블록 10년 누적 고용, S1 대비(%p)": cum[s] - cum["S1"], "외생 취업자 2035Q4, S1 대비(%)": 100 * (emp[s] / emp["S1"] - 1)} if levels_prev is not None: ep = levels_prev[levels_prev["분기"] == "2035Q4"].set_index("시나리오")["고용_취업자"] r["이전 판 외생 취업자 2035Q4, S1 대비(%)"] = 100 * (ep[s] / ep["S1"] - 1) rows.append(r) return pd.DataFrame(rows).set_index("시나리오") def calibrated_values(mod_path: str, names: list[str]) -> dict[str, float]: """.mod 파일의 캘리브레이션 값(이름 = 값;)을 읽는다.""" txt = open(mod_path, encoding="utf-8", errors="ignore").read() out = {} for n in names: m = re.search(rf"^\s*{n}\s*=\s*([-\d.eE+]+)\s*;", txt, flags=re.M) out[n] = float(m.group(1)) if m else float("nan") return out def refetch_labor(key_bok: str) -> pd.DataFrame: """ECOS 901Y027 경제활동인구(I61B)·취업자(I61BA), 분기, 천 명.""" return pd.DataFrame({"경활": fetch_ecos(key_bok, "901Y027", "I61B", "Q", "1999Q3"), "취업자": fetch_ecos(key_bok, "901Y027", "I61BA", "Q", "1999Q3")}) def refetch_inflation_2000(key_bok: str) -> list[float]: """ECOS 901Y009 소비자물가지수(월→분기 평균)로 2000년 전년동기 상승률(%, 소수 둘째 자리)을 다시 계산한다.""" c = _m_to_q(fetch_ecos(key_bok, "901Y009", "0", "M", "199901", "200012")) return [round(float(100 * (np.log(c.iloc[4 + k]) - np.log(c.iloc[k]))), 2) for k in range(4)] PARAM_LABEL = { "gamma_kor1": "한국 테일러 준칙 — 금리 평활", "gamma_kor2": "한국 테일러 준칙 — 기대물가 반응", "gamma_kor3": "한국 테일러 준칙 — 산출갭 반응", "beta_kor1": "한국 IS — 산출갭 지속성", "beta_kor2": "한국 IS — 산출갭 기대(선행항)", "beta_kor3": "한국 IS — 실질 회사채금리갭 반응", "lambda_kor1": "한국 필립스 — 기대물가 가중", "lambda_kor2": "한국 필립스 — 산출갭 기울기", "lambda_kor3": "한국 필립스 — 물가목표 가중", "alpha_kor1": "한국 오쿤 — 실업갭 지속성", "alpha_kor2": "한국 오쿤 — 산출갭 계수", "gamma_us1": "미국 테일러 준칙 — 금리 평활", "gamma_us2": "미국 테일러 준칙 — 기대물가 반응", "gamma_us3": "미국 테일러 준칙 — 산출갭 반응", "beta_us1": "미국 IS — 산출갭 지속성", "beta_us2": "미국 IS — 산출갭 기대(선행항)", "lambda_us1": "미국 필립스 — 기대물가 가중", "lambda_us2": "미국 필립스 — 산출갭 기울기", "lambda_us3": "미국 필립스 — 물가목표 가중", "rho_oil1": "유가 상승률 — 정상상태 회귀 가중", } def to_stationary(levels: pd.DataFrame, tf_map: dict | None = None) -> tuple[pd.DataFrame, pd.DataFrame]: """수준 외생 → 산업 VARX 투입형(수준↔수준, 차분↔차분). build_varx_exog.py와 같은 처리, 실질금리는 2026-10 정비로 1차 차분. tf_map: 변환표(기본 EXOG_TF). EXOG_TF_PREV를 주면 이전 판 투입형(실질금리 수준)을 만든다. 반환: (과거 2000Q1–2025Q4 공통=S1, 미래 2026Q1–2035Q4 시나리오별) """ tf_map = tf_map or EXOG_TF def tf(s, kind): s = s.astype(float) out = s if kind == "level" else (s.diff() if kind == "diff" else 100 * np.log(s).diff()) return out.bfill() # 첫 분기 결측은 다음 값으로 채운다(원본 그대로) rows = [] for s, g in levels.groupby("시나리오", sort=False): g = g.reset_index(drop=True) t = {c: tf(g[c], tf_map[c]) for c in EXOG_COLS} for i, qq in enumerate(g["분기"]): rows.append([s, qq] + [round(float(t[c][i]), 4) for c in EXOG_COLS]) out = pd.DataFrame(rows, columns=["시나리오", "분기"] + EXOG_COLS) hist = out[(out["시나리오"] == "S1") & (out["분기"] < "2026Q1")][["분기"] + EXOG_COLS].reset_index(drop=True) scen = out[out["분기"] >= "2026Q1"][["시나리오", "분기"] + EXOG_COLS].reset_index(drop=True) return hist, scen def write_levels_xlsx(levels: pd.DataFrame, path: str) -> None: """수준 외생표를 전달본과 같은 배치로 쓴다: '수준_전체'(과거는 '실적(공통)' 한 번, 미래는 시나리오별) + 외생별 '미래_*' 시트.""" hist = levels[(levels["시나리오"] == "S1") & (levels["분기"] < "2026Q1")].assign(시나리오="실적(공통)") fut = levels[levels["분기"] >= "2026Q1"] with pd.ExcelWriter(path, engine="openpyxl") as xw: pd.concat([hist, fut])[["시나리오", "분기"] + EXOG_COLS].to_excel(xw, sheet_name="수준_전체", index=False) for c in EXOG_COLS: fut.pivot(index="분기", columns="시나리오", values=c)[SCENARIOS].reset_index().to_excel( xw, sheet_name=f"미래_{c}", index=False) def compare_versions(new: pd.DataFrame, prev: pd.DataFrame, keys: list[str], cols: list[str]) -> pd.DataFrame: """새 판 − 이전 판(키로 맞춤): 열별 최대·평균 절대차와 최대차가 생긴 행.""" m = new.merge(prev, on=keys, suffixes=("", "_이전")) rows = [] for c in cols: d = m[c].astype(float) - m[f"{c}_이전"].astype(float) i = d.abs().idxmax() rows.append({"항목": c, "최대 절대차": float(d.abs().max()), "평균 절대차": float(d.abs().mean()), "최대차 위치": " ".join(str(m.loc[i, k]) for k in keys), "새 값": float(m.loc[i, c]), "이전 값": float(m.loc[i, f"{c}_이전"])}) return pd.DataFrame(rows).set_index("항목") def read_levels_xlsx(path: str) -> pd.DataFrame: """공표 수준 외생표(엑셀, '수준_전체' 시트) → 시나리오·분기·6종. 과거 '실적(공통)'은 S1–S5 모두에 붙인다.""" x = pd.read_excel(path, sheet_name="수준_전체") hist = x[x["시나리오"] == "실적(공통)"] fut = x[x["시나리오"] != "실적(공통)"] parts = [pd.concat([hist.assign(시나리오=s), fut[fut["시나리오"] == s]]) for s in SCENARIOS] return pd.concat(parts, ignore_index=True)[["시나리오", "분기"] + EXOG_COLS] def max_abs_diff(a: pd.DataFrame, b: pd.DataFrame, keys: list[str], cols: list[str]) -> pd.Series: """키로 맞춘 두 표의 열별 최대 절대차(행 불일치는 무한대로 표시).""" m = a.merge(b, on=keys, suffixes=("_a", "_b"), how="outer", indicator=True) out = {} for c in cols: d = (m[f"{c}_a"].astype(float) - m[f"{c}_b"].astype(float)).abs() out[c] = float("inf") if (m["_merge"] != "both").any() or d.isna().any() else float(d.max()) return pd.Series(out)