From b2f551beb2f5fa63ad02ab24029e88dca54faa9a Mon Sep 17 00:00:00 2001 From: R script <1695515+ms609@users.noreply.github.com> Date: Mon, 3 Aug 2026 16:23:25 +0100 Subject: [PATCH] fix(hsj): an inapplicable secondary contributes nothing; nesting validates Two red-team findings, both settled against Hopkins & St John (2021) rather than against the existing tests. T-396 -- an inapplicable secondary was scored as a state. Where a secondary was coded "-" but its controlling primary did not CERTAINLY code the structure absent (primary "?", or a matrix that codes the structure present and the secondary inapplicable anyway), fitch_label_char() kept "-" as an ordinary concrete state. Disjoint from every other state, it could be propagated inwards by the uppass and label an internal node in the middle of the PRESENT region, mismatching every secondary of the block at once, so score_hierarchy_block() charged that branch d = m -- the full alpha -- for a node that by construction has no inapplicable secondaries. The paper decides this, and against the shipped test: d counts "nonmatching secondary characters" among those that APPLY, inapplicable secondaries "have no influence on the estimated dissimilarity" (p.5), and treating them as a separate state "increases the dissimilarity of all pairwise comparisons ... overweighting the [controlling] primary character" (p.5) -- which is the FitchS behaviour HSJ exists to replace, not a stricter reading of it. inapp_bit is now stripped from every secondary's state set at every tip, and a cell left with no applicable state is non-constraining. test-ts-hsj.R:857 pinned the over-charge as expected (c(0, 1, 1, 0), with a comment rationalising score("-") == 1). That expectation is wrong and is now c(0, 1, 0, 0); the old value exceeded every concrete resolution of the cell AND the missing-data treatment, so no reading of the data produced it. Measured over-charge was exactly alpha/m per affected branch (witness: 1 vs 0.5 at alpha=0.5, 2 vs 1 at alpha=1). T-395 -- every nested hierarchy was rejected. .ParseOneBlock() records a sub-controller both as a dependent of its parent and as its own block's controlling character; that dual role IS nesting. .ValidateBlock() counted the second occurrence as a double claim, so the documented example in ?CharacterHierarchy could not be validated, scored, or searched -- and R CMD check could not see it, because CharacterHierarchy() never validates: the example constructs fine and is never scored. A child block no longer re-claims its controlling character. Genuine double claims, including a nested block's dependents colliding with another block, still error. HierarchyFromNames() nesting was dead for an independent reason: it tested intersect(subIdx, supIdx), which is necessarily empty since no one name carries both prefixes, and matched tags only as far as their first underscore, which collapsed every depth onto the outermost tag. Nesting is now expressed by tag extension (sup_tail_tip extends sup_tail) with longest-match attachment, so sub_tail_tip_gloss belongs to sup_tail_tip. A shared prefix without an underscore boundary does not nest. Documented and exampled. Verification - New tests fail 6/6 against a baseline build with both fixes reverted, and pass with them; test-ts-hsj.R fails 2 pre-fix, 150 pass post-fix. - Fix gates all green: hsj-paper-oracle 17/17, hsj-token-permutation invariant, hsj-totalwords-zero-noop searched. - Nested HSJ scores are rooting-invariant over all 6 rootings; xform reports its own "nested not supported" error instead of failing validation first. - 673 pass / 0 fail across hsj, hierarchy, xform, sankoff, collapse, recode, tree_length and resample suites. - check_man and spell_check_package(vignettes = TRUE) clean. Also fixes five tests in test-CharacterHierarchy.R that called internal functions without a TreeSearch::: prefix and so errored under a temp-library install, including ".HierarchyToBlocks flattens nested hierarchies" -- the one test of the nested path, which had never executed. It passes. Co-Authored-By: Claude Opus 5 --- NEWS.md | 46 +++++++ R/CharacterHierarchy.R | 145 ++++++++++++++--------- man/HierarchyFromNames.Rd | 19 ++- src/ts_hsj.cpp | 70 ++++++++--- tests/testthat/test-CharacterHierarchy.R | 134 ++++++++++++++++++++- tests/testthat/test-ts-hsj.R | 41 +++++-- 6 files changed, 361 insertions(+), 94 deletions(-) diff --git a/NEWS.md b/NEWS.md index 3492cf2d4..47724a8cb 100644 --- a/NEWS.md +++ b/NEWS.md @@ -1,5 +1,51 @@ # To integrate into 2.0.0 notes +- `inapplicable = "hsj"` no longer charges for an inapplicable secondary + character. Where a secondary was coded `"-"` but its controlling primary did + not certainly code the structure absent -- because the primary was `"?"`, or + because the matrix codes the structure present and the secondary inapplicable + anyway -- the `"-"` was treated as an ordinary state of that secondary. Being + disjoint from every other state, it could then be propagated inwards and used + to label an internal node in the middle of a region where the structure *is* + present, mismatching every secondary of the block at once and charging that + branch the full weight of the scaling parameter. + + This was the behaviour the method exists to avoid. Hopkins & St John (2021) + count only the secondary characters that apply, and note that treating + inapplicable cells as a separate state "increases the dissimilarity of all + pairwise comparisons", overweighting the controlling primary and favouring + clades that separate taxa possessing the structure from those lacking it. An + inapplicable secondary now contributes nothing to the dissimilarity, whatever + its controlling primary codes. + + **HSJ scores may therefore fall on matrices that code a secondary + inapplicable where its primary does not code absence**, by the scaling + parameter divided by the number of secondary characters in the block, for each + affected branch. Scores are unchanged wherever the matrix is coded + consistently, and remain independent of where the tree is rooted. + +- Nested hierarchies now validate, so a `CharacterHierarchy()` describing + tertiary characters can be scored under `inapplicable = "hsj"`. A + sub-controlling character is deliberately recorded both as a dependent of the + character above it and as the controlling character of the one below -- that + dual role is what nesting means -- but validation counted the second + occurrence as one character appearing in two blocks, and so rejected every + nested hierarchy that could be written, including the one documented in + `?CharacterHierarchy`. Genuine double claims are still rejected. + + `HierarchyFromNames()` now detects nesting as its documentation describes. A + controlling character whose tag extends another's, as `sup_tail_tip` extends + `sup_tail`, is nested beneath it, and each dependent attaches to the longest + tag it extends, so `sub_tail_tip_gloss` belongs to `sup_tail_tip` rather than + to `sup_tail`. Previously a tag was matched only as far as its first + underscore, which collapsed every depth onto the outermost tag; nesting was + silently dropped and could not be expressed at all. A shared prefix without + an underscore boundary does not nest, so `sup_tailfin` remains independent of + `sup_tail`. + + The x-transformation still does not implement nesting, and now says so + directly instead of failing validation first. + - `inapplicable = "xform"` scores are now reported at a canonical rooting, so a reported score is reproducible. The x-transformation's step matrix is asymmetric -- a gain costs one more than the number of secondary characters it diff --git a/R/CharacterHierarchy.R b/R/CharacterHierarchy.R index 9a6f07e9a..47fd16ee1 100644 --- a/R/CharacterHierarchy.R +++ b/R/CharacterHierarchy.R @@ -183,7 +183,15 @@ ValidateHierarchy <- function(hierarchy, dataset) { claimed <- integer(0) - .ValidateBlock <- function(node, depth = 1L) { + # `ctrlClaimedByParent` is TRUE for a nested block, whose controlling + # character is *deliberately* also a dependent of its parent: .ParseOneBlock() + # records a sub-controller in both places, because it is simultaneously a + # secondary of the character above it and the primary of the one below. That + # dual role is the whole content of "nested", so counting it as a second claim + # rejected every nested hierarchy that could be written -- including this + # file's own documented example. The double-claim check still applies in full + # to the block's own dependents, and to a controlling character at top level. + .ValidateBlock <- function(node, depth = 1L, ctrlClaimedByParent = FALSE) { ctrl <- node$controlling deps <- node$dependents @@ -198,14 +206,15 @@ ValidateHierarchy <- function(hierarchy, dataset) { } # Check no double-claiming - overlap <- intersect(allIdx, claimed) + newIdx <- if (ctrlClaimedByParent) deps else allIdx + overlap <- intersect(newIdx, claimed) if (length(overlap) > 0L) { stop(sprintf( "Character(s) %s appear in multiple hierarchy blocks.", paste(overlap, collapse = ", ") )) } - claimed <<- c(claimed, allIdx) + claimed <<- c(claimed, newIdx) # Check controlling character is binary (has exactly states "0" and "1", # possibly with inapplicable/missing) @@ -236,9 +245,10 @@ ValidateHierarchy <- function(hierarchy, dataset) { } } - # Recurse into children + # Recurse into children. A child's controlling character was already + # claimed just above, as one of this block's dependents. for (child in node$children) { - .ValidateBlock(child, depth + 1L) + .ValidateBlock(child, depth + 1L, ctrlClaimedByParent = TRUE) } } @@ -254,9 +264,16 @@ ValidateHierarchy <- function(hierarchy, dataset) { #' #' Parse character names following the TNT convention where controlling #' characters are named `sup_` and their dependent characters are -#' named `sub_[_suffix]`. Tags must match between a controlling -#' character and its dependents. Nested hierarchies are detected when a -#' `sub_` character is also a `sup_` for further characters. +#' named `sub_[_suffix]`. Each dependent is attached to the longest +#' `sup_` tag that its own tag extends, at an underscore boundary. +#' +#' Nested hierarchies are written by giving a controlling character a tag that +#' itself extends another controlling character's tag: `sup_tail_tip` controls +#' the tag `tail_tip`, and because `tail_tip` extends `tail` it is also a +#' dependent of `sup_tail`. Longest-match is what keeps `sub_tail_tip_gloss` +#' with `sup_tail_tip` rather than with `sup_tail`. A shared prefix alone does +#' not nest: `sup_tailfin` is independent of `sup_tail`, since `tailfin` does +#' not extend `tail` at an underscore boundary. #' #' @param charNames Character vector of names, one per original character. #' @@ -268,6 +285,11 @@ ValidateHierarchy <- function(hierarchy, dataset) { #' "sup_wing", "sub_wing_venation", "eyes") #' HierarchyFromNames(names) #' +#' # Nesting: `tail_tip` extends `tail`, so character 3 is both a dependent of +#' # character 1 and the controlling character of character 4. +#' HierarchyFromNames(c("sup_tail", "sub_tail_colour", +#' "sup_tail_tip", "sub_tail_tip_gloss")) +#' #' @family tree scoring #' @seealso [CharacterHierarchy()] #' @export @@ -284,72 +306,79 @@ HierarchyFromNames <- function(charNames) { return(NULL) } - # Extract tags + # Extract tags. A tag may itself contain underscores, and that is what makes + # nesting expressible: `sup_tail_tip` controls the tag "tail_tip", and because + # "tail_tip" extends the existing tag "tail" it is simultaneously a dependent + # of `sup_tail`. This is the only way one character can play both roles, since + # a single name cannot carry both the `sup_` and `sub_` prefix -- which is why + # the previous `intersect(subIdx, supIdx)` test could never fire, leaving the + # documented nesting support unreachable. + # + # Every sub_/sup_ character attaches to the LONGEST sup_ tag its own tag + # extends. Longest-wins is what keeps `sub_tail_tip_gloss` with `sup_tail_tip` + # instead of with `sup_tail`; matching only the first underscore-delimited + # component (as this function used to) collapsed every depth onto the outermost + # tag, so no dependent could ever reach a nested controller. supTags <- sub("^sup_", "", charNames[supIdx]) subTagsFull <- sub("^sub_", "", charNames[subIdx]) - # The tag is the first component before any additional underscore-suffix - # e.g. "sub_tail_colour" → tag = "tail" - subTags <- sub("_.*", "", subTagsFull) - # Build mapping: tag → controlling index, tag → dependent indices + # Build mapping: tag → controlling index tagToSup <- setNames(supIdx, supTags) - # Group sub characters by tag - tagToSubs <- split(subIdx, subTags) + # The longest sup_ tag that `tag` sits under: an exact match, or an extension + # at an underscore boundary (so "tailfin" does NOT sit under "tail"). + # `exclude` stops a sup_ tag being its own parent. + .ParentTag <- function(tag, exclude = "") { + eligible <- supTags[supTags != exclude] + if (!length(eligible)) return(NA_character_) + under <- tag == eligible | startsWith(tag, paste0(eligible, "_")) + if (!any(under)) return(NA_character_) + cand <- eligible[under] + cand[[which.max(nchar(cand))]] + } - # Check for sub_ characters referencing nonexistent sup_ tags - orphanTags <- setdiff(names(tagToSubs), supTags) - if (length(orphanTags) > 0L) { + subParent <- vapply(subTagsFull, .ParentTag, character(1), USE.NAMES = FALSE) + orphan <- is.na(subParent) + if (any(orphan)) { warning(sprintf( "sub_ characters reference tags with no corresponding sup_: %s", - paste(orphanTags, collapse = ", ") + paste(unique(subTagsFull[orphan]), collapse = ", ") )) } - # Detect nested hierarchies: a sub_ character that is also a sup_ - # Find sub_ chars that are also in supIdx - subAlsoSup <- intersect(subIdx, supIdx) - - # Build hierarchy - # First pass: create flat blocks for all sup_ tags - args <- list() - for (tag in supTags) { - ctrl <- tagToSup[[tag]] - subs <- tagToSubs[[tag]] - if (is.null(subs)) subs <- integer(0) - - # Check which subs are themselves controlling (nested hierarchy) - nestedSubs <- intersect(subs, supIdx) - flatSubs <- setdiff(subs, supIdx) - - if (length(nestedSubs) == 0L) { - # Simple block - args[[as.character(ctrl)]] <- as.integer(subs) - } else { - # Nested: build list with named sub-hierarchies - block <- as.list(as.integer(flatSubs)) - for (ns in nestedSubs) { - nsTag <- supTags[supIdx == ns] - nsSubs <- tagToSubs[[nsTag]] - if (is.null(nsSubs)) nsSubs <- integer(0) - block[[as.character(ns)]] <- as.integer(nsSubs) - } - args[[as.character(ctrl)]] <- block + # A sup_ tag that extends another sup_ tag is a nested controller. + supParent <- vapply(supTags, function(s) .ParentTag(s, exclude = s), + character(1), USE.NAMES = FALSE) + + # Dependents of one tag: its own sub_ characters, plus a named sub-hierarchy + # for each sup_ tag nested directly beneath it. Recursive, so nesting works + # to arbitrary depth. A nested controller is NOT added to the flat dependents + # here -- .ParseOneBlock() records a named sub-controller as a dependent of the + # enclosing block itself. + .BuildTag <- function(tag) { + flat <- as.integer(subIdx[!orphan & subParent == tag]) + kidTags <- supTags[!is.na(supParent) & supParent == tag] + if (!length(kidTags)) { + return(flat) + } + block <- as.list(flat) + for (k in kidTags) { + block[[as.character(tagToSup[[k]])]] <- .BuildTag(k) } + block } - # Filter out sup_ chars whose index also appears in subIdx - # (they'll be included as children of their parent) - topLevelSup <- setdiff(supIdx, subIdx) - if (length(topLevelSup) == 0L) { - # All sup_ characters are also sub_ — circular or all nested. - # Fall back to treating all as top-level with a warning. - warning("All sup_ characters are also sub_ characters. ", + topTags <- supTags[is.na(supParent)] + if (!length(topTags)) { + # Every sup_ tag extends another, which needs a cycle and so cannot arise + # from prefix matching; kept as a guard rather than a reachable branch. + warning("Every sup_ tag extends another sup_ tag. ", "Treating all as top-level.") - topLevelSup <- supIdx + topTags <- supTags } - topLevelCtrls <- as.character(topLevelSup) - args <- args[topLevelCtrls] + + args <- lapply(topTags, .BuildTag) + names(args) <- as.character(tagToSup[topTags]) do.call(CharacterHierarchy, args) } diff --git a/man/HierarchyFromNames.Rd b/man/HierarchyFromNames.Rd index 436e0f622..06456751a 100644 --- a/man/HierarchyFromNames.Rd +++ b/man/HierarchyFromNames.Rd @@ -16,15 +16,28 @@ detected. \description{ Parse character names following the TNT convention where controlling characters are named \verb{sup_} and their dependent characters are -named \verb{sub_[_suffix]}. Tags must match between a controlling -character and its dependents. Nested hierarchies are detected when a -\code{sub_} character is also a \code{sup_} for further characters. +named \verb{sub_[_suffix]}. Each dependent is attached to the longest +\code{sup_} tag that its own tag extends, at an underscore boundary. +} +\details{ +Nested hierarchies are written by giving a controlling character a tag that +itself extends another controlling character's tag: \code{sup_tail_tip} controls +the tag \code{tail_tip}, and because \code{tail_tip} extends \code{tail} it is also a +dependent of \code{sup_tail}. Longest-match is what keeps \code{sub_tail_tip_gloss} +with \code{sup_tail_tip} rather than with \code{sup_tail}. A shared prefix alone does +not nest: \code{sup_tailfin} is independent of \code{sup_tail}, since \code{tailfin} does +not extend \code{tail} at an underscore boundary. } \examples{ names <- c("sup_tail", "sub_tail_colour", "sub_tail_shape", "sup_wing", "sub_wing_venation", "eyes") HierarchyFromNames(names) +# Nesting: `tail_tip` extends `tail`, so character 3 is both a dependent of +# character 1 and the controlling character of character 4. +HierarchyFromNames(c("sup_tail", "sub_tail_colour", + "sup_tail_tip", "sub_tail_tip_gloss")) + } \seealso{ \code{\link[=CharacterHierarchy]{CharacterHierarchy()}} diff --git a/src/ts_hsj.cpp b/src/ts_hsj.cpp index 9c68a193f..436c2252e 100644 --- a/src/ts_hsj.cpp +++ b/src/ts_hsj.cpp @@ -128,21 +128,44 @@ static CanonOrder build_canon_order(const TreeState& tree) { // which is what state_sets must hold to make the Fitch downpass/uppass below // correct for ambiguous tokens. // -// `pri_free[t]` marks the tips at which this secondary carries no constraint: -// those whose controlling primary CANNOT code the structure present, so the -// character does not exist there and its "-" is not a state it takes (T-374). -// score_hierarchy_block() computes it and documents why the test is that -// strict one rather than "may be absent". Admitting "-" as an ordinary -// concrete state -- as this function formerly did, and as the comment here -// formerly asserted was deliberate -- let the uppass propagate it INWARDS and -// resolve a node in the middle of the PRESENT region to it, where it is -// disjoint from every present neighbour in every secondary at once, and -// score_hierarchy_block() charged that branch d = m, the full alpha, for a -// node that by construction has no inapplicable secondaries. That over-charge -// is wrong under any rooting (the paper's d counts "nonmatching secondary -// characters", p.5, among characters that APPLY), and because whether it fired -// depended on the DELTRAN direction it was also the dominant source of -// T-374's rooting-dependence. +// THE INAPPLICABLE TOKEN IS NOT A STATE OF A SECONDARY CHARACTER. This is the +// paper's central claim, not a convention we are free to pick: Hopkins & St John +// (2021) define d as "the number of nonmatching secondary characters dependent +// on that primary" (p.5) among the characters that APPLY, and state plainly that +// where secondaries are inapplicable to a taxon "they have no influence on the +// estimated dissimilarity" (p.5). The paper introduces HSJ precisely to avoid +// the alternative: "Treating inapplicable characters as a new, separate state +// will ... skew the analysis, because having a new separate state increases the +// dissimilarity of all pairwise comparisons ... This results in overweighting +// the [controlling] primary character and favors clades that separate taxa with +// secondary characters from those without" (p.5). Admitting "-" as a concrete +// state is therefore not a stricter reading of HSJ; it is the FitchS behaviour +// HSJ exists to replace. +// +// So `inapp_bit` is stripped from every secondary's observed set, at every tip, +// regardless of what the controlling primary codes. Two distinct routes reach a +// non-constraining tip and both must be handled here: +// * `pri_free[t]` -- the controlling primary CANNOT code the structure present, +// so the character does not exist at t (T-374). score_hierarchy_block() +// computes this and documents why the test is that strict one rather than +// the laxer "may be absent". +// * the secondary's own token carries no applicable state once "-" is removed +// -- i.e. the cell is coded "-" (or an ambiguity resolving only to "-") +// while the primary does NOT certainly code absence. That combination is +// contradictory coding -- ValidateHierarchy() enforces the converse +// direction (no applicable secondary where the primary codes absence) but +// not this one -- and whatever a validator decides to do about it, the score +// must not invent a state for it (T-396). Before this, such a tip +// kept a single concrete inapp state that contributed tie-break support and +// that the uppass could propagate INWARDS, resolving a node in the middle of +// the PRESENT region to it -- disjoint from every present neighbour in every +// secondary at once, so score_hierarchy_block() charged that branch d = m, +// the full alpha, for a node that by construction has no inapplicable +// secondaries. Measured over-charge was exactly alpha/m per affected +// branch, and it exceeded every concrete resolution as well as the missing +// treatment, so no reading of the data made it right. +// Because whether the old over-charge fired depended on the DELTRAN direction, +// it was also the dominant source of T-374's rooting-dependence. static int fitch_label_char( const TreeState& tree, const std::vector& tip_labels, @@ -150,6 +173,7 @@ static int fitch_label_char( int n_orig_chars, const std::vector& token_states, int n_levels, + uint32_t inapp_bit, const std::vector& pri_free, const CanonOrder& co, std::vector& state_sets) @@ -157,13 +181,20 @@ static int fitch_label_char( int n_tip = tree.n_tip; int n_node = tree.n_node; + // Applicable state set of this secondary at tip t: the states its token + // denotes, minus the inapplicable state, which is not one of them (above). + // Zero means the cell constrains nothing. + auto applicable_at = [&](int t) -> uint32_t { + return token_states[tip_labels[t * n_orig_chars + char_idx]] & ~inapp_bit; + }; + // The applicable domain: the states this character is observed in at tips // where it actually applies. Wildcarding to this rather than to all // n_levels bits keeps the tie-break arrays below as small as they were. uint32_t domain = 0; for (int t = 0; t < n_tip; ++t) { if (!pri_free[t]) { - domain |= token_states[tip_labels[t * n_orig_chars + char_idx]]; + domain |= applicable_at(t); } } // The character applies nowhere: it constrains nothing. Give every node one @@ -173,8 +204,8 @@ static int fitch_label_char( uint32_t used_mask = 0; std::vector observed(n_tip); for (int t = 0; t < n_tip; ++t) { - int label = tip_labels[t * n_orig_chars + char_idx]; - observed[t] = pri_free[t] ? domain : token_states[label]; + const uint32_t applicable = applicable_at(t); + observed[t] = (pri_free[t] || applicable == 0) ? domain : applicable; state_sets[t] = observed[t]; used_mask |= observed[t]; } @@ -401,7 +432,8 @@ static double score_hierarchy_block( std::vector buf(n_node); for (int j = 0; j < m; ++j) { fitch_label_char(tree, tip_labels, block.secondary_chars[j], - n_orig_chars, token_states, n_levels, pri_free, co, buf); + n_orig_chars, token_states, n_levels, inapp_bit, + pri_free, co, buf); for (int nd = 0; nd < n_node; ++nd) { sec_states[j * n_node + nd] = buf[nd]; } diff --git a/tests/testthat/test-CharacterHierarchy.R b/tests/testthat/test-CharacterHierarchy.R index 96a483232..17cd66d63 100644 --- a/tests/testthat/test-CharacterHierarchy.R +++ b/tests/testthat/test-CharacterHierarchy.R @@ -152,7 +152,7 @@ test_that(".NonHierarchyWeights subtracts hierarchy chars", { h <- CharacterHierarchy("1" = 2L) w_orig <- attr(ds, "weight") - w_adj <- .NonHierarchyWeights(ds, h) + w_adj <- TreeSearch:::.NonHierarchyWeights(ds, h) # Adjusted weights should be non-negative @@ -171,7 +171,7 @@ test_that(".BuildTipLabels creates correct matrix", { ds <- phangorn::phyDat(mat, type = "USER", levels = c("-", "0", "1"), ambiguity = "?") - tl <- .BuildTipLabels(ds) + tl <- TreeSearch:::.BuildTipLabels(ds) expect_equal(nrow(tl), 3L) expect_equal(ncol(tl), 3L) # Values should be 0-based token indices @@ -180,7 +180,7 @@ test_that(".BuildTipLabels creates correct matrix", { test_that(".HierarchyToBlocks converts to 0-based flat list", { h <- CharacterHierarchy("1" = 2:4, "5" = 6:7) - blocks <- .HierarchyToBlocks(h) + blocks <- TreeSearch:::.HierarchyToBlocks(h) expect_length(blocks, 2) expect_equal(blocks[[1]]$primary, 0L) expect_equal(blocks[[1]]$secondaries, 1:3) @@ -190,7 +190,7 @@ test_that(".HierarchyToBlocks converts to 0-based flat list", { test_that(".HierarchyToBlocks flattens nested hierarchies", { h <- CharacterHierarchy("1" = list(2, 4, "3" = 9:10)) - blocks <- .HierarchyToBlocks(h) + blocks <- TreeSearch:::.HierarchyToBlocks(h) expect_gte(length(blocks), 2) # First block: primary=0, secondaries should include 1 and 3 (chars 2 and 4) expect_equal(blocks[[1]]$primary, 0L) @@ -213,7 +213,7 @@ test_that(".NonHierarchyWeights preserves non-hierarchy patterns", { h <- CharacterHierarchy("1" = 2L) idx <- attr(ds, "index") w_orig <- attr(ds, "weight") - w_adj <- .NonHierarchyWeights(ds, h) + w_adj <- TreeSearch:::.NonHierarchyWeights(ds, h) # Character 3 is not in the hierarchy; its pattern should keep its weight # unless it shares a pattern with a hierarchy character @@ -225,3 +225,127 @@ test_that(".NonHierarchyWeights preserves non-hierarchy patterns", { expect_gte(w_adj[pat], 0L) } }) + + +# ========================================================================= +# Nested hierarchies: the documented example must actually validate (T-395) +# ========================================================================= +# .ParseOneBlock() deliberately records a sub-controller BOTH as a dependent of +# its parent and as the controlling character of its own block -- that dual role +# is the entire content of "nested". .ValidateBlock() used to count the second +# occurrence as a double claim, so every nested hierarchy that could be written +# was rejected, including this file's own roxygen example. R CMD check could not +# see it because CharacterHierarchy() itself never validates: the example +# constructs fine and is never scored. + +test_that("a nested hierarchy validates, and the documented example works", { + # Char 1 controls {2, 3, 4, 5}; char 3 additionally controls {9, 10}. + nested <- CharacterHierarchy("1" = list(2, 3, 4, 5, "3" = 9:10)) + + expect_equal(nested[[1]]$controlling, 1L) + expect_equal(sort(nested[[1]]$dependents), c(2L, 3L, 4L, 5L)) + expect_equal(length(nested[[1]]$children), 1L) + expect_equal(nested[[1]]$children[[1]]$controlling, 3L) + expect_equal(nested[[1]]$children[[1]]$dependents, 9:10) + + # Coding invariants: secondaries "-" where their controller codes absence; + # char 3 binary where it applies; chars 9-10 "-" wherever char 3 is not "1". + mat <- rbind( + t1 = c("1", "0", "1", "0", "1", "0", "1", "0", "0", "1"), + t2 = c("1", "1", "1", "1", "0", "1", "0", "1", "1", "1"), + t3 = c("1", "0", "0", "1", "1", "0", "1", "1", "-", "-"), + t4 = c("0", "-", "-", "-", "-", "1", "0", "0", "-", "-"), + t5 = c("0", "-", "-", "-", "-", "0", "1", "1", "-", "-"), + t6 = c("0", "-", "-", "-", "-", "1", "1", "0", "-", "-")) + ds <- phangorn::phyDat(mat, type = "USER", levels = c("-", "0", "1"), + ambiguity = "?") + + expect_silent(ValidateHierarchy(nested, ds)) + + # A genuine double claim must STILL be rejected -- the fix must not have + # simply disabled the check. Char 2 is claimed by both blocks here. + expect_error( + ValidateHierarchy(CharacterHierarchy("1" = 2:3, "5" = c(2L, 6L)), ds), + "multiple hierarchy blocks" + ) + # And a nested block's own dependents are still checked against other blocks. + expect_error( + ValidateHierarchy(CharacterHierarchy("1" = list(2, "3" = 9:10), "5" = 9L), + ds), + "multiple hierarchy blocks" + ) +}) + +test_that("a nested hierarchy scores, and stays rooting-invariant", { + nested <- CharacterHierarchy("1" = list(2, 3, 4, 5, "3" = 9:10)) + mat <- rbind( + t1 = c("1", "0", "1", "0", "1", "0", "1", "0", "0", "1"), + t2 = c("1", "1", "1", "1", "0", "1", "0", "1", "1", "1"), + t3 = c("1", "0", "0", "1", "1", "0", "1", "1", "-", "-"), + t4 = c("0", "-", "-", "-", "-", "1", "0", "0", "-", "-"), + t5 = c("0", "-", "-", "-", "-", "0", "1", "1", "-", "-"), + t6 = c("0", "-", "-", "-", "-", "1", "1", "0", "-", "-")) + ds <- phangorn::phyDat(mat, type = "USER", levels = c("-", "0", "1"), + ambiguity = "?") + tree <- Renumber(RenumberTips( + ape::read.tree(text = "(((t1,t2),t3),(t4,(t5,t6)));"), rownames(mat))) + + score <- TreeLength(tree, ds, hierarchy = nested, inapplicable = "hsj", + hsj_alpha = 1) + expect_true(is.finite(score)) + expect_gt(score, 0) + + # The HSJ score is a minimum over labellings of a sum of SYMMETRIC + # dissimilarities on an unrooted tree, so it cannot depend on the rooting. + rooted <- vapply(rownames(mat), function(tip) { + TreeLength(Renumber(RootTree(tree, tip)), ds, hierarchy = nested, + inapplicable = "hsj", hsj_alpha = 1) + }, double(1)) + expect_equal(diff(range(rooted)), 0) + + # The x-transformation genuinely does not implement nesting; its error must + # be the informative one, not a double-claim complaint from the validator. + expect_error(RecodeHierarchy(ds, nested), "[Nn]ested") +}) + +test_that("HierarchyFromNames detects nesting by tag extension", { + # Flat, documented case: unchanged. + flat <- HierarchyFromNames(c("sup_tail", "sub_tail_colour", "sub_tail_shape", + "sup_wing", "sub_wing_venation", "eyes")) + expect_equal(length(flat), 2L) + expect_equal(sort(flat[[1]]$dependents), c(2L, 3L)) + expect_equal(length(flat[[1]]$children), 0L) + + # Nested: `tail_tip` extends `tail`, so char 3 is both a dependent of char 1 + # and the controller of char 4. The old first-component tag match collapsed + # every depth onto the outermost tag, so this returned two flat blocks -- one + # of them empty -- and the nesting was silently lost. + nested <- HierarchyFromNames(c("sup_tail", "sub_tail_colour", + "sup_tail_tip", "sub_tail_tip_gloss")) + expect_equal(length(nested), 1L) + expect_equal(nested[[1]]$controlling, 1L) + expect_equal(sort(nested[[1]]$dependents), c(2L, 3L)) + expect_equal(length(nested[[1]]$children), 1L) + expect_equal(nested[[1]]$children[[1]]$controlling, 3L) + expect_equal(nested[[1]]$children[[1]]$dependents, 4L) + + # Longest-match, not first-match: the deeper dependent must not be captured + # by the shallower tag. + expect_false(4L %in% nested[[1]]$dependents) + + # A shared prefix without an underscore boundary is NOT nesting. + fin <- HierarchyFromNames(c("sup_tail", "sup_tailfin", "sub_tailfin_x")) + expect_equal(length(fin), 2L) + + # Three levels deep. + deep <- HierarchyFromNames(c("sup_a", "sub_a_x", "sup_a_b", "sub_a_b_y", + "sup_a_b_c", "sub_a_b_c_z")) + expect_equal(length(deep), 1L) + expect_equal(deep[[1]]$children[[1]]$controlling, 3L) + expect_equal(deep[[1]]$children[[1]]$children[[1]]$controlling, 5L) + expect_equal(deep[[1]]$children[[1]]$children[[1]]$dependents, 6L) + + # An orphan sub_ still warns. + expect_warning(HierarchyFromNames(c("sup_tail", "sub_nose_shape")), + "no corresponding sup_") +}) diff --git a/tests/testthat/test-ts-hsj.R b/tests/testthat/test-ts-hsj.R index c25c60470..a4d31a483 100644 --- a/tests/testthat/test-ts-hsj.R +++ b/tests/testthat/test-ts-hsj.R @@ -834,14 +834,32 @@ test_that("HSJ secondary '?' obeys the resolution invariant (T-375)", { # collapses the whole score to the secondary's ordinary Fitch step count. # # Tree ((t1,t2),(t3,t4)); primaries all "1"; secondary t1=t2=t3="0", t4 - # varies. Hand-derived: t4="0" ties all four -> 0 steps. t4="1" or t4="-" - # each disagree with the (t3,t4) clade's neighbour -> 1 step (the downpass - # intersect((t3=0),(t4=1 or -)) is empty, forcing a union). t4="?" must - # resolve to whichever concrete state is compatible AND cheapest -- here - # that's "0" (matching t1/t2/t3), giving 0 steps, so score("?") == 0 == - # min(0, 1, 1). Before the fix, fitch_label_char() bit-encoded the "?" - # TOKEN index as its own concrete state bit, indistinguishable from a - # genuine mismatch, and scored 1 -- violating the invariant (1 > 0). + # varies. Hand-derived: t4="0" ties all four -> 0 steps. t4="1" disagrees with + # the (t3,t4) clade's neighbour -> 1 step (the downpass intersect((t3=0), + # (t4=1)) is empty, forcing a union). t4="?" must resolve to whichever concrete + # state is compatible AND cheapest -- here that's "0" (matching t1/t2/t3), + # giving 0 steps. Before the fix, fitch_label_char() bit-encoded the "?" TOKEN + # index as its own concrete state bit, indistinguishable from a genuine + # mismatch, and scored 1 -- violating the invariant (1 > 0). + # + # t4="-" MUST ALSO SCORE 0, and this expectation CHANGED on 2026-08-03 (T-396): + # it previously asserted 1, reasoning that "-" disagrees with its neighbour + # just as "1" does. That reasoning treats the inapplicable token as an ordinary + # third state, which is precisely the behaviour Hopkins & St John (2021) + # introduce HSJ to avoid -- "Treating inapplicable characters as a new, + # separate state[] will ... skew the analysis, because having a new separate + # state increases the dissimilarity of all pairwise comparisons" (p.5). The + # paper's d counts "nonmatching secondary characters" among those that APPLY, + # and where secondaries are inapplicable to a taxon "they have no influence on + # the estimated dissimilarity" (p.5). An inapplicable secondary therefore + # contributes nothing to d whatever its controlling primary codes. + # + # Note every primary here is "1" (present), so this is the CONTRADICTORY coding + # case: "-" in a secondary whose controlling primary says the structure exists. + # The old value was not a defensible alternative reading -- it exceeded every + # concrete resolution AND the missing-data treatment, so no resolution of the + # cell produced it. Whether such a matrix should be rejected by + # ValidateHierarchy() is a separate question; the score must not invent a state. h <- CharacterHierarchy("1" = 2L) tree <- Renumber(RenumberTips( ape::read.tree(text = "((t1,t2),(t3,t4));"), paste0("t", 1:4))) @@ -854,8 +872,13 @@ test_that("HSJ secondary '?' obeys the resolution invariant (T-375)", { } scores <- vapply(c("0", "1", "-", "?"), score_for, double(1)) - expect_equal(unname(scores), c(0, 1, 1, 0)) + expect_equal(unname(scores), c(0, 1, 0, 0)) expect_lte(scores[["?"]], min(scores[c("0", "1", "-")])) + # An inapplicable secondary must cost no more than the cheapest concrete + # resolution of that cell -- the paper's "no influence" property. Asserted + # separately from the exact vector so a future change to the matrix still + # checks the property rather than only the number. + expect_lte(scores[["-"]], min(scores[c("0", "1")])) })