skbio
PCoA
Principal Coordinate Analysis
Python
data analysis

How to get skbio PCoA Principal Coordinate Analysis results?

ML System Design practice on Codemia

Design recommenders, ranking systems and training pipelines the way ML interviews actually ask for them, with worked solutions.

Practice ML system design

Introduction

Principal Coordinate Analysis (PCoA), also called Classical Multidimensional Scaling, is a dimensionality reduction technique that takes a distance matrix and produces a set of coordinates in a lower-dimensional space that best preserves the original distances. It is widely used in ecology and bioinformatics to visualize differences between microbial communities, genetic samples, or ecological sites. The scikit-bio (skbio) library provides skbio.stats.ordination.pcoa() to perform PCoA in Python.

Prerequisites

Install scikit-bio and its dependencies:

bash
pip install scikit-bio

scikit-bio requires NumPy, SciPy, and pandas. It works on Python 3.8+.

Creating a Distance Matrix

PCoA requires a skbio.DistanceMatrix as input. You can create one from a 2D array or compute it from sample data:

python
1import numpy as np
2from skbio import DistanceMatrix
3
4# From a symmetric distance array
5data = np.array([
6    [0.0, 0.5, 0.9, 0.8],
7    [0.5, 0.0, 0.6, 0.7],
8    [0.9, 0.6, 0.0, 0.3],
9    [0.8, 0.7, 0.3, 0.0]
10])
11
12sample_ids = ['SampleA', 'SampleB', 'SampleC', 'SampleD']
13dm = DistanceMatrix(data, ids=sample_ids)
14print(dm)

From ecological abundance data using a beta diversity metric:

python
1from skbio.diversity import beta_diversity
2
3# Rows = samples, columns = species counts
4counts = np.array([
5    [10, 20, 30, 0],
6    [15, 25, 5, 10],
7    [0, 5, 40, 20],
8    [30, 0, 10, 15]
9])
10
11sample_ids = ['Site1', 'Site2', 'Site3', 'Site4']
12
13# Bray-Curtis dissimilarity (common in ecology)
14dm = beta_diversity('braycurtis', counts, ids=sample_ids)

Running PCoA

python
1from skbio.stats.ordination import pcoa
2
3results = pcoa(dm)
4print(type(results))  # <class 'skbio.stats.ordination.OrdinationResults'>

Accessing PCoA Results

The OrdinationResults object contains several important attributes:

Coordinates (Sample Scores)

python
1# Principal coordinates for each sample
2coords = results.samples
3print(coords)
4#                  PC1       PC2       PC3
5# SampleA -0.382...  0.123...  0.045...
6# SampleB -0.115...  -0.298... -0.012...
7# SampleC  0.341...  0.087...  -0.089...
8# SampleD  0.156...  0.088...  0.056...
9
10# Access as NumPy array
11coords_array = results.samples.values
12print(coords_array.shape)  # (4, 3) — 4 samples, 3 PCs

Proportion Explained (Eigenvalues)

python
1# Proportion of variance explained by each axis
2print(results.proportion_explained)
3# PC1    0.652...
4# PC2    0.231...
5# PC3    0.117...
6
7# Total variance explained by first 2 axes
8total_2d = results.proportion_explained.iloc[:2].sum()
9print(f"First 2 axes explain {total_2d:.1%} of variance")
10
11# Raw eigenvalues
12print(results.eigvals)
13# PC1    0.198...
14# PC2    0.070...
15# PC3    0.036...

Plotting PCoA Results

python
1import matplotlib.pyplot as plt
2
3coords = results.samples
4
5fig, ax = plt.subplots(figsize=(8, 6))
6ax.scatter(coords['PC1'], coords['PC2'], s=100, edgecolors='black')
7
8# Label each point
9for idx, row in coords.iterrows():
10    ax.annotate(idx, (row['PC1'], row['PC2']),
11                textcoords='offset points', xytext=(5, 5))
12
13pct1 = results.proportion_explained['PC1'] * 100
14pct2 = results.proportion_explained['PC2'] * 100
15ax.set_xlabel(f'PC1 ({pct1:.1f}%)')
16ax.set_ylabel(f'PC2 ({pct2:.1f}%)')
17ax.set_title('PCoA of Sample Distances')
18ax.axhline(0, color='gray', linewidth=0.5)
19ax.axvline(0, color='gray', linewidth=0.5)
20plt.tight_layout()
21plt.savefig('pcoa_plot.png', dpi=150)
22plt.show()

