Skip to content

Secondary object#171

Open
akcd1 wants to merge 28 commits into
mainfrom
secondary_object
Open

Secondary object#171
akcd1 wants to merge 28 commits into
mainfrom
secondary_object

Conversation

@akcd1

@akcd1 akcd1 commented Dec 5, 2025

Copy link
Copy Markdown
Collaborator

Description

Thank you for your contribution to Brieflow!
Please succinctly summarize your proposed change.
What motivated you to make this change?

Added the option to segment ang phenotype additional objects contained within cells.

Please also link to any relevant issues that your code is associated with.

What is the nature of your change?

  • Bug fix (fixes an issue).
  • Enhancement (adds functionality).
  • Breaking change (fix or feature that would cause existing functionality to not work as expected).
  • This change requires a documentation update.

Checklist

Please ensure that all boxes are checked before indicating that a pull request is ready for review.

  • My code follows the conventions of this project.
  • I have updated the pyproject.toml to reflect the change as designated by semantic versioning.
  • I have checked linting and formatting with ruff check and ruff format.
  • I have performed a self-review of my own code.
  • I have commented my code, particularly in hard-to-understand areas.
  • I have made corresponding changes to the documentation.
  • My changes generate no new warnings.
  • I have deleted all non-relevant text in this pull request template.

@akcd1 akcd1 requested a review from mat10d December 5, 2025 20:09
@mat10d

mat10d commented Dec 7, 2025

Copy link
Copy Markdown
Collaborator

@akcd1 the tests are failing on something that suggests that you haven't properly rebased main. I would remerge asap.

@akcd1

akcd1 commented Dec 8, 2025

Copy link
Copy Markdown
Collaborator Author

re-based

@mat10d mat10d marked this pull request as draft March 30, 2026 17:07
acepedadiaz and others added 21 commits April 10, 2026 15:12
* conditional cytoplasm mask

* initial commit, refactoring second obj segmentation

* added thresholding options

* added visualization of thresholded output

* ruff check and format

* adds helper visualization function

* snakemake integration

* updated parameters and added scaffold for ML-based segmentation

* ruff check and format

* code consolidation

* improved documentation

* code consolidation and reformatting

* remove foci channel from second_obj pheno extraction, rename foci_channel_index elsewhere

* fixed bug foci_channel param passing

* fixed import of second obj features

* conditional cytoplasm mask

* Integrated secondary object segmentation with ml models

* ruff format

* disabled secondary object detection for test data

* consistent naming with config and notebook

* fixed ML parameter passing in script identify_second_objs.py

* fixed general config parameters

---------

Co-authored-by: edeny <edeny@fry.wi.mit.edu>
Co-authored-by: acepedadiaz <acepedadiaz@fry.wi.mit.edu>
* adds fine-tuning of Cellpose, code consolidation

* optimization of model performance calculations

* added image normalization for improved alignment performance

* updated to _cp nomenclature

* update to merge_phenotype_cp

* added intermediate tsv reporting alignment offsets and confidence

* added intermediate tsv reporting alignment offsets and confidence

* removed params for percentile normalization

* added alignment metrics

* removed percentile normalization, added reporting, normalization bug fixed

* minor formatting

* added alignment metrics

* minor formatting

* added alignment metrics

---------

Co-authored-by: acepedadiaz <acepedadiaz@fry.wi.mit.edu>
Co-authored-by: acepedadiaz <acepedadiaz@c13b1.wi.mit.edu>
Route per-tile offset_y/offset_x from alignment_metrics TSV through
extract_phenotype_cp so they flow naturally through merge_phenotype_cp
into final_merge with zero merge module changes.

- Add alignment metrics TSV as 5th input to extract_phenotype_cp rule
- Broadcast offset columns to each cell row in extract_phenotype.py
- Fix stitch_phenotype input to use align_phenotype[0] (image only)

Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
Introduce aggregate_cells_second_objs rule between final_merge and
split_datasets with 4 strategies for handling one-to-many cell-object
relationships: "none", "single", "all", "average".

