Regulatory Genomics Research Skill

Systematic investigation of gene regulation through transcription factor binding, chromatin state, and regulatory element annotation. Integrates JASPAR (TF motifs), ENCODE (functional genomics experiments), RegulomeDB (regulatory variant scoring), and UCSC cCREs.

Domain Reasoning

Regulatory element identification requires converging lines of evidence: sequence conservation alone is insufficient (many conserved sequences are not regulatory), chromatin accessibility is necessary but not sufficient (open chromatin can be structural), TF binding peaks require motif validation, and eQTL evidence ties the element to a transcriptional outcome. No single data type is sufficient. A high-confidence regulatory element requires at least two independent evidence types, and ideally all four.

LOOK UP DON'T GUESS

• TF binding motifs: retrieve from jaspar_search_matrices and jaspar_get_matrix; do not describe motifs from memory.
• Experimental ChIP-seq data: search ENCODE_search_experiments; do not assume a TF has been profiled in a given cell type.
• cCRE annotations for a genomic region: call UCSC_get_encode_cCREs with exact coordinates; do not guess element types.
• Regulatory impact of a variant: query RegulomeDB_query_variant; never estimate regulatory importance from position alone.

---

KEY PRINCIPLES:

1. English-first queries - Use English gene/TF names in all tool calls; respond in user's language
2. Evidence layering - Combine motif (JASPAR) + experimental (ENCODE ChIP-seq) + variant (RegulomeDB) evidence
3. Coordinate precision - Genome coordinates must specify assembly (GRCh38 preferred)
4. Negative results documented - Report when a TF has no ChIP-seq data in ENCODE

---

COMPUTE, DON'T DESCRIBE

When analysis requires computation (statistics, data processing, scoring, enrichment), write and run Python code via Bash. Don't describe what you would do — execute it and report actual results. Use ToolUniverse tools to retrieve data, then Python (pandas, scipy, statsmodels, matplotlib) to analyze it.

When to Use

• "What transcription factors bind near gene X?"
• "Does this SNP affect a regulatory element?"
• "Find CTCF binding sites in liver tissue"
• "What are the enhancers active in this cell type?"
• "Show me ChIP-seq experiments for H3K27ac in T cells"
• "Is rs1234567 in a regulatory region?"
• "What TF motifs overlap this genomic region?"
• "Find ENCODE experiments for ATAC-seq in cancer cell lines"
• "Predict the effect of this non-coding variant on expression / chromatin accessibility"
• "Predict regulatory activity (expression, accessibility, TF binding) directly from a DNA sequence"
• "Which of these enhancer variants is predicted to be most disruptive?"

---

Key Tools

