[eboyer221]

MDR: multidrug-resistance targets and design spec

This PR adds multidrug resistance (MDR) as a new unit, next to drug and drug class, with two targets, and adds design/mdr.md, the spec for both.

1. Yes/no MDR (target = "mdr"): a genome is MDR if it is resistant to 2 or more drug classes. It is an ordinary Resistant/Susceptible target, so only the profile changes.
2. Combination (target = "resistant_classes"): a port of old amRml's MDR model. Each combination of resistant classes (e.g. MAC_PEN_TMD) is one label of a multinomial glmnet model. On S. flexneri it reproduces old amRml exactly (see Testing).

Unlike old MDR, both targets go through the same stages and rules as drugs and classes, and run in the year and country modes, not only the whole dataset. cross_drug and lodo only take drug rows, so they leave MDR out.

Up to date with runner through de13896 (cross-drug fit reuse). The merge only needed conflicts in fit-models.R resolved: .fitTask() now keeps both the reused fit and the MDR label levels. Fit reuse only applies to cross_drug, which never includes MDR. After the merge, the S. flexneri comparison below is unchanged.

What's in this PR

file	what it does
design/mdr.md	MDR spec: what old amRml did, a data check, the design for each stage, and the comparison with old amRml
R/registry.R	.MULTICLASS_TARGETS and .isMulticlass(): the one place that declares the combination target multi-class, next to .MODES. Eligibility and fitting read it
R/orb-profile.R	profileORB() gains mdr_min_classes (default 2) and the $mdr table, and adds unit = "mdr" rows for both targets. n_resistant now counts anything not Susceptible, so a combination counts as resistant
R/eligibility.R	eligibleScopes() gains min_mdr_label_genomes (default 25, old min_n) and the per-scope label filter. Scopes gain n_labels, n_dropped, test_n_labels, test_n_dropped
R/fit-models.R	multinomial path: model, folds, metrics, per-label importance, new results columns
R/workflow.R	passes min_mdr_label_genomes through
man/	regenerated. fitModels.Rd and runModelingWorkflow.Rd also pick up the BPPARAM argument from #52, which hadn't been documented yet
tests/	profile and eligibility tests for both targets; new test-fit-models.R

How it works

profile <- profileORB(readORB("path/to/orb"))
profile$mdr   # genome_id, n_classes_tested, n_classes_resistant, mdr
eligibility <- eligibleScopes(profile)   # MDR scopes alongside drugs and classes

Yes/no MDR. Each genome gets one status, built from the class phenotypes the profile already computes:

• mdr: Resistant in at least mdr_min_classes classes. Phenotype Resistant.
• not_mdr: tested in at least mdr_min_classes classes, but Resistant in fewer. Phenotype Susceptible.
• undetermined: tested in fewer classes, so MDR can't be ruled out. Listed in $mdr with no phenotype, so it never enters the model. This matters on real data: 931 of 2,907 E. faecium genomes were tested in only one class.

Combination.

• Profile: each genome tested in a class gets its resistant classes in radix order, joined with _, or Susceptible. That matches amRdata's resistant_classes for every genome in S. flexneri, S. argenteus and E. faecium.
• Eligibility: each scope keeps the labels with at least min_mdr_label_genomes training genomes and drops every genome, training or test, with another label. On the whole dataset this gives old amRml's genome set exactly. In cross and leave-one-out modes it also drops test genomes whose label the model never trained on. The usual rules then apply, with "phenotype" meaning label: the scope needs 2 labels (single_phenotype), and its rarest label must fill the folds. Its test set needs only min_test_genomes.
• Fitting: a task whose target is in .MULTICLASS_TARGETS uses multinom_reg() with the same engine, grid and MCC selection. Folds are now assigned per label present, not only to Resistant and Susceptible rows; Resistant/Susceptible fold assignment is unchanged (tested).

Results columns (shared tables, no separate MDR tables):

table	new column	for drugs and classes	for the combination target
importance	label	NA	the label the coefficient belongs to; "POS" points to it
performance	n_labels	2	number of labels
predictions	prob_predicted	probability of the predicted phenotype	probability of the predicted label

For the combination target, performance has only mcc and nmcc, as old amRml reported; res_prop, res_prop_test and prob_resistant are NA.

Testing

