PaperlayerView exact revision on arXiv

Source-linked AI summary

SAUSS: Stochastic Approximation with Unbiased Simulated Scores for Limited Dependent Variable Models

Sokbae Lee, Yuan Liao, Myung Hwan Seo, Youngki Shin

arXiv:2608.25304v1stat.MEcs.LGecon.EM

TL;DR

Multinomial probit estimation combines flexible substitution patterns with costly high-dimensional probability calculations, while fixed-budget simulated likelihood can be biased. SAUSS addresses this with averaged mini-batch stochastic approximation using unbiased simulated scores, achieving broadly comparable accuracy with substantially less computation under the reported configurations.

  • Problem

    Multinomial probit models offer flexible covariance structures and richer substitution patterns, but high-dimensional region probabilities and repeated full-sample likelihood evaluations make estimation computationally demanding, while fixed-budget simulation can bias the log-likelihood.

  • Method

    SAUSS combines exact conditional simulation and unbiased simulated scores with mini-batch stochastic approximation and iterate averaging, while accounting for mini-batch and simulation variation in inference.

  • Results

    SAUSS produces broadly similar estimates and accuracy to GHK-based simulated maximum likelihood while requiring more than 100 times less computation time under the reported configurations.

  • Takeaways & Limitations

    SAUSS provides a computationally feasible likelihood-based approach as the number of alternatives increases and extends to limited dependent variable models with suitable score representations and exact conditional simulation.

  • Takeaways & Limitations

    Conditional on a fixed sample, the iterates’ path variation reflects randomized-algorithm variation around the sample target and does not account for sampling uncertainty.

Abstract

from arXiv · show

Multinomial choice models allow flexible substitution patterns but become computationally demanding with many alternatives or observations. With a fixed per-observation simulation budget, simulated maximum likelihood introduces simulation bias, while each optimization step requires a full-sample likelihood evaluation. We propose Stochastic Approximation with Unbiased Simulated Scores (SAUSS), an averaged stochastic approximation based on conditionally unbiased mini-batch score estimates. Each iteration uses a fixed mini-batch regardless of sample size. For multinomial probit, accept-reject sampling provides exact conditional draws and unbiased score estimates for any fixed number of accepted draws. Under local conditions, asymptotic theory for the averaged estimator and the partial-sum process of the SAUSS iterates incorporates mini-batch and simulation variability and supports random-scaling and plug-in inference. In simulations and an application, SAUSS gives comparable results in less than 1% of the computation time of simulated maximum likelihood. SAUSS extends to limited dependent variable models with conditional-expectation score representations and exact conditional sampling.

1 Introduction

Multinomial probit offers richer substitution patterns than multinomial logit but faces costly high-dimensional probability calculations and full-sample optimization. SAUSS addresses these computational and simulation issues with unbiased mini-batch scores, averaging, and inference that accounts for algorithmic variation.

  • Multinomial probit accommodates richer substitution patterns through a flexible covariance structure, but each likelihood contribution is a high-dimensional multivariate normal region probability.
  • As the number of alternatives or observations grows, higher-dimensional probabilities, more covariance parameters, and repeated full-sample evaluations increase computational burden.
  • Fixed-budget simulated likelihood can remain biased because taking logarithms of unbiased simulated probabilities generally does not preserve unbiasedness.
  • SAUSS embeds exact conditional accept–reject score simulation in mini-batch stochastic approximation with iterate averaging, allowing fixed-size batches and unbiased update directions.
  • GHK-based simulated maximum likelihood requires more than 100 times as much computation time as SAUSS under the reported benchmark and application configurations.
  • The theory separates mini-batch and score-simulation variation and supports random-scaling and plug-in inference for the averaged estimator.

2 The SAUSS Framework

The SAUSS framework replaces biased smooth simulated likelihood optimization with stochastic approximation driven by conditionally unbiased simulated score directions. Fixed mini-batches, accepted score draws, and iterate averaging yield an estimator whose uncertainty reflects both sampling mechanisms.

  • Motivation: Simulated maximum likelihood can be biased at a fixed simulation budget, while increasing that budget raises the cost of every mini-batch evaluation.
  • Unbiased simulated scores: MSS instead simulates likelihood scores directly, and exact conditional draws can make each simulated score unbiased for any fixed number of accepted draws.
  • SAUSS recursion: SAUSS uses the negative simulated score as a conditionally unbiased stochastic direction, without requiring the simulated direction to be the derivative of a smooth objective.
  • SAUSS recursion: Each iteration samples ordered mini-batch observations and fresh score-simulation inputs, with a decreasing learning rate and a recursion targeting the sample estimator conditionally on the data.
  • Iterate averaging: The reported estimator is a Polyak–Ruppert average of stored iterates after a chosen averaging start index.
  • Inference: For every fixed R ≥1, the ideal stochastic direction remains conditionally unbiased, although algorithmic uncertainty depends jointly on T, m, and R.

