The BLP Contraction Mapping

Summary

Estimation requires inverting the random-coefficients share system to recover the vector of mean utilities ${\mathbf{\delta }}_t$ that rationalizes the observed market shares ${\mathcal{S}}_t$ for a given guess of ${\theta }_2$. Berry (1994) / BLP (1995) show this inversion can be computed by iterating the fixed point ${\mathbf{\delta }}_t^{h+1}\leftarrow {\mathbf{\delta }}_t^h+\log {\mathcal{S}}_t-\log {\mathbf{s}}_t({\mathbf{\delta }}_t^h,{\tilde{\theta }}_2)$, which is a contraction mapping with a unique fixed point. This is the “inner loop” of the nested fixed-point algorithm. Modern best practice replaces plain iteration with acceleration (SQUAREM) or Jacobian-based (Levenberg-Marquardt) solvers, which are 3-12x faster.

Overview

For each candidate ${\theta }_2$, and separately and in parallel for each market $t$, we must solve the $J_t$ equations in $J_t$ unknowns ${\mathbf{\delta }}_t$ that set predicted shares equal to observed shares. Holding ${\theta }_2$ fixed turns one large $N$-dimensional nonlinear system into $T$ small $J_t$-dimensional systems, which is the source of BLP’s scalability. The mapping ${\mathbf{\delta }}_t\equiv D_t^{-1}({\mathcal{S}}_t,{\tilde{\theta }}_2)$ is the Berry inversion.

Main Content

The system to solve ^share-system

In market $t$ we seek the $J_t$-vector ${\mathbf{\delta }}_t$ satisfying

$${\mathcal{S}}_{jt}=s_{jt}({\mathbf{\delta }}_t\mid {\theta }_2)=\int \frac{\exp ({\delta }_{jt}+{\mu }_{ijt})}{{\sum }_{k\in J_t}\exp ({\delta }_{kt}+{\mu }_{ikt})}f({\mathbf{\mu }}_{it}\mid {\tilde{\theta }}_2)\mathrm{d}{\mathbf{\mu }}_{it}.$$

A unique solution exists mathematically, but it cannot be solved exactly numerically; instead we solve to a tolerance expressed in the log difference in shares:

$$\Vert \log {\mathcal{S}}_t-\log {\mathbf{s}}_t({\mathbf{\delta }}_t,{\theta }_2){\Vert }_{\infty }\le {ϵ}^{tol}.$$

Tolerance is a tradeoff: too loose and numerical error propagates to $\hat{\theta }$ (Dubé et al. 2012); too tight and it can never be met. The recommended ${ϵ}^{tol}$ is between 1E-14 and 1E-12 (machine epsilon $\approx$ 1E-16 in double precision).

The BLP fixed point is a contraction ^contraction

Berry et al. (1995) show the relation $f({\mathbf{\delta }}_t)={\mathbf{\delta }}_t$ given by

$$f:{\mathbf{\delta }}_t^{h+1}\leftarrow {\mathbf{\delta }}_t^h+\log {\mathcal{S}}_t-\log {\mathbf{s}}_t({\mathbf{\delta }}_t^h,{\tilde{\theta }}_2)$$

is a contraction mapping. Iteration is linearly convergent at a rate proportional to $L({\tilde{\theta }}_2)/[1-L({\tilde{\theta }}_2)]$, where the Lipschitz constant (Dubé et al. 2012) is

$$L({\tilde{\theta }}_2)=\underset{{\mathbf{\delta }}_t}{\max }\Vert I_{J_t}-\frac{\partial \log {\mathbf{s}}_t}{\partial {\mathbf{\delta }}_t}({\mathbf{\delta }}_t,{\tilde{\theta }}_2){\Vert }_{\infty }<1.$$

A smaller $L$ converges faster. A larger outside-good share generally implies a smaller Lipschitz constant; as the outside share shrinks, convergence takes increasingly many steps. (Empirically, shrinking $s_{0t}$ from 0.91 to 0.27 raised iteration counts by ~5x.)

