② 가구고정효과(within) 분석 + 2020 종량세 서술적 DiD

통합 가구파일 NaSTaBHW.csv(지속 가구ID hid × 17차 wide). 2원 고정효과(가구+연도)로 within 소득-주류 관계를 식별.

import pandas as pd, numpy as np, csv, os
import statsmodels.formula.api as smf
D=r"D:\Project\2026_2\주류세\data\nastab"; INT=os.path.join(D,"NaSTaBHW.csv")
MISS={-1,-9,-8,-99,-999,999999}
with open(INT,encoding="utf-8-sig") as f: hdr=set(next(csv.reader(f)))
need=["hid"]; percols={}
for w in range(1,18):
    ww=f"{w:02d}"
    d={"con":f"h{ww}ca002","inc":f"h{ww}ca004","alc1":f"h{ww}cd002","alc2":f"h{ww}cd004",
       "wt":f"h{ww}cwt","wine":f"h{ww}cd006","mak":f"h{ww}cd008"}
    d={k:v for k,v in d.items() if v in hdr}; percols[w]=d; need+=list(d.values())
df=pd.read_csv(INT,usecols=need,encoding="utf-8-sig",low_memory=False)
print("가구수:",len(df),"| 읽은 열:",len(need))

가구수: 15559 | 읽은 열: 107

long 패널 구성 (가구×차수)

def clean(s):
    s=pd.to_numeric(s,errors="coerce"); return s.where(~s.isin(MISS))
recs=[]
for w in range(1,18):
    c=percols[w]
    if "inc" not in c or "con" not in c: continue
    sub=pd.DataFrame({"hid":df["hid"]}); sub["wave"]=w; sub["year"]=2006+w
    sub["con"]=clean(df[c["con"]]).where(lambda x:x>0)
    sub["inc"]=clean(df[c["inc"]]).where(lambda x:x>0)
    sub["alc"]=clean(df[c["alc1"]]).fillna(0)+clean(df[c["alc2"]]).fillna(0)
    sub["wt"]=clean(df[c["wt"]]) if "wt" in c else 1.0
    sub["wine"]=clean(df[c["wine"]]).fillna(0) if "wine" in c else np.nan
    sub["mak"]=clean(df[c["mak"]]).fillna(0) if "mak" in c else np.nan
    recs.append(sub)
p=pd.concat(recs,ignore_index=True)
p=p[p["inc"].notna()&p["con"].notna()&p["wt"].notna()&(p["wt"]>0)]
p["share"]=p["alc"]/p["con"]*100; p["linc"]=np.log(p["inc"]); p["lincsq"]=p["linc"]**2
nobs=p.groupby("hid").size(); p=p[p["hid"].isin(nobs[nobs>=2].index)]
print(f"관측 {len(p):,} · 가구 {p['hid'].nunique():,} · 평균 {len(p)/p['hid'].nunique():.1f}회")

관측 66,307 · 가구 5,706 · 평균 11.6회

가구 within 디민 후 2원 고정효과 회귀 (가구클러스터 SE)

g=p.groupby("hid")
for c in ["share","linc","lincsq"]: p[c+"_w"]=p[c]-g[c].transform("mean")
m1=smf.ols("share_w ~ linc_w + C(year)",data=p).fit(cov_type="cluster",cov_kwds={"groups":p["hid"]})
m2=smf.ols("share_w ~ linc_w + lincsq_w + C(year)",data=p).fit(cov_type="cluster",cov_kwds={"groups":p["hid"]})
print(f"[선형] log소득 = {m1.params['linc_w']:+.4f} (p{m1.pvalues['linc_w']:.1e})")
print(f"[2차] log소득={m2.params['linc_w']:+.3f}(p{m2.pvalues['linc_w']:.2g}), (log소득)^2={m2.params['lincsq_w']:+.4f}(p{m2.pvalues['lincsq_w']:.2g})")

[선형] log소득 = -0.0902 (p6.2e-05)
[2차] log소득=+0.162(p0.24), (log소득)^2=-0.0233(p0.058)

해석: 같은 가구에서 소득이 오르면 주류 지출비중은 소폭 감소(−0.09/log소득, p<0.001). → 횡단면 역U자의 "저소득 최저" 구간은 가구 구성효과(저소득=고령·비음주)이고, 가구 내 행동은 소득탄력성<1. 중산층 최다부담은 가구 간 구조적 사실.

2020 종량세 서술적 DiD — 막걸리·전통주(처치) vs 와인(대조)

2024 기준판매비율은 아직 데이터 밖(최신 2023). 2020 종량세는 코로나 중첩·맥주 분리불가 → 서술적 참고치.

dd=p[p["wine"].notna()&p["mak"].notna()].copy(); dd=dd[(dd.year>=2016)&(dd.year<=2023)]
dd["post"]=(dd.year>=2020).astype(int)
A=dd[["hid","post","wt","mak"]].rename(columns={"mak":"exp"}); A["treat"]=1
B=dd[["hid","post","wt","wine"]].rename(columns={"wine":"exp"}); B["treat"]=0
dl=pd.concat([A,B],ignore_index=True)
m3=smf.ols("exp ~ treat*post",data=dl).fit(cov_type="cluster",cov_kwds={"groups":dl["hid"]})
print(f"DiD 상호작용(treat:post) = {m3.params['treat:post']:+.4f}  (p {m3.pvalues['treat:post']:.3g})")
print("→ 통계적 유의성 없음. 세제 인과 DiD는 18차(2025조사)·2024 기준판매비율로 재추정 권장.")

DiD 상호작용(treat:post) = -0.0207  (p 0.167)
→ 통계적 유의성 없음. 세제 인과 DiD는 18차(2025조사)·2024 기준판매비율로 재추정 권장.