CF · QJE Revise & Resubmit · Task 1

Rebuilding the value label on GWAS evidence

The referees call the DisGeNET value label circular — it's built from publication counts. We rebuild L/M/H on independent TIGA/GWAS evidence and rerun the main spec across every candidate definition. The behavioral herding result (column 2) is what has to survive.

Main results — disease-level regressions

The behavioral herding test is column 2, share_gene_after (genes explored per publication after 2000) — a measure independent of how value was labelled. Coefficients are on Medium and High diseases relative to Low; a negative sign means more valuable diseases are explored less broadly (herding). Table 1 reports it under the paper's controls and under each added GWAS-coverage control, with a Yes/No panel recording what enters each specification. Table 2 gives the two breakthrough outcomes.

Table 1. Herding (share_gene_after) on the value tier, by definition and control set

Notes. OLS estimates; standard errors clustered by disease class in parentheses. ^*** p<0.01, ^** p<0.05, ^* p<0.1. Dependent variable share_gene_after = distinct genes explored per post-2000 publication; a negative coefficient is herding. "—" = coefficient not identified: in the GWAS-covered subsample every disease is Medium or High, so there is no Low reference group. Controlling for # TIGA traits slightly strengthens the estimate; # GWAS-hit opportunities slightly weakens it (it is partly collinear with the value tier), but every full-sample estimate stays negative and significant — the herding is not a GWAS-coverage artifact.

Table 2. Breakthrough outcomes (TIGA-rebuilt), paper controls

Notes. Post2000_TIGA_H = 1 if the disease has ≥1 High-tier gene first discovered after 2000. Because the value tier (RHS) is fixed in the pre-2000 exploration window, this is a forward-looking breakthrough, not a mechanical restatement of the label. Delay_To_TIGA_H = years from 1980 to the first High-tier discovery (censored at 39). A positive Post2000 and a negative Delay both say higher-value diseases reach breakthroughs more, and sooner. Subsample (Variant B) High coefficients are similar and significant; Medium is again not identified. Stars as in Table 1.

How to read the rebuttal

Lead with Variant A (full sample, uncovered-as-Low). Present the metric rows — study count, p-value, GWAS composite, gene-level pleiotropy — as a robustness panel: the herding conclusion is invariant to which GWAS signal defines value. The correlation between the DisGeNET and TIGA labels is only ≈0.11, yet the behavioral result is the same — that independence is the point. The subsample (Variant B) stands as a deliberately conservative check whose only weakness is statistical power.

Appendix · Data pipeline & variables

Every row in the table above runs through the identical pipeline below — only the value metric and the zero-handling rule change between definitions. The crosswalk (step 3), which is the expensive part, is built once and reused. Source: src/run_definition_sweep.py.

Cleaning & build steps

1. Load the pair universe. 429,111 gene–disease pairs from disgenet_v1.dta (14,619 diseases × 17,181 genes), each carrying the DisGeNET score. Gene IDs mapped to HGNC symbols via the associations file so they can be joined to TIGA.
2. Load & normalize TIGA. tiga_gene-trait_stats.tsv read with "NA" treated as missing; EFO ids normalized (EFO_0004908 → EFO:0004908, upper-cased) so ontology ids match across files.
3. Build the CUI → GWAS-trait crosswalk (once). Each disease code (CUI) is linked to TIGA traits through disease_mappings.tsv, trying five bridges in order and keeping any that resolve to a trait TIGA actually studied: (i) EFO id, (ii) EFO label name, (iii) DisGeNET disease name, (iv) MONDO id, (v) ORDO id. Result: 1,388 / 14,619 diseases (9.5%) map to ≥1 GWAS-covered trait — genuine GWAS sparsity.
4. Score each pair. For a metric x, a pair's score is the max of x over the gene × the disease's mapped traits. Three states result: NaN if the disease is uncovered (unmeasurable), 0 if covered but this gene has no GWAS hit (a genuine Low), or a positive value. ~8,654 pairs score positive (6,233 for odds ratio).
5. Cut into L / M / H (zero-inflated). Because most pairs are 0, a naive percentile collapses to L/H with no Medium. Instead the positive pairs are percentile-ranked and split at the 60th / 90th percentiles into M / H; zeros are Low. The two variants differ only in uncovered diseases — A codes them Low (full sample); B drops them (covered subsample). Gene-level geneNtrait is ranked directly across all pairs at 60/90.
6. Assemble the regression sample. From gene_empirics.dta, keep disease_tag==1, marginal_disease!=1, and non-missing share_gene_after → 5,515 diseases.
7. Build disease-level value vars. Labels map to representative values (L=30, M=75, H=95); a disease's value tier is the max value among pairs discovered by the exploration cutoff (YearInitial ≤ 2000), re-cut at 60/90 — the pre-2000 RHS. Breakthrough outcomes: Post2000_TIGA_H = has ≥1 High gene first discovered after 2000 (forward-looking); Delay_To_TIGA_H = years from 1980 to its first High discovery.
8. Run the main spec. For each outcome: y ~ C(category) + count_pmid_after_fill + C(diseaseclass), OLS with standard errors clustered by diseaseclass — identical to the paper's specification.
9. Add GWAS-coverage controls (robustness). Re-run column 2 adding n_traits (# TIGA traits the disease maps to) and, separately, n_gwas_opp (# gene–disease pairs with any GWAS hit), to confirm the effect is not driven by some diseases being better covered by GWAS.

Variables included

★ = the recommended metric (median p-value) and the outcome that carries the rebuttal (share_gene_after). Every metric is drawn from GWAS evidence, not publication counts — which is what breaks the circularity the referees flagged.