Ordination

The ordination module provides various ordination methods commonly used in community ecology.

NMDS (Non-metric Multidimensional Scaling)

nuee.ordination.metaMDS(X: ndarray | DataFrame, k: int = 2, distance: str = 'bray', trymax: int = 20, maxit: int = 1000, trace: bool = False, autotransform: bool = True, wascores: bool = True, expand: bool = True, random_state: int | None = None, **kwargs) → OrdinationResult

Non-metric Multidimensional Scaling with automatic transformation.

This function provides a high-level interface for NMDS ordination, following the conventions of the R vegan package’s metaMDS function. It automatically handles data transformation and uses multiple random starts to find the best ordination solution.

Parameters:

• X (np.ndarray or pd.DataFrame) – Community data matrix with samples in rows and species in columns. Values should be non-negative abundances or counts.

• k (int, default=2) – Number of dimensions for the ordination. Common choices are 2 or 3.

• distance (str, default="bray") – Distance metric to use for calculating dissimilarities. See nuee.vegdist() for available options.

• trymax (int, default=20) – Maximum number of random starts to find the best solution. Higher values increase computation time but may find better solutions.

• maxit (int, default=1000) – Maximum number of iterations for each random start.

• trace (bool, default=False) – If True, print progress information including stress values.

• autotransform (bool, default=True) – If True, automatically apply square root transformation to abundance data (values > 1) to reduce the influence of dominant species.

• wascores (bool, default=True) – If True, calculate weighted average species scores based on site scores and species abundances.

• expand (bool, default=True) – If True, expand the result to include additional information.

• random_state (int, optional) – Seed used for reproducible random starts. If None, each run is initialised independently.

• **kwargs (dict) – Additional parameters passed to the NMDS class.

Returns:

Result object containing:

• pointspd.DataFrame

  Site (sample) scores in the ordination space

• speciespd.DataFrame, optional

  Species scores (if wascores=True)

• stressfloat

  Final stress value (lower is better)

• convergedbool

  Whether the solution converged

Return type:

OrdinationResult

Notes

NMDS is a rank-based ordination method that attempts to preserve the rank order of dissimilarities between samples. Unlike metric methods like PCA, NMDS makes no assumptions about the linearity of relationships.

Stress values provide a measure of fit: - < 0.05: excellent - 0.05 - 0.10: good - 0.10 - 0.20: acceptable - > 0.20: poor (consider increasing k or using a different method)

The function uses multiple random starts (trymax) because NMDS can get stuck in local optima. The solution with the lowest stress is returned.

Examples

Basic NMDS ordination:

>>> import nuee
>>> species = nuee.datasets.varespec()
>>> nmds = nuee.metaMDS(species, k=2, distance="bray")
>>> print(f"Stress: {nmds.stress:.3f}")
Stress: 0.133

With custom parameters:

>>> nmds = nuee.metaMDS(species, k=3, distance="bray", trymax=50, trace=True)
Applying square root transformation
Applying Wisconsin double standardization
NMDS stress: 0.0869
NMDS converged

Visualize the ordination:

>>> import matplotlib.pyplot as plt
>>> fig = nuee.plot_ordination(nmds, display="sites")
>>> plt.show()

See also

NMDS

Lower-level NMDS class

rda

Redundancy Analysis (constrained ordination)

cca

Canonical Correspondence Analysis

pca

Principal Component Analysis

References

[1]

Kruskal, J.B. (1964). Nonmetric multidimensional scaling: a numerical method. Psychometrika 29, 115-129.

[2]

Minchin, P.R. (1987). An evaluation of relative robustness of techniques for ecological ordinations. Vegetatio 69, 89-107.

[3]

Oksanen, J., et al. (2020). vegan: Community Ecology Package. https://CRAN.R-project.org/package=vegan

class nuee.ordination.NMDS(n_components: int = 2, max_iter: int = 1000, n_init: int = 20, eps: float = 1e-12, random_state: int | None = None, dissimilarity: str = 'bray', n_jobs: int | None = None)

Bases: OrdinationMethod

Non-metric Multidimensional Scaling for community ecology.

This implementation follows the approach used in vegan’s metaMDS function, including multiple random starts and stress evaluation. NMDS is a rank-based ordination method that attempts to represent ecological distances in a reduced dimensional space.

Parameters:

• n_components (int, default=2) – Number of dimensions for the embedding (ordination axes).

