# Bio Applied Regulatory Analysis

> Scan DNA for promoter/regulatory elements: TATA box regex search, CpG island detection via GC%/observed-over-expected sliding windows, and PFM to PWM (log-odds) construction/scanning for TFBS. Use when locating a TSS, calling CpG islands, building a position weight matrix from aligned binding sites, or scanning a promoter with a JASPAR/TRANSFAC-style motif.

- Skill: `pavel-kravchenko/bio-applied-regulatory-analysis` (Agent Skill)
- Install (CLI): `npx skillmds@latest add pavel-kravchenko/bio-applied-regulatory-analysis`
- Raw SKILL.md: https://api.skillmd.com/api/skills/pavel-kravchenko/bio-applied-regulatory-analysis/raw
- Safety review: pending
- Works with: Claude Code, Claude.ai, OpenAI Codex
- Category: Coding & Dev Tools
- Author: pavel-kravchenko (https://skillmd.com/u/pavel-kravchenko)
- Updated: 2026-09-17
- Page: https://skillmd.com/skills/pavel-kravchenko/bio-applied-regulatory-analysis

---


# Promoter and Regulatory Sequence Analysis

## When to Use

- Searching a promoter/upstream region for a TATA box, Inr, or other core promoter motif
- Calling CpG islands (or checking a region against the classic ≥200bp / GC≥50% / O:E≥0.6 rule)
- Building a PFM/PWM from a set of aligned transcription factor binding sites and scanning a sequence for hits
- Profiling GC content, dinucleotide frequency, or motif density relative to a known/candidate TSS
- Testing whether a TFBS is enriched in one gene set (e.g. stress-responsive) vs. another (e.g. housekeeping)

## Version Compatibility

- Python ≥ 3.10, NumPy ≥ 1.24, pandas ≥ 2.0, SciPy ≥ 1.11 (for `scipy.stats.mannwhitneyu`)
- No genome-specific databases required — all examples work on plain strings; swap in real sequences via BioPython (`Bio.SeqIO`) as needed

## Prerequisites

- `pip install numpy pandas scipy biopython`
- Familiarity with Python regex (`re` module) and basic sequence coordinates (TSS = +1, upstream = negative)
- Related skill: `bio-sequence-manipulation-motif-search` for generic motif scanning; `bio-chip-seq-motif-analysis` for ChIP-seq-derived motifs

## Key Regulatory Elements

| Element | Location | Function |
|---------|----------|----------|
| Promoter | −1 to −1000 bp from TSS | Recruits RNA Pol II |
| Enhancer | Distal (kb–Mb away) | Boosts transcription |
| Silencer | Variable | Represses transcription |
| Insulator | Between elements | Blocks enhancer–promoter crosstalk |

- **TSS** = position +1; upstream positions are negative
- **TATA box**: consensus `TATAAA` at ~−30; present in ~10–20% of human genes (TATA-less promoters use Inr, DPE instead)
- **CpG islands**: near ~70% of human gene promoters; classic criteria: length ≥200 bp, GC ≥50%, CpG observed/expected (O/E) ≥0.6

## TATA Box Detection

**Goal:** find TATA-box-like motifs in a promoter sequence and report their position relative to the TSS.
**Approach:** regex search for the exact consensus (`TATAAA`) or the relaxed IUPAC pattern (`TATAWAW`, W = A/T).

```python
import re


def find_tata_boxes(sequence, strict=True):
    """Find TATA box motifs in a DNA sequence.

    strict=True:  exact consensus TATAAA
    strict=False: relaxed IUPAC pattern TATA[AT]A[AT] (TATAWAW)
    Returns a list of (start_position, matched_motif) tuples.
    """
    sequence = sequence.upper()
    pattern = 'TATAAA' if strict else 'TATA[AT]A[AT]'
    return [(m.start(), m.group()) for m in re.finditer(pattern, sequence)]


# Example: TSS at position 1000, TATA box inserted at -30
tss = 1000
promoter_seq = "N" * 970 + "TATAAAG" + "N" * 1023
for pos, motif in find_tata_boxes(promoter_seq):
    print(f"Position {pos} (TSS{pos - tss:+d}): {motif}")
```

## CpG Island Detection

**Goal:** call CpG islands in a promoter and merge overlapping sliding-window hits into island coordinates.
**Approach:** slide a window across the sequence, compute GC% and CpG O/E per window, keep windows passing both thresholds, then merge adjacent/overlapping passing windows.

```python
def cpg_island_scanner(sequence, window=200, step=10, gc_thresh=0.5, oe_thresh=0.6):
    """Sliding-window CpG island detection.

    Returns a list of (start, end, gc_content, cpg_oe_ratio) for qualifying windows.
    Use step=1 for exact boundaries (slower); step=10-50 for a fast first pass.
    """
    sequence = sequence.upper()
    islands = []
    for i in range(0, len(sequence) - window + 1, step):
        win = sequence[i:i + window]
        n_c, n_g = win.count('C'), win.count('G')
        n_cpg = win.count('CG')
        gc = (n_c + n_g) / window
        oe = (n_cpg * window) / (n_c * n_g) if n_c > 0 and n_g > 0 else 0.0
        if gc >= gc_thresh and oe >= oe_thresh:
            islands.append((i, i + window, gc, oe))
    return islands


def merge_island_windows(windows, max_gap=0):
    """Merge overlapping/adjacent qualifying windows into contiguous CpG islands."""
    if not windows:
        return []
    sorted_wins = sorted(windows, key=lambda x: x[0])
    merged = [list(sorted_wins[0])]
    for start, end, gc, oe in sorted_wins[1:]:
        if start <= merged[-1][1] + max_gap:
            merged[-1][1] = max(merged[-1][1], end)
            merged[-1][2] = max(merged[-1][2], gc)
            merged[-1][3] = max(merged[-1][3], oe)
        else:
            merged.append([start, end, gc, oe])
    return merged


raw = cpg_island_scanner(promoter_seq, window=200, step=10)
for start, end, gc, oe in merge_island_windows(raw, max_gap=50):
    print(f"Island {start}-{end} ({end - start} bp)  GC={gc:.2f}  CpG O/E={oe:.2f}")
```

## PWM/PFM Construction and Scanning

**Goal:** build a transcription-factor motif model from aligned binding sites and scan a sequence for matches (JASPAR/TRANSFAC-style TFBS scan).
**Approach:** count bases per position (PFM), convert to log2-odds vs. background (PWM), then slide the motif across the target sequence and sum log-odds scores.

```python
import numpy as np


def build_pfm(sites):
    """Build a position frequency matrix from aligned binding sites (same length)."""
    length = len(sites[0])
    pfm = {b: [0] * length for b in 'ACGT'}
    for site in sites:
        for i, base in enumerate(site.upper()):
            pfm[base][i] += 1
    return pfm


def pfm_to_pwm(pfm, pseudocount=0.5, bg=None):
    """Convert a PFM to a log2-odds PWM (adds a pseudocount to avoid -inf scores)."""
    bg = bg or {'A': 0.25, 'C': 0.25, 'G': 0.25, 'T': 0.25}
    n = sum(pfm[b][0] for b in 'ACGT')
    length = len(pfm['A'])
    pwm = {b: [0.0] * length for b in 'ACGT'}
    for b in 'ACGT':
        for i in range(length):
            freq = (pfm[b][i] + pseudocount) / (n + 4 * pseudocount)
            pwm[b][i] = np.log2(freq / bg[b])
    return pwm


def score_sequence(subseq, pwm):
    """Score one sequence window against a PWM: sum of per-position log-odds."""
    return sum(pwm[b][i] for i, b in enumerate(subseq.upper()) if b in pwm)


def scan_with_pwm(sequence, pwm, threshold=0.0):
    """Slide a PWM across a sequence; return (pos, subseq, score) above threshold."""
    sequence = sequence.upper()
    motif_len = len(pwm['A'])
    hits = []
    for i in range(len(sequence) - motif_len + 1):
        sub = sequence[i:i + motif_len]
        score = score_sequence(sub, pwm)
        if score >= threshold:
            hits.append((i, sub, score))
    return hits


tata_sites = ['TATAAAG', 'TATAAAT', 'TATAAAA', 'TATATAG', 'TATAAAC',
              'TATAAAT', 'TATAAAG', 'TATAAAT', 'TATAAAA', 'TATAAAG']
pwm = pfm_to_pwm(build_pfm(tata_sites))
for pos, sub, score in sorted(scan_with_pwm(promoter_seq, pwm, threshold=5.0),
                               key=lambda x: -x[2])[:5]:
    print(f"Position {pos} (TSS{pos - tss:+d}): {sub}  score={score:.2f}")
```

## TSS Prediction Signals

| Signal | Peak location | Method |
|--------|--------------|--------|
| TATA box | −30 | Motif scan |
| CpG island | centered on TSS | O/E ratio |
| TFBS density | −200 upstream | PWM scan |
| CAGE signal | +1 | Experimental |

## TF Motif Databases

| Database | Description | URL |
|----------|-------------|-----|
| JASPAR | Open-access, curated | jaspar.elixir.no |
| HOCOMOCO | Human/mouse from ChIP-seq | hocomoco11.autosome.org |
| TRANSFAC | Comprehensive (commercial) | genexplain.com/transfac |

## Pitfalls

- **CpG underrepresentation**: vertebrate bulk genome has CpG O/E ~0.2 due to methylation-driven mutation; islands (O/E ≥0.6) are real anomalies, not noise
- **Coordinate systems**: BED = 0-based half-open; VCF/GFF = 1-based inclusive — mixing causes off-by-one errors when reporting TSS-relative positions
- **PWM threshold selection**: score ≥80% of the max possible PWM score is a common heuristic; too low a threshold floods results with false positives
- **Pseudocount matters**: without pseudocounts, a single zero count in the PFM gives a −∞ PWM score for any sequence with that base at that position
- **Multiple testing**: scanning a genome-wide promoter set for TFBS enrichment requires FDR correction (e.g. Benjamini-Hochberg), not raw p-values
- **Window step size**: `step=1` in `cpg_island_scanner` gives exact boundaries but is O(n); use a larger step for a fast first pass, then refine only around merged candidate islands

## See Also

- `bio-sequence-manipulation-motif-search` — generic sequence motif searching
- `bio-chip-seq-motif-analysis` — motif discovery/enrichment from ChIP-seq peaks
- `bio-chip-seq-peak-annotation` — annotating peaks relative to TSS/promoters
- `jaspar-database` — fetching real PFMs/PWMs from JASPAR