Algorithms

From subsamples
to one distance ordering.

Exact computation, shared work across scales, and the choices behind regression inference.

Compute the exact weighted average

  1. DistanceLocate neighbors
  2. OrderRank responses
  3. WeightsApply the scale
  4. SumReturn the estimate
Inputs: training X and y, query x, integer scale s. Output: μ̂s(x).
StepCalculationReuse or cost
FitStore training data and construct w(n,s).O(n) weight calculation, cached once per scale
DistanceCompute squared Euclidean distances to x.O(nd), shared across scales at this query
OrderStable-sort training indices by distance.O(n log n), ties retain original row order
PredictTake the weighted sum of ordered responses.O(n), no enumeration of subsets
w1 = s/n,
wi+1 = wi(n−i−s+1)/(n−i),   i = 1,…,n−s.

The recurrence evaluates the binomial ratio without constructing large binomial coefficients. The maintained implementation keeps the exact nonzero support, rather than cutting off small weights arbitrarily. One query costs O(nd + n log n).

Squared distances preserve Euclidean ordering. Direct coordinate differences avoid a cancellation problem that can occur when distances are formed by subtracting large squared norms.

One ordering serves both scales

Cache the rank weights for s₁ and s₂, or their combined vector. Use the same distance order for the two base estimates.

wcombinedi = a1wi(n,s1) + a2wi(n,s2),
a1 = −1 / [(s2/s1)2/d−1],   a2 = 1−a1.

The coefficients use the actual integer scales and distance dimension. Rounding a proposed second scale changes the ratio; recompute the coefficients for that ratio. A stable logarithmic/exponential formulation helps when ratios are close.

Bias-cancellation conditions · NumPy implementation.

Reuse work in resampling

NumPy inference with fixed scales
RoutineStatistical operationComputational operation
BootstrapDraw n paired observations with replacement. Recompute the estimator for each draw.Reuse distance ordering through counts and cumulative weights. Exact ties use a stable-sort fallback.
JackknifeDelete each observation and apply the paper's variance formula.Prefix/suffix sums avoid n refits and repeated sorts.
Multiple queriesUse the same resampled observation indices across queries.Preserve pairing so cross-query relationships are retained in the draws.

The base estimator's subsampling scheme and the inference bootstrap are distinct. BNN uses all fixed-size subsets without replacement; the bootstrap resamples the training dataset with replacement and then evaluates BNN on that resample.

The bootstrap sample variance uses denominator B−1. The jackknife is centered on the original estimate as in the paper's formula, and scales must be admissible for the n−1 sample.

Inference scope · Implementation notes.

Scale selection defines the procedure

The package's held-out helper minimizes global validation MSE over supplied pairs. Each validation query's distance ordering is shared across pairs, making the search cheaper than fitting and sorting anew for every candidate.

The paper replication includes its pointwise tuning procedures. Specify which rule is used when reporting results. Bootstrap inference with fixed selected scales does not automatically include the randomness of the selection step.

Read timings as measured workloads

The current package measurements concern 64 regression queries with n = 20,000, d = 4, scales 20/40, and float64. Fit and transfer are excluded. The records use warmed medians of five runs; CUDA measurement boundaries synchronize, and JAX completion is awaited.

Recorded regression prediction times
Path64 queriesContext
NumPy CPU0.0867 sIndependent reference formulation; separate generated sample
PyTorch CPU0.0295 sNumPy agreement and original scalar-call comparison
PyTorch CUDA · RTX 30900.0024 sSynchronized boundaries; NumPy and original comparison
JAX CPU0.0379 sWarmed execution; first call includes compilation

NumPy's generated sample differs from the backend-comparison sample. These records document measured configurations rather than a controlled ranking of all backends. JAX's first call was 0.132 s including compilation. Environment and raw records.

Saved paper computation table

The original paper's timing workload includes derivative-based estimation and first-stage calculations, rather than only the regression queries above. Scales are 7/14. Each cell gives point-estimate time / 100-bootstrap time in seconds.

Historical reported workload
nNumba CPUPyTorch CPUPyTorch GPU
5,0000.0066 / 0.04680.0404 / 1.15230.0032 / 0.2872
20,0000.0234 / 0.08140.0361 / 5.27780.0139 / 0.7940
80,0000.1050 / 0.18670.1290 / 15.95750.0305 / 2.5163
320,0000.3619 / 1.02410.3108 / 29.42740.0511 / 3.5945
1,280,0001.4909 / 5.71441.0065 / 116.95380.0817 / 6.3980

Retained from the saved source table. Original GPU timing code lacks explicit synchronization at measurement boundaries. This is historical reported evidence, separate from the synchronized package measurements.

Batching and hardware parallelism help when enough queries are available. For small workloads, device transfer and compilation can dominate. Warmed kernel time and end-to-end latency answer different practical questions.

The metric and workspace are part of computation

Euclidean distance on differently scaled features changes the neighborhood. Choose transformations and scaling intentionally, and apply them consistently. The dimension affects locality and the bias-correction coefficients.

Accelerator prediction works in bounded query batches. The direct-difference distance workspace scales with batch size, n, and d. Reusing fitted weights removes repeated construction but does not remove distance calculation or sorting for a new query.

Exact ordering is retained in the maintained implementation. A future approximate neighbor search would need its own accuracy analysis. Run the current implementation.