bio-crispr-screens-combinatorial-screens

What this skill does

## Version Compatibility

Reference examples tested with: MAGeCK 0.5.9+ (for MLE with interaction terms), Inzolia library annotation (Bayle 2024), pandas 2.2+, numpy 1.26+, scipy 1.12+, matplotlib 3.8+.

Before using code patterns, verify installed versions match. If versions differ:
- CLI: `mageck --version`; `mageck mle --help`
- For Cas12a libraries: verify against published Inzolia / in4mer / Big Papi annotations

If code throws ImportError, AttributeError, or TypeError, introspect the installed package and adapt the example to match the actual API rather than retrying.

## Combinatorial CRISPR Screen Analysis

**"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.

- CLI: `mageck mle` with explicit interaction terms for paired-Cas9 (Big Papi-style)
- Python: custom GI scoring for Cas12a multiplex (in4mer / Inzolia)
- Modality: enCas12a / LbCas12a single-array multiplex (preferred for paralog screens)

## Combinatorial Architecture Decision Tree

| Goal | Architecture | Library | Why |
|------|--------------|---------|-----|
| Paralog buffering, identify synthetic lethal paralog pairs | enCas12a single-array 4-guide multiplex | Inzolia (Bayle 2024) | Cas9 single-KO misses ~half of paralog-buffered essentials |
| Test specific pathway pair (e.g., DNA repair branches) | Big Papi (paired-Cas9 dual sgRNA cassette) | Custom | Mature methodology; SpCas9 well-characterized |
| 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:**
- Big Papi without proper paired-guide cloning: silent fusion of two guides into one cassette gives single perturbation
- Cas12a screens analyzed as Cas9 screens: MAGeCK normalization fails because Cas12a has different cut profile
- in4mer 4-guide arrays without all-singleton controls: GI scoring requires single-gene baselines

## Cas9 vs Cas12a for Multiplex

| Property | Cas9 paired (Big Papi) | Cas12a multiplex (in4mer / Inzolia) |
|----------|------------------------|--------------------------------------|
| Multiplex capacity per cassette | 2 sgRNAs (paired) | 4 (in4mer); 2 (standard Cas12a) |
| sgRNA processing | Two separate U6 promoters | 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 | ~2,000 paired cassettes | ~250 4-guide cassettes (in4mer) |
| Validated libraries | Limited (mostly custom) | Inzolia: 18k 4-guide arrays for 4k 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. The 30% library-size reduction vs paired-Cas9 makes it more cost-effective for genome-scale.

## The Paralog Buffering Phenomenon

**Dede et al 2020 *Genome Biol* 21:262; Thompson et al 2021** demonstrated that approximately half of constitutively-expressed essential genes are never detected in Cas9 single-KO screens. 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; 79% reproduce in ≥2 lines, 58% in all 3. These pairs were not findable by single-gene Cas9 screens, requiring combinatorial methodology.

**Examples:**
- **MAPK1 (ERK2) + MAPK3 (ERK1):** ERK family redundancy in proliferation
- **PIK3CA + PIK3CB:** PI3K alpha/beta redundancy
- **AKT1 + AKT2:** AKT family redundancy
- **HSP90AA1 + HSP90AB1:** HSP90 alpha/beta redundancy
- **STAG1 + STAG2:** Cohesion complex paralogs

Each is buffered: loss of one is tolerated; loss of both is lethal.

## Genetic Interaction (GI) Scoring

**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).

```python
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['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:**
- GI z-score < -2: Synthetic lethal (double-KO more lethal than expected) -- candidate drug target combinations
- GI z-score > 2: Synthetic rescue (double-KO less lethal than expected) -- compensatory pathway / paradoxical hit
- GI z-score -1 to 1: No interaction; effects are additive

## Run Combinatorial Screen Analysis (MAGeCK MLE with Interaction Indicator)

**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).

```bash
# 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` | Addition