3 The LDV Score Identity

The LDV score identity rewrites a likelihood score as a conditional expectation of a complete-data score by representing each observed outcome as a parameter-independent latent-region event. Exact conditional draws then provide unbiased score contributions for SAUSS.

  • Score identity: With the region fixed, the parameter enters through the conditional latent density, so differentiating the log probability yields a conditional expectation of a complete-data score.
  • Fixed-region representation: The identity fixes the observed outcome as a latent-region event whose region depends on data but not on the parameter.
  • Regularity conditions: Under positivity, differentiability, and an integrable-envelope condition, the Hajivassiliou–McFadden identity justifies this score representation.
  • Simulation: Independent exact latent draws make the simulated score unbiased for any fixed R ≥1, although accept–reject implementations can make the score map discontinuous in the parameter.
  • Scope: The identity covers canonical LDV settings including panel probit, correlated-error multinomial choice, and censored Tobit contributions.

4 SAUSS for Multinomial Probit

The MNP implementation represents choice probabilities as multivariate normal probabilities over observed-choice regions and estimates their loss gradients using unbiased accept–reject score simulation. SAUSS then applies these estimates in an averaged mini-batch stochastic approximation recursion with normalized covariance parameters.

  • Model and choice regions: MNP models utility differences relative to a base alternative, with observed choices represented by fixed regions in differenced latent-utility space.The base alternative is a labeling convention; it does not restrict substitution patterns.
  • Likelihood and descent target: Choice probability Pij(θ) is the integral of a multivariate normal density over the region associated with alternative j.The density has mean ∆Xiβ and covariance Ω.
  • Covariance normalization: The covariance is parameterized as Ω = LL⊤ with L11 = 1, positive remaining diagonal elements, and unrestricted strictly lower-triangular elements.Log-diagonal coordinates enforce positivity, and d = 1 reduces to the binary probit normalization Ω = 1.
  • Likelihood and descent target: SAUSS targets the exact likelihood score directly, using the negative simulated score as a descent direction rather than maximizing the simulated likelihood objective.This avoids making the simulated probability estimate part of the optimized criterion.
  • Unbiased score simulation: Accept–reject sampling retains latent normal draws that fall in the observed-choice region, producing an unbiased score estimate for any fixed R ≥1 accepted draws.The simulator runs until the required accepted draws are obtained, so proposal-draw runtime is random.
  • SAUSS recursion: Each iteration samples a mini-batch of m observations, computes fresh accepted-draw gradient contributions, updates the parameter with a decreasing learning rate, and averages iterates.The computational budget is governed by T, m, and R; m may be much smaller than n.

5 Asymptotic Theory

Under local regularity, SAUSS's averaged iterates are consistent and asymptotically normal, while its functional limit accounts for both mini-batch and simulation variability. The theory also covers fixed mini-batches, preconditioning, burn-in, and smooth parameter transformations.

  • Theorem 1: Theorem 1 establishes consistency, a linear representation, and a functional central limit theorem for the SAUSS recursion with one-observation updates.The result is stated under local identification, smooth exact loss-gradients, unbiased simulated loss-gradients, moment conditions, and variance regularity.
  • Functional limit: The asymptotic process is driven by a p-dimensional standard Wiener process and incorporates the covariance of the stochastic update noise.The conclusions hold on the event that the recursion remains in the local neighborhood, with an explicit ϵ-error for its complement.
  • Scope and limitation: Ordinary weak convergence without the ϵ-error is unavailable under the local assumptions because unbounded simulation noise can drive the recursion outside the neighborhood.Removing the ϵ-formulation would require stability conditions outside the local region, which are generally unsuitable for MNP because the population objective can flatten as latent variances increase.
  • Mini-batches: Fixed mini-batches preserve the theorem with V1,R replaced by Vm,R, because averaging m conditionally independent contributions divides conditional variance by m.The theorem is proved for m = 1, while the extension applies to fixed m under independent sampling with replacement and conditionally independent simulation draws.
  • Preconditioning: A fixed symmetric positive-definite preconditioner preserves the Polyak–Ruppert covariance H−1Vm,RH−1 after transforming back to the original coordinates.The preconditioner changes the centered-coordinate drift and innovation covariance and adjusts the small-step threshold.
  • Transformations: Smooth transformations of θ0 obey the same functional limit as linear functionals of the original iterates, enabling direct path-based inference for scalar transformations.The result requires twice continuous differentiability of the transformation and avoids a delta-method variance estimator for the transformed path.

