-- ----------------------------------------------------------------------
-- provsql 1.10.0 -> 1.11.0
--
-- New SQL surface since 1.10.0:
-- * Maintained provenance mappings. create_provenance_mapping(...,
-- maintained => true) registers the mapping in the new
-- provenance_mapping_registry; provenance_guard then appends each
-- genuine insert to it (keyed to the freshly minted input token), so
-- the mapping stays current AND survives the provsql rewrites that data
-- modification performs. This fixes the temporal validity reported for
-- a row after it is deleted/updated: a view-based mapping keyed validity
-- on the live (rewritten) provsql column and so lost the row's original
-- interval, while a maintained mapping keeps it keyed to the original
-- input token (the child a later monus/update gate wraps).
-- cleanup_table_info forgets a mapping when either table is dropped.
-- * create_provenance_mapping_view is removed (superseded by maintained
-- mapping tables). The base time_validity_view, a plain view over the
-- append-only update_provenance log, is unchanged by this upgrade.
-- * The continuous-distribution surface: the Gamma / Chi-squared /
-- Log-normal / Weibull / Pareto / Beta families, the discrete counts
-- (poisson / binomial / geometric / hypergeometric / negative_binomial
-- over categorical_from_log_pmf), rv_families(), the gmm and
-- empirical_samples / empirical_cdf constructors.
-- * Function application on random variables: ^ / pow / power / ln /
-- exp / sqrt, and the RV-typed CASE lowering targets provenance_case /
-- rv_case with the random_variable btree operator class that lets
-- GREATEST / LEAST parse (greatest / least order statistics, min / max
-- / rv_sum_or_null aggregates with per-aggregate identities).
-- * The aggregate-carrier CASE: agg_case, agg_gate_value /
-- agg_guard_holds (actual-world display values), agg_defined_event and
-- the conditional-on-defined agg_raw_moment (incl. the AVG arm over
-- agg_avg_moment_exact), and the extended agg_token_value_text.
-- * SQL-standard statistic aggregates over random_variable rows:
-- covar_pop / covar_samp / corr / stddev_pop / stddev_samp and the
-- ordered-set percentile_cont, with their rv_*_impl rewrite targets,
-- rv_stat_* circuit builders, and the rv_percentile_state transition
-- type.
-- * Readouts: probability (uuid alias + boolean-predicate overload with
-- the (A) | (B) conditioning operator), quantile / rv_quantile,
-- covariance / correlation / stddev, entropy / kl /
-- mutual_information with their rv_* bindings.
-- ----------------------------------------------------------------------
SET search_path TO provsql;
-- Guarded-selection gate for the RV / aggregate CASE lowering (rv_case /
-- agg_case below). ADD VALUE is upgrade-safe here: the value is only
-- referenced as a literal inside function bodies created by this script,
-- never materialised in this transaction (same pattern as 1.9.0 -> 1.10.0's
-- 'conditioned').
ALTER TYPE provenance_gate ADD VALUE IF NOT EXISTS 'case';
-- 'observe' gate: the likelihood-weighting observation leaf built by
-- observe(random_variable, datum). Only referenced at runtime (as the
-- create_gate type text), never in this script, so ADD VALUE is
-- transaction-safe on the PG12+ upgrade path exactly like 'case' above.
ALTER TYPE provenance_gate ADD VALUE IF NOT EXISTS 'observe';
-- Registry backing maintained mappings.
CREATE TABLE IF NOT EXISTS provsql.provenance_mapping_registry(
mapping oid PRIMARY KEY,
source oid NOT NULL,
attribute name NOT NULL
);
CREATE INDEX IF NOT EXISTS provenance_mapping_registry_source_idx
ON provsql.provenance_mapping_registry(source);
-- strip_annotations: identity-keyed consumers peel the transparent
-- annotation wrapper (dual of annotate).
CREATE OR REPLACE FUNCTION strip_annotations(token UUID) RETURNS UUID AS
$$
WITH RECURSIVE peel(g) AS (
SELECT token
UNION ALL
SELECT (provsql.get_children(p.g))[1] FROM peel p
WHERE provsql.get_gate_type(p.g) = 'annotation'
)
SELECT g FROM peel WHERE provsql.get_gate_type(g) <> 'annotation' LIMIT 1;
$$ LANGUAGE sql STABLE PARALLEL SAFE;
-- create_provenance_mapping gains the `maintained` argument: a signature
-- change, so drop the old four-argument form before recreating.
DROP FUNCTION IF EXISTS create_provenance_mapping(text, regclass, text, bool);
CREATE OR REPLACE FUNCTION create_provenance_mapping(
newtbl text,
oldtbl regclass,
att text,
preserve_case bool DEFAULT 'f',
maintained bool DEFAULT false
) RETURNS void AS
$$
DECLARE
BEGIN
-- Idempotence: when the mapping table already exists, leave it alone
-- with a NOTICE (re-runnable setup scripts / notebook cells). Drop it
-- first to rebuild a stale mapping.
IF (CASE WHEN preserve_case THEN to_regclass(format('%I', newtbl))
ELSE to_regclass(newtbl) END) IS NOT NULL THEN
RAISE NOTICE 'mapping table % already exists', newtbl;
RETURN;
END IF;
-- ON COMMIT DROP only fires at COMMIT: several mapping creations in
-- one transaction (a notebook cell, a setup script run via psql -1)
-- would otherwise collide on the leftover temp table. The to_regclass
-- probe (rather than DROP IF EXISTS) keeps the first call NOTICE-free.
IF to_regclass('pg_temp.tmp_provsql') IS NOT NULL THEN
DROP TABLE tmp_provsql;
END IF;
EXECUTE format('CREATE TEMP TABLE tmp_provsql ON COMMIT DROP AS TABLE %s', oldtbl);
ALTER TABLE tmp_provsql RENAME provsql TO provenance;
-- The mapping is keyed by gate identity (input-token UUIDs), so peel any
-- transparent annotation wrapper (e.g. the inversion-free certificate a
-- certified query attaches to its row roots) off the captured tokens.
UPDATE tmp_provsql SET provenance = provsql.strip_annotations(provenance)
WHERE provsql.get_gate_type(provenance) = 'annotation';
IF preserve_case THEN
EXECUTE format('CREATE TABLE %I AS SELECT %s AS value, provenance FROM tmp_provsql', newtbl, att);
EXECUTE format('CREATE INDEX ON %I(provenance)', newtbl);
ELSE
EXECUTE format('CREATE TABLE %s AS SELECT %s AS value, provenance FROM tmp_provsql', newtbl, att);
EXECUTE format('CREATE INDEX ON %s(provenance)', newtbl);
END IF;
IF maintained THEN
-- Register so genuine inserts into oldtbl keep the mapping current
-- (see provenance_guard); keyed to the input token, so it survives the
-- provsql rewrites that data modification performs.
INSERT INTO provsql.provenance_mapping_registry(mapping, source, attribute)
VALUES (
(CASE WHEN preserve_case THEN to_regclass(format('%I', newtbl))
ELSE to_regclass(newtbl) END)::oid,
oldtbl::oid, att)
ON CONFLICT (mapping)
DO UPDATE SET source = EXCLUDED.source, attribute = EXCLUDED.attribute;
END IF;
END
$$ LANGUAGE plpgsql;
-- provenance_guard appends to registered mappings on genuine inserts.
CREATE OR REPLACE FUNCTION provenance_guard()
RETURNS TRIGGER AS $$
DECLARE
_m RECORD;
BEGIN
IF TG_OP = 'INSERT' THEN
IF NEW.provsql IS NULL THEN
-- A genuine insert: mint a fresh atomic input variable. This is the
-- one place a new input token is born, so it is also where any
-- maintained mapping on this table is extended (keyed to that token).
-- Data-modification re-insertions (INSERT ... SELECT * FROM OLD_TABLE)
-- carry a supplied provsql and take the ELSE branch, so they are
-- correctly skipped: the validity stays keyed to the original input,
-- which is exactly the child a later monus/update gate wraps.
NEW.provsql := public.uuid_generate_v4();
FOR _m IN SELECT mapping, attribute
FROM provsql.provenance_mapping_registry WHERE source = TG_RELID
LOOP
EXECUTE format(
'INSERT INTO %s(value, provenance) SELECT ($1).%I, $2',
_m.mapping::regclass, _m.attribute)
USING NEW, NEW.provsql;
END LOOP;
ELSE
PERFORM provsql.set_table_info(TG_RELID, 'opaque');
END IF;
ELSIF TG_OP = 'UPDATE' THEN
IF NEW.provsql IS DISTINCT FROM OLD.provsql THEN
PERFORM provsql.set_table_info(TG_RELID, 'opaque');
END IF;
END IF;
RETURN NEW;
END;
$$ LANGUAGE plpgsql SET search_path=provsql,pg_temp,public
SECURITY DEFINER;
-- cleanup_table_info also forgets registry entries when a table is dropped.
CREATE OR REPLACE FUNCTION cleanup_table_info()
RETURNS event_trigger AS
$$
DECLARE
r RECORD;
BEGIN
FOR r IN
SELECT objid FROM pg_event_trigger_dropped_objects()
WHERE object_type IN ('table', 'foreign table', 'materialized view')
LOOP
PERFORM provsql.remove_table_info(r.objid);
-- Forget any maintained mapping whose source or mapping table is gone.
DELETE FROM provsql.provenance_mapping_registry
WHERE source = r.objid OR mapping = r.objid;
END LOOP;
END
$$ LANGUAGE plpgsql;
-- The view-based mapping helper is removed (superseded by maintained tables).
DROP FUNCTION IF EXISTS create_provenance_mapping_view(text, regclass, text, bool);
/**
* @brief Placeholder for @c "(predicate) | (predicate)" on two events.
*
* Conditions one comparison event on another when both operands are written
* as comparisons rather than pre-built tokens (e.g.
* @c "probability((x >= 2000) | (x >= 1000))"): an @c random_variable /
* @c agg_token comparison is statically @c boolean-typed, so neither the
* @c "uuid | uuid" (@c cond) nor the @c "uuid | boolean" (@c cond_predicate)
* operator resolves. Never executes: the ProvSQL planner hook lowers each
* Boolean operand to its event gate and emits @c cond(target, evidence), so
* the result carries the correlation-aware @c Pr(A ∧ B) / Pr(B). Returns
* @c uuid, so @c "A | B" is a first-class event token in every position
* (a @c probability(uuid) argument, a projected column, a further @c "|").
*/
CREATE OR REPLACE FUNCTION predicate_cond_predicate(target boolean, evidence boolean)
RETURNS UUID AS
$$
BEGIN
RAISE EXCEPTION '(predicate) | (predicate) must be rewritten by the ProvSQL '
'planner hook: both operands must be Boolean combinations of '
'random_variable / aggregate comparisons (is provsql.active off?)';
END
$$ LANGUAGE plpgsql IMMUTABLE STRICT PARALLEL SAFE;
DO $$ BEGIN
IF NOT EXISTS (
SELECT 1 FROM pg_operator o
JOIN pg_namespace n ON n.oid = o.oprnamespace
WHERE n.nspname = 'provsql' AND o.oprname = '|'
AND o.oprleft = 'pg_catalog.bool'::regtype::oid
AND o.oprright = 'pg_catalog.bool'::regtype::oid
AND o.oprcode <> 0
) THEN
CREATE OPERATOR | (LEFTARG=boolean, RIGHTARG=boolean, PROCEDURE=predicate_cond_predicate);
END IF;
END $$;
/**
* @brief Create a guarded-selection gate over scalar (RV) children.
*
* Builds a deterministic @c gate_case from the flattened wire list
* @c [guard_1, value_1, ..., guard_k, value_k, default] (odd length): the
* value of the first guard event that holds, else the default (first-match
* semantics). Each guard is a Boolean event token (a @c gate_cmp or Boolean
* combination); each value and the default are scalar-producing gates
* (@c gate_rv, @c gate_value, @c gate_arith, another @c gate_case, ...). The
* token UUID is derived deterministically from @p children so identical
* @c CASE expressions share their gate.
*
* @param children Flattened guard/value wires ending with the default
* (@c array_length must be odd and @c >= 1).
* @return UUID of the (possibly pre-existing) @c gate_case.
*/
CREATE OR REPLACE FUNCTION provenance_case(
children UUID[]
)
RETURNS UUID AS
$$
DECLARE
case_token UUID;
BEGIN
IF array_length(children, 1) IS NULL OR array_length(children, 1) % 2 = 0 THEN
RAISE EXCEPTION 'provenance_case expects an odd number of children '
'(guard/value pairs followed by a default), got %',
coalesce(array_length(children, 1), 0);
END IF;
case_token := public.uuid_generate_v5(
uuid_ns_provsql(),
concat('case', children::text)
);
PERFORM create_gate(case_token, 'case', children);
RETURN case_token;
END
$$ LANGUAGE plpgsql
SET search_path=provsql,pg_temp,public
SECURITY DEFINER
IMMUTABLE
PARALLEL SAFE
STRICT;
/**
* @brief Deterministic truth of a Boolean guard sub-circuit over aggregate
* comparisons, evaluated in the actual world (all input tuples present).
*
* Guards are the shapes @c having_Expr_to_provenance_cmp mints: @c cmp gates
* over aggregate-valued children (comparison-operator OID in @c info1),
* @c times / @c plus combinations (AND / OR, with negation pushed into the
* comparison operators), and the @c one / @c zero indicators of regular
* (aggregate-free) conditions. Uses Kleene three-valued logic: returns
* @c NULL on any other gate shape, or when an operand's deterministic value
* cannot be resolved.
*/
CREATE OR REPLACE FUNCTION agg_guard_holds(token UUID)
RETURNS boolean AS
$$
DECLARE
gt provenance_gate := get_gate_type(token);
ch uuid[];
opname text;
l numeric;
r numeric;
all_true boolean;
any_true boolean;
any_null boolean;
BEGIN
IF gt = 'one' THEN
RETURN true;
ELSIF gt = 'zero' THEN
RETURN false;
ELSIF gt IN ('times', 'plus') THEN
SELECT bool_and(h), bool_or(h), bool_or(h IS NULL)
INTO all_true, any_true, any_null
FROM (SELECT provsql.agg_guard_holds(c) AS h
FROM unnest(get_children(token)) AS c) AS s;
IF gt = 'times' THEN
-- AND: false dominates unknown (bool_and skips NULL inputs, so it is
-- false exactly when some child is false).
RETURN CASE WHEN NOT all_true THEN false
WHEN any_null THEN NULL
ELSE true END;
ELSE
-- OR: true dominates unknown.
RETURN CASE WHEN any_true THEN true
WHEN any_null THEN NULL
ELSE false END;
END IF;
ELSIF gt = 'cmp' THEN
ch := get_children(token);
l := agg_gate_value(ch[1]);
r := agg_gate_value(ch[2]);
IF l IS NULL OR r IS NULL THEN
RETURN NULL;
END IF;
SELECT oprname INTO opname
FROM pg_catalog.pg_operator WHERE oid = (get_infos(token)).info1;
RETURN CASE opname
WHEN '<' THEN l < r
WHEN '<=' THEN l <= r
WHEN '=' THEN l = r
WHEN '<>' THEN l <> r
WHEN '>=' THEN l >= r
WHEN '>' THEN l > r
END;
END IF;
RETURN NULL;
END
$$ LANGUAGE plpgsql STABLE STRICT PARALLEL SAFE
SET search_path=provsql,pg_temp,public;
/**
* @brief Deterministic (actual-world) scalar value of an aggregate-carrying
* gate.
*
* Resolves the value an aggregate expression takes on the actual data -- the
* value an @c agg_token display cell carries: @c agg / @c arith gates record
* it in @c extra (set by aggregate evaluation and @c agg_arith_make), a
* @c value gate carries its constant, a @c semimod wraps a value gate, a
* @c conditioned gate has its target's value, and a @c case gate selects the
* first branch whose guard holds in the actual world (per
* @c agg_guard_holds), else the default. Returns @c NULL when the gate is
* not aggregate-carrying or the value cannot be resolved (e.g. a
* non-numeric aggregate).
*/
CREATE OR REPLACE FUNCTION agg_gate_value(token UUID)
RETURNS numeric AS
$$
DECLARE
gt provenance_gate := get_gate_type(token);
ch uuid[];
n integer;
holds boolean;
BEGIN
IF gt IN ('agg', 'arith', 'value') THEN
BEGIN
RETURN get_extra(token)::numeric;
EXCEPTION WHEN others THEN
RETURN NULL; -- non-numeric aggregate (e.g. min over text)
END;
ELSIF gt = 'semimod' THEN
RETURN agg_gate_value((get_children(token))[2]);
ELSIF gt = 'conditioned' THEN
RETURN agg_gate_value((get_children(token))[1]);
ELSIF gt = 'case' THEN
ch := get_children(token);
n := array_length(ch, 1);
FOR i IN 1 .. (n - 1) / 2 LOOP
holds := agg_guard_holds(ch[2 * i - 1]);
IF holds IS NULL THEN
RETURN NULL;
ELSIF holds THEN
RETURN agg_gate_value(ch[2 * i]);
END IF;
END LOOP;
RETURN agg_gate_value(ch[n]);
END IF;
RETURN NULL;
END
$$ LANGUAGE plpgsql STABLE STRICT PARALLEL SAFE
SET search_path=provsql,pg_temp,public;
/**
* @brief Recover the @c "value (*)" display string for an aggregation gate
*
* Companion helper to the @c provsql.aggtoken_text_as_uuid GUC. With
* the GUC on, an @c agg_token cell prints as the underlying provenance
* UUID, which is convenient for tooling that wants to click through to
* the circuit but loses the human-readable aggregate value. This
* function takes such a UUID and returns the original @c "value (*)"
* string by reading the gate's @c extra (set by aggregate evaluation
* for @c agg gates, and by @c agg_arith_make for the @c arith gates
* that agg_token arithmetic mints); for the other aggregate-carrying
* gates (@c case, @c conditioned, @c semimod, @c value) the value is
* resolved through the circuit by @c agg_gate_value. Returns @c NULL
* if @p token does not resolve to an aggregate-carrying gate.
*
* @param token UUID of an @c agg gate (typically obtained from an
* @c agg_token cell when @c aggtoken_text_as_uuid is on,
* or via a manual UUID cast otherwise).
*/
CREATE OR REPLACE FUNCTION agg_token_value_text(token UUID)
RETURNS text AS
$$
SELECT CASE
-- agg gates: extra is set by aggregate evaluation; arith gates
-- (agg_token arithmetic): extra is recorded by agg_arith_make.
WHEN provsql.get_gate_type(token) IN ('agg', 'arith')
THEN provsql.get_extra(token) || ' (*)'
-- other aggregate-carrying gates: resolve the actual-world value
-- through the circuit.
WHEN provsql.get_gate_type(token) IN ('case', 'conditioned', 'semimod', 'value')
THEN provsql.agg_gate_value(token)::text || ' (*)'
ELSE NULL
END;
$$ LANGUAGE sql STABLE STRICT PARALLEL SAFE;
/**
* @brief Construct a gamma-distribution random variable with shape @p k
* (any positive real) and rate @p lambda
*
* The gamma distribution generalises Erlang to non-integer shape; its
* CDF is the regularised lower incomplete gamma, evaluated in closed
* form by the analytic passes. Sums of independent gammas with the
* same rate fold to a single gamma in the simplifier.
*
* Validation:
* - @p k must be finite and strictly positive. An integer @p k (in
* @c integer range) is silently routed through @c erlang -- the gamma
* with integer shape *is* Erlang -- so gamma(2, λ) shares
* its gate encoding and closure interplay with erlang(2, λ).
* - @p lambda must be finite and strictly positive.
*
* @warning VOLATILE is load-bearing; see the warning on
* @ref normal.
*
* @sa Wikipedia: Gamma distribution
*/
CREATE OR REPLACE FUNCTION gamma(k double precision, lambda double precision)
RETURNS random_variable AS
$$
DECLARE
token uuid;
BEGIN
IF NOT provsql.is_finite_float8(k) THEN
RAISE EXCEPTION 'provsql.gamma: k must be finite (got %)', k;
END IF;
IF k <= 0 THEN
RAISE EXCEPTION 'provsql.gamma: k must be strictly positive (got %)', k;
END IF;
IF NOT provsql.is_finite_float8(lambda) THEN
RAISE EXCEPTION 'provsql.gamma: lambda must be finite (got %)', lambda;
END IF;
IF lambda <= 0 THEN
RAISE EXCEPTION 'provsql.gamma: lambda must be strictly positive (got %)', lambda;
END IF;
IF k = floor(k) AND k <= 2147483647 THEN
RETURN provsql.erlang(k::integer, lambda);
END IF;
token := public.uuid_generate_v4();
PERFORM provsql.create_gate(token, 'rv');
PERFORM provsql.set_extra(token, 'gamma:' || k || ',' || lambda);
RETURN provsql.random_variable_make(token);
END
$$ LANGUAGE plpgsql STRICT VOLATILE PARALLEL SAFE;
/**
* @brief Construct a chi-squared random variable with @p k degrees of
* freedom: syntactic sugar for gamma(k/2, 1/2)
*
* @p k is accepted as @c double @c precision so fractional degrees of
* freedom work; it must be finite and strictly positive. Even degrees
* of freedom route through @c erlang via @c gamma's integer-shape rule.
*
* @warning VOLATILE is load-bearing; see the warning on
* @ref normal.
*
* @sa Wikipedia: Chi-squared distribution
*/
CREATE OR REPLACE FUNCTION chi_squared(k double precision)
RETURNS random_variable AS
$$
BEGIN
IF NOT provsql.is_finite_float8(k) THEN
RAISE EXCEPTION 'provsql.chi_squared: k must be finite (got %)', k;
END IF;
IF k <= 0 THEN
RAISE EXCEPTION 'provsql.chi_squared: k must be strictly positive (got %)', k;
END IF;
RETURN provsql.gamma(k / 2, 0.5);
END
$$ LANGUAGE plpgsql STRICT VOLATILE PARALLEL SAFE;
/**
* @brief Construct a log-normal random variable: @c exp of a
* Normal(@p mu, @p sigma), parameterised by the underlying
* normal (so its median is exp(mu) and its mean
* exp(mu + sigma^2/2))
*
* The multiplicative counterpart of @c normal: products of independent
* lognormals fold to a lognormal in the simplifier, and the
* exp(normal(...)) / ln(lognormal(...)) bridges fold
* in both directions, so log-scale models stay closed-form.
*
* Validation mirrors @c normal: both parameters must be finite,
* @p sigma non-negative; the degenerate @c sigma = 0 case is silently
* routed through @c as_random (a Dirac at exp(mu)).
*
* @warning VOLATILE is load-bearing; see the warning on
* @ref normal.
*
* @sa Wikipedia: Log-normal distribution
*/
CREATE OR REPLACE FUNCTION lognormal(mu double precision, sigma double precision)
RETURNS random_variable AS
$$
DECLARE
token uuid;
BEGIN
IF NOT provsql.is_finite_float8(mu) OR NOT provsql.is_finite_float8(sigma) THEN
RAISE EXCEPTION 'provsql.lognormal: parameters must be finite (got mu=%, sigma=%)', mu, sigma;
END IF;
IF sigma < 0 THEN
RAISE EXCEPTION 'provsql.lognormal: sigma must be non-negative (got %)', sigma;
END IF;
IF sigma = 0 THEN
RETURN provsql.as_random(exp(mu));
END IF;
token := public.uuid_generate_v4();
PERFORM provsql.create_gate(token, 'rv');
PERFORM provsql.set_extra(token, 'lognormal:' || mu || ',' || sigma);
RETURN provsql.random_variable_make(token);
END
$$ LANGUAGE plpgsql STRICT VOLATILE PARALLEL SAFE;
/**
* @brief Construct a Weibull random variable with shape @p k and
* scale @p lambda
*
* @p lambda is the SCALE (the 63.2% quantile), not a rate: @c k = 1 is
* the exponential with rate 1/lambda, and that case is
* silently routed through @c exponential to share its gate. The shape
* tunes the hazard: @c k < 1 infant mortality, @c k > 1 wear-out.
* Quantiles are exact, truncated moments are closed-form (via the
* regularised incomplete gamma), and the min of i.i.d. Weibulls has a
* closed-form mean (min-stability).
*
* Validation: both parameters must be finite and strictly positive.
*
* @warning VOLATILE is load-bearing; see the warning on
* @ref normal.
*
* @sa Wikipedia: Weibull distribution
*/
CREATE OR REPLACE FUNCTION weibull(k double precision, lambda double precision)
RETURNS random_variable AS
$$
DECLARE
token uuid;
BEGIN
IF NOT provsql.is_finite_float8(k) OR NOT provsql.is_finite_float8(lambda) THEN
RAISE EXCEPTION 'provsql.weibull: parameters must be finite (got k=%, lambda=%)', k, lambda;
END IF;
IF k <= 0 OR lambda <= 0 THEN
RAISE EXCEPTION 'provsql.weibull: parameters must be strictly positive (got k=%, lambda=%)', k, lambda;
END IF;
IF k = 1 THEN
RETURN provsql.exponential(1 / lambda);
END IF;
token := public.uuid_generate_v4();
PERFORM provsql.create_gate(token, 'rv');
PERFORM provsql.set_extra(token, 'weibull:' || k || ',' || lambda);
RETURN provsql.random_variable_make(token);
END
$$ LANGUAGE plpgsql STRICT VOLATILE PARALLEL SAFE;
/**
* @brief Construct a Pareto random variable with scale (minimum)
* @p xm and shape @p alpha
*
* The canonical heavy-tailed power law. Raw moments are @b infinite
* for alpha <= k and reported as Infinity (the mean
* for alpha <= 1, the variance for alpha <= 2)
* rather than estimated; quantiles, truncated moments, conditional
* sampling (self-similarity: X | X > a is Pareto(a, alpha)),
* and Pareto-vs-Pareto comparisons are all exact.
*
* Validation: both parameters must be finite and strictly positive.
*
* @warning VOLATILE is load-bearing; see the warning on
* @ref normal.
*
* @sa Wikipedia: Pareto distribution
*/
CREATE OR REPLACE FUNCTION pareto(xm double precision, alpha double precision)
RETURNS random_variable AS
$$
DECLARE
token uuid;
BEGIN
IF NOT provsql.is_finite_float8(xm) OR NOT provsql.is_finite_float8(alpha) THEN
RAISE EXCEPTION 'provsql.pareto: parameters must be finite (got xm=%, alpha=%)', xm, alpha;
END IF;
IF xm <= 0 OR alpha <= 0 THEN
RAISE EXCEPTION 'provsql.pareto: parameters must be strictly positive (got xm=%, alpha=%)', xm, alpha;
END IF;
token := public.uuid_generate_v4();
PERFORM provsql.create_gate(token, 'rv');
PERFORM provsql.set_extra(token, 'pareto:' || xm || ',' || alpha);
RETURN provsql.random_variable_make(token);
END
$$ LANGUAGE plpgsql STRICT VOLATILE PARALLEL SAFE;
CREATE OR REPLACE FUNCTION inverse_gamma(alpha double precision, beta double precision)
RETURNS random_variable AS
$$
DECLARE
token uuid;
BEGIN
IF NOT provsql.is_finite_float8(alpha) OR NOT provsql.is_finite_float8(beta) THEN
RAISE EXCEPTION 'provsql.inverse_gamma: parameters must be finite (got alpha=%, beta=%)', alpha, beta;
END IF;
IF alpha <= 0 OR beta <= 0 THEN
RAISE EXCEPTION 'provsql.inverse_gamma: parameters must be strictly positive (got alpha=%, beta=%)', alpha, beta;
END IF;
token := public.uuid_generate_v4();
PERFORM provsql.create_gate(token, 'rv');
PERFORM provsql.set_extra(token, 'inverse_gamma:' || alpha || ',' || beta);
RETURN provsql.random_variable_make(token);
END
$$ LANGUAGE plpgsql STRICT VOLATILE PARALLEL SAFE;
CREATE OR REPLACE FUNCTION inverse_gaussian(mu double precision, lambda double precision)
RETURNS random_variable AS
$$
DECLARE
token uuid;
BEGIN
IF NOT provsql.is_finite_float8(mu) OR NOT provsql.is_finite_float8(lambda) THEN
RAISE EXCEPTION 'provsql.inverse_gaussian: parameters must be finite (got mu=%, lambda=%)', mu, lambda;
END IF;
IF mu <= 0 OR lambda <= 0 THEN
RAISE EXCEPTION 'provsql.inverse_gaussian: parameters must be strictly positive (got mu=%, lambda=%)', mu, lambda;
END IF;
token := public.uuid_generate_v4();
PERFORM provsql.create_gate(token, 'rv');
PERFORM provsql.set_extra(token, 'inverse_gaussian:' || mu || ',' || lambda);
RETURN provsql.random_variable_make(token);
END
$$ LANGUAGE plpgsql STRICT VOLATILE PARALLEL SAFE;
CREATE OR REPLACE FUNCTION wald(mu double precision, lambda double precision)
RETURNS random_variable AS
$$
SELECT provsql.inverse_gaussian(mu, lambda);
$$ LANGUAGE sql VOLATILE PARALLEL SAFE;
/**
* @brief Build a discrete (categorical) random variable from outcomes
* and UNNORMALISED log-masses
*
* The shared back end of the discrete count constructors (@c poisson,
* @c binomial, @c geometric, @c hypergeometric,
* @c negative_binomial), and directly usable for any custom discrete
* pmf: the log-masses are shifted by their maximum (so only relative
* magnitudes matter and no @c exp underflows), outcomes whose relative
* mass is below 1e-15 are dropped, and the rest is
* renormalised before being handed to @c categorical. Working in log
* space keeps arbitrarily large parameters stable (e.g. a
* Poisson(1000) pmf whose linear-space recurrence would
* underflow at @c exp(-1000)).
*
* @param outcomes outcome values, same length as @p log_pmf
* @param log_pmf natural logs of the (unnormalised) masses
*/
CREATE OR REPLACE FUNCTION categorical_from_log_pmf(
outcomes double precision[], log_pmf double precision[])
RETURNS random_variable AS
$$
DECLARE
n int := array_length(outcomes, 1);
max_lp double precision := '-Infinity';
kept_o double precision[] := '{}';
kept_p double precision[] := '{}';
total double precision := 0;
v double precision;
i int;
BEGIN
IF n IS NULL OR n = 0 OR n <> coalesce(array_length(log_pmf, 1), 0) THEN
RAISE EXCEPTION 'provsql.categorical_from_log_pmf: outcomes and log_pmf must be non-empty arrays of the same length';
END IF;
FOR i IN 1..n LOOP
IF log_pmf[i] > max_lp THEN max_lp := log_pmf[i]; END IF;
END LOOP;
IF max_lp = '-Infinity' THEN
RAISE EXCEPTION 'provsql.categorical_from_log_pmf: all masses are zero';
END IF;
FOR i IN 1..n LOOP
v := exp(log_pmf[i] - max_lp);
IF v >= 1e-15 THEN
kept_o := array_append(kept_o, outcomes[i]);
kept_p := array_append(kept_p, v);
total := total + v;
END IF;
END LOOP;
FOR i IN 1..array_length(kept_p, 1) LOOP
kept_p[i] := kept_p[i] / total;
END LOOP;
RETURN provsql.categorical(kept_p, kept_o);
END
$$ LANGUAGE plpgsql STRICT VOLATILE PARALLEL SAFE;
/**
* @brief Construct a Poisson random variable with mean @p lambda, as a
* truncated categorical
*
* The pmf is enumerated over [max(0, λ-12√λ), λ+12√λ+30] (the
* omitted tails carry ~1e-30 of mass) by the log-space recurrence
* ln p(k+1) = ln p(k) + ln λ - ln(k+1) and handed to
* @c categorical_from_log_pmf, so moments, quantiles, and (in)equality
* comparisons are exact over the enumerated support. @c lambda = 0 is
* a Dirac at @c 0 (routed through @c as_random); supports up to 10000
* outcomes (λ up to ~170000), beyond which it raises -- approximate
* huge means by @c normal(λ, √λ) instead.
*
* @sa Wikipedia: Poisson distribution
*/
CREATE OR REPLACE FUNCTION poisson(lambda double precision)
RETURNS random_variable AS
$$
DECLARE
lo int;
hi int;
outcomes double precision[] := '{}';
lps double precision[] := '{}';
lp double precision := 0;
k int;
BEGIN
IF NOT provsql.is_finite_float8(lambda) OR lambda < 0 THEN
RAISE EXCEPTION 'provsql.poisson: lambda must be finite and non-negative (got %)', lambda;
END IF;
IF lambda = 0 THEN
RETURN provsql.as_random(0);
END IF;
lo := greatest(0, floor(lambda - 12 * sqrt(lambda)))::int;
hi := ceil(lambda + 12 * sqrt(lambda))::int + 30;
IF hi - lo + 1 > 10000 THEN
RAISE EXCEPTION 'provsql.poisson: support window of % outcomes exceeds 10000; approximate with normal(%, sqrt(%))', hi - lo + 1, lambda, lambda;
END IF;
-- ln p(0) = -λ; walk the recurrence, keeping only the window.
lp := -lambda;
FOR k IN 1..hi LOOP
lp := lp + ln(lambda) - ln(k::double precision);
IF k >= lo THEN
outcomes := array_append(outcomes, k::double precision);
lps := array_append(lps, lp);
END IF;
END LOOP;
IF lo = 0 THEN
outcomes := array_prepend(0::double precision, outcomes);
lps := array_prepend(-lambda, lps);
END IF;
RETURN provsql.categorical_from_log_pmf(outcomes, lps);
END
$$ LANGUAGE plpgsql STRICT VOLATILE PARALLEL SAFE;
/**
* @brief Construct a Beta(α, β) random variable on the unit interval
*
* The conjugate prior of Bernoulli / binomial success probabilities:
* closed-form moments, CDF through the regularised incomplete beta,
* quantiles through the generic CDF bisection over the finite
* @c [0, 1] support, and closed-form truncated moments (interval
* conditioning). Beta(1, 1) IS Uniform(0, 1) and is
* silently routed through @c uniform to share its richer closed forms.
*
* Validation: both shapes must be finite and strictly positive.
*
* @warning VOLATILE is load-bearing; see the warning on
* @ref normal.
*
* @sa Wikipedia: Beta distribution
*/
CREATE OR REPLACE FUNCTION beta(alpha double precision, beta double precision)
RETURNS random_variable AS
$$
DECLARE
token uuid;
BEGIN
IF NOT provsql.is_finite_float8(alpha) OR NOT provsql.is_finite_float8(beta) THEN
RAISE EXCEPTION 'provsql.beta: parameters must be finite (got alpha=%, beta=%)', alpha, beta;
END IF;
IF alpha <= 0 OR beta <= 0 THEN
RAISE EXCEPTION 'provsql.beta: parameters must be strictly positive (got alpha=%, beta=%)', alpha, beta;
END IF;
IF alpha = 1 AND beta = 1 THEN
RETURN provsql.uniform(0, 1);
END IF;
token := public.uuid_generate_v4();
PERFORM provsql.create_gate(token, 'rv');
PERFORM provsql.set_extra(token, 'beta:' || alpha || ',' || beta);
RETURN provsql.random_variable_make(token);
END
$$ LANGUAGE plpgsql STRICT VOLATILE PARALLEL SAFE;
/**
* @brief Construct a Binomial(n, p) random variable (number of
* successes in @p n independent trials), as a categorical
*
* Enumerated over {0..n} by the log-space recurrence
* ln p(k+1) = ln p(k) + ln((n-k)/(k+1)) + ln(p/(1-p))
* (outcomes below 1e-15 relative mass are dropped). @c p = 0 /
* @c p = 1 are Diracs at @c 0 / @c n; @c n is capped at 10000.
*
* @sa Wikipedia: Binomial distribution
*/
CREATE OR REPLACE FUNCTION binomial(n integer, p double precision)
RETURNS random_variable AS
$$
DECLARE
outcomes double precision[] := '{}';
lps double precision[] := '{}';
lp double precision;
k int;
BEGIN
IF n IS NULL OR n < 0 THEN
RAISE EXCEPTION 'provsql.binomial: n must be non-negative (got %)', n;
END IF;
IF NOT provsql.is_finite_float8(p) OR p < 0 OR p > 1 THEN
RAISE EXCEPTION 'provsql.binomial: p must be in [0, 1] (got %)', p;
END IF;
IF n > 10000 THEN
RAISE EXCEPTION 'provsql.binomial: n = % exceeds 10000; approximate with normal(n*p, sqrt(n*p*(1-p)))', n;
END IF;
IF n = 0 OR p = 0 THEN
RETURN provsql.as_random(0);
END IF;
IF p = 1 THEN
RETURN provsql.as_random(n);
END IF;
lp := n * ln(1 - p); -- ln p(0)
outcomes := array_append(outcomes, 0::double precision);
lps := array_append(lps, lp);
FOR k IN 0..(n - 1) LOOP
lp := lp + ln((n - k)::double precision / (k + 1)) + ln(p / (1 - p));
outcomes := array_append(outcomes, (k + 1)::double precision);
lps := array_append(lps, lp);
END LOOP;
RETURN provsql.categorical_from_log_pmf(outcomes, lps);
END
$$ LANGUAGE plpgsql STRICT VOLATILE PARALLEL SAFE;
/**
* @brief Construct a Geometric(p) random variable -- the number of
* TRIALS up to and including the first success (support
* starting at 1; subtract 1 for the failures convention)
*
* P(X = k) = (1-p)^{k-1} p, enumerated up to the 1e-15
* relative-mass tail and renormalised. @c p = 1 is a Dirac at @c 1.
*
* @sa Wikipedia: Geometric distribution
*/
CREATE OR REPLACE FUNCTION geometric(p double precision)
RETURNS random_variable AS
$$
DECLARE
k_max int;
outcomes double precision[] := '{}';
lps double precision[] := '{}';
k int;
BEGIN
IF NOT provsql.is_finite_float8(p) OR p <= 0 OR p > 1 THEN
RAISE EXCEPTION 'provsql.geometric: p must be in (0, 1] (got %)', p;
END IF;
IF p = 1 THEN
RETURN provsql.as_random(1);
END IF;
k_max := 1 + ceil(ln(1e-15) / ln(1 - p))::int;
IF k_max > 10000 THEN
RAISE EXCEPTION 'provsql.geometric: support window of % outcomes exceeds 10000 (p = % is too small); approximate with exponential(%)', k_max, p, p;
END IF;
FOR k IN 1..k_max LOOP
outcomes := array_append(outcomes, k::double precision);
lps := array_append(lps, (k - 1) * ln(1 - p) + ln(p));
END LOOP;
RETURN provsql.categorical_from_log_pmf(outcomes, lps);
END
$$ LANGUAGE plpgsql STRICT VOLATILE PARALLEL SAFE;
/**
* @brief Construct a Hypergeometric(N, K, n) random variable: the
* number of marked items among @p n draws WITHOUT replacement
* from a population of @p pop_n items of which @p k_marked are
* marked
*
* The exact finite support [max(0, n-(N-K)), min(n, K)] is
* enumerated by the pmf ratio recurrence (in log space, so large
* populations cannot overflow) and normalised -- exact "sampling
* without replacement" probabilities with no combinatorial functions
* needed.
*
* @sa Wikipedia: Hypergeometric distribution
*/
CREATE OR REPLACE FUNCTION hypergeometric(pop_n integer, k_marked integer, n integer)
RETURNS random_variable AS
$$
DECLARE
lo int;
hi int;
outcomes double precision[] := '{}';
lps double precision[] := '{}';
lp double precision := 0; -- relative log-mass; normalised later
k int;
BEGIN
IF pop_n IS NULL OR k_marked IS NULL OR n IS NULL
OR pop_n < 0 OR k_marked < 0 OR n < 0
OR k_marked > pop_n OR n > pop_n THEN
RAISE EXCEPTION 'provsql.hypergeometric: need 0 <= k_marked, n <= pop_n (got pop_n=%, k_marked=%, n=%)', pop_n, k_marked, n;
END IF;
lo := greatest(0, n - (pop_n - k_marked));
hi := least(n, k_marked);
IF hi - lo + 1 > 10000 THEN
RAISE EXCEPTION 'provsql.hypergeometric: support window of % outcomes exceeds 10000', hi - lo + 1;
END IF;
outcomes := array_append(outcomes, lo::double precision);
lps := array_append(lps, lp);
FOR k IN lo..(hi - 1) LOOP
-- pmf(k+1)/pmf(k) = (K-k)(n-k) / ((k+1)(N-K-n+k+1))
lp := lp + ln((k_marked - k)::double precision * (n - k))
- ln((k + 1)::double precision * (pop_n - k_marked - n + k + 1));
outcomes := array_append(outcomes, (k + 1)::double precision);
lps := array_append(lps, lp);
END LOOP;
RETURN provsql.categorical_from_log_pmf(outcomes, lps);
END
$$ LANGUAGE plpgsql STRICT VOLATILE PARALLEL SAFE;
/**
* @brief Construct a negative-binomial random variable: the number of
* FAILURES before the @p r-th success (support starting at 0),
* with real @p r > 0 allowed (the Polya / overdispersed-count
* parameterisation, the Poisson-Gamma mixture)
*
* P(X = k) = C(k+r-1, k) p^r (1-p)^k, enumerated by the
* log-space recurrence
* ln p(k+1) = ln p(k) + ln((k+r)/(k+1)) + ln(1-p) up to the
* 1e-15 relative-mass tail. @c p = 1 is a Dirac at @c 0.
*
* @sa Wikipedia: Negative binomial distribution
*/
CREATE OR REPLACE FUNCTION negative_binomial(r double precision, p double precision)
RETURNS random_variable AS
$$
DECLARE
outcomes double precision[] := '{}';
lps double precision[] := '{}';
lp double precision;
max_lp double precision;
mean double precision;
k int := 0;
BEGIN
IF NOT provsql.is_finite_float8(r) OR r <= 0 THEN
RAISE EXCEPTION 'provsql.negative_binomial: r must be finite and strictly positive (got %)', r;
END IF;
IF NOT provsql.is_finite_float8(p) OR p <= 0 OR p > 1 THEN
RAISE EXCEPTION 'provsql.negative_binomial: p must be in (0, 1] (got %)', p;
END IF;
IF p = 1 THEN
RETURN provsql.as_random(0);
END IF;
mean := r * (1 - p) / p;
lp := r * ln(p); -- ln p(0)
max_lp := lp;
outcomes := array_append(outcomes, 0::double precision);
lps := array_append(lps, lp);
LOOP
lp := lp + ln((k + r) / (k + 1)) + ln(1 - p);
k := k + 1;
IF lp > max_lp THEN max_lp := lp; END IF;
outcomes := array_append(outcomes, k::double precision);
lps := array_append(lps, lp);
EXIT WHEN k > mean AND lp < max_lp + ln(1e-15);
IF k >= 10000 THEN
RAISE EXCEPTION 'provsql.negative_binomial: support window exceeds 10000 outcomes (r=%, p=%)', r, p;
END IF;
END LOOP;
RETURN provsql.categorical_from_log_pmf(outcomes, lps);
END
$$ LANGUAGE plpgsql STRICT VOLATILE PARALLEL SAFE;
/**
* @brief Catalog of the registered continuous-distribution families.
*
* One row per @c gate_rv family known to this build of the extension:
* @c name is the on-disk token (the part before the colon in the gate's
* @c extra encoding), @c nparams the parameter count, @c param_names the
* conventional parameter symbols in @c extra order (e.g.
* {μ, σ}), and @c label a short display glyph (e.g. @c "N",
* @c "Γ"). UI clients (ProvSQL Studio's circuit inspector) read this to
* render families they were not hard-coded for, so a newly added family
* shows up without a client release.
*/
CREATE OR REPLACE FUNCTION rv_families()
RETURNS TABLE(name TEXT, nparams INT, param_names TEXT[], label TEXT) AS
'provsql','rv_families' LANGUAGE C STABLE PARALLEL SAFE;
/**
* @brief Gaussian-mixture-model (GMM) constructor.
*
* Packages the common fitted-density pattern -- a categorical choice
* among Normal components -- into one call:
*
* @code
* provsql.gmm(weights => ARRAY[0.3, 0.5, 0.2],
* means => ARRAY[120.0, 380.0, 1200.0],
* stddevs => ARRAY[40.0, 90.0, 250.0])
* @endcode
*
* No new gate: the mixture decomposes into a stick-breaking cascade of
* Bernoulli @c gate_mixture nodes over @c gate_rv Normal leaves
* (component @c i is selected with conditional probability
* @c w_i / (w_i + ... + w_n), so the joint selection probabilities are
* exactly @p weights), which every evaluator already handles: moments
* are closed-form through the mixture recursion, sampling is exact,
* and comparisons ride the existing mixture machinery. Zero-weight
* components are skipped; a single positive-weight component returns
* its Normal directly (no mixture node).
*
* Validation mirrors @c categorical: same-length non-empty arrays,
* weights finite in [0, 1] summing to 1 within @c 1e-9; the
* component parameters are validated by @c provsql.normal (finite
* @c mu, non-negative @c sigma; @c sigma @c = @c 0 degenerates to a
* Dirac component).
*
* @sa @c mixture, @c categorical, @c normal
* @sa Wikipedia: Mixture model
*/
CREATE OR REPLACE FUNCTION gmm(
weights double precision[],
means double precision[],
stddevs double precision[])
RETURNS random_variable AS
$$
DECLARE
n integer;
w_sum double precision := 0.0;
i integer;
acc random_variable := NULL;
remaining double precision := 0.0;
BEGIN
IF weights IS NULL OR means IS NULL OR stddevs IS NULL THEN
RAISE EXCEPTION 'provsql.gmm: weights, means, and stddevs must be non-null';
END IF;
n := array_length(weights, 1);
IF n IS NULL OR n < 1 THEN
RAISE EXCEPTION 'provsql.gmm: weights must be non-empty';
END IF;
IF array_length(means, 1) <> n OR array_length(stddevs, 1) <> n THEN
RAISE EXCEPTION 'provsql.gmm: weights, means, and stddevs must have the same length (got %, %, %)',
n, array_length(means, 1), array_length(stddevs, 1);
END IF;
FOR i IN 1..n LOOP
IF weights[i] IS NULL OR weights[i] = 'NaN'::float8
OR weights[i] < 0 OR weights[i] > 1 THEN
RAISE EXCEPTION 'provsql.gmm: weights[%] must be in [0,1] (got %)',
i, weights[i];
END IF;
w_sum := w_sum + weights[i];
END LOOP;
IF abs(w_sum - 1.0) > 1e-9 THEN
RAISE EXCEPTION 'provsql.gmm: weights must sum to 1 within 1e-9 (got %)', w_sum;
END IF;
-- Stick-breaking, built back to front: acc holds the mixture of
-- components i+1..n, and prepending component i selects it with
-- conditional probability w_i / (w_i + ... + w_n).
FOR i IN REVERSE n..1 LOOP
IF weights[i] <= 0.0 THEN
CONTINUE;
END IF;
IF acc IS NULL THEN
acc := provsql.normal(means[i], stddevs[i]);
remaining := weights[i];
ELSE
remaining := remaining + weights[i];
acc := provsql.mixture(least(1.0, weights[i] / remaining),
provsql.normal(means[i], stddevs[i]), acc);
END IF;
END LOOP;
RETURN acc;
END
$$ LANGUAGE plpgsql STRICT VOLATILE PARALLEL SAFE;
/**
* @brief Empirical-samples constructor: the ecdf of a sample bundle as
* a @c random_variable.
*
* Loads a Monte Carlo / MCMC / bootstrap sample array as the discrete
* distribution putting mass @c 1/n on each draw (duplicates merge, so a
* value drawn @c k times carries @c k/n) -- the standard empirical
* distribution. Reduces entirely to @ref categorical, so the whole
* exact discrete surface applies: moments are the sample moments,
* comparisons against constants are decided analytically ("fraction of
* samples below c"), and quantiles are the exact empirical quantiles.
*
* @code
* -- Bulk load via array_agg over a sample table
* INSERT INTO model_posteriors
* SELECT param, provsql.empirical_samples(array_agg(value))
* FROM mcmc_chain GROUP BY param;
* @endcode
*
* At most 10000 distinct values (the categorical block cap): thin the
* chain or bin the samples (e.g. with @c width_bucket) beyond that.
*
* @sa @ref categorical, @ref empirical_cdf
* @sa Wikipedia: Empirical distribution function
*/
CREATE OR REPLACE FUNCTION empirical_samples(samples double precision[])
RETURNS random_variable AS
$$
DECLARE
n integer;
sorted double precision[];
outcomes double precision[] := '{}';
probs double precision[] := '{}';
v double precision;
prev double precision;
run integer := 0;
started boolean := false;
BEGIN
n := array_length(samples, 1);
IF n IS NULL OR n < 1 THEN
RAISE EXCEPTION 'provsql.empirical_samples: samples must be non-empty';
END IF;
sorted := ARRAY(SELECT s FROM unnest(samples) AS s ORDER BY 1);
FOREACH v IN ARRAY sorted LOOP
IF v IS NULL OR v = 'NaN'::float8
OR v = 'Infinity'::float8 OR v = '-Infinity'::float8 THEN
RAISE EXCEPTION
'provsql.empirical_samples: samples must be finite (got %)', v;
END IF;
IF started AND v = prev THEN
run := run + 1;
ELSE
IF started THEN
outcomes := outcomes || prev;
probs := probs || (run::double precision / n);
END IF;
prev := v;
run := 1;
started := true;
END IF;
END LOOP;
outcomes := outcomes || prev;
probs := probs || (run::double precision / n);
IF array_length(outcomes, 1) > 10000 THEN
RAISE EXCEPTION
'provsql.empirical_samples: at most 10000 distinct values are '
'supported (got %); thin the chain or bin the samples (e.g. with '
'width_bucket)', array_length(outcomes, 1);
END IF;
RETURN provsql.categorical(probs, outcomes);
END
$$ LANGUAGE plpgsql STRICT VOLATILE PARALLEL SAFE;
/**
* @brief Empirical-CDF constructor: a piecewise-linear CDF table as a
* @c random_variable.
*
* Loads a tabulated CDF -- simulation output percentile tables, risk
* models, expert-elicited forecasts -- as the distribution whose CDF is
* @c cdf[i] at @c grid[i], linear in between: mass
* @c cdf[i+1] @c - @c cdf[i] spread uniformly over
* (grid[i], grid[i+1]), plus (when @c cdf[1] @c > @c 0) an
* atom of mass @c cdf[1] at @c grid[1] for the probability at or below
* the grid start. Packaged, like @ref gmm, as a stick-breaking cascade
* of Bernoulli @ref mixture nodes over @ref uniform components (and the
* optional @ref as_random atom), so moments and sampling are exact
* through the existing mixture machinery; comparisons ride Monte Carlo.
*
* @code
* provsql.empirical_cdf(
* grid => ARRAY[0.0, 0.5, 1.0, 2.0, 5.0, 10.0, 20.0],
* cdf => ARRAY[0.32, 0.51, 0.67, 0.82, 0.94, 0.99, 1.0])
* @endcode
*
* Validation: same-length arrays of at least two entries, @p grid
* strictly increasing and finite, @p cdf non-decreasing within
* [0, 1] and ending at @c 1 within @c 1e-9.
*
* @sa @ref gmm, @ref empirical_samples
* @sa Wikipedia: Cumulative distribution function
*/
CREATE OR REPLACE FUNCTION empirical_cdf(grid double precision[],
cdf double precision[])
RETURNS random_variable AS
$$
DECLARE
n integer;
i integer;
acc random_variable := NULL;
remaining double precision := 0.0;
w double precision;
comp random_variable;
BEGIN
n := array_length(grid, 1);
IF n IS NULL OR n < 2 THEN
RAISE EXCEPTION 'provsql.empirical_cdf: grid must have at least two entries';
END IF;
IF array_length(cdf, 1) <> n THEN
RAISE EXCEPTION 'provsql.empirical_cdf: grid and cdf must have the same length (got % and %)',
n, array_length(cdf, 1);
END IF;
IF n > 10000 THEN
RAISE EXCEPTION 'provsql.empirical_cdf: at most 10000 grid points are supported (got %)', n;
END IF;
FOR i IN 1..n LOOP
IF grid[i] IS NULL OR grid[i] = 'NaN'::float8
OR grid[i] = 'Infinity'::float8 OR grid[i] = '-Infinity'::float8 THEN
RAISE EXCEPTION 'provsql.empirical_cdf: grid[%] must be finite (got %)', i, grid[i];
END IF;
IF i > 1 AND NOT grid[i] > grid[i-1] THEN
RAISE EXCEPTION 'provsql.empirical_cdf: grid must be strictly increasing (grid[%] = %, grid[%] = %)',
i-1, grid[i-1], i, grid[i];
END IF;
IF cdf[i] IS NULL OR cdf[i] = 'NaN'::float8 OR cdf[i] < 0 OR cdf[i] > 1 THEN
RAISE EXCEPTION 'provsql.empirical_cdf: cdf[%] must be in [0,1] (got %)', i, cdf[i];
END IF;
IF i > 1 AND cdf[i] < cdf[i-1] THEN
RAISE EXCEPTION 'provsql.empirical_cdf: cdf must be non-decreasing (cdf[%] = %, cdf[%] = %)',
i-1, cdf[i-1], i, cdf[i];
END IF;
END LOOP;
IF abs(cdf[n] - 1.0) > 1e-9 THEN
RAISE EXCEPTION 'provsql.empirical_cdf: cdf must end at 1 within 1e-9 (got %)', cdf[n];
END IF;
-- Stick-breaking cascade, back to front: component i = 1 is the atom
-- at the grid start (mass cdf[1]); component i >= 2 is
-- uniform(grid[i-1], grid[i]) with mass cdf[i] - cdf[i-1].
FOR i IN REVERSE n..1 LOOP
w := CASE WHEN i = 1 THEN cdf[1] ELSE cdf[i] - cdf[i-1] END;
IF w <= 0.0 THEN
CONTINUE;
END IF;
comp := CASE WHEN i = 1 THEN provsql.as_random(grid[1])
ELSE provsql.uniform(grid[i-1], grid[i]) END;
IF acc IS NULL THEN
acc := comp;
remaining := w;
ELSE
remaining := remaining + w;
acc := provsql.mixture(least(1.0, w / remaining), comp, acc);
END IF;
END LOOP;
RETURN acc;
END
$$ LANGUAGE plpgsql STRICT VOLATILE PARALLEL SAFE;
/**
* @brief @c random_variable ^ @c random_variable (gate_arith POW).
*
* Real-valued branch only: evaluation raises if a negative base is
* drawn together with a non-integer exponent (write
* pow(greatest(x, 0), p) for the non-negative branch).
*/
CREATE OR REPLACE FUNCTION random_variable_pow(
a random_variable, b random_variable)
RETURNS random_variable AS
$$
SELECT provsql.random_variable_make(
provsql.provenance_arith(
7, -- PROVSQL_ARITH_POW
ARRAY[(a)::uuid,
(b)::uuid]));
$$ LANGUAGE sql IMMUTABLE STRICT PARALLEL SAFE;
/**
* @brief Natural logarithm of a @c random_variable (gate_arith LN).
*
* Defined on @c [0, +Infinity): evaluation raises if a negative value
* is drawn (restrict the argument's support); a draw of exactly @c 0
* yields @c -Infinity.
*/
CREATE OR REPLACE FUNCTION ln(a random_variable)
RETURNS random_variable AS
$$
SELECT provsql.random_variable_make(
provsql.provenance_arith(
8, -- PROVSQL_ARITH_LN
ARRAY[(a)::uuid]));
$$ LANGUAGE sql IMMUTABLE STRICT PARALLEL SAFE;
/** @brief @c e^x for a @c random_variable (gate_arith EXP). Total. */
CREATE OR REPLACE FUNCTION exp(a random_variable)
RETURNS random_variable AS
$$
SELECT provsql.random_variable_make(
provsql.provenance_arith(
9, -- PROVSQL_ARITH_EXP
ARRAY[(a)::uuid]));
$$ LANGUAGE sql IMMUTABLE STRICT PARALLEL SAFE;
/**
* @brief @c pow / @c power spellings of the @c ^ operator, mirroring
* PostgreSQL's numeric surface. Scalar exponents resolve
* through the implicit numeric-to-rv casts:
* pow(x, 0.5) is x ^ 0.5.
*/
CREATE OR REPLACE FUNCTION pow(a random_variable, b random_variable)
RETURNS random_variable AS
$$
SELECT provsql.random_variable_pow(a, b);
$$ LANGUAGE sql IMMUTABLE STRICT PARALLEL SAFE;
CREATE OR REPLACE FUNCTION power(a random_variable, b random_variable)
RETURNS random_variable AS
$$
SELECT provsql.random_variable_pow(a, b);
$$ LANGUAGE sql IMMUTABLE STRICT PARALLEL SAFE;
/**
* @brief Square root of a @c random_variable: sugar for
* x ^ 0.5 (no gate or opcode of its own). Evaluation
* raises on a negative draw, like any non-integer exponent.
*/
CREATE OR REPLACE FUNCTION sqrt(a random_variable)
RETURNS random_variable AS
$$
SELECT provsql.random_variable_pow(a, provsql.as_random(0.5));
$$ LANGUAGE sql IMMUTABLE STRICT PARALLEL SAFE;
DO $$ BEGIN
IF NOT EXISTS (
SELECT 1 FROM pg_operator o
JOIN pg_namespace n ON n.oid = o.oprnamespace
WHERE n.nspname = 'provsql' AND o.oprname = '^'
AND o.oprleft = 'provsql.random_variable'::regtype::oid
AND o.oprright = 'provsql.random_variable'::regtype::oid
AND o.oprcode <> 0
) THEN
CREATE OPERATOR ^ (
LEFTARG = random_variable,
RIGHTARG = random_variable,
PROCEDURE = random_variable_pow
);
END IF;
END $$;
/**
* @brief btree comparison support for @c random_variable -- always an error.
*
* A @c random_variable is a distribution, not a scalar, so it has no total
* order: sorting (@c ORDER @c BY), de-duplicating (@c DISTINCT), grouping, and
* the built-in @c GREATEST / @c LEAST all reduce to this btree comparison
* proc, which raises a clear diagnostic rather than a placeholder message.
*
* The proc exists only so a DEFAULT btree operator class can be declared for
* @c random_variable -- which is what lets PostgreSQL's @c GREATEST / @c LEAST
* grammar parse over random variables so the planner hook can lift it into a
* @c gate_arith @c MAX / @c MIN order statistic. When the hook is active the
* @c GREATEST / @c LEAST node is rewritten before it ever calls this proc.
*/
CREATE OR REPLACE FUNCTION random_variable_btree_cmp(
a random_variable, b random_variable) RETURNS integer AS
$$
BEGIN
RAISE EXCEPTION 'comparison or ordering of random_variable values is '
'meaningless: a random_variable is a distribution, not a scalar'
USING HINT =
'Compare them as a probabilistic event -- in a WHERE / JOIN clause or '
'with probability(x > y); take order statistics with provsql.greatest / '
'provsql.least (or the min / max aggregates); summarise numerically with '
'expected / variance / support.';
END
$$ LANGUAGE plpgsql IMMUTABLE STRICT PARALLEL SAFE;
DO $$ BEGIN
IF NOT EXISTS (
SELECT 1 FROM pg_opclass c
JOIN pg_namespace n ON n.oid = c.opcnamespace
JOIN pg_am am ON am.oid = c.opcmethod
WHERE n.nspname = 'provsql' AND c.opcname = 'random_variable_ops'
AND am.amname = 'btree'
) THEN
CREATE OPERATOR CLASS random_variable_ops
DEFAULT FOR TYPE random_variable USING btree AS
OPERATOR 1 <,
OPERATOR 2 <=,
OPERATOR 3 =,
OPERATOR 4 >=,
OPERATOR 5 >,
FUNCTION 1 random_variable_btree_cmp(random_variable, random_variable);
END IF;
END $$;
/**
* @name Order statistics over random_variable
*
* Same-row @c greatest / @c least over @c random_variable arguments: the
* order-statistic counterpart of the element-wise @c "+ - * /" operators.
* They lower to a single @c gate_arith with the @c MAX / @c MIN opcode over
* the argument circuits, the same n-ary shape the @c max / @c min aggregates
* build. Evaluation is Monte-Carlo-correct out of the box (@c std::max /
* @c std::min over the jointly-sampled children, so shared base RVs stay
* coupled); closed forms for i.i.d. families come from the analytic
* order-statistic pass.
*
* PostgreSQL's built-in @c GREATEST / @c LEAST are dedicated syntax (a
* @c MinMaxExpr requiring a btree comparison), not overloadable functions, so
* the surface is the schema-qualified @c provsql.greatest(...) /
* @c provsql.least(...). @c NULL arguments are ignored, matching the built-in
* (an all-@c NULL / empty call returns @c NULL).
* @{
*/
-- "greatest" / "least" are col_name keywords, so the CREATE FUNCTION name
-- must be quoted; callers reach them qualified as provsql.greatest(...).
-- Idempotence: max / min ignore repeats, so identical children (same gate)
-- are de-duplicated -- greatest(x, x, y) == greatest(x, y) -- and a single
-- surviving child collapses to itself -- greatest(x) == x. DISTINCT also sorts
-- the children, so the argument order does not matter for gate sharing. (Two
-- independent draws of the same distribution are distinct gates and are NOT
-- de-duplicated.)
CREATE OR REPLACE FUNCTION "greatest"(VARIADIC args random_variable[])
RETURNS random_variable AS
$$
DECLARE
children uuid[];
BEGIN
IF args IS NULL THEN
RETURN NULL;
END IF;
SELECT array_agg(DISTINCT (a)::uuid) INTO children
FROM unnest(args) a WHERE a IS NOT NULL;
IF children IS NULL OR array_length(children, 1) IS NULL THEN
RETURN NULL;
END IF;
IF array_length(children, 1) = 1 THEN
RETURN provsql.random_variable_make(children[1]);
END IF;
RETURN provsql.random_variable_make(
provsql.provenance_arith(5, children)); -- 5 = PROVSQL_ARITH_MAX
END
$$ LANGUAGE plpgsql IMMUTABLE PARALLEL SAFE;
CREATE OR REPLACE FUNCTION "least"(VARIADIC args random_variable[])
RETURNS random_variable AS
$$
DECLARE
children uuid[];
BEGIN
IF args IS NULL THEN
RETURN NULL;
END IF;
SELECT array_agg(DISTINCT (a)::uuid) INTO children
FROM unnest(args) a WHERE a IS NOT NULL;
IF children IS NULL OR array_length(children, 1) IS NULL THEN
RETURN NULL;
END IF;
IF array_length(children, 1) = 1 THEN
RETURN provsql.random_variable_make(children[1]);
END IF;
RETURN provsql.random_variable_make(
provsql.provenance_arith(6, children)); -- 6 = PROVSQL_ARITH_MIN
END
$$ LANGUAGE plpgsql IMMUTABLE PARALLEL SAFE;
/**
* @brief Build a @c random_variable from a guarded-selection @c gate_case.
*
* Thin @c random_variable wrapper over @c provenance_case (defined with the
* other gate builders, since it is uuid-only), the target of the planner-hook
* @c CASE-over-RV rewrite: the hook flattens the branches into
* @c [guard_1, value_1, ..., default] and emits this call so an RV-typed
* @c CASE surfaces as a first-class @c random_variable.
*/
CREATE OR REPLACE FUNCTION rv_case(
children UUID[]
)
RETURNS random_variable AS
$$
SELECT provsql.random_variable_make(provsql.provenance_case(children));
$$ LANGUAGE sql IMMUTABLE STRICT PARALLEL SAFE;
/**
* @brief Build an @c agg_token from a guarded-selection @c gate_case.
*
* The aggregate-carrier analogue of @c rv_case: a thin @c agg_token wrapper
* over the carrier-agnostic @c provenance_case, the target of the planner-hook
* lowering of a searched @c CASE whose guards are aggregate comparisons and
* whose branches are aggregates. The branches (and default) are already
* flattened into @c [guard_1, value_1, ..., default] UUIDs. The display cell
* carries the actual-world CASE value -- the branch selected on the actual
* data, resolved by @c agg_gate_value, exactly as a bare aggregate's cell
* carries its actual-world value. The probabilistic result is produced by
* the measure evaluators (``expected`` / ``probability`` / possible-worlds /
* Monte Carlo) from the gate, not the token's cell.
*/
CREATE OR REPLACE FUNCTION agg_case(
children UUID[]
)
RETURNS agg_token AS
$$
SELECT provsql.agg_token_make(t, coalesce(provsql.agg_gate_value(t), 0))
FROM (SELECT provsql.provenance_case(children) AS t) AS s;
$$ LANGUAGE sql IMMUTABLE STRICT PARALLEL SAFE;
/**
* @brief Identity-parameterised per-row wrap for an RV-returning aggregate.
*
* Generalises the two-argument @ref rv_aggregate_semimod: the else-branch
* (a row's contribution when its provenance is false) is
* @c as_random(@p identity) instead of the additive @c as_random(0). The
* planner-hook rewrite bakes each aggregate's own identity element into the
* wrap -- @c 1 for @c product, @f$-\infty@f$ / @f$+\infty@f$ for @c max /
* @c min -- so the aggregate's final function is a plain fold over the
* per-row mixtures with no gate inspection. @c sum keeps the two-argument
* form (@c identity @c = @c 0).
*/
CREATE OR REPLACE FUNCTION rv_aggregate_semimod(
prov uuid, rv random_variable, identity double precision)
RETURNS random_variable AS
$$
SELECT provsql.mixture(prov, rv, provsql.as_random(identity));
$$ LANGUAGE sql IMMUTABLE STRICT PARALLEL SAFE;
/**
* @brief Per-row denominator wrap for @c avg(random_variable): the
* provenance indicator @f$\mathbf{1}\{\varphi\}@f$.
*
* The row contributes @c 1 to the running count when present and @c 0 when
* absent, so @c sum over these wraps is the provenance-weighted count
* @f$\sum_i \mathbf{1}\{\varphi_i\}@f$. The planner-hook rewrites
* @c avg(x) into @c rv_sum_or_null(rv_aggregate_semimod(prov, x)) @c /
* @c sum(rv_aggregate_indicator(prov)) -- the "@c AVG @c = @c SUM @c /
* @c COUNT" identity lifted into the @c random_variable algebra -- so
* @c avg rides entirely on @c sum's fold and never inspects a gate.
*/
CREATE OR REPLACE FUNCTION rv_aggregate_indicator(prov uuid)
RETURNS random_variable AS
$$
SELECT provsql.rv_aggregate_semimod(prov, provsql.as_random(1::double precision));
$$ LANGUAGE sql IMMUTABLE STRICT PARALLEL SAFE;
/**
* @brief Numerator final function for the @c avg rewrite: @c sum, but
* @c NULL on an empty group.
*
* Identical to @ref sum_rv_ffunc except that an empty group returns
* @c NULL rather than the additive identity @c as_random(0). The
* planner-hook @c avg rewrite emits
* @c rv_sum_or_null(rv_aggregate_semimod(prov, x)) @c /
* @c sum(rv_aggregate_indicator(prov)); @c random_variable_div is
* @c STRICT, so an empty group propagates the numerator's @c NULL and
* @c avg is @c NULL -- the standard SQL @c AVG convention -- while a
* non-empty group behaves exactly like @c sum.
*/
CREATE OR REPLACE FUNCTION rv_sum_or_null_ffunc(state uuid[])
RETURNS random_variable AS
$$
BEGIN
IF state IS NULL OR array_length(state, 1) IS NULL THEN
RETURN NULL;
END IF;
IF array_length(state, 1) = 1 THEN
RETURN provsql.random_variable_make(state[1]);
END IF;
RETURN provsql.random_variable_make(
provsql.provenance_arith(0, state)); -- 0 = PROVSQL_ARITH_PLUS
END
$$ LANGUAGE plpgsql IMMUTABLE PARALLEL SAFE;
DO $$ BEGIN
IF NOT EXISTS (
SELECT 1 FROM pg_aggregate a
JOIN pg_proc p ON p.oid = a.aggfnoid
JOIN pg_namespace n ON n.oid = p.pronamespace
WHERE n.nspname = 'provsql' AND p.proname = 'rv_sum_or_null'
AND p.pronargs = 1
AND p.proargtypes[0] = 'provsql.random_variable'::regtype::oid
) THEN
CREATE AGGREGATE rv_sum_or_null(random_variable) (
SFUNC = sum_rv_sfunc,
STYPE = uuid[],
INITCOND = '{}',
FINALFUNC = rv_sum_or_null_ffunc
);
END IF;
END $$;
/**
* @brief Final function for @c avg(random_variable).
*
* @c avg lifts the "@c AVG @c = @c SUM @c / @c COUNT" identity into the
* @c random_variable algebra:
* @f[
* \mathrm{AVG}(x) \;=\; \frac{\sum_i \mathbf{1}\{\varphi_i\} \cdot X_i}
* {\sum_i \mathbf{1}\{\varphi_i\}}.
* @f]
* In a provenance-tracked query the planner-hook rewrites @c avg(x) into
* @c rv_sum_or_null(rv_aggregate_semimod(prov, x)) @c /
* @c sum(rv_aggregate_indicator(prov)) (see
* @c make_rv_aggregate_expression), so both the numerator and the
* provenance-weighted count denominator are built by @c sum's fold and no
* gate is inspected. This FFUNC is therefore reached only on an
* @em untracked call, where every row is unconditionally present: the
* numerator is @c sum over the raw per-row RVs and the denominator is the
* plain row count @c n (each row contributing @c as_random(1)).
*
* Empty group: returns @c NULL, matching standard SQL @c AVG (and unlike
* @c sum, whose empty group is the additive identity @c as_random(0)):
* the caller cannot otherwise disambiguate "0 rows" from "rows summing
* to 0".
*/
CREATE OR REPLACE FUNCTION avg_rv_ffunc(state uuid[])
RETURNS random_variable AS
$$
DECLARE
n integer;
i integer;
num_token uuid;
denom_token uuid;
denom_state uuid[] := '{}';
one_uuid uuid;
BEGIN
IF state IS NULL THEN
RETURN NULL;
END IF;
n := array_length(state, 1);
IF n IS NULL THEN
RETURN NULL;
END IF;
one_uuid := (provsql.as_random(1::double precision))::uuid;
FOR i IN 1..n LOOP
denom_state := array_append(denom_state, one_uuid);
END LOOP;
IF n = 1 THEN
num_token := state[1];
denom_token := denom_state[1];
ELSE
num_token := provsql.provenance_arith(0, state); -- 0 = PLUS
denom_token := provsql.provenance_arith(0, denom_state); -- 0 = PLUS
END IF;
RETURN provsql.random_variable_make(
provsql.provenance_arith(
3, -- 3 = PROVSQL_ARITH_DIV
ARRAY[num_token, denom_token]));
END
$$ LANGUAGE plpgsql IMMUTABLE PARALLEL SAFE;
/**
* @brief Final function for @c product(random_variable): fold a
* @c gate_arith TIMES root over the per-row contributions.
*
* Multiplicative analogue of @c sum(random_variable):
* @f[
* \mathrm{PRODUCT}(x) \;=\; \prod_i \big(\mathbf{1}\{\varphi_i\} \cdot X_i
* + \mathbf{1}\{\neg\varphi_i\} \cdot 1\big)
* \;=\; \prod_{i : \varphi_i} X_i.
* @f]
* Each per-row contribution already carries the multiplicative identity
* as its absent-row value: a provenance-tracked query wraps the argument
* as @c mixture(prov_i, X_i, as_random(1)) (identity baked in by the
* three-argument @ref rv_aggregate_semimod), and an untracked call passes
* the raw RV through. So the FFUNC is a plain fold with no gate
* inspection: @c gate_arith(TIMES, state).
*
* Reuses @c sum_rv_sfunc as the state-transition function. Empty group:
* the multiplicative identity @c as_random(1) -- the counterpart to
* @c sum's empty-group @c as_random(0). Singleton group: the single
* child directly, without a one-child TIMES root.
*/
CREATE OR REPLACE FUNCTION product_rv_ffunc(state uuid[])
RETURNS random_variable AS
$$
BEGIN
IF state IS NULL OR array_length(state, 1) IS NULL THEN
RETURN provsql.as_random(1::double precision);
END IF;
IF array_length(state, 1) = 1 THEN
RETURN provsql.random_variable_make(state[1]);
END IF;
RETURN provsql.random_variable_make(
provsql.provenance_arith(1, state)); -- 1 = PROVSQL_ARITH_TIMES
END
$$ LANGUAGE plpgsql IMMUTABLE PARALLEL SAFE;
/**
* @brief Final function for @c max(random_variable) / @c min(random_variable):
* fold a @c gate_arith MAX / MIN root over the per-row contributions.
*
* The order-statistic analogues of @c sum / @c product:
* @f[
* \mathrm{MAX}(x) = \max_{i : \varphi_i} X_i, \qquad
* \mathrm{MIN}(x) = \min_{i : \varphi_i} X_i.
* @f]
* A row absent in a world (its provenance @f$\varphi_i@f$ false) must not
* perturb the extremum, so it contributes the order-statistic identity
* @f$\mp\infty@f$. That identity is baked into each per-row contribution
* upstream: a provenance-tracked query wraps the argument as
* @c mixture(prov_i, X_i, as_random(∓∞)) (via the three-argument
* @ref rv_aggregate_semimod), and an untracked call passes the raw RV
* through. So the FFUNC is a plain fold with no gate inspection:
* @c gate_arith(@p op, state).
*
* Empty group: the identity @c as_random(@p identity) (@f$-\infty@f$ /
* @f$+\infty@f$), the extremum counterpart to @c sum's @c as_random(0).
* Singleton group: the single child directly.
*/
CREATE OR REPLACE FUNCTION extremum_rv_ffunc(
state uuid[], op integer, identity double precision)
RETURNS random_variable AS
$$
BEGIN
IF state IS NULL OR array_length(state, 1) IS NULL THEN
RETURN provsql.as_random(identity);
END IF;
IF array_length(state, 1) = 1 THEN
RETURN provsql.random_variable_make(state[1]);
END IF;
RETURN provsql.random_variable_make(
provsql.provenance_arith(op, state));
END
$$ LANGUAGE plpgsql IMMUTABLE PARALLEL SAFE;
CREATE OR REPLACE FUNCTION max_rv_ffunc(state uuid[])
RETURNS random_variable AS
$$
-- 5 = PROVSQL_ARITH_MAX; empty-group / row-absent identity -inf.
SELECT provsql.extremum_rv_ffunc(state, 5, '-Infinity'::double precision);
$$ LANGUAGE sql IMMUTABLE PARALLEL SAFE;
CREATE OR REPLACE FUNCTION min_rv_ffunc(state uuid[])
RETURNS random_variable AS
$$
-- 6 = PROVSQL_ARITH_MIN; empty-group / row-absent identity +inf.
SELECT provsql.extremum_rv_ffunc(state, 6, 'Infinity'::double precision);
$$ LANGUAGE sql IMMUTABLE PARALLEL SAFE;
DO $$ BEGIN
IF NOT EXISTS (
SELECT 1 FROM pg_aggregate a
JOIN pg_proc p ON p.oid = a.aggfnoid
JOIN pg_namespace n ON n.oid = p.pronamespace
WHERE n.nspname = 'provsql' AND p.proname = 'max'
AND p.pronargs = 1
AND p.proargtypes[0] = 'provsql.random_variable'::regtype::oid
) THEN
CREATE AGGREGATE max(random_variable) (
SFUNC = sum_rv_sfunc,
STYPE = uuid[],
INITCOND = '{}',
FINALFUNC = max_rv_ffunc
);
END IF;
END $$;
DO $$ BEGIN
IF NOT EXISTS (
SELECT 1 FROM pg_aggregate a
JOIN pg_proc p ON p.oid = a.aggfnoid
JOIN pg_namespace n ON n.oid = p.pronamespace
WHERE n.nspname = 'provsql' AND p.proname = 'min'
AND p.pronargs = 1
AND p.proargtypes[0] = 'provsql.random_variable'::regtype::oid
) THEN
CREATE AGGREGATE min(random_variable) (
SFUNC = sum_rv_sfunc,
STYPE = uuid[],
INITCOND = '{}',
FINALFUNC = min_rv_ffunc
);
END IF;
END $$;
-- ---------------------------------------------------------------------
-- SQL-standard statistic aggregates over random_variable rows:
-- covar_pop / covar_samp / corr (two-argument), stddev_pop / stddev_samp
-- (one-argument), and the ordered-set percentile_cont.
--
-- Row presence is carried by a per-row 0/1 indicator RV: the public
-- aggregates use the certain indicator as_random(1) (every row present),
-- and a provenance-tracked query is rewritten by the planner hook
-- (make_rv_aggregate_expression) to the rv_*_impl aggregates whose extra
-- leading argument is rv_aggregate_indicator(prov), so a row absent in a
-- world drops out of every sum, the count, and the percentile member set.
-- The moment statistics are built from indicator-weighted power sums with
-- existing gate_arith opcodes (e.g. covar_pop = SXY/N - (SX/N)(SY/N)); a
-- world where the statistic is undefined (N = 0, or N = 1 for the sample
-- forms) evaluates to NaN, the established undefined-world convention the
-- moment estimators skip. percentile_cont is the one gate the arithmetic
-- cannot express: it mints the PROVSQL_ARITH_PERCENTILE gate_arith
-- (interleaved [ind_1, x_1, ...] wires, fraction in extra) that the Monte
-- Carlo sampler evaluates by sorting each draw's present values and
-- interpolating.
-- ---------------------------------------------------------------------
/** @brief State transition for the one-argument RV statistic aggregates
* (@c stddev_pop / @c stddev_samp): append the certain indicator and the
* row's RV as a pair. NULL rows are skipped (standard SQL). */
CREATE OR REPLACE FUNCTION rv_stat1_sfunc(state uuid[], x random_variable)
RETURNS uuid[] AS
$$
SELECT CASE
WHEN x IS NULL THEN state
ELSE state || ARRAY[(provsql.as_random(1::double precision))::uuid,
(x)::uuid]
END;
$$ LANGUAGE sql IMMUTABLE PARALLEL SAFE;
/** @brief State transition for the two-argument RV statistic aggregates
* (@c covar_pop / @c covar_samp / @c corr): append the certain indicator
* and the row's RV pair as a triple. Rows with either side NULL are
* skipped (standard SQL covariance semantics). */
CREATE OR REPLACE FUNCTION rv_stat2_sfunc(
state uuid[], x random_variable, y random_variable)
RETURNS uuid[] AS
$$
SELECT CASE
WHEN x IS NULL OR y IS NULL THEN state
ELSE state || ARRAY[(provsql.as_random(1::double precision))::uuid,
(x)::uuid, (y)::uuid]
END;
$$ LANGUAGE sql IMMUTABLE PARALLEL SAFE;
/** @brief Indicator-carrying state transition for the one-argument
* @c rv_*_impl statistic aggregates: the planner-hook rewrite passes the
* row's provenance indicator @c rv_aggregate_indicator(prov) as @p ind. */
CREATE OR REPLACE FUNCTION rv_stat1_impl_sfunc(
state uuid[], ind random_variable, x random_variable)
RETURNS uuid[] AS
$$
SELECT CASE
WHEN x IS NULL THEN state
ELSE state || ARRAY[coalesce((ind)::uuid,
(provsql.as_random(1::double precision))::uuid),
(x)::uuid]
END;
$$ LANGUAGE sql IMMUTABLE PARALLEL SAFE;
/** @brief Indicator-carrying state transition for the two-argument
* @c rv_*_impl statistic aggregates. */
CREATE OR REPLACE FUNCTION rv_stat2_impl_sfunc(
state uuid[], ind random_variable, x random_variable, y random_variable)
RETURNS uuid[] AS
$$
SELECT CASE
WHEN x IS NULL OR y IS NULL THEN state
ELSE state || ARRAY[coalesce((ind)::uuid,
(provsql.as_random(1::double precision))::uuid),
(x)::uuid, (y)::uuid]
END;
$$ LANGUAGE sql IMMUTABLE PARALLEL SAFE;
/**
* @brief Mint the indicator-weighted power-sum gates shared by the
* covariance / stddev final functions.
*
* @p state is the flat interleaved aggregate state -- pairs
* @c [ind, x, ...] (@p stride 2) or triples @c [ind, x, y, ...]
* (@p stride 3). Emits @c gate_arith tokens for
* @f$N = \sum_i \mathbf{1}_i@f$, @f$SX = \sum_i \mathbf{1}_i x_i@f$,
* @f$SXX = \sum_i \mathbf{1}_i x_i^2@f$ and, at stride 3, @f$SY@f$,
* @f$SXY@f$, @f$SYY@f$. The per-row indicator gate is shared between
* @f$N@f$ and every product it weighs, so the Monte Carlo per-iteration
* cache keeps the row's presence coupled across all the sums (and a
* repeated child @c [ind, x, x] reuses the same draw of @c x, giving
* @f$x^2@f$, not two independent draws).
*/
CREATE OR REPLACE FUNCTION rv_stat_sum_tokens(
state uuid[], stride integer,
OUT n_tok uuid, OUT sx_tok uuid, OUT sxx_tok uuid,
OUT sy_tok uuid, OUT sxy_tok uuid, OUT syy_tok uuid)
AS
$$
DECLARE
nrows integer := coalesce(array_length(state, 1), 0) / stride;
inds uuid[] := '{}';
xs uuid[] := '{}';
xxs uuid[] := '{}';
ys uuid[] := '{}';
xys uuid[] := '{}';
yys uuid[] := '{}';
ind uuid;
x uuid;
y uuid;
BEGIN
FOR i IN 1..nrows LOOP
ind := state[(i-1) * stride + 1];
x := state[(i-1) * stride + 2];
inds := array_append(inds, ind);
xs := array_append(xs, provenance_arith(1, ARRAY[ind, x]));
xxs := array_append(xxs, provenance_arith(1, ARRAY[ind, x, x]));
IF stride = 3 THEN
y := state[(i-1) * stride + 3];
ys := array_append(ys, provenance_arith(1, ARRAY[ind, y]));
xys := array_append(xys, provenance_arith(1, ARRAY[ind, x, y]));
yys := array_append(yys, provenance_arith(1, ARRAY[ind, y, y]));
END IF;
END LOOP;
n_tok := provenance_arith(0, inds);
sx_tok := provenance_arith(0, xs);
sxx_tok := provenance_arith(0, xxs);
IF stride = 3 THEN
sy_tok := provenance_arith(0, ys);
sxy_tok := provenance_arith(0, xys);
syy_tok := provenance_arith(0, yys);
END IF;
END
$$ LANGUAGE plpgsql IMMUTABLE PARALLEL SAFE
SET search_path=provsql,pg_temp,public SECURITY DEFINER;
/** @brief Population-variance gate @f$SXX/N - (SX/N)^2@f$ from the
* power-sum tokens. */
CREATE OR REPLACE FUNCTION rv_stat_var_pop_token(
n_tok uuid, s_tok uuid, ss_tok uuid)
RETURNS uuid AS
$$
SELECT provsql.provenance_arith(2, ARRAY[
provsql.provenance_arith(3, ARRAY[ss_tok, n_tok]),
provsql.provenance_arith(1, ARRAY[
provsql.provenance_arith(3, ARRAY[s_tok, n_tok]),
provsql.provenance_arith(3, ARRAY[s_tok, n_tok])])]);
$$ LANGUAGE sql IMMUTABLE STRICT PARALLEL SAFE;
/** @brief Sample-variance gate @f$(SXX - SX^2/N) / (N - 1)@f$ from the
* power-sum tokens (NaN in a world with @f$N \le 1@f$, the undefined-world
* convention). */
CREATE OR REPLACE FUNCTION rv_stat_var_samp_token(
n_tok uuid, s_tok uuid, ss_tok uuid)
RETURNS uuid AS
$$
SELECT provsql.provenance_arith(3, ARRAY[
provsql.provenance_arith(2, ARRAY[
ss_tok,
provsql.provenance_arith(3, ARRAY[
provsql.provenance_arith(1, ARRAY[s_tok, s_tok]), n_tok])]),
provsql.provenance_arith(2, ARRAY[
n_tok, (provsql.as_random(1::double precision))::uuid])]);
$$ LANGUAGE sql IMMUTABLE STRICT PARALLEL SAFE;
/** @brief @f$\sqrt{\max(v, 0)}@f$ gate over a variance token: the max-clamp
* removes the tiny negative values float error can produce (variance is
* mathematically non-negative), so the POW domain guard never fires. */
CREATE OR REPLACE FUNCTION rv_stat_sqrt_token(v_tok uuid)
RETURNS uuid AS
$$
SELECT provsql.provenance_arith(7, ARRAY[
provsql.provenance_arith(5, ARRAY[
v_tok, (provsql.as_random(0::double precision))::uuid]),
(provsql.as_random(0.5::double precision))::uuid]);
$$ LANGUAGE sql IMMUTABLE STRICT PARALLEL SAFE;
/** @brief Population-covariance gate @f$SXY/N - (SX/N)(SY/N)@f$ from the
* power-sum tokens. */
CREATE OR REPLACE FUNCTION rv_stat_covar_pop_token(
n_tok uuid, sx_tok uuid, sy_tok uuid, sxy_tok uuid)
RETURNS uuid AS
$$
SELECT provsql.provenance_arith(2, ARRAY[
provsql.provenance_arith(3, ARRAY[sxy_tok, n_tok]),
provsql.provenance_arith(1, ARRAY[
provsql.provenance_arith(3, ARRAY[sx_tok, n_tok]),
provsql.provenance_arith(3, ARRAY[sy_tok, n_tok])])]);
$$ LANGUAGE sql IMMUTABLE STRICT PARALLEL SAFE;
/** @brief Final function for @c covar_pop(random_variable, random_variable). */
CREATE OR REPLACE FUNCTION covar_pop_rv_ffunc(state uuid[])
RETURNS random_variable AS
$$
DECLARE
t record;
BEGIN
IF state IS NULL OR array_length(state, 1) IS NULL THEN
RETURN NULL;
END IF;
SELECT * INTO t FROM rv_stat_sum_tokens(state, 3);
RETURN random_variable_make(
rv_stat_covar_pop_token(t.n_tok, t.sx_tok, t.sy_tok, t.sxy_tok));
END
$$ LANGUAGE plpgsql IMMUTABLE PARALLEL SAFE
SET search_path=provsql,pg_temp,public SECURITY DEFINER;
/** @brief Final function for @c covar_samp(random_variable, random_variable):
* @f$(SXY - SX\,SY/N) / (N-1)@f$. */
CREATE OR REPLACE FUNCTION covar_samp_rv_ffunc(state uuid[])
RETURNS random_variable AS
$$
DECLARE
t record;
BEGIN
IF state IS NULL OR array_length(state, 1) IS NULL THEN
RETURN NULL;
END IF;
SELECT * INTO t FROM rv_stat_sum_tokens(state, 3);
RETURN random_variable_make(
provenance_arith(3, ARRAY[
provenance_arith(2, ARRAY[
t.sxy_tok,
provenance_arith(3, ARRAY[
provenance_arith(1, ARRAY[t.sx_tok, t.sy_tok]), t.n_tok])]),
provenance_arith(2, ARRAY[
t.n_tok, (as_random(1::double precision))::uuid])]));
END
$$ LANGUAGE plpgsql IMMUTABLE PARALLEL SAFE
SET search_path=provsql,pg_temp,public SECURITY DEFINER;
/** @brief Final function for @c corr(random_variable, random_variable):
* @f$\mathrm{covar\_pop} / \sqrt{\max(v_x v_y, 0)}@f$ (a zero-variance
* world divides to @f$\pm\infty@f$ / NaN, the undefined-world convention,
* matching SQL's NULL for a zero-stddev input). */
CREATE OR REPLACE FUNCTION corr_rv_ffunc(state uuid[])
RETURNS random_variable AS
$$
DECLARE
t record;
vx uuid;
vy uuid;
BEGIN
IF state IS NULL OR array_length(state, 1) IS NULL THEN
RETURN NULL;
END IF;
SELECT * INTO t FROM rv_stat_sum_tokens(state, 3);
vx := rv_stat_var_pop_token(t.n_tok, t.sx_tok, t.sxx_tok);
vy := rv_stat_var_pop_token(t.n_tok, t.sy_tok, t.syy_tok);
RETURN random_variable_make(
provenance_arith(3, ARRAY[
rv_stat_covar_pop_token(t.n_tok, t.sx_tok, t.sy_tok, t.sxy_tok),
rv_stat_sqrt_token(provenance_arith(1, ARRAY[vx, vy]))]));
END
$$ LANGUAGE plpgsql IMMUTABLE PARALLEL SAFE
SET search_path=provsql,pg_temp,public SECURITY DEFINER;
/** @brief Final function for @c stddev_pop(random_variable). */
CREATE OR REPLACE FUNCTION stddev_pop_rv_ffunc(state uuid[])
RETURNS random_variable AS
$$
DECLARE
t record;
BEGIN
IF state IS NULL OR array_length(state, 1) IS NULL THEN
RETURN NULL;
END IF;
SELECT * INTO t FROM rv_stat_sum_tokens(state, 2);
RETURN random_variable_make(
rv_stat_sqrt_token(
rv_stat_var_pop_token(t.n_tok, t.sx_tok, t.sxx_tok)));
END
$$ LANGUAGE plpgsql IMMUTABLE PARALLEL SAFE
SET search_path=provsql,pg_temp,public SECURITY DEFINER;
/** @brief Final function for @c stddev_samp(random_variable). */
CREATE OR REPLACE FUNCTION stddev_samp_rv_ffunc(state uuid[])
RETURNS random_variable AS
$$
DECLARE
t record;
BEGIN
IF state IS NULL OR array_length(state, 1) IS NULL THEN
RETURN NULL;
END IF;
SELECT * INTO t FROM rv_stat_sum_tokens(state, 2);
RETURN random_variable_make(
rv_stat_sqrt_token(
rv_stat_var_samp_token(t.n_tok, t.sx_tok, t.sxx_tok)));
END
$$ LANGUAGE plpgsql IMMUTABLE PARALLEL SAFE
SET search_path=provsql,pg_temp,public SECURITY DEFINER;
DO $$ BEGIN
IF NOT EXISTS (
SELECT 1 FROM pg_aggregate a
JOIN pg_proc p ON p.oid = a.aggfnoid
JOIN pg_namespace n ON n.oid = p.pronamespace
WHERE n.nspname = 'provsql' AND p.proname = 'covar_pop'
AND p.pronargs = 2
AND p.proargtypes[0] = 'provsql.random_variable'::regtype::oid
AND p.proargtypes[1] = 'provsql.random_variable'::regtype::oid
) THEN
CREATE AGGREGATE covar_pop(random_variable, random_variable) (
SFUNC = rv_stat2_sfunc,
STYPE = uuid[],
INITCOND = '{}',
FINALFUNC = covar_pop_rv_ffunc
);
END IF;
END $$;
DO $$ BEGIN
IF NOT EXISTS (
SELECT 1 FROM pg_aggregate a
JOIN pg_proc p ON p.oid = a.aggfnoid
JOIN pg_namespace n ON n.oid = p.pronamespace
WHERE n.nspname = 'provsql' AND p.proname = 'covar_samp'
AND p.pronargs = 2
AND p.proargtypes[0] = 'provsql.random_variable'::regtype::oid
AND p.proargtypes[1] = 'provsql.random_variable'::regtype::oid
) THEN
CREATE AGGREGATE covar_samp(random_variable, random_variable) (
SFUNC = rv_stat2_sfunc,
STYPE = uuid[],
INITCOND = '{}',
FINALFUNC = covar_samp_rv_ffunc
);
END IF;
END $$;
DO $$ BEGIN
IF NOT EXISTS (
SELECT 1 FROM pg_aggregate a
JOIN pg_proc p ON p.oid = a.aggfnoid
JOIN pg_namespace n ON n.oid = p.pronamespace
WHERE n.nspname = 'provsql' AND p.proname = 'corr'
AND p.pronargs = 2
AND p.proargtypes[0] = 'provsql.random_variable'::regtype::oid
AND p.proargtypes[1] = 'provsql.random_variable'::regtype::oid
) THEN
CREATE AGGREGATE corr(random_variable, random_variable) (
SFUNC = rv_stat2_sfunc,
STYPE = uuid[],
INITCOND = '{}',
FINALFUNC = corr_rv_ffunc
);
END IF;
END $$;
DO $$ BEGIN
IF NOT EXISTS (
SELECT 1 FROM pg_aggregate a
JOIN pg_proc p ON p.oid = a.aggfnoid
JOIN pg_namespace n ON n.oid = p.pronamespace
WHERE n.nspname = 'provsql' AND p.proname = 'stddev_pop'
AND p.pronargs = 1
AND p.proargtypes[0] = 'provsql.random_variable'::regtype::oid
) THEN
CREATE AGGREGATE stddev_pop(random_variable) (
SFUNC = rv_stat1_sfunc,
STYPE = uuid[],
INITCOND = '{}',
FINALFUNC = stddev_pop_rv_ffunc
);
END IF;
END $$;
DO $$ BEGIN
IF NOT EXISTS (
SELECT 1 FROM pg_aggregate a
JOIN pg_proc p ON p.oid = a.aggfnoid
JOIN pg_namespace n ON n.oid = p.pronamespace
WHERE n.nspname = 'provsql' AND p.proname = 'stddev_samp'
AND p.pronargs = 1
AND p.proargtypes[0] = 'provsql.random_variable'::regtype::oid
) THEN
CREATE AGGREGATE stddev_samp(random_variable) (
SFUNC = rv_stat1_sfunc,
STYPE = uuid[],
INITCOND = '{}',
FINALFUNC = stddev_samp_rv_ffunc
);
END IF;
END $$;
DO $$ BEGIN
IF NOT EXISTS (
SELECT 1 FROM pg_aggregate a
JOIN pg_proc p ON p.oid = a.aggfnoid
JOIN pg_namespace n ON n.oid = p.pronamespace
WHERE n.nspname = 'provsql' AND p.proname = 'rv_covar_pop_impl'
AND p.pronargs = 3
AND p.proargtypes[0] = 'provsql.random_variable'::regtype::oid
AND p.proargtypes[1] = 'provsql.random_variable'::regtype::oid
AND p.proargtypes[2] = 'provsql.random_variable'::regtype::oid
) THEN
CREATE AGGREGATE rv_covar_pop_impl(
random_variable, random_variable, random_variable) (
SFUNC = rv_stat2_impl_sfunc,
STYPE = uuid[],
INITCOND = '{}',
FINALFUNC = covar_pop_rv_ffunc
);
END IF;
END $$;
DO $$ BEGIN
IF NOT EXISTS (
SELECT 1 FROM pg_aggregate a
JOIN pg_proc p ON p.oid = a.aggfnoid
JOIN pg_namespace n ON n.oid = p.pronamespace
WHERE n.nspname = 'provsql' AND p.proname = 'rv_covar_samp_impl'
AND p.pronargs = 3
AND p.proargtypes[0] = 'provsql.random_variable'::regtype::oid
AND p.proargtypes[1] = 'provsql.random_variable'::regtype::oid
AND p.proargtypes[2] = 'provsql.random_variable'::regtype::oid
) THEN
CREATE AGGREGATE rv_covar_samp_impl(
random_variable, random_variable, random_variable) (
SFUNC = rv_stat2_impl_sfunc,
STYPE = uuid[],
INITCOND = '{}',
FINALFUNC = covar_samp_rv_ffunc
);
END IF;
END $$;
DO $$ BEGIN
IF NOT EXISTS (
SELECT 1 FROM pg_aggregate a
JOIN pg_proc p ON p.oid = a.aggfnoid
JOIN pg_namespace n ON n.oid = p.pronamespace
WHERE n.nspname = 'provsql' AND p.proname = 'rv_corr_impl'
AND p.pronargs = 3
AND p.proargtypes[0] = 'provsql.random_variable'::regtype::oid
AND p.proargtypes[1] = 'provsql.random_variable'::regtype::oid
AND p.proargtypes[2] = 'provsql.random_variable'::regtype::oid
) THEN
CREATE AGGREGATE rv_corr_impl(
random_variable, random_variable, random_variable) (
SFUNC = rv_stat2_impl_sfunc,
STYPE = uuid[],
INITCOND = '{}',
FINALFUNC = corr_rv_ffunc
);
END IF;
END $$;
DO $$ BEGIN
IF NOT EXISTS (
SELECT 1 FROM pg_aggregate a
JOIN pg_proc p ON p.oid = a.aggfnoid
JOIN pg_namespace n ON n.oid = p.pronamespace
WHERE n.nspname = 'provsql' AND p.proname = 'rv_stddev_pop_impl'
AND p.pronargs = 2
AND p.proargtypes[0] = 'provsql.random_variable'::regtype::oid
AND p.proargtypes[1] = 'provsql.random_variable'::regtype::oid
) THEN
CREATE AGGREGATE rv_stddev_pop_impl(random_variable, random_variable) (
SFUNC = rv_stat1_impl_sfunc,
STYPE = uuid[],
INITCOND = '{}',
FINALFUNC = stddev_pop_rv_ffunc
);
END IF;
END $$;
DO $$ BEGIN
IF NOT EXISTS (
SELECT 1 FROM pg_aggregate a
JOIN pg_proc p ON p.oid = a.aggfnoid
JOIN pg_namespace n ON n.oid = p.pronamespace
WHERE n.nspname = 'provsql' AND p.proname = 'rv_stddev_samp_impl'
AND p.pronargs = 2
AND p.proargtypes[0] = 'provsql.random_variable'::regtype::oid
AND p.proargtypes[1] = 'provsql.random_variable'::regtype::oid
) THEN
CREATE AGGREGATE rv_stddev_samp_impl(random_variable, random_variable) (
SFUNC = rv_stat1_impl_sfunc,
STYPE = uuid[],
INITCOND = '{}',
FINALFUNC = stddev_samp_rv_ffunc
);
END IF;
END $$;
/**
* @brief Mint the @c PROVSQL_ARITH_PERCENTILE gate: the continuous
* percentile (SQL @c percentile_cont) over a group of RV rows.
*
* @p pairs is the interleaved wire list @c [ind_1, x_1, ..., ind_n, x_n]
* (each @p ind_i a 0/1 presence-indicator RV). The @p fraction is
* text-encoded in the gate's @c extra and participates in the token UUID
* (two percentiles of the same group at different fractions are distinct
* gates). Per Monte Carlo draw, the sampler collects the values whose
* indicator draws 1, sorts them, and linearly interpolates at the
* fraction; a draw with no present row is NaN (undefined world).
*/
CREATE OR REPLACE FUNCTION rv_percentile_make(fraction double precision,
pairs uuid[])
RETURNS random_variable AS
$$
DECLARE
token uuid;
BEGIN
IF fraction IS NULL THEN
RETURN NULL;
END IF;
IF fraction < 0 OR fraction > 1 THEN
RAISE EXCEPTION
'percentile_cont: fraction must be between 0 and 1 (got %)', fraction;
END IF;
token := public.uuid_generate_v5(
uuid_ns_provsql(),
concat('arith', '10', pairs::text, fraction::text));
PERFORM create_gate(token, 'arith', pairs);
PERFORM set_infos(token, 10); -- 10 = PROVSQL_ARITH_PERCENTILE
PERFORM set_extra(token, fraction::text);
RETURN random_variable_make(token);
END
$$ LANGUAGE plpgsql IMMUTABLE PARALLEL SAFE
SET search_path=provsql,pg_temp,public SECURITY DEFINER;
/** @brief State transition for the public ordered-set
* @c percentile_cont(float8) WITHIN GROUP (ORDER BY random_variable):
* append the certain indicator and the row's RV. Only reachable on
* untracked input (a provenance-tracked query is rewritten to
* @c rv_percentile_impl before planning), where the sort over
* @c random_variable raises the ordering-is-meaningless diagnostic
* first -- so in practice this runs only for empty input. */
CREATE OR REPLACE FUNCTION percentile_cont_rv_sfunc(
state uuid[], x random_variable)
RETURNS uuid[] AS
$$
SELECT provsql.rv_stat1_sfunc(state, x);
$$ LANGUAGE sql IMMUTABLE PARALLEL SAFE;
/** @brief Final function for the public ordered-set @c percentile_cont:
* receives the direct @p fraction argument after the state. */
CREATE OR REPLACE FUNCTION percentile_cont_rv_ffunc(
state uuid[], fraction double precision)
RETURNS random_variable AS
$$
SELECT CASE
WHEN state IS NULL OR array_length(state, 1) IS NULL THEN NULL
ELSE provsql.rv_percentile_make(fraction, state)
END;
$$ LANGUAGE sql IMMUTABLE PARALLEL SAFE;
DO $$ BEGIN
IF NOT EXISTS (
SELECT 1 FROM pg_aggregate a
JOIN pg_proc p ON p.oid = a.aggfnoid
JOIN pg_namespace n ON n.oid = p.pronamespace
WHERE n.nspname = 'provsql' AND p.proname = 'percentile_cont'
AND p.pronargs = 2
AND p.proargtypes[0] = 'pg_catalog.float8'::regtype::oid
AND p.proargtypes[1] = 'provsql.random_variable'::regtype::oid
) THEN
CREATE AGGREGATE percentile_cont(double precision ORDER BY random_variable) (
SFUNC = percentile_cont_rv_sfunc,
STYPE = uuid[],
INITCOND = '{}',
FINALFUNC = percentile_cont_rv_ffunc
);
END IF;
END $$;
DO $$ BEGIN
IF NOT EXISTS (
SELECT 1 FROM pg_type t
JOIN pg_namespace n ON n.oid = t.typnamespace
WHERE n.nspname = 'provsql' AND t.typname = 'rv_percentile_state'
) THEN
CREATE TYPE rv_percentile_state AS (
fraction double precision,
tokens uuid[]
);
END IF;
END $$;
/** @brief State transition for @c rv_percentile_impl, the planner-hook
* rewrite target of a provenance-tracked @c percentile_cont: stashes the
* (group-constant) fraction and appends the indicator/value pair. */
CREATE OR REPLACE FUNCTION rv_percentile_impl_sfunc(
state rv_percentile_state, fraction double precision,
ind random_variable, x random_variable)
RETURNS rv_percentile_state AS
$$
SELECT ROW(
coalesce((state).fraction, fraction),
CASE
WHEN x IS NULL THEN (state).tokens
ELSE (state).tokens ||
ARRAY[coalesce((ind)::uuid,
(provsql.as_random(1::double precision))::uuid),
(x)::uuid]
END)::provsql.rv_percentile_state;
$$ LANGUAGE sql IMMUTABLE PARALLEL SAFE;
/** @brief Final function for @c rv_percentile_impl. */
CREATE OR REPLACE FUNCTION rv_percentile_impl_ffunc(state rv_percentile_state)
RETURNS random_variable AS
$$
SELECT CASE
WHEN state IS NULL OR array_length((state).tokens, 1) IS NULL THEN NULL
ELSE provsql.rv_percentile_make((state).fraction, (state).tokens)
END;
$$ LANGUAGE sql IMMUTABLE PARALLEL SAFE;
DO $$ BEGIN
IF NOT EXISTS (
SELECT 1 FROM pg_aggregate a
JOIN pg_proc p ON p.oid = a.aggfnoid
JOIN pg_namespace n ON n.oid = p.pronamespace
WHERE n.nspname = 'provsql' AND p.proname = 'rv_percentile_impl'
AND p.pronargs = 3
AND p.proargtypes[0] = 'pg_catalog.float8'::regtype::oid
AND p.proargtypes[1] = 'provsql.random_variable'::regtype::oid
AND p.proargtypes[2] = 'provsql.random_variable'::regtype::oid
) THEN
CREATE AGGREGATE rv_percentile_impl(
double precision, random_variable, random_variable) (
SFUNC = rv_percentile_impl_sfunc,
STYPE = rv_percentile_state,
INITCOND = '(,"{}")',
FINALFUNC = rv_percentile_impl_ffunc
);
END IF;
END $$;
/**
* @brief Short alias of @c probability_evaluate.
*
* Bound to the same C symbol as @c probability_evaluate, so
* @c probability(token) is exactly @c probability_evaluate(token).
* Provided to match the concise polymorphic surface of @c expected,
* @c variance, and @c support: callers are not forced to spell out
* @c probability_evaluate.
*
* @param token provenance token to evaluate
* @param method knowledge compilation method (NULL for default)
* @param arguments additional arguments for the method
*/
CREATE OR REPLACE FUNCTION probability(
token UUID,
method text = NULL,
arguments text = NULL)
RETURNS DOUBLE PRECISION AS
'provsql','probability_evaluate' LANGUAGE C STABLE;
/**
* @brief Probability of a Boolean event over random variables.
*
* The @c (boolean) overload of @c probability lets a query ask for the
* probability of an event with the natural infix grammar, e.g.
* @c probability(x @c > @c y @c AND @c x @c < @c z). When the argument
* carries a probabilistic (random_variable / aggregate) comparison, the
* planner hook intercepts the call and rewrites it into
* @c probability_evaluate over the argument's event token (a @c gate_cmp /
* Boolean combination); the body below is then never reached.
*
* When the argument is a purely deterministic Boolean (no probabilistic
* comparison) the hook leaves the call alone and the body runs, so the
* probability of a definite event is simply @c 1 when it holds and @c 0 when
* it does not (@c NULL propagates). This makes @c probability total over
* Booleans -- @c probability(1 @c > @c 0) is @c 1, @c probability(region @c =
* @c 'north') is a per-row @c 0/1 -- and it works even with
* @c provsql.active off. @c NOT strict so a default-NULL @c method does not
* short-circuit the cast.
*
* The predicate surface deliberately lives only on the short @c probability
* name, not on @c probability_evaluate: a Boolean overload of the latter
* would make @c probability_evaluate('') ambiguous (an unknown
* literal matches both the @c uuid and the @c boolean overload), breaking
* existing string-literal callers. @c probability is new, so it carries the
* predicate overload without that hazard.
*/
CREATE OR REPLACE FUNCTION probability(
predicate boolean,
method text = NULL,
arguments text = NULL)
RETURNS DOUBLE PRECISION AS
$$
SELECT predicate::integer::double precision;
$$ LANGUAGE sql IMMUTABLE PARALLEL SAFE;
/** @brief Exact E[AVG^k | COUNT >= 1] over independent rows (the joint
* (sum, count) fold); NULL when the shape is out of scope (shared
* leaves, compound contributors), signalling @c agg_raw_moment's avg
* arm to fall back to the Monte-Carlo scalar path. */
CREATE OR REPLACE FUNCTION agg_avg_moment_exact(token uuid, k integer)
RETURNS double precision
AS 'provsql','agg_avg_moment_exact' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
/**
* @brief Boolean event "this aggregate-carrying gate's value is defined
* (non-NULL) in the world".
*
* Backs the conditional-on-defined convention of the aggregate moment
* readouts: @c sum / @c count (and constants) have a value in every
* world -- the empty group is the real value @c 0 -- so their defined
* event is @c gate_one(); @c min / @c max / @c avg (and any other
* aggregate) are @c NULL on an empty group, so their defined event is
* "some contributing row is present", the OR of the semimod children's
* row tokens; a @c case gate's value is defined iff its first-match
* selected branch's value is (the same region walk as the moment
* evaluator, conjoined per branch). Anything else (an @c arith
* composite, whose agg_token running value is total) counts as always
* defined.
*/
CREATE OR REPLACE FUNCTION agg_defined_event(token uuid)
RETURNS uuid AS $$
DECLARE
gt provenance_gate := get_gate_type(token);
fname varchar;
toks uuid[];
wires uuid[];
nw integer;
m integer;
i integer;
running_neg uuid := gate_one();
parts uuid[] := '{}';
BEGIN
IF gt = 'agg' THEN
SELECT proname INTO fname
FROM pg_proc WHERE oid = (get_infos(token)).info1;
IF fname IN ('sum', 'count') THEN
RETURN gate_one();
END IF;
SELECT array_agg((get_children(c))[1]) INTO toks
FROM unnest(get_children(token)) AS c;
IF toks IS NULL THEN
RETURN gate_zero(); -- structurally empty aggregate: never defined
END IF;
RETURN provenance_plus(toks);
ELSIF gt = 'case' THEN
wires := get_children(token);
nw := array_length(wires, 1);
m := (nw - 1) / 2;
FOR i IN 1..m LOOP
parts := parts || provenance_times(
running_neg, wires[2 * i - 1],
agg_defined_event(wires[2 * i]));
running_neg := provenance_times(running_neg,
provenance_not(wires[2 * i - 1]));
END LOOP;
parts := parts || provenance_times(running_neg,
agg_defined_event(wires[nw]));
RETURN provenance_plus(parts);
END IF;
-- value / arith / anything else: a value exists in every world.
RETURN gate_one();
END
$$ LANGUAGE plpgsql STABLE STRICT PARALLEL SAFE
SET search_path=provsql,pg_temp,public SECURITY DEFINER;
/**
* @brief Compute the raw moment E[X^k | prov] of an agg_token aggregate
*
* Sister of @c expected() for the agg_token side of the polymorphic
* @c moment / @c variance / @c central_moment dispatch. Supports the
* same aggregation functions as @c expected: SUM (which COUNT
* normalises to at the gate level via @c Aggregation.cpp:322), MIN,
* MAX, and AVG (exact over independent / laminar rows via the joint
* (sum, count) distribution, Monte-Carlo scalar fallback otherwise).
* MIN / MAX / AVG are NULL on an empty group, so their moments are
* CONDITIONAL on the aggregate being defined -- NULL only when it never
* is; SUM / COUNT treat the empty world as the real value 0.
*
* Strategy:
* - SUM: with X = Σᵢ Iᵢ·vᵢ (Iᵢ the per-row inclusion indicator,
* vᵢ the row's value), expanding X^k and taking expectation gives
* @f$E[X^k] = \sum_{(i_1,\ldots,i_k) \in \{1..n\}^k} v_{i_1}\cdots v_{i_k}
* \cdot P(\bigwedge_{i \in \text{distinct}(i_1..i_k)} I_i)@f$.
* We enumerate the @f$n^k@f$ tuples, conjoin the distinct inclusion
* tokens (and @p prov when conditioning), and evaluate the
* probability via @c probability_evaluate.
* - MIN / MAX: replace @c v with @c v^k in the rank-based
* enumeration that @c expected already uses; @c MAX is handled by
* sign-flipping per the existing trick (negate vs. rerank), with
* the outer multiplier becoming @f$(-1)^k@f$ instead of just @f$-1@f$.
*
* Cost: SUM is @f$O(n^k)@f$ probability evaluations -- tractable for
* small @p k or small @p n; for larger sizes, prefer reaching for the
* sampler. MIN / MAX stay linear in @p n.
*/
CREATE OR REPLACE FUNCTION agg_raw_moment(
token agg_token,
k integer,
prov UUID = gate_one(),
method text = NULL,
arguments text = NULL)
RETURNS DOUBLE PRECISION AS $$
DECLARE
aggregation_function VARCHAR;
child_pairs uuid[];
pair_children uuid[];
n integer;
i integer;
j integer;
vals float8[];
toks uuid[];
total float8;
total_probability float8;
tup integer[];
d integer;
prod_v float8;
distinct_tok uuid[];
conj_token uuid;
prob float8;
sign_max float8;
BEGIN
IF token IS NULL OR k IS NULL THEN
RETURN NULL;
END IF;
IF k < 0 THEN
RAISE EXCEPTION 'agg_raw_moment(): k must be non-negative (got %)', k;
END IF;
-- Aggregate-carrier CASE (a gate_case over aggregate branches): a first-match
-- guarded selection. The moment is CONDITIONAL on the CASE's value being
-- defined (NULL only when it never is, mirroring the MIN/MAX convention):
-- E[pick^k | defined ∧ prov]
-- = Σ_i P(region_i ∧ def_i) · E[value_i^k | region_i ∧ def_i]
-- / Σ_i P(region_i ∧ def_i),
-- where region_i = (¬g_1 ∧ … ∧ ¬g_{i-1}) ∧ g_i ∧ prov is the world set that
-- selects branch i (the default's region is "all guards false") and def_i is
-- the branch's defined event (agg_defined_event: gate_one for sum / count /
-- constants, "some row present" for min / max / avg, recursive for a nested
-- CASE). Both factors are exact: probability() over the region ∧ def event,
-- and the conditional aggregate moment (a recursive agg_raw_moment on the
-- branch aggregate, which conditions on its own definedness within the
-- region, so the two factors weigh the same worlds). The regions are
-- mutually exclusive, so the terms sum with no inclusion-exclusion, and
-- correlation between a guard and its branch (shared input tuples) is
-- carried by the conditioning, exactly as HAVING carries it. When every
-- branch is defined everywhere, the defined mass equals P(prov) and the
-- formula reduces to the plain region-weighted sum.
IF get_gate_type(token) = 'case' THEN
IF k = 0 THEN
RETURN 1;
END IF;
DECLARE
wires uuid[] := get_children(token);
nw integer := array_length(get_children(token), 1);
m integer := (array_length(get_children(token), 1) - 1) / 2;
running_neg uuid := gate_one();
region_full uuid;
prov_p float8;
p float8;
total float8 := 0;
def_mass float8 := 0;
ci integer;
vuid uuid;
bm float8;
BEGIN
prov_p := probability(prov);
IF prov_p IS NULL OR prov_p <= 0 THEN
RETURN NULL; -- impossible conditioning event
END IF;
-- Branches 1..m are the guarded WHENs; branch m+1 is the ELSE default,
-- whose region is "all guards false".
FOR ci IN 1 .. m + 1 LOOP
IF ci <= m THEN
region_full := provenance_times(running_neg, wires[2 * ci - 1], prov);
vuid := wires[2 * ci];
running_neg :=
provenance_times(running_neg, provenance_not(wires[2 * ci - 1]));
ELSE
region_full := provenance_times(running_neg, prov);
vuid := wires[nw];
END IF;
p := probability(provenance_times(region_full,
agg_defined_event(vuid)));
IF p > 0 THEN
-- E[value_i^k | region_i ∧ def_i]: a constant branch is a Dirac
-- (c^k, exact); a single aggregate or nested CASE is exact via
-- agg_raw_moment (whose MIN/MAX/CASE arms condition on their own
-- definedness within the region); an arithmetic / composite branch
-- takes the Monte-Carlo scalar path (which composes with the
-- aggregate leaves).
IF get_gate_type(vuid) = 'value' THEN
bm := power(CAST(get_extra(vuid) AS float8), k);
ELSIF get_gate_type(vuid) IN ('agg', 'case') THEN
bm := agg_raw_moment(agg_token_make(vuid, 0), k, region_full,
method, arguments);
ELSE
bm := rv_moment(vuid, k, false, region_full);
END IF;
total := total + p * bm;
def_mass := def_mass + p;
END IF;
END LOOP;
IF def_mass <= epsilon() THEN
RETURN NULL; -- the CASE's value is never defined under prov
END IF;
RETURN total / def_mass;
END;
END IF;
IF get_gate_type(token) <> 'agg' THEN
IF get_gate_type(token) IN ('arith', 'conditioned') THEN
RAISE EXCEPTION 'expected / variance / moment over an arithmetic '
'combination of aggregates (e.g. SUM(x) + SUM(y) or SUM(x) + 5), or a '
'conditioning of one, is not yet supported: a moment can be taken only '
'over a single aggregate (SUM / COUNT / MIN / MAX), optionally '
'conditioned (SUM(x) | C)'
USING HINT = 'Take the moment of each aggregate separately, or condition '
'the bare aggregate.';
ELSE
RAISE EXCEPTION USING MESSAGE='Wrong gate type for agg_raw_moment computation';
END IF;
END IF;
IF k = 0 THEN
RETURN 1;
END IF;
SELECT pp.proname::varchar FROM pg_proc pp
WHERE oid=(get_infos(token)).info1
INTO aggregation_function;
child_pairs := get_children(token);
n := COALESCE(array_length(child_pairs, 1), 0);
IF aggregation_function = 'sum' OR aggregation_function = 'count' THEN
-- count(col) keeps the COUNT identity at the gate level but its value is a
-- SUM of per-row 0/1 indicators, so its moments are computed exactly like
-- SUM (and its empty group is the real value 0, like SUM). count(*)
-- arrives here as 'sum' (it normalises to F_SUM_INT4); count(col) as 'count'.
-- Trivial empty aggregation: SUM = 0, so SUM^k = 0 for k >= 1.
-- Note: agg_token semantics treat the "no row included" world as
-- SUM = 0, so this stays consistent with k = 1 (= expected()).
IF n = 0 THEN
RETURN 0;
END IF;
-- Extract per-child token + value arrays.
vals := ARRAY[]::float8[];
toks := ARRAY[]::uuid[];
FOR i IN 1..n LOOP
pair_children := get_children(child_pairs[i]);
toks := toks || pair_children[1];
vals := vals || CAST(get_extra(pair_children[2]) AS float8);
END LOOP;
-- Enumerate all k-tuples (i_1, ..., i_k) in {1..n}^k. tup is the
-- current tuple; we step through them in lexicographic order.
total := 0;
tup := array_fill(1, ARRAY[k]);
LOOP
prod_v := 1;
FOR j IN 1..k LOOP
prod_v := prod_v * vals[tup[j]];
END LOOP;
SELECT array_agg(DISTINCT toks[idx]) INTO distinct_tok
FROM unnest(tup) AS idx;
IF prov <> gate_one() THEN
distinct_tok := distinct_tok || prov;
END IF;
conj_token := provenance_times(VARIADIC distinct_tok);
prob := probability_evaluate(conj_token, method, arguments);
total := total + prod_v * prob;
d := k;
WHILE d >= 1 AND tup[d] = n LOOP
tup[d] := 1;
d := d - 1;
END LOOP;
EXIT WHEN d = 0;
tup[d] := tup[d] + 1;
END LOOP;
ELSIF aggregation_function = 'min' OR aggregation_function = 'max' THEN
-- Rank enumeration: per distinct value v, P(MIN = v) is the
-- probability that some t_i with v_i=v is true and all t_j with
-- smaller v are false. For MAX we negate values so the same
-- "smaller-than" rank logic computes MIN-of-negated, then flip.
-- The outer multiplier picks up the right sign for the k-th moment
-- of MAX: E[MAX^k] = (-1)^k * E[MIN(-v)^k], so sign_max = (-1)^k.
sign_max := CASE
WHEN aggregation_function = 'max'
THEN power(-1::float8, k)
ELSE 1
END;
-- MIN/MAX over the empty input world are NULL (no elements), not ±Infinity:
-- SQL returns one row with a NULL value. The moment is therefore CONDITIONAL
-- on the aggregate being defined (non-empty) -- the empty world is excluded
-- and the result renormalised by P(prov AND non-empty). (count, whose empty
-- value 0 is a real value, keeps the empty world; sum keeps it too, as 0.)
IF n = 0 THEN
RETURN NULL; -- structurally empty: MIN/MAX undefined
END IF;
-- Numerator E[MIN^k . 1{prov AND non-empty}] (the rank sum naturally omits
-- the empty world, since every term requires a present token).
WITH tok_value AS (
SELECT (get_children(c))[1] AS tok,
(CASE WHEN aggregation_function='max' THEN -1 ELSE 1 END)
* CAST(get_extra((get_children(c))[2]) AS DOUBLE PRECISION) AS v
FROM UNNEST(child_pairs) AS c
) SELECT sign_max * COALESCE(SUM(p * power(v, k)), 0) FROM (
SELECT t1.v AS v,
probability_evaluate(
CASE WHEN prov = gate_one()
THEN provenance_monus(provenance_plus(ARRAY_AGG(t1.tok)),
provenance_plus(ARRAY_AGG(t2.tok)))
ELSE provenance_times(prov,
provenance_monus(provenance_plus(ARRAY_AGG(t1.tok)),
provenance_plus(ARRAY_AGG(t2.tok)))) END,
method, arguments) AS p
FROM tok_value t1 LEFT OUTER JOIN tok_value t2 ON t1.v > t2.v
GROUP BY t1.v) tmp
INTO total;
-- Denominator P(prov AND non-empty) = P(prov (x) (+) tokens).
SELECT probability_evaluate(
CASE WHEN prov = gate_one()
THEN provenance_plus(ARRAY_AGG(tok))
ELSE provenance_times(prov, provenance_plus(ARRAY_AGG(tok))) END,
method, arguments)
FROM (SELECT (get_children(c))[1] AS tok FROM UNNEST(child_pairs) AS c) s
INTO total_probability;
IF total_probability <= epsilon() THEN
RETURN NULL; -- never defined under prov: MIN/MAX undefined
END IF;
RETURN total / total_probability; -- already conditional; skip generic norm
ELSIF aggregation_function = 'avg' THEN
-- AVG = SUM/COUNT is a ratio of two correlated world-dependent
-- quantities, so the k-tuple expansion above does not apply. Like
-- MIN/MAX, AVG over the empty world is NULL, so its moment conditions
-- on the aggregate being defined (COUNT >= 1), NULL when it never is.
-- Two routes:
-- * EXACT (independent rows, unconditional): the joint (sum, count)
-- PMF folded in C by agg_avg_moment_exact --
-- E[AVG^k | COUNT>=1] = Σ_{(s,c), c>=1} (s/c)^k pmf(s,c) / P(c>=1).
-- * Monte-Carlo scalar fallback otherwise (an outer conditioning
-- event, shared leaves, compound contributors): rv_moment samples
-- the agg gate per world; its NaN-skip on empty draws implements
-- the same conditional-on-defined convention, at the
-- provsql.rv_mc_samples budget (0 raises, per convention).
IF n = 0 THEN
RETURN NULL; -- structurally empty: AVG undefined
END IF;
IF prov = gate_one() THEN
total := agg_avg_moment_exact((token)::uuid, k);
IF total IS NOT NULL THEN
RETURN total;
END IF;
END IF;
RETURN rv_moment((token)::uuid, k, false, prov);
ELSE
RAISE EXCEPTION USING MESSAGE=
'Cannot compute moment for aggregation function ' || aggregation_function;
END IF;
-- Conditional normalisation: E[X^k · 1_A] / P(A) = E[X^k | A].
IF prov <> gate_one()
AND total <> 0
AND total <> 'Infinity'::float8
AND total <> '-Infinity'::float8 THEN
total := total / probability_evaluate(prov, method, arguments);
END IF;
RETURN total;
END
$$ LANGUAGE plpgsql PARALLEL SAFE SET search_path=provsql SECURITY DEFINER;
/**
* @brief Internal: rv-side quantile computation.
*
* C entry point behind the polymorphic @c quantile dispatcher.
* Closed-form inverse CDF where the family has one (Normal via
* Beasley-Springer-Moro polished by Newton steps, Uniform and
* Exponential by algebraic inversion), generic monotone-CDF bisection
* otherwise (Erlang, Gamma), exact generalised inverse for categorical
* mixtures, and the empirical Monte Carlo quantile for compound scalar
* circuits. A non-trivial @p prov conditions (truncates) the
* distribution first, in closed form when the event reduces to an
* interval on a bare @c gate_rv.
*/
CREATE OR REPLACE FUNCTION rv_quantile(
token uuid, p double precision,
prov uuid DEFAULT gate_one())
RETURNS double precision
AS 'provsql','rv_quantile' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
/**
* @brief Compute the p-quantile (inverse CDF) of a probabilistic scalar
*
* @f$F^{-1}(p) = \min\{x : P(X \le x) \ge p\}@f$ for @f$p \in [0,1]@f$:
* medians (@c p = 0.5), percentiles, Value-at-Risk, and credible
* intervals. @c p = 0 / @c p = 1 return the (possibly infinite)
* support edges. Polymorphic dispatcher mirroring @c expected /
* @c moment: @c random_variable routes through @c rv_quantile
* (analytical inverse CDF / MC), plain numerics are their own quantile
* (a Dirac's inverse CDF is constant), and the optional @p prov
* argument conditions on a provenance event, e.g.
* quantile(x | (x > 0), 0.5) for the median of a truncated
* distribution.
*/
CREATE OR REPLACE FUNCTION quantile(
input ANYELEMENT,
p double precision,
prov UUID = gate_one(),
method text = NULL,
arguments text = NULL)
RETURNS DOUBLE PRECISION AS $$
BEGIN
IF p IS NULL THEN
RETURN NULL;
END IF;
IF p <> p OR p < 0 OR p > 1 THEN
RAISE EXCEPTION 'quantile: p must be in [0, 1] (got %)', p;
END IF;
IF pg_typeof(input) = 'random_variable'::regtype THEN
IF input IS NULL THEN
RETURN NULL;
END IF;
-- See variance(): rv_quantile handles the conditional/unconditional
-- dispatch internally based on the resolved prov gate type.
RETURN provsql.rv_quantile(
rv_conditioned_target((input::random_variable)::uuid), p,
rv_conditioned_prov((input::random_variable)::uuid, prov));
END IF;
IF pg_typeof(input) IN ('smallint'::regtype, 'integer'::regtype,
'bigint'::regtype, 'numeric'::regtype,
'real'::regtype, 'double precision'::regtype) THEN
-- A deterministic scalar is a Dirac: every quantile is the value.
RETURN input::double precision;
END IF;
RAISE EXCEPTION 'quantile() is not yet supported for input type %', pg_typeof(input);
END
$$ LANGUAGE plpgsql PARALLEL SAFE SET search_path=provsql SECURITY DEFINER;
/** @brief C entry point behind @ref covariance (uuid-level binding). */
CREATE OR REPLACE FUNCTION rv_covariance(x uuid, y uuid, prov uuid)
RETURNS double precision
AS 'provsql','rv_covariance' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
/**
* @brief Covariance Cov(X, Y) = E[XY] − E[X]·E[Y] of two random variables.
*
* The bivariate readout complementing the univariate moment surface
* (@ref expected / @ref variance / @ref moment / @ref central_moment).
* Exact tiers: an exact @c 0 when the two arguments' stochastic-leaf
* footprints are structurally independent (given @p prov), a variance
* readout when the two arguments coincide, and the closed-form
* @c E[XY] − E[X]·E[Y] whenever every factor decomposes analytically.
* When some factor has no closed form, a SINGLE coupled Monte-Carlo pass
* over the joint circuit draws @c (x, y) pairs (shared leaves produce one
* draw both observe) and returns the sample covariance -- the estimator's
* noise then scales with the covariance signal itself, not with the
* product of the means as the naive three-run E[XY] − E[X]·E[Y]
* subtraction would.
*
* @param prov optional conditioning event (a provenance @c uuid); the
* default @c gate_one() is the unconditional covariance. Conditioning
* is applied jointly: the Monte-Carlo pass rejection-samples the pair on
* @p prov, giving @c Cov(X, Y | prov).
*/
CREATE OR REPLACE FUNCTION covariance(
x random_variable, y random_variable, prov uuid DEFAULT gate_one())
RETURNS double precision AS $$
SELECT provsql.rv_covariance((x)::uuid, (y)::uuid, prov);
$$ LANGUAGE sql PARALLEL SAFE STABLE SET search_path=provsql SECURITY DEFINER;
/**
* @brief Standard deviation σ(X) = √Var(X) of a random variable.
*
* A thin numeric readout over @ref variance: the square root is taken on
* the scalar @c double result, so no RV-level @c sqrt is involved and this
* carries no dependency on RV function application (@c pow / @c sqrt).
* @c NULL propagates from a @c NULL input; the order-2 central moment is
* non-negative by construction, so the root is always real.
*
* @param prov optional conditioning event; default @c gate_one()
* (unconditional).
*/
CREATE OR REPLACE FUNCTION stddev(
x random_variable, prov uuid DEFAULT gate_one())
RETURNS double precision AS $$
SELECT sqrt(provsql.variance(x, prov));
$$ LANGUAGE sql PARALLEL SAFE STABLE SET search_path=provsql SECURITY DEFINER;
/** @brief C entry point behind @ref correlation (uuid-level binding). */
CREATE OR REPLACE FUNCTION rv_correlation(x uuid, y uuid, prov uuid)
RETURNS double precision
AS 'provsql','rv_correlation' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
/**
* @brief Pearson correlation ρ(X, Y) = Cov(X, Y) / (σ(X)·σ(Y)).
*
* Same exact tiers as @ref covariance; on the Monte-Carlo path the
* covariance and BOTH standard deviations are read off the same coupled
* pass, instead of stacking five independent estimates (three for the
* covariance, one per standard deviation). Returns @c NULL when either
* standard deviation is @c 0 (a degenerate / constant variable, for which
* correlation is undefined) rather than raising a division-by-zero.
*
* @param prov optional conditioning event; default @c gate_one()
* (unconditional).
*/
CREATE OR REPLACE FUNCTION correlation(
x random_variable, y random_variable, prov uuid DEFAULT gate_one())
RETURNS double precision AS $$
SELECT provsql.rv_correlation((x)::uuid, (y)::uuid, prov);
$$ LANGUAGE sql PARALLEL SAFE STABLE SET search_path=provsql SECURITY DEFINER;
/** @brief C entry point behind @ref entropy (uuid-level binding). */
CREATE OR REPLACE FUNCTION rv_entropy(token uuid, prov uuid)
RETURNS double precision
AS 'provsql','rv_entropy' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
/**
* @brief Entropy H(X) of a random variable, in nats.
*
* Shannon entropy for a discrete distribution (a categorical / discrete
* count / constant -- a point mass has entropy @c 0), differential
* entropy for a continuous one (quadrature of @c -f ln f over the
* family's integration range; also exact through independent-arm
* Bernoulli mixture trees such as @ref gmm's). Shapes with no
* closed density (arithmetic composites) and the conditional form fall
* back to a Monte Carlo histogram plug-in estimate at the
* @c provsql.rv_mc_samples budget.
*
* @param prov optional conditioning event; default @c gate_one()
* (unconditional).
*/
CREATE OR REPLACE FUNCTION entropy(
x random_variable, prov uuid DEFAULT gate_one())
RETURNS double precision AS $$
SELECT provsql.rv_entropy((x)::uuid, prov);
$$ LANGUAGE sql PARALLEL SAFE STABLE SET search_path=provsql SECURITY DEFINER;
/** @brief C entry point behind @ref kl (uuid-level binding). */
CREATE OR REPLACE FUNCTION rv_kl(p uuid, q uuid)
RETURNS double precision
AS 'provsql','rv_kl' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
/**
* @brief Kullback-Leibler divergence KL(P || Q), in nats.
*
* Exact: the defining sum for two discrete distributions (matching
* outcomes by value) and the defining integral (quadrature over P's
* integration window) for two continuous ones, including
* independent-arm mixture trees. Returns @c Infinity when P is not
* absolutely continuous with respect to Q -- an outcome of P that Q
* gives zero mass, mismatched kinds (discrete vs continuous), or a
* region of P's support where Q's density (under)flows to zero. Both
* arguments must resolve to closed-form densities; arithmetic
* composites and conditioned variables raise.
*/
CREATE OR REPLACE FUNCTION kl(p random_variable, q random_variable)
RETURNS double precision AS $$
SELECT provsql.rv_kl((p)::uuid, (q)::uuid);
$$ LANGUAGE sql PARALLEL SAFE STABLE SET search_path=provsql SECURITY DEFINER;
/** @brief C entry point behind @ref mutual_information (uuid-level
* binding). */
CREATE OR REPLACE FUNCTION rv_mutual_information(x uuid, y uuid)
RETURNS double precision
AS 'provsql','rv_mutual_information' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
/**
* @brief Mutual information I(X; Y), in nats.
*
* Exactly @c 0 for structurally independent variables (disjoint
* stochastic-leaf footprints, the same test the moment evaluators use);
* @c H(X) for a discrete variable paired with itself and @c Infinity
* for a continuous one (I(X;X) diverges). A genuinely correlated pair
* (shared leaves) is estimated by a 2-D histogram plug-in over coupled
* joint Monte Carlo draws -- both roots evaluated against the same
* per-iteration cache, so shared leaves keep their joint law -- at the
* @c provsql.rv_mc_samples budget.
*/
CREATE OR REPLACE FUNCTION mutual_information(
x random_variable, y random_variable)
RETURNS double precision AS $$
SELECT provsql.rv_mutual_information((x)::uuid, (y)::uuid);
$$ LANGUAGE sql PARALLEL SAFE STABLE SET search_path=provsql SECURITY DEFINER;
-- ---------------------------------------------------------------------
-- NULL-semantics surface: the deprecated aggregation_evaluate driver is
-- retired, the core combinators gain three-valued-logic NULL handling,
-- avg gains its value-aware presence indicator, and the explicit
-- zero-filtering predicates are added. Bodies are transcribed verbatim
-- from provsql.common.sql so md5(prosrc) matches a fresh install.
-- ---------------------------------------------------------------------
-- Retire both aggregation_evaluate overloads (a NULL-token crash vector).
-- The C symbol stays in provsql.so as an always-NULL stub for the 1.0.0
-- fixture; only the SQL surface is dropped.
DROP FUNCTION IF EXISTS aggregation_evaluate(uuid, regclass, regproc, regproc, regproc, anyelement, regproc, regproc, regproc, regproc);
DROP FUNCTION IF EXISTS aggregation_evaluate(uuid, regclass, regproc, regproc, regproc, anyelement, regtype, regproc, regproc, regproc, regproc);
CREATE OR REPLACE FUNCTION provenance_times(VARIADIC tokens uuid[])
RETURNS UUID AS
$$
DECLARE
times_token uuid;
filtered_tokens uuid[];
canonical uuid;
BEGIN
-- A NULL element reads as the ⊗-neutral 1: it is the token slot of an
-- untracked source (a join against an untracked table), which is
-- certain. Contrast provenance_plus / provenance_monus, where NULL
-- reads as the ⊕- / ⊖-right-neutral 0: each combinator maps NULL to
-- its own neutral element. Nothing may therefore hand a NULL to ⊗
-- meaning "false"; a comparison with a NULL operand goes through
-- provenance_cmp, which returns gate_zero for it.
SELECT array_agg(t) FROM unnest(tokens) t WHERE t IS NOT NULL AND t <> gate_one() INTO filtered_tokens;
-- Dispatch on the FILTERED count: a single survivor short-circuits
-- to that token directly (no useless single-child times gate); zero
-- survivors collapse to the identity. Using array_length(tokens, 1)
-- here would miss the [one, cmp] → [cmp] case, leaving the cmp wrapped
-- in a one-child times when its only sibling was gate_one().
CASE coalesce(array_length(filtered_tokens, 1), 0)
WHEN 0 THEN
times_token:=gate_one();
WHEN 1 THEN
times_token:=filtered_tokens[1];
ELSE
-- Computed separately from the filtering aggregate above: an
-- ORDER BY aggregate there would make the planner feed *both*
-- aggregates sorted input, scrambling the stored children order.
SELECT uuid_generate_v5(uuid_ns_provsql(),
concat('times-canonical', array_agg(t ORDER BY t)))
FROM unnest(filtered_tokens) t
INTO canonical;
IF get_gate_type(canonical) = 'times' THEN
-- A deliberate pre-creation at the canonical address: same
-- children, same product.
times_token := canonical;
ELSE
times_token := uuid_generate_v5(uuid_ns_provsql(),concat('times',filtered_tokens));
PERFORM create_gate(times_token, 'times', ARRAY_AGG(t)) FROM UNNEST(filtered_tokens) AS t WHERE t IS NOT NULL;
END IF;
END CASE;
RETURN times_token;
END
$$ LANGUAGE plpgsql SET search_path=provsql,pg_temp,public SECURITY DEFINER PARALLEL SAFE;
CREATE OR REPLACE FUNCTION provenance_plus(tokens uuid[])
RETURNS UUID AS
$$
DECLARE
c INTEGER;
plus_token uuid;
filtered_tokens uuid[];
canonical uuid;
BEGIN
-- A NULL element reads as the ⊕-neutral 0: it stands for a row absent
-- from the disjunction (a null-padded antijoin row whose token array
-- slot is NULL), not for an untracked source. Contrast provenance_times,
-- where NULL reads as the ⊗-neutral 1 (untracked source): each
-- combinator maps NULL to its own neutral element.
SELECT array_agg(t) FROM unnest(tokens) t
WHERE t IS NOT NULL AND t <> gate_zero()
INTO filtered_tokens;
c:=array_length(filtered_tokens, 1);
IF c = 0 THEN
plus_token := gate_zero();
ELSIF c = 1 THEN
plus_token := filtered_tokens[1];
ELSE
-- Computed separately from the filtering aggregate above: an ORDER
-- BY aggregate there would make the planner feed *both* aggregates
-- sorted input, scrambling the stored (aggregation-order) children.
SELECT uuid_generate_v5(uuid_ns_provsql(),
concat('plus-canonical', array_agg(t ORDER BY t)))
FROM unnest(filtered_tokens) t
INTO canonical;
IF get_gate_type(canonical) = 'plus' THEN
-- A deliberate pre-creation at the canonical address: same
-- children, same sum.
plus_token := canonical;
ELSE
plus_token := uuid_generate_v5(
uuid_ns_provsql(),
concat('plus', filtered_tokens));
PERFORM create_gate(plus_token, 'plus', filtered_tokens);
END IF;
END IF;
RETURN plus_token;
END
$$ LANGUAGE plpgsql STRICT SET search_path=provsql,pg_temp,public SECURITY DEFINER PARALLEL SAFE IMMUTABLE;
CREATE OR REPLACE FUNCTION provenance_monus(token1 UUID, token2 UUID)
RETURNS UUID AS
$$
DECLARE
monus_token uuid;
BEGIN
IF token1 IS NULL THEN
RAISE EXCEPTION USING MESSAGE='provenance_monus is called with first argument NULL';
END IF;
IF token2 IS NULL THEN
-- The ⊖-right-neutral 0: a NULL second argument is the no-match case
-- of the difference operator's LEFT OUTER JOIN (nothing to subtract),
-- so X ⊖ NULL = X ⊖ 0 = X. Note this is NOT the NULL ≡ 1 reading of
-- provenance_times; each combinator maps NULL to its own neutral.
RETURN token1;
END IF;
IF token1 = token2 THEN
-- X-X=0
monus_token:=gate_zero();
ELSIF token1 = gate_zero() THEN
-- 0-X=0
monus_token:=gate_zero();
ELSIF token2 = gate_zero() THEN
-- X-0=X
monus_token:=token1;
ELSE
monus_token:=uuid_generate_v5(uuid_ns_provsql(),concat('monus',token1,token2));
PERFORM create_gate(monus_token, 'monus', ARRAY[token1::uuid, token2::uuid]);
END IF;
RETURN monus_token;
END
$$ LANGUAGE plpgsql SET search_path=provsql,pg_temp,public SECURITY DEFINER PARALLEL SAFE IMMUTABLE;
CREATE OR REPLACE FUNCTION provenance_cmp(
left_token UUID,
comparison_op OID,
right_token UUID
)
RETURNS UUID AS
$$
DECLARE
cmp_token UUID;
BEGIN
-- A comparison with a NULL operand (a NULL random_variable cell, or an
-- aggregate that is NULL on the instance) is unknown under SQL's 3VL in
-- every possible world: the row is annotated zero. The function must
-- not be STRICT: a NULL result would read as the neutral token
-- (provenance_times drops it), silently turning "unknown" into
-- "certainly true".
IF left_token IS NULL OR right_token IS NULL OR comparison_op IS NULL THEN
RETURN gate_zero();
END IF;
-- deterministic v5 namespace id
cmp_token := public.uuid_generate_v5(
uuid_ns_provsql(),
concat('cmp', left_token::text, comparison_op::text, right_token::text)
);
-- wire it up in the circuit
PERFORM create_gate(cmp_token, 'cmp', ARRAY[left_token, right_token]);
PERFORM set_infos(cmp_token, comparison_op::integer);
RETURN cmp_token;
END
$$ LANGUAGE plpgsql
SET search_path=provsql,pg_temp,public
SECURITY DEFINER
IMMUTABLE
PARALLEL SAFE;
CREATE OR REPLACE FUNCTION provenance_delta
(token UUID)
RETURNS UUID AS
$$
DECLARE
delta_token uuid;
BEGIN
-- NULL token ≡ 1 (untracked source), and δ(1) = 1. Tested first: the
-- equality comparisons below are not NULL-safe.
IF token IS NULL THEN
return gate_one();
END IF;
IF token = gate_zero() OR token = gate_one() THEN
return token;
END IF;
delta_token:=uuid_generate_v5(uuid_ns_provsql(),concat('delta',token));
PERFORM create_gate(delta_token,'delta',ARRAY[token::uuid]);
RETURN delta_token;
END
$$ LANGUAGE plpgsql SET search_path=provsql,pg_temp,public SECURITY DEFINER PARALLEL SAFE;
CREATE OR REPLACE FUNCTION rv_aggregate_indicator(prov uuid, rv random_variable)
RETURNS random_variable AS
$$
SELECT CASE WHEN rv IS NULL THEN NULL
ELSE provsql.rv_aggregate_indicator(prov) END;
$$ LANGUAGE sql IMMUTABLE PARALLEL SAFE;
CREATE OR REPLACE FUNCTION true_nonzero(token uuid)
RETURNS boolean AS
'provsql', 'true_nonzero' LANGUAGE C PARALLEL SAFE STABLE;
CREATE OR REPLACE FUNCTION nonzero(token uuid,
semiring text DEFAULT NULL,
mapping regclass DEFAULT NULL)
RETURNS boolean AS
$$
BEGIN
IF token IS NULL THEN
RETURN true;
END IF;
IF semiring IS NULL THEN
RETURN provsql.true_nonzero(token);
ELSIF semiring = 'boolean' THEN
RETURN provsql.provenance_evaluate_compiled(token, mapping, 'boolean', TRUE);
ELSIF semiring = 'counting' THEN
RETURN provsql.provenance_evaluate_compiled(token, mapping, 'counting', 1) <> 0;
ELSE
RAISE EXCEPTION 'nonzero: unsupported semiring "%" (supported: boolean, counting; NULL for the universal zero test)', semiring;
END IF;
END
$$ LANGUAGE plpgsql PARALLEL SAFE STABLE;
CREATE OR REPLACE FUNCTION present(token uuid)
RETURNS boolean AS
$$
SELECT provsql.nonzero(token, 'boolean');
$$ LANGUAGE sql PARALLEL SAFE STABLE;
-- A backend warmed under 1.10.0 caches InvalidOid for the new 'case' enum
-- value; force a fresh lookup on the next get_constants() call.
-- ----------------------------------------------------------------------
-- Latent random variables + likelihood-weighting posterior inference.
--
-- This surface (distribution constructors taking random_variable
-- parameters, rv_parametric1/2, the evidence / observe / given
-- conditioning constructors and the and_agg fold, shapley_observe, the
-- agg_token -> random_variable bridge, and the collapsed-moment / variance
-- readouts) was added to the base SQL after 1.10.0 but had not been
-- replicated here. Functions use CREATE OR REPLACE (idempotent); the
-- non-replaceable objects (the 'observe' enum value near the top, the
-- agg_token -> random_variable cast, the and_agg aggregate, and the unary
-- | given operators) are each guarded so the script stays idempotent and
-- PostgreSQL 10/11-safe (no CREATE OR REPLACE AGGREGATE, no bare CREATE).
-- ----------------------------------------------------------------------
-- given_predicate(boolean) was folded into the given(boolean) overload;
-- drop the stale prefix operator (it points at given_predicate) and the
-- function so the upgraded catalog matches a fresh install. The
-- given(boolean) form and its | operator are (re)created further down.
DROP OPERATOR IF EXISTS | (NONE, boolean);
DROP FUNCTION IF EXISTS given_predicate(boolean);
CREATE OR REPLACE FUNCTION provenance_times(VARIADIC tokens uuid[])
RETURNS UUID AS
$$
DECLARE
times_token uuid;
filtered_tokens uuid[];
canonical uuid;
BEGIN
-- A NULL element reads as the ⊗-neutral 1: it is the token slot of an
-- untracked source (a join against an untracked table), which is
-- certain. Contrast provenance_plus / provenance_monus, where NULL
-- reads as the ⊕- / ⊖-right-neutral 0: each combinator maps NULL to
-- its own neutral element. Nothing may therefore hand a NULL to ⊗
-- meaning "false"; a comparison with a NULL operand goes through
-- provenance_cmp, which returns gate_zero for it.
SELECT array_agg(t) FROM unnest(tokens) t WHERE t IS NOT NULL AND t <> gate_one() INTO filtered_tokens;
-- Dispatch on the FILTERED count: a single survivor short-circuits
-- to that token directly (no useless single-child times gate); zero
-- survivors collapse to the identity. Using array_length(tokens, 1)
-- here would miss the [one, cmp] → [cmp] case, leaving the cmp wrapped
-- in a one-child times when its only sibling was gate_one().
CASE coalesce(array_length(filtered_tokens, 1), 0)
WHEN 0 THEN
times_token:=gate_one();
WHEN 1 THEN
times_token:=filtered_tokens[1];
ELSE
-- Computed separately from the filtering aggregate above: an
-- ORDER BY aggregate there would make the planner feed *both*
-- aggregates sorted input, scrambling the stored children order.
SELECT uuid_generate_v5(uuid_ns_provsql(),
concat('times-canonical', array_agg(t ORDER BY t)))
FROM unnest(filtered_tokens) t
INTO canonical;
IF get_gate_type(canonical) = 'times' THEN
-- A deliberate pre-creation at the canonical address: same
-- children, same product.
times_token := canonical;
ELSE
times_token := uuid_generate_v5(uuid_ns_provsql(),concat('times',filtered_tokens));
PERFORM create_gate(times_token, 'times', ARRAY_AGG(t)) FROM UNNEST(filtered_tokens) AS t WHERE t IS NOT NULL;
END IF;
END CASE;
RETURN times_token;
END
$$ LANGUAGE plpgsql SET search_path=provsql,pg_temp,public SECURITY DEFINER PARALLEL SAFE IMMUTABLE;
CREATE OR REPLACE FUNCTION provenance_delta
(token UUID)
RETURNS UUID AS
$$
DECLARE
delta_token uuid;
BEGIN
-- NULL token ≡ 1 (untracked source), and δ(1) = 1. Tested first: the
-- equality comparisons below are not NULL-safe.
IF token IS NULL THEN
return gate_one();
END IF;
IF token = gate_zero() OR token = gate_one() THEN
return token;
END IF;
delta_token:=uuid_generate_v5(uuid_ns_provsql(),concat('delta',token));
PERFORM create_gate(delta_token,'delta',ARRAY[token::uuid]);
RETURN delta_token;
END
$$ LANGUAGE plpgsql SET search_path=provsql,pg_temp,public SECURITY DEFINER PARALLEL SAFE IMMUTABLE;
CREATE OR REPLACE FUNCTION provenance_aggregate(
aggfnoid integer,
aggtype integer,
val anyelement,
tokens uuid[],
is_scalar boolean DEFAULT false)
RETURNS agg_token AS
$$
DECLARE
c INTEGER;
agg_tok uuid;
agg_val varchar;
BEGIN
-- Drop the NULL placeholders array_agg keeps for rows that did not produce a
-- semimod gate (provenance_semimod returns NULL for a NULL aggregated value),
-- so a NULL input never participates in the aggregate.
tokens := array_remove(tokens, NULL);
c:=COALESCE(array_length(tokens, 1), 0);
agg_val = CAST(val as VARCHAR);
IF c = 0 THEN
agg_tok := gate_zero();
ELSE
-- aggfnoid must be part of the UUID: SUM(id) and AVG(id) over the
-- same children would otherwise collapse to a single gate, and
-- their concurrent set_infos calls would overwrite each other's
-- aggregation operator (resulting in the wrong agg_kind being
-- read by provsql_having under cross-backend contention). The
-- scalar-aggregation flag must likewise be hashed: a scalar and a
-- grouped aggregate over identical children carry different info2 and
-- must stay distinct gates, else the concurrent set_infos calls would
-- clobber the flag. The flag is stored in the high bit of info2 (the
-- low 31 bits keep the result-type OID); aggtype itself is passed clean
-- so the agg_token->scalar cast still finds a valid type.
agg_tok := uuid_generate_v5(
uuid_ns_provsql(),
concat('agg',aggfnoid,tokens,CASE WHEN is_scalar THEN 'S' ELSE '' END));
PERFORM create_gate(agg_tok, 'agg', tokens);
PERFORM set_infos(agg_tok, aggfnoid,
CASE WHEN is_scalar THEN aggtype | (-2147483648) ELSE aggtype END);
PERFORM set_extra(agg_tok, agg_val);
END IF;
RETURN '( '||agg_tok||' , '||agg_val||' )';
END
$$ LANGUAGE plpgsql PARALLEL SAFE SET search_path=provsql,pg_temp,public SECURITY DEFINER IMMUTABLE;
CREATE OR REPLACE FUNCTION provenance_semimod(val anyelement, token UUID)
RETURNS UUID AS
$$
DECLARE
semimod_token uuid;
value_token uuid;
BEGIN
-- A NULL value means this row does not participate in the aggregate (SQL
-- aggregates ignore NULL inputs; only count(*) counts rows unconditionally,
-- and it passes a constant 1 here). Produce no semimod gate so the row is
-- skipped when provenance_aggregate builds the agg gate.
IF val IS NULL THEN
RETURN NULL;
END IF;
SELECT uuid_generate_v5(uuid_ns_provsql(),concat('value',CAST(val AS VARCHAR)))
INTO value_token;
SELECT uuid_generate_v5(uuid_ns_provsql(),concat('semimod',value_token,token))
INTO semimod_token;
--create value gates
PERFORM create_gate(value_token,'value');
PERFORM set_extra(value_token, CAST(val AS VARCHAR));
--create semimod gate
PERFORM create_gate(semimod_token,'semimod',ARRAY[token::uuid,value_token]);
RETURN semimod_token;
END
$$ LANGUAGE plpgsql PARALLEL SAFE SET search_path=provsql,pg_temp,public SECURITY DEFINER IMMUTABLE;
CREATE OR REPLACE FUNCTION agg_token_to_random_variable(a agg_token)
RETURNS random_variable AS
$$ SELECT provsql.random_variable_make(provsql.agg_token_uuid($1)); $$
LANGUAGE sql IMMUTABLE STRICT PARALLEL SAFE;
DO $$ BEGIN
IF NOT EXISTS (
SELECT 1 FROM pg_cast c
JOIN pg_type s ON s.oid = c.castsource
JOIN pg_type t ON t.oid = c.casttarget
WHERE s.typname = 'agg_token' AND t.typname = 'random_variable'
) THEN
CREATE CAST (agg_token AS random_variable)
WITH FUNCTION agg_token_to_random_variable(agg_token) AS IMPLICIT;
END IF;
END $$;
CREATE OR REPLACE FUNCTION rv_parametric2(
family text,
p1_tok uuid, p1_lit double precision,
p2_tok uuid, p2_lit double precision)
RETURNS random_variable AS
$$
DECLARE
token uuid;
wires uuid[] := ARRAY[]::uuid[];
s1 text;
s2 text;
BEGIN
IF p1_tok IS NOT NULL THEN
wires := wires || p1_tok;
s1 := '$' || (array_length(wires, 1) - 1);
ELSE
IF NOT provsql.is_finite_float8(p1_lit) THEN
RAISE EXCEPTION 'provsql.%: literal parameter must be finite (got %)',
family, p1_lit;
END IF;
s1 := p1_lit::text;
END IF;
IF p2_tok IS NOT NULL THEN
wires := wires || p2_tok;
s2 := '$' || (array_length(wires, 1) - 1);
ELSE
IF NOT provsql.is_finite_float8(p2_lit) THEN
RAISE EXCEPTION 'provsql.%: literal parameter must be finite (got %)',
family, p2_lit;
END IF;
s2 := p2_lit::text;
END IF;
token := public.uuid_generate_v4();
PERFORM provsql.create_gate(token, 'rv', wires);
PERFORM provsql.set_extra(token, family || ':' || s1 || ',' || s2);
RETURN provsql.random_variable_make(token);
END
$$ LANGUAGE plpgsql VOLATILE PARALLEL SAFE;
CREATE OR REPLACE FUNCTION rv_parametric1(family text, p_tok uuid)
RETURNS random_variable AS
$$
DECLARE
token uuid;
BEGIN
token := public.uuid_generate_v4();
PERFORM provsql.create_gate(token, 'rv', ARRAY[p_tok]);
PERFORM provsql.set_extra(token, family || ':$0');
RETURN provsql.random_variable_make(token);
END
$$ LANGUAGE plpgsql STRICT VOLATILE PARALLEL SAFE;
CREATE OR REPLACE FUNCTION normal(mu random_variable, sigma double precision)
RETURNS random_variable AS $$ SELECT provsql.rv_parametric2('normal', ($1)::uuid, NULL, NULL, $2); $$
LANGUAGE sql STRICT VOLATILE PARALLEL SAFE;
CREATE OR REPLACE FUNCTION normal(mu double precision, sigma random_variable)
RETURNS random_variable AS $$ SELECT provsql.rv_parametric2('normal', NULL, $1, ($2)::uuid, NULL); $$
LANGUAGE sql STRICT VOLATILE PARALLEL SAFE;
CREATE OR REPLACE FUNCTION normal(mu random_variable, sigma random_variable)
RETURNS random_variable AS $$ SELECT provsql.rv_parametric2('normal', ($1)::uuid, NULL, ($2)::uuid, NULL); $$
LANGUAGE sql STRICT VOLATILE PARALLEL SAFE;
CREATE OR REPLACE FUNCTION logistic(mu random_variable, s double precision)
RETURNS random_variable AS $$ SELECT provsql.rv_parametric2('logistic', ($1)::uuid, NULL, NULL, $2); $$
LANGUAGE sql STRICT VOLATILE PARALLEL SAFE;
CREATE OR REPLACE FUNCTION logistic(mu double precision, s random_variable)
RETURNS random_variable AS $$ SELECT provsql.rv_parametric2('logistic', NULL, $1, ($2)::uuid, NULL); $$
LANGUAGE sql STRICT VOLATILE PARALLEL SAFE;
CREATE OR REPLACE FUNCTION logistic(mu random_variable, s random_variable)
RETURNS random_variable AS $$ SELECT provsql.rv_parametric2('logistic', ($1)::uuid, NULL, ($2)::uuid, NULL); $$
LANGUAGE sql STRICT VOLATILE PARALLEL SAFE;
CREATE OR REPLACE FUNCTION uniform(a random_variable, b double precision)
RETURNS random_variable AS $$ SELECT provsql.rv_parametric2('uniform', ($1)::uuid, NULL, NULL, $2); $$
LANGUAGE sql STRICT VOLATILE PARALLEL SAFE;
CREATE OR REPLACE FUNCTION uniform(a double precision, b random_variable)
RETURNS random_variable AS $$ SELECT provsql.rv_parametric2('uniform', NULL, $1, ($2)::uuid, NULL); $$
LANGUAGE sql STRICT VOLATILE PARALLEL SAFE;
CREATE OR REPLACE FUNCTION uniform(a random_variable, b random_variable)
RETURNS random_variable AS $$ SELECT provsql.rv_parametric2('uniform', ($1)::uuid, NULL, ($2)::uuid, NULL); $$
LANGUAGE sql STRICT VOLATILE PARALLEL SAFE;
CREATE OR REPLACE FUNCTION exponential(lambda random_variable)
RETURNS random_variable AS $$ SELECT provsql.rv_parametric1('exponential', ($1)::uuid); $$
LANGUAGE sql STRICT VOLATILE PARALLEL SAFE;
CREATE OR REPLACE FUNCTION gamma(k random_variable, lambda double precision)
RETURNS random_variable AS $$ SELECT provsql.rv_parametric2('gamma', ($1)::uuid, NULL, NULL, $2); $$
LANGUAGE sql STRICT VOLATILE PARALLEL SAFE;
CREATE OR REPLACE FUNCTION gamma(k double precision, lambda random_variable)
RETURNS random_variable AS $$ SELECT provsql.rv_parametric2('gamma', NULL, $1, ($2)::uuid, NULL); $$
LANGUAGE sql STRICT VOLATILE PARALLEL SAFE;
CREATE OR REPLACE FUNCTION gamma(k random_variable, lambda random_variable)
RETURNS random_variable AS $$ SELECT provsql.rv_parametric2('gamma', ($1)::uuid, NULL, ($2)::uuid, NULL); $$
LANGUAGE sql STRICT VOLATILE PARALLEL SAFE;
CREATE OR REPLACE FUNCTION lognormal(mu random_variable, sigma double precision)
RETURNS random_variable AS $$ SELECT provsql.rv_parametric2('lognormal', ($1)::uuid, NULL, NULL, $2); $$
LANGUAGE sql STRICT VOLATILE PARALLEL SAFE;
CREATE OR REPLACE FUNCTION lognormal(mu double precision, sigma random_variable)
RETURNS random_variable AS $$ SELECT provsql.rv_parametric2('lognormal', NULL, $1, ($2)::uuid, NULL); $$
LANGUAGE sql STRICT VOLATILE PARALLEL SAFE;
CREATE OR REPLACE FUNCTION lognormal(mu random_variable, sigma random_variable)
RETURNS random_variable AS $$ SELECT provsql.rv_parametric2('lognormal', ($1)::uuid, NULL, ($2)::uuid, NULL); $$
LANGUAGE sql STRICT VOLATILE PARALLEL SAFE;
CREATE OR REPLACE FUNCTION weibull(k random_variable, lambda double precision)
RETURNS random_variable AS $$ SELECT provsql.rv_parametric2('weibull', ($1)::uuid, NULL, NULL, $2); $$
LANGUAGE sql STRICT VOLATILE PARALLEL SAFE;
CREATE OR REPLACE FUNCTION weibull(k double precision, lambda random_variable)
RETURNS random_variable AS $$ SELECT provsql.rv_parametric2('weibull', NULL, $1, ($2)::uuid, NULL); $$
LANGUAGE sql STRICT VOLATILE PARALLEL SAFE;
CREATE OR REPLACE FUNCTION weibull(k random_variable, lambda random_variable)
RETURNS random_variable AS $$ SELECT provsql.rv_parametric2('weibull', ($1)::uuid, NULL, ($2)::uuid, NULL); $$
LANGUAGE sql STRICT VOLATILE PARALLEL SAFE;
CREATE OR REPLACE FUNCTION pareto(xm random_variable, alpha double precision)
RETURNS random_variable AS $$ SELECT provsql.rv_parametric2('pareto', ($1)::uuid, NULL, NULL, $2); $$
LANGUAGE sql STRICT VOLATILE PARALLEL SAFE;
CREATE OR REPLACE FUNCTION pareto(xm double precision, alpha random_variable)
RETURNS random_variable AS $$ SELECT provsql.rv_parametric2('pareto', NULL, $1, ($2)::uuid, NULL); $$
LANGUAGE sql STRICT VOLATILE PARALLEL SAFE;
CREATE OR REPLACE FUNCTION pareto(xm random_variable, alpha random_variable)
RETURNS random_variable AS $$ SELECT provsql.rv_parametric2('pareto', ($1)::uuid, NULL, ($2)::uuid, NULL); $$
LANGUAGE sql STRICT VOLATILE PARALLEL SAFE;
CREATE OR REPLACE FUNCTION beta(alpha random_variable, beta double precision)
RETURNS random_variable AS $$ SELECT provsql.rv_parametric2('beta', ($1)::uuid, NULL, NULL, $2); $$
LANGUAGE sql STRICT VOLATILE PARALLEL SAFE;
CREATE OR REPLACE FUNCTION beta(alpha double precision, beta random_variable)
RETURNS random_variable AS $$ SELECT provsql.rv_parametric2('beta', NULL, $1, ($2)::uuid, NULL); $$
LANGUAGE sql STRICT VOLATILE PARALLEL SAFE;
CREATE OR REPLACE FUNCTION beta(alpha random_variable, beta random_variable)
RETURNS random_variable AS $$ SELECT provsql.rv_parametric2('beta', ($1)::uuid, NULL, ($2)::uuid, NULL); $$
LANGUAGE sql STRICT VOLATILE PARALLEL SAFE;
CREATE OR REPLACE FUNCTION inverse_gamma(alpha random_variable, beta double precision)
RETURNS random_variable AS $$ SELECT provsql.rv_parametric2('inverse_gamma', ($1)::uuid, NULL, NULL, $2); $$
LANGUAGE sql STRICT VOLATILE PARALLEL SAFE;
CREATE OR REPLACE FUNCTION inverse_gamma(alpha double precision, beta random_variable)
RETURNS random_variable AS $$ SELECT provsql.rv_parametric2('inverse_gamma', NULL, $1, ($2)::uuid, NULL); $$
LANGUAGE sql STRICT VOLATILE PARALLEL SAFE;
CREATE OR REPLACE FUNCTION inverse_gamma(alpha random_variable, beta random_variable)
RETURNS random_variable AS $$ SELECT provsql.rv_parametric2('inverse_gamma', ($1)::uuid, NULL, ($2)::uuid, NULL); $$
LANGUAGE sql STRICT VOLATILE PARALLEL SAFE;
CREATE OR REPLACE FUNCTION inverse_gaussian(mu random_variable, lambda double precision)
RETURNS random_variable AS $$ SELECT provsql.rv_parametric2('inverse_gaussian', ($1)::uuid, NULL, NULL, $2); $$
LANGUAGE sql STRICT VOLATILE PARALLEL SAFE;
CREATE OR REPLACE FUNCTION inverse_gaussian(mu double precision, lambda random_variable)
RETURNS random_variable AS $$ SELECT provsql.rv_parametric2('inverse_gaussian', NULL, $1, ($2)::uuid, NULL); $$
LANGUAGE sql STRICT VOLATILE PARALLEL SAFE;
CREATE OR REPLACE FUNCTION inverse_gaussian(mu random_variable, lambda random_variable)
RETURNS random_variable AS $$ SELECT provsql.rv_parametric2('inverse_gaussian', ($1)::uuid, NULL, ($2)::uuid, NULL); $$
LANGUAGE sql STRICT VOLATILE PARALLEL SAFE;
CREATE OR REPLACE FUNCTION logistic(mu double precision, s double precision)
RETURNS random_variable AS
$$
DECLARE
token uuid;
BEGIN
IF NOT provsql.is_finite_float8(mu) OR NOT provsql.is_finite_float8(s) THEN
RAISE EXCEPTION 'provsql.logistic: parameters must be finite (got mu=%, s=%)', mu, s;
END IF;
IF s < 0 THEN
RAISE EXCEPTION 'provsql.logistic: scale s must be non-negative (got %)', s;
END IF;
IF s = 0 THEN
RETURN provsql.as_random(mu);
END IF;
token := public.uuid_generate_v4();
PERFORM provsql.create_gate(token, 'rv');
PERFORM provsql.set_extra(token, 'logistic:' || mu || ',' || s);
RETURN provsql.random_variable_make(token);
END
$$ LANGUAGE plpgsql STRICT VOLATILE PARALLEL SAFE;
CREATE OR REPLACE FUNCTION poisson(lambda random_variable)
RETURNS random_variable AS
$$ SELECT provsql.rv_parametric1('poisson', ($1)::uuid); $$
LANGUAGE sql STRICT VOLATILE PARALLEL SAFE;
CREATE OR REPLACE FUNCTION binomial(n integer, p random_variable)
RETURNS random_variable AS
$$ SELECT provsql.rv_parametric2('binomial', NULL, $1::double precision,
($2)::uuid, NULL); $$
LANGUAGE sql STRICT VOLATILE PARALLEL SAFE;
CREATE OR REPLACE FUNCTION geometric(p random_variable)
RETURNS random_variable AS
$$ SELECT provsql.rv_parametric1('geometric', ($1)::uuid); $$
LANGUAGE sql STRICT VOLATILE PARALLEL SAFE;
CREATE OR REPLACE FUNCTION negative_binomial(r double precision, p random_variable)
RETURNS random_variable AS
$$ SELECT provsql.rv_parametric2('negative_binomial', NULL, $1, ($2)::uuid, NULL); $$
LANGUAGE sql STRICT VOLATILE PARALLEL SAFE;
CREATE OR REPLACE FUNCTION negative_binomial(r random_variable, p double precision)
RETURNS random_variable AS
$$ SELECT provsql.rv_parametric2('negative_binomial', ($1)::uuid, NULL, NULL, $2); $$
LANGUAGE sql STRICT VOLATILE PARALLEL SAFE;
CREATE OR REPLACE FUNCTION negative_binomial(r random_variable, p random_variable)
RETURNS random_variable AS
$$ SELECT provsql.rv_parametric2('negative_binomial', ($1)::uuid, NULL, ($2)::uuid, NULL); $$
LANGUAGE sql STRICT VOLATILE PARALLEL SAFE;
CREATE OR REPLACE FUNCTION random_variable_cond(rv random_variable, cond uuid)
RETURNS random_variable AS
$$
DECLARE
tgt uuid;
ev uuid;
result uuid;
ch uuid[];
BEGIN
IF cond IS NULL OR cond = gate_one() THEN
RETURN rv;
END IF;
-- A point-equality "Y = c" on a bare random-variable leaf is an
-- OBSERVATION, not a truncation: rewrite it to the internal likelihood-
-- weighting evidence (its density / mass at c). This is what lets
-- "X | (normal(mu,1) = 8)" (a continuous point event, measure-zero as a
-- Boolean selection) condition as the disintegration rather than fold to
-- an infeasible event.
cond := provsql.evidence_as_observation(cond);
tgt := (rv)::uuid;
IF get_gate_type(tgt) = 'conditioned'
AND array_length(get_children(tgt), 1) = 2 THEN
-- Fold (X|A)|B = X|(A∧B): the rv-carrier conditioned gate is the
-- two-child [target, condition] shape; accumulate the new event.
ch := get_children(tgt);
tgt := ch[1];
ev := provenance_times(ch[2], cond);
ELSE
ev := cond;
END IF;
result := public.uuid_generate_v5(uuid_ns_provsql(),
concat('conditioned', tgt, ev));
PERFORM create_gate(result, 'conditioned', ARRAY[tgt, ev]);
RETURN (result)::random_variable;
END
$$ LANGUAGE plpgsql SET search_path=provsql,pg_temp,public
SECURITY DEFINER PARALLEL SAFE;
CREATE OR REPLACE FUNCTION given(evidence UUID) RETURNS UUID AS
$$
BEGIN
RETURN provsql.evidence_as_observation(evidence);
END
$$ LANGUAGE plpgsql VOLATILE PARALLEL SAFE
SET search_path=provsql,pg_temp,public SECURITY DEFINER;
DO $$ BEGIN
IF NOT EXISTS (
SELECT 1 FROM pg_operator o
JOIN pg_namespace n ON n.oid = o.oprnamespace
WHERE n.nspname = 'provsql' AND o.oprname = '|'
AND o.oprleft = 0
AND o.oprright = 'pg_catalog.uuid'::regtype::oid
AND o.oprcode <> 0
) THEN
CREATE OPERATOR | (RIGHTARG=UUID, PROCEDURE=given);
END IF;
END $$;
CREATE OR REPLACE FUNCTION given(predicate boolean) RETURNS UUID AS
$$
BEGIN
RAISE EXCEPTION 'given(predicate) / prefix | (predicate) must be rewritten '
'by the ProvSQL planner hook: the operand must be a Boolean combination '
'of random_variable / aggregate comparisons (is provsql.active off?)';
END
$$ LANGUAGE plpgsql IMMUTABLE STRICT PARALLEL SAFE;
DO $$ BEGIN
IF NOT EXISTS (
SELECT 1 FROM pg_operator o
JOIN pg_namespace n ON n.oid = o.oprnamespace
WHERE n.nspname = 'provsql' AND o.oprname = '|'
AND o.oprleft = 0
AND o.oprright = 'pg_catalog.bool'::regtype::oid
AND o.oprcode <> 0
) THEN
CREATE OPERATOR | (RIGHTARG=boolean, PROCEDURE=given);
END IF;
END $$;
CREATE OR REPLACE FUNCTION evidence_as_observation(ev uuid) RETURNS uuid AS
$$
DECLARE
ch uuid[];
i1 integer;
leaf uuid;
datum_gate uuid;
BEGIN
IF ev IS NULL OR provsql.get_gate_type(ev) <> 'cmp' THEN
RETURN ev;
END IF;
ch := provsql.get_children(ev);
IF array_length(ch, 1) <> 2 THEN
RETURN ev;
END IF;
-- The cmp stores the comparison OPERATOR's OID in info1; match on its name
-- '=' the same way the C-side cmpOpFromOid does (get_opname), rather than a
-- fixed operator OID (which varies per install / carrier type).
SELECT info1 INTO i1 FROM provsql.get_infos(ev);
IF (SELECT oprname FROM pg_catalog.pg_operator WHERE oid = i1) IS DISTINCT FROM '=' THEN
RETURN ev; -- not an equality
END IF;
IF provsql.get_gate_type(ch[1]) = 'rv'
AND provsql.get_gate_type(ch[2]) = 'value' THEN
leaf := ch[1]; datum_gate := ch[2];
ELSIF provsql.get_gate_type(ch[2]) = 'rv'
AND provsql.get_gate_type(ch[1]) = 'value' THEN
leaf := ch[2]; datum_gate := ch[1];
ELSE
RETURN ev; -- not a bare-leaf-vs-constant point event
END IF;
RETURN provsql.observe((leaf)::random_variable,
provsql.get_extra(datum_gate)::double precision);
END
$$ LANGUAGE plpgsql VOLATILE
SET search_path=provsql,pg_temp,public SECURITY DEFINER PARALLEL SAFE;
CREATE OR REPLACE FUNCTION observe(x random_variable, datum double precision)
RETURNS uuid AS
$$
DECLARE
leaf uuid := (x)::uuid;
result uuid;
BEGIN
IF provsql.get_gate_type(leaf) <> 'rv' THEN
RAISE EXCEPTION 'provsql.observe: the argument must be a bare '
'random-variable leaf (a gate_rv), got a % gate', provsql.get_gate_type(leaf)
USING HINT = 'observe binds a datum to a single distribution leaf; '
'observing a derived quantity (a sum, product, or comparison) needs '
'a change-of-variables density and is out of scope.';
END IF;
IF NOT provsql.is_finite_float8(datum) THEN
RAISE EXCEPTION 'provsql.observe: datum must be finite (got %)', datum;
END IF;
result := public.uuid_generate_v4();
PERFORM provsql.create_gate(result, 'observe', ARRAY[leaf]);
PERFORM provsql.set_extra(result, datum::text);
RETURN result;
END
$$ LANGUAGE plpgsql VOLATILE
SET search_path=provsql,pg_temp,public SECURITY DEFINER PARALLEL SAFE;
CREATE OR REPLACE FUNCTION and_agg_sfunc(state uuid, ev uuid)
RETURNS uuid AS
$$
SELECT provsql.provenance_times(state, ev);
$$ LANGUAGE sql PARALLEL SAFE;
DO $$ BEGIN
IF NOT EXISTS (
SELECT 1 FROM pg_aggregate a
JOIN pg_proc p ON p.oid = a.aggfnoid
JOIN pg_namespace n ON n.oid = p.pronamespace
WHERE n.nspname = 'provsql' AND p.proname = 'and_agg'
AND p.pronargs = 1
AND p.proargtypes[0] = 'pg_catalog.uuid'::regtype::oid
) THEN
CREATE AGGREGATE and_agg(uuid) (
SFUNC = and_agg_sfunc,
STYPE = uuid
);
END IF;
END $$;
CREATE OR REPLACE FUNCTION evidence(evidence uuid)
RETURNS double precision
AS 'provsql','rv_evidence' LANGUAGE C STRICT PARALLEL SAFE;
CREATE OR REPLACE FUNCTION observe_atoms(evidence uuid)
RETURNS uuid[] AS
$$
WITH RECURSIVE walk(tok) AS (
SELECT evidence
UNION
SELECT c
FROM walk, LATERAL unnest(provsql.get_children(walk.tok)) AS c
WHERE provsql.get_gate_type(walk.tok) = 'times'
)
SELECT array_agg(tok ORDER BY tok)
FROM walk
WHERE provsql.get_gate_type(tok) = 'observe';
$$ LANGUAGE sql STABLE PARALLEL SAFE SET search_path=provsql,pg_temp,public;
CREATE OR REPLACE FUNCTION shapley_observe(
target uuid, evidence uuid, payoff text DEFAULT 'expected')
RETURNS TABLE(observation uuid, value double precision) AS
$$
DECLARE
atoms uuid[];
n int;
nmasks int;
pv double precision[];
popc int[];
fact double precision[];
mask int;
i int;
j int;
cnt int;
subset uuid[];
ev_s uuid;
sh double precision;
bit int;
s_size int;
BEGIN
IF payoff NOT IN ('expected', 'variance') THEN
RAISE EXCEPTION 'provsql.shapley_observe: payoff must be ''expected'' or '
'''variance'' (got %)', payoff;
END IF;
atoms := provsql.observe_atoms(evidence);
n := coalesce(array_length(atoms, 1), 0);
IF n = 0 THEN
RAISE EXCEPTION 'provsql.shapley_observe: evidence contains no observe() '
'atoms (got a % gate)', provsql.get_gate_type(evidence);
END IF;
IF n > 12 THEN
RAISE EXCEPTION 'provsql.shapley_observe: exact attribution over % '
'observations is exponential; capped at 12 (sampling-based '
'attribution is future work)', n;
END IF;
-- factorials 0!..n! (fact[k+1] = k!)
fact := ARRAY[1::double precision];
FOR i IN 1..n LOOP fact := fact || (fact[i] * i); END LOOP;
nmasks := (1 << n);
pv := array_fill(NULL::double precision, ARRAY[nmasks]);
popc := array_fill(0, ARRAY[nmasks]);
-- Payoff value function for every subset of observations.
FOR mask IN 0 .. nmasks - 1 LOOP
subset := ARRAY[]::uuid[];
cnt := 0;
FOR i IN 0 .. n - 1 LOOP
IF (mask >> i) & 1 = 1 THEN
subset := subset || atoms[i + 1];
cnt := cnt + 1;
END IF;
END LOOP;
popc[mask + 1] := cnt;
IF cnt = 0 THEN
ev_s := provsql.gate_one(); -- prior (no evidence)
ELSE
ev_s := provsql.provenance_times(VARIADIC subset);
END IF;
IF payoff = 'expected' THEN
pv[mask + 1] := provsql.rv_moment(target, 1, false, ev_s);
ELSE
pv[mask + 1] := provsql.rv_moment(target, 2, true, ev_s);
END IF;
END LOOP;
-- Shapley value of each observation atom.
FOR i IN 0 .. n - 1 LOOP
sh := 0;
bit := (1 << i);
FOR mask IN 0 .. nmasks - 1 LOOP
IF (mask >> i) & 1 = 0 THEN -- subsets S not containing i
s_size := popc[mask + 1];
-- weight |S|! (n-|S|-1)! / n!
sh := sh + (fact[s_size + 1] * fact[n - s_size] / fact[n + 1])
* (pv[(mask | bit) + 1] - pv[mask + 1]);
END IF;
END LOOP;
observation := atoms[i + 1];
value := sh;
RETURN NEXT;
END LOOP;
END
$$ LANGUAGE plpgsql VOLATILE
SET search_path=provsql,pg_temp,public SECURITY DEFINER;
CREATE OR REPLACE FUNCTION agg_collapsed_moment(token uuid, k integer)
RETURNS double precision
AS 'provsql','agg_collapsed_moment' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
CREATE OR REPLACE FUNCTION agg_collapsed_moments(token uuid)
RETURNS double precision[]
AS 'provsql','agg_collapsed_moments' LANGUAGE C IMMUTABLE STRICT PARALLEL SAFE;
CREATE OR REPLACE FUNCTION agg_raw_moment(
token agg_token,
k integer,
prov UUID = gate_one(),
method text = NULL,
arguments text = NULL)
RETURNS DOUBLE PRECISION AS $$
DECLARE
aggregation_function VARCHAR;
child_pairs uuid[];
pair_children uuid[];
n integer;
i integer;
j integer;
vals float8[];
toks uuid[];
total float8;
total_probability float8;
tup integer[];
d integer;
prod_v float8;
distinct_tok uuid[];
conj_token uuid;
prob float8;
sign_max float8;
BEGIN
IF token IS NULL OR k IS NULL THEN
RETURN NULL;
END IF;
IF k < 0 THEN
RAISE EXCEPTION 'agg_raw_moment(): k must be non-negative (got %)', k;
END IF;
-- Aggregate-carrier CASE (a gate_case over aggregate branches): a first-match
-- guarded selection. The moment is CONDITIONAL on the CASE's value being
-- defined (NULL only when it never is, mirroring the MIN/MAX convention):
-- E[pick^k | defined ∧ prov]
-- = Σ_i P(region_i ∧ def_i) · E[value_i^k | region_i ∧ def_i]
-- / Σ_i P(region_i ∧ def_i),
-- where region_i = (¬g_1 ∧ … ∧ ¬g_{i-1}) ∧ g_i ∧ prov is the world set that
-- selects branch i (the default's region is "all guards false") and def_i is
-- the branch's defined event (agg_defined_event: gate_one for sum / count /
-- constants, "some row present" for min / max / avg, recursive for a nested
-- CASE). Both factors are exact: probability() over the region ∧ def event,
-- and the conditional aggregate moment (a recursive agg_raw_moment on the
-- branch aggregate, which conditions on its own definedness within the
-- region, so the two factors weigh the same worlds). The regions are
-- mutually exclusive, so the terms sum with no inclusion-exclusion, and
-- correlation between a guard and its branch (shared input tuples) is
-- carried by the conditioning, exactly as HAVING carries it. When every
-- branch is defined everywhere, the defined mass equals P(prov) and the
-- formula reduces to the plain region-weighted sum.
IF get_gate_type(token) = 'case' THEN
IF k = 0 THEN
RETURN 1;
END IF;
DECLARE
wires uuid[] := get_children(token);
nw integer := array_length(get_children(token), 1);
m integer := (array_length(get_children(token), 1) - 1) / 2;
running_neg uuid := gate_one();
region_full uuid;
prov_p float8;
p float8;
total float8 := 0;
def_mass float8 := 0;
ci integer;
vuid uuid;
bm float8;
BEGIN
prov_p := probability(prov);
IF prov_p IS NULL OR prov_p <= 0 THEN
RETURN NULL; -- impossible conditioning event
END IF;
-- Branches 1..m are the guarded WHENs; branch m+1 is the ELSE default,
-- whose region is "all guards false".
FOR ci IN 1 .. m + 1 LOOP
IF ci <= m THEN
region_full := provenance_times(running_neg, wires[2 * ci - 1], prov);
vuid := wires[2 * ci];
running_neg :=
provenance_times(running_neg, provenance_not(wires[2 * ci - 1]));
ELSE
region_full := provenance_times(running_neg, prov);
vuid := wires[nw];
END IF;
p := probability(provenance_times(region_full,
agg_defined_event(vuid)));
IF p > 0 THEN
-- E[value_i^k | region_i ∧ def_i]: a constant branch is a Dirac
-- (c^k, exact); a single aggregate or nested CASE is exact via
-- agg_raw_moment (whose MIN/MAX/CASE arms condition on their own
-- definedness within the region); an arithmetic / composite branch
-- takes the Monte-Carlo scalar path (which composes with the
-- aggregate leaves).
IF get_gate_type(vuid) = 'value' THEN
bm := power(CAST(get_extra(vuid) AS float8), k);
ELSIF get_gate_type(vuid) IN ('agg', 'case') THEN
bm := agg_raw_moment(agg_token_make(vuid, 0), k, region_full,
method, arguments);
ELSE
bm := rv_moment(vuid, k, false, region_full);
END IF;
total := total + p * bm;
def_mass := def_mass + p;
END IF;
END LOOP;
IF def_mass <= epsilon() THEN
RETURN NULL; -- the CASE's value is never defined under prov
END IF;
RETURN total / def_mass;
END;
END IF;
IF get_gate_type(token) <> 'agg' THEN
IF get_gate_type(token) IN ('arith', 'conditioned') THEN
RAISE EXCEPTION 'expected / variance / moment over an arithmetic '
'combination of aggregates (e.g. SUM(x) + SUM(y) or SUM(x) + 5), or a '
'conditioning of one, is not yet supported: a moment can be taken only '
'over a single aggregate (SUM / COUNT / MIN / MAX), optionally '
'conditioned (SUM(x) | C)'
USING HINT = 'Take the moment of each aggregate separately, or condition '
'the bare aggregate.';
ELSE
RAISE EXCEPTION USING MESSAGE='Wrong gate type for agg_raw_moment computation';
END IF;
END IF;
IF k = 0 THEN
RETURN 1;
END IF;
SELECT pp.proname::varchar FROM pg_proc pp
WHERE oid=(get_infos(token)).info1
INTO aggregation_function;
child_pairs := get_children(token);
n := COALESCE(array_length(child_pairs, 1), 0);
IF aggregation_function = 'sum' OR aggregation_function = 'count' THEN
-- count(col) keeps the COUNT identity at the gate level but its value is a
-- SUM of per-row 0/1 indicators, so its moments are computed exactly like
-- SUM (and its empty group is the real value 0, like SUM). count(*)
-- arrives here as 'sum' (it normalises to F_SUM_INT4); count(col) as 'count'.
-- Trivial empty aggregation: SUM = 0, so SUM^k = 0 for k >= 1.
-- Note: agg_token semantics treat the "no row included" world as
-- SUM = 0, so this stays consistent with k = 1 (= expected()).
IF n = 0 THEN
RETURN 0;
END IF;
-- Collapsed fast path: a correlated COUNT / SUM whose per-row selection
-- events share a single continuous latent has an O(G·n) 1-D quadrature,
-- vastly cheaper than the O(n^k) tuple enumeration below (which is the
-- O(n^2) pair-probability bottleneck for the variance). Only fires
-- unconditionally (prov = one) and for k in {1, 2}; agg_collapsed_moment
-- returns NULL when the shared-latent pattern does not match, and we
-- fall through to the exact enumeration.
IF prov = gate_one() AND k <= 2 THEN
total := agg_collapsed_moment((token)::uuid, k);
IF total IS NOT NULL THEN
RETURN total;
END IF;
END IF;
-- Extract per-child token + value arrays.
vals := ARRAY[]::float8[];
toks := ARRAY[]::uuid[];
FOR i IN 1..n LOOP
pair_children := get_children(child_pairs[i]);
toks := toks || pair_children[1];
vals := vals || CAST(get_extra(pair_children[2]) AS float8);
END LOOP;
-- Enumerate all k-tuples (i_1, ..., i_k) in {1..n}^k. tup is the
-- current tuple; we step through them in lexicographic order.
total := 0;
tup := array_fill(1, ARRAY[k]);
LOOP
prod_v := 1;
FOR j IN 1..k LOOP
prod_v := prod_v * vals[tup[j]];
END LOOP;
SELECT array_agg(DISTINCT toks[idx]) INTO distinct_tok
FROM unnest(tup) AS idx;
IF prov <> gate_one() THEN
distinct_tok := distinct_tok || prov;
END IF;
conj_token := provenance_times(VARIADIC distinct_tok);
prob := probability_evaluate(conj_token, method, arguments);
total := total + prod_v * prob;
d := k;
WHILE d >= 1 AND tup[d] = n LOOP
tup[d] := 1;
d := d - 1;
END LOOP;
EXIT WHEN d = 0;
tup[d] := tup[d] + 1;
END LOOP;
ELSIF aggregation_function = 'min' OR aggregation_function = 'max' THEN
-- Rank enumeration: per distinct value v, P(MIN = v) is the
-- probability that some t_i with v_i=v is true and all t_j with
-- smaller v are false. For MAX we negate values so the same
-- "smaller-than" rank logic computes MIN-of-negated, then flip.
-- The outer multiplier picks up the right sign for the k-th moment
-- of MAX: E[MAX^k] = (-1)^k * E[MIN(-v)^k], so sign_max = (-1)^k.
sign_max := CASE
WHEN aggregation_function = 'max'
THEN power(-1::float8, k)
ELSE 1
END;
-- MIN/MAX over the empty input world are NULL (no elements), not ±Infinity:
-- SQL returns one row with a NULL value. The moment is therefore CONDITIONAL
-- on the aggregate being defined (non-empty) -- the empty world is excluded
-- and the result renormalised by P(prov AND non-empty). (count, whose empty
-- value 0 is a real value, keeps the empty world; sum keeps it too, as 0.)
IF n = 0 THEN
RETURN NULL; -- structurally empty: MIN/MAX undefined
END IF;
-- Numerator E[MIN^k . 1{prov AND non-empty}] (the rank sum naturally omits
-- the empty world, since every term requires a present token).
WITH tok_value AS (
SELECT (get_children(c))[1] AS tok,
(CASE WHEN aggregation_function='max' THEN -1 ELSE 1 END)
* CAST(get_extra((get_children(c))[2]) AS DOUBLE PRECISION) AS v
FROM UNNEST(child_pairs) AS c
) SELECT sign_max * COALESCE(SUM(p * power(v, k)), 0) FROM (
SELECT t1.v AS v,
probability_evaluate(
CASE WHEN prov = gate_one()
THEN provenance_monus(provenance_plus(ARRAY_AGG(t1.tok)),
provenance_plus(ARRAY_AGG(t2.tok)))
ELSE provenance_times(prov,
provenance_monus(provenance_plus(ARRAY_AGG(t1.tok)),
provenance_plus(ARRAY_AGG(t2.tok)))) END,
method, arguments) AS p
FROM tok_value t1 LEFT OUTER JOIN tok_value t2 ON t1.v > t2.v
GROUP BY t1.v) tmp
INTO total;
-- Denominator P(prov AND non-empty) = P(prov (x) (+) tokens).
SELECT probability_evaluate(
CASE WHEN prov = gate_one()
THEN provenance_plus(ARRAY_AGG(tok))
ELSE provenance_times(prov, provenance_plus(ARRAY_AGG(tok))) END,
method, arguments)
FROM (SELECT (get_children(c))[1] AS tok FROM UNNEST(child_pairs) AS c) s
INTO total_probability;
IF total_probability <= epsilon() THEN
RETURN NULL; -- never defined under prov: MIN/MAX undefined
END IF;
RETURN total / total_probability; -- already conditional; skip generic norm
ELSIF aggregation_function = 'avg' THEN
-- AVG = SUM/COUNT is a ratio of two correlated world-dependent
-- quantities, so the k-tuple expansion above does not apply. Like
-- MIN/MAX, AVG over the empty world is NULL, so its moment conditions
-- on the aggregate being defined (COUNT >= 1), NULL when it never is.
-- Two routes:
-- * EXACT (independent rows, unconditional): the joint (sum, count)
-- PMF folded in C by agg_avg_moment_exact --
-- E[AVG^k | COUNT>=1] = Σ_{(s,c), c>=1} (s/c)^k pmf(s,c) / P(c>=1).
-- * Monte-Carlo scalar fallback otherwise (an outer conditioning
-- event, shared leaves, compound contributors): rv_moment samples
-- the agg gate per world; its NaN-skip on empty draws implements
-- the same conditional-on-defined convention, at the
-- provsql.rv_mc_samples budget (0 raises, per convention).
IF n = 0 THEN
RETURN NULL; -- structurally empty: AVG undefined
END IF;
IF prov = gate_one() THEN
total := agg_avg_moment_exact((token)::uuid, k);
IF total IS NOT NULL THEN
RETURN total;
END IF;
END IF;
RETURN rv_moment((token)::uuid, k, false, prov);
ELSE
RAISE EXCEPTION USING MESSAGE=
'Cannot compute moment for aggregation function ' || aggregation_function;
END IF;
-- Conditional normalisation: E[X^k · 1_A] / P(A) = E[X^k | A].
IF prov <> gate_one()
AND total <> 0
AND total <> 'Infinity'::float8
AND total <> '-Infinity'::float8 THEN
total := total / probability_evaluate(prov, method, arguments);
END IF;
RETURN total;
END
$$ LANGUAGE plpgsql PARALLEL SAFE SET search_path=provsql SECURITY DEFINER;
CREATE OR REPLACE FUNCTION variance(
input ANYELEMENT,
prov UUID = gate_one(),
method text = NULL,
arguments text = NULL)
RETURNS DOUBLE PRECISION AS $$
DECLARE
m1 float8;
m2 float8;
BEGIN
IF pg_typeof(input) = 'random_variable'::regtype THEN
IF input IS NULL THEN
RETURN NULL;
END IF;
-- Conditioning on prov is handled inside rv_moment: when prov
-- resolves to gate_one() (the default, or load-time
-- simplification of any always-true sub-circuit) the
-- unconditional analytical path runs unchanged; otherwise the
-- joint-circuit loader unifies shared gate_rv leaves between
-- input and prov, and the conditional path runs either
-- truncated-distribution closed form or MC rejection.
RETURN provsql.rv_moment(
rv_conditioned_target((input::random_variable)::uuid), 2, true,
rv_conditioned_prov((input::random_variable)::uuid, prov));
END IF;
IF pg_typeof(input) = 'agg_token'::regtype THEN
IF input IS NULL THEN
RETURN NULL;
END IF;
-- Collapsed fast path: E[C] and E[C^2] from a single circuit load and plan
-- build, instead of two agg_raw_moment() calls that each reload. Mirrors
-- the guard in agg_raw_moment (unconditional only, prov = one); on any
-- mismatch agg_collapsed_moments returns NULL and we fall through to the
-- generic per-order path (which handles conditioning, SUM enumeration, ...).
IF rv_conditioned_prov(input::uuid, prov) = gate_one() THEN
DECLARE ms float8[];
BEGIN
ms := agg_collapsed_moments(
(agg_conditioned_target(input::agg_token))::uuid);
IF ms IS NOT NULL THEN
RETURN ms[2] - ms[1] * ms[1];
END IF;
END;
END IF;
m1 := agg_raw_moment(agg_conditioned_target(input::agg_token), 1,
rv_conditioned_prov(input::uuid, prov), method, arguments);
m2 := agg_raw_moment(agg_conditioned_target(input::agg_token), 2,
rv_conditioned_prov(input::uuid, prov), method, arguments);
IF m1 IS NULL OR m2 IS NULL THEN
RETURN NULL;
END IF;
RETURN m2 - m1 * m1;
END IF;
-- Bernoulli event token (see moment()): Var[X] = p(1 - p).
IF pg_typeof(input) = 'uuid'::regtype THEN
IF input IS NULL THEN
RETURN NULL;
END IF;
m1 := provsql.probability_evaluate(provsql.cond(input::uuid, prov),
method, arguments);
RETURN m1 * (1 - m1);
END IF;
RAISE EXCEPTION 'variance() is not yet supported for input type %', pg_typeof(input);
END
$$ LANGUAGE plpgsql PARALLEL SAFE SET search_path=provsql SECURITY DEFINER;
CREATE OR REPLACE FUNCTION moment(
input ANYELEMENT,
k integer,
prov UUID = gate_one(),
method text = NULL,
arguments text = NULL)
RETURNS DOUBLE PRECISION AS $$
BEGIN
IF pg_typeof(input) = 'random_variable'::regtype THEN
IF input IS NULL OR k IS NULL THEN
RETURN NULL;
END IF;
-- See variance() above: rv_moment handles the conditional/unconditional
-- dispatch internally based on the resolved prov gate type.
RETURN provsql.rv_moment(
rv_conditioned_target((input::random_variable)::uuid), k, false,
rv_conditioned_prov((input::random_variable)::uuid, prov));
END IF;
IF pg_typeof(input) = 'agg_token'::regtype THEN
RETURN agg_raw_moment(agg_conditioned_target(input::agg_token), k,
rv_conditioned_prov(input::uuid, prov), method, arguments);
END IF;
-- A bare provenance event token (a gate_cmp lifted from an RV comparison,
-- e.g. expected(x <= c)) is a Bernoulli indicator: X in {0,1}, so every raw
-- moment E[X^k] with k >= 1 equals P(event), and E[X^0] = 1. cond() applies
-- the optional conditioning prov (a no-op for the default gate_one()).
IF pg_typeof(input) = 'uuid'::regtype THEN
IF input IS NULL OR k IS NULL THEN
RETURN NULL;
END IF;
IF k = 0 THEN
RETURN 1;
END IF;
RETURN provsql.probability_evaluate(provsql.cond(input::uuid, prov),
method, arguments);
END IF;
RAISE EXCEPTION 'moment() is not yet supported for input type %', pg_typeof(input);
END
$$ LANGUAGE plpgsql PARALLEL SAFE SET search_path=provsql SECURITY DEFINER;
CREATE OR REPLACE FUNCTION central_moment(
input ANYELEMENT,
k integer,
prov UUID = gate_one(),
method text = NULL,
arguments text = NULL)
RETURNS DOUBLE PRECISION AS $$
DECLARE
mu float8;
total float8;
i integer;
raw_i float8;
binom float8;
-- iterative binomial coefficient C(k, i)
k_double float8;
BEGIN
IF pg_typeof(input) = 'random_variable'::regtype THEN
IF input IS NULL OR k IS NULL THEN
RETURN NULL;
END IF;
-- See variance() above: rv_moment handles the conditional/unconditional
-- dispatch internally based on the resolved prov gate type.
RETURN provsql.rv_moment(
rv_conditioned_target((input::random_variable)::uuid), k, true,
rv_conditioned_prov((input::random_variable)::uuid, prov));
END IF;
IF pg_typeof(input) = 'agg_token'::regtype THEN
IF input IS NULL OR k IS NULL THEN
RETURN NULL;
END IF;
IF k < 0 THEN
RAISE EXCEPTION 'central_moment(): k must be non-negative (got %)', k;
END IF;
IF k = 0 THEN RETURN 1; END IF;
IF k = 1 THEN RETURN 0; END IF;
mu := agg_raw_moment(agg_conditioned_target(input::agg_token), 1,
rv_conditioned_prov(input::uuid, prov), method, arguments);
IF mu IS NULL THEN RETURN NULL; END IF;
-- mu may be ±Infinity for empty MIN / MAX with positive empty
-- probability; central_moment is undefined in that case.
IF mu = 'Infinity'::float8 OR mu = '-Infinity'::float8 THEN
RETURN mu;
END IF;
total := 0;
binom := 1; -- C(k, 0)
k_double := k;
FOR i IN 0..k LOOP
raw_i := agg_raw_moment(agg_conditioned_target(input::agg_token), i,
rv_conditioned_prov(input::uuid, prov), method, arguments);
IF raw_i IS NULL THEN RETURN NULL; END IF;
total := total + binom * power(-mu, k - i) * raw_i;
-- C(k, i+1) = C(k, i) * (k - i) / (i + 1)
IF i < k THEN
binom := binom * (k_double - i) / (i + 1);
END IF;
END LOOP;
RETURN total;
END IF;
-- Bernoulli event token (see moment()): with p = P(event),
-- E[(X-p)^k] = (1-p)(-p)^k + p(1-p)^k; k = 0 -> 1, k = 1 -> 0.
IF pg_typeof(input) = 'uuid'::regtype THEN
IF input IS NULL OR k IS NULL THEN
RETURN NULL;
END IF;
IF k < 0 THEN
RAISE EXCEPTION 'central_moment(): k must be non-negative (got %)', k;
END IF;
IF k = 0 THEN RETURN 1; END IF;
IF k = 1 THEN RETURN 0; END IF;
mu := provsql.probability_evaluate(provsql.cond(input::uuid, prov),
method, arguments);
RETURN (1 - mu) * power(-mu, k) + mu * power(1 - mu, k);
END IF;
RAISE EXCEPTION 'central_moment() is not yet supported for input type %', pg_typeof(input);
END
$$ LANGUAGE plpgsql PARALLEL SAFE SET search_path=provsql SECURITY DEFINER;
-- Annotation transparency for identity-sensitive consumers: the
-- conjunctive-leaves walker and the reachability edge gatherer see
-- through the transparent gate_annotation wrapper.
CREATE OR REPLACE FUNCTION token_conjunctive_leaves(token uuid)
RETURNS uuid[] AS
$$
WITH RECURSIVE walk(g) AS (
SELECT token
UNION
SELECT c FROM walk w, unnest(provsql.get_children(w.g)) AS c
WHERE provsql.get_gate_type(w.g) IN ('times', 'project', 'eq', 'annotation')
)
SELECT CASE WHEN bool_and(provsql.get_gate_type(g)
IN ('times', 'project', 'eq', 'annotation', 'input'))
THEN array_agg(DISTINCT g)
FILTER (WHERE provsql.get_gate_type(g) = 'input')
ELSE NULL END
FROM walk;
$$ LANGUAGE sql STABLE;
CREATE OR REPLACE FUNCTION gather_reachability_edges(
IN rel regclass,
IN source_attribute TEXT,
IN destination_attribute TEXT,
IN extra_vertices TEXT[],
IN edge_quals TEXT DEFAULT NULL,
IN rel_sql TEXT DEFAULT NULL,
OUT sources INT[],
OUT destinations INT[],
OUT tokens UUID[],
OUT probabilities DOUBLE PRECISION[],
OUT block_keys UUID[],
OUT block_indices INT[],
OUT extra_ids INT[],
OUT vertices TEXT[])
AS
$$
DECLARE
tkind text;
bkey_expr text;
sel_probs text;
sel_bkeys text;
sel_bidx text;
verbosity int := coalesce(current_setting('provsql.verbose_level', true)::int, 0);
BEGIN
-- Consult the per-table characterisation registry (TID / BID / OPAQUE,
-- maintained by add_provenance / repair_key and the CTAS lineage hook):
-- a TID relation is certified all-independent-inputs, a BID relation
-- holds input or mulinput rows with the block structure given by the
-- registry's key columns. Derived (OPAQUE), unregistered, or
-- subquery-defined edges take the fully dynamic per-token path.
IF rel IS NOT NULL AND rel_sql IS NULL THEN
tkind := (provsql.get_table_info(rel::oid)).kind;
END IF;
IF tkind NOT IN ('tid', 'bid') THEN
tkind := NULL;
END IF;
IF tkind = 'bid' THEN
SELECT string_agg(quote_ident(a.attname) || '::text', ' || '','' || '
ORDER BY k.ord)
INTO bkey_expr
FROM unnest((provsql.get_table_info(rel::oid)).block_key)
WITH ORDINALITY AS k(attnum, ord)
JOIN pg_attribute a ON a.attrelid = rel AND a.attnum = k.attnum;
-- An empty registry key means the whole table is one block.
bkey_expr := coalesce(bkey_expr, quote_literal(''));
END IF;
IF tkind IS NOT NULL AND verbosity >= 20 THEN
-- The function-level client_min_messages = warning (which silences
-- the CTAS / DROP TABLE chatter) would also swallow this notice;
-- lift it for the one RAISE. The function-level SET restores the
-- caller's value at exit regardless.
PERFORM set_config('client_min_messages', 'notice', true);
RAISE NOTICE 'ProvSQL: catalog characterises % as %', rel, upper(tkind);
PERFORM set_config('client_min_messages', 'warning', true);
END IF;
-- Materialize the edges with their tokens; the planner hook resolves
-- provenance() over the tracked relation, and remove_provenance strips
-- the automatic provsql column so the later aggregation is plain SQL.
-- For a BID relation the synthetic per-block key (a v5 UUID over the
-- registry key columns' values) is computed here, while the columns
-- are in scope.
DROP TABLE IF EXISTS provsql_reachability_edges_tmp;
EXECUTE format(
'CREATE TEMP TABLE provsql_reachability_edges_tmp AS '
|| 'SELECT %1$I::text AS u, %2$I::text AS v, '
|| 'provsql.strip_annotations(provsql.provenance()) AS token%5$s '
|| 'FROM %3$s WHERE %1$I IS NOT NULL AND %2$I IS NOT NULL%4$s',
source_attribute, destination_attribute,
CASE WHEN rel_sql IS NULL THEN rel::text
ELSE '(' || rel_sql || ') AS provsql_edge_subquery' END,
CASE WHEN edge_quals IS NULL THEN ''
ELSE ' AND (' || edge_quals || ')' END,
CASE WHEN tkind = 'bid'
THEN ', public.uuid_generate_v5(provsql.uuid_ns_provsql(), '
|| quote_literal('bidblock' || rel::text || ':')
|| ' || ' || bkey_expr || ') AS bkey'
ELSE ', NULL::uuid AS bkey' END);
PERFORM provsql.remove_provenance('provsql_reachability_edges_tmp');
DROP TABLE IF EXISTS provsql_reachability_support_tmp;
IF tkind IS NULL THEN
-- Dynamic path: validate the token shapes and, for conjunction-shaped
-- (join-defined) tokens, the pairwise disjointness of their supports.
IF EXISTS (SELECT 1 FROM provsql_reachability_edges_tmp
WHERE provsql.get_gate_type(token) NOT IN ('input', 'mulinput', 'times',
'project', 'eq')) THEN
DROP TABLE provsql_reachability_edges_tmp;
RAISE EXCEPTION 'reachability: the provenance of % must consist of base input, repair_key, or conjunctive join tokens', coalesce(rel::text, 'the edge query');
END IF;
CREATE TEMP TABLE provsql_reachability_support_tmp AS
SELECT t.token, l.leaf
FROM (SELECT DISTINCT token FROM provsql_reachability_edges_tmp
WHERE provsql.get_gate_type(token) IN ('times', 'project', 'eq')) t,
LATERAL unnest(provsql.token_conjunctive_leaves(t.token)) AS l(leaf);
IF EXISTS (SELECT 1
FROM (SELECT DISTINCT token FROM provsql_reachability_edges_tmp) t
WHERE provsql.get_gate_type(t.token) IN ('times', 'project', 'eq')
AND provsql.token_conjunctive_leaves(t.token) IS NULL) THEN
DROP TABLE provsql_reachability_support_tmp;
DROP TABLE provsql_reachability_edges_tmp;
RAISE EXCEPTION 'reachability: a join-defined edge token is not a pure conjunction of base tuples';
END IF;
IF EXISTS (SELECT 1 FROM (
SELECT leaf FROM provsql_reachability_support_tmp
UNION ALL
SELECT DISTINCT token FROM provsql_reachability_edges_tmp
WHERE provsql.get_gate_type(token) = 'input'
) all_leaves
GROUP BY leaf HAVING count(*) > 1) THEN
DROP TABLE provsql_reachability_support_tmp;
DROP TABLE provsql_reachability_edges_tmp;
RAISE EXCEPTION 'reachability: join-defined edges share base tuples (their supports overlap), so they are not independent';
END IF;
END IF;
-- Per-kind classification expressions for the final aggregation: a TID
-- relation needs no per-row gate introspection at all; a BID relation
-- one get_gate_type per row (the input/mulinput split), block keys from
-- the precomputed column-derived key and indices by numbering within
-- the block; the dynamic path reads the gates.
IF tkind = 'tid' THEN
sel_probs := 'coalesce(provsql.get_prob(e.token), 1.0)';
sel_bkeys := $sql$'00000000-0000-0000-0000-000000000000'::uuid$sql$;
sel_bidx := '0';
ELSIF tkind = 'bid' THEN
sel_probs := 'coalesce(provsql.get_prob(e.token), 1.0)';
sel_bkeys := $sql$CASE WHEN provsql.get_gate_type(e.token) = 'mulinput'
THEN e.bkey
ELSE '00000000-0000-0000-0000-000000000000'::uuid END$sql$;
sel_bidx := 'e.bidx';
ELSE
sel_probs := $sql$CASE WHEN provsql.get_gate_type(e.token) IN ('times','project','eq')
THEN (SELECT CASE WHEN bool_or(coalesce(provsql.get_prob(s.leaf),1.0) = 0)
THEN 0.0
ELSE exp(sum(ln(coalesce(provsql.get_prob(s.leaf),1.0)))) END
FROM provsql_reachability_support_tmp s
WHERE s.token = e.token)
ELSE coalesce(provsql.get_prob(e.token), 1.0) END$sql$;
sel_bkeys := $sql$CASE WHEN provsql.get_gate_type(e.token) = 'mulinput'
THEN (provsql.get_children(e.token))[1]
ELSE '00000000-0000-0000-0000-000000000000'::uuid END$sql$;
sel_bidx := $sql$CASE WHEN provsql.get_gate_type(e.token) = 'mulinput'
THEN (provsql.get_infos(e.token)).info1 ELSE 0 END$sql$;
END IF;
EXECUTE format(
$sql$
WITH verts AS (
SELECT u AS x FROM provsql_reachability_edges_tmp
UNION SELECT v FROM provsql_reachability_edges_tmp
UNION SELECT unnest($1)),
ids AS (
SELECT x, (row_number() OVER (ORDER BY x))::int AS id FROM verts)
SELECT array_agg(iu.id), array_agg(iv.id),
array_agg(e.token),
array_agg(%s),
array_agg(%s),
array_agg(%s),
(SELECT array_agg(i.id ORDER BY ev.ord)
FROM unnest($1) WITH ORDINALITY AS ev(x, ord)
JOIN ids i ON i.x = ev.x),
(SELECT array_agg(x ORDER BY id) FROM ids)
FROM (SELECT t.*,
(row_number() OVER (PARTITION BY t.bkey))::int AS bidx
FROM provsql_reachability_edges_tmp t) e
JOIN ids iu ON iu.x = e.u
JOIN ids iv ON iv.x = e.v
$sql$, sel_probs, sel_bkeys, sel_bidx)
INTO sources, destinations, tokens, probabilities, block_keys,
block_indices, extra_ids, vertices
USING extra_vertices;
DROP TABLE provsql_reachability_edges_tmp;
DROP TABLE IF EXISTS provsql_reachability_support_tmp;
END
-- No SET search_path: the deparsed edge subquery (and the regclass
-- rendering) must resolve against the caller's search_path; the ProvSQL
-- calls above are schema-qualified instead.
$$ LANGUAGE plpgsql SET client_min_messages = warning;
SELECT reset_constants_cache();