• max_iter (int, default=1000) – Maximum number of iterations for the optimization algorithm.

• n_init (int, default=20) – Number of random initializations. The solution with minimum stress is returned.

• eps (float, default=1e-12) – Convergence tolerance for the stress value.

• random_state (int, optional) – Random seed for reproducibility. If None, the random state is not fixed.

• dissimilarity (str, default="bray") – Distance metric to use. See nuee.vegdist() for available options.

• n_jobs (int, optional) – Number of parallel jobs for computation. If None, uses a single core.

n_components

Number of dimensions.

Type:

int

max_iter

Maximum iterations.

Type:

int

n_init

Number of random starts.

Type:

int

eps

Convergence tolerance.

Type:

float

See also

metaMDS

High-level interface with automatic transformations

MDS

Scikit-learn’s MDS implementation

Notes

NMDS is particularly useful when: - Relationships between samples are non-linear - You want to use a specific distance metric - You have presence/absence or abundance data

The stress value indicates goodness-of-fit: - Values < 0.05 indicate excellent fit - Values 0.05-0.10 indicate good fit - Values > 0.20 indicate poor fit

References

[1]

Kruskal, J.B. (1964). Nonmetric multidimensional scaling. Psychometrika 29, 115-129.

Examples

>>> from nuee.ordination import NMDS
>>> import nuee
>>> species = nuee.datasets.varespec()
>>> nmds = NMDS(n_components=2, n_init=20)
>>> result = nmds.fit(species)
>>> print(f"Stress: {result.stress:.3f}")
Stress: 0.073

__init__(n_components: int = 2, max_iter: int = 1000, n_init: int = 20, eps: float = 1e-12, random_state: int | None = None, dissimilarity: str = 'bray', n_jobs: int | None = None)

Initialize NMDS with specified parameters.

fit(X: ndarray | DataFrame, distance: ndarray | None = None, **kwargs) → OrdinationResult

Fit NMDS to the data.

Parameters:

• X – Community data matrix (samples x species) or distance matrix

• distance – Precomputed distance matrix (optional)

• **kwargs – Additional parameters

Returns:

OrdinationResult with NMDS coordinates and stress

RDA (Redundancy Analysis)

nuee.ordination.rda(X: ndarray | DataFrame, Y: ndarray | DataFrame | None = None, Z: ndarray | DataFrame | None = None, formula: str | None = None, data: DataFrame | None = None, scale: bool = False, center: bool = True, **kwargs) → ConstrainedOrdinationResult

Redundancy Analysis (RDA).

RDA is a constrained ordination method that finds linear combinations of explanatory variables that best explain the variation in the response matrix.

Supplying Z fits a partial RDA: the effect of the conditioning variables is removed from both the response and the explanatory matrix before the constrained decomposition, and the inertia it accounts for is reported separately as partial_chi.

Parameters:

• X – Response matrix (samples x species)

• Y – Explanatory matrix (samples x variables)

• Z – Conditioning matrix for partial RDA (optional)

• formula – Formula string (e.g., “~ var1 + var2”)

• data – DataFrame containing variables for formula

• scale – Whether to scale species to unit variance

• center – Whether to centre species column-wise

• **kwargs – Additional parameters

Returns:

ConstrainedOrdinationResult with RDA results. Inertia is expressed as variance, and every component shares that scale, so partial_chi + sum(constrained_eig) + sum(unconstrained_eig) equals tot_chi (partial_chi being None for a non-partial fit). See inertia_partition() for the decomposition as a table.

Examples

# Simple RDA result = rda(species_data, environmental_data)

# RDA with formula result = rda(species_data, formula=”~ pH + temperature”, data=env_data)

# Partial RDA, conditioning on covariates result = rda(species_data, environmental_data, Z=conditioning_data) result.inertia_partition()

class nuee.ordination.RDA(formula: str | None = None, scale: bool = False, center: bool = True)

Bases: OrdinationMethod

Redundancy Analysis for constrained ordination.

This implementation follows the classical algorithm:

1. Centre/scale the response matrix as requested, then divide by sqrt(n - 1) so that the total sum of squares is a variance. This is the inertia scale vegan reports, and every eigenvalue inherits it.

2. Centre the constraining variables.

3. Solve the weighted least-squares regression analytically.

4. Perform SVD on fitted (constrained) and residual (unconstrained) matrices.

