bio-crispr-screens-jacks-analysis

What this skill does

## Version Compatibility

Reference examples tested with: JACKS 0.2.0+ (felicityallen/JACKS), pandas 2.2+, numpy 1.26+, scipy 1.12+, matplotlib 3.8+.

Before using code patterns, verify installed versions match. If versions differ:
- CLI: `python run_JACKS.py --help` (run_JACKS.py at the JACKS repo root after clone)
- Python: from jacks.jacks_io import runJACKS; help(runJACKS)
- GitHub: install via `git clone https://github.com/felicityallen/JACKS && cd JACKS && pip install .`

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

## JACKS CRISPR Screen Analysis

**"Analyze CRISPR screens with guide-level efficacy modeling"** -> Jointly model per-sgRNA log-fold-change across one or more screens as the product of gene essentiality and guide efficacy, sharing efficacy across screens with the same library so that low-quality guides are down-weighted automatically.

- CLI: `python run_JACKS.py countfile replicatefile guidemappingfile [options]` (script at JACKS repo root)
- Python: `from jacks.jacks_io import runJACKS` for programmatic use; lower-level `from jacks.infer import inferJACKS`
- Output: per-gene essentiality (`X1`), per-sgRNA efficacy (`X1`), log-likelihood ratio per gene

## The JACKS Model (under the hood)

**Why this matters for postdoc-level use:** JACKS decomposes the observed per-sgRNA log-fold-change as:

```
LFC[i, c] = gene_effect[g(i), c] * guide_efficacy[i] + noise
```

where `i` is sgRNA index, `c` is screen condition, `g(i)` is the gene targeted by sgRNA i. Gene effect varies by condition (different cell lines, different treatments) but guide efficacy is intrinsic to the sgRNA sequence and is treated as constant across screens. The model fits both parameters via variational Bayes with hierarchical priors:

- `guide_efficacy[i] ~ Beta(alpha, beta)` shared across all sgRNAs (hyperparameters fit empirically)
- `gene_effect[g, c] ~ Normal(0, sigma_c^2)` per condition

The variational posterior gives expected guide efficacy and gene effect; log-likelihood-ratio tests against a null (zero gene effect) provide gene-level significance.

**Critical assumption:** Guide efficacy is treated as cell-line independent within the same chemistry. Allen 2019 reports per-sgRNA Cas9 KO efficacy correlates ~0.7 across cell lines (within-chemistry), supporting library-shared efficacy. **However**, efficacy is NOT shareable across chemistries: Cas9 KO efficacy != CRISPRi knockdown efficiency != CRISPRa activation efficiency. JACKS must be run separately per chemistry; use only within the same chemistry on the same library.

## When JACKS Outperforms MAGeCK and BAGEL2

| Scenario | Advantage | Quantified gain (Allen 2019) |
|----------|-----------|------------------------------|
| Multi-screen joint analysis (>=3 screens with same library) | Efficacy shared; noise averaged | 12-21% lower error vs MAGeCK; 9% vs BAGEL2; 97-99% of cell lines improved |
| Reusing public reference screens (DepMap, Sanger Score) as efficacy prior | Transfer learning | New screens can be smaller; 2.5x sample-size reduction |
| Libraries with broad efficacy variance (e.g. older GeCKOv2) | Down-weights known weak guides | Larger gain than on Brunello (already efficacy-filtered) |
| Heterogeneous quality (mixed plasmid quality across screens) | Per-screen noise estimation | Cleaner per-condition gene effects |

## When JACKS Is Not the Right Tool

- **Single screen, no prior efficacy:** JACKS has nothing to leverage; MAGeCK or BAGEL2 work as well.
- **Single timepoint / two-condition essentiality:** RRA or BAGEL2 simpler and equivalent.
- **Heavy-selection drug screens:** drugZ explicit for chemogenomic; JACKS less sensitive.
- **Cancer-cell-line copy-number screens:** Chronos preferred; jointly models CN bias + screen quality; JACKS does neither.
- **Cross-chemistry sharing (e.g. CRISPRi + Cas9):** Efficacy is chemistry-specific; do not share.

## Run JACKS Joint Analysis

**Goal:** Jointly analyze multiple CRISPR screens performed with the same library and chemistry.

