A2-MN Calibrated-Moment Newton Solver

Seeks Gamma(a, b) hyperprior parameters that match target moments for the number of clusters K_J under a Dirichlet process prior, and reports convergence only after an independent higher-order check.

Usage

DPprior_a2_newton(
  J,
  mu_K,
  var_K,
  a0 = NULL,
  b0 = NULL,
  tol_F = .TOL_NEWTON,
  tol_rel = 1e-08,
  tol_step = 1e-10,
  max_iter = 20L,
  damping = TRUE,
  use_fallback = TRUE,
  M = .QUAD_NODES_DEFAULT,
  M_verify = NULL,
  verification_abs_tol = 1e-10,
  verification_rel_tol = 1e-08,
  verbose = FALSE
)

Arguments

J: Integer; sample size (number of observations/sites). Must be >= 2.
mu_K: Numeric; target prior mean \(E[K_J]\). Must satisfy \(1 < \mu_K < J\).
var_K: Numeric; target prior variance \(\mathrm{Var}(K_J)\). Must be positive.
a0: Numeric or NULL; initial shape parameter. If NULL, computed via DPprior_a1.
b0: Numeric or NULL; initial rate parameter. If NULL, computed via DPprior_a1.
tol_F: Positive numeric; absolute part of the componentwise moment tolerance. Retained under its historical name for compatibility.
tol_rel: Non-negative numeric; relative part of the componentwise tolerance. Each component uses tol_F + tol_rel * max(abs(target), 1).
tol_step: Numeric; stopping tolerance for Newton step size. Default: 1e-10.
max_iter: Integer; maximum Newton iterations. Default: 20.
damping: Logical; if TRUE, use backtracking line search for damped Newton updates. Default: TRUE.
use_fallback: Logical; if TRUE, use Nelder-Mead fallback when Newton fails to converge. Default: TRUE.
M: Integer; number of quadrature nodes for moment computation. Default: 80.
M_verify: Optional integer; independent verification order. It must be at least max(2*M, M+40) and no greater than 512. The required order is selected automatically when available. If the required order exceeds 512, the function fails closed with a typed no-candidate condition. The numerical source candidate is then available only in the condition's non-authoritative compatibility view, with guidance to refit using M <= 256.
verification_abs_tol, verification_rel_tol: Non-negative absolute and relative tolerances for selected-versus-verification order stability. At least one must be positive.
verbose: Logical; if TRUE, print iteration progress. Default: FALSE.

Value

A canonical dpprior.result/1 DPprior_fit object with authoritative nested components:

parameters: Gamma shape, rate, and parameterization.
J: Integer; sample size
target: Canonical immutable K target.
achieved,residuals: Selected-order public moment evidence.
method: Character; "A2-MN" or "A2-MN+NM" if fallback was used
status: One of "converged", "boundary", "approximate", or "failed".
usable,verified: Logical status fields derived from the numerical and independent-verification gates.
computation: Orders, fixed scaling, typed attempts, candidate ledger, fallback, termination, and trace.
verification: Independent higher-order recomputation, componentwise adequacy, and order-stability evidence.
provenance,compatibility: Canonical provenance plus a quarantined non-authoritative legacy view.

Finite results also retain exact compatibility aliases a, b, converged, iterations, termination, fit, diagnostics, trace, and attempts. If no independent verifier order is available, the function signals a dpprior_a2_no_candidate condition containing a valid parameterless failed result in its result field.

Details

This implements TSMM Stage 2 (A2-MN) from Lee (2026). The A2-MN algorithm uses Newton's method in log-scale to ensure positivity of the Gamma parameters. The Jacobian uses score-function identities and their order-refinement diagnostics. Residual components are scaled by max(abs(target), 1) and checked individually under an explicit absolute-plus-relative budget; a small aggregate norm cannot hide a failed component.

Algorithm Steps:

Initialize: \((a_0, b_0)\) from A1 closed-form or user-provided
Log-parameterize: \(\eta = (\log a, \log b)\)
For each iteration:
- Compute moments \((M_1, V)\) and Jacobian \(J_F\)
- Compute natural-unit and scaled component residuals
- Transform Jacobian to log-scale: \(J_{\log} = J_F \cdot \text{diag}(a, b)\)
- Newton step: \(\Delta = -J_{\log}^{-1} F\)
- Backtracking line search (if damping enabled)
- Update: \(\eta \leftarrow \eta + \lambda \Delta\)
If needed, run the predeclared Nelder–Mead fallback on \((\log a,\log b)\), never on raw positive parameters
Independently recompute the selected candidate at a distinct higher quadrature order; retain the fit-order candidate unchanged

Termination and status:

Optimizer exit, a small step, or line-search exhaustion is never by itself convergence.
converged requires the selected- and higher-order residual gates plus the quadrature-stability gate.
approximate retains a finite but unusable candidate when a solver or available independent-verification gate is unmet.
An unavailable required verifier order produces a typed failed/no-candidate result; it never exposes the numerical source candidate as public scientific output.
Singular/ill-conditioned Jacobians, line-search failure, stagnation, fallback outcome, and verification failure have stable reason codes in computation and the compatibility audit view.

