Skip to content

TemporalHazard 1.2.8

Choose a tag to compare

@ehrlinger ehrlinger released this 02 Sep 20:11
· 1619 commits to main since this release
f790c73

New features

  • New hzr_theta_names() returns the names of a multiphase theta
    vector, in the order theta requires, before any fit runs. Use it to check
    a hand-written starting vector against the specification it belongs to:

    phases <- list(early = hzr_phase("cdf"), late = hzr_phase("g3"))
    stopifnot(length(theta0) == length(hzr_theta_names(phases)))
    setNames(theta0, hzr_theta_names(phases))

    theta is positional and its entries are not on a common scale -- the late
    phase logs mu and tau but carries gamma, alpha and eta naturally --
    and wrapping the wrong element in log() produces a fit, not an error.
    A comment describing the order is therefore not enough for a template that
    ships to many studies, which is what prompted this: the order is a property
    of the phase specification and changes the moment an author adds, removes or
    retypes a phase.

    The function is a thin wrapper over the naming the optimizer itself uses,
    not a second implementation of it. hazard()'s optimizer, the score test's
    re-expansion and this function now all go through one internal helper, so
    the order documented here cannot drift from the order a fit produces --- a
    test asserts the two are identical for three different phase
    specifications. Phase validation is shared too, so unnamed phases get the
    same phase_1 / phase_2 labels a fit will give them.

  • hzr_stepwise()'s $steps frame gains a stat_type column, saying
    what the stat on that row is and so which reference distribution
    recomputes its p-value: "score_q" (chi-square on df), "wald_z"
    (standard normal) or "wald_chisq" (chi-square on df).

    df could not tell these apart. A scalar Wald is reported as a z, not as
    its square, so it and a score Q are both recorded at df = 1 while calling
    for different distributions. The gap was widest on a candidate rescued by
    the Wald fallback added in 1.2.7: that row carries a Wald z under
    criterion = "score", where every neighbouring row carries a Q. The
    selection was right and p_value was right, but a reader recomputing a
    p-value from stat the way the neighbouring rows permit got an answer wrong
    by dozens of orders of magnitude. Under criterion = "score" the entry
    rows reading "wald_z" are the ones the fallback rescued --- filter on
    action == "enter" & stat_type == "wald_z". Drop rows are always
    Wald-tested under that criterion, following SAS, so they read "wald_z"
    whether or not the fallback fired.

  • hzr_bootstrap() reports Wald fallbacks in select mode, through
    $n_wald_fallback_replicates and $n_wald_fallbacks, with a warning when
    either is non-zero. Every replicate runs under suppressWarnings(), so a
    run in which the fallback fired throughout previously reported nothing at
    all -- and the bootstrap is where a wholesale substitution matters most,
    since those entries drive the pooled selection frequencies.

  • The nomogram parser accepts lines = as well as path =. A multi-fit
    listing needs each nomogram attributed to the fit whose block contains it,
    not to the file. Previously a caller who had already split the listing by
    fit had to write each block back out to a temporary file to parse it, which
    is a workaround rather than an interface. Passing the lines directly now
    works, and is what makes the multiple-nomogram warning actionable.

