End-to-end CLIP-seq pipeline from FASTQ to ENCODE-compliant binding sites, single-nucleotide crosslink maps, annotation, motifs, and (optionally) differential binding. Use when running the full Yeo lab eCLIP / iCLIP / iCLIP2 / iCLIP3 / irCLIP / PAR-CLIP analysis with SMInput control, protocol-specific UMI extraction, ENCODE STAR parameters, CLIPper or Skipper peak calling with stringent log2 FC and -log10 p thresholds, IDR rescue and self-consistency QC, and downstream motif registration with mCross or PEKA.
End-to-end CLIP-seq pipeline from FASTQ to ENCODE-compliant binding sites, single-nucleotide crosslink maps, annotation, motifs, and (optionally) differential binding. Use when running the full Yeo lab eCLIP / iCLIP / iCLIP2 / iCLIP3 / irCLIP / PAR-CLIP analysis with SMInput control, protocol-specific UMI extraction, ENCODE STAR parameters, CLIPper or Skipper peak calling with stringent log2 FC and -log10 p thresholds, IDR rescue and self-consistency QC, and downstream motif registration with mCross or PEKA.
Before using code patterns, verify installed versions match. If versions differ:
CLI: <tool> --version then <tool> --help to confirm flags
Python: pip show <package> then help(module.function) to check signatures
R: packageVersion('<pkg>') then ?function_name to verify parameters
If code throws unexpected errors, introspect the installed tool and adapt the example rather than retrying.
CLIP-seq End-to-End Pipeline
"Analyze my CLIP-seq data from raw FASTQ to ENCODE-compliant binding sites" -> Orchestrate protocol-specific UMI extraction, 3'-only adapter trimming (preserving the R2 5' truncation = crosslink site -1), ENCODE STAR alignment, UMI-based deduplication, library complexity QC, peak calling against SMInput with stringent thresholds (log2 FC >= 3 AND -log10 p >= 3), single-nucleotide crosslink-site detection, ChIPseeker annotation with CLIP-appropriate tssRegion, motif discovery with GC-matched background and CL-position registration, and optional differential binding between conditions.
Pipeline Overview
FASTQ + SMInput
-> [clip-preprocessing] UMI extract + 3' adapter trim (-q 6 -m 18) + two-pass for eCLIP
-> [clip-alignment] STAR ENCODE block (alignEndsType EndToEnd, mismatch 0.04 or 0.07 for PAR-CLIP) + UMI dedup
-> [clip-qc] preseq, FRiP, IDR rescue + self-consistency, read distribution
-> [clip-peak-calling] CLIPper + SMInput log2 norm (stringent: log2 FC >= 3, -log10 p >= 3) OR Skipper (210-320% more sites)
-> [crosslink-site-detection] PureCLIP or CTK CITS for single-nt CL positions
-> [binding-site-annotation] ChIPseeker (tssRegion=c(-100,100), level=transcript) + RBP-Maps for splicing factors
-> [clip-motif-analysis] HOMER + mCross (registered) + RBNS Kd cross-check
-> [differential-clip] DEWSeq window-level NB with type:condition interaction (optional)
CLIP Variant Selection
Variant
When to use
UMI pattern
STAR mismatch ceiling
Detection signal
eCLIP (Van Nostrand 2016)
ENCODE comparability; SMInput available
10 nt R1
0.04
R2 5' truncation
iCLIP / iCLIP2 / iCLIP3
Single-end; high motif specificity
NNNXXXXNN (3+4+2; demux first)
0.04
R1 5' truncation
irCLIP / FLASH
Non-radioactive; fast
Protocol-specific
0.04
Truncation
PAR-CLIP
Photoactivatable nucleoside (4SU); HEK293/K562
4 nt typical
0.07 (raised for T->C)
T->C transitions
miCLIP / miCLIP2
m6A modification
iCLIP-style
0.04
Truncation + C->T at m6A
STAMP / scSTAMP
Antibody-free; in vivo or single-cell
NA (no UV)
0.04 (RNA-seq mode)
C->U editing (RBP-APOBEC1 fusion)
chimeric eCLIP / miR-eCLIP
Direct miRNA-target pairs
10 nt R1
0.04
Chimeric reads
Step 1: Quality Control of Raw FASTQ
# Initial QC
fastqc raw_R1.fq.gz raw_R2.fq.gz -o qc/raw/
# Inspect first 12 bases of 100 reads to verify UMI pattern matches the prep
zcat raw_R1.fq.gz | awk 'NR%4==2' | head -100 | cut -c1-12 | sort | uniq -c | sort -rn | head# Random barcode positions show ~25% per base; library barcodes are fixed
Step 2: Preprocessing (Protocol-Specific)
Goal: Convert raw CLIP FASTQ into UMI-deduplicated, alignment-ready FASTQ while preserving the R2 5' end (= crosslink site -1) that drives single-nucleotide resolution downstream.
Approach: Use the protocol-matched UMI pattern (10 nt eCLIP, NNNXXXXNN iCLIP, 4 nt PAR-CLIP), run umi_tools extract to move random barcodes to read names, then apply cutadapt with 3'-only adapter trimming at -q 6 -m 18 (permissive 5' to protect the truncation base). eCLIP uses two-pass trimming to remove read-through inline adapters from R2 5' only; iCLIP and PAR-CLIP use single-pass.
# eCLIP: 10 nt UMI on R1; two-pass adapter trim for read-through# See clip-seq/clip-preprocessing for protocol-specific patterns
umi_tools extract \
--bc-pattern=NNNNNNNNNN \
--stdin=raw_R1.fq.gz --read2-in=raw_R2.fq.gz \
--stdout=R1.umi.fq.gz --read2-out=R2.umi.fq.gz \
--log=qc/umi_extract.log
# Pass 1: 3' adapter on both reads# -q 6 is intentionally permissive; aggressive trimming destroys R2 5' = CL site -1
cutadapt \
-a AGATCGGAAGAGCACACGTCT \
-A AGATCGGAAGAGCGTCGTGTAGGGAAAGAGTGT \
--quality-base 33 -q 6 -m 18 \
-j 8 \
-o R1.p1.fq.gz -p R2.p1.fq.gz \
R1.umi.fq.gz R2.umi.fq.gz \
> qc/cutadapt_pass1.log 2>&1
# Pass 2: strip read-through 5' adapter from R2 only (NEVER -g on R1)
cutadapt \
-G GATCGTCGGACTGTAGAACTCTGAAC \
--quality-base 33 -q 6 -m 18 \
-j 8 \
-o R1.trim.fq.gz -p R2.trim.fq.gz \
R1.p1.fq.gz R2.p1.fq.gz \
>> qc/cutadapt_pass2.log 2>&1
For PAR-CLIP: same UMI extraction but downstream alignment raises --outFilterMismatchNoverReadLmax from 0.04 to 0.07 (the T->C signature would otherwise be filtered as sequencing error). See clip-seq/clip-preprocessing for full per-protocol guidance.
For PAR-CLIP: change --outFilterMismatchNoverReadLmax 0.04 to 0.07. For repeat-binding RBPs (MATR3, ZFP36, FUS at LINE-1, HNRNPK at SINEs): change --outFilterMultimapNmax 1 to 100 and add --outSAMmultNmax -1, then run CLAM downstream for EM-based multi-mapper assignment. See clip-seq/clip-alignment for full guidance.
CLIP libraries have 40-70% PCR duplication BY DESIGN (the IP enriches a small molecule pool). Low duplication usually means failed IP, not a good library. The unique-fragment count after UMI dedup is the actual quality metric. See clip-seq/clip-qc for full five-gate diagnostic.
For maximum sensitivity (210-320% more sites than CLIPper for mRNA-binding RBPs), use Skipper Snakemake workflow with the same SMInput control. Mandatory for FASTKD2 / mt-RBPs which CLIPper misses on chrM. See clip-seq/clip-peak-calling for the full caller taxonomy.
Single-nt CL sites feed mCross motif registration and allele-specific binding analyses. They are NOT a replacement for the broad peak list; complementary outputs. See clip-seq/crosslink-site-detection.
Step 7: IDR Across Replicates
# Sort each replicate's compressed BED by signal (log2 FC, column 5)sort -k5,5gr peaks/rep1.compressed.bed > peaks/rep1.sorted.bed
sort -k5,5gr peaks/rep2.compressed.bed > peaks/rep2.sorted.bed
# True replicates threshold 0.05
idr --samples peaks/rep1.sorted.bed peaks/rep2.sorted.bed \
--input-file-type bed --rank 5 \
--output-file qc/idr.true.out \
--idr-threshold 0.05 \
--plot --log-output-file qc/idr.log
# ENCODE rule: rescue + self-consistency ratios both < 2 to pass# Pseudo-replicate IDR (split BAM in half) at threshold 0.10
Default ChIPseeker tssRegion=c(-3000, 3000) over-extends for CLIP (would label 30-50% peaks as "Promoter"). Splicing factors additionally need RBP-Maps (Yeo lab) for the 1400 nt cassette-exon regulatory metagene. See clip-seq/binding-site-annotation.
UV254 crosslinking has a strong U bias (~60-80% CL events at U); naive logos centered on CL positions are U-enriched even for non-U-binding RBPs. mCross corrects this by registering motif relative to the CL offset. See clip-seq/clip-motif-analysis.
Step 10: Differential Binding (Optional, Across Conditions)
# DEWSeq window-level NB with the interaction-term design# The interaction `~ type + condition + type:condition` tests whether IP/SMInput ratio shifts;# naive `~ condition` confounds binding with expression changes.
library(DEWSeq)
counts <- read.table('counts/merged.tsv', sep='\t', header=TRUE, row.names=1)
colData <- data.frame(
type =c('ip','ip','ip','ip','sminput','sminput','sminput','sminput'),
condition =c('treat','treat','ctrl','ctrl','treat','treat','ctrl','ctrl'))
dds <- DESeqDataSetFromSlidingWindows(
countData=counts, colData=colData,
annotObj='annotation_windows.bed',
design =~ type + condition + type:condition
)
dds <- DESeq(dds)
res <- results(dds, name='typeip.conditiontreat')
See clip-seq/differential-clip for full DEWSeq workflow and the htseq-clip preprocessing required upstream.
