Mixalot
The content on this page was written by AI under human supervision.
Mixalot is a Python package that computes the Bayesian evidence of finite mixture models for count data exactly, as a rational number, together with the exact posterior probability of each possible number of components. You give it a vector of integer counts (how many of N observations fell in each of k categories) and a number of components g. It returns a Python Fraction; for very large problems you print its logarithm rather than its digits. The package collects the evaluators written for the paper Exact Bayesian evidence for mixture models behind one import, with a self-test suite and a checksum verification of every bundled file.
What it does
The evidence (or marginal likelihood) of a model is the probability of the data averaged over the prior on its parameters. The ratio of two evidences is a Bayes factor, the Bayesian way to decide between "one population" and "two". For a mixture of g categorical distributions on k states it is an integral over every component's probabilities and the mixing weights, in practice almost always estimated by Monte Carlo. Mixalot evaluates it exactly: with uniform Dirichlet priors the integral collapses to a finite sum over the ways the N observations can be allocated among components, computed by dynamic programming in big-integer arithmetic. The default convention is a fixed observation sequence (no multinomial coefficient) with Dirichlet(1) priors; measure="lebesgue" divides out the Dirichlet normalizing constants $(g-1)!\,((k-1)!)^g$ instead. A separate bundled evaluator written in the conventions of Lin, Sturmfels and Xu (arXiv:0805.3602) reproduces their printed benchmark values (core.lsx_example).
mixalot.Z, mixalot.gstar and mixalot.bayes_factor cover small and medium problems. Larger ones go through the bundled evaluators in mixalot.engines, each specialized in one direction. There are closed forms in harmonic numbers for one binary variable, a route for thousands of states, one for any number of components, and the exact infinite-component (Dirichlet-process) limit. Two-way contingency tables are computed modulo many primes and reassembled by the Chinese remainder theorem. A fourth route, the boxwalk subpackage (from mixalot import boxwalk), evaluates the underlying unit-box moment integrals $\int_{[0,1]^n}\prod_k P_k(x)^{u_k}\,dx$ exactly for integer polynomials $P_k$, by a recurrence walk modulo primes followed by rational reconstruction. mixalot.advise() reports the measured reach of the exact routes. With two components, N around $10^5$ is easy; with three, N = 3,000 needs about 1 GB of memory and N = 10,000 about 25 to 30 GB; with four or more, N up to roughly 500 to 700.
Two evaluators handle structured data. seg_v1 takes an ordered sequence of count vectors (chapters of a book, say) and computes the evidence that it came from g contiguous segments, each with its own profile, plus the posterior on where the boundaries fall. frozen_comp_v1 takes components with known profiles and unknown Dirichlet weights and returns the evidence and the Bayes factor for whether a given component is present. Each has an independently written twin (seg_blind, frozen_comp_blind); the self-tests require the two frozen-component implementations to agree exactly and cross-check seg_v1's exact route against its floating-point route. Three small modules extend this to authorship-style questions on an ordered token sequence. seam gives the exact evidence that one known profile produced a prefix and another the suffix, with the posterior on the changeover point. dcm replaces the multinomial with a Dirichlet-compound-multinomial, so word rates may vary from stretch to stretch. nullcal reads a mixed-versus-pure Bayes factor against the same statistic computed on objects known to be pure, because a mixture's extra freedom also absorbs ordinary rate variation. A Monte Carlo comparison suite is bundled too, for checking a sampler against exact values where both run, along with evaluators for two related likelihoods: Etienne's neutral-biodiversity sampling formula (with its supporting exact, interval-arithmetic and certified maximum-likelihood engines) and the telegraph model of gene expression.
The evidence routes return exact rationals; the exceptions say so in their names or modes (Z_float, lnZ_g_float, the interval mode of the Etienne evaluator, the Monte Carlo suite). A non-integer count raises NonIntegerCountError instead of being rounded. mixalot.verify() re-hashes every bundled file against recorded checksums and raises VendorTamperError, refusing to compute, if anything has changed, is missing or has been added. One limit is statistical: for a bare vector of counts from single draws, a mixture of categoricals is itself a categorical, so the data say little about g and at small N the posterior mostly reflects the prior. The exact number reports that faithfully but cannot remove it. Structure (ordered units, known profiles) is what gives power, and plant, recover_blind and mic_power measure detection power on synthetic data shaped like yours. The package handles categorical count data only.
Examples
Exact evidence and the component-count posterior. The manual's quickstart, against a checkout of the repository:
import sys; sys.path.insert(0, "<your-checkout>/tools/mixalot") import mixalot mixalot.verify() mixalot.Z([4,1,3,2], 2) mixalot.gstar([4,1,3,2], gmax=3) from mixalot.engines import zseries, seg_v1, frozen_comp_v1, w4_blind_gf
verify() returns a report whose vendor_ok list names every bundled file that passed; its source_missing entries refer to the original source trees, which are not distributed, and are expected. Z returns Fraction(1010921, 2330808480000), the two-component evidence for ten observations on four states. gstar returns a dictionary with posterior (one Fraction per g from 1 to gmax, summing to one), map_g, p_ge2, evidence, prior and gmax. Here the posterior is roughly 0.22, 0.35, 0.43 for g = 1, 2, 3, so P(g ≥ 2) is near 0.78 with no value of g singled out, which is the weak-identifiability caution above: ten bare counts barely constrain g.
Command line. Put the counts in a JSON file, either a bare list [3, 1, 4, 1, 5] or an object {"U": [3,1,4,1,5], "gmax": 3}, and run
python3 -m mixalot.cli counts.json [--gmax 4]
The output is one JSON block: U, k, N, gmax, the posterior for each g as an exact fraction with a float beside it, map_g, p_ge2, the evidence for each g, and the advise() text for that N. With no file it prints the input schema.
Exact versus Monte Carlo at one point. The self-test suite and examples/worked_examples.py both run this comparison, two components of one binary variable with counts (50, 50):
from mixalot.engines import load
cf = load("closed_form_1var")
z = cf.Z_closed(50, 50)
est = load("estimators")
r = est.run("m1", 2, [50, 50], "nested", seed=1)
z is the exact evidence as a Fraction; its logarithm is −71.0794. r is a dictionary with logZ_hat, err_est, diagnostics, settings and wall_s; the recorded test log has logZ_hat = −71.0204, 0.6 standard errors from the exact value, so nested sampling is calibrated at this size. python3 examples/worked_examples.py prints this and four more demonstrations: the exact log-evidence at N = 100,000 in roughly ten seconds, 300 states, the Dirichlet-process limit $Z_{\rm DPM}=5/48$ for counts (2, 1), and the exact finite-g correction $g\,(Z_g - Z_{\rm DPM}) = -1/48$. python3 examples/large_kgn_table.py [--quick] writes the full exact-versus-sampler table with timings.
Routines
Top level (import mixalot)
verify(quiet=False)— checksum check of every bundled file; raisesVendorTamperErroron any changed, missing or added file; returns the report.Z(U, g, measure="dirichlet")— exact evidence as a Fraction; stops with an error above N = 2,000,000 and names the scale routes to use.gstar(U, gmax=3, prior=None)— exact posterior over g = 1..gmax; returnsposterior,map_g,p_ge2,evidence,prior,gmax.bayes_factor(U, g1, g2)— Z(U, g1) / Z(U, g2), exact.advise(N=None, g=None)— measured memory and time reach of the exact routes.plant(k, g_true, N, seed, separation="strong"),recover_blind(k=4, N=60, gmax=3)— synthetic counts from a known mixture, and a plant-and-recover loop reporting how oftengstarfinds the true g.core.Zdpm(U, alpha)— exact evidence in the Dirichlet-process limit.core.lsx_example(name)— reproduces a Lin–Sturmfels–Xu benchmark ('swiss','coin10','coin242'); True on agreement.core.p_g1(X, gmax=6)— probability that an ordered list of unit count vectors is a single segment, viamic_power.engines.load(name),engines.<name>— loads a bundled evaluator after re-checking its checksum; the source is executed directly and Python's bytecode cache is never used.boxwalk— the box-integral subpackage,from mixalot import boxwalk(see Boxwalk below).
Authorship modules (from mixalot import nullcal, seam, dcm)
nullcal.calibrate(claim_bf, known_bfs),nullcal.blend_vs_pure_bf(U, gmix=2, Z=None),nullcal.calibrate_objects(U_claim, known_Us, gmix=2, Z=None)— place a mixed-versus-pure log10 Bayes factor within its distribution over known-pure objects (counts, quantiles and a verdictINSIDE-NULL,TAILorABOVE-NULL); the last two compute the factors exactly, throughbigginside the package or through aZcallable you pass (raisesEngineUnavailableotherwise).seam.z_seam(seq, pA, pB),seam.z_seam_brute(seq, pA, pB)— exact evidence that profile A produced a prefix and profile B the suffix of an ordered sequence of category indices, changeover point uniform, plus the posterior mode and quantiles of the changeover; raisesProfileErrorif a profile does not sum to one; the brute version enumerates directly.dcm.profile_scaled,dcm.dcm_pure,dcm.dcm_mixture,dcm.dcm_seam,dcm.fit_kappa— exact Dirichlet-compound-multinomial evidences (pure, token mixture with weight f ~ Beta(1,1), and change-point) with one shared concentration kappa fitted on known-pure objects over a rational grid.
Command line and scripts
python3 -m mixalot.cli counts.json [--gmax G]— posterior over the number of components, as JSON.python3 battery/battery.py [--full]— the self-test suite, including a deliberately corrupted copy that must refuse to run and the boxwalk self-test; writesbattery/BATTERY.txt; exit code 1 on any failure.python3 examples/worked_examples.py,python3 examples/large_kgn_table.py [--quick]— the demonstrations in Examples.mixalot/selftest.py [--planted fixture]— tests fornullcal,seamanddcm, runnable from any directory: fixture checksums, the Federalist known-author reference counts,z_seamagainstz_seam_brute, anddcmagainst direct enumeration over all token assignments; ends with[SELFTEST] ALL PASS; with--planted fixturea deliberately altered fixture must be refused with exit code 3.mixalot/acceptance_no55.py [--planted]— reproduces an independent implementation's outputs for Federalist No. 55 under thedcmmodel from the bundled fixtures (ten values, each to $10^{-30}$) and checks the bundled blend-model values for the same paper;--plantedflips one reference digit and must fail with exit code 1.scripts/stamp_pins.py— maintainer script that rewrites the checksums after a deliberate update of the bundled files.
Bundled evaluators (mixalot.engines.<name>)
bigg.Z_bigg— general-g exact evidence; the routine behindmixalot.Z.bigk2.Z_fast2— two components, thousands of states (print the logarithm; the integers exceed Python's default print limit).closed_form_1var.Z_closed,formula_emitter.emit_formula,formula_emitter_dirichlet.emit_dirichlet— the harmonic-number closed form for two components of one variable, evaluated or written out term by term.zseries.z_ray_1var,z_ray_table_modp,Z_table_modp— the series Z(n·U0), exactly or modulo a prime below $2^{25}$, for one variable or a two-way table.lsx55_exact— script for the 3×3, N = 132 table of Lin, Sturmfels and Xu (their Example 5.5) by many primes and rational reconstruction, about 6.5 minutes and 3.6 GB per prime. It writes into the current directory, so run it from a scratch directory.lsx_direct.run_example,Z_phi,selftest— independent evaluator in the Lin–Sturmfels–Xu conventions; as a script,--selftest,--example coin10or--table '4,2;2,4'.w4_blind_gf.Zg_gf,Zdpm_gf— finite-g evidence with Dirichlet(alpha/g) weights and its exact Dirichlet-process limit.biggen,bigk,w1_brute,w1_collapsed,w4_dpm_limit,f_nu,production_sweep— reference implementations and drivers used in the cross-checks (production_sweep, likemic_power, needsMIXALOT_PILOT_DIR).swap_route.lnZ_g_float,lnZ_dpm_float— float64 log-evidence for k = 2 up to very many components (FFT convolution); a demonstration checked against the exact value at N = 100; it is not a general evaluator.estimators.run(family, k, counts, estimator, seed)— Monte Carlo estimators (hm,bridge,ss,chib,nested) at default settings for one binary variable (m1) or a four-flip binomial (m4); same seed, same answer.annihilator.pf_from_series,find_recurrence_fast— a copy of the Annihilator recurrence finder, for testing an evidence series for a linear recurrence.seg_v1.Z_exact,Z_float,boundary_posterior_exact,boundary_posterior_float,selfcheck— contiguous-segment evidence and boundary posteriors for an ordered list of count vectors, exact or log-domain float;seg_blind.evidence,posterior_g,boundary_posteriorare the independent twin.mic_power.p_g1,sample_units,arm— plants a boundary between two profiles at given unit sizes and counts detections; refuses to load unless the environment variableMIXALOT_PILOT_DIRnames the directory holdingseg_v1.py.frozen_comp_v1.evidence,presence_bayes_factor,frozen_comp_blind.evidence_blind,presence_bayes_factor_blind— exact evidence for a blend of known profiles with Dirichlet weights, and the Bayes factor for including one more profile; the second pair is the independent twin.etienne_evaluate.P_exact,logP_ball— Etienne's sampling formula of neutral biodiversity theory for a species-abundance vector, exact or as a certified interval; as a script,--abund "1,1,2,3,5,8" --theta 7.047958 --m 0.22635923(add--mode exactfor the rational value) or--selftest.ball_engine.P_exact,ball_engine.logP_ball,etienne_oracle.etienne_P,phase2_certify.interval_newton,rf_engine.logP_rf,multisample.multisample_P,multisample_ball,phase2_multisample_eqI,hier_kron.hier_lnP_kron,hierarchical3,hier_ball,hier_eval,gate_engines— the engines behindetienne_evaluate: exact product-tree and interval ("ball") evaluation of the Etienne formula, an independent urn-process implementation used as the reference value, certified maximum likelihood by interval Newton, the equal-immigration multi-sample factorization, and a three-level hierarchical evaluator; loading the last of these runs the exact-versus-reference agreement check (a fraction of a second).phase2_multisample_eqIrun as a script writes a JSON file into the current directory, so run it from a scratch directory.telegraph_evaluate.pmf_certified,loglik_certified,pmf_float,loglik_float,selftest— transcript-count distribution and log-likelihood of the two-state gene-expression model; as a script,pmforloglikwith--mode certified|mp|float.
Boxwalk (from mixalot import boxwalk with tools/mixalot on sys.path; the subpackage mixalot/boxwalk/)
boxwalk.load_spec(path_or_dict)— a problem specification:variables,polynomials(label → {exponent string: integer coefficient}),target_u(label → count),orderandoptions; examplemixalot/boxwalk/examples/jc_quartet.json(fifteen Jukes–Cantor quartet kernels;examples/make_jc_quartet.pybeside it regenerates the file).boxwalk.plan(spec)— the walk plan: class order, mechanism per segment, memory per prime and a program checksum.boxwalk.walk(spec, plan, primes)— the moment window and the value modulo a batch of primes.boxwalk.produce(spec, plan=None, nprimes=8, batch=4, procs=None, outdir=None)— the exact Fraction with aMANIFEST.jsonrecording the checks it passed (walk against dense evaluation at a truncated target, a corrupted program detected, two disjoint prime sets reconstructing the same value); stops with "increase nprimes" rather than return a wrong value.boxwalk.emit_fiber_recurrence(spec, ray_class=None, depth=80),boxwalk.gcrd_reduce— the integer-coefficient recurrence of Z along one exponent direction, verified on a prime not used in the fit, and the shift-operator calculus (fiber.gcrd,lclm,right_divides,op_mul,annihilates_mod) for reducing and combining such operators.boxwalk.verify_manifest(spec, manifest, nfresh=2)— re-checks a recorded manifest by an independent replay: spec hash, every recorded residue, a local reconstruction of the value, fresh primes through dense evaluation (cost estimated first; refuses overmax_cells), and the pivot-program hash when a plan is supplied.boxwalk.home()— the subpackage's directory, wherecli.py,selftest.py,examples/and the referenceSELFTEST.jsonlive.python3 -m mixalot.boxwalk plan SPEC.json [--slab-only]/produce SPEC.json OUTDIR [--nprimes N] [--batch K] [--procs P]/fiber SPEC.json [--ray CLASS] [--depth D]/verify SPEC.json MANIFEST.json [--nfresh N] [--replan] [--max-cells C]— the same from the shell, run fromtools/mixalot/(mixalot/boxwalk/cli.pyalso runs as a plain script from any directory).python3 -m mixalot.boxwalk selftest— boxwalk's own tests: writesSELFTEST.jsoninto the subpackage's own directory unlessBOXWALK_SELFTEST_OUTnames another file; the copy in the repository is the reference record, and the package's main self-test suite runs these tests too.
Used on this site
- Mixture models — the package collects that paper's exact evaluators; the downloadable bundle on the page is this package.
- Disputed authorship — evaluated the blend and change-of-source evidence integrals for the Federalist, Shakespeare and Bible analyses.
- The Voynich Manuscript — the exact mixture evidence behind the page-by-page dialect classification.
Requirements and source
Python 3 with the BootLoops common dependencies (mpmath, sympy, numpy, python-flint) plus scipy and dynesty (pip install scipy dynesty); the HMC-based comparison estimators also import jax. The boxwalk subpackage additionally needs Singular on the PATH (without it the Singular part of its self-test is reported as skipped) and fits recurrences with the package's bundled copy of the Annihilator (BOXWALK_ANNIHILATOR=/path/to/annihilator.py selects another copy, such as the stand-alone tools/annihilator/). Self-tests: python3 battery/battery.py [--full] from the package directory, about two minutes (set MIXALOT_BATTERY_OUT to write the BATTERY.txt record elsewhere); --full adds multi-minute checks, some of which skip unless the environment variables MIXALOT_REFDATA_PILOT, MIXALOT_REFDATA_SWEEP and MIXALOT_PILOT_DIR point at reference data that is not part of the repository. The authorship modules have their own script, mixalot/selftest.py, and boxwalk has python3 -m mixalot.boxwalk selftest. Code: tools/mixalot/ in BootLoops' bootloops-dev repository (GitHub organization BootLoops-ai), released under the MIT license; boxwalk is the subpackage mixalot/boxwalk/ inside it, with its own README. A packaged copy with the license files, and the manual, are also hosted here: mixalot-package.tar.gz, MANUAL.md; that copy uses an earlier layout in which boxwalk is a separate boxwalk/ directory beside mixalot/, with its own guide and the command line python3 cli.py ... run from that directory.