Color by Group Metadata

python
1import pandas as pd
2
3metadata = pd.DataFrame({
4    'group': ['Control', 'Control', 'Treatment', 'Treatment']
5}, index=sample_ids)
6
7colors = {'Control': 'blue', 'Treatment': 'red'}
8
9fig, ax = plt.subplots(figsize=(8, 6))
10for group, color in colors.items():
11    mask = metadata['group'] == group
12    samples = coords.loc[mask]
13    ax.scatter(samples['PC1'], samples['PC2'],
14               c=color, label=group, s=100, edgecolors='black')
15
16ax.legend()
17ax.set_xlabel(f'PC1 ({pct1:.1f}%)')
18ax.set_ylabel(f'PC2 ({pct2:.1f}%)')
19plt.tight_layout()
20plt.show()

Using the Built-in Plot Method

OrdinationResults has a built-in plot() method for quick visualization:

python
1# Simple 2D plot
2fig = results.plot(
3    df=metadata,              # DataFrame with sample metadata
4    column='group',           # Column to color by
5    title='PCoA',
6    cmap='Set1',
7    s=80
8)

For 3D plots (requires mpl_toolkits):

python
1fig = results.plot(
2    df=metadata,
3    column='group',
4    axes=(0, 1, 2),  # Plot PC1 vs PC2 vs PC3
5    title='PCoA 3D'
6)

Statistical Testing with PERMANOVA

After PCoA visualization, test whether groups are statistically different:

python
1from skbio.stats.distance import permanova
2
3result = permanova(dm, metadata['group'], permutations=999)
4print(result)
5# method name    PERMANOVA
6# test statistic name    pseudo-F
7# sample size    4
8# number of groups    2
9# test statistic    ...
10# p-value    ...

Common Pitfalls

  • Non-symmetric or negative diagonal: DistanceMatrix requires a symmetric matrix with zeros on the diagonal. Even floating-point rounding can cause DistanceMatrix(data) to fail. Use (data + data.T) / 2 to force symmetry and np.fill_diagonal(data, 0) to fix the diagonal.
  • Negative eigenvalues: Some distance metrics (like UniFrac or custom metrics) can produce distance matrices that are not Euclidean, resulting in negative eigenvalues. PCoA will still run, but the proportion explained may not sum to 1. Use pcoa(dm, method='fsvd', number_of_dimensions=3) for faster computation that avoids this issue.
  • Too few samples: PCoA with fewer than 4 samples produces at most 3 axes and may not give meaningful results. Ensure adequate sample size for your study.
  • Interpreting distances: PCoA preserves distances, not correlations. Points close together on the plot are similar according to your distance metric. The choice of metric (Bray-Curtis, UniFrac, Euclidean) strongly affects interpretation.
  • Missing sample IDs: If your distance matrix IDs do not match your metadata index, coloring by group will fail silently or raise an error. Always verify alignment with dm.ids and metadata.index.

Summary

  • Use skbio.stats.ordination.pcoa(distance_matrix) to run PCoA
  • Access coordinates with results.samples (DataFrame with PC1, PC2, ..., columns)
  • Check variance explained with results.proportion_explained — report this in publications
  • Create distance matrices with skbio.diversity.beta_diversity() or from a NumPy array via DistanceMatrix()
  • Plot with results.samples and matplotlib, or use the built-in results.plot() method
  • Follow up with PERMANOVA (skbio.stats.distance.permanova) to test group significance

Related reading
Free course
Beginner
7 lessons
2 hours
Tackling System Design Interview Problems

A short course that equips you with the skills to approach system design interviews methodically.

Start the free course
Track what you have practised

A free account saves your progress, solutions and study plan across every problem on Codemia.

ML System Design practice on Codemia

Design recommenders, ranking systems and training pipelines the way ML interviews actually ask for them, with worked solutions.

Practice ML system design

All Rights Reserved.