crispr-screens/combinatorial-screens/SKILL.md
Designs and analyzes combinatorial CRISPR screens covering paired-Cas9 (Big Papi, Najm 2018), enhanced AsCas12a multiplex (enCas12a, DeWeirdt 2021), in4mer 4-guide-array Cas12a (Esmaeili Anvar N et al 2024 Nat Commun 15:3577) and the Inzolia paralog-pair library, paralog-buffering detection (Dede 2020 Genome Biol; Thompson 2021 Nat Commun 12:1302), genetic-interaction (GI) scoring as observed_double_LFC minus expected_additive_double_LFC, synthetic-lethal and synthetic-rescue interaction interpretation, the half-of-essentiality buffered by paralogs phenomenon, multiplex screen statistical analysis with MAGeCK MLE interaction terms, and the relationship to single-cell combinatorial Perturb-seq. Use when designing a paralog or pathway-pair screen, choosing between paired-Cas9 (Big Papi) and Cas12a multiplex (Inzolia), interpreting genetic interaction scores, identifying synthetic-lethal targets for drug development, or scaling beyond single-gene CRISPR screens.
npx skillsauth add GPTomics/bioSkills bio-crispr-screens-combinatorial-screensInstall this skill globally with one command. Works with Claude Code, Cursor, and Windsurf.
3 of 9 scanners reported clean
Some scanners were skipped, did not run, or reported a non-clean status. Review each row below.
Reference examples tested with: MAGeCK 0.5.9+ (for MLE with interaction terms), Inzolia library annotation (Esmaeili Anvar 2024), pandas 2.2+, numpy 1.26+, scipy 1.12+, matplotlib 3.8+.
Before using code patterns, verify installed versions match. If versions differ:
mageck --version; mageck mle --helpIf code throws ImportError, AttributeError, or TypeError, introspect the installed package and adapt the example to match the actual API rather than retrying.
"Run a combinatorial CRISPR screen to find synthetic-lethal interactions" -> Design a paired or multiplex library, screen for double-knockout fitness, score per-pair genetic interaction (GI = observed_double - expected_additive), and identify synthetic-lethal (negative GI) and synthetic-rescue (positive GI) interactions.
mageck mle with explicit interaction terms for paired-Cas9 (Big Papi-style)| Goal | Architecture | Library | Why | |------|--------------|---------|-----| | Paralog buffering, identify synthetic lethal paralog pairs | enCas12a single-array 4-guide multiplex | Inzolia (Esmaeili Anvar 2024) | Cas9 single-KO misses paralog-buffered essentials (42% of constitutively expressed genes never score, Dede 2020) | | Test specific pathway pair (e.g., DNA repair branches) | Big Papi (orthologous SaCas9 + SpCas9; two sgRNAs from U6 and H1 in pPapi) | Custom | Mature methodology; orthologous enzymes avoid repeated-element recombination | | Combinatorial 3-way / 4-way knockout | in4mer (4-guide single Cas12a array) | Custom (in4mer) | Single transcript processed by Cas12a; multi-gene | | Single-cell Perturb-seq with multi-pert per cell | Combinatorial Perturb-seq + Cas9 multiplex | Custom | Single-cell readout of multi-perturbation effects | | Drug-modifier + KO interaction | Cas9 KO + drug treatment | Standard libraries | Drug as second "perturbation" |
Fails when:
| Property | Cas9 paired (Big Papi) | Cas12a multiplex (in4mer / Inzolia) | |----------|------------------------|--------------------------------------| | Multiplex capacity per cassette | 2 sgRNAs (paired) | 4 (in4mer); 2 (standard Cas12a) | | sgRNA processing | U6 and H1 promoters driving sgRNAs for two orthologous Cas9s | Single transcript processed by Cas12a itself | | sgRNA inhibition with multiple targets | None | None (Cas12a's intrinsic processing handles all) | | Library size for 1,000 pairs | ~4,000-6,000 paired cassettes (4-6 per pair) | ~2,000 arrays (2 per pair) plus singleton controls | | Validated libraries | Limited (mostly custom) | Inzolia: ~49k 4-guide arrays covering 19,687 genes plus ~4,435 paralog pairs | | Per-perturbation editing efficiency | High (each sgRNA independently) | Variable (Cas12a less efficient on some targets) | | Best for | Pairwise GI of specific interest | Genome-scale paralog buffering; multi-gene perturbation |
Recommendation: For modern paralog screens, use Cas12a multiplex with the Inzolia library. It is ~30% smaller than a typical monogenic Cas9 library while additionally covering ~4,000 paralog pairs (Esmaeili Anvar 2024), which makes it more cost-effective at genome scale.
Dede et al 2020 Genome Biol 21:262 showed that a large share of constitutively expressed genes are never scored as essential in any Cas9 single-KO fitness screen (3,032 of 7,282; 42%), and that these never-essentials are strongly enriched for paralogs. The reason: gene paralogs perform redundant essential functions. Loss of one paralog is buffered by the other; only loss of both creates the essentiality phenotype.
Quantified impact: 24 synthetic-lethal paralog pairs identified in Dede 2020 across 3 cell lines; 19 of 24 (79%) reproduce in >=2 lines, 14 of 24 (58%) in all 3. These pairs were not findable by single-gene Cas9 screens, requiring combinatorial methodology.
Examples:
Each is buffered: loss of one is tolerated; loss of both is lethal.
Goal: Identify pairs where the double-knockout fitness differs from the additive expectation.
Approach: From per-pair and per-singleton fitness data, compute GI = observed_double_LFC - (single_A_LFC + single_B_LFC). Synthetic lethal: GI < threshold (more depleted than additive). Synthetic rescue: GI > threshold (less depleted than additive).
import pandas as pd
import numpy as np
from scipy.stats import zscore
def gi_score(paired_lfc_df, single_lfc_df):
'''Score genetic interactions from paired vs single LFCs.
paired_lfc_df: rows = paired-KO; columns = ['gene_A', 'gene_B', 'paired_lfc']
single_lfc_df: rows = single-KO; columns = ['gene', 'single_lfc']
'''
single = dict(zip(single_lfc_df['gene'], single_lfc_df['single_lfc']))
df = paired_lfc_df.copy()
df['single_A_lfc'] = df['gene_A'].map(single)
df['single_B_lfc'] = df['gene_B'].map(single)
df['expected_additive'] = df['single_A_lfc'] + df['single_B_lfc']
df['gi_score'] = df['paired_lfc'] - df['expected_additive']
df = df.dropna(subset=['gi_score']) # a single missing singleton would NaN every z-score
df['gi_z'] = zscore(df['gi_score'])
df['gi_class'] = np.where(df['gi_z'] < -2, 'synthetic_lethal',
np.where(df['gi_z'] > 2, 'synthetic_rescue', 'no_interaction'))
return df.sort_values('gi_z')
Interpretation:
Goal: Use MAGeCK MLE to estimate the effect of each gene independently and the additional effect when both genes are simultaneously perturbed.
Approach: Design matrix encodes single-A, single-B, double-AB conditions; the interaction column is set to 1 only for double-KO samples. The resulting beta for that column captures the extra effect beyond the sum of single-gene betas. Note: MAGeCK MLE does not natively perform a formal interaction-significance test, but the interaction|beta and |fdr columns serve as the GI estimate; for formal interaction testing, compute GI = observed_double_lfc - (single_A_lfc + single_B_lfc) explicitly (see GI scoring section below).
# Design matrix encoding double-KO as a separate "interaction" indicator
# Conditions: NT (control), A_KO, B_KO, A_B_KO
cat > combo_design.txt <<EOF
Samples baseline geneA geneB interaction
NT_r1 1 0 0 0
NT_r2 1 0 0 0
A_r1 1 1 0 0
A_r2 1 1 0 0
B_r1 1 0 1 0
B_r2 1 0 1 0
AB_r1 1 1 1 1
AB_r2 1 1 1 1
EOF
mageck mle \
--count-table combo_counts.txt \
--design-matrix combo_design.txt \
--output-prefix combo_mle
# Output: per-gene beta scores per design column
# The "interaction" column beta captures additional joint effect beyond additive
Interpretation of MAGeCK MLE output:
| Column | Meaning |
|--------|---------|
| geneA|beta | Single-A effect |
| geneB|beta | Single-B effect |
| interaction|beta | Additional effect under joint perturbation beyond sum of singles |
| interaction|p-value, |fdr | Significance vs zero |
A significantly negative interaction|beta is synthetic lethal; positive is synthetic rescue. For formal GI hypothesis testing, prefer the explicit GI scoring approach (next section) over MAGeCK MLE interpretation, since MAGeCK MLE does not validate the additive null.
Esmaeili Anvar 2024 Nat Commun 15:3577 introduced in4mer, a Cas12a multiplex architecture where each array contains 4 guides processed by Cas12a's intrinsic crRNA-processing activity. The Inzolia library is the canonical implementation, covering the protein-coding genome plus ~4,435 paralog pairs.
Library design:
# Per-pair analysis from in4mer screen
def in4mer_pair_analysis(paired_counts_df, gene_pairs, value_cols):
'''Aggregate cassette-level counts to per-pair statistics.
paired_counts_df: rows = cassettes, with a cassette_id COLUMN (reset_index first if it is the index).
gene_pairs: DataFrame with cassette_id and gene_A, gene_B columns.
value_cols: the numeric sample/LFC columns to aggregate.
'''
merged = paired_counts_df.merge(gene_pairs, on='cassette_id')
return merged.groupby(['gene_A', 'gene_B'])[value_cols].agg(['mean', 'std', 'count'])
Trigger: A dual-sgRNA construct built from repeated elements -- two copies of the U6 promoter, or two copies of the SpCas9 tracrRNA scaffold. Mechanism: Najm 2018 reports that repetitive elements in lentiviral vectors, including the U6 promoter and multiple copies of the tracrRNA sequence, drive high levels of recombination and reduce combinatorial screen efficiency. Big Papi avoids this by pairing two orthologous enzymes (SaCas9 + SpCas9), whose scaffolds differ, and expressing the two sgRNAs from distinct U6 and H1 promoters. Symptom: Constructs collapse to a single perturbation; measured GI scores are diluted toward zero. Fix: Use the pPapi architecture (orthologous Cas9s, U6 + H1) rather than duplicated U6/tracr elements; verify construct integrity by amplicon sequencing of clones.
Trigger: Cas12a less efficient than Cas9 at some loci; some guides in the 4-guide array don't cut. Mechanism: Cas12a editing rate varies by sequence context; some loci edit at <30%. Symptom: Specific pairs missing expected effects despite cassette presence. Fix: Pilot Cas12a efficiency at the loci before full screen; use enCas12a (enhanced) variant; for known low-efficiency loci, supplement with Cas9.
Trigger: Library lacks single-gene controls (only paired knockouts). Mechanism: GI = paired - expected_additive requires single-gene LFC; without them, expected cannot be computed. Symptom: Cannot score GI; only paired LFCs available. Fix: Design library to include singletons (place gene A with 3 placeholder guides; gene B with 3 placeholders); re-run with full design.
Trigger: Using public single-gene LFCs (e.g., DepMap) as the baseline for paired-screen GI scoring. Mechanism: Single-gene effects are cell-line specific; using HCT116 single-gene LFCs to score K562 paired-screen GIs is invalid. Symptom: GI scores look noisy; many false positives. Fix: Include singleton controls in the screen; or use cell-line-matched DepMap data.
Trigger: Paired KO of two cell-cycle-impacting genes; the double-effect saturates cell cycle. Mechanism: If A_KO causes 50% growth arrest and B_KO causes 50%, the combined 75% arrest is already saturating proliferation; additive expectation overestimates double-effect, generating false "synthetic-rescue." Symptom: GI scores positive for pairs of essential cell-cycle genes; biologically unexpected. Fix: Use log-space (LFC) GI scoring rather than linear; saturation is less severe in log-space. Alternative: model with logistic / saturable response curve.
Trigger: Inzolia library has uneven cassette representation; some pairs at 10x lower coverage than others. Mechanism: Standard library QC (Gini, skew) applies; low-coverage cassettes yield noisier LFCs. Symptom: GI z-scores vary 2-3x across cassettes targeting the same pair. Fix: Standard library QC; for low-coverage pairs, aggregate fewer cassettes but with more sequencing depth; or drop low-coverage pairs from analysis.
For high-stakes synthetic-lethal hits (drug-target nomination), validate by:
| Threshold | Value | Source / Rationale | |-----------|-------|--------------------| | Synthetic lethal GI z-score | <-2 | Standard convention | | Synthetic rescue GI z-score | >2 | Standard convention | | No interaction | -1 to +1 | Within additive expectation | | Cas9 paired-screen cassette count per pair | 4-6 | Standard library convention | | Cas12a 4-guide arrays per paralog pair (Inzolia) | 2 (2 guides per gene, order swapped) | Esmaeili Anvar 2024 | | Singletons in combinatorial library | At least 4-6 per single gene | For stable expected_additive | | Cells per cassette for stable GI | 500+ at infection | Standard pooled-screen coverage | | Cas12a editing efficiency for inclusion | >50% | Below = unreliable signal |
| Error / symptom | Cause | Solution | |-----------------|-------|----------| | Dual-sgRNA construct acts as single | Recombination between repeated U6/tracr elements | Use pPapi (orthologous SaCas9 + SpCas9, U6 + H1) | | Cas12a low editing | Locus-specific inefficiency | Pilot loci first; use enCas12a | | Cannot compute GI | No singletons in library | Re-design to include all-singletons | | GI scores noisy | Library skew | Standard library QC; aggregate cassettes | | Many false "rescue" GIs | Saturation in linear-space | Use log-space (LFC) GI scoring | | Drug-target paralog shows no GI in screen | Cell-line-specific buffering | Cross-validate with multiple lines |
development
Installs 425 bioinformatics skills covering sequence analysis, RNA-seq, single-cell, variant calling, metagenomics, structural biology, and 56 more categories. Use when setting up bioinformatics capabilities or when a bioinformatics task requires specialized skills not yet installed.
testing
Chains a somatic (tumor-normal) SNV/indel and structural-variant pipeline end to end with GATK Mutect2 (or Strelka2), wiring the somatic-specific machinery - panel-of-normals and gnomAD germline-resource priors, GetPileupSummaries/CalculateContamination, and LearnReadOrientationModel FFPE/oxoG orientation-bias filtering fed into FilterMutectCalls. Use when calling somatic mutations from a tumor-normal pair (or tumor-only with PoN caveats), deciding which artifact filter removes which class of false positive, reasoning about VAF/purity/ploidy and clonal-vs-subclonal detection, adding somatic SV/CNV or TMB/MSI/signatures, or routing variants to AMP/ASCO/CAP tier and oncogenicity interpretation (never germline ACMG).
development
End-to-end pooled and single-cell CRISPR screen analysis from FASTQ to hit genes. Orchestrates library design QC, guide counting, six-stage screen QC (plasmid Gini, replicate Pearson, CEGv2 PR-AUC, copy-number artifact), method-appropriate hit calling across MAGeCK RRA/MLE, BAGEL2, drugZ, JACKS, and Chronos, cancer-cell-line copy-number correction (CRISPRcleanR / Chronos), batch correction for multi-batch screens, and the specialized branches for combinatorial paralog screens, single-cell Perturb-seq, base-editor variant-function screens, prime-editor screens, and in vivo bottleneck-aware screens. Use when analyzing any pooled CRISPR screen end-to-end, matching the hit-calling method to the experimental design, integrating copy-number correction into the pipeline, or branching the workflow for single-cell, combinatorial, base-editor, prime-editor, or in vivo variants.
development
Transcribe DNA to RNA and translate to protein using Biopython, with NCBI codon-table selection, CDS validation, and six-frame ORF finding. Use when converting a CDS or ORF to its amino-acid sequence, selecting a non-standard (mitochondrial, bacterial, ciliate) genetic code, validating a coding sequence, or scanning all reading frames.