Tool	Purpose	Key Params
jaspar_search_matrices	Find TF binding motifs by TF name or organism	name, species, collection, tax_id
jaspar_get_matrix	Get full PWM/PFM for a specific JASPAR matrix	matrix_id (e.g., "MA0139.1")
JASPAR_get_transcription_factors	List all TF matrices (paginated)	page, page_size
UniBind_search_datasets	Find curated, experimentally-derived direct TF-DNA binding-site datasets	tf_name, species, cell_line, collection (Robust/Permissive)
UniBind_get_dataset	Get binding-site detail for one UniBind dataset (JASPAR motifs, score/CentriMo thresholds, BED/FASTA URLs)	dataset_id
UniBind_list_tfs	List/filter TFs profiled in UniBind	search (substring), limit
ENCODE_search_experiments	Search ENCODE ChIP-seq/ATAC-seq/WGBS experiments	assay_title, target, biosample_term_name, limit
ENCODE_search_histone_experiments	Search histone mark ChIP-seq specifically	histone_mark, biosample_term_name, limit
ENCODE_search_chromatin_accessibility	Search ATAC-seq/DNase-seq experiments	biosample_term_name, limit
ENCODE_get_experiment	Get full metadata for a specific ENCODE experiment	accession (e.g., "ENCSR000EGM")
ENCODE_search_annotations	Search ENCODE cCRE and chromatin state annotations	annotation_type, biosample_term_name, limit
ENCODE_get_chromatin_state	Search ChromHMM segmentation data	biosample_term_name, limit
UCSC_get_encode_cCREs	Get cCREs overlapping a genomic region	chrom, start, end
RegulomeDB_query_variant	Score regulatory impact of a variant	rsid (e.g., "rs4994")
ENCODE_search_biosamples	Find available cell lines/tissues in ENCODE	term_name, biosample_type, limit
EnsemblReg_get_motif_features	TF binding motif instances in a region (Ensembl Regulatory Build)	region, species
EnsemblReg_get_constrained_elements	Evolutionarily constrained (purifying-selection) elements in a region — the conservation evidence type called out below	region, species
EnsemblReg_get_binding_matrix	Full PWM/PFM for an Ensembl binding-matrix stable ID	binding_matrix_id (e.g. ENSPFM0320)
MEME_fimo_scan / MEME_discover_motifs / MEME_tomtom_compare / MEME_list_databases	Scan sequences for known motifs / discover novel motifs de novo / compare a motif against JASPAR-HOCOMOCO-CIS-BP / list available comparison databases	sequences, motif_text (MEME format)
HOCOMOCO_search_motifs / HOCOMOCO_get_motif	Curated ChIP-seq-derived TF motif by gene name / full motif detail + optional PWM	query (gene name) / motif_id
Dfam_search_families / Dfam_get_family / Dfam_get_annotations	Find transposable-element/repeat families / get one family's detail / get repeat annotations for a genomic region	name_prefix/clade/repeat_type / accession / chrom,start,end

