Extracted main text — title through conclusion, appendix excluded. This is what our citation measures are computed over, published so the extraction can be checked by eye.
68,817 characters · 12 sections · 56 citation commands
Faster estimation of dynamic discrete choice models using index invertibility
In dynamic discrete choice modeling, estimation of the structural parameters that underlie economic decisions is often computationally challenging. Many available estimators for the structural parameter of interest $\theta_{0}\in\Theta$ are extremum estimators:
For instance, the criterion function $\hat{Q}$ may be the log-likelihood function john1988maximum, a pseudo log-likelihood function hotz1993conditional,arcidiacono2011conditional or a minimum distance function pesendorfer2008asymptotic. While these estimators offer appealing theoretical properties, they often impose substantial computational demands for multiple reasons. First, evaluating the criterion function may involve solving the model through costly fixed-point iteration or by simulation. Second, the criterion function's global concavity is not always guaranteed, often necessitating the use of global optimization methods or initializing a local optimization algorithm at various starting values. A relevant case is finite mixture models whose likelihood function may lack global concavity \parencites[e.g.,][ p. 182]{robert1999monte}[]{arcidiacono2011conditional}.
In this paper we harness the index restrictions inherent in many structural models to introduce an estimator for $\theta_0$ that offers substantial computational advantages and is asymptotically equivalent (of arbitrarily high order) to $\hat{\theta}^*$. Our focus is on models satisfying a condition we term `index invertibility'. Drawing from the semiparametric index regression literature, we describe a model as index invertible if its reduced form parameters are an invertible function of a vector of linear indices ahn2018simple. We establish that index invertibility implies a set of linear equality constraints which constrain $\theta_0$ to belong in a subspace of $\Theta$,\footnote{The subspace may be a strict subspace of $\Theta$ when there is at least one continuous covariate. We conjecture that it is possible to extend our method with inequality constraints when there is no continuous covariate khan2018discussion.} thereby reducing the dimension of the optimization problem presented in equation (ref). The main contribution of our paper is to propose an estimator which implements the constraints implied by index invertibility, and prove its asymptotic equivalence to the computationally intensive estimator $\hat{\theta}^*$.
Arguably, the class of index invertible structural econometric models is very broad. First, we prove that a broad class of dynamic discrete choice models with persistent unobserved heterogeneity (i.e., unobserved state variables) and partially linear flow utility satisfy index invertibility. In this leading example of index invertibility, the reduced form parameters are the conditional choice probabilities (defined as the probability of each choice conditional upon the observed state variables) which depend on multiple indices which govern the flow utility and transition of the state variables. Second, we do not restrict nor require specification of the number of indices required to attain index invertibility. Of course, as we show formally, the computational gains of our approach may diminish as the number of indices required to achieve index invertibility grows. Finally, the condition encompasses many invertible index models in the literature (see, for example, ahn2018simple and references therein).
Our approach is based on the observation that index invertibility implies a set of equality constraints on the structural parameter. Namely, we show that under index invertibility, the true parameter value $\theta_0$ satisfies $$\Sigma_{0}\boldsymbol{\gamma}(\theta_0)=0$$ for a known linear function $\boldsymbol{\gamma}(\cdot)$ and a nonparametrically identified matrix ${\Sigma}_0$. If $\Sigma_0$ were known, to solve the population version of equation (ref) it would be sufficient to search among $\boldsymbol{\gamma}(\theta)$ in the nullspace of $\Sigma_0$. Our estimator builds upon this idea and is defined by the following two steps: first, given an estimator $\hat{\Sigma}$ for $\Sigma_{0}$ (e.g., kernel smoothing in Section (ref)) we compute
Solving the optimization problem in equation (ref) is computationally simpler than the unconstrained problem in equation (ref) as it is only necessary to search over parameter values in $\{\theta\in\Theta\colon\hat{\Sigma}\boldsymbol{\gamma}(\theta)=0\}\subseteq\Theta$.\footnote{Since the constraints are linear in parameters, once estimated, they can be imposed generically at negligible computational cost. See Section (ref) for a general description.} The second step is to apply Newton-Raphson updates from $\tilde{\theta}$ in the direction of the target estimator $\hat\theta^*=\arg\max_{\theta\in\Theta}\hat{Q}(\theta)$ robinson1988stochastic. The resulting estimator is asymptotically equivalent to the more computationally intensive $\hat{\theta}^*$. Notably, given the typical statistical justification for $\hat{\theta}^*$ relies on asymptotic approximations (of a certain order), our proposed estimator inherits the favorable statistical properties of $\hat{\theta}^*$ but is computationally more efficient. To illustrate, if $\hat{\theta}^*$ stands for the parametric maximum likelihood estimator, then, under standard conditions, our method can achieve the Cram\'{e}r-Rao bound at lower computational cost by leveraging semiparametric index restrictions. Regardless, we emphasize that, under the conditions elaborated in the following sections, our method can be used to target any extremum estimator $\hat\theta^*$ within a model satisfying index invertibility---for example, if $\hat\theta^*$ is motivated by computational considerations, our method may be used to further reduce the computational burden.
As computational efficiency motivates our estimator, it is natural to explore the magnitude of possible computational benefits. Section (ref) provides some theoretical insights on this question. Recall that the computational gains arise from imposing the constraints $\Sigma_0\boldsymbol{\gamma}(\theta_0)=0$. Thus a key determinant of the computational benefits of our estimator is the rank of $\Sigma_0$: the larger the rank of $\Sigma_0$, the more restrictions $\Sigma_0\boldsymbol{\gamma}(\theta_0)=0$ places on $\theta_0$. Using the definition of $\Sigma_0$ (equation (ref)), we develop a series of results on the rank of $\Sigma_0$. Our results suggest two situations where the computational gains of our method will be large: either in the presence of many continuous covariates, or in the presence of at least one continuous covariate that satisfies a particular rectangular support condition. Our results also suggest that the rank of $\Sigma_0$ may decrease with the number of indices required to attain index invertibility.
To illustrate the advantages of our method, we consider a Monte Carlo experiment based on the utility function specification of toivanen2005market.\footnote{We consider a single-agent model instead of the two-agent model in toivanen2005market.} They estimate a dynamic model of firm entry into the U.K. fast food market between 1991 and 1995. In their model, firm profits from entry are determined by market size, which is modeled as depending on a long vector of socio-economic variables bresnahan1991entry,toivanen2005market,aguirregabiria2020identification. Due to the availability of these continuous socio-economic variables, our method is able to feasibly apply 8 restrictions to the 11-dimensional payoff parameter vector. Specifically, we consider two standard estimators as target estimators for our method: one based on the approach of bajari2011simple, the other on arcidiacono2011conditional. By simulating data from this model, we demonstrate that our estimator is, on average, around 20 times computationally more efficient than the standard approaches to estimating the model,\footnote{Our simulations demonstrate two benefits of effectively reducing the dimension of an optimization problem that requires repeated initializations: first, fewer initializations are required to reasonably cover the lower dimensional space; second, the lower dimensional problem is easier to solve per initialization. In our simulations, roughly, we use three times fewer initializations for the lower dimensional problem, and our estimator is around six times faster per initialization. Of course, while the number of required initializations may increase exponentially with the dimension of the optimization problem, it remains at the user's discretion and may be parallelized.} and provide empirical validation of our main theoretical result.
Our proposed method aims to contribute to a large literature on the computational aspects of structural modeling and, in particular, dynamic discrete choice \parencites[e.g., ][]{hotz1994simulation,arcidiacono2011conditional,bajari2011simple,su2012constrained,arcidiacono2013approximating,kristensen2021solving}. Rather than proposing an alternative to computationally advantageous estimators in the literature, our method can be used to improve the computation time for any estimator that can be expressed as the maximizer of a smooth sample criterion function. For instance, in our Monte Carlo experiment, we apply our method to estimators based on \textcites{arcidiacono2011conditional,bajari2011simple}, the latter of which is especially known to be computationally attractive. Parts of this paper are closely related to ahn2018simple, who develop a computationally simple estimator for a class of invertible index models. Whereas their paper focuses on identification and estimation of the index parameter, we allow the index parameter to be one part of a broader structural model and harness the semiparametric index restrictions for computational purposes.
The remainder of the paper is structured as follows. To fix ideas on our leading example of index invertibility, Section (ref) considers the specific context of dynamic discrete choice models. We then formally introduce the general model and index invertibility, and derive the equality constraints implied by index invertibility (Section (ref)). Section (ref) derives bounds on the rank of $\Sigma_0$, an important determinant of the number of independent restrictions in $\Sigma_0\boldsymbol{\gamma}(\theta_0)=0$. Section (ref) outlines the estimator and derives its equivalence to the computationally intensive estimator, our main result. Finally, Section (ref) presents our Monte Carlo experiment. In the Appendices we gather proofs, a proposed consistent estimator for $\Sigma_0$, and additional details on the Monte Carlo exercise.
In this section, we introduce the idea of `index invertibility' within the context of dynamic discrete choice (DDC) models. We have two goals for this section: first, to illustrate the idea of index invertibility through a simple and classical DDC model (Section (ref)); second, to introduce a broad, empirically relevant class of DDC models and prove that they satisfy index invertibility (Section (ref)). In Section (ref), we discuss index invertibility in a general model. We emphasize that the model in this section is just one example of sufficient conditions for index invertibility within a DDC model: the conditions specified in this section do not preclude other DDC models from possessing the index invertibility property (cf. Section 5).
In each period $t=1,2,\dots,T=\infty$, an agent observes a state variable $S_{t}$ and chooses an action $A_{t}\in\mathcal{A}=\{0,1\}$ to maximize their expected discounted utility.\footnote{The result of this section applies to $T<\infty$ (i.e., a non-stationary problem). We present only the $T=\infty$ case for notational ease.} The state variable is composed of two subvectors, $Z_{t}$ and $\epsilon_{t}$ which are observed and unobserved to the econometrician, respectively. We suppose the distribution of $\epsilon_t\in\mathbb{R}$ is known and has full support. The agent has time-separable utility and discounts future payoffs by the known rate $\beta_0\in[0,1)$.\footnote{Here and throughout the paper we use the subscript $0$ to indicate the true parameter value.} The period $t$ payoff for action $1$ is given by $Z_t^\intercal\gamma_0+\epsilon_t$, where $\gamma_0\in\mathbb{R}^{\dim(Z)}$ is unknown. The flow payoff to action $A_t=0$ is 0.
We assume $\{Z_t,\epsilon_t,A_t\}$ follows a stationary first-order Markov process and satisfies the following conditional independence assumption:
We further impose $$F_Z(z'|z,a)=G(z',\delta_0^\intercal z,a)$$ for some function $G$ and $\delta_0\in\mathbb{R}^{\dim(Z)\times J_2}$. This assumption is substantive when $J_2<\dim(Z)$,\footnote{The existence of $(\delta_0,G)$ is without loss of generality since one can always set $\delta_0$ equal to the identity matrix and $G=F_Z$ (in this case, $J_2=\dim(Z)$). As we explain below, the computational gain from our method comes from the rank of the matrix $\Sigma_0$ in equation (ref). In the trivial case of $\delta_0$ equaling identity, the rank of $\Sigma_0$ is zero, and our method does not provide a computational gain.} and is testable from observed data since $F_Z$ is nonparametrically identifiable.\footnote{For Lemma (ref), it is sufficient that the expected future utility $E[{v}(Z_{t+1})\mid Z_t=z,A_t=a]$ depends on $z$ only through $\delta_0^\intercal z$.} Since $G$ and $\delta_0$ are nonparametrically identified, they can be consistently estimated in a computationally feasible manner. Notice that the above conditions allow the special case that the transition of the state variable is deterministic (e.g., the lagged choice is a state variable), but in general the state transition may be unknown but nonparametrically identified.
Now define the conditional choice probability function as $$ \Pi_0(z)=\Pr(A_{t}=1\mid{Z}_{t}=z)=\Pr\left(v(1,z)+\epsilon_{t}>v(0,z)\right), $$ where $v(a,z) = \gamma_0^\intercal z 1\{a=1\}+\beta E[{v}(Z_{t+1})\mid Z_t=z,A_t=a]$ and $v(z)$ is the equilibrium ex-ante value function.\footnote{The ex-ante value function is defined as the discounted sum of future payoffs from optimal behavior given $Z_t=z$ but before the agent observes $\epsilon_t$ and chooses $A_t$, which we assume is well defined as a function of $z$ (john1988maximum. See, e.g., aguirregabiria2007sequential or bugni2021iterated.} In this simple model, the unknown structural parameter includes the payoff parameter $\gamma$, and a nuisance parameter $\delta$ which governs the state transition probabilities and is identified directly from the data. We now show that, if the state transition index $\delta_0^\intercal z$ is held fixed, the conditional choice probabilities $\Pi_0$ are an invertible function of the payoff index $\gamma_0^\intercal z$, an example of a property we term {`index invertibility'} (Assumption (ref)):
The condition in equation (ref) provides the basis of the computational savings of our method. To conclude the current section, we offer some intuition for this connection. To explain, suppose $\delta_0$ and $\Pi_0$ are known (note that they are identified directly from the data and can be estimated in a computationally attractive manner). If there are two values $z=z_1$ and $z=z_2$ that yield the same value of $(\delta_0^\intercal z,\Pi_{0}(z))$, then equation (ref) implies that the true payoff parameter $\gamma_0$ satisfies $\gamma_0^\intercal (z_1-z_2)=0$. Thus, purely from knowledge of this $z_1,z_2$, we learn that $\gamma_0$ belongs in the strict subspace of $\mathbb{R}^{\dim(\gamma)}$, $\{\gamma:\gamma^\intercal (z_1-z_2)=0\}$, and thus the effective dimension of $\gamma_0$ has decreased. In Section (ref) we formally introduce our estimator which builds on this idea, broadening the number of pairs of $Z$ which are used to find the lower dimensional subspace (i.e., Theorem (ref)), and taking into account that $\delta_0$ and $\Pi_0$ may be estimated.
In this section we introduce a broad class of dynamic discrete choice models that satisfy index invertibility. The choice set is $\mathcal{A}=\{0,1,\ldots,J_1\}$, and the state variable consists of three subvectors, $Z_{t}$, $\lambda_{t}$ and $\epsilon_{t}$ where $Z_{t}$ and $(\lambda_{t},\epsilon_{t})$ are observed and unobserved to the econometrician, respectively. The unobserved components $\epsilon_{t}$ are action specific, i.e., $\epsilon_{t}\in\mathbb{R}^{J_1+1}$, and we assume $\epsilon_{t}$ has full support. The period $t$ payoff for action $a$ is given by
where $\gamma_0(a)\in\mathbb{R}^{\dim(Z)}$, $\delta_{u,0}\in\mathbb{R}^{\dim(Z)\times J_{2,u}}$, and $f$ is a (possibly) nonlinear function.\footnote{We do not claim that the flow payoff is identified nonparametrically without, e.g., restrictions on $\delta_{u,0}$ and $f$. In applications of interest, $\delta_{u,0}$ may be a known matrix that governs which elements of $Z_{t}$ enter $f$, and $f$ may be a parametric function. As explained in Section (ref), we assume that some structural parameter of interest $\theta_0$ is point identified, and focus on harnessing the semiparametric structure to efficiently estimate the point identified parameter.} Note that, as a special case, some parts of the index parameter $\gamma(a)$ may be zero (and known), so that the corresponding elements of the observed state $Z_{t}$ affect payoffs only via the nonlinear function $f$. In other words, the flow utility function may be `partially linear' in its arguments---it allows for the case that only a subset of the state variables enter linearly.
We collect the parameter $\gamma_0(a)$ over $a\in\mathcal{A}$ in a matrix $\gamma_0=\left[\gamma_0(1),\ldots,\gamma_0(J_1)\right]\in\mathbb{R}^{\dim(Z)\times J_1}$, and impose the outside good assumption $ \gamma_0(0)=0$ and $f(\delta_{1,0}^\intercal Z_t,0,\lambda_{t})=0$. As usual, suppose $\{Z_t,\epsilon_t,\lambda_t,A_t\}$ follows a stationary first-order Markov process and satisfies
for some $(\delta_{F,0},G)$ such that $$F_Z(z'|z,a)=G(z',\delta_{F,0}^\intercal z, a).$$ The condition in (ref) strengthens the standard conditional independence assumption of rust1994structural by imposing a stronger independence condition on $\epsilon_t$, and a type of conditional independence between $\lambda_t$ and $Z_t$. Within this model, we define $\Pi_{0}(z)=\{\Pi_0(a,z)\colon a = 0, 1, \ldots, J_{1}\}$ as a mixture of the $\lambda_{t}$-specific conditional choice probability function, that is $$ \Pi_0(a,z)=\int\Pr\left(a=\arg\max_{\tilde{a}\in\mathcal{A}}\left\{v(\tilde{a},Z_t,\lambda_t)+\epsilon_{t}(\tilde{a})\right\}\mid Z_t=z,\lambda_{t}=l\right)dF_{\lambda_t}(l), $$ where $v(a,Z_t,\lambda_t) =\gamma_0(a)^\intercal Z_t+f(\delta_{u,0}^\intercal Z_t,a,\lambda_t)+\beta E[{v}(Z_{t+1},\lambda_{t+1})\mid Z_t,\lambda_t,A_t=a]$ and $v(Z_{t+1},\lambda_{t+1})$ is the equilibrium ex-ante value function. In Appendix (ref), we invoke the results of kasahara2009nonparametric to show that $\Pi_0$ can estimated directly from the data, and prove the next result that states that the above model is index invertible.
Relative to the illustrative model of Section (ref), the model for which we establish index invertibility in Theorem (ref) is substantially richer. First, it allows for multinomial choice (e.g., occupational choice, number of operating plants, etc.). Second, the model includes persistent unobserved heterogeneity $\lambda_{t}$. In the general case that $\lambda_{t}$ is time-varying, this variable may be referred to as the (vector of) `unobserved states'; in the special case that $\lambda_{t}$ is fixed across time (i.e., $\lambda_{t}=\lambda$), it is often called the agent/market's `type' or permanent unobserved heterogeneity. Furthermore, the persistent unobserved heterogeneity may enter the flow payoff function non-linearly with states through $f$. Fourth, the flow utility function is allowed to be non-linear in states. Fifth, as noted above, as a special case the flow utility may be `partially linear'---that is, only some part of the flow utility may be linear---which appears to be a common functional form in empirical work (see, among others, \textcites{arcidiacono2005affirmative,todd2006assessing,kennan2011effect,beffy2012choosing}).\footnote{At some notational cost, we could generalize the flow utility function in equation (ref) to $h(\gamma_0(a)^\intercal Z_t)+f(\delta_{u,0}^\intercal Z_t,a,\lambda_{t})+\epsilon_t(a)$ where $h$ is an invertible function. For example, $h$ may be the exponential function.} Thus, anticipating the results of Section (ref), our method has the potential to generate computational gains in DDC models with `partially linear' flow utility in the presence of continuous observed state variables.
In this paper we are interested in learning a finite-dimensional parameter vector $\theta_0\in\Theta$ that is identified as the unique maximum of a population criterion function ${Q}(\theta)$: $$ {\theta}_0=\arg\max_{\theta\in\Theta}{Q}(\theta). $$ For instance, $\theta_0$ may be the point-identified structural parameter (sub-)vector in a dynamic discrete choice model. If $\hat{Q}$ is an estimator for $Q$, one may estimate $\theta_0$ by $$ \hat{\theta}^*=\arg\max_{\theta\in\Theta}\hat{Q}(\theta). $$ However, in many cases, finding the maximum of $\hat{Q}(\theta)$ over $\Theta$ may be computationally challenging. For instance, if $\hat{Q}(\theta)$ represents the sample log-likelihood function of a DDC model rust1994structural,aguirregabiria2002swapping, $\hat{Q}(\theta)$ may lack a known closed form, requiring iterative or simulation methods to compute. Moreover, in the presence of unobserved types and/or states, $\hat{Q}(\theta)$ may lack global concavity arcidiacono2011conditional, necessitating repeated initialization of the optimization algorithm.
Our goal is to obtain an asymptotically equivalent estimator to $\hat{\theta}^*$ in a computationally feasible way. We achieve this by incorporating the following restriction into the optimization.
We refer to Assumption (ref) as index invertibility. It states that for a known linear function $\boldsymbol{\gamma}(\theta)$ of the parameter of interest, the random variable $Z$ can be used to construct a vector of indices $[\gamma_0,\delta_0]^\intercal Z$ for which $\Pi_0(Z)$ is an invertible function of $\gamma_0^\intercal Z$, while $\delta_0^\intercal z$ is held fixed. It is worth noting that the qualifier $\delta_0^\intercal z_1=\delta_0^\intercal z_2$ is included to make Assumption (ref) apply more generally: we allow for the case that $\delta_0$ is the $ \dim (Z)\times 1$ zero vector. We can interpret $J_1+J_2$ as the number of indices required to achieve index invertibility. For example, in the DDC example of Section (ref), $J_1$ is the number of index parameters that enter the flow utility (i.e., the size of choice set less the outside option), and $J_2$ is the number of indices that enter the state transition. As we show in Section (ref), the computational benefits of our approach depend importantly on $J_1$ and $J_2$.
In practice, one's choice of $(\delta_0,\Pi_0)$ is best guided by application-specific knowledge and the chosen structural model. For example, in the model of Section 5, one part of the flow profits depends on a set of time-invariant socio-economic variables. To be more specific, although profits depend on both the (time-varying) number of operating stores and the (time-invariant) socio-economic variables, we assume that the available market size depends only on the latter. Then, as we show, the probability of operating at least one store in a given time period is increasing in the available market size. In this case, we can use $\delta_0$ to hold fixed the time period, and $\Pi_0$ is naturally related to the probability that the firm operates at least one store.
In order to exploit index invertibility, we define the matrix
where $Z_1$ and $Z_2$ are independent random variables with the same marginal distribution as $Z$. The next result shows that $\Sigma_{0}$ characterizes the equality constraints that are implied by index invertibility.
Theorem (ref) shows that index invertibility (Assumption (ref)) implies that $\theta_0$ satisfies $\dim (Z)\times J_1$ equality constraints, namely that $\theta_0\in\{\theta\in\Theta\colon{\Sigma}_0\gamma_0=0\}\subseteq\Theta$. If $\Sigma_0$ were known, then imposing the equality constraints in the optimization problem (equation (ref)) necessarily eases the computational burden, since the search is limited to a smaller set of possible parameter values. Our estimator (described in Section (ref)) builds on these ideas.
Equation (ref) suggests there are $\dim (Z)\times J_1$ restrictions on $\theta_0$, however, in practice, these restrictions may be linearly dependent. From equation (ref), a key determinant of the number of linearly independent restrictions is the rank of $\Sigma_0$: in the extreme case that $\Sigma_0=0$, there are no restrictions on $\theta_0$ from $\Sigma_0\gamma_0=0$; in the other extreme case that the rank of $\Sigma_0$ is $\dim (Z)-1$, then there are $(\dim (Z)-1)\times J_1$ restrictions on $\boldsymbol{\gamma}(\theta)$ ahn2018simple. Given the importance of the number of linearly independent restrictions to the benefits of imposing the equality constraints, we now provide some results on the rank of $\Sigma_0$ (Section (ref)).
In this section, we consider the rank of $\Sigma_{0}\in\mathbb{R}^{\dim(Z)\times \dim(Z)}$, which determines the strength of restrictions implied by index invertibility. Under index invertibility, each column of the structural parameter $\gamma_0\in\mathbb{R}^{\dim(Z)\times J_1}$ belongs in the nullspace of $\Sigma_{0}$, which has dimension $\dim(Z)-\mathrm{rank}(\Sigma_0)$ by the rank-nullity theorem. Ergo, the effective number of restrictions on $\gamma_0$ implied by index invertibility is rank$(\Sigma_0)\times J_1$. That is, the larger the rank of $\Sigma_{0}$, the greater the computational advantage of imposing the equality constraints $\Sigma_0\gamma_0=0$. In this section, we provide some sufficient conditions under which a lower bound on $\mathrm{rank}(\Sigma_0)$ can be characterized.
To summarize the findings in broad terms, the results of this section provide two routes to achieving a high $\mathrm{rank}(\Sigma_0)$: either by having many continuous components of $Z$ (Theorem (ref)), or by having one continuous component of $Z$ that satisfies a particular support condition (Theorem (ref)). In terms of practical guidance for applied work, the results of this section can be viewed as suggestive of the type of empirical settings where our results may generate large effective dimension reduction. Namely, in an index invertible model with either many continuous covariates, or one continuous covariate that satisfies a rectangular support condition, our method may generate large computational savings. For example, in the dynamic discrete choice model studied in Section (ref), $Z$ includes 9 continuous components, and we can apply Theorem (ref) to show $\Sigma_0$ implies 8 restrictions on the structural parameter.
The first theorem provides a lower bound on the rank of $\Sigma_0$ which depends on the number of continuous components of $Z$, but may be lower when the number of indices $J_1+J_2$ is larger.
By the interior-point assumption, the variable $Z_A$ is continuously distributed (given $Z_B$). The term $\mathrm{rank}(Var([\gamma_0,\delta_0]^\intercal[Z_{A}^\intercal,0^\intercal]^\intercal))$ represents how many components in $[\Pi_0(Z),Z^\intercal\delta_0]^\intercal$ are continuously distributed. It is naturally bounded above by $J_1+J_2$, the number of indices required to achieve index invertibility---Theorem (ref) states that is preferable for this number to be small relative to the number of continuous components of $Z$. In particular, the second part of Theorem (ref) states it is desirable for the non-structural index $\delta_0^\intercal Z$ to depend only on discrete components of $Z$. In this case $J_1$ is an upper bound for $\mathrm{rank}(Var([\gamma_0,\delta_0]^\intercal[Z_{A}^\intercal,0^\intercal]^\intercal))$. To provide a concrete example, in a dynamic discrete choice problem this would occur if the state transition depended only upon lagged actions and discrete state variables.
The second theorem states that if one component of $Z$ satisfies an additional condition, then the lower bound on $\mathrm{rank}(\Sigma_0)$ does not depend on the {number} of continuous components of $Z$. To show this result, we modify the arguments of horowitz1996direct to the current framework.
Relative to Theorem (ref), Theorem (ref) provides an improved lower bound by depending on the length of $Z$ instead of the number of continuous components in $Z$. This improved bound is available when $(\gamma_0,\delta_0)^\intercal Z$ satisfies a rectangular support assumption.
In this section, we introduce our estimator for $\theta_0$. Our method is motivated by the computational difficulty of an available estimator $\hat\theta^*=\arg\max_{\theta\in\Theta}\hat{Q}(\theta)$, where $\Theta$ is a subset of a Euclidean space and $\hat{Q}:\Theta\rightarrow\mathbb{R}$ is a sample criterion function. As discussed previously, in many cases the estimator $\hat{\theta}^*$ is computationally heavy, or may even be computationally infeasible in practice. For example, maximum likelihood estimation of finite-mixture dynamic discrete models is considered extremely computationally costly, to such a degree that alternative estimators are often preferred arcidiacono2011conditional.
Our estimator $\hat{\theta}$ for $\theta_0$ is constructed in the following two steps:
In the first step, we form a preliminary estimator $\tilde{\theta}$ by maximizing the sample criterion function subject to the estimated constraints $\hat{\Sigma}\boldsymbol{\gamma}(\theta)=0$. Since these are (low dimensional) linear equality constraints, they can be easily imposed in practice with negligible computational cost (see Section (ref) for one approach). The second step consists of $L$ Newton-Raphson iterates from the preliminary estimator $\tilde{\theta}$. The main result of this section (Theorem (ref)) states that the number of Newton-Raphson iterates controls the rate at which $\hat{\theta}-\hat{\theta}^*$ converges to zero as sample size $n$ diverges (i.e., $n\rightarrow\infty$). The remainder of this section is dedicated to showing this result, which will use two additional assumptions.
The first step of our estimator solves a maximization problem subject to the estimated constraint $\hat{\Sigma}\gamma=0$. Naturally, we require that the estimated constraint $\hat{\Sigma}\gamma=0$ provides a good approximation to ${\Sigma}_0\gamma=0$, which we formalize in Assumption (ref).
The first part of Assumption (ref) states that $\hat{\Sigma}$ is consistent for ${\Sigma_0}$. Notably, the rate of convergence need not be known by the econometrician. In particular, we allow the rate of convergence to be arbitrarily slow: Theorem (ref) implies that even if the convergence rate is slow, only moderate increases in $L$ are required to attain fast convergence between our estimator and the computationally intensive estimator. Many nonparametric methods can achieve consistent estimation (e.g., kernel smoothing, nearest neighbor, splines, or series estimators). In Section (ref), we provide conditions for consistent estimation using kernel smoothing ahn2018simple.
The second part of Assumption (ref) states that $\hat{\Sigma}$ is rank-consistent, which ensures that $\hat{\Sigma}\gamma=0$ imposes the same number of linearly independent constraints as $\Sigma_0\gamma=0$ with probability approaching one. Given a consistent estimator $\tilde\Sigma$, which may or may not have the same rank as $\Sigma_0$, one may construct a rank-consistent estimator by a low rank approximation.\footnote{Instead of this approach of $\hat\Sigma$, we may be able to apply a rank estimator, e.g., in chen2019improved. Since our results rely only on the convergence rate of $\tilde\theta$ and the rank is correctly estimated with probability approaching one, we conjecture that estimating the rank does not change our main result.} To explain, let $\hat\lambda_1\geq\cdots\geq\hat\lambda_K$ be the eigenvalues of $\tilde{\Sigma}$, and $\hat\nu_1,\cdots,\hat\nu_K$ be the corresponding eigenvectors. Define the low-rank approximation
where $\kappa$ is a threshold value. The following result (Lemma (ref)) states that the low-rank approximation $\hat{\Sigma}$ satisfies Assumption (ref) as long as $\kappa$ converges to zero slowly.
It is worthwhile noting that any estimator for $\Sigma_0$ satisfying Assumption (ref) will be constructed using estimators of $\Pi_0$ and $\delta_0$. Typically, both the reduced form parameter $\Pi_0$ and the nuisance parameter $\delta_0$ are identified directly from the data, and can thus be estimated without constructing the function $\hat{Q}$. For example, in the Monte Carlo exercise of Section (ref), $\delta_0$ is known a priori, and $\Pi_0$ are conditional probabilities which we estimate via kernel regression.
The computationally intensive estimator $\hat\theta^*$ is an example of an extremum estimator. Assumption (ref) imposes mild regularity conditions that are typical in extremum estimation problems.
We now state the main theoretical result of this paper.
Theorem (ref) states that our estimator $\hat{\theta}$ is asymptotically equivalent to the computationally more intensive estimator $\hat{\theta}^*$. In particular, that the difference $\hat\theta -\hat\theta^*$ converges to zero at the rate $\max\{\|\hat{\Sigma}-{\Sigma_0}\|,n^{-1/2}\}^{2^L}$. Let us now provide some intuition for the rate of convergence. First, the term $\max\{\|\hat{\Sigma}-{\Sigma_0}\|,n^{-1/2}\}$ represents the convergence rate of $\tilde\theta-\hat\theta^*$, i.e., the difference between the start-up and target estimators for the Newton-Raphson iterations. The convergence rate can be understood as follows. Because $\hat\Sigma\tilde\theta=0$, the difference $\tilde\theta-\hat\theta^*$ is proportional to $\hat{\Sigma}\hat{\theta}^*$ whose convergence rate depends on $\hat{\Sigma}-{\Sigma}_0$ and $n^{-1/2}$ (from $\hat{\theta}^*-\theta_0$). Second, the exponent $2^L$ represents the effect of $L$ Newton-Raphson iterations from the start-up estimator $\tilde{\theta}$. As in robinson1988stochastic, the rate of convergence of $\hat\theta$ to $\hat{\theta}^*$ increases exponentially in the number of Newton-Raphson updates $L$.
A practical consideration for our estimator is how to choose the number of Newton-Raphson iterations $L$. Our main theoretical result (Theorem (ref)) suggests that $L$ should be chosen to achieve the desired rate of convergence between $\hat\theta$ and $\hat\theta^*$. For example, if $\hat\theta^*$ is justified by first-order asymptotics, then $L$ can be chosen to achieve first-order asymptotic equivalence between $\hat\theta$ and $\hat\theta^*$, which is attained with one Newton-Raphson update when $\hat{\Sigma}-{\Sigma_0}=o_p(n^{-1/4})$. If $\hat\theta^*$ has desirable higher-order asymptotic properties, then $L$ can be set to a larger number. Importantly, because $L$ impacts the rate of convergence through the exponent $2^L$, fast convergence of $\hat\theta-\hat\theta^*$ can be attained for moderate $L$. Of course, extra Newton-Raphson iterations impose additional computation costs. However, our experience in simulations suggests that the computational cost of Newton-Raphson updates (i.e., Step 2 of our estimator) may be small relative to solving the constrained optimization problem (i.e., Step 1). Overall, consideration of theoretical and empirical aspects suggests choosing $L$ as small as possible to achieve the desired degree of asymptotic equivalence.
This section investigates the performance of our proposed estimator in a Monte Carlo simulation. The two primary objectives in this section are, first, to quantify the computational benefits of our approach and, second, to provide empirical support for our main asymptotic equivalence result (Theorem (ref)). Section (ref) contains further details on the Monte Carlo design and implementation.
Our Monte Carlo design is based on the empirical setting of toivanen2005market, which analyzes firm entry into the U.K. fast food market between 1991 and 1995. Restricting attention to the largest two firms, their analysis divides the U.K. into 422 local markets and records information about each market and the firms' decisions on how many stores to operate in each market. To fit it within the scope of this paper, we model a single firm's decision as a dynamic discrete choice problem with the profit function in the spirit of bresnahan1991entry,toivanen2005market,aguirregabiria2020identification.
In each period and geographic market, a firm decides whether to open an additional store, upon observation of the state variables. The firm's decision in market $i$ and time $t$ is $A_{i,t}\in\{0,1\}$, which takes value $1$ if the firm opens a store in market $i$ at time $t$, and $0$ otherwise. In each period $t$, the vector of state variables known by the firm in market $i$ is $(N_{i,t},W_i^\intercal,\lambda_{i},\epsilon_{i,t})^\intercal$ where $N_{i,t}$ is the number of incumbent stores (that is, prior to the realization of $A_{i,t}$) that we assume is bounded above by 3, $W_i$ are variables that affect the size of market $i$, $\lambda_{i}$ is the market type, and $\epsilon_{i,t}$ is an idiosyncratic shock. Firms are assumed to be forward looking, choosing $A_{i,t}$ to maximize expected discounted profits. We assume that $(N_{i,t},W_i^\intercal)^\intercal$ and $(\lambda_i,\epsilon_{i,t})^\intercal$ are observed and unobserved to the econometrician, respectively.
Denoting the observed state variable $(N_{i,t},W_i^\intercal)^\intercal$, period profits from opening an additional store in market $i$ at time $t$ (i.e., $A_{i,t}=1$) are equal to $$ \left(\lambda_i + \theta_W^\intercal W_i\right) - \left(\theta_{FC} N_{i,t} + \theta_{EC} \mathsf{1}(N_{i,t}=0) + \epsilon_{i,t}\right). $$ The first component $(\lambda_i + \theta_W^\intercal W_i)$ of the period profits is the marginal revenue from opening an additional store, i.e., the market size. Following the literature toivanen2005market,aguirregabiria2020identification, we allow market size to depend on an unobserved market type ($\lambda_i$) and a long vector of continuous demographic and socioeconomic variables. For example, in toivanen2005market, market size depends on total population, youth population, and pensioner-age population. In aguirregabiria2020identification market size depends additionally on population density, the local unemployment rate and local GDP per capita. Inspired by these papers, we set $W_i\in\mathbb{R}^{\dim(W_i)}$ with $\dim(W_i)=9$. Also, we impose that $W_i$ and $\lambda_i$ are statistically independent. The second component $\left(\theta_{FC} N_{i,t} + \theta_{EC} \mathsf{1}(N_{i,t}=0) + \epsilon_{i,t}\right)$ represents the marginal cost of opening an additional store, which depends on the firm's local experience. Also following toivanen2005market,aguirregabiria2020identification, we assume $\epsilon_{i,t}$ is an opening cost shock that is known to follow the standard normal distribution. The period payoff from not opening an additional store (i.e., $A_{i,t}=0$) is normalized to zero. Firms discount future payoffs at a discount factor $\beta=0.95$, which we assume to be known. Finally, observe that $N_{i,t+1}$ is uniquely determined by $(N_{i,t},A_{i,t})$, i.e., $N_{i,t+1}=\max\{N_{i,t}+A_{i,t},3\}$, so the transition of the state variable is deterministic. We collect the part of the structural parameter that governs the flow payoff as $\theta=(\theta_W^\intercal,\theta_{FC},\theta_{EC})^\intercal\in\mathbb{R}^{9+1+1}$. The remaining component of the structural parameter is the distribution of types $F_{\lambda}$.
We now discuss Assumption (ref) and the dimension reduction of our proposed method. Let $Z_{i,t}=(W_i^\intercal,t)^\intercal$, $\gamma(\theta)=(\theta_W^\intercal, 0)^\intercal$, $\delta_0=(0_{\dim(W)}^\intercal,1)^\intercal$, and $\Pi_0(z)=Pr(N_{i,t+1}\ge 1\mid W_i=w)$. With these definitions, Assumption (ref) holds if, for every pair $z$ and $\tilde{z}$ with $t=\tilde{t}$,
That is, the probability of having a store operating in a given period is higher if and only if the payoff index $\theta_W^\intercal W_i$ is higher. In Appendix (ref) we show that, under some conditions, equation (ref) holds. Thus $\Sigma_0\boldsymbol{\gamma}(\theta_0)=0$ where $\Sigma_0=E[(Z-\tilde{Z})(Z-\tilde{Z})^\intercal\mid \Pi_0(Z)=\Pi_0(\tilde{Z}),\delta_0^\intercal Z=\delta_0^\intercal\tilde{Z}]$ for $Z$, $\tilde{Z}$ independent draws of $Z_{i,t}$. Given that $\mathrm{Supp}(W_i)$ has a non-empty interior and $\delta_{0}^\intercal Z_{i,t}=t$ is discrete, the second part of Theorem (ref) applies and $\mathrm{rank}(\Sigma_0)\ge 9-1=8$, so that the identity $\Sigma_0\boldsymbol{\gamma}(\theta_0)=0$ may reduce the dimension of the optimization problem by $\dim(W_i)-1=8$.
Finally, we describe our choice of simulation parameters. In our design, we observe the set $\{N_{i,t},W_i,A_{i,t}:t=1,2\ldots,8\}_{i=1}^n$ for $n=100$, $200$, $350$, and $500$ i.i.d. markets (recall in the dataset of toivanen2005market, $n=422$ and $T=5$). Our simulation results are based on 100 draws from this data generating process with the following value of the structural parameter: $\mathrm{Supp}(\lambda_i)=\{0.1, 1\}$, $\Pr(\lambda_i=1)=0.63$, and $\theta_0=((-0.3,-0.2,-0.1,0.1,,0.2,0.3,0.4,0.5,-0.6),0.5,0.5)^\intercal$. We initialize $N_{i,1}=0$ and draw $W_i$ from the uniform distribution over $[0,1]^{\dim(W_i)}$.
We consider two target estimators for our method: one based upon arcidiacono2011conditional, the other upon bajari2011simple. Both of our target estimators are extremum estimators that use the sample log-likelihood function as the criterion, that is
where the likelihood contribution $\ell_{i,t}(v,\theta)$ is the model-implied probability of $A_{i,t}$ conditional upon $N_{i,t},W_i$ and $\lambda_i=v$ evaluated at $\theta$. The estimators differ in their approach to modeling the distribution of $\lambda_i$. Specifically, while both approaches model ${\lambda_i}$ as a discrete random variable, arcidiacono2011conditional's estimator fixes the number of support points but not their locations, whereas bajari2011simple's estimator specifies a (possibly large) number of fixed points that are known to contain the support.\footnote{The estimator of bajari2011simple approach has also been analyzed as a nonparametric estimator in certain discrete choice problems \parencites[e.g.,][]{fox2016simple,bunting2020}.} Section (ref) provides details on our definition of these estimators. For the estimator of arcidiacono2011conditional, we correctly specify that $\lambda_i$ has two unknown points of support. For the estimator based on bajari2011simple, we allow the support to contain any or all of 21 uniformly spaced points from $-0.5$ to $1.5$. Henceforth, we refer to our formulation of arcidiacono2011conditional's estimator as $\hat\theta^*_{EM}$, and the estimator based upon bajari2011simple as $\hat\theta^*_{H}$.
For each target estimator, our estimator is constructed in the following two steps:
For the target estimators $\hat\theta^*_{H}$ and $\hat\theta^*_{EM}$, we denote our corresponding estimator as $\hat\theta_{H}$ and $\hat\theta_{EM}$, and the first step estimators as $\tilde\theta_{H}$ and $\tilde\theta_{EM}$, respectively. To illustrate the computational comparison with two standard approaches to dynamic discrete choice estimation, we also compute $\hat\theta^*_{H}$ and $\hat\theta^*_{EM}$ for each simulated dataset.
Estimation of this model may be computationally intensive for two main reasons. First, the choice probabilities $\ell_{i,t}(v,\theta)$ are generated as the solution to a dynamic programming problem, which is solved by iterating a contraction mapping until convergence separately for each candidtate $(v,\theta)$.\footnote{To solve the model, we iterate on the conditional choice probability mapping bugni2021iterated.} Second, the presence of unobserved types $\lambda_i$ means that the observed data is an unknown mixture of type-specific choice models. Even if $\lambda_i$ is assumed to have finite support, its presence means that the log-likelihood function may not be globally concave, which creates the risk that the optimization algorithm converges to a non-global optimum. To mitigate this risk, one may rerun the optimization algorithm a number of times, each run starting from a different initial value robert1999monte. One ideal approach would be to initialize the estimator at each element of $\times_{d=1}^D \{i_{d,1},i_{d,2},\ldots,i_{d,P}\}$ where $D$ is the dimension of the optimization problem, and $i_{d,p},p=1,\dots,P$ are the vertices of a one dimensional grid. We choose $P=11$ in this section. For moderate $D$, this ideal approach quickly becomes infeasible. To see this, the ideal approach described above would require ${P}^3$ initializations to form $\tilde\theta$ with $D=3$, and $P^{11}$ initializations to form $\hat\theta^*$ with $D=11$. Instead, we proceed by randomly sampling $2\times D$ times from $\times_{d=1}^D \{i_{d,1},\ldots,i_{d,P}\}$ without replacement, where we choose $i_{d,p},p=1,\dots,P$ to be evenly spaced points around an initial guess based on the (parametric) pseudo maximum likelihood estimator. See Section (ref) for a formal specification of the initial points. We then define $\tilde\theta$ and $\hat\theta^*$ as the maximizer of the sample criterion over the set of initial points.\footnote{More precisely, for $j=H,EM$, if $\mathcal{I}$ is the set of initial points and $\tilde\theta_{j}^{(i)}$ (resp. $\hat\theta_{j}^{*,(i)}$) is the solution to the constrained (resp. unconstrained) optimization problem initialized at $i\in\mathcal{I}$, then $\tilde\theta_{j}=\operatorname*{arg\,max}_{\{\tilde\theta_{j}^{(i)}\colon i\in\mathcal{I} \}}\hat{Q}(\theta)$ (resp. $\hat\theta_{j}^*=\operatorname*{arg\,max}_{\{\hat\theta_{j}^{*,(i)}\colon i\in\mathcal{I} \}}\hat{Q}(\theta)$).}
Table (ref) displays the mean computation time in minutes for our first step estimator $\tilde\theta$, our final estimator $\hat\theta$, and the target estimator $\hat\theta^*$, for each sample size and estimation method.\footnote{Simulations were run in the Julia programming language, and replication files are available at \href{https://doi.org/10.5281/zenodo.15116658}{https://doi.org/10.5281/zenodo.15116658}. Results for $\hat\theta_H$ were generated using Julia 1.11 on a standard desktop computer, and results for $\hat\theta_{AM}$ were generated using Julia 1.10 on the JuliaHub cloud-based server. The replication files contain a complete description of the two environments.} A number of observations can be made. First, our estimator is substantially faster than the standard estimator. For example, for $n=500$, our method for targeting the estimator based on bajari2011simple is over 20 times faster than the standard approach for $n=500$, which translates to reducing the computation time from around 9 hours to just over 25 minutes. For $n=500$, our method for targeting arcidiacono2011conditional's estimator is over 16 times faster than the standard approach, which translates to reducing the computation time from over 9 hours to around one half hour. Second, the computational cost of the Newton-Raphson iterates (i.e., Step 2 of our estimator) is minor relative to solving the constrained optimization problem (i.e., Step 1 of our estimator). Namely, the difference between the computational times for obtaining $\tilde\theta$ and for obtaining $\hat\theta$ is small. Third, we observe that the computational burden increases approximately linearly with sample size. This means the computational time savings are substantially bigger when $n$ is large---that is, the gain from our method is larger for the `harder' computational problem.
The computational gains from our method described in Table (ref) can be decomposed as the product of two effects. Namely, our method reduces the number of grid points (and thus initializations) designed to cover the parameter space, and also reduces the computational burden per initialization of the optimization algorithm.\footnote{As described above, we choose 7 ($=2\times3+1$) and 23 ($=2\times 11+1$) random starting values for the constrained $\hat\theta$ and unconstrained $\hat\theta^*$, respectively. However, clearly, the computational gain from reducing the dimension of the grid can be made exponentially larger by requiring more points for each dimension of the grid. We also point out that parallelization can help mitigate this part of the computational burden.} We explore this decomposition in Table (ref) which presents the average runtime per initialization, for each sample size and estimation method. We observe that, per initialization, our method is around 6 times faster than the standard, unconstrained approach on average. For example, for $n=500$, our method of targeting the estimator based upon bajari2011simple is around 6 times faster per iteration than the standard approach, which translates to reducing the computation time from 24 minutes to around 4 minutes per run. Similarly, for $n=500$, the constrained EM algorithm is around 4.5 times faster per iteration than the unconstrained approach, which translates to reducing the computation time from around 25 minutes to 6 minutes per run.
We now turn to results on the convergence between our estimator and the corresponding target estimators. For each target estimator, Table (ref) shows that the $\sqrt{n}$-scaled root mean squared error of $\hat\theta-\hat\theta^*$ tends to decrease with sample size. Thus it appears that the difference of the estimators is converging to zero at faster than $\sqrt{n}$ rate, which provides empirical support for our main theorem (Theorem (ref)).\footnote{For sake of comparison, we report the $\sqrt{n}$-scaled root mean squared error of $\tilde\theta-\hat\theta^*$ in Table (ref) of Appendix (ref). In general, the $\sqrt{n}$-scaled root mean squared error of $\hat\theta-\hat\theta^*$ is much smaller than that of $\tilde\theta-\hat\theta^*$.}
To conclude, as established in Tables (ref) and (ref), our estimator can be computed substantially faster than the target estimator. Table (ref) provides evidence that the two estimators are asymptotically equivalent, consistent with our main result (Theorem (ref)). Taken together, the two findings demonstrate that our method can be used to complement both the estimator of arcidiacono2011conditional and bajari2011simple: regardless of which approach to modeling $F_\lambda$ is preferred, our method can be used to form an asymptotically equivalent estimator, with a much lighter computational burden.
In this paper we provide a method to simplify estimation of structural economic models whose reduced form parameters are an invertible function of a number of linear indices. A leading example of such `index invertibility' are dynamic discrete choice models with persistent unobserved heterogeneity and partially linear flow payoffs. Index invertibility implies a set of equality constraints which restrict the structural parameter of interest to belong in a subspace of the parameter space. We propose an estimator that imposes the equality constraints, and show it is asymptotically equivalent to the unconstrained estimator. The proposed constrained estimator may be computationally advantageous due to the effective reduction in the dimension of the optimization problem. Furthermore, we provide a number of results on the extent of effective dimension reduction, and demonstrate our method in Monte Carlo simulations.
\paragraph{Acknowledgments}
We thank Aureo de Paula and the review panel for excellent comments and suggestions that have improved the manuscript. We also thank Yanqin Fan and seminar participants at Boston College, Georgetown, Duke, Simon Fraser and the University of Melbourne for helpful comments. The usual disclaimer applies.
\addcontentsline{toc}{section}{\refname} \printbibliography