ipynb/chris_conlon-blp.ipynb

{
    "cells": [
        {
            "cell_type": "markdown",
            "metadata": {},
            "source": "# Random Coefficients Logit Tutorial with the Automobile Data"
        },
        {
            "cell_type": "code",
            "execution_count": 1,
            "metadata": {},
            "outputs": [
                {
                    "data": {
                        "text/plain": [
                            "'0.7.0'"
                        ]
                    },
                    "execution_count": 1,
                    "metadata": {},
                    "output_type": "execute_result"
                }
            ],
            "source": [
                "import pyblp\n",
                "import numpy as np\n",
                "import pandas as pd \n",
                "\n",
                "pyblp.options.digits = 2\n",
                "pyblp.options.verbose = False\n",
                "pyblp.__version__"
            ]
        },
        {
            "cell_type": "markdown",
            "metadata": {},
            "source": "In this tutorial, we'll use data from [Berry, Levinsohn, and Pakes (1995)](https://pyblp.readthedocs.io/en/stable/references.html#berry-levinsohn-and-pakes-1995) to solve the paper's automobile problem.\n\n\n## Application of Random Coefficients Logit with the Automobile Data\n\nThis tutorial is similar to the [fake cereal tutorial](nevo.ipynb), but exhibits some other features of pyblp:\n\n- Incorporating a supply side into demand estimation.\n- Allowing for simple price-income demographic effects.\n- Calculating clustered standard errors.\n\n\n### Loading the Data\n\nWe'll use [pandas](https://pandas.pydata.org/) to load two sets of data:\n\n1. `product_data`, which contains prices, shares, and other product characteristics.\n2. `agent_data`, which contains draws from the distribution of heterogeneity."
        },
        {
            "cell_type": "code",
            "execution_count": 1,
            "metadata": {},
            "outputs": [
                {
                    "data": {
                        "text/html": [
                            "<div>\n",
                            "<style scoped>\n",
                            "    .dataframe tbody tr th:only-of-type {\n",
                            "        vertical-align: middle;\n",
                            "    }\n",
                            "\n",
                            "    .dataframe tbody tr th {\n",
                            "        vertical-align: top;\n",
                            "    }\n",
                            "\n",
                            "    .dataframe thead th {\n",
                            "        text-align: right;\n",
                            "    }\n",
                            "</style>\n",
                            "<table border=\"1\" class=\"dataframe\">\n",
                            "  <thead>\n",
                            "    <tr style=\"text-align: right;\">\n",
                            "      <th></th>\n",
                            "      <th>market_ids</th>\n",
                            "      <th>clustering_ids</th>\n",
                            "      <th>car_ids</th>\n",
                            "      <th>firm_ids</th>\n",
                            "      <th>region</th>\n",
                            "      <th>shares</th>\n",
                            "      <th>prices</th>\n",
                            "      <th>hpwt</th>\n",
                            "      <th>air</th>\n",
                            "      <th>mpd</th>\n",
                            "      <th>...</th>\n",
                            "      <th>demand_instruments2</th>\n",
                            "      <th>demand_instruments3</th>\n",
                            "      <th>demand_instruments4</th>\n",
                            "      <th>demand_instruments5</th>\n",
                            "      <th>supply_instruments0</th>\n",
                            "      <th>supply_instruments1</th>\n",
                            "      <th>supply_instruments2</th>\n",
                            "      <th>supply_instruments3</th>\n",
                            "      <th>supply_instruments4</th>\n",
                            "      <th>supply_instruments5</th>\n",
                            "    </tr>\n",
                            "  </thead>\n",
                            "  <tbody>\n",
                            "    <tr>\n",
                            "      <th>0</th>\n",
                            "      <td>1971</td>\n",
                            "      <td>AMGREM71</td>\n",
                            "      <td>129</td>\n",
                            "      <td>15</td>\n",
                            "      <td>US</td>\n",
                            "      <td>0.001051</td>\n",
                            "      <td>4.935802</td>\n",
                            "      <td>0.528997</td>\n",
                            "      <td>0</td>\n",
                            "      <td>1.888146</td>\n",
                            "      <td>...</td>\n",
                            "      <td>0.566217</td>\n",
                            "      <td>0.365328</td>\n",
                            "      <td>0.659480</td>\n",
                            "      <td>0.141017</td>\n",
                            "      <td>-0.011161</td>\n",
                            "      <td>1.478879</td>\n",
                            "      <td>-0.546875</td>\n",
                            "      <td>-0.163302</td>\n",
                            "      <td>-0.833091</td>\n",
                            "      <td>0.301411</td>\n",
                            "    </tr>\n",
                            "    <tr>\n",
                            "      <th>1</th>\n",
                            "      <td>1971</td>\n",
                            "      <td>AMHORN71</td>\n",
                            "      <td>130</td>\n",
                            "      <td>15</td>\n",
                            "      <td>US</td>\n",
                            "      <td>0.000670</td>\n",
                            "      <td>5.516049</td>\n",
                            "      <td>0.494324</td>\n",
                            "      <td>0</td>\n",
                            "      <td>1.935989</td>\n",
                            "      <td>...</td>\n",
                            "      <td>0.566217</td>\n",
                            "      <td>0.290959</td>\n",
                            "      <td>0.173552</td>\n",
                            "      <td>0.128205</td>\n",
                            "      <td>-0.079317</td>\n",
                            "      <td>1.088327</td>\n",
                            "      <td>-0.546875</td>\n",
                            "      <td>-0.095609</td>\n",
                            "      <td>-0.390314</td>\n",
                            "      <td>0.289947</td>\n",
                            "    </tr>\n",
                            "    <tr>\n",
                            "      <th>2</th>\n",
                            "      <td>1971</td>\n",
                            "      <td>AMJAVL71</td>\n",
                            "      <td>132</td>\n",
                            "      <td>15</td>\n",
                            "      <td>US</td>\n",
                            "      <td>0.000341</td>\n",
                            "      <td>7.108642</td>\n",
                            "      <td>0.467613</td>\n",
                            "      <td>0</td>\n",
                            "      <td>1.716799</td>\n",
                            "      <td>...</td>\n",
                            "      <td>0.566217</td>\n",
                            "      <td>0.599771</td>\n",
                            "      <td>-0.546387</td>\n",
                            "      <td>0.002634</td>\n",
                            "      <td>0.021034</td>\n",
                            "      <td>0.609213</td>\n",
                            "      <td>-0.546875</td>\n",
                            "      <td>-0.449818</td>\n",
                            "      <td>0.400461</td>\n",
                            "      <td>0.434632</td>\n",
                            "    </tr>\n",
                            "    <tr>\n",
                            "      <th>3</th>\n",
                            "      <td>1971</td>\n",
                            "      <td>AMMATA71</td>\n",
                            "      <td>134</td>\n",
                            "      <td>15</td>\n",
                            "      <td>US</td>\n",
                            "      <td>0.000522</td>\n",
                            "      <td>6.839506</td>\n",
                            "      <td>0.426540</td>\n",
                            "      <td>0</td>\n",
                            "      <td>1.687871</td>\n",
                            "      <td>...</td>\n",
                            "      <td>0.566217</td>\n",
                            "      <td>0.620544</td>\n",
                            "      <td>-1.122968</td>\n",
                            "      <td>0.089023</td>\n",
                            "      <td>-0.090014</td>\n",
                            "      <td>0.207461</td>\n",
                            "      <td>-0.546875</td>\n",
                            "      <td>-0.454159</td>\n",
                            "      <td>0.934641</td>\n",
                            "      <td>0.331099</td>\n",
                            "    </tr>\n",
                            "    <tr>\n",
                            "      <th>4</th>\n",
                            "      <td>1971</td>\n",
                            "      <td>AMAMBS71</td>\n",
                            "      <td>136</td>\n",
                            "      <td>15</td>\n",
                            "      <td>US</td>\n",
                            "      <td>0.000442</td>\n",
                            "      <td>8.928395</td>\n",
                            "      <td>0.452489</td>\n",
                            "      <td>0</td>\n",
                            "      <td>1.504286</td>\n",
                            "      <td>...</td>\n",
                            "      <td>0.566217</td>\n",
                            "      <td>0.877198</td>\n",
                            "      <td>-1.258575</td>\n",
                            "      <td>-0.153840</td>\n",
                            "      <td>0.038013</td>\n",
                            "      <td>0.385211</td>\n",
                            "      <td>-0.546875</td>\n",
                            "      <td>-0.728959</td>\n",
                            "      <td>1.146654</td>\n",
                            "      <td>0.520555</td>\n",
                            "    </tr>\n",
                            "  </tbody>\n",
                            "</table>\n",
                            "<p>5 rows × 25 columns</p>\n",
                            "</div>"
                        ],
                        "text/plain": [
                            "   market_ids clustering_ids  car_ids  firm_ids region    shares    prices  \\\n",
                            "0        1971       AMGREM71      129        15     US  0.001051  4.935802   \n",
                            "1        1971       AMHORN71      130        15     US  0.000670  5.516049   \n",
                            "2        1971       AMJAVL71      132        15     US  0.000341  7.108642   \n",
                            "3        1971       AMMATA71      134        15     US  0.000522  6.839506   \n",
                            "4        1971       AMAMBS71      136        15     US  0.000442  8.928395   \n",
                            "\n",
                            "       hpwt  air       mpd         ...           demand_instruments2  \\\n",
                            "0  0.528997    0  1.888146         ...                      0.566217   \n",
                            "1  0.494324    0  1.935989         ...                      0.566217   \n",
                            "2  0.467613    0  1.716799         ...                      0.566217   \n",
                            "3  0.426540    0  1.687871         ...                      0.566217   \n",
                            "4  0.452489    0  1.504286         ...                      0.566217   \n",
                            "\n",
                            "   demand_instruments3  demand_instruments4  demand_instruments5  \\\n",
                            "0             0.365328             0.659480             0.141017   \n",
                            "1             0.290959             0.173552             0.128205   \n",
                            "2             0.599771            -0.546387             0.002634   \n",
                            "3             0.620544            -1.122968             0.089023   \n",
                            "4             0.877198            -1.258575            -0.153840   \n",
                            "\n",
                            "   supply_instruments0  supply_instruments1  supply_instruments2  \\\n",
                            "0            -0.011161             1.478879            -0.546875   \n",
                            "1            -0.079317             1.088327            -0.546875   \n",
                            "2             0.021034             0.609213            -0.546875   \n",
                            "3            -0.090014             0.207461            -0.546875   \n",
                            "4             0.038013             0.385211            -0.546875   \n",
                            "\n",
                            "   supply_instruments3  supply_instruments4  supply_instruments5  \n",
                            "0            -0.163302            -0.833091             0.301411  \n",
                            "1            -0.095609            -0.390314             0.289947  \n",
                            "2            -0.449818             0.400461             0.434632  \n",
                            "3            -0.454159             0.934641             0.331099  \n",
                            "4            -0.728959             1.146654             0.520555  \n",
                            "\n",
                            "[5 rows x 25 columns]"
                        ]
                    },
                    "execution_count": 1,
                    "metadata": {},
                    "output_type": "execute_result"
                }
            ],
            "source": [
                "product_data = pd.read_csv(pyblp.data.BLP_PRODUCTS_LOCATION)\n",
                "product_data.head()"
            ]
        },
        {
            "cell_type": "markdown",
            "metadata": {},
            "source": "The `product_data` contains market IDs, product IDs, firm IDs, shares, prices, a number of product characteristics, and some pre-computed excluded instruments. The product IDs are called `clustering_ids` because they will be used to compute clustered standard errors. For more information about the instruments and the example data as a whole, refer to the [`data`](https://pyblp.readthedocs.io/en/stable/_api/pyblp.data.html#module-pyblp.data) module.\n\nThe `agent_data` contains market IDs, integration weights $w_{ti}$, integration nodes $\\nu_{ti}$, and demographics $d_{ti}$. Here we use $I_{t} = 200$ equally weighted Monte Carlo draws per market.\n\nIn non-example problems, it is usually a better idea to use many more draws, or a more sophisticated [`Integration`](https://pyblp.readthedocs.io/en/stable/_api/pyblp.Integration.html#pyblp.Integration) configuration such as sparse grid quadrature."
        },
        {
            "cell_type": "code",
            "execution_count": 1,
            "metadata": {},
            "outputs": [
                {
                    "data": {
                        "text/html": [
                            "<div>\n",
                            "<style scoped>\n",
                            "    .dataframe tbody tr th:only-of-type {\n",
                            "        vertical-align: middle;\n",
                            "    }\n",
                            "\n",
                            "    .dataframe tbody tr th {\n",
                            "        vertical-align: top;\n",
                            "    }\n",
                            "\n",
                            "    .dataframe thead th {\n",
                            "        text-align: right;\n",
                            "    }\n",
                            "</style>\n",
                            "<table border=\"1\" class=\"dataframe\">\n",
                            "  <thead>\n",
                            "    <tr style=\"text-align: right;\">\n",
                            "      <th></th>\n",
                            "      <th>market_ids</th>\n",
                            "      <th>weights</th>\n",
                            "      <th>nodes0</th>\n",
                            "      <th>nodes1</th>\n",
                            "      <th>nodes2</th>\n",
                            "      <th>nodes3</th>\n",
                            "      <th>nodes4</th>\n",
                            "      <th>income</th>\n",
                            "    </tr>\n",
                            "  </thead>\n",
                            "  <tbody>\n",
                            "    <tr>\n",
                            "      <th>0</th>\n",
                            "      <td>1971</td>\n",
                            "      <td>0.005</td>\n",
                            "      <td>0.548814</td>\n",
                            "      <td>0.457760</td>\n",
                            "      <td>0.564690</td>\n",
                            "      <td>0.395537</td>\n",
                            "      <td>0.392173</td>\n",
                            "      <td>9.728478</td>\n",
                            "    </tr>\n",
                            "    <tr>\n",
                            "      <th>1</th>\n",
                            "      <td>1971</td>\n",
                            "      <td>0.005</td>\n",
                            "      <td>0.715189</td>\n",
                            "      <td>0.376918</td>\n",
                            "      <td>0.839746</td>\n",
                            "      <td>0.844017</td>\n",
                            "      <td>0.041157</td>\n",
                            "      <td>7.908957</td>\n",
                            "    </tr>\n",
                            "    <tr>\n",
                            "      <th>2</th>\n",
                            "      <td>1971</td>\n",
                            "      <td>0.005</td>\n",
                            "      <td>0.602763</td>\n",
                            "      <td>0.702335</td>\n",
                            "      <td>0.376884</td>\n",
                            "      <td>0.150442</td>\n",
                            "      <td>0.923301</td>\n",
                            "      <td>11.079404</td>\n",
                            "    </tr>\n",
                            "    <tr>\n",
                            "      <th>3</th>\n",
                            "      <td>1971</td>\n",
                            "      <td>0.005</td>\n",
                            "      <td>0.544883</td>\n",
                            "      <td>0.207324</td>\n",
                            "      <td>0.499676</td>\n",
                            "      <td>0.306309</td>\n",
                            "      <td>0.406235</td>\n",
                            "      <td>17.641671</td>\n",
                            "    </tr>\n",
                            "    <tr>\n",
                            "      <th>4</th>\n",
                            "      <td>1971</td>\n",
                            "      <td>0.005</td>\n",
                            "      <td>0.423655</td>\n",
                            "      <td>0.074280</td>\n",
                            "      <td>0.081302</td>\n",
                            "      <td>0.094570</td>\n",
                            "      <td>0.944282</td>\n",
                            "      <td>12.423995</td>\n",
                            "    </tr>\n",
                            "  </tbody>\n",
                            "</table>\n",
                            "</div>"
                        ],
                        "text/plain": [
                            "   market_ids  weights    nodes0    nodes1    nodes2    nodes3    nodes4  \\\n",
                            "0        1971    0.005  0.548814  0.457760  0.564690  0.395537  0.392173   \n",
                            "1        1971    0.005  0.715189  0.376918  0.839746  0.844017  0.041157   \n",
                            "2        1971    0.005  0.602763  0.702335  0.376884  0.150442  0.923301   \n",
                            "3        1971    0.005  0.544883  0.207324  0.499676  0.306309  0.406235   \n",
                            "4        1971    0.005  0.423655  0.074280  0.081302  0.094570  0.944282   \n",
                            "\n",
                            "      income  \n",
                            "0   9.728478  \n",
                            "1   7.908957  \n",
                            "2  11.079404  \n",
                            "3  17.641671  \n",
                            "4  12.423995  "
                        ]
                    },
                    "execution_count": 1,
                    "metadata": {},
                    "output_type": "execute_result"
                }
            ],
            "source": [
                "agent_data = pd.read_csv(pyblp.data.BLP_AGENTS_LOCATION)\n",
                "agent_data.head()"
            ]
        },
        {
            "cell_type": "markdown",
            "metadata": {},
            "source": "### Setting up the Problem\n\nUnlike the fake cereal problem, we won't absorb any fixed effects in the automobile problem, so the linear part of demand $X_1$ has more components. We also need to specify a formula for the random coefficients $X_2$, including a random coefficient on the constant, which captures correlation among all inside goods.\n\nThe primary new addition to the model relative to the fake cereal problem is that we add a supply side formula for product characteristics that contribute to marginal costs, $X_3$. The [patsy](https://patsy.readthedocs.io/en/stable/)-style formulas support functions of regressors such as the `log` function used below.\n\nWe stack the three product formulations in order: $X_1$, $X_2$, and $X_3$."
        },
        {
            "cell_type": "code",
            "execution_count": 1,
            "metadata": {},
            "outputs": [
                {
                    "data": {
                        "text/plain": [
                            "(1 + hpwt + air + mpd + space,\n",
                            " 1 + prices + hpwt + air + mpd + space,\n",
                            " 1 + log(hpwt) + air + log(mpg) + log(space) + trend)"
                        ]
                    },
                    "execution_count": 1,
                    "metadata": {},
                    "output_type": "execute_result"
                }
            ],
            "source": [
                "product_formulations = (\n",
                "   pyblp.Formulation('1 + hpwt + air + mpd + space'),\n",
                "   pyblp.Formulation('1 + prices + hpwt + air + mpd + space'),\n",
                "   pyblp.Formulation('1 + log(hpwt) + air + log(mpg) + log(space) + trend')\n",
                ")\n",
                "product_formulations"
            ]
        },
        {
            "cell_type": "markdown",
            "metadata": {},
            "source": "The original specification for the automobile problem includes the term $\\log(y_i - p_j)$, in which $y$ is income and $p$ are prices. Instead of including this term, which gives rise to a host of numerical problems, we'll follow [Berry, Levinsohn, and Pakes (1999)](https://pyblp.readthedocs.io/en/stable/references.html#berry-levinsohn-and-pakes-1999) and use its first-order linear approximation, $p_j / y_i$. \n\nThe agent formulation for demographics, $d$, includes a column of $1 / y_i$ values, which we'll interact with $p_j$. To do this, we will treat draws of $y_i$ as demographic variables."
        },
        {
            "cell_type": "code",
            "execution_count": 1,
            "metadata": {},
            "outputs": [
                {
                    "data": {
                        "text/plain": [
                            "I(1 / income)"
                        ]
                    },
                    "execution_count": 1,
                    "metadata": {},
                    "output_type": "execute_result"
                }
            ],
            "source": [
                "agent_formulation = pyblp.Formulation('0 + I(1 / income)')\n",
                "agent_formulation"
            ]
        },
        {
            "cell_type": "markdown",
            "metadata": {},
            "source": "As in the cereal example, the [`Problem`](https://pyblp.readthedocs.io/en/stable/_api/pyblp.Problem.html#pyblp.Problem) can be constructed by combining the `product_formulations`, `product_data`, `agent_formulation`, and `agent_data`. "
        },
        {
            "cell_type": "code",
            "execution_count": 1,
            "metadata": {},
            "outputs": [
                {
                    "data": {
                        "text/plain": [
                            "Dimensions:\n",
                            "=======================================================\n",
                            " T    N     F    I     K1    K2    K3    D    MD    MS \n",
                            "---  ----  ---  ----  ----  ----  ----  ---  ----  ----\n",
                            "20   2217  26   4000   5     6     6     1    11    12 \n",
                            "=======================================================\n",
                            "\n",
                            "Formulations:\n",
                            "======================================================================================\n",
                            "       Column Indices:            0          1        2       3          4         5  \n",
                            "-----------------------------  --------  ---------  -----  --------  ----------  -----\n",
                            " X1: Linear Characteristics       1        hpwt      air     mpd       space          \n",
                            "X2: Nonlinear Characteristics     1       prices    hpwt     air        mpd      space\n",
                            "  X3: Cost Characteristics        1      log(hpwt)   air   log(mpg)  log(space)  trend\n",
                            "       d: Demographics         1/income                                               \n",
                            "======================================================================================"
                        ]
                    },
                    "execution_count": 1,
                    "metadata": {},
                    "output_type": "execute_result"
                }
            ],
            "source": [
                "problem = pyblp.Problem(product_formulations, product_data, agent_formulation, agent_data)\n",
                "problem"
            ]
        },
        {
            "cell_type": "markdown",
            "metadata": {},
            "source": "The problem outputs a table of dimensions:\n\n- $T$ denotes the number of markets.\n- $N$ is the length of the dataset (the number of products across all markets).\n- $F$ denotes the number of firms.\n- $I = \\sum_t I_t$ is the total number of agents across all markets (200 draws per market times 20 markets).\n- $K_1$ is the number of linear demand characteristics.\n- $K_2$ is the number of nonlinear demand characteristics.\n- $K_3$ is the number of linear supply characteristics.\n- $D$ is the number of demographic variables.\n- $M_D$ is the number of demand instruments, including exogenous regressors.\n- $M_S$ is the number of supply instruments, including exogenous regressors.\n\nThe formulations table describes all four formulas for linear characteristics, nonlinear characteristics, cost characteristics, and demographics."
        },
        {
            "cell_type": "markdown",
            "metadata": {},
            "source": "### Solving the Problem"
        },
        {
            "cell_type": "markdown",
            "metadata": {},
            "source": "The only remaining decisions are:\n\n- Choosing $\\Sigma$ and $\\Pi$ starting values, $\\Sigma_0$ and $\\Pi_0$.\n- Potentially choosing bounds for $\\Sigma$ and $\\Pi$.\n- Choosing a form for marginal costs, $c_{jt}$: either a linear or log-linear functional form.\n\nThe decisions we will use are:\n\n- Use published estimates as our starting values in $\\Sigma_0$.\n- Interact the inverse of income, $1 / y_i$, only with prices, and use the published estimate on $\\log(y_i - p_j)$ as our starting value for $\\alpha$ in $\\Pi_0$.\n- Bound $\\Sigma_0$ to be positive since it is a diagonal matrix where the diagonal consists of standard deviations.\n- Constrain the $p_j / y_i$ coefficient to be negative. Specifically, we'll use a bound that's slightly smaller than zero because when the parameter is exactly zero, there are matrix inversion problems with estimating marginal costs.\n\nWhen using a routine that supports bounds, [`Problem`](https://pyblp.readthedocs.io/en/stable/_api/pyblp.Problem.html#pyblp.Problem) chooses some default bounds. These bounds are often very wide, so it's usually a good idea to set your own more restrictive bounds."
        },
        {
            "cell_type": "code",
            "execution_count": 1,
            "metadata": {},
            "outputs": [],
            "source": [
                "initial_sigma = np.diag([3.612, 0, 4.628, 1.818, 1.050, 2.056])\n",
                "initial_pi = np.c_[[0, -43.501, 0, 0, 0, 0]]\n",
                "sigma_bounds = (\n",
                "   np.zeros_like(initial_sigma),\n",
                "   np.diag([100, 0, 100, 100, 100, 100])\n",
                ")\n",
                "pi_bounds = (\n",
                "   np.c_[[0, -100, 0, 0, 0, 0]],\n",
                "   np.c_[[0, -0.1, 0, 0, 0, 0]]\n",
                ")"
            ]
        },
        {
            "cell_type": "markdown",
            "metadata": {},
            "source": "Note that there are only 5 nonzeros on the diagonal of $\\Sigma$, which means that we only need 5 columns of integration nodes to integrate over these 5 dimensions of unobserved heterogeneity. Indeed, `agent_data` contains exactly 5 columns of nodes. If we were to ignore the $\\log(y_i - p_j)$ term (by not configuring $\\Pi$) and include a term on prices in $\\Sigma$ instead, we would have needed 6 columns of integration nodes in our `agent_data`.\n\nA linear marginal cost specification is the default setting, so we'll need to use the `costs_type` argument of [`Problem.solve`](https://pyblp.readthedocs.io/en/stable/_api/pyblp.Problem.solve.html#pyblp.Problem.solve) to employ the log-linear specification used by [Berry, Levinsohn, and Pakes (1995)](https://pyblp.readthedocs.io/en/stable/references.html#berry-levinsohn-and-pakes-1995). A downside of this specification is that nonpositive estimated marginal costs can create problems for the optimization routine when computing $\\log c(\\hat{\\theta})$. We'll use the `costs_bounds` argument to bound marginal costs from below by a small number. \n\nFinally, as in the original paper, we'll use `W_type` and `se_type` to cluster by product IDs, which were specified as `clustering_ids` in `product_data`."
        },
        {
            "cell_type": "code",
            "execution_count": 1,
            "metadata": {},
            "outputs": [
                {
                    "data": {
                        "text/plain": [
                            "Problem Results Summary:\n",
                            "==================================================================================================================================\n",
                            "                                                                                                   Smallest    Largest    Clipped \n",
                            "Computation  GMM   Optimization   Objective   Fixed Point  Contraction  Objective    Gradient      Hessian     Hessian    Marginal\n",
                            "   Time      Step   Iterations   Evaluations  Iterations   Evaluations    Value    Infinity Norm  Eigenvalue  Eigenvalue   Costs  \n",
                            "-----------  ----  ------------  -----------  -----------  -----------  ---------  -------------  ----------  ----------  --------\n",
                            " 00:31:03     2         39           49          72322       217903     +2.1E+05     +2.8E+03      +1.8E+02    +1.8E+04      0    \n",
                            "==================================================================================================================================\n",
                            "\n",
                            "Nonlinear Coefficient Estimates (Robust SEs Adjusted for 999 Clusters in Parentheses):\n",
                            "=====================================================================================================\n",
                            "Sigma:      1         prices       hpwt        air         mpd        space     |   Pi:     1/income \n",
                            "------  ----------  ----------  ----------  ----------  ----------  ----------  |  ------  ----------\n",
                            "  1      +1.2E+00    +0.0E+00    +0.0E+00    +0.0E+00    +0.0E+00    +0.0E+00   |    1      +0.0E+00 \n",
                            "        (+2.1E+00)                                                              |                    \n",
                            "                                                                                |                    \n",
                            "prices               +0.0E+00    +0.0E+00    +0.0E+00    +0.0E+00    +0.0E+00   |  prices   -1.2E+01 \n",
                            "                                                                                |          (+4.8E+00)\n",
                            "                                                                                |                    \n",
                            " hpwt                            +0.0E+00    +0.0E+00    +0.0E+00    +0.0E+00   |   hpwt    +0.0E+00 \n",
                            "                                (+1.0E+01)                                      |                    \n",
                            "                                                                                |                    \n",
                            " air                                         +0.0E+00    +0.0E+00    +0.0E+00   |   air     +0.0E+00 \n",
                            "                                            (+2.6E+00)                          |                    \n",
                            "                                                                                |                    \n",
                            " mpd                                                     +2.6E+00    +0.0E+00   |   mpd     +0.0E+00 \n",
                            "                                                        (+1.4E+00)              |                    \n",
                            "                                                                                |                    \n",
                            "space                                                                +1.1E+00   |  space    +0.0E+00 \n",
                            "                                                                    (+1.9E+00)  |                    \n",
                            "=====================================================================================================\n",
                            "\n",
                            "Beta Estimates (Robust SEs Adjusted for 999 Clusters in Parentheses):\n",
                            "==========================================================\n",
                            "    1          hpwt        air         mpd        space   \n",
                            "----------  ----------  ----------  ----------  ----------\n",
                            " -7.6E+00    +4.8E+00    +2.6E+00    -2.2E+00    +1.8E+00 \n",
                            "(+1.2E+00)  (+2.7E+00)  (+1.4E+00)  (+1.5E+00)  (+1.3E+00)\n",
                            "==========================================================\n",
                            "\n",
                            "Gamma Estimates (Robust SEs Adjusted for 999 Clusters in Parentheses):\n",
                            "======================================================================\n",
                            "    1       log(hpwt)      air       log(mpg)   log(space)    trend   \n",
                            "----------  ----------  ----------  ----------  ----------  ----------\n",
                            " +2.3E+00    +6.4E-01    +8.8E-01    -3.9E-01    +2.8E-01    +1.4E-02 \n",
                            "(+1.4E-01)  (+1.2E-01)  (+7.1E-02)  (+1.2E-01)  (+1.2E-01)  (+3.5E-03)\n",
                            "======================================================================"
                        ]
                    },
                    "execution_count": 1,
                    "metadata": {},
                    "output_type": "execute_result"
                }
            ],
            "source": [
                "results = problem.solve(\n",
                "    initial_sigma,\n",
                "    initial_pi,\n",
                "    sigma_bounds=sigma_bounds,\n",
                "    pi_bounds=pi_bounds,\n",
                "    costs_type='log',\n",
                "    costs_bounds=(0.001, None),\n",
                "    W_type='clustered',\n",
                "    se_type='clustered'\n",
                ")\n",
                "results"
            ]
        },
        {
            "cell_type": "markdown",
            "metadata": {},
            "source": "There are some discrepancies between our results and the original paper. The instruments we constructed to are meant to mimic the original instruments, which we were unable to re-construct perfectly. We also use different agent data and the first-order linear approximation to the $\\ln(y_i - p_j)$ term."
        }
    ],
    "metadata": {
        "kernelspec": {
            "display_name": "Python 3",
            "language": "python",
            "name": "python3"
        },
        "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.6.6"
        }
    },
    "nbformat": 4,
    "nbformat_minor": 2
}