A1 projection may be requested explicitly to obtain an initialization only. Its projection record is retained in the non-authoritative initialization audit view, while the canonical requested A2 target is immutable and is never silently projected.

References

Lee, J. (2026). Design-Conditional Prior Elicitation for Dirichlet Process Mixtures. arXiv preprint arXiv:2602.06301.

Examples

# Basic usage
fit <- DPprior_a2_newton(J = 50, mu_K = 5, var_K = 8)
print(fit)
#> DPprior Prior Elicitation Result
#> ============================================= 
#> 
#> Schema: dpprior.result/1
#> Method: A2-MN (mode: a2_moment)
#> Status: converged; usable: yes; verified: yes
#> 
#> Target (J = 50):
#>   E[K_J]   = 5.0000
#>   Var(K_J) = 8.0000
#> 
#> Canonical candidate:
#>   alpha ~ Gamma(a = 2.0361, b = 1.6051)
#>   Achieved E[K_J] = 5.000000; Var(K_J) = 8.000000
#>   Maximum absolute moment residual = 1.89e-12
#> 
#> Canonical guidance: scaled component residuals and independent higher-order verification passed

# Inspect the selected-order match and independent verification
achieved <- exact_K_moments(
  50, fit$parameters$a, fit$parameters$b
)
cat(sprintf("Target E[K]=5, Achieved E[K]=%.10f\n", achieved$mean))
#> Target E[K]=5, Achieved E[K]=5.0000000000
cat(sprintf("Target Var=8, Achieved Var=%.10f\n", achieved$var))
#> Target Var=8, Achieved Var=8.0000000000

# Compare A1 vs A2 accuracy
a1 <- DPprior_a1(J = 50, mu_K = 5, var_K = 8)
a1_mom <- exact_K_moments(
  50, a1$parameters$a, a1$parameters$b
)
a2_mom <- exact_K_moments(
  50, fit$parameters$a, fit$parameters$b
)
cat(sprintf("A1 mean error: %.6f\n", abs(a1_mom$mean - 5)))
#> A1 mean error: 0.538649
cat(sprintf("A2 mean error: %.2e\n", abs(a2_mom$mean - 5)))
#> A2 mean error: 1.16e-12

# View iteration trace (includes step size and Jacobian determinant)
head(fit$computation$trace)
#>   iter        a         b       M1         V     residual residual_mean
#> 1    1 4.000000 3.9120230 4.461351  4.783136 3.261649e+00 -5.386492e-01
#> 2    2 1.178648 0.9119679 4.909045 10.854547 2.855995e+00 -9.095498e-02
#> 3    3 1.849987 1.4593445 4.976321  8.388598 3.893184e-01 -2.367904e-02
#> 4    4 2.029715 1.6000585 4.999257  8.012326 1.234807e-02 -7.433076e-04
#> 5    5 2.036085 1.6050485 4.999999  8.000014 1.381980e-05 -8.210875e-07
#> 6    6 2.036093 1.6050541 5.000000  8.000000 2.222027e-12  1.162626e-12
#>   residual_variance   scaled_mean scaled_variance standardized_norm
#> 1     -3.216864e+00 -1.077298e-01   -4.021080e-01      2.605909e+07
#> 2      2.854547e+00 -1.819100e-02    3.568183e-01      2.245304e+07
#> 3      3.885976e-01 -4.735808e-03    4.857470e-02      3.065838e+06
#> 4      1.232568e-02 -1.486615e-04    1.540710e-03      9.723509e+04
#> 5      1.379539e-05 -1.642175e-07    1.724424e-06      1.088179e+02
#> 6      1.893596e-12  2.325251e-13    2.366995e-13      2.022562e-05
#>   max_budget_ratio step    step_norm line_search_iterations accepted   det_Jlog
#> 1     3.574293e+07    1 1.900957e+00                      1     TRUE  -5.300961
#> 2     3.171718e+07    1 6.513523e-01                      1     TRUE -21.408588
#> 3     4.317751e+06    1 1.306524e-01                      1     TRUE -15.275642
#> 4     1.369520e+05    1 4.417739e-03                      1     TRUE -14.290104
#> 5     1.532821e+02    1 4.948559e-06                      1     TRUE -14.257384
#> 6     2.103996e-05   NA           NA                      0    FALSE -14.257347
#>   reciprocal_condition  jacobian_status derivative_status
#> 1            0.1193782 well_conditioned         converged
#> 2            0.1404766 well_conditioned         converged
#> 3            0.1554857 well_conditioned         converged
#> 4            0.1562930 well_conditioned         converged
#> 5            0.1563039 well_conditioned         converged
#> 6            0.1563039 well_conditioned         converged
#>              reason_code
#> 1          step_accepted
#> 2          step_accepted
#> 3          step_accepted
#> 4          step_accepted
#> 5          step_accepted
#> 6 residual_tolerance_met