- New second_obj_utils.py with strategy dispatch and merge logic
- New aggregate_cells_second_objs.py snakemake script
- Conditional DAG: uses aggregate output when second_obj_detection is
  enabled, falls back to final_merge otherwise
- Add cell_summary columns to metadata cols so they are excluded from
  feature aggregation

Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
@akcd1 akcd1 force-pushed the secondary_object branch from 66c3463 to 38272f4 Compare April 10, 2026 20:27
acepedadiaz and others added 2 commits April 10, 2026 16:43
Test images have no secondary objects. Set second_obj_detection: false
in test config and remove secondary object metadata columns from
cell_data_metadata_cols.tsv since they won't exist in test data.

Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
acepedadiaz and others added 5 commits April 10, 2026 17:10
- merge_phenotype rule now falls back to extract_phenotype_cp when
  second_obj_detection is disabled, fixing MissingInputException in CI
- Add r-prefix to docstring containing backslashes in segment_cellpose.py

Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
Filter aggregate_cells_second_objs from targets/mappings when
second_obj_detection is disabled, matching the pattern in phenotype
targets. Fixes MissingInputException in CI tests. Also add period
to cellpose_training.py docstring to satisfy ruff D415.

Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
calculate_offsets was changed to return (offsets, errors) but three
call sites in align_cycles.py were not updated to destructure the
tuple, causing ValueError in apply_offsets when offset elements
were arrays instead of scalars.

Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
The script was updated to load alignment_metrics from snakemake.input[1],
but the SBS extract_sbs_info rule only provides one input (nuclei tiff).
The phenotype extract_phenotype_info rule provides both. Guard the metrics
loading with a length check so the script works for both callers.

Co-Authored-By: Claude Opus 4.6 (1M context) <noreply@anthropic.com>
…ary_object (joint scrubbed) (#222)

* feat(aggregate): add compartment_combo_subset to filter feature columns by compartment

* feat(aggregate): add resolve_aggregate_combos validation helper

* feat(shared): register compartment_combo filename prefix (CmCo-)

* feat(aggregate): add compartment_combo wildcard to aggregate output paths

* feat(aggregate): add compartment_combo wildcard to aggregate rules

* feat(shared): add compartment_combo to bootstrap rule_utils helpers

* feat(cluster): add compartment_combo wildcard to cluster paths and rules

Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com>

* feat(shared): propagate compartment_combo through metrics, targets, cluster_eval

* feat(aggregate): apply compartment_combo_subset in split_datasets and read wildcard in bootstrap_gene

* test(small_test_analysis): add compartment_combo to aggregate fixtures

* test(small_test_analysis): add compartment_combo to cluster_combo fixture

* fix(aggregate): preserve int64 plate dtype in aggregate_cells_second_objs

The script cast second_objs_df's plate/well to str, then _aggregate_single
merged that against cells_df (plate as int64) and failed with:
  ValueError: You are trying to merge on int64 and object columns for key 'plate'

Both merge_final.parquet and phenotype_second_objs.parquet have plate as
int64; the str coercion was unnecessary. Compare plate as int and leave
well as its natural string form.

Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com>

* fix(aggregate): use known compartment prefixes as drop universe

all_compartments_run was derived from the compartments present in the
TSV. Compartments never requested in any combo (e.g. cytoplasm when all
combos only use cell/nucleus/second_obj) were absent from the universe,
so compartment_combo_subset never added them to `excluded` and their
feature columns leaked through into every filtered parquet. Source the
universe from COMPARTMENT_PREFIXES so undeclared compartments get
dropped instead of passed through.

Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com>

* feat(aggregate): make montage generation opt-out via generate_montages

Read config["aggregate"]["generate_montages"] (default True) in
targets/aggregate.smk and evaluate MONTAGE_TARGETS_ALL to an empty list
when False. This lets users skip prepare_montage_data, generate_montage,
and initiate_montage when montages aren't wanted, without needing to
invoke a different snakemake target. Default preserves existing behavior
for configs that don't set the flag.

Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com>

* fix(aggregate): harmonize pool schema at align + generate_feature_table

Per-well missing_values_filter drops different columns in different wells
when class composition varies (e.g. Uninfected-only wells vs mixed
Bystander+Infected wells). Downstream pooling via pyarrow's dataset
union then fills those cols with NaN across whole wells, which breaks
centerscale_by_batch (zero-variance batch) and PCA (NaN input).

Add harmonize_pool_schema(paths, metadata_cols, drop_cols_threshold) in
lib/aggregate/filter.py that:
  1. Intersects per-file schemas. Cols present in some but not all files
     are dropped up front and logged per-file so the operator can see
     which wells lacked which columns (the "missing class in a well"
     case the user explicitly wants to tolerate).
  2. Applies drop_cols_threshold at the pool level on the intersected
     col set, matching the convention from missing_values_filter.

Wire into align.py and generate_feature_table.py: call the helper once
to produce a canonical metadata+feature col set, pass it as `columns=`
to every ds.dataset scan (PCA sample, class read, batched transform,
feature-table batch loop). Remove the per-batch .dropna(axis=1) in
align.py's batch loop (redundant; could previously drop inconsistent
col sets across batches and break PCA transform). Add defensive
per-subset row-drop on residual NaN.

