From 0f7241d0e5c85542cda948a26d40112a8accea21 Mon Sep 17 00:00:00 2001 From: hyperpolymath <6759885+hyperpolymath@users.noreply.github.com> Date: Fri, 23 Oct 2026 07:56:08 +0000 Subject: [PATCH 1/2] Merge PR #224 validated-statistics + cladistics fixes onto main, resolve conflicts, fix review findings MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Bring PR #224's resolved #221 changes (BigInt validated run statistics, cladistics/CladeCumulus repairs, CI failure digest) onto current main, handling the four overlapping files: - validated_statistics.jl: keep SQL projection + ordering and BigInt accumulation; adopt #221's pre-allocation MAX_CELLS cell check moved inside the accumulation so an oversized result is refused before it is allocated. - clade_cumulus.jl route: keep main's AbstractString boundary (accepts SubString request fields) and duckdb_tables() view filter; the PR's stricter ::String signature would MethodError on SubString study/run. - clade_cumulus.jl analysis: keep main's correct QueryCandidate field names (distal_len/pendant_len/like_weight); fix the weight-sum check to be asymmetric per Codex P2 โ€” tolerate a lower shortfall down to 1 - PLACEMENT_WEIGHT_SUM_TOL but reject any overweight beyond a tight PLACEMENT_WEIGHT_OVER_TOL bound. - ci.yml: publish a digest of every failing assertion (testset, expression, evaluated values) before the log tail so the Checks API shows the full failure list. - phylogeny.jl: import ..Tools so PH.Tools._tools resolves. - trees.jl: _tree_format/_valid_tree_file accept AbstractString for SubString safety. - server.jl: include cladistics routes once. - tests: zero count in a sample with reads is a defined 0 proportion; fractional JSON counts refused; helper view unregistered after use. Files intentionally not taken from the PR head: src/core/doi.jl (PR drops the `module DOI` wrapper, which is invalid top-level Julia) and src/server/routes/clade_cumulus.jl (main's version is the stricter fix). Co-authored-by: arena-agent <297053741+arena-agent@users.noreply.github.com> --- .github/workflows/ci.yml | 21 +++++++-- src/analysis/clade_cumulus.jl | 15 ++++--- src/pipeline/phylogeny.jl | 1 + src/server/routes/trees.jl | 4 +- src/server/routes/validated_statistics.jl | 53 +++++++++++++---------- src/server/server.jl | 1 - test/unit/test_clade_cumulus.jl | 12 +++-- test/unit/test_exact_summaries.jl | 10 ++--- test/unit/test_validated_statistics.jl | 4 ++ 9 files changed, 79 insertions(+), 42 deletions(-) diff --git a/.github/workflows/ci.yml b/.github/workflows/ci.yml index d51aa441..5fa55493 100644 --- a/.github/workflows/ci.yml +++ b/.github/workflows/ci.yml @@ -400,11 +400,26 @@ jobs: shell: bash run: | set +e - julia --project=. -t2 --compiled-modules=no test/runtests.jl numeric_policy exact_summaries high_precision derived_analyses validated_statistics 2>&1 | tee /tmp/julia-validated-tests.log + julia --project=. -t2 --compiled-modules=no test/runtests.jl numeric_policy exact_summaries high_precision derived_analyses validated_statistics 2>&1 | tee /tmp/validated-statistics-tests.log status="${PIPESTATUS[0]}" if [ "$status" -ne 0 ]; then - enc() { head -c 3500 | sed -e 's/%/%25/g' -e 's/\r/%0D/g' | awk '{printf "%s%%0A", $0}'; } - echo "::error title=Validated statistics failed (tail)::$(tail -c 3500 /tmp/julia-validated-tests.log | enc)" + enc() { head -c 3000 | sed -e 's/%/%25/g' -e 's/\r/%0D/g' | awk '{printf "%s%%0A", $0}'; } + # One annotation cannot hold a failing run, and a tail hides every + # failure but the last one. Publish a digest of all failing + # assertions first - the testset, the expression and the evaluated + # values are what identify the bug - then the end of the log. + grep -E "Test Failed at |Error During Test at |^ Expression: |^ Evaluated: |^ Test threw exception" \ + /tmp/validated-statistics-tests.log > /tmp/validated-statistics-failures.log || true + rm -f /tmp/validated-statistics-digest-* + split -b 2500 /tmp/validated-statistics-failures.log /tmp/validated-statistics-digest- || true + index=0 + for chunk in /tmp/validated-statistics-digest-*; do + [ -e "$chunk" ] || continue + index=$((index + 1)) + [ "$index" -le 8 ] || break + echo "::error title=Validated-statistics failure digest $index::$(enc < "$chunk")" + done + echo "::error title=Validated-statistics log tail::$(tail -c 3000 /tmp/validated-statistics-tests.log | enc)" fi exit "$status" diff --git a/src/analysis/clade_cumulus.jl b/src/analysis/clade_cumulus.jl index 63cbdf24..9ac3f4f0 100644 --- a/src/analysis/clade_cumulus.jl +++ b/src/analysis/clade_cumulus.jl @@ -28,11 +28,15 @@ const DEFAULT_TAXONOMIC_RANKS = [ "Class", "Order", "Family", "Genus", "Species", ] -# Placement tools may truncate a query's candidate list, so retained -# like-weight ratios can legitimately sum to slightly less than 1 (a single -# retained placement at 0.98 is complete evidence). Sums far from 1 indicate -# a malformed or partial record and stay undecided. +# Placement tools truncate a query's candidate list at a likelihood-weight +# cutoff, so retained like-weight ratios can legitimately sum to a little +# under 1 (a single retained placement at 0.98 is complete evidence). A shortfall +# beyond the tolerance is a missing or discarded placement and stays undecided. +# A sum over 1 is never legitimate: weights are mutually exclusive placement +# ratios, so overweight evidence is malformed and is rejected by a tight upper +# bound while the lower shortfall is tolerated. const PLACEMENT_WEIGHT_SUM_TOL = 0.05 +const PLACEMENT_WEIGHT_OVER_TOL = 1e-9 ## Phylogenetic tree structures & Newick / jplace parsers mutable struct PhyloNode @@ -573,7 +577,8 @@ function parse_jplace(raw_json::AbstractString)::JplaceDocument if !isempty(bad_reason) push!(queries, QueryPlacement(q_ord, q_names, QueryCandidate[], :undecided, bad_reason, nothing)) - elseif abs(weight_sum - 1.0) > PLACEMENT_WEIGHT_SUM_TOL + elseif weight_sum > 1.0 + PLACEMENT_WEIGHT_OVER_TOL || + weight_sum < 1.0 - PLACEMENT_WEIGHT_SUM_TOL push!(queries, QueryPlacement(q_ord, q_names, cands, :undecided, "weight_sum_out_of_tolerance ($(round(weight_sum; digits=6)))", nothing)) else diff --git a/src/pipeline/phylogeny.jl b/src/pipeline/phylogeny.jl index b394d9dc..b978ef02 100644 --- a/src/pipeline/phylogeny.jl +++ b/src/pipeline/phylogeny.jl @@ -34,6 +34,7 @@ export REFERENCE_STEPS, PLACEMENT_STEPS, Step, placement_settings, run_workflow, using ..PipelineLog using ..Validation using ..Tools: tool_bin, _sq, _run_logged, _safe_optional_args, _num + import ..Tools using ..RemoteExec: run_remote, checked_target import ..Provenance diff --git a/src/server/routes/trees.jl b/src/server/routes/trees.jl index f3cc8e07..beec713a 100644 --- a/src/server/routes/trees.jl +++ b/src/server/routes/trees.jl @@ -17,9 +17,9 @@ const _tree_view_lock = ReentrantLock() _trees_dir(study::String) = joinpath(ServerState.projects_dir(), study, "trees") _tree_view_path(study::String, file::String) = joinpath(_trees_dir(study), ".views", file * ".json") -_tree_format(file::String) = endswith(lowercase(file), ".jplace") ? "jplace" : "newick" +_tree_format(file::AbstractString) = endswith(lowercase(file), ".jplace") ? "jplace" : "newick" -function _valid_tree_file(file::String)::Bool +function _valid_tree_file(file::AbstractString)::Bool _valid_name(file) && any(ext -> endswith(lowercase(file), ext), _TREE_EXTENSIONS) end diff --git a/src/server/routes/validated_statistics.jl b/src/server/routes/validated_statistics.jl index fb14d584..e6cf4ffd 100644 --- a/src/server/routes/validated_statistics.jl +++ b/src/server/routes/validated_statistics.jl @@ -310,8 +310,6 @@ function _validated_integer_sample_columns(con, table::String, columns::Vector{S samples end -"""Read and aggregate an integer-safe run table. Only GROUP BY is done in SQL; -all count addition is Julia BigInt. No SQL SUM and no Float64 matrix is used.""" function _reject_float_count_filters(sample_columns::Vector{String}, params::Dict{String,String}) for sample in sample_columns, prefix in ("col_min.", "col_max.") key = prefix * sample @@ -322,6 +320,15 @@ function _reject_float_count_filters(sample_columns::Vector{String}, params::Dic nothing end +"""Read and aggregate an integer-safe run table without SQL arithmetic. + +SQL only projects rows; every count is added in Julia `BigInt`. A `SUM` of +`BIGINT`/`HUGEINT` columns is evaluated in `Int64`/`Int128` and can overflow +before the exact total is known - two sample columns of `2^126` already sum to +`2^127`, one past `typemax(Int128)` - so no `SUM` and no `Float64` matrix is +used. Rows that share a feature label are therefore summed here, after SQL has +ordered them by label. +""" function _validated_duckdb_counts(con, table::String, columns::Vector{String}; prefix=nothing, taxonomic_level=nothing, where_clause::String="", where_params=Any[]) @@ -348,45 +355,47 @@ function _validated_duckdb_counts(con, table::String, columns::Vector{String}; "COALESCE(NULLIF(CAST($(_sql_identifier(taxon)) AS VARCHAR), ''), 'Unassigned')" end count_columns = join([_sql_identifier(c) for c in samples], ", ") - # Rows sharing a feature (e.g. several sequences in one genus) are summed - # here, in BigInt, not in SQL: DuckDB refuses bare count columns beside - # GROUP BY, and SUM over HUGEINT counts would overflow exactly the values - # this path exists to keep exact. ORDER BY keeps features sorted and each - # feature's rows adjacent. + # Project the rows and let SQL sort them by feature; no addition happens + # here. `GROUP BY feature` alone is a DuckDB binder error, because the count + # columns are neither grouped nor aggregated, and `SUM` cannot take their + # place without overflowing. Sorting by label makes equal labels adjacent + # and leaves the accumulator independent of row order. sql = "SELECT $feature_expression AS feature, $count_columns FROM $(_sql_identifier(table)) " * "$where_clause ORDER BY feature LIMIT $(MAX_CELLS + 1)" frame = DataFrame(DBInterface.execute(con, sql, where_params)) - nrow(frame) <= MAX_CELLS || throw(ResourceLimitError("rows", MAX_CELLS, nrow(frame))) nrow(frame) > 0 || throw(ArgumentError("no rows remain after filtering")) + nrow(frame) <= MAX_CELLS || throw(ResourceLimitError("table rows", MAX_CELLS, nrow(frame))) features = String[] - feature_rows = Vector{BigInt}[] + totals = Dict{String,Vector{BigInt}}() for (i, row) in enumerate(eachrow(frame)) ismissing(row.feature) && throw(UnsupportedRepresentationError( "a non-missing feature label", "missing", "validated run statistics row $i")) - feature = string(row.feature) - if isempty(features) || features[end] != feature - push!(features, feature) - push!(feature_rows, zeros(BigInt, length(samples))) + label = string(row.feature) + sums = get(totals, label, nothing) + if isnothing(sums) + length(features) + 1 <= MAX_CELLS รท max(length(samples), 1) || + throw(ResourceLimitError("table cells", MAX_CELLS, + (length(features) + 1) * length(samples))) + sums = zeros(BigInt, length(samples)) + totals[label] = sums + push!(features, label) end - target = feature_rows[end] for (j, sample) in enumerate(samples) value = row[Symbol(sample)] ismissing(value) && throw(UnsupportedRepresentationError( - "a non-missing integer count", "NULL", "feature '$feature', sample '$sample'")) + "a non-missing integer count", "NULL", "feature '$label', sample '$sample'")) (value isa Integer && !(value isa Bool)) || throw(UnsupportedRepresentationError( - "an integer count", repr(value), "feature '$feature', sample '$sample'")) + "an integer count", repr(value), "feature '$label', sample '$sample'")) value >= 0 || throw(ArgumentError( - "negative count $value at feature '$feature', sample '$sample'")) - target[j] += BigInt(value) # every addition is exact BigInt + "negative count $value at feature '$label', sample '$sample'")) + sums[j] += BigInt(value) # every addition is exact, so no overflow end end counts = Matrix{BigInt}(undef, length(features), length(samples)) - for (i, feature_row) in enumerate(feature_rows) - counts[i, :] = feature_row + for (i, label) in enumerate(features) + counts[i, :] = totals[label] end - length(features) * length(samples) <= MAX_CELLS || throw(ResourceLimitError( - "table cells", MAX_CELLS, length(features) * length(samples))) ExactSummaries.validate_coordinates(counts, samples, features) counts, samples, features end diff --git a/src/server/server.jl b/src/server/server.jl index ebf5abb6..eca368a0 100644 --- a/src/server/server.jl +++ b/src/server/server.jl @@ -118,7 +118,6 @@ module Server include(joinpath(@__DIR__, "routes", "phylogeny.jl")) include(joinpath(@__DIR__, "routes", "cladistics.jl")) include(joinpath(@__DIR__, "routes", "clade_cumulus.jl")) - include(joinpath(@__DIR__, "routes", "cladistics.jl")) ## Bad request middleware # Helpers deep in a route (e.g. `_req_group`) throw this to reject a request. diff --git a/test/unit/test_clade_cumulus.jl b/test/unit/test_clade_cumulus.jl index 6e12a0c9..51d9d041 100644 --- a/test/unit/test_clade_cumulus.jl +++ b/test/unit/test_clade_cumulus.jl @@ -252,8 +252,8 @@ done echo "(MMTIP_000001:0.05,(MMTIP_000002:0.03,MMTIP_000003:0.04)0.910:0.12,MMTIP_000004:0.08);" > "$OUT" """) chmod(mock_ft, 0o755) - old_ft = get(MetaManifold.Tools._tools, "fasttree", nothing) - MetaManifold.Tools._tools["fasttree"] = mock_ft + old_ft = get(PH.Tools._tools, "fasttree", nothing) + PH.Tools._tools["fasttree"] = mock_ft try s4 = [ "ASV_1" => "ACGTACGTACGTACGTACGT", @@ -267,9 +267,9 @@ echo "(MMTIP_000001:0.05,(MMTIP_000002:0.03,MMTIP_000003:0.04)0.910:0.12,MMTIP_0 r["sha256"] == bytes2hex(open(sha256, mock_ft)), ft_res.tool_records) finally if old_ft === nothing - delete!(MetaManifold.Tools._tools, "fasttree") + delete!(PH.Tools._tools, "fasttree") else - MetaManifold.Tools._tools["fasttree"] = old_ft + PH.Tools._tools["fasttree"] = old_ft end end end @@ -310,6 +310,10 @@ echo "(MMTIP_000001:0.05,(MMTIP_000002:0.03,MMTIP_000003:0.04)0.910:0.12,MMTIP_0 con = DBInterface.connect(db) DuckDB.register_data_frame(con, df, "src_df") DBInterface.execute(con, "CREATE TABLE merged AS SELECT * FROM src_df") + # The registration is a view backed by a replacement scan that only + # exists on this connection, so it is removed again: left in the file + # it makes the view unreadable for every later connection. + DuckDB.unregister_data_frame(con, "src_df") DBInterface.close!(con) close(db) before_db_sha = bytes2hex(open(sha256, db_path)) diff --git a/test/unit/test_exact_summaries.jl b/test/unit/test_exact_summaries.jl index 098667c8..d78f23e3 100644 --- a/test/unit/test_exact_summaries.jl +++ b/test/unit/test_exact_summaries.jl @@ -14,10 +14,10 @@ const EXACT_FEATURES = ["f1", "f2", "f3"] @test !s.approximate @test s.samples[1].total == 3 @test s.samples[1].proportions == fill(big(1)//3, 3) - # A zero count in a sample with a positive total has an exact proportion - # of 0 (the independent Python fractions reference below agrees); only a - # zero total leaves proportions undefined. - @test s.samples[2].proportions == Union{Rational{BigInt},Nothing}[big(1)//1, big(0)//1, big(0)//1] + # Sample B has [2, 0, 0] reads: a zero count in a sample that has reads is a + # defined proportion of zero. Only a zero *total* leaves proportions + # undefined, and this sample's total is 2. + @test s.samples[2].proportions == Union{Rational{BigInt},Nothing}[big(1), big(0), big(0)] @test length(s.groups) == 2 @test s.groups[1].label == "g1" @test s.groups[2].counts == BigInt[2, 0, 0] @@ -31,7 +31,7 @@ const EXACT_FEATURES = ["f1", "f2", "f3"] @test storage["input_status"] == "exact" @test storage["samples"][1]["proportions"] == ["1/3", "1/3", "1/3"] @test storage["samples"][1]["total"] == 3 - @test storage["samples"][2]["proportions"][2] == "0/1" + @test storage["samples"][2]["proportions"] == ["1/1", "0/1", "0/1"] @test length(jcs_sha256(storage)) == 64 display = summary_to_display(s; digits=8) @test occursin("exact", display) diff --git a/test/unit/test_validated_statistics.jl b/test/unit/test_validated_statistics.jl index 7a7465d2..bea38a57 100644 --- a/test/unit/test_validated_statistics.jl +++ b/test/unit/test_validated_statistics.jl @@ -106,6 +106,10 @@ end raw = """{"numeric_policy":{"mode":"exact_counts"},"sample_labels":["A"],"feature_labels":["f"],"counts":[[$token]]}""" @test _vs_post("/api/v1/statistics/summaries", raw; raw=true).status == 422 end + + # Fractional input stays invalid as well, even when it is nearly 1. + fractional = """{"numeric_policy":{"mode":"exact_counts"},"sample_labels":["A"],"feature_labels":["f"],"counts":[[1.0000000001]]}""" + @test _vs_post("/api/v1/statistics/summaries", fractional; raw=true).status == 422 negative_string = merge(_VS_EXACT_BODY, Dict("persist" => false, "counts" => [["-1", 1], [1, 0], [1, 0]])) @test _vs_post("/api/v1/statistics/summaries", negative_string).status == 422 From f3b0d2c68050e21c95e9db6c9b84ba51cf62957e Mon Sep 17 00:00:00 2001 From: hyperpolymath <6759885+hyperpolymath@users.noreply.github.com> Date: Fri, 23 Oct 2026 16:46:34 +0000 Subject: [PATCH 2/2] test(cladistics): pin the zero-read half of the vacated-clade tip rule `_reaggregate_after_move` treats a childless clade as a tip when it holds reads *or* when its rank is a sequence terminal. The existing emptied-wrapper test pins the first half (a clade vacated by a move is not a tip). The second half -- a terminal holding zero reads in this table is still a tip -- had no test at all, so nothing would have caught it being dropped. The new case moves a GenusB terminal into a clade holding a read-bearing GenusA tip, a zero-read GenusA tip and an unassigned zero-read tip, and runs both terminal ranks (`"tip"` and `"ASV"`). Two monophyly assertions are the regression guard, and they fail in the right direction if the rule is weakened: with the rule intact GenusA has two tips, so neither reads `true`; drop the rank half and the zero-read tip leaves the global tip count, GenusA looks like a singleton, the read-bearing tip flips to `true` and the zero-read tip gets no verdict (`nothing`) instead of `false`. Taxon purity is asserted too, but as a plain invariant rather than the guard: it is 0.5 either way, because the tie-break between two equally sized genera lands on the same number by a different route. Readings for both the intact and the weakened rule were derived by porting `_reaggregate_after_move` and running the fixture, not by guessing, so the expected values are the algorithm's. Also documents why the `.jplace` test unregisters its `src_df` registration: it is a connection-local replacement scan, and left in the file it makes the view unreadable for every later connection. Co-authored-by: arena-agent <297053741+arena-agent@users.noreply.github.com> --- test/unit/test_clade_cumulus.jl | 63 +++++++++++++++++++++++++++++++++ 1 file changed, 63 insertions(+) diff --git a/test/unit/test_clade_cumulus.jl b/test/unit/test_clade_cumulus.jl index 4d943b97..b2d09497 100644 --- a/test/unit/test_clade_cumulus.jl +++ b/test/unit/test_clade_cumulus.jl @@ -857,4 +857,67 @@ echo "(MMTIP_000001:0.05,(MMTIP_000002:0.03,MMTIP_000003:0.04)0.910:0.12,MMTIP_0 end + @testset "A zero-read terminal still counts as a tip after a move" begin + # The tip rule in `_reaggregate_after_move` has two halves: a clade with + # no children is a tip if it holds reads *or* if it is a sequence + # terminal. The emptied-wrapper test above pins the first half (a + # vacated internal clade is not a tip). This pins the second: a + # terminal that happens to hold zero reads in this table is still a + # tip, so it still counts toward the global per-taxon tip total behind + # the monophyly verdict. + # + # These two monophyly assertions are what move. With the rule intact + # GenusA has two tips, so neither is a singleton and both read `false`. + # Drop the `|| rank in ("tip", "ASV")` half and the zero-read tip leaves + # the count, GenusA looks like a singleton, `positive` flips to `true` + # and `zero` gets no verdict at all (`nothing`) instead of `false`. + # Taxon purity happens not to move here -- it is 0.5 either way, since + # the tie-break between equally sized genera lands on the same number + # by a different route -- so it is asserted as a plain invariant, not + # as the regression guard. + _zero_node(id, rank, count, taxon, children = CC.CladeNode[]) = + CC.CladeNode( + id, id, rank, count, + count + sum(ch.cumulative_count for ch in children; init = 0.0), + 0.0, # cumulative_frequency + 0.0, # residual_count + 5.0, # cloud_size + 0.1, # branch_length + 1.0, # support + 0.0, # subtree_pd + taxon, + taxon === nothing ? nothing : 1.0, + nothing, # monophyletic -- recomputed + :present_in_every, + "#2e7d32", + children, + ) + + for rank in ("tip", "ASV") + positive = _zero_node("positive", rank, 10.0, "GenusA") + zero = _zero_node("zero", rank, 0.0, "GenusA") + unassigned = _zero_node("unassigned", rank, 0.0, nothing) + moving = _zero_node("moving", rank, 20.0, "GenusB") + target = _zero_node("target", "clade", 0.0, "GenusA", [positive, zero, unassigned]) + root = _zero_node("root", "root", 0.0, "GenusA", [target, moving]) + + proposed = CC._reaggregate_after_move( + root, moving, target.id, EP.Case("zero_read_tips", EP.WorldCandidate[]), 30.0, + ) + + # Four tips under `target` once `moving` joins it: two GenusA (one + # with reads, one without), one unassigned, one GenusB. + @test proposed.taxon_purity โ‰ˆ 0.5 atol = 1e-12 + @test CC._find_clade(proposed, target.id).taxon_purity โ‰ˆ 0.5 atol = 1e-12 + # GenusA has two tips, so neither of them is a monophyletic singleton. + @test CC._find_clade(proposed, positive.id).monophyletic == false + @test CC._find_clade(proposed, zero.id).monophyletic == false + # A terminal with no reads stays inactive, so it contributes no PD. + @test CC._find_clade(proposed, zero.id).subtree_pd == 0.0 + # Reads are relocated, not created or lost. + @test proposed.cumulative_count == 30.0 + end + end + + end