Sequence-based deep-learning models (predict, don't just annotate)

The tools above tell you what is known to be at a locus (databases). These models instead predict regulatory activity directly from the DNA sequence, and — by scoring a reference vs. alternate window — predict what a non-coding variant does. RegulomeDB ranks a variant by overlap with existing annotations; these give a quantitative, tissue-aware effect size even for novel variants with no annotation. Reach for them when annotation is silent or when the question is "how much does this allele change regulation".

Tool	Op	Predicts	Context	Access
AlphaGenome_atlas_lookup_variant / AlphaGenome_atlas_scan_interval	precomputed lookup: one SNV / every SNV in a region	unified AVI_SCORE (AlphaGenome + AlphaMissense) + per-track scores, no live model run	up to 10 kb per scan	hosted API — ALPHA_GENOME_API_KEY (much higher query rate than the live ops below)
AlphaGenome_predict_interval / AlphaGenome_score_variant / AlphaGenome_score_ism_variants	profile region / score variant / saturation-mutagenesis scan	RNA-seq, ATAC, CAGE, splice tracks (frontier accuracy, single-base); ISM ranks every substitution in a ≤500 bp window by effect — good for "which base in this element matters"	up to 1 Mb (ISM window ≤500 bp)	hosted API — ALPHA_GENOME_API_KEY
run_enformer_predict / run_enformer_variant_effect	profile / score	5,313 human (+1,643 mouse) tracks: expression, chromatin, TF binding	196 kb	remote MCP server
run_borzoi_predict / run_borzoi_variant_effect	profile / score	RNA-seq coverage (expression / polyA / splicing emphasis), 7,611 tracks	524 kb	remote MCP server
run_chrombpnet_predict / run_chrombpnet_variant_effect	profile / score	chromatin accessibility (ATAC / DNase), base-resolution profile + counts	~2 kb	remote MCP server
Evo2_score_variant	score	genome-foundation-model delta log-likelihood; coding and non-coding	up to 1 Mb	hosted NIM — NVIDIA_API_KEY
gi_* (opt-in, GENOMIC_INTELLIGENCE_MCP_URL)	profile (6 tasks) + composite find-genes-then-predict	promoter, splice sites, enhancer activity, chromatin state, expression, gene annotation	—	public demo server, shared rate-limited quota, no key

Picking one: for a known single-nucleotide variant, check AlphaGenome_atlas_lookup_variant (or atlas_scan_interval for a whole region) first — it's a precomputed database read covering all ~9B possible human SNVs, far cheaper than a live call, and its AVI_SCORE already fuses AlphaGenome's regulatory signal with AlphaMissense's coding-impact model into one number. Fall back to the live AlphaGenome_* operations for indels, custom/synthetic sequences, or full per-track detail (predict_variant/predict_sequence/score_interval are also available for raw-track and no-variant use cases); run_enformer_* / run_borzoi_* are the published, self-hostable equivalents (Enformer for general regulation, Borzoi when expression/splicing is the question); run_chrombpnet_* when the question is specifically chromatin accessibility; Evo2_score_variant as a sequence-only check that also covers coding variants. Outputs are Δ (alt − ref) effect sizes, not calibrated probabilities (except AVI_SCORE, which is a unified impact score) — rank/calibrate against known variants. If no key/server is provisioned, fall back to the annotation tools above and say so.

---

Workflow

Phase 1: TF Motif Discovery (JASPAR)

When asked about TF binding motifs or what TFs might regulate a gene:

1. jaspar_search_matrices(name="TF_NAME", species="Homo sapiens")
   -> Returns list of matrices with matrix_id, collection, version

2. jaspar_get_matrix(matrix_id="MA0139.1")
   -> Returns full PFM/PWM matrix, sequence logo URL, binding sites URL

3. For broad TF family search:
   jaspar_search_matrices(species="Homo sapiens", collection="CORE")
   -> Filter by TF family name in results

JASPAR Collections:

• CORE: High-quality, non-redundant matrices (best for most use cases)
• CNE: Conserved non-coding elements
• POLII: RNA Pol II binding sites

Key Response Fields:

• matrix_id: Versioned ID (e.g., "MA0139.1") — use for jaspar_get_matrix
• name: TF gene symbol
• sequence_logo: URL to binding site logo PNG/SVG
• collection: Which JASPAR collection

Phase 1b: Direct TF-DNA Binding Sites (UniBind)

JASPAR gives the motif model; UniBind gives curated, experimentally-derived direct binding sites (motif-anchored, base-pair resolution, per ChIP-seq experiment via the DAMO/ChIP-eat pipeline). Use it to find which experiments support binding for a TF and to get downloadable BED/FASTA of the sites — a middle layer between JASPAR motifs and raw ChIP-seq peaks.

# Find curated binding-site datasets for a TF (filters compose)
tu.run_tool("UniBind_search_datasets",
            {"tf_name": "CTCF", "species": "Homo sapiens",
             "collection": "Robust", "page_size": 5})
#   -> list of {tf_name, total_peaks, dataset_id, dataset_url}

# Full binding-site detail for one dataset
tu.run_tool("UniBind_get_dataset",
            {"dataset_id": "EXP030726.neural_stem_cells.SMAD3"})
#   -> tf_name, cell_line, jaspar_id[], tfbs_models[] each with
#      jaspar_id, total_tfbs, score_threshold, adj_centrimo_pvalue,
#      bed_url, fasta_url, summary_plot_url

# Discover valid tf_name values (search is client-side substring)
tu.run_tool("UniBind_list_tfs", {"search": "SMAD"})   # -> [SMAD2, SMAD3, SMAD4]

Notes: species is the scientific name ('Homo sapiens', not a taxid); collection is 'Robust' (high-confidence) or 'Permissive'; public, no API key.

Phase 1c: Ensembl Regulatory Build (motifs + the missing conservation evidence type)

The Domain Reasoning above lists four converging evidence types for a high-confidence regulatory element — conservation, accessibility, TF binding, eQTL — but conservation had no tool until now. EnsemblReg_get_constrained_elements fills that gap directly: it returns evolutionarily constrained (purifying-selection) elements in a region, independent of JASPAR/UniBind/ENCODE's binding-evidence tools. Combine it with EnsemblReg_get_motif_features (a second, database-scale source of TF motif instances, complementary to JASPAR/UniBind) to build the two-or-more-evidence-types case this skill's reasoning already calls for.

# TF binding motif instances in a region (species defaults to homo_sapiens)
tu.run_tool("EnsemblReg_get_motif_features", {"region": "7:140424943-140524564"})
#   -> motif_count, motif_features[] each with stable_id,
#      transcription_factor_complex, binding_matrix_stable_id, score,
#      start/end/strand

# Evolutionary conservation evidence for the same region
tu.run_tool("EnsemblReg_get_constrained_elements", {"region": "17:7661779-7687538"})
#   -> element_count, constrained_elements[] each with score (higher =
#      stronger conservation), start/end

# Full PWM for a binding matrix found above (e.g. from motif_features)
tu.run_tool("EnsemblReg_get_binding_matrix", {"binding_matrix_id": "ENSPFM0320"})
#   -> associated_tfs[], nucleotide-frequency matrix, threshold, consensus

Operational notes (live-verified, not assumed): these three endpoints are slow — get_motif_features typically takes 10-25s, get_binding_matrix 50-60s, get_constrained_elements up to 25s for a small region. Use a generous timeout (60-90s) rather than treating a slow response as a hang. get_constrained_elements is also unreliable on larger regions (e.g. the ~100kb BRAF locus reproducibly failed or timed out on repeated attempts while the smaller TP53 region succeeded) — prefer regions well under 50kb for this specific endpoint, and retry once before concluding the call failed.

Phase 1d: TF Motif Discovery and Comparison (MEME Suite) + Curated Reference Motifs (HOCOMOCO)

JASPAR gives a curated motif library; the MEME Suite tools instead let you discover a motif de novo from a set of related sequences and compare any motif (discovered or known) against reference databases. HOCOMOCO is a second curated reference library (ChIP-seq-derived, human/mouse only) — complementary to JASPAR, not a replacement, since the two are built from different underlying peak sets and sometimes disagree on model quality for the same TF.

# 1. Scan sequences for known motif occurrences (submits a job to
#    meme-suite.org, polls ~10-30s)
tu.run_tool("MEME_fimo_scan", {
    "sequences": ">promoter\nGGATCCGCGCGCTATAAAAGGATCC...",
    "motif_text": "MEME version 5\n\nALPHABET= ACGT\n...",  # full MEME-format motif
    "pvalue_threshold": 0.0001
})
#   -> hits[] with start/stop/strand/score/pvalue/qvalue/matched_sequence

# 2. De novo motif discovery across >=2 related sequences (~30-120s)
tu.run_tool("MEME_discover_motifs", {
    "sequences": ">seq1\n...\n>seq2\n...",  # co-regulated gene promoters, etc.
    "nmotifs": 3
})
#   -> motifs[] with consensus, width, sites, evalue, probability_matrix

# 3. Identify which known TF a discovered (or any) motif resembles
#    (~5-20s)
tu.run_tool("MEME_tomtom_compare", {
    "query_motif": "MEME version 5\n...",
    "target_db": "JASPAR2026_vertebrates"  # or "HOCOMOCO_v12"; see list_databases
})
#   -> matches[] ranked by evalue, with target_id (e.g. JASPAR MA0106.3)

# 4. See available target databases for tomtom_compare (no remote call,
#    cached metadata)
tu.run_tool("MEME_list_databases", {"category_filter": "JASPAR"})

HOCOMOCO as a curated cross-check — given a TF name, get its curated consensus/PWM/quality grade directly, no discovery step needed:

tu.run_tool("HOCOMOCO_search_motifs", {"query": "CTCF"})
#   -> [{motif_id: "CTCF.H14CORE.0.P.B", gene_name_human, quality, consensus, ...}]
tu.run_tool("HOCOMOCO_get_motif", {"motif_id": "CTCF.H14CORE.0.P.B", "include_pwm": true})
#   -> real verified response: consensus "hbRCCRShAGRKGGCGShvn", quality "B",
#      tfclass (superclass/class/family/subfamily), uniprot_ac_human "P49711",
#      pwm (if include_pwm=true)

Quality grades run A (best) to D (lowest); data_sources "P" = ChIP-Seq peaks, "S" = HT-SELEX. Use HOCOMOCO's quality grade to sanity-check a JASPAR/MEME-discovered motif for the same TF — a low grade or missing entry is itself informative (that TF's binding may be less well characterized).