Bug fixes

  • objective is now recorded on the fit, so refits keep the estimand.
    hazard(objective = "sas") was accepted and acted on by the optimizer, but
    the choice was never stored on the object. Everything that rebuilds a
    hazard() call from fit$spec therefore reverted to
    objective = "likelihood": every hzr_stepwise() candidate refit, the
    Wald-fallback refit and each accepted move, plus the score test's numeric
    Hessian and gradient.

    Selecting on a "sas" base fit consequently differenced delta_logLik,
    aic and delta_aic across two different estimands, and the score
    statistic was computed against a likelihood the fit had not been fitted to.
    On the esophagectomy reference the two objectives differ by about 22
    log-likelihood units -- larger than most single-variable effects -- and the
    run produced a full $steps table with no error and no warning.

    fit$spec$objective now records it, and one accessor feeds every consumer
    so they cannot drift apart. A fit that carries no objective predates the
    argument and is read as "likelihood", which is what it was.

    Note the asymmetry this had created: hzr_bootstrap() refit mode
    re-evaluates the stored call, which did carry objective = "sas", so it was
    already consistent; select mode goes through the scope refit and was not.

    Following from this, objective cannot be overridden per refit. Passing it
    through hzr_stepwise()'s ... is accepted when it restates the base fit's
    own objective and refused when it conflicts, with a message saying why:
    a candidate refit under a different objective would make delta_logLik,
    aic and delta_aic differences between two estimands rather than between
    two models. It is refused rather than quietly ignored, because a full
    $steps table that silently disregarded an explicit argument is the failure
    this entry is about.

  • A candidate that neither criterion could test is now reported as such.
    criterion = "score" declines a candidate whose observed information is
    indefinite at beta = 0 -- which happens when the effect is large -- and
    then refits it and tests it by Wald instead. That rescue can converge,
    producing a perfectly good point estimate, while its Hessian is singular:
    no standard error, so the Wald test cannot be computed either.

    The fallback dropped that case silently. The row kept the score's reason,
    information_indefinite, which describes the first of two independent
    failures and says nothing about the second, and no counter moved. A
    strongly predictive variable vanished and the screen rendered as an honest
    "nothing met slentry" -- the same failure this package has twice fixed
    elsewhere, where a clean-looking screen and an honest null result cannot be
    told apart.

    Such rows now report fallback_no_variance, distinct from
    information_indefinite (the rescue errored or did not converge, and is
    listed in $criteria$refit_failures). hzr_stepwise() and
    hzr_bootstrap() both warn on either, saying the candidate was tested by
    neither criterion and which mechanism applied. The
    information_indefinite prose claimed "that refit also failed", which was
    wrong for the new case and sent readers to an empty refit_failures; both
    it and the bootstrap warning now describe the two mechanisms separately.

    No behaviour changed in what gets selected: a candidate that could not be
    tested is still not entered. What changed is that the run says so.

  • fit$fit$weak now distinguishes "no ridge" from "not checked". The
    weak-identification diagnostic introduced in 1.2.7 returned NULL both when
    a fit had been examined and found well identified and when it could not be
    examined at all -- most importantly when no Hessian was available, which
    happens on an install without the suggested numDeriv package and on fits whose
    rows are left- or interval-censored, where the analytic Hessian declines by
    design. NEWS offered fit$fit$weak as the programmatic check, so on those
    installs it certified as well identified a fit nothing had looked at.

    The field now takes three values: a list when a ridge was found, NULL when
    the fit was examined and is well identified, and NA when the check could
    not run. Test it with is.list(fit$fit$weak) rather than !is.null().
    summary() prints a note in the NA case saying the check did not run, and
    a fit imported with hzr_read_outhaz() -- which has no R Hessian to examine
    -- now reports NA rather than reading as certified clean.

  • A fit with more than one flat direction says so. The detector reported a
    single direction, which invited reading every parameter it did not name as
    identified. It now reports n_directions, the number of distinct
    near-flat directions found, and the warning says when there is more than
    one. Distinctness is counted over the parameters spanning each direction,
    not over eigenvectors: a ridge between two parameters clears the gate twice,
    once on the flat direction and once on its stiff partner, so counting
    eigenvectors would report one ridge as two.

  • A rescued candidate that goes on to win is no longer fitted twice. The
    Wald fallback refits each candidate it rescues, and the acceptance step then
    refit the winner again with identical arguments. The two fits were
    bit-identical, so this was cost rather than incorrectness -- but it doubled
    the price of every accepted fallback entry, against a criterion whose whole
    advantage is that it does not refit per candidate. The rescuing fit is now
    kept and reused, as the Wald path already did with its candidate fits.

  • DELTA is not implemented, and two comments said it was absorbed. The
    headers of R/decomposition.R and R/argument_mapping.R both stated that
    the C DELTA parameter's time transformation
    B(t) = (exp(delta * t) - 1) / delta is "absorbed by decompos()". It is
    not absorbed; it is unimplemented, and delta = 0 is assumed. DELTA
    enters the C reference in three separate places -- it builds rho from
    B(t_half) rather than t_half, it replaces the time argument with
    B(t), and delta * t enters the log-density additively so the density
    carries a factor of exp(delta * t) -- and R computes the delta = 0
    branch of all three.

    The comment was the harmful part. It made the omission look deliberate and
    safe, so a reader looking for exactly this discrepancy was told to stop
    looking, while a PROC HAZARD job with DELTA != 0 was reproduced against a
    different function with no error. Both comments now say what is true.

    The SAS-facing paths now distinguish the two cases rather than treating
    DELTA as one unmapped keyword. hzr_read_outhaz() already stopped on a
    non-zero DELTA; hzr_translate_sas() now records PARMS DELTA = <nonzero>
    as untranslated with a reason saying the emitted call fits a different
    model, and the .lst natural-estimates parser warns when a listing carries
    one. DELTA = 0 and a bare FIXDELTA are treated as faithful translations,
    because that is the branch R implements -- previously both the safe and the
    unsafe case produced the identical generic note "PARMS keyword has no phase
    target", which distinguished nothing.

  • A multiphase fit now records whether Conservation of Events was actually
    applied.
    CoE counts exact events, so it is disabled whenever any status
    falls outside {0, 1} -- which interval or left censoring guarantees -- and
    whenever the model has fewer than two phases. That is deliberate and
    correct. What was missing is that nothing on the returned object said it
    had happened
    : a caller who passed control = list(conserve = TRUE) got a
    fit carrying conserve = TRUE over a computation that did not run.

    It is not an edge case. ICENSOR appears on 42 to 74 blocks per production
    study, and ICENSOR guarantees status leaves {0, 1}, so the auto-disable
    fires constantly. Production also writes NOCONSERVE on 16 to 68 blocks per
    study, so R and SAS usually agree on the outcome -- but for different
    reasons, and a job carrying both CONSERVE and ICENSOR is exactly where
    the two could diverge unobserved.

    A fitted multiphase object now carries two new fields, both alongside the
    requested conserve under fit$spec$control:

    • fit$spec$control$conserve_applied -- logical, whether CoE was actually
      applied;
    • fit$spec$control$conserve_disabled_reason -- one of "not_requested",
      "unsupported_censoring", "single_phase", "no_events" or
      "setup_failed", and NA when CoE was applied.

    Read fit$spec$control$conserve_applied, not fit$spec$control$conserve:
    the latter says only what you asked for. conserve is a
    dist = "multiphase" control; the single-distribution fits do not use it.

    The reason is recorded rather than a bare logical because the causes want
    different responses -- and a bare FALSE reads as "you turned it off" to a
    user who did the opposite.

  • The SAS .lst nomogram parser no longer returns the first of several
    tables silently.
    .hzr_parse_sas_nomogram() matched every nomogram header
    in a listing and read only the first, with no warning and nothing in the
    return value to say a second existed. That is the same shape as the three
    layout defects fixed in 1.2.1 -- silent, and discoverable only by pointing
    the parser at a second study.

    It now warns when a listing holds more than one, and attaches n_found to
    the returned frame whatever the count, so a caller can distinguish "one
    nomogram" from "the first of several" without reading the source. The
    behaviour is otherwise unchanged: the first table is still what comes back.

    Every file in the corpus this was found against happens to print exactly one
    nomogram, so the parser got the right answer there -- by luck of the corpus
    rather than by construction.