6 Inference

SAUSS inference accounts for both mini-batch and score-simulation variability in averaged iterates. It offers path-based random scaling and plug-in sandwich inference, with distinct interpretations under population versus fixed-sample sampling.

  • Inference consequences: SAUSS has two algorithmic variation sources—mini-batch selection and score simulation—whose magnitude relative to sampling uncertainty depends on n, T, m, and R.These sources are incorporated into the inference theory for the averaged estimator and iterate process.
  • Random scaling: Random-scaling inference uses the solution path and does not require estimates of H or the simulation variance.It can also be applied directly to a smooth scalar transformation tracked along the recursion.
  • Finite-run variance: The finite-run variance approximation adds algorithmic uncertainty from terminating after T iterations to sampling uncertainty in the exact sample estimator.It also reflects reductions from larger mini-batches and additional accepted draws.
  • Sandwich inference: Sandwich inference estimates the first-order algorithmic covariance when H and V_m,R are consistently estimated.It is closer to conventional MLE inference when computational uncertainty is intended to be small.
  • Inference target: Under population sampling, path variation targets θ0; conditional on a fixed sample, it describes only randomized-algorithm variation around the sample target θ̂_n.The random-scaling self-normalizer is not a consistent estimator of the algorithmic covariance matrix.

7 Monte Carlo Experiments

Monte Carlo experiments assess SAUSS in multinomial probit models as the number of alternatives grows and against GHK simulated maximum likelihood. SAUSS remains computationally feasible and has similar small-choice-set accuracy under the reported configurations.

  • Design and implementation: The baseline design fixes n = 2000, uses 1000 replications, and varies the number of alternatives over J ∈ {4, 8, 16, 32, 64}.The covariance dimension ranges from 3 × 3 to 63 × 63.
  • Design and implementation: SAUSS uses mini-batches of m = 10, R = 5 accepted draws per observation, and a trial cap of Kmax = 20,000.Each iteration samples observations with replacement and omits an observation’s contribution when no draw is accepted within the cap.
  • Scaling with alternatives: As J increases from 4 to 64, the coefficient criterion improves from 0.0205 to 0.0140, while acceptance rates decline from 0.281 to 0.033.At J = 64, median computation time is 80.5 seconds per replication and no-acceptance events remain below 0.02 percent.
  • Main results: At J = 4, SAUSS and GHK-SML differ by at most 0.0058 across the three accuracy criteria.Both estimators use the same 1,000 simulated data sets and normalization.
  • Main results: At J = 4, SAUSS has a median runtime of approximately 0.5 seconds versus 71.8 seconds for GHK-SML, a factor of about 140.The reported configurations produce similar aggregate accuracy in this small-choice-set setting.
  • Stress designs: Under stress designs with unequal covariance structure and a rare alternative, SAUSS remains effective despite heterogeneous acceptance rates.In D3, rare-alternative observations have acceptance rate 0.26 versus 0.66 for the most frequently chosen alternative.

8 Application: Maternal Labor-Supply Intentions

The application re-estimates a maternal labor-supply multinomial probit model with SAUSS and compares it with GHK simulated maximum likelihood. Estimates and significance classifications are broadly similar, while SAUSS is much faster under the reported implementations.

  • Application setup: The application analyzes stated choices among not working, part-time, and full-time employment for N = 2,873 respondents.The model includes five alternative-specific belief and norm measures and six demographic controls.
  • Application setup: SAUSS uses 1,600 epochs, mini-batches of m = 10, five accepted simulator draws per observation, averaging, and t^-0.501 step-size decay.Published estimates are rescaled to the common L11 = 1 normalization.
  • Estimates: Each SAUSS coefficient lies within one GHK standard error of its counterpart, with matching signs and 5% significance classifications.The SAUSS coefficients are somewhat larger, with fitted Ω22 = 0.58 versus 0.45 for GHK.
  • Uncertainty: The estimated algorithmic variance component is 0.16% of total variance, so the reported variance is dominated by sampling uncertainty.Across ten seeds, coefficient estimates vary by only three to five percent of one reported standard error.
  • Uncertainty: Using an outer product of simulated scores in place of H understates standard errors when R is fixed.The resulting standard errors are a median 0.56 times the simulation-adjusted values at R = 5 and 0.96 times them at R = 200.
  • Computation: The reported fit takes 33 seconds with SAUSS versus 4,180 seconds with the GHK likelihood.At 200 GHK draws, the GHK fit takes 163 seconds and changes by at most 0.16 standard errors.