Old vs new, S. flexneri (amRdata main, October 2026; old amRml from #42):

• Matrices: all 13 old MDR matrices (.parquet2MDRMatrix(), min_n = 25) are identical to the new whole-dataset combination matrices: the same 288 genomes and 4 labels, the same features, and the same non-zero values.
• Fit, AMRFinder binary, seed 123: old runMLPipeline() and new fitModels() choose the same penalty (1e-4) and mixture (0.2), test on the same 59 genomes, predict the same label for every one, and score MCC 0.85. All 521 non-zero coefficients match exactly.

Other data:

• S. flexneri, yes/no: 300 MDR, 75 not MDR, 4 undetermined, matching an independent count from metadata.parquet. The whole-dataset scope is eligible, and it built and fitted end to end on gene features with no changes to matrix building or fitting.
• S. argenteus: old amRml fits a 1-label combination model (PEN, 65 genomes). The new code records single_phenotype. The yes/no target (5 MDR, 68 not MDR) fails the fold check.
• Year and country modes: S. flexneri's kept genomes fall in one year bin and one country, so MDR is rejected there for the expected reasons. E. faecium is the only dataset where they would run, and its ORB's amRml_input is still marked as building, so it hasn't been tested.

R CMD check: 0 errors. The one warning (tidyselect not in Imports, from .buildRecipe() on runner) and the notes don't come from this PR.

Unit tests: 92 tests; all new ones pass. 7 tests in test-orb-read.R and test-orb-profile.R already fail on runner (7ba0cca), without this PR. Three of them (two on missing_features, one on moved ORBs) test checks that the ORB-structure changes commented out; I haven't looked into the other four. Four existing eligibility and profile tests were updated because the new target adds scopes or rows. Their old expectations no longer apply; nothing else changed in them.

Questions for reviewers

• Undetermined genomes in the combination target. The yes/no target leaves out genomes tested in fewer than 2 classes (undetermined), but the combination target keeps them, labelled from the one class they were tested in. That's deliberate, to match old amRml exactly. Should the combination target leave them out too? It would change results compared with old amRml.

• The yes/no gene-feature fits scored MCC = 1.00 on 75 test genomes. I expect that's because the classes behind MDR here (penicillins, trimethoprim, macrolides) have well-known resistance genes in the gene features, but I haven't checked the top features yet.

• New scope columns (n_labels, n_dropped, test_n_labels, test_n_dropped) and results columns (label, n_labels, prob_predicted) should be included in the global rename.

• Old runMDRmodels() defaulted to split = c(0.8, 0); I compared with split = c(1, 0), as in the earlier Fitting and running models through a single entry point #52 comparison. Worth confirming that's the split we want to match.

[eboyer221]

@amcim @epbrenner #54 (MDR) is ready for review. It's up to date with runner and has no conflicts. The commits are split by stage, so it may be easiest to go commit by commit. On S. flexneri it reproduces old amRml's MDR matrices and fit exactly; details are in the description.

Two things that touch your PRs:

1. Feature rescoring (#53) and the MDR combination target. The combination target has importance per label (a label column). #53's exportMergedResults() drops label, and rescoring sums importance per drug, so the 4 labels' coefficients for a feature would get mixed together. Old amRml avoided this by keeping MDR in separate MDR_* files.

I have a fix ready if you want it: branch mdr-rescoring-exclude (off feature-rescoring, not pushed yet). It leaves the combination target out of the merged export, with a test. The yes/no mdr target still exports like a drug. It does nothing until #54 merges. Happy to open it as a PR into feature-rescoring, or you can handle it differently, for example per-label rescoring, or after the rename.

2. Merge order. #52, #53 and #54 all change fit-models.R. I'd suggest #52, then #54, then the rename, then #53 once #17 is done, so #53 only has to be updated once against everything else. Does that work for you?

Also, #54 adds .MULTICLASS_TARGETS / .isMulticlass() to registry.R, next to .MODES, so stages check a declared property instead of target names. It might be worth folding into the rename.

[AbhirupaGhosh]

I have done my first pass and things look ok to me. No major error or logical problem detected.

[amcim]

Overall this looks pretty good

@eboyer221 Answer to what you raised

Old runMDRmodels() defaulted to split = c(0.8, 0); I compared with split = c(1, 0), as in the earlier #52 comparison. Worth confirming that's the split we want to match.

• Basically there was something that was very misleading at best in the original code, where if its c(1,0) it would still hold out 20%, and c(.8,0) never actually run, it would throw an error. So c(1,0) is the right setting to match

Your point about PR #53

• What you propose makes sense, we just need to see how PR Feature rescoring module and outputting merged results and model fits #53 needs to change to account for how the output is now stored and whatever changes to magic were made. So lets keep note of this for followup

merge order

• I am glad you brought this up as I was going to, and I fully agree. Once PR Fitting and running models through a single entry point #52 is merged, we merge amrml-rewrite branch into here and that should deal with the merge conflicts

• As for the rename, I would like to do this during code review week. Now that @epbrenner has gotten the output more how he wants it, making things consistent is for sure doable now

[eboyer221]

Merging. Noting what's still open so we can come back to it:

Added after review: @AbhirupaGhosh suggested dropping MDR when a dataset has fewer than 2 classes. The yes/no target already had no model there (every genome is undetermined), but the combination target was still built and repeated the class model. New eligibility rule single_class_labels (85a1860) skips any combination scope whose kept labels involve fewer than 2 classes. S. flexneri and S. argenteus results are unchanged.

Not tested yet:
MDR in the year and country modes. Only E. faecium has enough years and countries among MDR genomes to run them, and its ORB's amRml_input is still marked as building.

Follow-ups in other PRs:

• Feature rescoring module and outputting merged results and model fits #53 feature rescoring. Its export drops label, so rescoring would mix up the combination model's per-label scores. The fix (leave the combination target out of the merged export) is on mdr-rescoring-exclude. When Feature rescoring module and outputting merged results and model fits #53 is updated, also check Evan's new network.R, which reads importance too.
• Rename (code review week). New MDR names to include: label, n_labels, prob_predicted in results; n_labels, n_dropped, test_n_labels, test_n_dropped in eligibility; min_mdr_label_genomes, mdr_min_classes; .MULTICLASS_TARGETS / .isMulticlass().

Not from this PR:
-test-fit-resume.R ("dense and sparse fits agree across serial, sockets and reuse settings") fails on amrml-rewrite with could not find function ".saveFitCheckpoint" when run with devtools::test(). Likely the socket workers not seeing internal functions when the package isn't installed.