5. Store the full decomposition so downstream tools can reconstruct scaled scores mirroring vegan’s behaviour.

For a partial RDA the conditioning step happens between 1 and 3, and the already-prepared residual must not be rescaled a second time; see _partial_rda() and the prepared argument of _simple_rda().

__init__(formula: str | None = None, scale: bool = False, center: bool = True)

Parameters:

• formula – R-style model formula describing the constraints (optional).

• scale – Whether to standardise species to unit variance before analysis.

• center – Whether to centre species (response) data column-wise.

fit(X: ndarray | DataFrame, Y: ndarray | DataFrame | None = None, Z: ndarray | DataFrame | None = None, **kwargs) → ConstrainedOrdinationResult

Fit RDA to the data.

Parameters:

• X – Response matrix (samples x species).

• Y – Constraining variables (samples x predictors) or DataFrame used with a formula.

• Z – Conditioning matrix for partial RDA. When supplied, its effect is removed from both X and Y, and the inertia it accounts for is reported as partial_chi rather than folded into the constrained or unconstrained components.

Returns:

All inertia components share the scale of tot_chi; see ConstrainedOrdinationResult.inertia_partition().

Return type:

ConstrainedOrdinationResult

CCA (Canonical Correspondence Analysis)

nuee.ordination.cca(Y: ndarray | DataFrame, X: ndarray | DataFrame | None = None, formula: str | None = None, **kwargs) → ConstrainedOrdinationResult | OrdinationResult

Canonical Correspondence Analysis (or CA when X is None).

Parameters:

• Y – Species data matrix (sites x species).

• X – Environmental data matrix (sites x variables) or DataFrame for formula evaluation. When None and no formula is given, an unconstrained Correspondence Analysis (CA) is performed.

• formula – R-style formula string referencing columns in X.

Return type:

ConstrainedOrdinationResult (CCA) or OrdinationResult (CA).

class nuee.ordination.CCA(formula: str | None = None)

Bases: OrdinationMethod

Canonical Correspondence Analysis.

__init__(formula: str | None = None)

fit(Y: ndarray | DataFrame, X: ndarray | DataFrame | None = None, **kwargs) → ConstrainedOrdinationResult | OrdinationResult

Fit CCA (or CA when X is None) to the data.

PCA (Principal Component Analysis)

nuee.ordination.pca(X: ndarray | DataFrame, n_components: int | None = None, scale: bool = True, center: bool = True, **kwargs) → OrdinationResult

Principal Component Analysis.

Parameters:

• X – Data matrix (samples x variables).

• n_components – Number of components to keep.

• scale – Whether to scale variables (unit variance) prior to analysis.

• center – Whether to subtract column means prior to analysis.

Return type:

OrdinationResult with PCA results.

class nuee.ordination.PCA(n_components: int | None = None, scale: bool = True, center: bool = True)

Bases: OrdinationMethod

Principal Component Analysis for ordination.

__init__(n_components: int | None = None, scale: bool = True, center: bool = True)

fit(X: ndarray | DataFrame, **kwargs) → OrdinationResult

Fit PCA to the data using an explicit SVD.

Environmental Fitting

nuee.ordination.envfit(ordination: OrdinationResult, env: ndarray | DataFrame, permutations: int = 999, scaling: int | str | None = 2, choices: Sequence[int] | None = None, random_state: int | None = None, **kwargs) → Dict[str, Any]

Fit environmental vectors to an ordination configuration.

Parameters:

• ordination – OrdinationResult providing site scores.

• env – Environmental data matrix (matching site order).

• permutations – Number of permutations used for significance testing.

• scaling – Scaling option passed to ordination.get_scores.

• choices – Optional 1-based axis indices to include; defaults to all axes.

• random_state – Seed for the permutation generator.

Returns:

Dictionary mirroring vegan’s envfit output layout with a vectors entry containing scores, r, r², and p-values.

Return type:

dict

Procrustes Analysis

nuee.ordination.procrustes(X: ndarray | DataFrame, Y: ndarray | DataFrame, scale: bool = True) → dict

Procrustes analysis to compare two ordinations.

Parameters:

• X – First ordination matrix

• Y – Second ordination matrix

• scale – Whether to scale configurations

Returns:

Dictionary with Procrustes results

Variable Selection

