{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "e71d712c",
   "metadata": {},
   "source": [
    "# ② 가구고정효과(within) 분석 + 2020 종량세 서술적 DiD\n",
    "통합 가구파일 `NaSTaBHW.csv`(지속 가구ID `hid` × 17차 wide). 2원 고정효과(가구+연도)로 within 소득-주류 관계를 식별."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "65d269ab",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-07-08T19:07:09.779546Z",
     "iopub.status.busy": "2026-07-08T19:07:09.779546Z",
     "iopub.status.idle": "2026-07-08T19:07:26.702220Z",
     "shell.execute_reply": "2026-07-08T19:07:26.698716Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "가구수: 15559 | 읽은 열: 107\n"
     ]
    }
   ],
   "source": [
    "import pandas as pd, numpy as np, csv, os\n",
    "import statsmodels.formula.api as smf\n",
    "D=r\"D:\\Project\\2026_2\\주류세\\data\\nastab\"; INT=os.path.join(D,\"NaSTaBHW.csv\")\n",
    "MISS={-1,-9,-8,-99,-999,999999}\n",
    "with open(INT,encoding=\"utf-8-sig\") as f: hdr=set(next(csv.reader(f)))\n",
    "need=[\"hid\"]; percols={}\n",
    "for w in range(1,18):\n",
    "    ww=f\"{w:02d}\"\n",
    "    d={\"con\":f\"h{ww}ca002\",\"inc\":f\"h{ww}ca004\",\"alc1\":f\"h{ww}cd002\",\"alc2\":f\"h{ww}cd004\",\n",
    "       \"wt\":f\"h{ww}cwt\",\"wine\":f\"h{ww}cd006\",\"mak\":f\"h{ww}cd008\"}\n",
    "    d={k:v for k,v in d.items() if v in hdr}; percols[w]=d; need+=list(d.values())\n",
    "df=pd.read_csv(INT,usecols=need,encoding=\"utf-8-sig\",low_memory=False)\n",
    "print(\"가구수:\",len(df),\"| 읽은 열:\",len(need))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "3dfa3312",
   "metadata": {},
   "source": [
    "## long 패널 구성 (가구×차수)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "151906df",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-07-08T19:07:26.707218Z",
     "iopub.status.busy": "2026-07-08T19:07:26.707218Z",
     "iopub.status.idle": "2026-07-08T19:07:27.421729Z",
     "shell.execute_reply": "2026-07-08T19:07:27.419729Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "관측 66,307 · 가구 5,706 · 평균 11.6회\n"
     ]
    }
   ],
   "source": [
    "def clean(s):\n",
    "    s=pd.to_numeric(s,errors=\"coerce\"); return s.where(~s.isin(MISS))\n",
    "recs=[]\n",
    "for w in range(1,18):\n",
    "    c=percols[w]\n",
    "    if \"inc\" not in c or \"con\" not in c: continue\n",
    "    sub=pd.DataFrame({\"hid\":df[\"hid\"]}); sub[\"wave\"]=w; sub[\"year\"]=2006+w\n",
    "    sub[\"con\"]=clean(df[c[\"con\"]]).where(lambda x:x>0)\n",
    "    sub[\"inc\"]=clean(df[c[\"inc\"]]).where(lambda x:x>0)\n",
    "    sub[\"alc\"]=clean(df[c[\"alc1\"]]).fillna(0)+clean(df[c[\"alc2\"]]).fillna(0)\n",
    "    sub[\"wt\"]=clean(df[c[\"wt\"]]) if \"wt\" in c else 1.0\n",
    "    sub[\"wine\"]=clean(df[c[\"wine\"]]).fillna(0) if \"wine\" in c else np.nan\n",
    "    sub[\"mak\"]=clean(df[c[\"mak\"]]).fillna(0) if \"mak\" in c else np.nan\n",
    "    recs.append(sub)\n",
    "p=pd.concat(recs,ignore_index=True)\n",
    "p=p[p[\"inc\"].notna()&p[\"con\"].notna()&p[\"wt\"].notna()&(p[\"wt\"]>0)]\n",
    "p[\"share\"]=p[\"alc\"]/p[\"con\"]*100; p[\"linc\"]=np.log(p[\"inc\"]); p[\"lincsq\"]=p[\"linc\"]**2\n",
    "nobs=p.groupby(\"hid\").size(); p=p[p[\"hid\"].isin(nobs[nobs>=2].index)]\n",
    "print(f\"관측 {len(p):,} · 가구 {p['hid'].nunique():,} · 평균 {len(p)/p['hid'].nunique():.1f}회\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "96dcb47c",
   "metadata": {},
   "source": [
    "## 가구 within 디민 후 2원 고정효과 회귀 (가구클러스터 SE)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "2dabd83a",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-07-08T19:07:27.426715Z",
     "iopub.status.busy": "2026-07-08T19:07:27.425718Z",
     "iopub.status.idle": "2026-07-08T19:07:28.523905Z",
     "shell.execute_reply": "2026-07-08T19:07:28.522901Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[선형] log소득 = -0.0902 (p6.2e-05)\n",
      "[2차] log소득=+0.162(p0.24), (log소득)^2=-0.0233(p0.058)\n"
     ]
    }
   ],
   "source": [
    "g=p.groupby(\"hid\")\n",
    "for c in [\"share\",\"linc\",\"lincsq\"]: p[c+\"_w\"]=p[c]-g[c].transform(\"mean\")\n",
    "m1=smf.ols(\"share_w ~ linc_w + C(year)\",data=p).fit(cov_type=\"cluster\",cov_kwds={\"groups\":p[\"hid\"]})\n",
    "m2=smf.ols(\"share_w ~ linc_w + lincsq_w + C(year)\",data=p).fit(cov_type=\"cluster\",cov_kwds={\"groups\":p[\"hid\"]})\n",
    "print(f\"[선형] log소득 = {m1.params['linc_w']:+.4f} (p{m1.pvalues['linc_w']:.1e})\")\n",
    "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})\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "5351268d",
   "metadata": {},
   "source": [
    "**해석:** 같은 가구에서 소득이 오르면 주류 지출비중은 **소폭 감소**(−0.09/log소득, p<0.001). → 횡단면 역U자의 \"저소득 최저\" 구간은 **가구 구성효과**(저소득=고령·비음주)이고, 가구 내 행동은 소득탄력성<1. **중산층 최다부담**은 가구 간 구조적 사실."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2de2247b",
   "metadata": {},
   "source": [
    "## 2020 종량세 서술적 DiD — 막걸리·전통주(처치) vs 와인(대조)\n",
    "2024 기준판매비율은 아직 데이터 밖(최신 2023). 2020 종량세는 코로나 중첩·맥주 분리불가 → **서술적 참고치**."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "4c4fa214",
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-07-08T19:07:28.528541Z",
     "iopub.status.busy": "2026-07-08T19:07:28.527182Z",
     "iopub.status.idle": "2026-07-08T19:07:28.622183Z",
     "shell.execute_reply": "2026-07-08T19:07:28.620675Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "DiD 상호작용(treat:post) = -0.0207  (p 0.167)\n",
      "→ 통계적 유의성 없음. 세제 인과 DiD는 18차(2025조사)·2024 기준판매비율로 재추정 권장.\n"
     ]
    }
   ],
   "source": [
    "dd=p[p[\"wine\"].notna()&p[\"mak\"].notna()].copy(); dd=dd[(dd.year>=2016)&(dd.year<=2023)]\n",
    "dd[\"post\"]=(dd.year>=2020).astype(int)\n",
    "A=dd[[\"hid\",\"post\",\"wt\",\"mak\"]].rename(columns={\"mak\":\"exp\"}); A[\"treat\"]=1\n",
    "B=dd[[\"hid\",\"post\",\"wt\",\"wine\"]].rename(columns={\"wine\":\"exp\"}); B[\"treat\"]=0\n",
    "dl=pd.concat([A,B],ignore_index=True)\n",
    "m3=smf.ols(\"exp ~ treat*post\",data=dl).fit(cov_type=\"cluster\",cov_kwds={\"groups\":dl[\"hid\"]})\n",
    "print(f\"DiD 상호작용(treat:post) = {m3.params['treat:post']:+.4f}  (p {m3.pvalues['treat:post']:.3g})\")\n",
    "print(\"→ 통계적 유의성 없음. 세제 인과 DiD는 18차(2025조사)·2024 기준판매비율로 재추정 권장.\")"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "name": "py3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.12.10"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}