Productivity

tooluniverse-epigenomics

mims-harvard/tooluniverse · updated Apr 8, 2026

$npx skills add https://github.com/mims-harvard/tooluniverse --skill tooluniverse-epigenomics
summary

Production-ready skill combining Python computation (pandas, scipy, numpy, pysam, statsmodels) with ToolUniverse annotation tools for epigenomics analysis.

skill.md

Genomics and Epigenomics Data Processing

Production-ready skill combining Python computation (pandas, scipy, numpy, pysam, statsmodels) with ToolUniverse annotation tools for epigenomics analysis.

LOOK UP, DON'T GUESS

When uncertain about any scientific fact, SEARCH databases first.

When to Use

Methylation data, ChIP-seq peaks, ATAC-seq, multi-omics integration, genome-wide epigenomic statistics. Keywords: methylation, CpG, ChIP-seq, ATAC-seq, histone, chromatin, epigenetic.

NOT for: RNA-seq DEG, variant calling, gene enrichment, protein structure.


Key Principles

  1. Data-first - Load/inspect before analysis
  2. Question-driven - Extract specific numeric answer
  3. Coordinate system awareness - Track genome build (hg19/hg38/mm10), chr prefix
  4. Statistical rigor - FDR correction, effect size filtering
  5. CpG identification - Parse Illumina probe IDs, genomic coordinates

Workflow

Phase 0: Question Parsing

Identify data files, specific statistic, thresholds, genome build. Categorize by keywords. See ANALYSIS_PROCEDURES.md for decision tree.

Phase 1: Methylation Processing

  • Load beta/M-value matrix (CSV/TSV/parquet/HDF5)
  • Filter by variance, missing rate, probe type, chromosome, CpG island relation
  • Differential methylation: T-test/Wilcoxon between groups + FDR
  • Age-related CpG: Pearson/Spearman correlation + FDR
  • Chromosome density: CpG count / chromosome length

Phase 2: ChIP-seq Peak Analysis

  • Load BED/narrowPeak/broadPeak, normalize chromosomes
  • Peak stats, annotation to genes, overlap analysis (Jaccard)

Phase 3: ATAC-seq

  • NFR detection (<150bp peaks), region classification

Phase 4: Multi-Omics Integration

  • Methylation-expression correlation per probe-gene (Pearson/Spearman + FDR)
  • ChIP-seq + expression: promoter peaks vs expression levels

Phase 5: Clinical Data

  • Missing data analysis across modalities, complete case identification

Phase 6: ToolUniverse Annotation

ENCODE tools:

  • ENCODE_search_rnaseq_experiments: assay_type ("total RNA-seq" default; fall back to "polyA plus RNA-seq"), biosample, limit
  • ENCODE_search_histone_experiments: target (e.g., "H3K27ac"), cell_type/tissue/biosample, limit

GEO tools: GEO_search_rnaseq_datasets, GEO_search_atacseq_datasets -- both accept limit or max_results

GTEx tools:

  • GTEx_get_median_gene_expression: gene_symbol (NOT Ensembl ID)
  • GTEx_query_eqtl: gene_symbol, tissue_id (case-sensitive exact, e.g., "Whole_Blood")

Other: ensembl_lookup_gene (requires species='homo_sapiens'), ensembl_get_regulatory_features (NO "chr" prefix), SCREEN_get_regulatory_elements, ChIPAtlas_* (requires operation param), SRA_search_experiments (library_strategy: "ChIP-Seq"/"Bisulfite-Seq"/"ATAC-seq")

Phase 7: Genome-Wide Statistics

Global mean/median beta, probe variance, chromosome density, DMP counts.

See CODE_REFERENCE.md for full implementations.


Common Patterns

Pattern Key Steps
Differential methylation Filter probes → groups → t-test → FDR → threshold
Age-related CpG density Correlate with age → FDR → map to chr → density ratio
Multi-omics missing data Extract IDs → intersect → check NaN → complete case count
ChIP-seq annotation Load peaks → annotate genes → classify regions
Methylation-expression Align samples → correlate → FDR → anti-correlations

GTEx Tissue IDs

Whole_Blood, Liver, Lung, Breast_Mammary_Tissue, Brain_Cortex, Heart_Left_Ventricle, Kidney_Cortex, Thyroid, Adipose_Subcutaneous, Muscle_Skeletal


Evidence Grading

Grade Criteria
Strong padj < 0.01 AND abs(delta-beta) >= 0.2, replicated
Moderate padj < 0.05 AND abs(delta-beta) >= 0.1
Weak padj < 0.05 but delta-beta < 0.1
Insufficient padj >= 0.05 or no replication

Delta-beta >= 0.2 = strong effect. ChIP-seq: q < 0.01, FE >= 2 for confidence. ATAC-seq NFR < 150bp = active regulatory. Always apply BH FDR. Verify genome build consistency.


Limitations

  • No pybedtools/pyBigWig: pure Python intervals
  • Illumina-centric (450K/EPIC); uses t-test/Wilcoxon (not limma)
  • No peak calling (assumes pre-called)
  • API rate limits: ~20 genes per batch

Reference Files

CODE_REFERENCE.md, TOOLS_REFERENCE.md, ANALYSIS_PROCEDURES.md, QUICK_START.md