nuee.ordination.ordistep(ordination: OrdinationResult, env: ndarray | DataFrame, direction: str = 'both', permutations: int = 199, p_enter: float = 0.05, p_remove: float = 0.1, max_steps: int = 50, random_state: int | None = None, **kwargs) → Dict[str, Any]

Stepwise model selection for constrained ordination.

Terms are added and/or dropped one at a time using permutation tests, in the manner of vegan’s ordistep. A candidate is tested as a constraint conditioned on the terms already retained, so each test is marginal.

Parameters:

• ordination – A fitted constrained ordination (e.g. from nuee.rda()), used as the source of the response matrix.

• env – Candidate explanatory variables, one per column.

• direction – "forward" only adds terms, "backward" starts from the full set and only drops them, "both" alternates.

• permutations – Permutations per test.

• p_enter – A term is added when its p-value is at or below this.

• p_remove – A retained term is dropped when its p-value rises above this.

• max_steps – Safety bound on the number of add/drop rounds.

• random_state – Seed, for reproducible selection.

Returns:

selected_variables (in the order chosen), p_values for those terms, the direction used, and the steps taken.

Return type:

dict

Raises:

ValueError – If direction is not one of the three accepted values, or the ordination carries no response matrix to select against.

Base Classes

class nuee.ordination.OrdinationResult(points: ndarray | None, species: ndarray | None = None, eigenvalues: ndarray | None = None, stress: float | None = None, converged: bool = True, call: Dict[str, Any] | None = None, site_u: ndarray | None = None, species_v: ndarray | None = None, singular_values: ndarray | None = None, site_u_unconstrained: ndarray | None = None, species_v_unconstrained: ndarray | None = None, singular_values_unconstrained: ndarray | None = None, row_weights: ndarray | None = None, column_weights: ndarray | None = None, scaling_backend: Callable[[int], Tuple[ndarray | None, ndarray | None]] | None = None)

Bases: object

Base class for ordination results.

points

Sample coordinates in ordination space

species

Species coordinates (if applicable)

eigenvalues

Eigenvalues for ordination axes

stress

Stress value (for NMDS)

converged

Whether the ordination converged

nobj

Number of objects (samples)

ndim

Number of dimensions

call

Dictionary containing method call parameters

__init__(points: ndarray | None, species: ndarray | None = None, eigenvalues: ndarray | None = None, stress: float | None = None, converged: bool = True, call: Dict[str, Any] | None = None, site_u: ndarray | None = None, species_v: ndarray | None = None, singular_values: ndarray | None = None, site_u_unconstrained: ndarray | None = None, species_v_unconstrained: ndarray | None = None, singular_values_unconstrained: ndarray | None = None, row_weights: ndarray | None = None, column_weights: ndarray | None = None, scaling_backend: Callable[[int], Tuple[ndarray | None, ndarray | None]] | None = None)

property proportion_explained: ndarray

Proportion of variance/inertia explained by each axis.

property cumulative_proportion: ndarray

Cumulative proportion of variance/inertia explained.

biplot(scaling: int | str | None = None, **kwargs)

Create a biplot.

For PCA the species loadings are shown as arrows from the origin. For constrained ordination (RDA / CCA), environmental variables are shown as arrows while species are shown as points.

Parameters:

• scaling – Scaling mode for site/species scores

• **kwargs – Additional plotting arguments

Returns:

matplotlib Figure object

get_scores(display: str = 'sites', scaling: int | str | None = None) → ndarray | Tuple[ndarray, ndarray]

Retrieve (optionally scaled) ordination scores.

Parameters:

• display – “sites”, “species”, or “both”

• scaling – None, 1/2/3 or “sites”/”species”/”symmetric”

Returns:

Requested scores as numpy arrays.

scores(display: str = 'sites', scaling: int | str | None = None) → ndarray | Tuple[ndarray, ndarray]

Retrieve (optionally scaled) ordination scores.

Parameters:

• display – “sites”, “species”, or “both”

• scaling – None, 1/2/3 or “sites”/”species”/”symmetric”

Returns:

Requested scores as numpy arrays.

class nuee.ordination.OrdinationMethod

Bases: ABC

Abstract base class for ordination methods.

abstractmethod fit(X: ndarray | DataFrame, **kwargs) → OrdinationResult

Fit the ordination method to the data.

Parameters:

• X – Community data matrix (samples x species)

• **kwargs – Method-specific parameters

Returns:

OrdinationResult object