Repeat Elements (Dfam)

Transposable elements and other repeats can create, destroy, or shift regulatory elements (SVA/Alu insertions carrying their own TF binding sites is a well-known mechanism) — a fifth context worth checking alongside the four evidence types in Domain Reasoning above when a regulatory element sits inside or near a repeat-masked region.

tu.run_tool("Dfam_search_families", {"name_prefix": "AluS"})
#   -> [{accession: "DF000000003", name: "AluSc", repeat_type: "SINE", ...}]
tu.run_tool("Dfam_get_family", {"accession": "DF000000003"})
#   -> consensus_sequence, classification, clades[], citations[]

# Real example: TP53 locus (chr17:7,661,779-7,687,538, GRCh38) — the same
# region used in the Ensembl Regulatory Build section above — carries two
# real SVA retroposon insertions (SVA_E, SVA_A) and an AluSc SINE:
tu.run_tool("Dfam_get_annotations", {
    "assembly": "hg38", "chrom": "chr17", "start": 7661779, "end": 7687538
})
#   -> [{accession: "DF000001071", query: "SVA_E", type: "Retroposon",
#        bit_score: 180, e_value: "1.3e-50", seq_start: 7661463, ...}, ...]

nrph=true (default) returns non-redundant profile hits only; set false for every raw hit including overlapping/lower-scoring alignments.