Pass config["aggregate"]["drop_cols_threshold"] through as a rule
param to both rules.

Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com>

* fix(aggregate): plumb perturbation columns through perturbation_score

perturbation_score.py hardcoded "gene_symbol_0", "nontargeting", and
"cell_barcode_0" across prepare_alignment_data, centerscale_on_controls,
and calculate_perturbation_scores — so any config where
perturbation_name_col != "gene_symbol_0" ran the scorer on mismatched
columns and either produced silently-wrong scores or blew up. Under
pandas nullable-string dtypes, (gene_symbol_0 == <sgRNA barcode>)
returned a BooleanArray with pd.NA, and .astype(int) collapsed it to
{0, 1, <NA>} — which LogisticRegression(solver="liblinear") then
rejected as n_classes >= 3.

Plumb perturbation_name_col, perturbation_id_col, control_key,
control_name_col, and batch_cols through perturbation_score and
_process_gene_subset from the caller. calculate_perturbation_scores
now requires perturbation_col (no default) and builds y via
.eq(gene).fillna(False).astype(bool).astype(int) so NaN rows always
map to 0 and the target is strictly binary.

align.py passes the additional config values into perturbation_score in
both the joint and non-joint paths.

Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com>

* fix(aggregate): dedupe construct_columns when pert_id_col == pert_col

When perturbation_id_col and perturbation_name_col point at the same
column (a valid config when aggregation is at the per-construct level,
e.g. both "cell_barcode_0"), building construct_columns as
[pert_id_col, pert_col, "cell_count"] + feature_cols listed that name
twice. construct_table[construct_columns] then produced duplicate
columns, and construct_table[pert_col] returned a 2-col DataFrame,
breaking the subsequent .str.contains(control_key) call with
AttributeError: 'DataFrame' object has no attribute 'str'.

Use dict.fromkeys to dedupe while preserving order.

Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com>

* fix(aggregate): use .loc instead of .iloc for keep_idx in perturbation_score