**Approach:** Provide a count matrix with all samples across all screens, a replicate map identifying which samples belong to which screen and condition, and a sgRNA-to-gene map. JACKS learns guide efficacy shared across screens and gene effects per screen.

```python
# Programmatic invocation
from jacks.jacks_io import runJACKS

# Input file paths
counts_path = 'counts.txt'                    # rows=sgRNA; first cols 'sgRNA' (or custom), then sample counts
replicate_map_path = 'replicatemap.txt'       # tab-separated with header: Replicate, Sample, Control
guide_map_path = 'guidemap.txt'               # tab-separated with header: sgRNA, Gene

# Replicate map format (tab-separated WITH header; column names match flags below)
# Replicate                Sample          Control
# Screen1_T1               Screen1_T       Screen1_C
# Screen1_T2               Screen1_T       Screen1_C
# Screen1_C1               Screen1_C       Screen1_C
# Screen2_T1               Screen2_T       Screen2_C
# Screen2_T2               Screen2_T       Screen2_C
# Screen2_C1               Screen2_C       Screen2_C

runJACKS(
    countfile=counts_path,
    replicatefile=replicate_map_path,
    guidemappingfile=guide_map_path,
    rep_hdr='Replicate',
    sample_hdr='Sample',
    ctrl_sample_hdr='Control',                # per-sample control specification
    sgrna_hdr='sgRNA',
    gene_hdr='Gene',
    outprefix='jacks_out',
    apply_w_hp=True,                          # hierarchical prior on efficacy
)
```

```bash
# Equivalent CLI run (run_JACKS.py is at the JACKS repo root after clone)
python run_JACKS.py \
    counts.txt \
    replicatemap.txt \
    guidemap.txt \
    --rep_hdr Replicate \
    --sample_hdr Sample \
    --ctrl_sample_hdr Control \              # per-sample control (or --common_ctrl_sample <name>)
    --sgrna_hdr sgRNA \
    --gene_hdr Gene \
    --outprefix jacks_out \
    --apply_w_hp                              # hierarchical prior on efficacy
# Outputs:
#   jacks_out_gene_JACKS_results.txt    gene-level: X1 (effect), X2 (std), p_neg, p_pos
#   jacks_out_grna_JACKS_results.txt    sgRNA-level: X1 (efficacy 0-1), X2 (std)
#   jacks_out_JACKS_results_full.pickle  full posterior for downstream
```

## Output Interpretation

| Column | Meaning | Direction |
|--------|---------|-----------|
| `X1` (gene file) | Posterior mean of gene_effect | Negative = essential (depleted); positive = enriched |
| `X2` (gene file) | Posterior std of gene_effect | Lower = more confident; combine with X1/X2 ratio for z-like scoring |
| `p_neg` / `p_pos` | One-sided posterior probability | Hits at p_neg <0.05 or p_pos <0.05 |
| `fdr_log10` | log10(FDR) of gene effect | Use <-1 for FDR<0.1; <-2 for FDR<0.01 |
| `X1` (sgRNA file) | Posterior mean of guide efficacy | 0 = ineffective, 1 = maximum efficacy |
| `X2` (sgRNA file) | Posterior std of efficacy | Confidence in efficacy estimate |

**Interpretation rule:** A gene is essential if X1 < 0 AND fdr_log10 < -1. The X1/X2 ratio gives a z-like statistic; |X1/X2| > 2 corresponds to ~95% credible deviation from zero. Sort by X1 (most negative first) for essentiality rank.

## Build Library-Wide Efficacy Prior from Reference Screens

**Goal:** Transfer learned efficacy from a large public screen panel to a new small screen.

**Approach:** Run JACKS on the reference panel (e.g. DepMap CRISPR screens with TKOv3 or Brunello), extract per-sgRNA efficacy posterior, and supply it as the prior for a new screen.

```python
def extract_efficacy_prior(reference_jacks_results):
    '''Build per-sgRNA efficacy prior (mean + std) from a large reference screen.'''
    df = pd.read_csv(reference_jacks_results, sep='\t')
    prior = df[['sgRNA', 'X1', 'X2']].rename(columns={'X1': 'eff_mean', 'X2': 'eff_std'})
    return prior

# Use in new JACKS run via --reffile <path>
# Reference: Allen 2019 Genome Research