Binding Internals: MEstimator

Python callbacks, Rust optimization, and fit-time inference

MEstimator is the callback-driven alternative to a built-in likelihood. Python defines the objective and observation-level scores; Rust owns L-BFGS, convergence checks, covariance construction, and bootstrap orchestration. The class reference gives the full signature. The implementation lives in src/estimators/mle.rs.

Callback contracts

• objective_fn(theta, data) returns a finite scalar to minimize and its one-dimensional NumPy gradient, of length len(theta).
• score_fn(theta, data) returns a finite (n, p) NumPy array of individual estimating-function contributions. It must have a stable shape at nearby parameter values.
• Neither callback should mutate its inputs. Both must refer to the same statistical model and support the resampling indices when bootstrap is used.

Cost and gradient requests cross the Python boundary separately, so expensive callbacks can dominate runtime. Built-in estimators avoid this callback overhead.

Fit and inference

fit(data, theta0) clears prior state, validates inputs, and runs L-BFGS with a More-Thuente line search. Exhausting the iteration budget raises an error. Before a successful fit returns, Rust evaluates the scores and constructs

\[ A=\left.\frac{\partial\,n^{-1}\sum_i\psi_i(\theta)} {\partial\theta'}\right|_{\hat\theta},\qquad B=n^{-1}\sum_i\psi_i(\hat\theta)\psi_i(\hat\theta)',\qquad \widehat V=\frac1n A^{-1}BA^{-T}. \]

The bread \(A\) is the central-difference derivative of the mean scores, not a second copy of \(B\). The parameter-specific finite-difference step is derivative_step * max(abs(theta[j]), 1). A singular bread or invalid callback output causes the fit to fail and clears the partial result.

summary() reuses the fit-time covariance. Mutating callback data afterward does not silently change that summary. compute_vcov() explicitly recomputes it using the retained data; bootstrap also reuses those data, so keep them immutable for the fitted object’s lifetime.

A complete callback example

This least-squares objective includes its own intercept column. MEstimator does not add one or reorder the parameters.

import numpy as np
import crabbymetrics as cm

rng = np.random.default_rng(91)
x = np.column_stack([np.ones(300), rng.normal(size=(300, 2))])
y = x @ np.array([0.4, 1.0, -0.5]) + rng.normal(size=300)
data = {"X": x, "y": y, "n": len(y)}

def rows(data):
    idx = data.get("indices", np.arange(data["n"]))
    return data["X"][idx], data["y"][idx]

def objective(theta, data):
    design, outcome = rows(data)
    residual = design @ theta - outcome
    return 0.5 * np.mean(residual**2), design.T @ residual / len(outcome)

def scores(theta, data):
    design, outcome = rows(data)
    return design * (design @ theta - outcome)[:, None]

model = cm.MEstimator(objective, scores, max_iterations=300, tolerance=1e-8)
model.fit(data, np.zeros(x.shape[1]))
result = model.summary()
print("Coefficients:", result["coef"])
print("Sandwich SE:", result["se"])
print("Converged:", result["converged"])
draws = model.bootstrap(20, seed=91)
print("Bootstrap draw shape:", draws.shape)
np.testing.assert_allclose(result["coef"], np.linalg.lstsq(x, y, rcond=None)[0], atol=1e-5)

Coefficients: [ 0.48858993  0.96647527 -0.54561321]
Sandwich SE: [0.05351078 0.05131308 0.05336422]
Converged: True
Bootstrap draw shape: (20, 3)

The bootstrap requires a dictionary containing n. Each draw shallow-copies that dictionary, supplies fresh row indices, and refits the objective. The example’s 20 draws demonstrate the interface, not a precise tail interval. Use GMM when the starting point is a vector of moment conditions and you need explicit weighting or dependence-robust inference.