Phase 2: ENCODE Experiment Search

When looking for ChIP-seq, ATAC-seq, or other functional genomics data:

For TF ChIP-seq:

ENCODE_search_experiments(
    assay_title="TF ChIP-seq",
    target="CTCF",              # TF gene name
    biosample_term_name="HepG2", # Cell line or tissue
    limit=10
)

For histone marks:

ENCODE_search_histone_experiments(
    histone_mark="H3K27ac",         # or H3K4me3, H3K27me3, H3K36me3
    biosample_term_name="liver",
    limit=10
)

For chromatin accessibility:

ENCODE_search_chromatin_accessibility(
    biosample_term_name="T cell",
    limit=10
)

For regulatory annotations (cCREs, ChromHMM):

ENCODE_search_annotations(
    annotation_type="candidate Cis-Regulatory Elements",
    biosample_term_name="K562",
    limit=10
)

Common assay_title values:

• "TF ChIP-seq" - Transcription factor binding
• "Histone ChIP-seq" - Histone modification
• "ATAC-seq" - Chromatin accessibility
• "DNase-seq" - Open chromatin (older method)
• "WGBS" - DNA methylation

Note: ENCODE_search_experiments returns experiment metadata only (accession, biosample, status). Use ENCODE_get_experiment(accession) to get file download links and detailed metadata.

Phase 3: cCRE Annotation (UCSC + ENCODE)

When annotating a specific genomic region:

UCSC_get_encode_cCREs(
    chrom="chr8",       # Chromosome (GRCh38)
    start=37966000,     # Start coordinate
    end=37967000        # End coordinate
)
# Returns cCREs with type: pELS (proximal enhancer), dELS (distal enhancer),
# PLS (promoter-like), CTCF-only, DNase-H3K4me3

cCRE Types:

