diff --git a/src/locusview/repository.py b/src/locusview/repository.py index 43df1a1..b0653ad 100644 --- a/src/locusview/repository.py +++ b/src/locusview/repository.py @@ -80,6 +80,20 @@ def eqtls_for_gene( """Return eQTL associations for ``gene_id`` within the given datasets.""" ... + def variant_chrom(self, rs_id: int) -> int | None: + """Resolve a variant's chromosome (needed to use the ``(chrom, rs_id)`` index), or None.""" + ... + + def eqtls_for_variant( + self, chrom: int, rs_id: int, dataset_ids: Sequence[int] + ) -> list[EqtlAssociation]: + """Reverse lookup: every association for one variant across the given datasets.""" + ... + + def gene_by_id(self, gene_id: int) -> Gene | None: + """Resolve a :class:`Gene` from its integer id (reverse of :meth:`resolve_gene`).""" + ... + # ── Shared pure helpers ────────────────────────────────────────────────────── @@ -191,6 +205,25 @@ def eqtls_for_gene( hits = [a for a in self._associations if a.gene_id == gene_id and a.dataset_id in wanted] return hits[:limit] + def variant_chrom(self, rs_id: int) -> int | None: + for a in self._associations: + if a.rs_id == rs_id: + return a.chrom + return None + + def eqtls_for_variant( + self, chrom: int, rs_id: int, dataset_ids: Sequence[int] + ) -> list[EqtlAssociation]: + wanted = set(dataset_ids) + return [ + a + for a in self._associations + if a.chrom == chrom and a.rs_id == rs_id and a.dataset_id in wanted + ] + + def gene_by_id(self, gene_id: int) -> Gene | None: + return next((g for g in self._genes if g.gene_id == gene_id), None) + # ── Real (locuscompare2 MySQL) ────────────────────────────────────────────── @@ -253,6 +286,48 @@ def eqtls_for_gene( finally: conn.close() + def variant_chrom(self, rs_id: int) -> int | None: + # The shards have no standalone rs_id index, but the 1000G LD tables are indexed on SNP_A. + # Scan autosomes (the eQTL shards are autosomal) for the first chromosome carrying the rsID. + snp = f"rs{rs_id}" + for chrom in range(1, 23): + if self._query( + f"SELECT 1 FROM tkg_p3v5a_ld_chr{chrom}_EUR WHERE SNP_A = %s LIMIT 1", (snp,) + ): + return chrom + return None + + def eqtls_for_variant( + self, chrom: int, rs_id: int, dataset_ids: Sequence[int] + ) -> list[EqtlAssociation]: + if not dataset_ids: + return [] + conn = self._connect() # one connection for the whole fan-out + try: + out: list[EqtlAssociation] = [] + for dataset_id in dataset_ids: + table = _shard_table(dataset_id) + # (chrom, rs_id) hits idx_eqtl_snp_join — a precise seek, ~100 rows not a full scan. + rows = self._run( + conn, + f"SELECT gene_id, rs_id, chrom, position, pvalue, beta, se " + f"FROM {table} WHERE chrom = %s AND rs_id = %s", + (chrom, rs_id), + ) + out.extend(_row_to_eqtl(dataset_id, r) for r in rows) + return out + finally: + conn.close() + + def gene_by_id(self, gene_id: int) -> Gene | None: + ensg = f"ENSG{gene_id:011d}" # 141510 -> "ENSG00000141510" + rows = self._query( + f"SELECT gene_name, gene_id, chr, start, end, strand FROM {_GENCODE_TABLE} " + f"WHERE gene_id LIKE %s LIMIT 1", + (f"{ensg}%",), + ) + return _row_to_gene(rows[0]) if rows else None + def pymysql_connection_factory() -> ConnectionFactory: """Build a factory that opens a read-only pymysql connection from ``LOCUSCOMPARE2_DB_*`` diff --git a/src/locusview/templates/variant.html b/src/locusview/templates/variant.html new file mode 100644 index 0000000..a165260 --- /dev/null +++ b/src/locusview/templates/variant.html @@ -0,0 +1,56 @@ + + + + + + {{ rsid }} — locusview + + + +

← locusview

+

{{ rsid }}

+

chromosome {{ chrom }} · eQTL reverse lookup — the genes this variant associates with

+ +

eQTLs

+

+ Showing associations from {{ n_tissues }} of {{ total_tissues }} datasets (bounded for the MVP). + Effect direction (the sign of β) is not interpretable — the source shards store + no effect allele (issue #18). +

+ + {% if rows %} +
+ + + + + + {% for r in rows %} + + + + + + + + + + {% endfor %} + +
TissueGeneEnsemblPositionp-valueβSE
{{ r.tissue }}{{ r.gene }}{{ r.ensembl_id }}{{ r.position }}{{ '%.2e'|format(r.pvalue) if r.pvalue is not none else '' }}{{ r.beta if r.beta is not none else '' }}{{ r.se if r.se is not none else '' }}
+
+ {% else %} +

No eQTLs found for this variant in the queried datasets.

+ {% endif %} + +

locusview v{{ version }}

+ + diff --git a/src/locusview/web.py b/src/locusview/web.py index 7c47226..1a11263 100644 --- a/src/locusview/web.py +++ b/src/locusview/web.py @@ -23,6 +23,7 @@ from locusview.repository import ( Dataset, FakeQtlRepository, + Gene, LocuscompareRepository, QtlRepository, pymysql_connection_factory, @@ -39,6 +40,12 @@ _MAX_TISSUES = 25 +def _rs_to_int(rsid: str) -> int | None: + """``"rs12345"`` -> ``12345``; anything else -> ``None``.""" + s = rsid.strip().lower() + return int(s[2:]) if s.startswith("rs") and s[2:].isdigit() else None + + def _default_repository() -> QtlRepository: """Pick a repository from config: the real DB if configured, else an empty fake.""" if get_db_settings().host: # pragma: no cover - requires DB config + network @@ -73,11 +80,13 @@ def search(q: str) -> Response: return RedirectResponse(f"/gene/{parsed.gene_symbol}", status_code=303) if parsed.kind is QueryKind.ENSEMBL_GENE and parsed.ensembl_id: return RedirectResponse(f"/gene/{parsed.ensembl_id}", status_code=303) + if parsed.kind is QueryKind.RSID and parsed.rsid: + return RedirectResponse(f"/variant/{parsed.rsid}", status_code=303) return _render( "not_found.html", status_code=404, query=q, - reason="Only gene search (symbol or Ensembl id) is available so far.", + reason="Search supports genes (symbol or Ensembl id) and variants (rsID) so far.", ) @app.get("/gene/{name}", response_class=HTMLResponse) @@ -115,4 +124,53 @@ def gene_page(name: str) -> HTMLResponse: total_tissues=len(all_datasets), ) + @app.get("/variant/{rsid}", response_class=HTMLResponse) + def variant_page(rsid: str) -> HTMLResponse: + rs_int = _rs_to_int(rsid) + if rs_int is None: + return _render( + "not_found.html", + status_code=404, + query=rsid, + reason="Not a valid rsID (expected e.g. rs12345).", + ) + chrom = repo.variant_chrom(rs_int) + if chrom is None: + return _render( + "not_found.html", + status_code=404, + query=f"rs{rs_int}", + reason="Could not locate this variant on an autosome in the reference panel.", + ) + all_datasets = repo.datasets() + datasets = all_datasets[:_MAX_TISSUES] + by_id: dict[int, Dataset] = {d.id: d for d in datasets} + hits = repo.eqtls_for_variant(chrom, rs_int, [d.id for d in datasets]) + hits = sorted(hits, key=lambda a: (a.pvalue is None, a.pvalue or 0.0)) # significant first + genes: dict[int, Gene | None] = {} + rows: list[dict[str, object]] = [] + for a in hits: + if a.gene_id not in genes: + genes[a.gene_id] = repo.gene_by_id(a.gene_id) + gene = genes[a.gene_id] + rows.append( + { + "tissue": by_id[a.dataset_id].tissue, + "gene": gene.symbol if gene else str(a.gene_id), + "ensembl_id": gene.ensembl_id if gene else "", + "position": a.position, + "pvalue": a.pvalue, + "beta": a.beta, + "se": a.se, + } + ) + return _render( + "variant.html", + rsid=f"rs{rs_int}", + chrom=chrom, + rows=rows, + n_tissues=len(datasets), + total_tissues=len(all_datasets), + ) + return app diff --git a/tests/test_repository.py b/tests/test_repository.py index cd9c711..a259fe1 100644 --- a/tests/test_repository.py +++ b/tests/test_repository.py @@ -212,3 +212,71 @@ def test_locuscompare_resolve_gene_by_ensembl_uses_like() -> None: def test_locuscompare_resolve_gene_not_found() -> None: factory, _ = _factory([]) assert LocuscompareRepository(factory).resolve_gene("NOPE") is None + + +# ── variant reverse-lookup ────────────────────────────────────────────────── + + +def test_fake_variant_chrom() -> None: + repo = FakeQtlRepository(associations=[EqtlAssociation(1, 7, 111, 17, 100, 0.01, 0.2, 0.05)]) + assert repo.variant_chrom(111) == 17 + assert repo.variant_chrom(999) is None + + +def test_fake_eqtls_for_variant() -> None: + repo = FakeQtlRepository( + associations=[ + EqtlAssociation(1, 7, 111, 17, 100, 0.01, 0.2, 0.05), + EqtlAssociation(2, 8, 111, 17, 100, 0.02, 0.1, 0.05), + EqtlAssociation(1, 9, 222, 17, 200, 0.5, 0.3, 0.05), # different variant + EqtlAssociation(3, 7, 111, 17, 100, 0.03, 0.1, 0.05), # dataset not requested + ] + ) + hits = repo.eqtls_for_variant(17, 111, [1, 2]) + assert sorted(a.gene_id for a in hits) == [7, 8] + + +def test_fake_gene_by_id() -> None: + gene = Gene(141510, "TP53", "ENSG00000141510.16", "17", 1, 2, "-") + repo = FakeQtlRepository(genes=[gene]) + assert repo.gene_by_id(141510) == gene + assert repo.gene_by_id(999) is None + + +def test_locuscompare_variant_chrom_finds_chromosome() -> None: + factory, log = _factory([(1,)]) # fake returns a hit on the first LD table queried + assert LocuscompareRepository(factory).variant_chrom(62062621) == 1 + assert "tkg_p3v5a_ld_chr1_EUR" in log[0][0] + assert log[0][1] == ("rs62062621",) + + +def test_locuscompare_variant_chrom_none_scans_all_autosomes() -> None: + factory, log = _factory([]) # no hit anywhere + assert LocuscompareRepository(factory).variant_chrom(1) is None + assert len(log) == 22 # chr1..chr22 + + +def test_locuscompare_eqtls_for_variant_uses_chrom_and_rsid() -> None: + factory, log = _factory([(141510, 62062621, 17, 7757304, "3e-4", "0.03", "0.05")]) + hits = LocuscompareRepository(factory).eqtls_for_variant(17, 62062621, [1, 2]) + assert [a.dataset_id for a in hits] == [1, 2] + assert all(params == (17, 62062621) for _, params in log) + assert "chrom = %s AND rs_id = %s" in log[0][0] + + +def test_locuscompare_eqtls_for_variant_empty_datasets() -> None: + factory, log = _factory([]) + assert LocuscompareRepository(factory).eqtls_for_variant(1, 1, []) == [] + assert log == [] + + +def test_locuscompare_gene_by_id_reconstructs_ensg() -> None: + factory, log = _factory([_GENCODE]) + gene = LocuscompareRepository(factory).gene_by_id(141510) + assert gene is not None and gene.symbol == "TP53" + assert log[0][1] == ("ENSG00000141510%",) # 141510 -> zero-padded ENSG + + +def test_locuscompare_gene_by_id_not_found() -> None: + factory, _ = _factory([]) + assert LocuscompareRepository(factory).gene_by_id(1) is None diff --git a/tests/test_web.py b/tests/test_web.py index 61df1f3..ad44829 100644 --- a/tests/test_web.py +++ b/tests/test_web.py @@ -66,6 +66,50 @@ def test_search_redirects_an_ensembl_id() -> None: def test_search_unsupported_query_is_404() -> None: - response = _gene_client().get("/search", params={"q": "rs12345"}, follow_redirects=False) + response = _gene_client().get("/search", params={"q": "1:1000-2000"}, follow_redirects=False) assert response.status_code == 404 - assert "gene search" in response.text.lower() + assert "supports" in response.text.lower() + + +# ── variant page (reverse lookup) ─────────────────────────────────────────── + + +def _variant_client() -> TestClient: + repo = FakeQtlRepository( + datasets=[Dataset(1, "Whole_Blood", "gtex-v8"), Dataset(2, "Liver", "gtex-v8")], + genes=[Gene(141510, "TP53", "ENSG00000141510.16", "17", 7661779, 7687550, "-")], + associations=[ + EqtlAssociation( + 1, 141510, 62062621, 17, 7757304, 3e-4, 0.03, 0.05 + ), # TP53 / Whole_Blood + EqtlAssociation(2, 141510, 62062621, 17, 7757304, 1e-3, 0.02, 0.05), # TP53 / Liver + EqtlAssociation( + 1, 99999, 62062621, 17, 7757304, 0.01, 0.10, 0.05 + ), # gene not in catalog + ], + ) + return TestClient(create_app(repository=repo)) + + +def test_variant_page_reverse_lookup() -> None: + response = _variant_client().get("/variant/rs62062621") + assert response.status_code == 200 + assert "rs62062621" in response.text + assert "TP53" in response.text # gene resolved from its integer id + assert "Whole_Blood" in response.text and "Liver" in response.text # both tissues + assert "99999" in response.text # unresolved gene falls back to its id + + +def test_variant_page_unknown_variant_is_404() -> None: + response = _variant_client().get("/variant/rs999999") # no association -> chrom unresolved + assert response.status_code == 404 + + +def test_variant_page_bad_rsid_is_404() -> None: + assert _variant_client().get("/variant/notarsid").status_code == 404 + + +def test_search_redirects_an_rsid() -> None: + response = _variant_client().get("/search", params={"q": "rs62062621"}, follow_redirects=False) + assert response.status_code == 303 + assert response.headers["location"] == "/variant/rs62062621"