// std/math — Extended math utilities
//
// Import with: import "std/math"
fn _is_number(value) {
let t = type_of(value)
return t == "int" || t == "float"
}
fn _require_number(value, label) {
require _is_number(value), "${label} must be numeric"
return value
}
fn _require_numeric_list(items, label = "items") {
require type_of(items) == "list", "${label} must be a list"
for item in items {
require _is_number(item), "${label} must contain only numbers"
}
return items
}
fn _require_non_empty_numeric_list(items, label = "items") {
_require_numeric_list(items, label)
require len(items) > 0, "${label} must not be empty"
return items
}
fn _require_probability(prob, label = "probability") {
_require_number(prob, label)
require prob > 0 && prob < 1, "${label} must be between 0 and 1"
return prob
}
fn _require_positive_number(value, label) {
_require_number(value, label)
require value > 0, "${label} must be greater than 0"
return value
}
fn _require_positive_int(value, label) {
_require_number(value, label)
require floor(value) == value, "${label} must be an integer"
require value > 0, "${label} must be at least 1"
return value
}
fn _require_non_negative_int(value, label) {
_require_number(value, label)
require floor(value) == value, "${label} must be an integer"
require value >= 0, "${label} must be non-negative"
return value
}
fn _require_vector(vec, label = "vector") {
_require_numeric_list(vec, label)
return vec
}
fn _require_same_vector_dimensions(a, b) {
_require_vector(a, "a")
_require_vector(b, "b")
require len(a) == len(b), "vectors must have the same length"
}
fn _sum_squared_deltas(a, b) {
_require_same_vector_dimensions(a, b)
var total = 0.0
var idx = 0
while idx < len(a) {
let delta = a[idx] - b[idx]
total = total + delta * delta
idx = idx + 1
}
return total
}
fn _require_points(points) {
require type_of(points) == "list", "points must be a list"
require len(points) > 0, "points must not be empty"
let dims = len(_require_vector(points[0], "point"))
require dims > 0, "points must have at least one dimension"
for point in points {
_require_vector(point, "point")
require len(point) == dims, "all points must have the same dimension"
}
return dims
}
fn _copy_list(items) {
var out = []
for item in items {
out = out.push(item)
}
return out
}
fn _weight_list(weights, label = "weights") {
let ws = _require_non_empty_numeric_list(weights, label)
for weight in ws {
require weight >= 0, "${label} must contain only non-negative numbers"
}
return ws
}
fn _require_parallel_numeric_lists(xs, ys, x_label = "xs", y_label = "ys") {
_require_non_empty_numeric_list(xs, x_label)
_require_non_empty_numeric_list(ys, y_label)
require len(xs) == len(ys), "${x_label} and ${y_label} must have the same length"
}
fn _score_item(item, score_fn = nil) {
if score_fn == nil {
return item
}
return score_fn(item)
}
fn _insert_ranked(ranked, entry, descending = false) {
var out = []
var inserted = false
for existing in ranked {
if !inserted {
let better = if descending {
entry.score > existing.score
} else {
entry.score < existing.score
}
if better {
out = out.push(entry)
inserted = true
}
}
out = out.push(existing)
}
if !inserted {
out = out.push(entry)
}
return out
}
fn _rank_items(items, score_fn = nil, descending = false) {
require type_of(items) == "list", "items must be a list"
var ranked = []
var idx = 0
while idx < len(items) {
let entry = {
index: idx,
item: items[idx],
score: _score_item(items[idx], score_fn)
}
ranked = _insert_ranked(ranked, entry, descending)
idx = idx + 1
}
return ranked
}
fn _zero_vector(size) {
var out = []
var idx = 0
while idx < size {
out = out.push(0.0)
idx = idx + 1
}
return out
}
fn _nearest_centroid(point, centroids) {
var best_index = 0
var best_distance = _sum_squared_deltas(point, centroids[0])
var idx = 1
while idx < len(centroids) {
let distance = _sum_squared_deltas(point, centroids[idx])
if distance < best_distance {
best_distance = distance
best_index = idx
}
idx = idx + 1
}
return {index: best_index, distance: best_distance}
}
fn _farthest_point(points, centroids) {
var best_index = 0
var best_distance = _nearest_centroid(points[0], centroids).distance
var idx = 1
while idx < len(points) {
let distance = _nearest_centroid(points[idx], centroids).distance
if distance > best_distance {
best_distance = distance
best_index = idx
}
idx = idx + 1
}
return best_index
}
fn _init_kmeans_centroids(points, k) {
var centroids = [_copy_list(points[0])]
while len(centroids) < k {
let next_index = _farthest_point(points, centroids)
centroids = centroids.push(_copy_list(points[next_index]))
}
return centroids
}
fn _assign_points(points, centroids) {
var assignments = []
var errors = []
for point in points {
let nearest = _nearest_centroid(point, centroids)
assignments = assignments.push(nearest.index)
errors = errors.push(nearest.distance)
}
return {assignments: assignments, errors: errors}
}
fn _largest_error_point(errors) {
var best_index = 0
var best_error = errors[0]
var idx = 1
while idx < len(errors) {
if errors[idx] > best_error {
best_error = errors[idx]
best_index = idx
}
idx = idx + 1
}
return best_index
}
fn _recompute_centroids(points, assignments, k, dims, errors) {
var sums = []
var counts = []
var idx = 0
while idx < k {
sums = sums.push(_zero_vector(dims))
counts = counts.push(0)
idx = idx + 1
}
idx = 0
while idx < len(points) {
let cluster = assignments[idx]
counts[cluster] = counts[cluster] + 1
var dim = 0
while dim < dims {
var row = sums[cluster]
row[dim] = row[dim] + points[idx][dim]
sums[cluster] = row
dim = dim + 1
}
idx = idx + 1
}
var centroids = []
idx = 0
while idx < k {
if counts[idx] == 0 {
let reseed = _largest_error_point(errors)
centroids = centroids.push(_copy_list(points[reseed]))
} else {
var centroid = []
var dim = 0
while dim < dims {
centroid = centroid.push(sums[idx][dim] / counts[idx])
dim = dim + 1
}
centroids = centroids.push(centroid)
}
idx = idx + 1
}
return {centroids: centroids, counts: counts}
}
/** Clamp a value between min and max. */
fn clamp(value, lo, hi) {
if value < lo { return lo }
if value > hi { return hi }
return value
}
/** Linear interpolation between a and b by t (0..1). */
fn lerp(a, b, t) {
return a + (b - a) * t
}
/** Map a value from one range to another. */
fn map_range(value, in_lo, in_hi, out_lo, out_hi) {
let t = (value - in_lo) / (in_hi - in_lo)
return lerp(out_lo, out_hi, t)
}
/** Convert degrees to radians. */
fn deg_to_rad(degrees) {
return degrees * pi / 180
}
/** Convert radians to degrees. */
fn rad_to_deg(radians) {
return radians * 180 / pi
}
/** Sum a list of numbers. */
fn sum(items) {
_require_numeric_list(items)
return items.reduce(0, { acc, x -> acc + x })
}
/** Average of a list of numbers. */
fn avg(items) {
_require_numeric_list(items)
if len(items) == 0 { return 0 }
return sum(items) / (len(items) * 1.0)
}
/** Alias for avg(items). */
fn mean(items) {
return avg(items)
}
/** Median of a list of numbers. */
fn median(items) {
let sorted = _require_non_empty_numeric_list(items).sort()
let n = len(sorted)
let mid = floor(n / 2)
if n % 2 == 1 {
return sorted[mid]
}
return (sorted[mid - 1] + sorted[mid]) / 2.0
}
/** R-7 percentile interpolation with p in [0, 100]. */
fn percentile(items, p) {
let sorted = _require_non_empty_numeric_list(items).sort()
_require_number(p, "percentile")
require p >= 0 && p <= 100, "percentile must be between 0 and 100"
if len(sorted) == 1 { return sorted[0] }
let n = len(sorted)
let h = 1 + (n - 1) * (p / 100.0)
let lower = floor(h)
let upper = ceil(h)
if lower == upper {
return sorted[lower - 1]
}
let weight = h - lower
return sorted[lower - 1] + weight * (sorted[upper - 1] - sorted[lower - 1])
}
/** Return indices that would sort items ascending. */
fn argsort(items, score_fn = nil) {
let ranked = _rank_items(items, score_fn)
return ranked.map({ entry -> entry.index })
}
/** Return the top k items by score, descending. */
fn top_k(items, k, score_fn = nil) {
let ranked = _rank_items(items, score_fn, true)
let limit = _require_non_negative_int(k, "k")
return ranked.take(limit).map({ entry -> entry.item })
}
/** Population or sample variance. */
fn variance(items, sample = false) {
let xs = _require_non_empty_numeric_list(items)
if sample {
require len(xs) >= 2, "sample variance requires at least 2 values"
}
let mu = mean(xs)
var total = 0.0
for x in xs {
let delta = x - mu
total = total + delta * delta
}
let denom = if sample { len(xs) - 1 } else { len(xs) }
return total / denom
}
/** Population or sample standard deviation. */
fn stddev(items, sample = false) {
return sqrt(variance(items, sample))
}
/** Scale a numeric list to [0, 1]. */
fn minmax_scale(items) {
let xs = _require_non_empty_numeric_list(items)
let lo = xs.min()
let hi = xs.max()
if hi == lo {
return xs.map({ _x -> 0.0 })
}
let span = (hi - lo) * 1.0
return xs.map({ x -> (x - lo) / span })
}
/** Standardize a numeric list to zero mean and unit variance. */
fn zscore(items, sample = false) {
let xs = _require_non_empty_numeric_list(items)
let mu = mean(xs)
let sigma = stddev(xs, sample)
if sigma == 0 {
return xs.map({ _x -> 0.0 })
}
return xs.map({ x -> (x - mu) / sigma })
}
/** Weighted arithmetic mean. */
fn weighted_mean(items, weights) {
let xs = _require_non_empty_numeric_list(items)
let ws = _weight_list(weights)
require len(xs) == len(ws), "items and weights must have the same length"
let total_weight = sum(ws)
require total_weight > 0, "weights must sum to a positive value"
var total = 0.0
var idx = 0
while idx < len(xs) {
total = total + xs[idx] * ws[idx]
idx = idx + 1
}
return total / total_weight
}
/** Draw one item according to non-negative weights. */
fn weighted_choice(items, weights = nil) {
require type_of(items) == "list", "items must be a list"
require len(items) > 0, "items must not be empty"
let ws = if weights == nil {
items.map({ _item -> 1.0 })
} else {
_weight_list(weights)
}
require len(items) == len(ws), "items and weights must have the same length"
let total_weight = sum(ws)
require total_weight > 0, "weights must sum to a positive value"
let target = random() * total_weight
var cumulative = 0.0
var idx = 0
while idx < len(items) {
cumulative = cumulative + ws[idx]
if target < cumulative {
return items[idx]
}
idx = idx + 1
}
return items[len(items) - 1]
}
/** Convert scores into probabilities. */
fn softmax(items, temperature = 1.0) {
let xs = _require_non_empty_numeric_list(items)
_require_positive_number(temperature, "temperature")
let max_x = xs.max()
let exps = xs.map({ x -> exp((x - max_x) / temperature) })
let total = sum(exps)
require total > 0, "softmax total must be positive"
return exps.map({ x -> x / total })
}
/** Standard normal probability density. */
fn normal_pdf(x, mu = 0, sigma = 1) {
_require_number(x, "x")
_require_number(mu, "mean")
_require_positive_number(sigma, "stddev")
let z = (x - mu) / sigma
return exp(-0.5 * z * z) / (sigma * sqrt(2 * pi))
}
/** Standard normal cumulative distribution. */
fn normal_cdf(x, mu = 0, sigma = 1) {
_require_number(x, "x")
_require_number(mu, "mean")
_require_positive_number(sigma, "stddev")
let z = (x - mu) / sigma
let sign = if z < 0 { -1 } else { 1 }
let az = abs(z) / sqrt(2)
let t = 1 / (1 + 0.3275911 * az)
let a1 = 0.254829592
let a2 = -0.284496736
let a3 = 1.421413741
let a4 = -1.453152027
let a5 = 1.061405429
let erf = sign * (1 - (((((a5 * t + a4) * t + a3) * t + a2) * t + a1) * t * exp(-az * az)))
return 0.5 * (1 + erf)
}
/** Inverse normal CDF using Acklam's approximation. */
fn normal_quantile(prob, mu = 0, sigma = 1) {
let p = _require_probability(prob)
_require_number(mu, "mean")
_require_positive_number(sigma, "stddev")
let a1 = -39.69683028665376
let a2 = 220.9460984245205
let a3 = -275.9285104469687
let a4 = 138.357751867269
let a5 = -30.66479806614716
let a6 = 2.506628277459239
let b1 = -54.47609879822406
let b2 = 161.5858368580409
let b3 = -155.6989798598866
let b4 = 66.80131188771972
let b5 = -13.28068155288572
let c1 = -0.007784894002430293
let c2 = -0.3223964580411365
let c3 = -2.400758277161838
let c4 = -2.549732539343734
let c5 = 4.374664141464968
let c6 = 2.938163982698783
let d1 = 0.007784695709041462
let d2 = 0.3224671290700398
let d3 = 2.445134137142996
let d4 = 3.754408661907416
let plow = 0.02425
let phigh = 1 - plow
var x = 0.0
if p < plow {
let q = sqrt(-2 * ln(p))
let num = (((((c1 * q + c2) * q + c3) * q + c4) * q + c5) * q + c6)
let den = ((((d1 * q + d2) * q + d3) * q + d4) * q + 1)
x = num / den
} else if p <= phigh {
let q = p - 0.5
let r = q * q
let num = (((((a1 * r + a2) * r + a3) * r + a4) * r + a5) * r + a6) * q
let den = (((((b1 * r + b2) * r + b3) * r + b4) * r + b5) * r + 1)
x = num / den
} else {
let q = sqrt(-2 * ln(1 - p))
let num = (((((c1 * q + c2) * q + c3) * q + c4) * q + c5) * q + c6)
let den = ((((d1 * q + d2) * q + d3) * q + d4) * q + 1)
x = -(num / den)
}
return mu + sigma * x
}
/** Dot product of two vectors. */
fn dot(a, b) {
_require_same_vector_dimensions(a, b)
var total = 0.0
var idx = 0
while idx < len(a) {
total = total + a[idx] * b[idx]
idx = idx + 1
}
return total
}
/** Euclidean vector norm. */
fn vector_norm(v) {
_require_vector(v)
return sqrt(dot(v, v))
}
/** Normalize a vector to unit length. */
fn vector_normalize(v) {
_require_vector(v)
let norm = vector_norm(v)
require norm > 0, "cannot normalize a zero vector"
return v.map({ x -> x / norm })
}
/** Cosine similarity between two vectors. */
fn cosine_similarity(a, b) {
let an = vector_norm(a)
let bn = vector_norm(b)
require an > 0 && bn > 0, "cosine similarity requires non-zero vectors"
return dot(a, b) / (an * bn)
}
/** Euclidean distance between two vectors. */
fn euclidean_distance(a, b) {
return sqrt(_sum_squared_deltas(a, b))
}
/** Manhattan distance between two vectors. */
fn manhattan_distance(a, b) {
_require_same_vector_dimensions(a, b)
var total = 0.0
var idx = 0
while idx < len(a) {
total = total + abs(a[idx] - b[idx])
idx = idx + 1
}
return total
}
/** Chebyshev distance between two vectors. */
fn chebyshev_distance(a, b) {
_require_same_vector_dimensions(a, b)
var best = 0.0
var idx = 0
while idx < len(a) {
let delta = abs(a[idx] - b[idx])
if delta > best {
best = delta
}
idx = idx + 1
}
return best
}
/** Population or sample covariance between two numeric lists. */
fn covariance(xs, ys, sample = false) {
_require_parallel_numeric_lists(xs, ys)
if sample {
require len(xs) >= 2, "sample covariance requires at least 2 values"
}
let mu_x = mean(xs)
let mu_y = mean(ys)
var total = 0.0
var idx = 0
while idx < len(xs) {
total = total + (xs[idx] - mu_x) * (ys[idx] - mu_y)
idx = idx + 1
}
let denom = if sample { len(xs) - 1 } else { len(xs) }
return total / denom
}
/** Pearson correlation between two numeric lists. */
fn correlation(xs, ys, sample = false) {
_require_parallel_numeric_lists(xs, ys)
if sample {
require len(xs) >= 2, "sample correlation requires at least 2 values"
}
let sigma_x = stddev(xs, sample)
let sigma_y = stddev(ys, sample)
require sigma_x > 0 && sigma_y > 0, "correlation requires non-constant inputs"
return covariance(xs, ys, sample) / (sigma_x * sigma_y)
}
/** Sliding-window moving average. */
fn moving_avg(items, window) {
let xs = _require_non_empty_numeric_list(items)
let size = _require_positive_int(window, "window")
require size <= len(xs), "window must not exceed the number of items"
var out = []
var idx = 0
while idx + size <= len(xs) {
out = out.push(mean(xs.slice(idx, idx + size)))
idx = idx + 1
}
return out
}
/** Exponential moving average, returning one value per input. */
fn ema(items, alpha) {
let xs = _require_non_empty_numeric_list(items)
let a = _require_probability(alpha, "alpha")
var out = [xs[0] * 1.0]
var prev = xs[0] * 1.0
var idx = 1
while idx < len(xs) {
prev = a * xs[idx] + (1 - a) * prev
out = out.push(prev)
idx = idx + 1
}
return out
}
/** Deterministic k-means clustering. */
fn kmeans(points, k, options = nil) {
let dims = _require_points(points)
_require_non_negative_int(k, "k")
require k >= 1, "k must be at least 1"
require k <= len(points), "k must not exceed the number of points"
let max_iterations = if options != nil && options.max_iterations != nil {
options.max_iterations
} else {
100
}
_require_positive_int(max_iterations, "max_iterations")
var centroids = _init_kmeans_centroids(points, k)
var assignments = []
var counts = []
var iterations = 0
var converged = false
var inertia = 0.0
while iterations < max_iterations {
let assigned = _assign_points(points, centroids)
let next_assignments = assigned.assignments
inertia = sum(assigned.errors)
let updated = _recompute_centroids(points, next_assignments, k, dims, assigned.errors)
centroids = updated.centroids
counts = updated.counts
iterations = iterations + 1
if next_assignments == assignments {
assignments = next_assignments
converged = true
break
}
assignments = next_assignments
}
if !converged {
let assigned = _assign_points(points, centroids)
assignments = assigned.assignments
inertia = sum(assigned.errors)
var final_counts = []
var idx = 0
while idx < k {
final_counts = final_counts.push(0)
idx = idx + 1
}
for cluster in assignments {
final_counts[cluster] = final_counts[cluster] + 1
}
counts = final_counts
}
return {
centroids: centroids,
assignments: assignments,
counts: counts,
iterations: iterations,
converged: converged,
inertia: inertia
}
}