From 62bb56df776c6abaf43ee85ee58af99e51cce4a8 Mon Sep 17 00:00:00 2001 From: Samuel Jenness Date: Fri, 21 Aug 2026 20:46:46 -0400 Subject: [PATCH] Fix active.sex and young.prop in the sexual-cessation configuration Two defects in build_netstats() when age.limits[2] sits above age.sexual.cessation, so the population keeps aging past the point where it stops forming partnerships. active.sex was zeroed on max(attr_age.grp) as realized rather than on the age group the breaks declare. Whenever the post-cessation group came out empty that marked the oldest sexually active group inactive and zeroed its degree on all three layers, silently. It hit the sampling path and the target_pop data frame path alike, so a projected population whose oldest member sits below the cessation age had the same problem. young.prop = 1 empties the post-cessation group by construction, which triggered the first defect and then failed much later in netest() on a target.stats length mismatch, since an age group with no members carries no ergm term. Both endpoints are now rejected with a message pointing at a usable value, and a band that comes out empty for other reasons warns instead. Neither change touches the default configuration, where the population upper limit and the cessation age coincide and there is no post-cessation band. --- DESCRIPTION | 4 +- NEWS.md | 25 ++++ R/NetStats.R | 55 ++++++++- man/build_netstats.Rd | 6 +- tests/testthat/test-sex-cessation.R | 182 ++++++++++++++++++++++++++++ 5 files changed, 267 insertions(+), 5 deletions(-) create mode 100644 tests/testthat/test-sex-cessation.R diff --git a/DESCRIPTION b/DESCRIPTION index eac6171..0ba7828 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,6 +1,6 @@ Package: ARTnet -Version: 2.9.1 -Date: 2026-06-08 +Version: 2.9.2 +Date: 2026-08-21 Title: Epidemic and Network Model Parameterization with the ARTnet Dataset Description: Epidemic and Network Model Parameterization with the ARTnet Dataset. Maintainer: Samuel Jenness diff --git a/NEWS.md b/NEWS.md index 4e025d7..07bd818 100644 --- a/NEWS.md +++ b/NEWS.md @@ -1,3 +1,28 @@ +## ARTnet v2.9.2 + +Two fixes to the sexual-cessation configuration of `build_netstats()`, where +`age.limits[2]` is set above `age.sexual.cessation` so the population keeps +aging past the point where it stops forming partnerships. + +* `active.sex` is now zeroed on the age group the breaks declare rather than on + the largest age group realized in the sampled population. When the + post-cessation group came out empty, the old behavior marked the oldest + *sexually active* group inactive and zeroed its degree on all three layers, + with no error and no warning. This affected both the sampling path and the + `target_pop` data frame path, where a projected population whose oldest + member sits below the cessation age would have hit the same problem. +* `young.prop` must now be strictly between 0 and 1 when a cessation band + exists, and is rejected with a message pointing at a usable value. Either + endpoint empties one side of the age split by construction, and an age group + with no members carries no `ergm` term, so the failure used to surface much + later as a `target.stats` length mismatch in `netest()`. A post-cessation + group that comes out empty for other reasons (a `network.size` too small for + the seeded share, or a `target_pop` with no members above the cessation age) + now warns for the same reason. + +Neither change affects the default configuration, where the population upper +limit and the cessation age coincide and there is no post-cessation band. + ## ARTnet v2.9.0 * `ARTnetData` is now a suggested package and example `EpiStats` and `Netstats` diff --git a/R/NetStats.R b/R/NetStats.R index db06bea..92ccc4d 100644 --- a/R/NetStats.R +++ b/R/NetStats.R @@ -191,6 +191,10 @@ #' @param young.prop The proportion of the population that should be below the age of sexual cessation. #' Default is NULL (meaning no re-weighting of the `age.pyramid` parameter is performed). #' This parameter is only used if the age of sexual cessation is less than the upper age bound. +#' Must be strictly between 0 and 1: either endpoint empties one side of the split, and an +#' age group with no members cannot carry an `ergm` term. To start a model with effectively +#' nobody above the cessation age, use a value just below 1 (for example 0.999) rather than +#' 1 itself, and let the post-cessation group fill over the burn-in. #' @param target_pop Optional specification of the synthetic target population. Defaults to NULL, #' which uses the legacy patchwork of reference sources (NCHS age pyramid + #' `ARTnetData::race.dist` + ARTnet's own degree / role / risk-quintile distributions). @@ -373,6 +377,30 @@ build_netstats <- function(epistats, netparams, time.unit <- epistats$time.unit + # `young.prop` splits the age pyramid at the sexual cessation age, giving + # `young.prop` of the sampling mass to the sexually active ages and the + # remainder to the post-cessation ages. The two endpoints are degenerate: + # each empties one side of the split entirely, and an age group with no + # members cannot carry an `ergm` term, because levels are derived from the + # values actually present on the network. `netest()` would then fail on a + # target.stats length mismatch, a long way from the cause. Reject them here. + if (sex.cess.mod == TRUE && !is.null(young.prop)) { + if (young.prop >= 1) { + stop("`young.prop` must be less than 1 when the age of sexual cessation ", + "is below the upper age limit. At 1 no node is placed in the ", + "post-cessation age group, so that group's `nodefactor` and ", + "`nodematch` terms would be missing from any model fitted on the ", + "resulting population. Use a value just below 1 (for example ", + "0.999, which seeds about 0.1% of `network.size` above the ", + "cessation age) and let the group fill over the model burn-in.") + } + if (young.prop <= 0) { + stop("`young.prop` must be greater than 0. At 0 the entire population ", + "is placed above the age of sexual cessation, leaving no sexually ", + "active nodes to fit a network model to.") + } + } + # Demographic Initialization ---------------------------------------------- @@ -510,6 +538,14 @@ build_netstats <- function(epistats, netparams, age.breaks <- out$demog$age.breaks <- epistats$age.breaks nquants <- length(netparams$inst$nf.risk.grp) + # The post-cessation age group, as declared by the age breaks rather than as + # realized in the sampled population. `build_epistats` puts + # `age.sexual.cessation` into `age.breaks`, so under `sex.cess.mod` the last + # group is exactly the post-cessation band. Reading the group off the data + # instead (`max(attr_age.grp)`) silently marks the oldest *sexually active* + # group inactive whenever the post-cessation band happens to be empty. + top.age.grp <- length(age.breaks) - 1L + # List-form distribution overrides for the sampling path. NULL means # "use the existing default source" (current behavior). .ov <- if (.tp$form == "list") .tp$overrides else list() @@ -531,7 +567,7 @@ build_netstats <- function(epistats, netparams, } else { as_active <- rep(1L, num) if (sex.cess.mod == TRUE) { - as_active[attr_age.grp == max(attr_age.grp, na.rm = TRUE)] <- 0L + as_active[attr_age.grp %in% top.age.grp] <- 0L } as_active } @@ -598,7 +634,7 @@ build_netstats <- function(epistats, netparams, # sexually active attribute attr_active.sex <- rep(1L, num) if (sex.cess.mod == TRUE) { - attr_active.sex[attr_age.grp == max(attr_age.grp)] <- 0L + attr_active.sex[attr_age.grp %in% top.age.grp] <- 0L } # race attribute @@ -640,6 +676,21 @@ build_netstats <- function(epistats, netparams, attr_role.class <- apportion_lr(num, 0:2, .dist_role.class, shuffled = TRUE) } + # A post-cessation band can still come out empty for reasons other than + # `young.prop`: a `network.size` too small for the seeded share, or a + # `target_pop` data frame with no members above the cessation age. The + # netstats object is self-consistent either way (the band's target statistics + # are zero), but no model with an `age.grp` term can be fitted on the + # population it describes, so say so here rather than leaving it to `netest`. + if (sex.cess.mod == TRUE && !any(attr_age.grp %in% top.age.grp)) { + warning("No node was placed in the post-cessation age group (age.grp ", + top.age.grp, ", ages ", age.breaks[top.age.grp], " and over). ", + "`nodefactor` and `nodematch` terms on `age.grp` will be missing ", + "that level, so fitting against these target statistics will fail ", + "on a length mismatch. Raise `network.size`, lower `young.prop`, ", + "or include post-cessation members in `target_pop`.") + } + # Common attr assignments (both paths) ----------------------------------- out$attr$age <- attr_age out$attr$sqrt.age <- attr_sqrt.age diff --git a/man/build_netstats.Rd b/man/build_netstats.Rd index f6c1e61..8bf008f 100644 --- a/man/build_netstats.Rd +++ b/man/build_netstats.Rd @@ -41,7 +41,11 @@ Hispanic, and White/Other).} \item{young.prop}{The proportion of the population that should be below the age of sexual cessation. Default is NULL (meaning no re-weighting of the \code{age.pyramid} parameter is performed). -This parameter is only used if the age of sexual cessation is less than the upper age bound.} +This parameter is only used if the age of sexual cessation is less than the upper age bound. +Must be strictly between 0 and 1: either endpoint empties one side of the split, and an +age group with no members cannot carry an \code{ergm} term. To start a model with effectively +nobody above the cessation age, use a value just below 1 (for example 0.999) rather than +1 itself, and let the post-cessation group fill over the burn-in.} \item{method}{Character. Either \code{"existing"} (default) or \code{"joint"}. \code{"existing"} reproduces the pre-refactor behavior byte-for-byte: target statistics for edges, nodefactor, and diff --git a/tests/testthat/test-sex-cessation.R b/tests/testthat/test-sex-cessation.R new file mode 100644 index 0000000..3b69fb4 --- /dev/null +++ b/tests/testthat/test-sex-cessation.R @@ -0,0 +1,182 @@ +# Tests for the sexual-cessation configuration of build_netstats(), where +# `age.limits[2]` sits above `age.sexual.cessation` so the population keeps +# aging past the point where it stops forming partnerships. +# +# Two defects motivated these: +# 1. `active.sex` was zeroed on `max(attr_age.grp)` as realized rather than on +# the age group the breaks declare. Whenever the post-cessation group came +# out empty, that marked the oldest *sexually active* group inactive and +# zeroed its degree, silently and without any error. +# 2. `young.prop = 1` empties the post-cessation group by construction, which +# triggered (1) and then failed much later in `netest()` on a target.stats +# length mismatch, because an age group with no members carries no `ergm` +# term. + +skip_without_artnetdata <- function() { + testthat::skip_if(system.file(package = "ARTnetData") == "", + "ARTnetData not installed") +} + +have_data <- system.file(package = "ARTnetData") != "" +if (have_data) { + # Sexual cessation at 65, population to 100: six age groups, the sixth being + # the post-cessation band. + set.seed(20260821L) + .ep_cess <- build_epistats(geog.lvl = "city", geog.cat = "Atlanta", + init.hiv.prev = c(0.33, 0.137, 0.084), + race = TRUE, time.unit = 7, + age.limits = c(15, 100), + age.sexual.cessation = 65) + .np_cess <- build_netparams(.ep_cess, smooth.main.dur = TRUE) + + # The default configuration, where the population limit and the cessation age + # coincide and there is no post-cessation band at all. + set.seed(20260821L) + .ep_plain <- build_epistats(geog.lvl = "city", geog.cat = "Atlanta", + init.hiv.prev = c(0.33, 0.137, 0.084), + race = TRUE, time.unit = 7) + .np_plain <- build_netparams(.ep_plain, smooth.main.dur = TRUE) +} + +mk_pop <- function(n, age_lo, age_hi) { + data.frame( + age = stats::runif(n, age_lo, age_hi), + race = sample(1:3, n, replace = TRUE, prob = c(0.5, 0.05, 0.45)), + deg.casl = sample(0:3, n, replace = TRUE), + deg.main = sample(0:2, n, replace = TRUE), + role.class = sample(0:2, n, replace = TRUE), + risk.grp = sample(1:5, n, replace = TRUE) + ) +} + + +# ---- The configuration itself ------------------------------------------------ + +test_that("sexual cessation adds a post-cessation age group and real mortality", { + skip_without_artnetdata() + expect_true(.ep_cess$sex.cess.mod) + expect_equal(.ep_cess$age.breaks, c(15, 25, 35, 45, 55, 65, 100)) + expect_equal(.ep_cess$age.grps, 6) + + set.seed(101L) + ns <- build_netstats(.ep_cess, .np_cess, expect.mort = 0.000478213, + network.size = 3000) + asmr <- ns$demog$asmr + rates <- as.matrix(asmr[asmr$age >= 65 & asmr$age <= 99, -1]) + # Under the default 15-to-65 configuration mortality is forced to 1 at 65. + # Here the cut has to move to the upper age limit instead. + expect_true(all(rates > 0 & rates < 1)) + expect_true(all(asmr[asmr$age == 100, -1] == 1)) + expect_equal(min(asmr$age[apply(asmr[, -1] >= 1, 1, any)]), 100) +}) + + +# ---- active.sex keys on the declared band, not the realized maximum ---------- + +test_that("active.sex is 0 exactly for nodes at or above the cessation age", { + skip_without_artnetdata() + set.seed(102L) + ns <- build_netstats(.ep_cess, .np_cess, expect.mort = 0.000478213, + network.size = 3000) + expect_gt(sum(ns$attr$age >= 65), 0) + expect_equal(ns$attr$active.sex == 0, ns$attr$age >= 65) + expect_equal(ns$attr$active.sex == 0, ns$attr$age.grp == 6) + # Inactive nodes cannot carry degree on any layer. + inactive <- ns$attr$active.sex == 0 + expect_true(all(ns$attr$deg.main[inactive] == 0)) + expect_true(all(ns$attr$deg.casl[inactive] == 0)) + expect_true(all(ns$attr$deg.tot[inactive] == 0)) +}) + +test_that("an empty post-cessation band does not retire the 55-64 group", { + skip_without_artnetdata() + # This is the regression test for defect (1). Every member of this target + # population is sexually active, so the post-cessation group is empty and the + # realized maximum age group is 5. Reading the band off the data marked all + # 55-to-64 year olds inactive and zeroed their degree. + set.seed(103L) + df <- mk_pop(2000, 20, 64) + expect_warning( + ns <- build_netstats(.ep_cess, .np_cess, expect.mort = 0.000478213, + target_pop = df), + "post-cessation age group" + ) + expect_equal(sum(ns$attr$age.grp == 6), 0) + expect_gt(sum(ns$attr$age.grp == 5), 0) + expect_true(all(ns$attr$active.sex == 1)) + expect_true(any(ns$attr$deg.main > 0)) +}) + +test_that("a populated band is unaffected by the fix", { + skip_without_artnetdata() + set.seed(104L) + df <- mk_pop(2000, 20, 90) + ns <- build_netstats(.ep_cess, .np_cess, expect.mort = 0.000478213, + target_pop = df) + expect_gt(sum(ns$attr$age.grp == 6), 0) + expect_equal(ns$attr$active.sex == 0, ns$attr$age >= 65) +}) + + +# ---- young.prop endpoints ----------------------------------------------------- + +test_that("young.prop = 1 is rejected with an actionable message", { + skip_without_artnetdata() + expect_error( + build_netstats(.ep_cess, .np_cess, expect.mort = 0.000478213, + network.size = 3000, young.prop = 1), + "must be less than 1" + ) + expect_error( + build_netstats(.ep_cess, .np_cess, expect.mort = 0.000478213, + network.size = 3000, young.prop = 1), + "0.999" + ) +}) + +test_that("young.prop = 0 is rejected", { + skip_without_artnetdata() + expect_error( + build_netstats(.ep_cess, .np_cess, expect.mort = 0.000478213, + network.size = 3000, young.prop = 0), + "must be greater than 0" + ) +}) + +test_that("young.prop just below 1 seeds the band and nothing else changes", { + skip_without_artnetdata() + set.seed(105L) + ns <- build_netstats(.ep_cess, .np_cess, expect.mort = 0.000478213, + network.size = 25000, young.prop = 0.999) + n_retired <- sum(ns$attr$age.grp == 6) + expect_gt(n_retired, 0) + expect_lt(n_retired / ns$demog$num, 0.01) + expect_equal(ns$attr$active.sex == 0, ns$attr$age >= 65) + # The band's target statistics are zero on every layer, which is what makes + # `ergm` pin those terms off and keep retired nodes tie-free. + for (layer in c("main", "casl", "inst")) { + nf <- ns[[layer]]$nodefactor_age.grp + expect_length(nf, 6) + expect_equal(unname(nf[6]), 0) + expect_gt(sum(nf[1:5]), 0) + } +}) + + +# ---- The default configuration is untouched ---------------------------------- + +test_that("young.prop is inert without a sexual-cessation band", { + skip_without_artnetdata() + # With `age.limits[2]` equal to the cessation age there is no split to + # reweight, so no value of young.prop should be rejected and none should + # change the result. + set.seed(106L) + a <- build_netstats(.ep_plain, .np_plain, expect.mort = 0.000478213, + network.size = 3000, young.prop = 1) + set.seed(106L) + b <- build_netstats(.ep_plain, .np_plain, expect.mort = 0.000478213, + network.size = 3000) + expect_equal(a$attr, b$attr) + expect_true(all(a$attr$active.sex == 1)) + expect_equal(length(.ep_plain$age.breaks) - 1L, 5L) +})