Skip to content

Commit afbedc5

Browse files
committed
Fix 4 bugs: ZIP fallback, Emax bound, NCA terminal min-points, log-rank tolerance
drug_synergy_tool.py — Bug #7 (CRITICAL): ZIP fallback was Bliss, not ZIP - ZIP requires Hill curve fitting; without it a Bliss independence surface was silently returned labelled as a ZIP delta score — scientifically wrong - Now returns a clear error naming the drug(s) that failed Hill fitting and directing the user to calculate_bliss for a non-parametric alternative drug_synergy_tool.py — Bug #8: Emax upper bound 1.5 in Loewe/CI Hill fitter - _fit_hill_for_loewe uses fractional inhibition (0-1); Emax > 1 is physically impossible, yet the scipy curve_fit upper bound allowed 1.5 - Corrected to 1.0 (Loewe & Muischnek 1926; Chou & Talalay 1984) nca_tool.py — Bug #4: terminal slope accepted 2 instead of 3 post-Tmax points - FDA NCA guidance requires 3+ points for lambda_z: a 2-point linear regression is a perfect fit by construction (R2=1 trivially), giving no quality signal - Now requires 3+ post-Tmax positive concentrations before using the full set; falls back to last n_points otherwise survival_tool.py — Bug #6: log-rank zero-variance check used absolute tolerance - Threshold 1e-10 too tight for large studies: float rounding errors over many event times scale as N x 2.2e-16, so 1e-10 would erroneously reject valid proportional distributions in datasets with thousands of events - Replaced with relative tolerance: 1e-8 x total_events (safe up to 10^7 events)
1 parent 888e2e3 commit afbedc5

3 files changed

Lines changed: 34 additions & 9 deletions

File tree

src/tooluniverse/drug_synergy_tool.py

Lines changed: 19 additions & 7 deletions
Original file line numberDiff line numberDiff line change
@@ -289,12 +289,24 @@ def fit_hill(doses, effects):
289289
params_b = fit_hill(db, effects_b_marginal)
290290

291291
if params_a is None or params_b is None:
292-
# Fallback: simplified ZIP using means
293-
expected_zip = (
294-
np.outer(effects_a_marginal, np.ones(len(db)))
295-
+ np.outer(np.ones(len(da)), effects_b_marginal)
296-
- np.outer(effects_a_marginal, effects_b_marginal)
297-
)
292+
# ZIP fundamentally requires Hill curve fits for each drug to compute
293+
# the expected interaction surface. Without fitted Hill parameters we
294+
# cannot compute a valid ZIP score — returning a Bliss-independence
295+
# surface here would silently report the wrong model.
296+
failed = []
297+
if params_a is None:
298+
failed.append("drug A")
299+
if params_b is None:
300+
failed.append("drug B")
301+
return {
302+
"status": "error",
303+
"error": (
304+
f"ZIP model requires Hill curve fitting for each drug, but fitting "
305+
f"failed for {' and '.join(failed)}. "
306+
"Ensure each drug has ≥3 non-zero dose points with measurable inhibition. "
307+
"Use calculate_bliss for a simpler non-parametric synergy score."
308+
),
309+
}
298310
else:
299311
# ZIP expected using Hill fits
300312
pred_a = np.array([hill_curve(d, *params_a) for d in da])
@@ -348,7 +360,7 @@ def _fit_hill_for_loewe(self, doses, effects):
348360
dm_init = float(np.median(d))
349361
m_init = 1.0
350362
p0 = [dm_init, m_init, emax_init]
351-
bounds = ([1e-15, 0.1, 1e-6], [np.inf, 10.0, 1.5])
363+
bounds = ([1e-15, 0.1, 1e-6], [np.inf, 10.0, 1.0])
352364

353365
def hill(x, dm, m, emax):
354366
x = np.maximum(x, 1e-15)

src/tooluniverse/nca_tool.py

Lines changed: 5 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -153,15 +153,19 @@ def _estimate_terminal_slope(
153153
return None, None, None
154154

155155
# Select terminal phase points
156+
# FDA NCA guidance requires ≥3 points for reliable λz: a 2-point regression
157+
# is a perfect fit by construction (R²=1 always), giving no quality signal.
156158
if tmax is not None:
157159
post_tmax = t_valid > tmax
158-
if np.sum(post_tmax) >= 2:
160+
if np.sum(post_tmax) >= 3:
159161
# FDA NCA: use all post-Tmax points for λz estimation
160162
t_term = t_valid[post_tmax]
161163
c_term = c_valid[post_tmax]
162164
else:
163165
# Fallback: not enough post-Tmax points, use last n_points
164166
n_use = min(n_points, len(t_valid))
167+
if n_use < 2:
168+
return None, None, None
165169
t_term = t_valid[-n_use:]
166170
c_term = c_valid[-n_use:]
167171
else:

src/tooluniverse/survival_tool.py

Lines changed: 10 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -355,7 +355,16 @@ def _log_rank_test(self, arguments: Dict[str, Any]) -> Dict[str, Any]:
355355
# split is proportional (O==E for every stratum). In this degenerate case
356356
# chi2 = 0 and p = 1.0 is the correct, well-defined answer: there is no
357357
# evidence of a difference between the two groups.
358-
if abs(O1_total - E1_total) < 1e-10 and abs(O2_total - E2_total) < 1e-10:
358+
# Use relative tolerance: floating-point sums over many events accumulate
359+
# rounding errors proportional to total_events × machine epsilon (~2.2e-16),
360+
# so an absolute threshold of 1e-10 is too tight for studies with thousands
361+
# of events. Relative tolerance 1e-8 is safe up to ~10^7 events.
362+
total_events = O1_total + O2_total
363+
rel_tol = 1e-8 * max(total_events, 1.0)
364+
if (
365+
abs(O1_total - E1_total) <= rel_tol
366+
and abs(O2_total - E2_total) <= rel_tol
367+
):
359368
chi2_stat = 0.0
360369
p_value = 1.0
361370
else:

0 commit comments

Comments
 (0)