Add native Hill dose-response fitting - #1097
Conversation
Codecov Report✅ All modified and coverable lines are covered by tests. Additional details and impacted files@@ Coverage Diff @@
## main #1097 +/- ##
==========================================
+ Coverage 79.95% 80.09% +0.14%
==========================================
Files 55 55
Lines 7533 7587 +54
==========================================
+ Hits 6023 6077 +54
Misses 1510 1510
🚀 New features to boost your workflow:
|
|
Great! Let me know when you want me to have a look, please. |
|
Thanks! It's ready for review whenever you have time. |
Zethson
left a comment
There was a problem hiding this comment.
Thanks! I think we're getting there.
| similar = ds.nearest_perturbations(ds_adata, "IFNGR2", target_col="gene_target") | ||
| ``` | ||
|
|
||
| See [perturbation space tutorial](https://pertpy.readthedocs.io/en/latest/tutorials/notebooks/perturbation_space.html). |
There was a problem hiding this comment.
We're removing it here and then stating this again after the Dose-response curve fitting section where we're however only stating that we're explaining the dose-response stuff there. I'd just keep this as the more general statement without specifying what can be found there.
| See [perturbation space tutorial](https://pertpy.readthedocs.io/en/latest/tutorials/notebooks/perturbation_space.html). | ||
| ### Dose-response curve fitting | ||
|
|
||
| `PerturbationSpace.dose_response` calculates a scalar distance from control for every perturbation and dose. |
| ### Dose-response curve fitting | ||
|
|
||
| `PerturbationSpace.dose_response` calculates a scalar distance from control for every perturbation and dose. | ||
| `PerturbationSpace.fit_dose_response` fits a four-parameter Hill curve to these values and reports the EC50 for this transcriptomic-distance response: |
|
|
||
| adata = pt.dt.srivatsan_2020_sciplex2() | ||
| adata = adata[adata.obs["dose_value"].notna()].copy() | ||
| adata.obs["dose_value"] = adata.obs["dose_value"].astype(float) |
There was a problem hiding this comment.
I wonder whether these types of changes should just be made to the dataloader?
| adata = pt.dt.srivatsan_2020_sciplex2() | ||
| adata = adata[adata.obs["dose_value"].notna()].copy() | ||
| adata.obs["dose_value"] = adata.obs["dose_value"].astype(float) | ||
| adata.obs["perturbation"] = adata.obs["perturbation"].astype(str) |
There was a problem hiding this comment.
I wonder whether these types of changes should just be made to the dataloader?
But why are we doing this in the first place? Is the problem the categorical dtype? I think we should keep it
| def _four_parameter_logistic( | ||
| dose: np.ndarray, e0: float, emax: float, log_midpoint: float, hill_coefficient: float | ||
| ) -> np.ndarray: | ||
| """Evaluate a four-parameter Hill curve with a positive Hill coefficient.""" |
There was a problem hiding this comment.
What does "evaluate" mean here? Maybe just dummy me doesn't immediately know what you mean here.
| Responses are rescaled internally for numerical stability; ``e0`` and ``emax`` retain the input units. | ||
| The standard error uses a local linear approximation. It is NaN, with a warning, if the parameter | ||
| covariance is non-finite or numerically rank deficient, or there are no residual degrees of freedom. | ||
| Parameters are retained for inspection. A small standard error, high R-squared or an in-range midpoint |
There was a problem hiding this comment.
| Parameters are retained for inspection. A small standard error, high R-squared or an in-range midpoint | |
| Parameters are retained for inspection. | |
| A small standard error, high R-squared or an in-range midpoint |
Nit: the second sentence should be on a new line for better diffs
| Parameters are retained for inspection. A small standard error, high R-squared or an in-range midpoint | ||
| does not establish that the doses capture both plateaus or that the Hill model is appropriate. | ||
|
|
||
| Examples: |
There was a problem hiding this comment.
See above. There must be a way to keep this example simpler
| } | ||
| ) | ||
|
|
||
| return pd.DataFrame.from_records(records) |
There was a problem hiding this comment.
Should this return a DF or can we write something useful into the AnnData object? The whole idea of AnnData is to to always write results into the object. Ideally not uns but probably better than return a DF.
|
|
||
|
|
||
| def test_dose_response(rng): | ||
| def test_dose_response(): |
PR Checklist
docsis updatedDescription of changes
Adds native
PerturbationSpace.fit_dose_responsefor the Hill-fitting part of #1036. It accepts the table returned bydose_response()or another scalar assay response and returns per-perturbation curve parameters, relative EC50/IC50, approximate standard error, R-squared and a dose-range indicator.Technical details
Uses a four-parameter Hill curve with SciPy and adds no new dependency.
response_typeselects EC50 or IC50 terminology; biological/control normalization remains the user's responsibility. The documentation explains the limitations of midpoint and uncertainty estimates. Bliss/Loewe and other combination analyses remain outside this PR.Additional context
Validation: 23 tests passed in
test_perturbation_space_extras.py; Ruff and format checks passed for the changed Python files.