Espartaco

“Is that to say we are against Free Trade? No, we are for Free Trade, because by Free Trade all economical laws, with their most astounding contradictions, will act upon a larger scale, upon the territory of the whole earth; and because from the uniting of all these contradictions in a single group, where they will stand face to face, will result the struggle which will itself eventuate in the emancipation of the proletariat.”

Karl Heinrich Marx · Marx-Engels Collected Works, Vol. VI, p. 290

EnglishEspañol

Tag: triangular arrays

  • Marx, Adam Smith, and the Law of Large Numbers

    Marx, Adam Smith, and the Law of Large Numbers

    A new research uses probability theory — and sixty years of U.S. economic data — to test one of the most consequential (and most overlooked) assumptions in political economy.


    The Assumption Hiding in Plain Sight

    If you’ve ever read a Marxist analysis of how profits equalize across industries, you’ve probably encountered something called the average rate of profit. The idea is straightforward: competition between capitalists drives different rates of profit in different sectors toward a common, system-wide average. This is one of the pillars of Marx’s theory of value in Capital, Volume III.

    But there’s a quieter assumption underneath this one — so quiet that most discussions never mention it explicitly. To arrive at a uniform profit rate, Marx first assumes a uniform rate of surplus value across all productive sectors. In plainer terms: he assumes that the degree to which workers are exploited — the ratio of unpaid labor to paid labor — is roughly the same everywhere, whether you work in steel manufacturing, food processing, or textiles.

    Adam Smith proposed this idea before Marx. Smith argued that if one job were obviously more exploitative (in the sense of yielding far more unpaid surplus per dollar of wages paid), workers and capital would flow toward or away from it until the differences vanished. Marx adopted this observation and, as scholar Jonathan Cogliano notes, elevated it to “the status of a central economic law” within his framework.

    Yet the assumption has been challenged from multiple directions — Marxist and non-Marxist alike. Is it actually justified? Or is it a convenient simplification that distorts our understanding of how capitalism works?

    José Mauricio Gómez Julián, of the Universidad Latina de Costa Rica, decided to approach the question from an unexpected angle: probability theory. His paper, published in Ciencia Económica (2022), asks whether the mathematical law that would need to hold for this assumption to be valid actually does hold — and then checks the answer against six decades of real-world data from the United States.


    The Mathematical Backbone: The Law of Large Numbers

    If you’ve taken any statistics course, you’ve likely met the Law of Large Numbers (LLN). It tells us that as you observe more and more instances of something random — coin flips, dice rolls, stock returns — the average of those observations settles down toward the true expected value.

    There are two versions:

    • The Weak Law (WLLN): With enough observations, the sample average is probably close to the expected value.
    • The Strong Law (SLLN): With enough observations, the sample average is almost certainly equal to the expected value — a much stronger guarantee.

    Gómez Julián’s insight is this: if you think of each productive sector of the economy as a random variable representing that sector’s rate of surplus value, then the LLN tells you what happens to the average across sectors as the number of sectors grows large. In mathematical language:

    • Strong Law: The probability that the average surplus-value rate across sectors equals the global expected value, in the limit, is exactly 1.
    • Weak Law: The probability that the average deviates from the global expected value by more than any tiny amount shrinks to zero.

    If either version holds, you get the result Marx needs: across a sufficiently large number of sectors, the rates of surplus value converge to a common value — uniformity, or at least a powerful tendency toward it.


    The Catch: Independence and Identical Distributions

    Here’s where things get interesting — and where the classical LLN hits a wall.

    The textbook version of the LLN requires two conditions:

    1. Independence: The random variables (sectoral surplus-value rates) must be statistically independent of each other.
    2. Identical distribution: Each variable must follow the same probability distribution.

    Neither condition holds for the real economy. And Gómez Julián is admirably upfront about this. Sectors are deeply intertwined — the steel industry depends on mining, manufacturing depends on steel, services depend on consumer spending powered by manufacturing wages. The idea that one sector’s surplus-value rate has no relationship to another’s is economically unrealistic. Furthermore, different industries have different cost structures, different labor intensities, and different technologies. There is no reason their surplus-value rates should follow the same statistical distribution.

    So does this kill the argument? Not at all. In fact, it’s the most intellectually interesting part of the paper.


    Non-Classical Varieties: When the Rules Relax

    Over the past several decades, mathematicians and econometricians have developed non-classical versions of the LLN that weaken or entirely drop the independence and identical-distribution requirements. Gómez Julián surveys several of these:

    • Li, Rao, and Wang (1995) showed the LLN holds for random variables arranged on a lattice structure under certain conditions — a structure that, as it happens, economic data naturally exhibits.
    • Adler and Rosalsky (1987) proved the law for weighted sums of independent, identically distributed random variables belonging to a normalized sum, generalizing the classical case.
    • Chen and Sung (2016) extended those results further: the variables no longer need to be identically distributed. They only need to be “stochastically dominated” by a single random variable, with certain weighting conditions.
    • Sung (2011) showed that the strong law can hold even when variables are dependent on each other, provided their probability moments (roughly, their averages and variability) satisfy certain finiteness conditions.

    The crucial point: these results collectively tell us that the LLN’s convergence conclusion can survive even when the classical assumptions are substantially violated — which is exactly the situation with sectoral surplus-value rates.

    Gómez Julián argues that the economic dynamics described by Smith — workers and capital moving between sectors in response to unequal advantages — are precisely the kind of compensatory dependence mechanism that these non-classical versions accommodate. The variables aren’t independent, but their dependence is structured in a way that still drives convergence.


    What the Data Actually Shows

    The theoretical argument is compelling, but Gómez Julián doesn’t stop there. He turns to sixty years of U.S. data (1960–2020), sourced from the Bureau of Economic Analysis (BEA), to see what the empirical evidence says.

    He calculates sectoral surplus-value rates using macroeconomic data on gross operating surplus (representing surplus labor time) and employee compensation (representing necessary labor time), following a standard operationalization of Marx’s categories. After carefully determining which sectors qualify as “productive” in the Marxist sense — a nontrivial task, since the service sector includes activities with very different relationships to surplus-value production — he arrives at 36 productive sectors.

    Here’s what the statistical analysis found:

    Finding 1: No Identical Distributions

    A probability distribution fitting exercise (using the Bayesian Information Criterion) revealed that the 36 sectors’ surplus-value rates follow a patchwork of different distributions — Log-Normal, Cauchy, Uniform, Weibull, and Logistic — with none following a normal distribution. The identical-distribution requirement of the classical LLN is not met.

    Finding 2: No Statistical Independence

    A Pearson correlation analysis across all 630 possible sector pairs yielded a mean correlation of about 0.08 and a median of about 0.14. While these may look small, a deeper cut reveals that roughly 40% of sector pairs have correlations of 0.3 or above — a level that’s practically meaningful. The sectors are not independent. This makes intuitive sense: industries are connected through supply chains, labor markets, and shared macroeconomic conditions.

    Finding 3: Differences Tend Toward Zero

    This is the key finding. When Gómez Julián computed the differences between each sector’s surplus-value rate and the global average (both the mean and the median), he found that these differences exhibit a strong tendency toward reciprocal nullification — positive differences roughly cancel out negative ones. The sum of all differences relative to the global mean was essentially zero (on the order of 10⁻¹⁴). The mean of differences relative to the global median was 0.0012 — vanishingly small.

    Distributional fitting on these differences revealed they follow a Cauchy distribution (when measured against the global mean) or a uniform distribution (against the global median), with the medians of these distributions sitting very close to zero.

    In plain language: sectors deviate from the average in different directions, and those deviations largely cancel each other out.


    Why This Matters

    Gómez Julián’s paper makes two types of contributions that are worth distinguishing:

    For Marxist political economy: If the uniformity assumption holds — even approximately, even as a tendency rather than an iron law — then a large body of research on the long-run behavior of the average rate of profit, both within countries and across the global economy, is on sounder footing than critics have suggested. Researchers studying the tendency of the rate of profit to fall (or not) can continue to work without needing to explicitly model sector-by-sector differences in exploitation rates, at least for aggregate, long-run analyses.

    For probability theory and economics: The paper demonstrates a productive intersection between a specific question in political economy and the deep mathematics of convergence theorems. It shows that the non-classical LLN theorems aren’t just abstract curiosities — they have direct relevance to understanding real economic phenomena. The structured dependence between economic sectors isn’t a bug that invalidates the mathematical framework; it’s a feature that the right version of the framework already accounts for.


    A Few Honest Caveats — And Why They No Longer Apply

    The original 2022 paper was unusually transparent about its limitations, and that transparency is one of its strengths. Rather than forcing the data into inappropriate statistical procedures, it openly acknowledged where the available inferential tools broke down.

    At the time, three important caveats remained.

    First, formal hypothesis testing had to be abandoned.

    The reason was purely statistical rather than economic. Classical inferential procedures—Student’s t tests, Wilcoxon tests, and most conventional non-parametric alternatives—are built on assumptions that the data simply did not satisfy. Sectoral surplus-value rates are neither independent nor identically distributed. They are linked through supply chains, technological change, capital mobility, and macroeconomic shocks. Even bootstrap procedures could not fully solve the problem because ordinary resampling may weaken dependence between resamples while leaving the internal dependence structure fundamentally unchanged. Consequently, the 2022 paper relied primarily on descriptive statistics together with probability-theoretic arguments instead of formal significance testing.

    Second, the classification of productive sectors inevitably involved theoretical judgment.

    Although the paper carefully justified the inclusion and exclusion of economic activities using Marxian categories and modern national accounting, reasonable scholars could still debate where certain services belong within the circuit of capital.

    Third, the empirical evidence came exclusively from the United States.

    The descriptive regularities were remarkably strong, but demonstrating that the same convergence mechanism operates under different institutional settings naturally remained an empirical question.

    Those were genuine limitations in 2022.

    Today, however, the first—and arguably the most important—of them has largely been overcome.

    A much more comprehensive methodological paper (Gómez Julián, 2026; SSRN 5172185) develops an entirely new inferential framework specifically designed for exactly the type of data that made the original analysis difficult: dependent, heterogeneous, and unbalanced observations. Instead of trying to force classical statistical tests to work outside the assumptions under which they were derived, the newer paper constructs hypothesis testing from the ground up for this class of problems.

    The key innovation is recognizing that the convergence of sectoral surplus-value rates is fundamentally a law-of-large-numbers problem under dependence, not an independent-samples problem. The framework therefore combines three complementary asymptotic structures—triangular arrays (TAC), correlation-weighted sums (WSC), and mixingale processes (MPC)—which respectively model hierarchical dependence, contemporaneous intersectoral dependence, and temporal dependence. Rather than treating these as competing approaches, the paper proves conditions under which they become metrically equivalent and therefore support the same inferential conclusions.

    The inferential consequences are substantial.

    Instead of abandoning significance testing because dependence invalidates classical procedures, the new framework explicitly extends the Neyman-Pearson paradigm to dependent observations, derives dependence-aware confidence regions, establishes rigorous Type I error control under strong-mixing assumptions, and integrates Bayesian and frequentist inference into a single coherent architecture. Robust procedures—including fixed-b heteroskedasticity-and-autocorrelation-robust inference, block bootstrap techniques that preserve dependence, adaptive conformal inference, composite and Whittle likelihoods, and hierarchical Bayesian estimation—serve as mutually reinforcing validation mechanisms rather than isolated alternatives.

    In other words, what had been acknowledged as a methodological limitation in the 2022 paper became the central research question of the later work.

    Rather than concluding that inference was impossible under dependence, the subsequent research asks a more fundamental question: what should hypothesis testing look like when dependence is the normal state of the data rather than an exception? The result is a unified inferential framework specifically intended for datasets that violate the assumptions of classical statistics—precisely the situation encountered with sectoral surplus-value rates.

    The other caveats remain, although they are considerably less problematic than before. The classification of productive sectors continues to depend on theoretical interpretation, because that issue belongs to political economy rather than statistics. Likewise, expanding the empirical analysis to additional countries remains a desirable avenue for future research. Yet the principal statistical objection—that no valid inferential procedure existed for dependent sectoral data—has now been directly addressed through a purpose-built mathematical framework.

    Looking back, the 2022 paper can therefore be read as identifying an important statistical obstacle, while the later work attempts to remove it. Together, the two papers form a coherent research program: first demonstrating that the convergence hypothesis is theoretically plausible and descriptively supported, and then developing the inferential machinery required to test that hypothesis rigorously without relying on unrealistic assumptions of independence or identical distributions.


    Gómez Julián, J.M. (2022). Sobre la validez del supuesto de uniformidad en las tasas de plusvalía sectorial desde la teoría de las probabilidades. Ciencia Económica, 11(17). DOI: 10.22201/fe.24484962e.2022.11.17.2

    Gómez Julián, J.M. (2026). Hypothesis Testing for Dependent Variables with Unbalanced Data: A Unified Framework: Theory, Robustness, and Software. SSRN Electronic Journal. DOI: 10.2139/ssrn.5172185.

  • When the Textbook Test Fails: How HTDV Brings Rigor to Dependent, Unbalanced Data

    When the Textbook Test Fails: How HTDV Brings Rigor to Dependent, Unbalanced Data

    You can also find this library at CRAN and download it directly from R and RStudio.

    The problem hiding in your data

    Picture a straightforward question: Is the average inflation rate in the United States meaningfully different from zero? You have monthly data going back decades. A classical t-test would seem like the natural tool — and it would be quietly, systematically wrong.

    The reason is that inflation figures do not bounce around independently. January’s number carries information about February’s. This autocorrelation corrupts the standard error that the t-test relies on, inflating the false-positive rate well beyond the nominal 5% you think you are signing up for. The same problem afflicts yield spreads, stock returns, sectoral profitability, regional employment — virtually any real-world time series you might want to compare.

    Now make it harder. The two groups you are comparing have different sample sizes — one sector has twenty years of data, another only five. The data may have heavy tails you cannot rule out. And your sample is finite, which means the asymptotic guarantees printed in your econometrics textbook are promises that may not have been kept yet.

    This is the terrain HTDV was built for. Short for Hypothesis Testing for Dependent Variables with Unbalanced Data, HTDV is an R package that answers a deceptively simple question — do these dependent, possibly unequally-sized samples come from the same population? — under the worst combination of conditions an applied statistician is likely to encounter.

    The central idea: triangulation, not trust

    The most common approach to dependent data is to reach for a single robust method — a heteroskedasticity-and-autocorrelation-consistent (HAC) standard error, say, or a block bootstrap — and hope it is calibrated. HTDV takes a structurally different stance: run three independent inferential methods in parallel and expose the disagreement between them as a signal.

    The three layers are:

    1. A hierarchical Bayesian fit via Hamiltonian Monte Carlo (HMC), implemented in Stan. This layer builds a full probability model of the data-generating process, places weakly informative priors on the dependence parameters, and produces a posterior distribution for the quantity of interest.
    2. A fixed-bandwidth HAR Wald test in the frequentist tradition of Kiefer and Vogelsang (2005). Rather than letting the bandwidth grow with the sample in the usual way, it holds the bandwidth at a fixed fraction of the sample size. This produces a non-standard asymptotic distribution that is better calibrated in finite samples than the conventional chi-square approximation.
    3. A stationary block bootstrap (Politis and Romano, 1994) with automatic block-length selection (Patton, Politis, and White, 2009). This resamples the data in blocks long enough to preserve the dependence structure, then constructs confidence intervals from the resampled distribution.

    A fourth, distribution-free layer — adaptive conformal inference (Gibbs and Candès, 2021) — is available for online prediction settings where no parametric assumption is palatable.

    The logic is forensic. Where all three layers agree, your conclusion is robust. Where they disagree, the pattern of disagreement tells you something specific about your data. If the Bayesian interval is dramatically wider than the HAR or bootstrap interval, your series likely has strong temporal persistence, and the asymptotic critical values that HAR and bootstrap rely on are losing their reliability. That gap is not a bug — it is the most informative thing the framework can show you.

    Why a single method is not enough

    It is fair to ask: if the Bayesian layer is the most reliable, why not just use it and discard the others? The answer is that each layer has a regime where it is the appropriate tool, and the framework’s job is to make the regime visible.

    HAR inference is computationally cheap — sub-second on typical data — and well-calibrated when persistence is low to moderate and sample sizes are large enough for asymptotics to bite. The block bootstrap shares those advantages while making fewer distributional assumptions. The Bayesian layer is the most computationally expensive (each fit can take tens of seconds) but is the only one that maintains nominal calibration under strong persistence at finite sample sizes, because it models the dependence explicitly rather than relying on asymptotic corrections.

    The package ships with a pre-registered factorial Monte Carlo study — 1,024 cells crossing sample size, autocorrelation, tail heaviness, imbalance ratio, and location shift, with 500 replications per cell across all three inferential layers — and the results are unambiguous. The Bayesian layer holds nominal size (mean 0.056 against a target of 0.05) and nominal coverage (mean 0.944 against a target of 0.95) across the entire grid. HAR and bootstrap, by contrast, inflate dramatically in the worst corners: under strong persistence and small samples, HAR’s empirical rejection rate under the null reaches 0.60, and its coverage drops to 0.29.

    The narrowness of the HAR and bootstrap intervals in those corners is not a sign of precision. It is a sign of miscalibration — the intervals are confidently wrong.

    The theory that holds it together

    Running three different methods on the same data and comparing the answers is sound practice, but it raises a mathematical question: under what conditions are the three methods even addressing the same inferential target? A Bayesian posterior on a triangular-array likelihood and a HAR-Wald statistic on a mixingale process are, on their face, different objects.

    HTDV’s theoretical backbone is a metric equivalence theorem that resolves this concern. The framework identifies three structurally distinct ways real-world data can violate the independence assumption — each corresponding to a different law-of-large-numbers regime:

    • Triangular Arrays Convergence (TAC): information accumulates through hierarchical aggregation. Think of input-output tables disaggregated into ever-finer sectors, where each “row” of the array adds more observations.
    • Weighted Sums with Correlation (WSC): the observations share a cross-sectional covariance structure. Regional markets that move together, trade flows between linked economies.
    • Mixingale Process Convergence (MPC): temporal memory that decays smoothly over time. Forecast errors, model residuals, prediction intervals that gradually lose contact with the past.

    The theorem proves that, under α-mixing with polynomial decay rate γ > 1 and finite moment conditions, these three regimes induce strictly equivalent metrics on the space of hypothesis-testing problems. The equivalence comes with explicit, computable finite-sample constants — exposed by the function htdv_equivalence_constants() — that tell you the maximum slack when translating a conclusion from one regime to another. For typical parameter values (γ = 2, q = 6, n = 500), the conversion slack is about 18%, a margin that is usually irrelevant for a hypothesis-testing decision.

    This is what makes the three-layer architecture mathematically well-defined rather than merely pragmatic. Without the equivalence theorem, comparing a Bayesian result on a TAC dataset with a HAR result on a WSC dataset would be comparing apples and oranges. The theorem certifies that the metrics are coercible to one another with computable error.

    The dependence assumption, plainly

    The framework assumes that the data are α-mixing with polynomial decay — meaning that the statistical dependence between observations dies off as they get farther apart in time, and it does so fast enough (at a rate faster than 1/k) that the long-range correlations are summable. This is a mild condition satisfied by most stationary time series in econometrics and finance, including ARMA processes, GARCH models, and a broad class of Markov chains.

    It is not satisfied by long-memory processes (where dependence decays more slowly than 1/k) or by unit-root processes (where dependence does not decay at all). The framework is honest about these limitations: it will fit near-unit-root data, but the posterior will widen correspondingly — which is the correct answer. For explicit unit-root testing, the standard ADF or Phillips-Perron tools remain the right choice.

    The Bayesian engine

    The hierarchical Bayesian core fits Stan models via the No-U-Turn Sampler (NUTS), the state-of-the-art Hamiltonian Monte Carlo variant. The models are parameterized around an AR(1) structure — the mean θ, the autocorrelation φ, and the innovation scale σ — with hierarchical priors on the dependence nuisance parameters that are weakly informative enough to respect admissible ranges without overwhelming the data.

    Five likelihood backends are available, corresponding to the three convergence regimes plus two parametric likelihood families: the Whittle likelihood (which works in the frequency domain, comparing the observed periodogram to a theoretical spectral density) and the composite likelihood (which works in the time domain, combining conditional densities over short blocks). Both are well-established in the time-series literature; the choice between them depends on whether you have more confidence in your spectral model or your conditional density model.

    A distinctive feature is the Berger-robust envelope — a method for combining posteriors across multiple fitted models into a single, wider posterior that hedges against the worst-likelihood-specification scenario. If you are unsure whether the Whittle or composite likelihood better describes your data, the envelope gives you an inferential answer that is honest about that uncertainty rather than forcing an arbitrary choice.

    After sampling, every fit must pass a five-gate diagnostic check before its posterior is admissible: split-R̂ below 1.01, bulk and tail effective sample sizes above 400, zero post-warmup divergences, and energy Bayesian fraction of missing information (E-BFMI) above 0.3. These are the standard HMC convergence diagnostics from the Stan ecosystem, enforced as a gate rather than offered as a suggestion.

    The validation: visible in the data

    The most compelling aspect of HTDV is that it does not merely claim to be well-calibrated — it ships the evidence. Two validation datasets are bundled with the package.

    The first is the factorial simulation described above, with its 3,069-row summary table accessible as a package dataset. The headline finding — that the Bayesian layer is the only one maintaining nominal calibration across the full design — is not an assertion but a reproducible fact. The full study took 31 hours on a 16-core workstation; the scripts to regenerate it from scratch are shipped in the package repository.

    The second is a set of three external benchmarks against published references on public-source data:

    • Post-1984 US CPI inflation, compared against Stock and Watson (2007).
    • Shiller’s log-CAPE ratio, compared against Campbell and Shiller (1998).
    • The US–Canada 10-year yield differential, compared against the naive iid Welch baseline.

    All three layers reproduce all three references with agreement in every case. But the width of the agreement tells the real story. The interval widths scale monotonically with the persistence of the underlying series. At moderate persistence (φ ≈ 0.45, the inflation series), the Bayesian interval is actually narrower than HAR — 0.81 times its width. At high persistence (φ ≈ 0.97, the CAPE series), the Bayesian interval is 2.8 times wider. At near-unit-root persistence (φ ≈ 0.99, the yield differential), it is 15 times wider.

    This gradient is the framework’s central empirical finding. Both layers are technically asymptotically valid. Only the Bayesian layer accounts honestly for the finite-sample uncertainty inflation that occurs as φ approaches 1. The HAR and bootstrap intervals do not widen because they know more — they fail to widen because their asymptotic critical values have not yet caught up with the data.

    When to use it — and when not to

    HTDV is the right tool when your data are time-dependent or spatially dependent, when your samples are of unequal size, when you suspect heavy tails but cannot rule them out, and when you need an inferential answer (a test or an interval) rather than a prediction. It is particularly valuable when the stakes are high enough that you want your conclusion to survive methodological scrutiny — the framework ships its own validation evidence precisely so that a reviewer can interrogate the calibration claims rather than taking them on faith.

    It is the wrong tool when your data are genuinely independent with finite variance — classical methods are simpler, equivalent, and faster. It is also not designed for long-memory processes, explicit unit-root testing, structural breaks (unless you segment the sample first), or forecasting. The framework is built for hypothesis testing and parameter estimation under uncertainty, not for predictive accuracy.

    An open architecture

    The package exposes its full infrastructure: the simulation engine (htdv_simstudy()), the equivalence constants calculator, the diagnostic suite, the posterior-predictive checks on dependence statistics, and the decision tools — ROPE-based decisions (Kruschke, 2018), bridge-sampling Bayes factors, WAIC and leave-future-out cross-validation, and predictive stacking (Yao, Vehtari, Simpson, and Gelman, 2018). Every function is documented with its underlying reference, so the user can trace any method back to its source.

    The complete function reference, mathematical foundations, tutorial walkthroughs (oriented toward novices, applied statisticians, and mathematicians respectively), and the full validation narrative are in the HTDV Wiki on GitHub. The package is installed with a single command — remotes::install_github("IsadoreNabi/HTDV") — and requires rstan as its only hard dependency.

    The larger point

    HTDV embodies a methodological philosophy worth stating explicitly: when no single inferential method is universally valid in the finite-sample regime, the honest response is not to pick the best one and hide its limitations, but to run several and make the disagreement visible. The framework’s value is not that it always gives you a narrower interval or a more powerful test. Its value is that it shows you — concretely, quantitatively — where your inference is on solid ground and where it is standing on asymptotic ice.

    The validation evidence makes this concrete. In 98% of the simulation cells, the Bayesian layer alone passes the calibration benchmarks. The HAR and bootstrap layers pass in the regime where asymptotics have bitten — low persistence, large samples — and fail predictably outside it. The framework does not hide that failure. It turns it into a signal.

    That signal is the product.


    HTDV is released under the MIT license. The companion paper, full validation vignette, and reproducibility scripts are available at github.com/IsadoreNabi/HTDV. For the complete mathematical foundations, function reference, and tutorials, see the project wiki.