Methylome Segmentation Workflow

This implementation-level companion to MethylSeg methodology maps the public workflow methods to the components that perform the work. For a runnable end-to-end example, see the full pipeline notebook.

Workflow layers

  1. Prepare an input table as SampleInfo with methylseg.MethylSegPathway.prepare_sample_info().

  2. Fit state-assignment components with methylseg.MethylSegPathway.fit_pathway().

  3. Segment one chromosome with methylseg.MethylSegPathway.generate_regions() or many with methylseg.MethylSegPathway.run_on_all_chroms().

  4. Clean raw regions with methylseg.MethylSegPathway.get_clean_regions().

  5. Use methylseg.MethylSegPathway.run_pathway() to perform steps 2 through 4 in one call.

Sample preparation

prepare_sample_info() delegates input parsing and filtering to methylseg.MethylDataPrep. It returns the canonical sample_info object plus removed_df, a table of CpGs removed during preprocessing.

from pathlib import Path
from methylseg import MethylSegPathway

sample_info, removed_df = MethylSegPathway.prepare_sample_info(
    sample_name="sample-a",
    sample_file=Path("sample-a.tsv.gz"),
    resolution="wgbs",  # or "450k"
    remove_low_coverage_like_cpgs=True,
)

Keep removed_df if it should appear as a background layer in label plots.

Fitting state assignment

Create a pathway with a training sample, then call fit_pathway(). This trains the assigner through methylseg.MethylStateAssigner.train_kmeans_for_sample(). For a rule-based configuration, it also optimizes analyzer cutoffs through methylseg.MethylStateAnalyzer.optimize_rule_params_random().

pathway = MethylSegPathway(
    train_sample_info=sample_info,
    out_dir=Path("out") / sample_info.sample_id,
)
pathway.fit_pathway()

Use train_chroms and max_cpg_per_chrom in the constructor to limit fitting data. See the KMeans notebook and the rule-based notebook for model-specific configuration.

Segmentation and cleaning

generate_regions() segments one chromosome, creates contiguous regions, writes raw BED files, and applies the configured minimum-length filter to the returned table.

chr1_regions = pathway.generate_regions(
    sample_info=sample_info,
    chrom="chr1",
    min_probes=3,
)

run_on_all_chroms() is the genome-scale counterpart. It accepts a chromosome subset, writes raw summaries, and by default cleans the regions and writes cleaned summaries.

summary_paths = pathway.run_on_all_chroms(
    sample_info=sample_info,
    chroms=["chr1", "chr2"],
    clean_regions=True,
)

For the Sticky HMM, selected chromosomes are segmented jointly. Other backends call generate_regions once per chromosome. Both paths use methylseg.MethylSegmentor.segment_sample() and methylseg.MethylSegmentor.create_regions().

get_clean_regions() merges and filters raw regions using pathway settings. It writes chromosome-local artifacts in clean_regions/. When called genome-wide, it can also write cleaned summaries; its metadata TSV files are used for cleaned-region plotting.

state_paths, metadata_path = pathway.get_clean_regions(
    regions_df=chr1_regions,
    sample_id=sample_info.sample_id,
    chrom="chr1",
)

What run_pathway() calls

run_pathway() is the convenience entry point. The table shows its orchestration path; summary writing is an implementation detail included to explain the resulting artifacts.

End-to-end call chain

Entry point

Delegates to

Result

run_pathway()

fit_pathway()

Trains KMeans and optionally optimizes rule cutoffs.

run_pathway()

run_on_all_chroms()

Resolves selected chromosomes and coordinates segmentation.

run_on_all_chroms()

generate_regions() per chromosome, or joint Sticky-HMM processing

Produces raw region tables and chromosome-level BED files.

Region generation

segmentor.segment_sample() then segmentor.create_regions()

Smooths CpG labels and groups contiguous states into regions.

run_on_all_chroms() with clean_regions=True

get_clean_regions()

Writes cleaned BED and metadata artifacts.

run_on_all_chroms()

Summary-file writing

Writes genome-wide raw summaries and, when requested, cleaned summaries.

Alternative usage patterns

One-call workflow

summary_paths = pathway.run_pathway(
    sample_info=sample_info,
    chroms=["chr1", "chr2"],
    clean_regions=True,
)

Explicit fitting and per-chromosome segmentation

pathway.fit_pathway()
regions_by_chrom = {
    chrom: pathway.generate_regions(sample_info=sample_info, chrom=chrom)
    for chrom in ["chr1", "chr2"]
}

Manual postprocessing

raw_regions = pathway.generate_regions(sample_info=sample_info, chrom="chr1")
pathway.get_clean_regions(
    regions_df=raw_regions,
    sample_id=sample_info.sample_id,
    chrom="chr1",
    merge_gap_bp=50_000,
    min_region_length=10_000,
)

Component-oriented usage

The pathway exposes assigner, analyzer, and segmentor for custom analysis after fitting.

pathway.fit_pathway()
meth_data, emissions = pathway.assigner.prepare_sample_for_clustering(
    sample_info=sample_info, chrom="chr1",
)
kmeans_result = pathway.assigner.apply_kmeans_to_sample(
    sample_info=sample_info, chrom="chr1",
)
rule_labels = pathway.analyzer.define_states_by_rules(
    sample_info=sample_info, chrom="chr1",
)
segmented_data, hmm_model = pathway.segmentor.segment_sample(
    sample_info=sample_info, chrom="chr1",
)
regions = pathway.segmentor.create_regions(region_min_probes=3)

For a complete component walkthrough, see the components example. Component calls have more state and ordering requirements than pathway methods.