Skip to content
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
75 changes: 75 additions & 0 deletions src/locusview/repository.py
Original file line number Diff line number Diff line change
Expand Up @@ -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 ──────────────────────────────────────────────────────

Expand Down Expand Up @@ -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) ──────────────────────────────────────────────

Expand Down Expand Up @@ -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_*``
Expand Down
56 changes: 56 additions & 0 deletions src/locusview/templates/variant.html
Original file line number Diff line number Diff line change
@@ -0,0 +1,56 @@
<!DOCTYPE html>
<html lang="en">
<head>
<meta charset="utf-8" />
<meta name="viewport" content="width=device-width, initial-scale=1" />
<title>{{ rsid }} — locusview</title>
<style>
body { font-family: system-ui, sans-serif; max-width: 60rem; margin: 3rem auto; padding: 0 1rem; line-height: 1.5; }
.muted { color: #666; }
table { border-collapse: collapse; width: 100%; margin-top: 1rem; font-size: 0.92rem; }
th, td { text-align: left; padding: 0.35rem 0.6rem; border-bottom: 1px solid #eee; }
th { border-bottom: 2px solid #ccc; }
td.num { text-align: right; font-variant-numeric: tabular-nums; }
.wrap { overflow-x: auto; }
</style>
</head>
<body>
<p><a href="/">&larr; locusview</a></p>
<h1>{{ rsid }}</h1>
<p class="muted">chromosome {{ chrom }} &middot; eQTL reverse lookup &mdash; the genes this variant associates with</p>

<h2>eQTLs</h2>
<p class="muted">
Showing associations from {{ n_tissues }} of {{ total_tissues }} datasets (bounded for the MVP).
Effect <em>direction</em> (the sign of &beta;) is not interpretable &mdash; the source shards store
no effect allele (<a href="https://github.com/boxiangliulab/locusview/issues/18">issue&nbsp;#18</a>).
</p>

{% if rows %}
<div class="wrap">
<table>
<thead>
<tr><th>Tissue</th><th>Gene</th><th>Ensembl</th><th>Position</th><th>p-value</th><th>&beta;</th><th>SE</th></tr>
</thead>
<tbody>
{% for r in rows %}
<tr>
<td>{{ r.tissue }}</td>
<td><a href="/gene/{{ r.gene }}">{{ r.gene }}</a></td>
<td class="muted">{{ r.ensembl_id }}</td>
<td class="num">{{ r.position }}</td>
<td class="num">{{ '%.2e'|format(r.pvalue) if r.pvalue is not none else '' }}</td>
<td class="num">{{ r.beta if r.beta is not none else '' }}</td>
<td class="num">{{ r.se if r.se is not none else '' }}</td>
</tr>
{% endfor %}
</tbody>
</table>
</div>
{% else %}
<p>No eQTLs found for this variant in the queried datasets.</p>
{% endif %}

<p class="muted">locusview v{{ version }}</p>
</body>
</html>
60 changes: 59 additions & 1 deletion src/locusview/web.py
Original file line number Diff line number Diff line change
Expand Up @@ -23,6 +23,7 @@
from locusview.repository import (
Dataset,
FakeQtlRepository,
Gene,
LocuscompareRepository,
QtlRepository,
pymysql_connection_factory,
Expand All @@ -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
Expand Down Expand Up @@ -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)
Expand Down Expand Up @@ -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
68 changes: 68 additions & 0 deletions tests/test_repository.py
Original file line number Diff line number Diff line change
Expand Up @@ -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
48 changes: 46 additions & 2 deletions tests/test_web.py
Original file line number Diff line number Diff line change
Expand Up @@ -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"
Loading