Disproportionality Analysis (PRR, ROR, IC, EBGM)
The statistical machinery of quantitative signal detection - the frequentist proportional reporting ratio (PRR) and reporting odds ratio (ROR), and the Bayesian shrinkage estimators from the BCPNN (information component, IC) and the MGPS (empirical Bayes geometric mean, EBGM) - whose observed and expected counts come from a 2x2 contingency table of report counts (EBGM's empirical prior is fit once across the full drug-event matrix, not from one pair's table alone) and which flag drug-event pairs reported more often than a specified expected baseline in a spontaneous-report database.
On this page
When a safety database shows a drug and a side effect appearing together a lot, analysts compute a small set of ratios to ask "more than expected, by how much, and how sure are we?" PRR and ROR are simple, hand-calculable ratios built straight from report counts. IC and EBGM are Bayesian versions that deliberately pull the ratio back toward "nothing unusual" when the count behind it is tiny - so two case reports do not produce a wild, unstable-looking signal. PRR, ROR, and IC's own-pair inputs all come from the same four-cell counting table; EBGM's do too, but EBGM also needs a "typical shrinkage" pattern learned once from every other drug-event pair in the database, so it cannot be computed from one pair's table alone. All four answer a reporting-pattern question, never a question about how many patients were actually harmed.
This entry is the formula-level companion to `signal-detection`: it works through exactly how PRR, ROR, and IC are computed from the same 2x2 table (and how EBGM's pair-specific observed/expected inputs come from that same table even though its shrinkage prior does not), why the four numbers are correlated but not identical, and precisely what each does and does not measure. Read `signal-detection` first for the screening workflow and regulatory context; read this entry for the arithmetic, the confidence and credibility intervals, and the shrinkage mechanics that make Bayesian disproportionality behave differently from its frequentist cousins at small counts.
The 2x2 contingency frame
PRR, ROR, and the unshrunk information component are each a function of one table, collapsed across the entire spontaneous-reporting database for one drug-event pair: `a` = reports naming both the drug and the event; `b` = reports naming the drug with any other event; `c` = reports naming the event with any other drug; `d` = all remaining reports. Row and column margins matter as much as `a` itself: `D = a+b` (all reports on the drug), `Etot = a+c` (all reports on the event), and `N = a+b+c+d` (all reports in the database). Nothing in this table is a person or a person-time unit - it is counts of reports - so every statistic below is a reporting-proportion contrast, never an incidence, a rate, or a risk.
PRR and ROR: the frequentist ratios
The proportional reporting ratio contrasts the drug's own event-reporting proportion against the proportion in all other drugs combined, deliberately excluding `a` and `b` from the comparator: PRR = [a/(a+b)] / [c/(c+d)]. The reporting odds ratio is the cross-product odds ratio on the same table: ROR = (a*d)/(b*c), with SE(ln ROR) = sqrt(1/a+1/b+1/c+1/d) giving the log-scale 95% CI = exp(ln ROR +/- 1.96*SE). PRR and ROR are nearly identical when the event's reporting proportion is small in both the drug's own reports and the comparator's (so a/(a+b) approximately a/b and c/(c+d) approximately c/d) - "rare in the database overall" is not sufficient on its own, since a common event concentrated within one drug's reports would not satisfy this - and diverge as the event becomes more common in either group, because the odds ratio does not share the ratio-of-proportions' ceiling effect; EudraVigilance and the Netherlands Pharmacovigilance Centre prefer ROR for this reason, while the MHRA's PRR>=2 & chi-square>=4 & n>=3 rule remains the most widely taught frequentist gate. Both are undefined or explosive when a cell is zero or near-zero; a 0.5 continuity correction on all four cells yields a finite estimate but can materially shift the point estimate and interval coverage in sparse or imbalanced tables - the direction of the shift depends on the specific cell configuration, not uniformly toward the null - so report both the corrected and uncorrected cells whenever a is small.
Bayesian shrinkage: BCPNN information component and MGPS EBGM
The Bayesian family starts from a different baseline than PRR: the expected count under total independence across the whole matrix, E[a] = D*Etot/N, which folds `a` itself back into both margins (unlike PRR's "all other drugs" comparator). The relative reporting ratio RRR = a/E[a] is the observed/expected ratio both the BCPNN and the MGPS shrink. The information component is IC = log2(RRR) = log2(a/E[a]); a positive IC means the pair is reported more than the whole-database independence model expects, and a common WHO-UMC screening convention is IC025 (the lower bound of the 95% credibility interval) > 0. The unshrunk IC computed directly from a and E[a] is unstable at small a for exactly the same reason PRR is - it is a ratio of small counts. The BCPNN's actual innovation is a weakly informative Bayesian prior on the underlying cell probabilities (Bate et al. 1998's original formulation used Beta priors with a normal approximation to the posterior; Noren et al.'s 2006 extension uses a Dirichlet-multinomial posterior evaluated by Monte Carlo) that pulls the posterior IC toward zero in proportion to how little data support it. A simple, commonly taught way to see the shrinkage mechanic by hand is to add a small constant gamma to both the observed and expected counts before taking the log: IC_shrunk = log2[(a+gamma)/(E[a]+gamma)], gamma often taken as 0.5. This approximation captures the qualitative behavior correctly (see the worked example) but production BCPNN implementations (the PhViD R package, WHO-UMC's VigiBase pipeline) compute the full posterior rather than this shortcut.
The MGPS is a fully empirical Bayes method: DuMouchel's 1999 Gamma-Poisson Shrinker fits a two-component Gamma mixture prior to the entire matrix of relative-reporting-ratios across every drug-event pair in the database (by maximum likelihood / EM), then computes EBGM for each pair as the geometric mean of the posterior distribution of the relative-reporting-ratio parameter lambda under that fitted mixture - EBGM = exp(E[ln(lambda) | data]), not the arithmetic posterior mean - with EB05 as the lower 5th percentile of that same posterior. Because the prior's two Gamma components are estimated from the database itself, EBGM cannot be hand-computed from a single 2x2 table the way PRR, ROR, or the simplified IC can - it requires the full-matrix fit. This is precisely why shrinkage matters most for sparse cells: in a database screening millions of drug-event pairs, the great majority of cells have a=0, 1, or 2, and an unshrunk RRR or IC computed from two case reports against a tiny expected count can be enormous and almost entirely noise. FDA materials commonly apply an EB05>=2 screening convention (occasionally stated as EB05>2), and WHO-UMC commonly applies IC025>0; neither is a formal multiple-testing correction, and crossing either boundary produces a statistical alert for review, not a validated safety signal on its own. Both trade some interpretability and computational simplicity for shrinkage that measurably reduces sparse-cell instability under the massive multiplicity of a full database scan, though how IC and EBGM compare to each other in calibration and false-discovery behavior is dataset- and reference-standard-dependent rather than settled by a single universal ranking.
Stratification and subgrouping
Every formula above is computed on the whole-database margins by default, which silently pools heterogeneous reporting sources (spontaneous, solicited, literature), regions, and time periods. Mantel-Haenszel-stratified PRR/ROR, or a logistic-regression-adjusted ROR with drug, age, sex, reporting year, and report-source as covariates, recovers a within-stratum disproportionality estimate and is the standard remedy when pooling is suspected of confounding the reporting association (masking, notoriety spikes concentrated in one region or one report type). None of the shrinkage estimators substitute for stratification - shrinkage stabilizes a sparse cell; it does not remove confounding by reporting source.
What disproportionality is NOT
None of PRR, ROR, IC, or EBGM has a person-time or person-count denominator; none can be converted into an incidence rate, an absolute risk, or a risk ratio without an entirely separate denominator-based study. A drug used by ten million patients and a drug used by ten thousand can show an identical PRR if their reporting proportions happen to match, even though their absolute event counts and true incidence differ by three orders of magnitude. Masking (a very large, dominant event or drug-class signal inflating a comparator margin so a real smaller signal is statistically hidden) and competition bias can distort every measure here, frequentist or Bayesian alike, because they distort the margins the ratio is built from - but the size and sometimes the direction of the distortion depends on which margin is inflated and which comparator or shrinkage the affected measure uses, so PRR, ROR, IC, and EBGM are not guaranteed to be biased by the same amount. Shrinkage does not fix masking either way; it only tempers instability from small counts. See `spontaneous-reporting-biases-rwe` for the full treatment of masking, competition bias, and notoriety/stimulated reporting.
Pros, cons, and trade-offs
- PRR/ROR vs Bayesian shrinkage (IC/EBGM): PRR and ROR are transparent, hand-computable, and the basis of explicit regulator-facing rules (MHRA's PRR>=2 & chi-square>=4 & n>=3). They over-fire badly at small a. IC and EBGM trade a fully closed-form calculation for shrinkage-stabilized behavior under massive multiplicity (the magnitude of the stabilization is dataset-dependent, not a fixed number), at the cost of needing either a fitted prior (EBGM, whole-database EM) or a posterior approximation (IC, per-pair but still Bayesian). Prefer PRR/ROR for a small, auditable, single-product screen; prefer IC/EBGM for routine data-mining across the full MedDRA grid.
- BCPNN IC vs MGPS EBGM: IC is computed per drug-event pair with a comparatively simple, fast-to-evaluate prior and is the WHO-UMC/VigiBase standard; EBGM requires an EM fit of a two-component Gamma mixture across the entire database and is the FDA/Empirica standard. EBGM's empirically-fitted, two-component prior is often reported to handle the extreme sparse tail (a=1, a=2) differently than IC's fixed weakly-informative prior - benchmark studies do not uniformly favor one over the other, and comparative calibration is dataset- and reference-standard-dependent - at materially higher computational cost and the operational requirement of re-fitting as the database grows.
- Unadjusted vs stratified/regression-adjusted disproportionality: The unadjusted whole-database 2x2 is fast and standard for a first pass, but pools heterogeneous report sources, regions, and time periods, hiding confounding-by-report-source and masking. Mantel-Haenszel stratification or a logistic-regression-adjusted ROR costs analyst time and covariate data but is the only remedy that addresses confounded margins rather than merely sparse-cell noise.
When NOT to use - and when it is actively misleading
- As a stand-in for incidence, risk, or a rate ratio. None of these four statistics has a denominator of exposed persons or person-time; do not report a PRR or an EBGM as if it were a risk ratio, and do not back-calculate an incidence from any of them.
- At very low counts without a shrinkage estimator. With a<3, PRR and ROR are numerically unstable and their confidence intervals are close to uninformative; this is precisely the regime EBGM and IC were built for, and reporting an unshrunk PRR of 40 from a=1 without disclosing the count is misleading regardless of how the number looks.
- Pooled across drug-event pairs of very different report volume as if the statistics were comparable in strength. An IC of 3.0 on a pair with a=5,000 is a far more stable finding than an IC of 3.0 on a pair with a=5; always report `a` and the interval width alongside the point estimate.
- Without stratifying for a plausible masking or notoriety pattern. A dominant class-level signal, a stimulated-reporting spike after a label change, or heterogeneous report sources pooled together can distort every measure here, though not necessarily by the same magnitude or direction across PRR/ROR and IC/EBGM since they are built from different comparators; see `spontaneous-reporting-biases-rwe` before trusting an unstratified result.
- As the terminal analysis for a regulatory or clinical decision. A signal of disproportionate reporting from any of these four statistics is an input to `signal-validation-prioritization-evaluation-rwe`, not a conclusion; confirmation requires a denominator-based design.
Data-source operational depth
- Spontaneous-reporting systems (FAERS, VigiBase, EudraVigilance): These formulas assume one row per deduplicated report-drug-event mention; failing to deduplicate before tabulating inflates `a` and, through it, every downstream statistic - though PRR/ROR and IC/EBGM propagate that inflation through different comparators, so the numeric impact is not identical across the four measures. Pre-specify the MedDRA level (Preferred Term vs Standardised MedDRA Query) before computing any of the four statistics, because grouping materially changes both `a` and `c`.
- Claims/EHR: These formulas are not the natural tool once a person-time denominator exists; running a disproportionality-style 2x2 over claims data for exploratory screening is sometimes done, but the resulting PRR/ROR/IC/EBGM must not be presented as equivalent to the SRS-derived versions - the underlying report-generation process (structured claims vs voluntary reports) differs completely, and sequential methods (`maxsprt-sequential-safety-surveillance-rwe`, TreeScan) are the denominator-aware alternative.
- Registry/literature/solicited reports: Different reporting dynamics than a spontaneous ICSR mean these sources should be stratified, not pooled, before computing any of the four statistics on a combined table.
Worked example (toy 2x2, computing PRR, ROR, and IC by hand)
A drug-event pair, Drug Y and QT-interval prolongation, is screened in a snapshot of 500,000 deduplicated reports. The collapsed table is a=25 (Drug Y and QT prolongation), b=4,975 (Drug Y with any other event), c=500 (any other drug and QT prolongation), d=494,500 (all remaining reports); D=a+b=5,000, Etot=a+c=525, N=500,000.
- PRR = (25/5,000) / (500/495,000) = 0.005 / 0.0010101 = 4.95.
- ROR = (25*494,500)/(4,975*500) = 12,362,500/2,487,500 = 4.97; SE(ln ROR) = sqrt(1/25+1/4,975+1/500+1/494,500) = 0.2054; 95% CI = exp(ln 4.97 +/- 1.96*0.2054) = (3.32, 7.43).
- Yates chi-square on this table = 71.4 (far exceeds 4).
- MHRA check: PRR=4.95>=2, chi-square=71.4>=4, a=25>=3 - all met, signal of disproportionate reporting by the frequentist rule.
- Expected count under whole-database independence: E[a] = D*Etot/N = 5,000*525/500,000 = 5.25 (note this differs from PRR's comparator, which never uses E[a]).
- RRR = a/E[a] = 25/5.25 = 4.76; unshrunk IC = log2(4.76) = 2.25.
- Shrunk IC with gamma=0.5: IC_shrunk = log2[(25+0.5)/(5.25+0.5)] = log2(4.43) = 2.15 - a modest pull toward zero because a=25 already carries reasonable information.
- Contrast with a sparse cell: if a=2 and E[a]=0.08 (a very rare background event), the unshrunk RRR=25 and IC=4.64 look dramatic off two case reports; the same gamma=0.5 shrinkage gives IC_shrunk= log2(2.5/0.58)=2.11 - cutting the apparent signal size by more than half and reflecting the real uncertainty in an estimate built from two reports. This is why EB05 and IC025, not raw PRR, are the operational gates for full-database mining.
Interpreting the output
PRR=4.95 and ROR=4.97 (95% CI 3.32-7.43) both mean the same qualitative thing on this table: QT prolongation is reported roughly five times more often with Drug Y than the comparator baseline, and all three MHRA criteria are met. IC=2.25 (unshrunk) says the same direction of effect on a different baseline (whole-database independence, not "all other drugs"), and a log2 scale of 2.25 corresponds to roughly 4.76-fold observed-over-expected. None of these numbers is an incidence, a risk, or a causal estimate - there is no exposed person-time and no unexposed comparator group. The correct downstream action for a = 25 with all three frequentist criteria met and a comfortably positive shrunk IC is clinical and epidemiologic triage under `signal-validation-prioritization-evaluation-rwe`, checking first for notoriety/stimulated reporting, duplicate cases, and masking, then - if the signal survives triage - designing a denominator-based confirmatory study.
Decision diagram
flowchart TD
Tab["2x2 table: a=25, b=4975, c=500, d=494500<br/>D=5000, Etot=525, N=500000"]
Tab --> PRR["PRR = a/(a+b) over c/(c+d)<br/>= 0.005 / 0.0010101 = 4.95"]
Tab --> ROR["ROR = a*d / b*c<br/>= 4.97, 95% CI 3.32-7.43"]
Tab --> Ea["E[a] = D*Etot/N = 5.25<br/>(whole-database independence baseline)"]
Ea --> IC["Unshrunk IC = log2(a/E[a])<br/>= log2(4.76) = 2.25"]
IC --> ICs["Shrunk IC (gamma=0.5)<br/>= log2(25.5/5.75) = 2.15"]
Ea --> EBGM["EBGM/EB05<br/>needs the FULL database's<br/>Gamma mixture prior, not just this table"]
PRR --> Rule1["MHRA: PRR>=2 AND chi2>=4 AND a>=3"]
ROR --> Rule1
ICs --> Rule2["WHO-UMC: IC025 > 0"]
EBGM --> Rule3["FDA: EB05 >= 2"]
Rule1 --> SDR{Signal of disproportionate reporting?}
Rule2 --> SDR
Rule3 --> SDR
SDR -->|yes| Next["signal-validation-prioritization-evaluation-rwe"]flowchart LR
subgraph PRRcomp["PRR / ROR comparator"]
direction TB
P1["numerator: a/(a+b)<br/>drug-specific reporting proportion"]
P2["denominator: c/(c+d)<br/>ALL OTHER DRUGS ONLY<br/>(a, b excluded)"]
end
subgraph ICcomp["IC / EBGM comparator"]
direction TB
I1["observed: a"]
I2["expected: E[a] = D*Etot/N<br/>WHOLE-DATABASE independence<br/>(a folded back into both margins)"]
end
PRRcomp -->|"4.95 on the worked example"| Out1["PRR"]
ICcomp -->|"4.76 on the worked example"| Out2["RRR feeding IC"]
Out2 --> Shrink["Bayesian prior pulls RRR/IC<br/>toward the null when a is small<br/>PRR/ROR have no such mechanism"]Worked example
Scenario
A pharmacovigilance statistician has a deduplicated 2x2 table for Drug Y and QT-interval prolongation from a 500,000-report database snapshot. She computes PRR, ROR with its 95% CI, the Yates chi-square, the MHRA signal flags, and both the unshrunk and shrunk (Bayesian-style) information component, to see how the frequentist and Bayesian measures agree and where they diverge.
Dataset
Collapsed 2x2 contingency table for Drug Y and QT-interval prolongation, 500,000 deduplicated reports.
| cell | label | count |
|---|---|---|
| a | Drug Y AND QT prolongation | 25 |
| b | Drug Y AND any OTHER event | 4975 |
| c | Any OTHER drug AND QT prolongation | 500 |
| d | Any OTHER drug AND any OTHER event | 494500 |
Steps
Result
PRR = 4.95; ROR = 4.97 (95% CI 3.32-7.43); Yates chi-square = 71.4; all three MHRA criteria met. Unshrunk IC = 2.25 and shrunk IC (gamma=0.5) = 2.15. All four numbers point the same direction - Drug Y and QT prolongation are reported roughly five times more than expected - but none of them is an incidence or a risk; the pair should move to signal validation and clinical review, not a causal conclusion.
Trade-offs
Runnable example
PRR, ROR with 95% CI, Yates chi-square, MHRA signal flag, and both the unshrunk and shrunk information component (configurable shrinkage constant gamma) computed from a deduplicated, MedDRA-coded report table (one row per report-drug-event mention, keyed by report_id so a resubmitted case cannot double-count a pair).
import pandas as pd
import numpy as np
from scipy.stats import chi2_contingency
def disproportionality(reports: pd.DataFrame,
drug_col: str = "drug_name",
event_col: str = "event_pt",
report_col: str = "report_id",
gamma: float = 0.5) -> pd.DataFrame:
# One row per (report, drug, event) mention, so a report counted twice (e.g. a
# follow-up submission) cannot inflate a, b, or c.
mentions = reports[[report_col, drug_col, event_col]].drop_duplicates()
N = reports[report_col].nunique() # a+b+c+d
n_drug = mentions.groupby(drug_col)[report_col].nunique() # a+b : reports naming the drug
n_event = mentions.groupby(event_col)[report_col].nunique() # a+c : reports naming the event
a_tab = mentions.groupby([drug_col, event_col])[report_col].nunique() # a : reports naming both
rows = []
for (drug, event), a in a_tab.items():
a = int(a)
b = int(n_drug[drug] - a) # drug, other events
c = int(n_event[event] - a) # event, other drugs
d = int(N - a - b - c) # all other reports
if min(a, b, c, d) < 0:
continue
# PRR and ROR: same 0.5 continuity correction on all four cells if any cell is
# zero, applied to BOTH ratios (correcting ROR only leaves PRR dividing by a
# zero comparator proportion when c=0).
pa, pb, pc, pd_ = (a, b, c, d) if min(a, b, c, d) > 0 else (a + .5, b + .5, c + .5, d + .5)
prr = (pa / (pa + pb)) / (pc / (pc + pd_))
ror = (pa * pd_) / (pb * pc)
se = np.sqrt(1/pa + 1/pb + 1/pc + 1/pd_)
ror_lo, ror_hi = np.exp(np.log(ror) - 1.96 * se), np.exp(np.log(ror) + 1.96 * se)
chi2 = chi2_contingency([[a, b], [c, d]], correction=True)[0] # Yates-corrected
# Whole-database independence baseline (different from PRR's comparator); uses the
# raw uncorrected cells, not the PRR/ROR continuity correction above.
expected_a = (a + b) * (a + c) / N # E[a]
rrr = a / expected_a if expected_a else np.nan
ic_unshrunk = np.log2(rrr) if rrr and rrr > 0 else np.nan
ic_shrunk = np.log2((a + gamma) / (expected_a + gamma)) if expected_a else np.nan
signal_freq = (prr >= 2) and (chi2 >= 4) and (a >= 3) # MHRA operational criteria
rows.append(dict(drug=drug, event=event, a=a, b=b, c=c, d=d,
PRR=prr, ROR=ror, ROR_lo=ror_lo, ROR_hi=ror_hi, chi2_yates=chi2,
E_a=expected_a, RRR=rrr, IC=ic_unshrunk, IC_shrunk=ic_shrunk,
MHRA_signal=signal_freq))
return pd.DataFrame(rows).sort_values("PRR", ascending=False)
def _worked_example_fixture() -> pd.DataFrame:
"""Report-level fixture reproducing a=25, b=4975, c=500, d=494500 for Drug Y / QT prolongation."""
pair_rows = ([("DrugY", "QT_prolongation")] * 25 + [("DrugY", "Other_Event_B")] * 4975 +
[("OtherDrug_C", "QT_prolongation")] * 500 + [("OtherDrug_D", "Other_Event_D")] * 494500)
return pd.DataFrame({"report_id": range(1, len(pair_rows) + 1),
"drug_name": [r[0] for r in pair_rows],
"event_pt": [r[1] for r in pair_rows]})
# Worked-example check (run this file directly): reproduces a=25, b=4975, c=500, d=494500,
# PRR=4.95, ROR=4.97, IC=2.25, IC_shrunk=2.15.
if __name__ == "__main__":
out = disproportionality(_worked_example_fixture())
row = out[(out.drug == "DrugY") & (out.event == "QT_prolongation")].iloc[0]
assert (row.a, row.b, row.c, row.d) == (25, 4975, 500, 494500)
assert round(row.PRR, 2) == 4.95 and round(row.ROR, 2) == 4.97
assert round(row.IC, 2) == 2.25 and round(row.IC_shrunk, 2) == 2.15PRR, ROR with 95% CI, Yates chi-square, MHRA flag, and unshrunk/shrunk information component in R. The canonical route is PhViD (PRR, ROR, BCPNN, and GPS/EBGM in one call, computing the full posterior rather than the gamma-constant shortcut);
library(data.table)
## ---- Canonical route: PhViD (all four estimators incl. full Bayesian posteriors) ----
## library(PhViD)
## ## as.PhViD() needs one row per (drug, event) pair with its AGGREGATED report count (a),
## ## not a report-level table -- pre-aggregate first:
## agg <- reports[, .(n11 = uniqueN(report_id)), by = .(drug_name, event_pt)]
## dm <- as.PhViD(as.data.frame(agg))
## PRR(dm, RR0 = 1, MIN.n11 = 3) # PRR with a>=3 gate
## ROR(dm, OR0 = 1, MIN.n11 = 3) # ROR() takes OR0, not RR0 -- different signature than PRR()
## BCPNN(dm, RR0 = 1, MIN.n11 = 3) # information component
## GPS(dm, RR0 = 1, MIN.n11 = 3) # EBGM/EB05 (whole-database EM fit)
## ## PhViD's default DECISION rule is an FDR-controlled signal list, NOT the fixed MHRA /
## ## IC025>0 / EB05>=2 screening conventions used elsewhere in this entry -- do not treat
## ## $SIGNALS as already applying those cutoffs. Read the point estimate and credibility/
## ## confidence bound off $ALLSIGNALS and apply the convention you intend explicitly (or set
## ## DECISION = 3 with RANKSTAT/DECISION.THRES matched to that convention).
## ---- Manual fallback: PRR, ROR + 95% CI, Yates chi-square, unshrunk/shrunk IC ----
disproportionality <- function(reports, gamma = 0.5) {
dt <- unique(as.data.table(reports)[, .(report_id, drug_name, event_pt)]) # one mention per report-drug-event triple
N <- dt[, uniqueN(report_id)] # a+b+c+d
n_drug <- dt[, .(nd = uniqueN(report_id)), by = drug_name]
n_event <- dt[, .(ne = uniqueN(report_id)), by = event_pt]
tab <- dt[, .(a = uniqueN(report_id)), by = .(drug_name, event_pt)]
tab <- merge(tab, n_drug, by = "drug_name")
tab <- merge(tab, n_event, by = "event_pt")
tab[, `:=`(b = nd - a, c = ne - a)]
tab[, d := N - a - b - c]
## 0.5 continuity correction on all four cells (applied to BOTH PRR and ROR) if any cell is
## zero -- correcting ROR only leaves PRR dividing by a zero comparator proportion when c=0.
tab[, `:=`(pa = ifelse(pmin(a,b,c,d) > 0, a, a + .5),
pb = ifelse(pmin(a,b,c,d) > 0, b, b + .5),
pc = ifelse(pmin(a,b,c,d) > 0, c, c + .5),
pd = ifelse(pmin(a,b,c,d) > 0, d, d + .5))]
tab[, PRR := (pa / (pa + pb)) / (pc / (pc + pd))]
tab[, ROR := (pa * pd) / (pb * pc)]
tab[, se := sqrt(1/pa + 1/pb + 1/pc + 1/pd)]
tab[, `:=`(ROR_lo = exp(log(ROR) - 1.96 * se),
ROR_hi = exp(log(ROR) + 1.96 * se))]
## Yates-corrected chi-square on each 2x2. Margins are cast to double before multiplying --
## (c+d)*(b+d) on a real FAERS-scale table (tens of millions of reports) overflows R's 32-bit
## integer type and silently returns NA. The continuity term is clamped at 0 so it cannot go
## negative and inflate chi2 near independence (matches scipy's chi2_contingency behavior).
tab[, `:=`(E11 = as.double(a + b) * (a + c) / N,
E12 = as.double(a + b) * (b + d) / N,
E21 = as.double(c + d) * (a + c) / N,
E22 = as.double(c + d) * (b + d) / N)]
tab[, chi2 := pmax(abs(a - E11) - .5, 0)^2 / E11 + pmax(abs(b - E12) - .5, 0)^2 / E12 +
pmax(abs(c - E21) - .5, 0)^2 / E21 + pmax(abs(d - E22) - .5, 0)^2 / E22]
## Whole-database independence baseline: RRR, unshrunk IC, gamma-shrunk IC (raw a, E_a --
## not the PRR/ROR continuity-corrected cells).
tab[, E_a := E11]
tab[, RRR := a / E_a]
tab[, IC := log2(RRR)]
tab[, IC_shrunk := log2((a + gamma) / (E_a + gamma))]
tab[, MHRA_signal := PRR >= 2 & chi2 >= 4 & a >= 3]
tab[order(-PRR), .(drug_name, event_pt, a, b, c, d, PRR, ROR, ROR_lo, ROR_hi, chi2,
E_a, RRR, IC, IC_shrunk, MHRA_signal)]
}
## Worked-example check: a report-level fixture with 25/4975/500/494500 rows for
## DrugY-QT_prolongation / DrugY-other / OtherDrug-QT_prolongation / OtherDrug-other
## reproduces PRR=4.95, ROR=4.97, chi2=71.4, IC=2.25, IC_shrunk=2.15 exactly.PRR, ROR, Yates chi-square, MHRA flag, and unshrunk/shrunk information component in SAS. PROC SQL assembles the per-pair 2x2 margins; a DATA step computes every statistic including the whole-database independence baseline E[a] that the IC shrinkage formulas use.
/* One drug-event mention per report so an ICSR cannot double-count the same pair. */
proc sql;
create table pairs as
select distinct report_id, drug_name, event_pt from work.reports;
select count(distinct report_id) into :N trimmed from pairs; /* a+b+c+d */
quit;
proc sql;
create table cells as
select p.drug_name, p.event_pt,
count(distinct p.report_id) as a, /* both */
(select count(distinct x.report_id) from pairs x
where x.drug_name = p.drug_name) as n_drug, /* a+b */
(select count(distinct y.report_id) from pairs y
where y.event_pt = p.event_pt) as n_event /* a+c */
from pairs p
group by p.drug_name, p.event_pt;
quit;
%let gamma = 0.5;
data dispro;
set cells;
b = n_drug - a;
c = n_event - a;
d = &N - a - b - c;
/* 0.5 continuity correction on all four cells (applied to BOTH PRR and ROR) if any cell
is zero -- correcting ROR only leaves PRR dividing by a zero comparator proportion
when c=0. */
ca = ifn(min(a,b,c,d) > 0, a, a + 0.5);
cb = ifn(min(a,b,c,d) > 0, b, b + 0.5);
cc = ifn(min(a,b,c,d) > 0, c, c + 0.5);
cd = ifn(min(a,b,c,d) > 0, d, d + 0.5);
/* PRR: comparator excludes a and b entirely. */
PRR = (ca / (ca + cb)) / (cc / (cc + cd));
/* ROR + 95% CI on the log scale. */
ROR = (ca * cd) / (cb * cc);
se = sqrt(1/ca + 1/cb + 1/cc + 1/cd);
ROR_lo = exp(log(ROR) - 1.96 * se);
ROR_hi = exp(log(ROR) + 1.96 * se);
/* Yates-corrected chi-square for the 2x2, on the raw (uncorrected) cells; the continuity
term is clamped at 0 so it cannot go negative and inflate chi2 near independence. */
E11 = (a + b) * (a + c) / &N; E12 = (a + b) * (b + d) / &N;
E21 = (c + d) * (a + c) / &N; E22 = (c + d) * (b + d) / &N;
chi2 = (max(abs(a - E11) - 0.5, 0))**2 / E11 + (max(abs(b - E12) - 0.5, 0))**2 / E12
+ (max(abs(c - E21) - 0.5, 0))**2 / E21 + (max(abs(d - E22) - 0.5, 0))**2 / E22;
/* Whole-database independence baseline (E11 above = E[a]); unshrunk and shrunk IC use the
raw a, E_a -- not the PRR/ROR continuity-corrected cells. */
E_a = E11;
RRR = a / E_a;
IC = log2(RRR);
IC_shrunk = log2((a + &gamma) / (E_a + &gamma));
MHRA_signal = (PRR >= 2 and chi2 >= 4 and a >= 3);
run;
proc sort data=dispro; by descending PRR; run;Citations
- [1]Evans SJW, Waller PC, Davis S. Use of proportional reporting ratios (PRRs) for signal generation from spontaneous adverse drug reaction reports. Pharmacoepidemiology and Drug Safety. 2001;10(6):483-486.
- [2]Rothman KJ, Lanes S, Sacks ST. The reporting odds ratio and its advantages over the proportional reporting ratio. Pharmacoepidemiology and Drug Safety. 2004;13(8):519-523.
- [3]Bate A, Lindquist M, Edwards IR, et al. A Bayesian neural network method for adverse drug reaction signal generation. European Journal of Clinical Pharmacology. 1998;54(4):315-321.
- [4]DuMouchel W. Bayesian data mining in large frequency tables, with an application to the FDA spontaneous reporting system. The American Statistician. 1999;53(3):177-190.
- [5]van Puijenbroek EP, Bate A, Leufkens HGM, Lindquist M, Orre R, Egberts ACG. A comparison of measures of disproportionality for signal detection in spontaneous reporting systems for adverse drug reactions. Pharmacoepidemiology and Drug Safety. 2002;11(1):3-10.