{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "c79e977d",
   "metadata": {},
   "source": [
    "\n",
    "File: 04-wooldridge.py\n",
    "\n",
    "\n",
    "Michel Bierlaire\n",
    "\n",
    "Sun Aug 25 18:45:36 2024\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "90261fc5",
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "import biogeme.biogeme_logging as blog\n",
    "import pandas as pd\n",
    "from IPython.core.display_functions import display\n",
    "from biogeme.biogeme import BIOGEME\n",
    "from biogeme.database import Database\n",
    "from biogeme.expressions import (\n",
    "    Beta,\n",
    "    Draws,\n",
    "    Expression,\n",
    "    MonteCarlo,\n",
    "    PanelLikelihoodTrajectory,\n",
    "    Variable,\n",
    "    log,\n",
    ")\n",
    "from biogeme.models import logit\n",
    "from biogeme.results_processing import get_pandas_estimated_parameters\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "cfc2da45",
   "metadata": {},
   "source": [
    "Consider the estimation of the dynamic choice model with panel effects that was performed in the previous laboratory,\n",
    "and reported below for reference. The objective of this laboratory is to address the \"initial condition problem\"\n",
    "using Wooldridge's method."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "db47a0b0",
   "metadata": {},
   "source": [
    "As the estimation time may be long, we ask Biogeme to report the details of the iterations."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "b4275b9a",
   "metadata": {},
   "outputs": [],
   "source": [
    "logger = blog.get_screen_logger(level=blog.INFO)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "4ed3991f",
   "metadata": {},
   "source": [
    "**Tip:**<div class=\"alert alert-block alert-info\">It is advised to start working with a low number of draws, until\n",
    "the script is working well. Then, increase the number of draws to 10000, say. Then, execute the script overnight.\n",
    "</div>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "3624f216",
   "metadata": {},
   "outputs": [],
   "source": [
    "number_of_draws = 10\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a566511a",
   "metadata": {},
   "source": [
    "# Dynamic Choice Models with Panel Effects"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2868b653",
   "metadata": {},
   "source": [
    "We analyze again the smoking behavior of individuals, as a function of their age and the price of tobacco using\n",
    "synthetic data. We develop a model that predicts, for every year, the probability to smoke or not."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ba4076bb",
   "metadata": {},
   "source": [
    "## Data"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "bf93b435",
   "metadata": {},
   "source": [
    "We postulate a true model for the data generation process. It is a mixture of logit models where the utility\n",
    "associated with \"not smoking\" is\n",
    "\\begin{equation}\n",
    "U_{0nt}= \\varepsilon_{0nt}\n",
    "\\end{equation}\n",
    "and the utility associated with \"smoking\" is\n",
    "\\begin{equation}\n",
    "U_{1nt}= \\beta_{nt} y_{n,t-1} + \\beta^p_{nt} P_{t} + c_n + \\varepsilon_{1nt},\n",
    "\\end{equation}\n",
    "where\n",
    "\n",
    "- $\\beta_{nt} = 10$,\n",
    "\n",
    "- $y_{n,t-1}=1$ if $n$ is smoking at time $t-1$, $0$ otherwise,\n",
    "\n",
    "- $\\beta^p_{nt} = -0.1$,\n",
    "\n",
    "- $P_t$ is the price of cigarets at time $t$,\n",
    "\n",
    "- $c_n$ is an individual specific constant that captures the a priori, intrinsic attraction of each individual\n",
    "towards smoking. It is assumed to be normally distributed in the population, with zero mean and standard deviation\n",
    "50: $N(0, 50^2)$,"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "511d4cee",
   "metadata": {},
   "source": [
    "## True value of the parameters"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "24e7235a",
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "true_parameters = pd.DataFrame(\n",
    "    {'Value': [-0.1, 10, 0, 50]},\n",
    "    index=['coef_price', 'beta_last_year', 'cte_mean', 'cte_std'],\n",
    ")\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "15de03c5",
   "metadata": {},
   "source": [
    "## Data"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "75644c37",
   "metadata": {},
   "source": [
    "We observe every individual only from the age of 45 and the age of 55."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "1dd5acd2",
   "metadata": {},
   "outputs": [],
   "source": [
    "df = pd.read_table('smoking55.dat', sep=',')\n",
    "display(df)\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e5019051",
   "metadata": {},
   "source": [
    "Different values for age"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "ae1ef656",
   "metadata": {},
   "outputs": [],
   "source": [
    "print(df['Age'].unique())\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "dc212365",
   "metadata": {},
   "source": [
    "The data contains the following columns:\n",
    "\n",
    "- the age of the individual,\n",
    "- the price of the cigarettes,\n",
    "- a variable that is 1 if the individual is smoking, 0 otherwise,\n",
    "- a variable that is 1 if the individual was smoking last year, 0 otherwise,\n",
    "- a unique id for each individual,\n",
    "- a variable that is 1 if the individual was smoking at the age of 45, in the beginning of the observation period."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "337e63f8",
   "metadata": {},
   "outputs": [],
   "source": [
    "database = Database('smoking55', df)\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "029f280f",
   "metadata": {},
   "source": [
    "Variables"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "05c6dd36",
   "metadata": {},
   "outputs": [],
   "source": [
    "Price = Variable('Price')\n",
    "Smoking = Variable('Smoking')\n",
    "LastYear = Variable('LastYear')\n",
    "Smoking45 = Variable('Smoking45')\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "3b3092ef",
   "metadata": {},
   "source": [
    "We declare that the data set contains panel data."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "244ad5b3",
   "metadata": {},
   "outputs": [],
   "source": [
    "database.panel('Id')\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "8fd60ebc",
   "metadata": {},
   "source": [
    "## Estimation procedure"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "924175e2",
   "metadata": {},
   "source": [
    "The following procedure estimates the choice model, and returns the estimated parameters in a Pandas format. If the\n",
    "model happens to have been already estimated, the estimation results are read from the pickle file and reported."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "446fdfb0",
   "metadata": {},
   "outputs": [],
   "source": [
    "def estimate(\n",
    "    the_logprob: Expression, the_name: str, the_database: Database\n",
    ") -> pd.DataFrame:\n",
    "    \"\"\"Estimates the choice model (or read the estimation results from\n",
    "        file if recycle is True), and returns the estimated parameters\n",
    "        in a Pandas format.\n",
    "\n",
    "    :param the_logprob: formula for the log likelihood function\n",
    "    :param the_name: name of the model (important for the output files)\n",
    "    :param the_database: database\n",
    "    :return: estimated values of the parameters and statistics.\n",
    "    \"\"\"\n",
    "    biogeme = BIOGEME(\n",
    "        the_database,\n",
    "        the_logprob,\n",
    "        number_of_draws=number_of_draws,\n",
    "    )\n",
    "    biogeme.model_name = the_name\n",
    "    results = biogeme.estimate()\n",
    "    print(results.short_summary())\n",
    "    pandas_results = get_pandas_estimated_parameters(estimation_results=results)\n",
    "    return pandas_results.set_index('Name')\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "43e35279",
   "metadata": {},
   "source": [
    "## Dynamic model with serial correlation"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c1fbdb30",
   "metadata": {},
   "source": [
    "In the previous quiz, we have estimated a dynamic model with panel effects to account for serial correlation."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "189cd376",
   "metadata": {},
   "outputs": [],
   "source": [
    "cte_mean = Beta('cte_mean', 0, None, None, 0)\n",
    "cte_std = Beta('cte_std', 1, None, None, 0)\n",
    "cte = cte_mean + cte_std * Draws('agent', 'NORMAL_ANTI')\n",
    "coef_price = Beta('coef_price', 0, None, None, 0)\n",
    "beta_last_year = Beta('beta_last_year', 0, None, None, 0)\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b11da4e6",
   "metadata": {},
   "source": [
    "Model specification"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "6e6fd17a",
   "metadata": {},
   "outputs": [],
   "source": [
    "v_smoking = beta_last_year * LastYear + coef_price * Price + cte\n",
    "v_non_smoking = 0\n",
    "v = {0: v_non_smoking, 1: v_smoking}\n",
    "conditional_observation_probability = logit(v, None, Smoking)\n",
    "conditional_individual_probability = PanelLikelihoodTrajectory(\n",
    "    conditional_observation_probability\n",
    ")\n",
    "log_probability = log(MonteCarlo(conditional_individual_probability))\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ccf1e7e4",
   "metadata": {},
   "source": [
    "Estimation"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "ef2c3a39",
   "metadata": {},
   "outputs": [],
   "source": [
    "r_serial_dynamic = estimate(log_probability, 'dynamic_model_serial', database)\n",
    "display(r_serial_dynamic)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "5d5484df",
   "metadata": {},
   "source": [
    "### Comparison of the estimates"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "f6c3719e",
   "metadata": {},
   "outputs": [],
   "source": [
    "summary = pd.concat(\n",
    "    [\n",
    "        true_parameters['Value'],\n",
    "        r_serial_dynamic['Value'],\n",
    "    ],\n",
    "    axis=1,\n",
    ")\n",
    "summary.columns = ['True', 'Dynamic + serial']\n",
    "summary.fillna('')\n",
    "display(summary)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f1e7d93f",
   "metadata": {},
   "source": [
    "We observe here the issue of the \"initial condition problem\". Although the model specification is correct (it is the\n",
    "same model as the data generation process), the values of the parameters are not correctly recovered. It is because\n",
    "the first observed choice, that is, the fact that an individual is smoking at the age of 45, is strongly correlated\n",
    "with the agent effect. This creates endogeneity. One visible consequence is the positive price coefficient."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e19148bb",
   "metadata": {},
   "source": [
    "Estimate the parameters using Wooldridge's method to address the endogeneity issue."
   ]
  }
 ],
 "metadata": {},
 "nbformat": 4,
 "nbformat_minor": 5
}