Documentation

  • Corrected the documented paths for two fit fields. The weak-identification
    result is at fit$fit$weak and the phase shares at fit$fit$phase_share;
    NEWS.md and, for the shares, ?hazard had both named them one level too
    high. Code following the documented recipe got FALSE from
    is.list(fit$weak) for every fit, ridge or not -- the same wrong answer the
    1.2.8 three-value change was made to prevent, reached by a different route.

  • stat_type no longer over-claims which rows the Wald fallback rescued.
    Under criterion = "score" drop rows are always Wald-tested, following SAS,
    so they read "wald_z" whether or not the fallback fired. The rescued rows
    are the entry rows: filter on
    action == "enter" & stat_type == "wald_z".

  • predict() and hzr_nelson() now say that SAS draws narrower bands.
    SAS %KAPLAN, %NELSONT and PROC HAZPRED all take their band width from
    CLEVEL, whose default is 0.68268948 -- documented in the macro source as
    "(1 sd)". That makes the multiplier 1 to seven decimals -- the literal is
    truncated -- so the band is one standard error, 68.3%, not 95%.

    Nothing here computed the wrong thing: parity is tested and passing, and
    hzr_kaplan() already documented the convention. The gap was that the other
    three entry points did not, and they are the ones a reader meets when
    checking an R fit against an existing SAS figure. At the R default of
    level = 0.95 the reproduced band is about 1.96 times wider than the one
    being checked against, with no error on either side -- so the two look like
    they disagree numerically when they do not.

    The defaults are unchanged. 0.95 is the right R-side default, and adopting
    SAS's silently would make predict() disagree with every other R modelling
    function. The help pages now carry the level to pass instead:
    level = 2 * stats::pnorm(1) - 1.

  • summary()'s documentation and the Inference and diagnostics vignette
    both listed the notes the method prints and had not been updated for the
    ridge note. Both now include it, and the vignette says plainly that a flat
    direction means the point estimates along it are unreliable, not only their
    standard errors.

  • hazard() now states that SAS's STEEPEST has no equivalent.
    PROC HAZARD jobs write STEEPEST QUASI together -- steepest descent, then
    quasi-Newton -- and STEEPEST appears 14 to 109 times per study across the
    corpus. QUASI/QUASINEWTON maps to method = "bfgs"; there is no
    steepest-descent option and no two-stage strategy. Since the multiphase
    likelihood is multimodal, a different descent path can land on a different
    optimum, so a fit translated from such a job may not reproduce SAS's
    estimates. hzr_translate_sas() already recorded the keyword as
    untranslated rather than dropping it; the control$method documentation now
    says why.

Internal

  • Two tests in test-score-wald-fallback.R did not catch the mutations their
    comments named. Widening .hzr_score_fallback_reasons to include constant
    and collinear left both green, because a degenerate candidate still fails
    to enter -- it merely costs a refit on the way out -- and the "noise stays
    out" assertion is guarded by slentry rather than by how narrow the
    fallback is (the fixture's own Wald p-values are 0.0997 and 0.149, so a
    fallback that refit everything would still decline both). Both now assert
    n_wald_fallbacks, which is the quantity that moves.

  • A roxygen block in R/score-test.R bound to the character vector declared
    after it rather than to the function it documents. @noRd, so no Rd was
    affected; source readability only.