keep_idx is built from perturbation_col.index[...] and nt_idx which hold
pandas label indices, but line 102 used .iloc (positional). When
_compute_perturbation_score_joint in align.py iterates groupby(\"class\")
slices without reset_index, cell_data's labels exceed its length and
.iloc raises IndexError: positional indexers are out-of-bounds. Lines
130-131 already use .loc correctly with the same gene_idx; make 102
consistent so the function works with either RangeIndex or preserved
label index inputs.

Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com>

* fix(aggregate): make eval_aggregate compartment-aware and NaN-safe

eval_aggregate.py hardcoded a cell→nucleus fallback chain when picking
*_mean feature columns for the distribution plot. For combos whose
compartment_combo is second_obj or cytoplasm, neither 'cell_*_mean' nor
'nucleus_*_mean' columns survive compartment_combo_subset, so
merge_feature_cols collapsed to [], and the downstream
plt.subplots(nrows=0, ...) raised ValueError.

Two layers:
- Script: derive the prefix priority from snakemake.wildcards.compartment_combo
  and iterate all known compartments (cell, nucleus, cytoplasm, second_obj)
  as a fallback chain, using the first non-empty set of *_mean columns.
- Lib plot_feature_distributions: guard empty inputs by returning a small
  placeholder figure with an explanatory label. eval_aggregate's other
  outputs (na_stats tsv/png) still write; the DAG proceeds.

Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com>

* feat(aggregate): carry secondary name column through aggregation

When aggregate is grouped by a construct ID (perturbation_name_col =
cell_barcode_0), the aggregated.tsv previously had no gene_symbol_0
column. Downstream consumers break: cluster benchmarking
(evaluate_resolution -> filter_complexes) tries to intersect
aggregated perturbation names with CORUM/KEGG gene_name fields and
finds nothing; uniprot merge in notebook 10 cell 20 fails on missing
gene_symbol_0. Both assume gene-level naming.

Add a carry_cols parameter to lib/aggregate/aggregate.py:aggregate()
that preserves columns functionally determined by group_cols in the
output — taking the single unique value per group and raising
ValueError if the invariant is violated. Wire the rule aggregate
to pass carry_cols=[control_name_col] when control_name_col is
different from perturbation_name_col (the normal setup for
barcode-level aggregation + gene-level control naming).

Net effect: when pert=cell_barcode_0 and control=gene_symbol_0,
the aggregated.tsv gains a gene_symbol_0 column and downstream
clustering/benchmarking/annotation all get the name they expect.
When pert == control (gene-level aggregation), carry_cols stays
None and behavior is unchanged.

Co-Authored-By: Claude Opus 4.7 (1M context) <noreply@anthropic.com>

* fix(aggregate): use positional indexing in missing_values_filter imputation

The imputation block read/wrote with .loc[batch_na_idx, cols], which fails
on inputs whose index has duplicate labels. Production callers always feed
unique-indexed frames (scripts/aggregate/filter.py uses ignore_index=True),
but any new caller that doesn't can trip a 'cannot reindex on an axis with
duplicate labels' ValueError. Switch to positional iloc indexing — same
behavior on clean inputs, robust to duplicate-indexed inputs.

* fix(aggregate): plumb control_name_col through generate_feature_table

centerscale_on_controls and the construct/gene control-table masks identified
controls via perturbation_name_col. When perturbation_name_col != control_name_col
(construct-level aggregation, e.g. cell_barcode_0 perturbations with gene_symbol_0
controls) no controls were found: centering-on-controls was silently skipped and
the control gene table was empty. Carry control_name_col through the construct
accumulation and use it for control identification, mirroring align.py. Falls back
to perturbation_name_col when unset, so gene-level runs are unchanged.

Co-Authored-By: Claude Opus 4.8 (1M context) <noreply@anthropic.com>

* refactor(aggregate): add type hints + pin schema in harmonize_pool_schema

- Pass schema=pq.read_schema(paths[0]) to ds.dataset() in harmonize_pool_schema so
  files are read against one reference schema instead of pyarrow auto-unifying
  (which raises ArrowInvalid on dtype mismatch across pooled per-well parquets).
- Add type hints to harmonize_pool_schema, compartment_combo_subset, and
  resolve_aggregate_combos for consistency with other functions in the PR.

Co-Authored-By: Claude Opus 4.8 (1M context) <noreply@anthropic.com>

---------

Co-authored-by: acepedadiaz <acepedadiaz@cheeserind.wi.mit.edu>
Co-authored-by: Claude Opus 4.7 (1M context) <noreply@anthropic.com>
Co-authored-by: acepedadiaz <acepedadiaz@wi.mit.edu>
@akcd1 akcd1 marked this pull request as ready for review June 9, 2026 15:46
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

3 participants