# -*- coding: utf-8 -*-
"""한치록·김덕파(2018) 맥주 수요체계(AIDS) 축소재현 — 국산 vs 수입 2재화, 종가세기 2006~2019
 자료: 국세청 주세신고(맥주) 국산/수입 출고량·과세표준(종가세=신고 출고금액)
 방법: LA-AIDS 지출몫 방정식(Stone 물가지수), 동차성 부과 → Marshallian 자체/교차/지출 탄력성
 주의: 연 14obs·집계 → 정밀도 낮음. 논문의 제품·용기별(수입프리미엄/국산코어…) 큰 탄력성은 스캐너 필요.
"""
import sys,io,json; sys.stdout=io.TextIOWrapper(sys.stdout.buffer,encoding='utf-8')
import numpy as np, pandas as pd, statsmodels.formula.api as smf
B='C:/Temp/주류세/공개데이터'
g=pd.read_csv(B+'/국세청_국세통계/주세신고_주종별_제조수입_2006_2025.csv')
g['값']=pd.to_numeric(g['값'].astype(str).str.replace(',',''),errors='coerce')
def piv(ind):
    p=g[(g.주종=='맥주')&(g.지표==ind)].pivot_table(index='연도',columns='구분',values='값')
    return p
vol=piv('출고량'); tax=piv('과세표준_종가세')   # tax=과세표준(≈신고 출고금액, 백만원 단위 가정)
d=pd.DataFrame(index=vol.index)
d['q_imp']=vol['수입분']; d['q_dom']=vol['합계']-vol['수입분']
d['x_imp']=tax['수입분']; d['x_dom']=tax['합계']-tax['수입분']     # 매출(지출)
d=d.dropna()
d=d[(d.index>=2006)&(d.index<=2019)]
d['p_imp']=d.x_imp/d.q_imp; d['p_dom']=d.x_dom/d.q_dom            # 단가=매출/물량
d['X']=d.x_imp+d.x_dom
d['w_imp']=d.x_imp/d.X; d['w_dom']=d.x_dom/d.X
# CPI 맥주로 실질화(선택) — 지출 실질화용
cpi=pd.read_csv(B+'/한국은행_ECOS/CPI_품목별_소주맥주_200601_202512.csv')
cb=cpi[cpi.품목=='맥주'].copy(); cb['yr']=cb.시점//100; cyr=cb.groupby('yr').값.mean()
d['cpi']=d.index.map(cyr)
# Stone 물가지수 lnP* = Σ w_i ln p_i
d['lnP']=d.w_imp*np.log(d.p_imp)+d.w_dom*np.log(d.p_dom)
d['ln_XP']=np.log(d.X/d.cpi)                     # 실질 총지출(맥주)
d['ln_rel']=np.log(d.p_imp/d.p_dom)              # 상대가격(동차성 부과)
d['t']=d.index-2006
print('===== 기초 (맥주 국산/수입, 2006~2019) =====')
print(d[['q_imp','q_dom','p_imp','p_dom','w_imp']].round(3).to_string())
print(f'\n수입 지출몫 w_imp: {d.w_imp.iloc[0]:.3f}(2006) → {d.w_imp.iloc[-1]:.3f}(2019)')
print(f'상대가격 p_imp/p_dom: {(d.p_imp/d.p_dom).iloc[0]:.2f} → {(d.p_imp/d.p_dom).iloc[-1]:.2f}')

def aids(df, use_trend):
    f='w_imp ~ ln_rel + ln_XP' + (' + t' if use_trend else '')
    m=smf.ols(f,data=df).fit(cov_type='HC1')
    g_ii=m.params['ln_rel']          # γ_imp,imp (동차성: γ_imp,dom=-γ)
    b_i =m.params['ln_XP']           # β_imp
    wi=df.w_imp.mean(); wd=df.w_dom.mean()
    # Marshallian (LA-AIDS): e_ij=-δ_ij+(γ_ij-β_i w_j)/w_i
    e_ii=-1+(g_ii - b_i*wi)/wi
    e_id= 0+(-g_ii - b_i*wd)/wi
    # 국산(adding-up): β_dom=-β_imp, γ_dom,dom=γ, γ_dom,imp=-γ
    e_dd=-1+(g_ii - (-b_i)*wd)/wd
    e_di= 0+(-g_ii - (-b_i)*wi)/wd
    eta_i=1+b_i/wi; eta_d=1+(-b_i)/wd
    return m,dict(gamma=g_ii,beta=b_i,e_imp_imp=e_ii,e_imp_dom=e_id,e_dom_dom=e_dd,e_dom_imp=e_di,eta_imp=eta_i,eta_dom=eta_d,wi=wi,wd=wd)

for lab,tr in [('추세 미통제',False),('추세(연도) 통제',True)]:
    m,r=aids(d,tr)
    print(f'\n===== AIDS [{lab}] (R²={m.rsquared:.2f}) =====')
    print(f'  γ(ln상대가)={r["gamma"]:+.3f}(t={m.tvalues["ln_rel"]:.2f})  β(ln실질지출)={r["beta"]:+.3f}(t={m.tvalues["ln_XP"]:.2f})')
    print(f'  [자체가격탄력성] 수입 e={r["e_imp_imp"]:+.2f} | 국산 e={r["e_dom_dom"]:+.2f}')
    print(f'  [교차가격탄력성] 수입수요의 국산가격탄력성(국산가↑→수입) {r["e_imp_dom"]:+.2f} · 국산수요의 수입가격탄력성(수입가↑→국산) {r["e_dom_imp"]:+.2f}')
    print(f'  [지출탄력성] 수입 {r["eta_imp"]:+.2f} | 국산 {r["eta_dom"]:+.2f}   (평균몫 수입 {r["wi"]:.2f}/국산 {r["wd"]:.2f})')
    if tr: OUT=r; OUT['R2']=m.rsquared; OUT['t_gamma']=float(m.tvalues['ln_rel']); OUT['t_beta']=float(m.tvalues['ln_XP'])

print('\n[원논문(한치록·김덕파 2018) 대조] 자체가격탄력성: 수입프리미엄 -4.8~-6.0, 국산코어 -1.8~-1.9')
print(' → 본 축소재현은 수입<->국산 2군 집계·연간이라 절댓값은 작으나, 수입이 국산보다 더 탄력적이면 정성 일치')
json.dump({k:(round(v,3) if isinstance(v,float) else v) for k,v in OUT.items()},
          open(B+'/문헌/han2018_aids_2good.json','w',encoding='utf-8'),ensure_ascii=False,indent=2)
d.round(4).to_csv(B+'/문헌/han2018_aids_data.csv',encoding='utf-8-sig')
print('\n저장: 문헌/han2018_aids_2good.json, han2018_aids_data.csv')