How parity is verified¶
R's RNG cannot be reproduced by NumPy value-for-value, so every quantity is split in two.
Deterministic quantities — densities, distribution functions, generators and
their derivatives, Pickands functions, τ, ρ, β, λ, iTau/iRho, Rosenblatt
transforms, the Kendall function, and every goodness-of-fit statistic — are
compared against fixtures generated by running the real R package, to ~1e-10.
Everything RNG-dependent uses R's own draws, exported to CSV, so that fitted parameters, standard errors and test statistics all become deterministic and comparable. Only bootstrap p-values remain stochastic, and those are checked statistically.
Fixtures are committed, so CI never needs R. Regenerate them with:
brew install r
Rscript -e 'install.packages("copula", repos="https://cloud.r-project.org")'
bash tools/rgolden/run_all.sh
pytest -q -m golden
Standard errors are calibrated, not merely compared¶
The only honest test of a standard error is whether it predicts the actual spread of the estimator. Measured against Monte-Carlo sampling distributions:
| estimator | mean SE ÷ empirical SD |
|---|---|
mpl (Clayton, n=1000) |
1.004 |
itau (Gaussian, n=1500) |
1.030 |
irho (Gaussian, n=1500) |
0.962 |
itau (Frank, n=800) |
1.030 |
itau (Clayton, d=4, n=600) |
0.967 |
Where R is the less accurate side¶
Each of these is pinned by a regression test, so a future fix in R is noticed rather than silently diverged from.
| R | rcopula |
|
|---|---|---|
| Spearman ρ for Clayton/Gumbel | ~1e-3 | ~1e-11 |
| Elliptical CDF, d = 4–5 | ~1e-4 (Miwa) |
~5e-9 |
Frank iTau |
self-inconsistent at the 10th digit | consistent |
| Plackett τ at independence | −0.000257 | 0 exactly |
| t-EV τ at negative ρ | 0.1225 | 0.0211 (true) |
var.icor for Frank |
0.128 | 0.289 (truth ≈ 0.280) |
pmvt with non-integer df |
refused | supported |
rho() for the t copula |
NA |
computed by quadrature |
The Frank one is the most consequential: R reports less than half the true sampling standard deviation, so confidence intervals built from it are half as wide as they should be.
Things R declines that this computes¶
- Spearman's ρ for the t copula. τ is shared by all elliptical copulas; ρ is not — \((6/\pi)\arcsin(\rho/2)\) is specific to the Gaussian. At ρ = 0.5, df = 4 the true value is 0.4690, not 0.4826.
- A variance for
joeCopulafitted byitau. - Standard errors for
itau/irhoabove two dimensions.
Validation that does not mention R at all¶
Agreeing with R is a useful check and an insufficient one: it would also be
satisfied by faithfully reproducing an R bug. tests/test_literature.py asks a
different question — does this agree with the definition? — and consults no
fixture at all. 240 checks, each computing a quantity from its published
definition by a completely different route than the package uses:
| quantity | definition evaluated | agreement |
|---|---|---|
| Spearman's ρ | \(12\iint C\,du\,dv - 3\), tensor Gauss–Legendre | ~1e-11 |
| Kendall's τ | \(1 - 4\iint \partial_u C\,\partial_v C\,du\,dv\) | ~1e-9 |
| Blomqvist's β | \(4C(\tfrac12,\tfrac12) - 1\) | exact |
| λ_L, λ_U | the limits along the diagonal, as sequences | converging |
| \(W \le C \le M\) | on a 1600-point grid | exact |
Every bivariate family is run through all of them, including the structural constructions (rotated, Khoudraji, mixture), because a definition that holds for six families and fails for the seventh is exactly what a hand-picked example misses.
Alongside those, the identities that connect different parts of the theory and would break loudly if either side were wrong:
- Gamma frailty gives Clayton and positive-stable frailty gives Gumbel — built by sampling the frailty and dividing exponentials by it, then measuring τ. The generator is checked directly as the frailty's Laplace transform.
- \(\tau = \frac{2}{\pi}\arcsin\rho\) for every elliptical copula at any degrees of freedom, while ρ_S demonstrably is not shared — the two measures are often described as interchangeable and here they are not.
- The Nataf transform is a Gaussian copula, so structural reliability inherits zero tail dependence.
- Mutual information is copula entropy, against \(-\tfrac12\log(1-\rho^2)\).
- The nonparanormal identity \(\Sigma = \sin(\pi\tau/2)\) under two brutal but increasing marginal transforms.
Two constants are pinned against mpmath at 25 digits rather than against
anything in this package:
| constant | value | note |
|---|---|---|
| ρ_S of Frank, τ at θ = 5 | 0.45670095816011690 | Debye integral at 30 digits |
| ρ_S of Clayton at θ = 2 | 0.68223383328065629 | \(12\iint C - 3\) at 25 digits |
The second doubles as a cross-family regression test: Clayton and Gumbel have the same Spearman ρ at τ = ½, verified equal to 20 digits in arbitrary precision despite different CDFs, different tails and different generators. Their ρ curves cross exactly there, and the test checks the crossing is transversal — unequal at τ = 0.45 and 0.55 — so it cannot be passed by two implementations that have collapsed into each other.