• PLS (Promoter-like): High DNase + H3K4me3 + H3K27ac signal near TSS
• pELS (Proximal Enhancer): High DNase + H3K27ac, within 2kb of TSS
• dELS (Distal Enhancer): High DNase + H3K27ac, >2kb from TSS
• CTCF-only: CTCF binding without enhancer marks
• DNase-H3K4me3: Unclassified accessible region

Phase 4: Regulatory Variant Scoring (RegulomeDB)

When assessing regulatory impact of a variant:

RegulomeDB_query_variant(rsid="rs4994")
# Returns:
#   regulome_score.ranking: "1a"-"7" (1a = highest regulatory evidence)
#   regulome_score.probability: 0-1 continuous score
#   tissue_specific_scores: dict of tissue -> score
#   overlapping features: eQTLs, TF binding, DNase peaks, motifs

RegulomeDB Score Interpretation:

Rank	Meaning
1a	eQTL + TF binding + matched TF motif + DNase peak
1b	eQTL + TF binding + DNase peak
1c	eQTL + TF binding or DNase peak
1d	eQTL + motif or protein binding
1e	eQTL + motif hit
1f	eQTL only
2a	TF binding + motif match + DNase
2b	TF binding + matched motif
2c	TF binding with/without motif
3a	DNase peak + motif
3b	DNase peak only
4	Motif hit only
5	Proximity to Footprint
6	Proximity to Footprint + TF
7	No evidence

Variants with rank 1a-2b are most likely to affect gene regulation.

---

Tool Parameter Reference

Tool	Required Params	Optional Params	Notes
jaspar_search_matrices	(none — returns all if empty)	name, species, collection, tax_id, page, page_size	Use name for TF name search
jaspar_get_matrix	matrix_id	—	Full version required: "MA0139.1" not "MA0139"
JASPAR_get_transcription_factors	(none)	page, page_size	Paginated; default page_size=10
jaspar_get_matrix_versions	base_id	—	base_id is unversioned (e.g., "MA0139")
ENCODE_search_experiments	(none — returns all if empty)	assay_title, target, biosample_term_name, limit	assay_title must match ENCODE vocabulary exactly
ENCODE_search_histone_experiments	(none)	histone_mark, biosample_term_name, limit	histone_mark: "H3K27ac", "H3K4me3", etc.
ENCODE_search_chromatin_accessibility	(none)	biosample_term_name, limit	Returns ATAC-seq and DNase-seq
ENCODE_get_experiment	accession	—	accession: "ENCSR..." format
ENCODE_search_annotations	(none)	annotation_type, biosample_term_name, limit	annotation_type: "candidate Cis-Regulatory Elements"
ENCODE_get_chromatin_state	(none)	biosample_term_name, limit	Returns ChromHMM segmentation
ENCODE_search_biosamples	(none)	term_name, biosample_type, limit	biosample_type: "cell line", "tissue", "primary cell"
UCSC_get_encode_cCREs	chrom, start, end	—	Coordinates in GRCh38; chrom format: "chr1"
RegulomeDB_query_variant	rsid	—	rsid format: "rs4994" (with rs prefix)

---

Common Patterns

Pattern 1: TF Binding Site Investigation

Goal: Find where TF X binds and what motif it recognizes
Flow:
  1. jaspar_search_matrices(name="CTCF") -> get matrix_id
  2. jaspar_get_matrix(matrix_id) -> get full PWM, logo URL
  3. ENCODE_search_experiments(assay_title="TF ChIP-seq", target="CTCF") -> experimental binding data
  4. For specific tissue: add biosample_term_name="HepG2"
Output: Motif logo + experimental binding evidence

Pattern 2: Regulatory Variant Interpretation

Goal: Assess if variant rs1234567 affects gene regulation
Flow:
  1. RegulomeDB_query_variant(rsid="rs1234567") -> score + overlapping features
  2. If score <= 2b: ENCODE_search_experiments(target=overlapping_TF) -> experimental evidence
  3. UCSC_get_encode_cCREs(chrom, start, end) -> check if variant in known cCRE
Output: Regulatory score + supporting evidence + cCRE context