9 Conclusion

The paper presents SAUSS as mini-batch stochastic approximation built around unbiased simulated scores for likelihood-based estimation. Results show comparable accuracy to GHK-based simulated maximum likelihood with substantially lower computation time, while retaining extensions to suitable limited dependent variable models.

  • Conclusion: SAUSS combines the method-of-simulated-scores identity with mini-batch stochastic approximation to handle discontinuous but unbiased score simulation.The MNP implementation uses differencing relative to alternative 1 and L11 = 1 normalization.
  • Conclusion: SAUSS achieves accuracy comparable to GHK-based simulated maximum likelihood in the small-choice-set benchmark and remains computationally feasible as alternatives increase.The application likewise produces broadly similar estimates with substantially lower computation time.

A.1 Setup and notation

The proofs specialize the stochastic approximation setup to mini-batch size m = 1, defining observation and simulation conditioning, drift, covariance, and asymptotic notation. They also specify the assumptions and run-length asymptotics used throughout.

  • The proof setup takes mini-batch size m = 1 and uses observations drawn independently across iterations by streaming or with-replacement sampling.
  • At iteration t, Dt is the drawn observation and Y ∗t,r(θ) is the r-th accepted latent draw from the exact simulator.
  • Gt(θ) averages one-draw loss gradients over simulation draws, while R(θ) is its conditional expectation and H(θ) = ∇θR(θ), with R(θ0) = 0.
  • The conditional simulation covariance is At(θ, Dt), and S = V1,R denotes the limiting covariance from Assumption 7.
  • Asymptotic statements let T and n0(T) diverge without a relative-rate restriction, while p, R, γ0, and a remain fixed.

A.2 Proofs of the main results

The main proofs decompose SAUSS errors into observation noise, simulation noise, and nonlinear drift remainders, then combine martingale limits with high-probability localization. This yields weak convergence, consistency, and inference results, including fixed-m extensions.

  • The innovation decomposes into simulation noise ζt and observation-sampling noise ξt, while ηt captures the nonlinear remainder from linearizing the mean drift.
  • The proofs compare actual iterates with continued innovations on the event E that the process remains in the local region, controlling the complement through P(Ec) ≤ ϵ + o(1).
  • Theorem 1 follows by applying a martingale functional central limit theorem to the innovation process and transferring the limit from continued to actual innovations.
  • Theorem 2 uses local Taylor expansion, a stopped remainder bound, and Lemmas 3–5 to establish weak convergence for the averaged estimator in the same ϵ-sense.
  • Corollary 1 combines Theorem 1 with the fixed-m extension in Remark 5, yielding positive-definite covariance and self-normalized or sandwich convergence results.

A.3 Supporting lemmas

The supporting lemmas establish local stability, high-probability containment, linear representations, and martingale functional limits for decreasing step sizes. They also identify the scope of the resulting consistency statement.

  • Lemma 2 establishes local drift and moment bounds, then shows the iterates remain in Bρ with probability at least 1 − ϵ − o(1) under a sufficiently small γ0.
  • The small-step requirement can alternatively be achieved by shifting the learning-rate sequence with a sufficiently large t0.
  • For likelihood problems, H is the population information matrix and is symmetric; the stability arguments otherwise require only the stated drift condition.
  • The recursion separates linear drift, martingale innovations, and a remainder bounded by C∥∆t−1∥2 within Bρ, enabling finite-horizon transition-matrix representations.
  • With γt = γ0(t − 1)^−a and a ∈ (1/2, 1), the lemmas control transient, early, and recent terms and establish stopped consistency and approximation remainders.
  • Fresh innovations at θ0 form a martingale-difference sequence whose functional central limit theorem leads, after premultiplication by −H−1, to the stated process limit.
  • The resulting consistency is only an ϵ-form statement and does not imply unconditional oP(1) for a fixed step-size scale.
Loading 2608.25304v1…