RCNL: the contraction must be dampened ^rcnl-dampen

For random-coefficients nested logit, the plain update is no longer a contraction; it must be dampened by $(1-\rho )$:

$${\mathbf{\delta }}_t\leftarrow {\mathbf{\delta }}_t+(1-\rho )[\log {\mathcal{S}}_t-\log {\mathbf{s}}_t({\mathbf{\delta }}_t,{\theta }_2)].$$

Convergence becomes arbitrarily slow as $\rho \to 1$ (more within-nest substitution), making RCNL harder to estimate. Without random coefficients (${\mu }_{ijt}=0$) the inversion has the closed form ${\delta }_{jt}=\log {\mathcal{S}}_{jt}-\log {\mathcal{S}}_{0t}-\rho \log {\mathcal{S}}_{j\mid ht}$ (Berry 1994).

Faster solvers: Newton / Levenberg-Marquardt and SQUAREM ^accelerated

Two families improve on plain iteration:

• Jacobian-based. Newton-Raphson updates ${\mathbf{\delta }}_t^{h+1}\leftarrow {\mathbf{\delta }}_t^h-\lambda {\Psi }_t^{-1}{\mathbf{s}}_t$, where ${\Psi }_t=\partial {\mathbf{s}}_t/\partial {\mathbf{\delta }}_t$. The Levenberg-Marquardt (LM) least-squares solver ${\min }_{{\mathbf{\delta }}_t}{\sum }_j[{\mathcal{S}}_{jt}-s_{jt}{]}^2$ is the fastest, most reliable Jacobian method; its update solves $[{\Psi }_t^{\prime }{\Psi }_t+\lambda \mathrm{diag}({\Psi }_t^{\prime }{\Psi }_t)]{\mathbf{x}}_t={\Psi }_t^{\prime }[{\mathcal{S}}_t-{\mathbf{s}}_t]$, interpolating between Gauss-Newton ($\lambda =0$) and gradient descent (large $\lambda$). Cost per iteration is high (computing $J_t\times J_t$ numerical-integral Jacobians).
• Accelerated fixed points. SQUAREM (Varadhan & Roland 2008) uses the residual $r^h=f({\mathbf{\delta }}_t^h)-{\mathbf{\delta }}_t^h$ and curvature to take an approximate Newton step without forming the Jacobian, at essentially the cost of plain iterations. It is 3-6x faster than direct iteration (3-8x fewer iterations in simulations) and is the recommended default. (No formal convergence guarantee, as the accelerated iteration is no longer a contraction; DF-SANE is a similar but slightly slower/less robust alternative.)

Examples

Inner loop for one guess of ${\theta }_2$ (NFXP step (a)):

1. For each market $t$, start ${\mathbf{\delta }}_t^0$ at the plain-logit solution $\log {\mathcal{S}}_t-\log {\mathcal{S}}_{0t}$ (or the previous ${\theta }_2$‘s solution, a “hot start” saving 10-20% of iterations).
2. Run SQUAREM (or LM) until $\Vert \log {\mathcal{S}}_t-\log {\mathbf{s}}_t{\Vert }_{\infty }\le$ 1E-14.
3. Return ${\hat{\mathbf{\delta }}}_t({\theta }_2)$ to the outer GMM loop. Markets run in parallel across processors.

Connections

• Random Coefficients Logit Model — defines the share integrals being inverted.
• GMM Estimation and Instruments for Price Endogeneity — the recovered ${\mathbf{\delta }}_t$ feeds the linear index and moment conditions (this is the “inner loop” of NFXP).
• Numerical Integration and Optimization in PyBLP — each share evaluation requires numerical integration; the contraction is the “inner loop,” optimization the “outer loop.”
• Method of Simulated Moments — the simulated shares make this a simulated fixed point.

See Also

• BLP Demand Estimation - Overview
• Supply Side and Markups
• _Index