Big Data Trajectories#

Processing large trajectory datasets with constant memory using H5MSM and iterators.

High-throughput molecular dynamics trajectories easily reach tens to hundreds of gigabytes. Loading entire trajectories into memory causes high memory pressure and crashes analysis scripts. MolSysMT addresses this challenge through the H5MSM format, an HDF5-based hierarchical storage architecture designed for out-of-core streaming and sub-millisecond random frame access.

In this recipe, we demonstrate how to inspect, selectively slice, and stream trajectories with constant \(O(1)\) memory usage using chunked iterators.

Added in version 1.0.0.

Inspecting Metadata#

MolSysMT queries metadata (atom counts, frame counts, box parameters, and topology) directly from the HDF5 header without loading coordinate datasets into RAM:

import molsysmt as msm

# Reference bundled H5MSM system
h5_file = msm.systems['chicken villin HP35']['chicken_villin_HP35_solvated.h5msm']

# Inspect file metadata at 0 MB memory cost
msm.info(h5_file)
n_atoms, n_structures = msm.get(h5_file, n_atoms=True, n_structures=True)
print(f"File contains {n_structures} structures with {n_atoms} atoms.")

File contains 1 structures with 4369 atoms.

Selective Slicing#

You can load only the necessary atoms (e.g. protein alpha-carbons) across specific frame ranges, skipping solvent and unneeded frames entirely at the disk read level:

# Load topology and trajectory bundle
top = msm.systems['chicken villin HP35']['chicken_villin_HP35_solvated.h5msm']
traj = msm.systems['chicken villin HP35']['traj_chicken_villin_HP35_solvated.dcd']
molsys = msm.convert([top, traj], to_form='molsysmt.MolSys')

# Read only C-alpha atoms for every second structure
ca_coords = msm.get(
    molsys,
    selection='atom_name=="CA"',
    structure_indices=range(0, n_structures, 2),
    coordinates=True
)

print(f"Loaded subset coordinates shape: {ca_coords.shape}")

WARNING: Structural attributes were dropped because only an item outside the structure axis of the molecular system provides them: atom_index, structure_chemical_state_index, structure_id, time. (Hint: A file holding a single reference conformation cannot supply a series for a whole trajectory. Take the attribute from the trajectory item, or convert first. Docs: https://www.uibcdf.org/MolSysMT)

Loaded subset coordinates shape: (1, 36, 3)

Streaming Iterators#

When processing massive multi-gigabyte trajectories, molsysmt.basic.Iterator enables chunked streaming, processing fixed-size batches in constant memory:

# Create a chunked iterator streaming 10 frames per chunk
iterator = msm.Iterator(molsys, selection='atom_name=="CA"', chunk=10, coordinates=True)

chunks_read = []
for chunk in iterator:
    chunks_read.append(chunk)

print(f"Processed {len(chunks_read)} chunks in streaming mode, chunk shape: {chunks_read[0].shape}.")

Processed 2 chunks in streaming mode, chunk shape: (10, 36, 3).

Writing H5MSM Files#

MolSysMT allows converting heterogeneous trajectory files into a single unified, ultra-fast H5MSM archive:

import os

# Convert subset to a new standalone H5MSM file
output_h5 = 'subset_traj.h5msm'
msm.convert(molsys, to_form='file:h5msm', output_filename=output_h5, selection='atom_name=="CA"')

print(f"Compact H5MSM created: {os.path.exists(output_h5)} (size: {os.path.getsize(output_h5)} bytes)")
if os.path.exists(output_h5):
    os.remove(output_h5)

Compact H5MSM created: True (size: 1060097 bytes)