biplotpy

Biplots for Python, with their graphics, ported from R and validated numerically against the R packages they come from: MultBiplotR (Vicente-Villardón, Vicente-González & Frutos-Bernal), GGEBiplots, agricolae, smacof, vegan, ade4, FactoMineR, ca, pls, psych, GPArotation, DistatisR, biplotbootGUI, MANOVA.RM and heplots.

import biplotpy as bp
iris = bp.load_iris()
hj = bp.hj_biplot(iris.iloc[:, :4], dim=3)
hj.summary()["col_qualities"]
hj.plot(groups=iris["Species"])
cb = bp.canonical_biplot(iris.iloc[:, :4], iris["Species"])
cb.manova; cb.plot(radius="bonferroni")

Validation

tests/test_vs_multbiplotr.py (fixtures generated by tests/fixtures/make_fixtures.R) checks, up to the sign of each axis: coordinates, eigenvalues, inertia, scale factor, structure correlations and contributions of 12 classical biplots (α = 0, 0.5, 1, HJ; five transformations); 11 initial transformations; the canonical biplot (group, variable, individual markers, weights, eigenvalues, Wilks F, ANOVAs, four radii, contributions, qualities); supplementary variables; calibrated scales (port of R's pretty()); MANOVA tables against R's summary.manova; LDA accuracy against MASS::lda.

Phase 2 (tests/test_phase2.py, fixtures from fixtures/make_phase2.R):

Phase 3 (tests/test_phase3.py, fixtures from fixtures/make_phase3.R):

Note on DISTATIS with the nuclear norm: DistatisR sums the square roots of the eigenvalues larger than 2 eps; a numerically zero eigenvalue of a double-centred matrix can pass that threshold in one library and not in the other, which moves the results in the 8th decimal.

Simulated recovery: STATIS recovers a common configuration and down-weights a pure-noise table; unfolding recovers an exact configuration (stress < 1e-3); logistic biplots recover the true probabilities, with correlation increasing with n and p (0.81 at 100 x 10, 0.97 at 400 x 80 for the binary model).

Caveat: the deviance "tests" of the logistic biplots use df = dim as in MultBiplotR, but the row scores are estimated from the same data, so the p-values are optimistic; use them as descriptive fit indices. Joint ML with a very small penalization overfits (near separation); use penalization >= 0.5 for small tables.

Phase 4 (tests/test_phase4.py, fixtures from fixtures/make_phase4.R):

Phase 5, version 0.4.0 (tests/test_phase5.py, fixtures from fixtures/make_phase5.R):

Notes. (1) AddOrdVars2Biplot converts the factors to integers and RidgeOrdinalLogistic then stops ("The dependent variable must be ordinal"); the reference is computed with that conversion disabled. (2) In R >= 4.2 PLSR1Bin stops inside RidgeBinaryLogistic (class(y) == "formula" on a matrix); the reference uses class(y)[1]. Its algorithm 2 refers to undefined objects and is not ported. (3) BinaryPLSR alternates conjugate-gradient fits from random starts; biplotpy draws the same starts (R's generator) and runs R's CG algorithm, but the path amplifies differences in the last digits of the arithmetic (BLAS summation order), so results agree to about 1e-4. (4) BootstrapSmacof stops (it reads Sol$Dhat, which SMACOF does not return); bootstrap_mds resamples the residuals between disparities and configuration distances, truncates negative dissimilarities at 0, refits from the original configuration and aligns by orthogonal Procrustes; the R reference runs the same scheme on smacof::mds. MultBiplotR::SimpleProcrustes, used by BootstrapDistance, applies the transpose of the optimal rotation (it fixes reflections but turns small rotations the wrong way); bootstrap_pcoa uses the optimal rotation by default and procrustes="MultBiplotR" reproduces R. In the constrained binary logistic biplot the deviance tests use dim degrees of freedom (R uses 1) and the percentages of correct classification threshold the probability at 0.5 (R thresholds the linear predictor).

Graphics

Every method has a plot (examples/gallery.py renders all 49 figures into examples/gallery/): classical biplots with arrows or calibrated axes, group concentration ellipses, convex hulls, alpha-bags and densities (ggbiplot, factoextra, biplotEZ), quality colouring, correlation circle, scree and contribution plots, 3-D biplot (bpca), bootstrap regions; canonical/CVA biplot with confidence circles, classification regions and calibrated axes (biplotEZ); HE plots (heplots); binary/ordinal logistic biplots and nominal prediction regions; STATIS interstructure, compromise and trajectories; dynamic biplot; unfolding; DISTATIS; CDPCA; the ten GGE views of GGEBiplots (which-won-where, mean vs stability, discriminativeness vs representativeness, rankings...); AMMI1/AMMI2; CA/MCA maps, CCA and RDA triplots; PCoA/MDS with supplementary variables; external logistic biplot; canonical analysis of distances; co-inertia; FA, NIPALS and PLS biplots; polychoric biplot; constrained and external ordinal logistic biplots; PLSR1Bin and BinaryPLSR biplots; log-frequency biplot; bootstrap clouds of PCoA and MDS; HE plots of several factors. Their geometry is checked by tests (calibrated axes read the CVA group means exactly; ellipses and alpha-bags have their nominal coverage) and every figure is rendered in the test suite.

Differences from MultBiplotR (corrections)

See biplotpy.FIXES: within-groups standardization divides by the standard deviation (MultBiplotR divides by the variance); the initial logit is applied (MultBiplotR leaves the data untransformed); supplementary rows are projected after transforming them (MultBiplotR projects them untransformed); STATIS honours the requested transformation; Dual STATIS is complete; the binary logistic objective and gradient are consistent; NIPALS, Coinertia, PLSR, FA scores, PolyBiplot, BinaryDistances and weighted PCoA are corrected (16 entries in total).

Not ported

Already available in Python, so not duplicated: CA, MCA, MFA, FAMD (prince); CCA, RDA, PCoA, PERMDISP/betadisper (scikit-bio); PPML with high-dimensional fixed effects (pyfixest); Ledoit-Wolf shrinkage and LDA (scikit-learn); nearest PSD matrix (statsmodels.stats.correlation_tools.cov_nearest); orthogonal Procrustes and generalized eigenproblems (scipy).

Not included, and why:

A tutorial in Spanish with an example of each family of methods is in examples/tutorial_es.ipynb.