Replogle atlas: query-specific perturbation nomination

Start with a response gene and obtain candidate perturbations for follow-up. This offline tutorial uses authentic archived sufficient statistics from Replogle et al. (2022). The reliability statistics retain information from the full original gene axis, while the effect matrix stores the queried gene coordinates. All candidates in this archived task are retained.

Install the accompanying wheel or run python -m pip install . from the package directory. No atlas download is required for this tutorial. Its computations recompute weights and nominations from those statistics. The fixed scientific parameters are five control directions and attenuation 0.9, as in manuscript v25.2. See scientific design for how the reliability calculation and gene readout are separated.

import numpy as np
import pandas as pd
from pathlib import Path
import pertrank as pr
from pertrank.datasets import load_example

data, provenance = load_example("replogle")
print(provenance["dataset"])
print("Candidates:", len(data["target_ids"]))
print("Weight genes:", int(data["weight_gene_count"]))
print("Stored effect columns:", len(data["gene_ids"]))
queries = data["selected_query_ids"].tolist()
print("Example queries:", queries)

Replogle et al. 2022 K562 genome-wide Perturb-seq
Candidates: 2000
Weight genes: 8202
Stored effect columns: 300
Example queries: ['ENSG00000004478', 'ENSG00000005249', 'ENSG00000016490']

Generate candidate lists

The example queries are the first three archived IDs in lexicographic order. They were chosen without consulting rankings or evaluation outcomes. The effect_gene_subset declaration preserves the distinction between the full reliability geometry and the displayed gene effects.

result = pr.from_statistics(
    data, queries, top_k=10, effect_gene_subset=True,
    weight_gene_count=int(data["weight_gene_count"]),
)
for method in ["Euclidean", "PertRank"]:
    np.testing.assert_allclose(
        result.weights[method].to_numpy(), data[("archived_weight_" if "archived_weight_" + method in data else "frozen_formula_weight_") + method],
        rtol=0, atol=1e-12,
    )
print("Both weight vectors match the saved implementation within 1e-12.")
query = queries[0]
shortlist = result.nominations.query(
    "method == 'PertRank' and query_id == @query and direction == 'up'"
)
display(shortlist[["rank", "target_id", "effect", "weight", "signed_score"]])

Both weight vectors match the saved implementation within 1e-12.

	rank	target_id	effect	weight	signed_score
20	1	ENSG00000185049	0.201208	0.271581	0.054644
21	2	ENSG00000232838	0.400243	0.132111	0.052876
22	3	ENSG00000196396	0.189119	0.259603	0.049096
23	4	ENSG00000120533	0.126978	0.370043	0.046987
24	5	ENSG00000155858	0.217506	0.176427	0.038374
25	6	ENSG00000141198	0.368378	0.082767	0.030490
26	7	ENSG00000156482	0.082666	0.367504	0.030380
27	8	ENSG00000151247	0.064372	0.398538	0.025655
28	9	ENSG00000014123	0.250634	0.091600	0.022958
29	10	ENSG00000099991	0.116327	0.195976	0.022797

Inspect how reliability changes the shortlist

The scalar changes each candidate's competitive strength while the readout uses its original gene-specific effect. Compare PertRank with Unweighted to inspect the complete weighting workflow, and with Euclidean to inspect the added control geometry. Candidate movement describes this archived list; accuracy requires a separately defined reference.

raw = result.nominations.query(
    "method == 'Unweighted' and query_id == @query and direction == 'up'"
)
changes = pd.DataFrame({
    "added_by_PertRank": pd.Series(sorted(set(shortlist.target_id) - set(raw.target_id))),
    "removed_from_raw": pd.Series(sorted(set(raw.target_id) - set(shortlist.target_id))),
})
display(changes)
display(result.weights[["n_input_cells", "Euclidean", "PertRank"]].head())

	added_by_PertRank	removed_from_raw
0	ENSG00000014123	ENSG00000021762
1	ENSG00000099991	ENSG00000087903
2	ENSG00000120533	ENSG00000102897
3	ENSG00000151247	ENSG00000136156
4	ENSG00000155858	ENSG00000136536
5	ENSG00000156482	ENSG00000141068
6	ENSG00000185049	ENSG00000148426
7	ENSG00000196396	ENSG00000184489

	n_input_cells	Euclidean	PertRank
target_id
ENSG00000189403	20.0	0.006677	0.003824
ENSG00000064961	20.0	0.000000	0.000000
ENSG00000111906	20.0	0.003330	0.000000
ENSG00000107614	20.0	0.093870	0.051334
ENSG00000111704	20.0	0.013092	0.024704

Export for review

Up and down are signed expression directions. Boundary-tie inclusion is the chance of inclusion under uniform tie breaking, not a biological confidence probability. The complete score table retains eligible candidates with zero weight as well as candidates outside the shortlist.

out = Path("outputs/replogle")
out.mkdir(parents=True, exist_ok=True)
result.nominations.to_csv(out / "nominations.csv", index=False)
result.weights.to_csv(out / "weights.csv", index=False)
result.all_candidate_scores.to_csv(out / "all_candidate_scores.csv", index=False)
print("Exported", len(result.nominations), "nomination rows to", out)

Exported 180 nomination rows to outputs/replogle

Use your own AnnData

pr.rank(adata, queries, input_scale="counts", full_gene_axis=True, strata_keys=["batch"], top_k=25) is the cell-input entry point. Its input must contain the full measured raw gene axis, selected perturbation cells and matched controls in one declared stratum. Supply the metadata columns shown in the README. Already log-normalized data use input_scale="log" with an explicitly declared analysis axis.

The executable 00_anndata_quickstart.py demonstrates the complete cell-input route on a synthetic software fixture. It is separate from these authentic archived atlas statistics. Consult input contracts before using a new screen, and the full manuscript archive for scientific evaluation. The provenance JSON packaged with each fixture records its source files, their hashes, units and selection rules.