cme-methylation 0.1.0

Chemical-master-equation model fitter for CpG methylation density data (Rust-accelerated)

CME Methylation Model

A chemical master equation (CME) kinetic model of CpG DNA methylation, fit to observed methylation distributions binned by local CpG density. The model captures cooperative methylation and demethylation dynamics and was used in the analyses reported in [paper title and citation here].

This code builds on prior work by Bonsu et al. (Ultrasens_DNAMethylation).

Repository layout

.
├── pyproject.toml                  # Packaging: builds the `cme-methylation` wheel
│                                   # (scripts as CLIs + bundled Rust backend)
├── cpg_density/                    # Rust CLI: genome FASTA → per-CpG density BED
│                                   # (preprocessing step 1; build with cargo)
├── python/
│   └── cme_methylation/            # Installed Python package
│       ├── fit_cme.py              # Main fitting / profile-likelihood / bootstrap script
│       │                           #   → installed as the `fit-cme` command
│       ├── build_density_distributions.py
│       │                           # Preprocess WGBS bedGraph data into
│       │                           # density-binned methylation histograms
│       │                           #   → installed as `build-density-distributions`
│       └── cpg_densities.py        # Pure-Python CpG-density helper
├── modeling/
│   └── cme_rust/                   # Rust backend (PyO3 + maturin), bundled into the wheel
│
├── notebooks/
│   ├── Plot Density and Methylation.ipynb
│   │                               # Plots of raw WGBS methylation summary
│   │                               # statistics by density bin
│   ├── Plot Model Output.ipynb
│   │                               # Plots fit-vs-data, StE/StW, CoE/CoW etc.
│   │                               # from fit_cme.py outputs
│   └── Plot Profile Likelihood Simple Model.ipynb
│                                   # Plots profile-likelihood curves
│
├── density_data/                   # Example .distribution.npz files (PBMC, ages 0/26/103)
├── requirements.txt
├── LICENSE                         # MIT
└── notebook_cell_tool.py           # Helper CLI for editing notebooks cell-by-cell

Installation

From PyPI

pip install cme-methylation

This installs the Rust backend (prebuilt wheel — no compiler needed) and two command-line tools, fit-cme and build-density-distributions.

From source (development)

The package is a single maturin project that compiles the Rust crate (modeling/cme_rust/) and installs the scripts as CLIs. From the repo root:

pip install maturin
maturin develop --release     # builds the Rust backend + installs fit-cme / build-density-distributions

The fitter falls back to a slower pure-Python path if the Rust backend is unavailable, but the prebuilt wheel / maturin develop always includes it.

CpG density tool (Rust CLI)

cpg_density/ is a standalone Rust CLI that computes per-CpG density from a genome FASTA — the first step of the from-scratch pipeline (see Usage §3). It is a separate, plain cargo build (no Python, no maturin):

cd cpg_density
cargo build --release
# binary: cpg_density/target/release/cpg_density

Usage

1. Fit the model on the example data

The included density_data/ directory contains three example .distribution.npz files (PBMC samples at ages 0, 26, and 103). To fit the collaborative model on these:

fit-cme \
    --input-dir ./density_data \
    --output-dir ./figures/example_fit \
    --name example_fit \
    --normalize StW \
    --shared-st --weight-by-cpg \
    --n-u 1.25 --n-m 1.25 \
    --bootstrap --n-jobs 8

Outputs land in ./figures/example_fit/:

• example_fit_params.csv — fitted kinetic parameters + fit statistics
• example_fit_{age}_profile.csv — profile-likelihood curves (with --profile)
• example_fit_bootstrap_* — bootstrap replicates (with --bootstrap)

2. Plot the results

Open the notebooks under notebooks/. The paths in Plot Model Output.ipynb and Plot Profile Likelihood Simple Model.ipynb are relative and assume you run Jupyter with the repo root as cwd.

Run one of the "CONFIG" cells in Plot Model Output.ipynb to choose a dataset. The default All CpGs only config points at modeling/figures/allcpg_stw_5meth_weight/ and the example density_data/ shipped with the repo. Other CONFIG cells reference dataset-specific paths (/path/to/your/...) that you need to edit to point at your own data.

Notes on Plot Density and Methylation.ipynb

This notebook reads raw WGBS bedGraph data (via polars_bio) and is included for reference. It is not runnable on the example .distribution.npz files — you need your own WGBS dataset and absolute paths pointing at it.

3. (Optional) Run the full pipeline on your own data

Starting from a genome FASTA and per-CpG WGBS methylation, two preprocessing steps produce the .distribution.npz that the fitter consumes.

3a. Compute CpG densities from the genome. cpg_density scans a genome FASTA and writes, for every CpG, the number of CpGs within ±window bp divided by window, as a 4-column BED (chrom start end density). The FASTA needs a .fai index (samtools faidx genome.fa):

samtools faidx genome.fa                 # once, if no .fai exists yet
cpg_density/target/release/cpg_density genome.fa --window 50 --threads 8
# -> genome_50_cpg_densities.bed

Works with any genome/assembly (hg19, hg38, mm10, …); CpG matching is case-insensitive, so soft-masked references are fine. Output is identical to the reference Python implementation but ~75× faster (full hg19 in ~8 s). Only --threads (default: all cores) and --window (default: 50) are tunable; --window sets the density definition and must match across a study.

3b. Build .distribution.npz. build_density_distributions.py combines the density BED from 3a with raw per-CpG methylation bedGraphs into the .distribution.npz format fit_cme.py expects:

build-density-distributions \
    --input-dir path/to/wgbs_bed_input \
    --density-file genome_50_cpg_densities.bed \
    --output-dir path/to/density_data_5meth \
    --n-meth-bins 5 \
    --n-bootstrap 100

Then fit as in §1, pointing --input-dir at your new output directory.

Input data format

.distribution.npz files contain:

key	shape	description
density_bins	(D,)	CpG density bin centers
meth_bins	(M+1,)	Methylation ratio bin edges
histograms	(D, M)	Count of CpGs per (density bin, methylation bin)
spacings	per-density-bin sequences	CpG spacings used in the cooperative-interaction kernel
total_cpg	(D,)	Total CpGs per density bin

Editing notebooks without opening them

Jupyter notebooks here are large. notebook_cell_tool.py is a small CLI for listing, viewing, replacing, inserting and deleting individual cells without loading the whole file:

python notebook_cell_tool.py list notebooks/<file>.ipynb
python notebook_cell_tool.py show notebooks/<file>.ipynb <cell-id-prefix>
python notebook_cell_tool.py replace notebooks/<file>.ipynb <cell-id> --source-file new.py

License

MIT — see LICENSE.