Nothing
tests/testthat/test_survmean.R
In popEpi: Functions for Epidemiological Analysis using Population Data
testthat::context("mean survival testing")
library("data.table")
library("Epi")
testthat::test_that("survmean() agrees with old results", {
popEpi:::skip_normally()
sr <- data.table(popEpi::sire)[dg_date < ex_date, ]
sr$agegr <- cut(sr$dg_age, c(0,45,60,Inf), right=FALSE)
x <- Lexis(entry = list(FUT = 0, AGE = dg_age, CAL = get.yrs(dg_date)),
exit = list(CAL = get.yrs(ex_date)),
data = sr,
exit.status = factor(status, levels = 0:2,
labels = c("alive", "canD", "othD")),
entry.status = factor(0, levels = 0:2,
labels = c("alive", "canD", "othD")),
merge = TRUE)
## observed survival
pm <- data.table(popEpi::popmort)
names(pm) <- c("sex", "CAL", "AGE", "haz")
sm <- survmean(Surv(time = FUT, event = lex.Xst != "alive") ~ agegr,
pophaz = pm, data = x,
breaks = list(FUT = seq(0, 10, 1/12)),
e1.breaks = list(FUT = c(seq(0, 10, 1/12), 11:100)))
## values to test against computed on 2016-03-04;
## git ref: 5077677
testthat::expect_equal(sm$est, c(33.951439, 21.611419, 7.604318), tol = 0.005, scale = 1)
testthat::expect_equal(sm$exp, c(45.25686, 31.22712, 13.06725), tol = 0.005, scale = 1)
})
# test_that("survmean() agrees with results computed using pkg survival", {
# popEpi:::skip_normally()
#
#
# BL <- list(fot= seq(0,15,1/24))
# eBL <- list(fot = unique(c(BL$fot, seq(15, 115,0.5))))
#
# sire2 <- data.table(popEpi::sire)[dg_date<ex_date, ]
# sire2$statusf <- factor(sire2$status, levels = 0:2,
# labels = c("alive", "canD", "othD"))
#
# x <- lexpand(sire2,
# birth = bi_date, entry = dg_date, exit = ex_date,
# status = statusf,
# breaks = NULL)
#
# popmort_sm <- setDT(copy(popEpi::popmort))
# setnames(popmort_sm, c("agegroup", "year"), c("age", "per"))
#
# sm <- survmean(Surv(fot, event = lex.Xst) ~ 1,
# breaks = BL, e1.breaks = eBL, r = 1,
# pophaz = popmort_sm, data = x)
#
# st <- survtab(Surv(fot, event = lex.Xst) ~ 1,
# breaks = BL,
# pophaz = popmort_sm, data = x)
#
# sm_curve <- copy(attr(sm, "survmean.meta")$curve)
# td01 <- readRDS("tests/testthat/survmean_test_data_01.rds")
#
# expect_equal(sm_curve[survmean_type == "est", surv],
# td01$est,
# scale = 1L,
# tolerance = 0.001)
# })
testthat::test_that("survmean expected survival curve corresponds to full Ederer I", {
popEpi:::skip_normally()
sr <- data.table(sire)[dg_date < ex_date, ]
sr$agegr <- cut(sr$dg_age, c(0,45,60,Inf), right=FALSE)
x <- Lexis(entry = list(FUT = 0, AGE = dg_age, CAL = get.yrs(dg_date)),
exit = list(CAL = get.yrs(ex_date)),
data = sr,
exit.status = factor(status, levels = 0:2,
labels = c("alive", "canD", "othD")),
entry.status = factor(0, levels = 0:2,
labels = c("alive", "canD", "othD")),
merge = TRUE)
pm <- copy(popEpi::popmort)
names(pm) <- c("sex", "CAL", "AGE", "haz")
BL <- list(FUT = seq(0, 10, 1/12))
eBL <- list(FUT = c(BL$FUT, seq(11,110,1/2)))
sm <- survmean(Surv(time = FUT, event = lex.Xst != "alive") ~ 1,
pophaz = pm, data = x,
breaks = BL,
e1.breaks = eBL)
## pure Ederer I curve
setDT(x)
x[, lex.dur := 110]
setattr(x, "class", c("Lexis", "data.table", "data.frame"))
e1 <- comp_e1(x = x, breaks = eBL,
pophaz = pm, survScale = "FUT")
setkeyv(e1, "FUT")
e1[, "delta" := diff(eBL$FUT)]
e1[, "l1" := c(1,surv.exp[-.N])]
sm.e1 <- e1[, sum((l1+surv.exp)/2L*delta)]
testthat::expect_equal(sm$exp, sm.e1, tol = 0.0005, scale = 1)
})
testthat::test_that("survmean period method is useful", {
popEpi:::skip_normally()
sr <- data.table(popEpi::sire)[dg_date < ex_date, ]
sr$agegr <- cut(sr$dg_age, c(0,45,60,Inf), right=FALSE)
x <- Lexis(entry = list(FUT = 0, AGE = dg_age, CAL = get.yrs(dg_date)),
exit = list(CAL = get.yrs(ex_date)),
data = sr,
exit.status = factor(status, levels = 0:2,
labels = c("alive", "canD", "othD")),
entry.status = factor(0, levels = 0:2,
labels = c("alive", "canD", "othD")),
merge = TRUE)
pm <- data.table(popEpi::popmort)
names(pm) <- c("sex", "CAL", "AGE", "haz")
BL <- list(FUT = seq(0, 10, 1/12))
eBL <- list(FUT = c(BL$FUT, seq(max(BL$FUT),110,1/2)))
eBL$FUT <- sort(unique(eBL$FUT))
sm <- survmean(Surv(time = FUT, event = lex.Xst != "alive") ~ 1,
subset = dg_date >= "1998-01-01" & dg_date < "2003-01-01",
pophaz = pm, data = x,
r = 1,
breaks = BL,
e1.breaks = eBL)
BL <- list(FUT = seq(0, 5, 1/12), CAL = c(1998,2003))
eBL <- list(FUT = c(BL$FUT, seq(max(BL$FUT), 10, 1/12), seq(10,110,1/2)))
eBL$FUT <- sort(unique(eBL$FUT))
smp <- survmean(Surv(time = FUT, event = lex.Xst != "alive") ~ 1,
pophaz = pm, data = x,
breaks = BL,
r = 1,
e1.breaks = eBL)
testthat::expect_equal(sm$obs, smp$obs)
testthat::expect_equal(sm$exp, smp$exp)
testthat::expect_equal(smp$est, 10.01542, tol = 0.0005, scale = 1)
testthat::expect_equal(sm$est, 10.0216, tol = 0.0005, scale = 1)
})
testthat::test_that("Dates and frac. yrs produce congruent results", {
popEpi:::skip_normally()
x <- data.table(popEpi::sire)
x <- x[dg_date<ex_date]
## phony group variable
set.seed(1L)
x$group <- rbinom(nrow(x), 1, 0.5)
## yrs
xy <- Lexis(entry = list(FUT = 0, AGE = dg_age, CAL = get.yrs(dg_date)),
exit = list(CAL = get.yrs(ex_date)),
data = x,
exit.status = factor(status, levels = 0:2,
labels = c("alive", "canD", "othD")),
entry.status = factor(0, levels = 0:2,
labels = c("alive", "canD", "othD")),
merge = TRUE)
## dates
xd <- Lexis(entry = list(FUT = 0L, AGE = dg_date-bi_date, CAL = dg_date),
exit = list(CAL = ex_date),
data = x,
exit.status = factor(status, levels = 0:2,
labels = c("alive", "canD", "othD")),
entry.status = factor(0, levels = 0:2,
labels = c("alive", "canD", "othD")),
merge = TRUE)
yd <- 365.242199
BLy <- list(FUT = seq(0, 9, 1/4))
BLd <- lapply(BLy, function(el) el * yd)
eBLy <- list(FUT = c(seq(0, 9, 1/4), 9.5, 10:75))
eBLd <- lapply(eBLy, function(el) el * yd)
pmy <- data.table(popEpi::popmort)
setnames(pmy, c("year", "agegroup"), c("CAL", "AGE"))
pmd <- data.table(pmy)
pmd[, CAL := as.Date(paste0(CAL, "-01-01"))]
pmd[, AGE := AGE * yd]
pmd[, haz := haz/yd]
#### hazard method
## observed survival & Ederer II
sty <- survmean(Surv(FUT, lex.Xst) ~ group, data = xy,
surv.method = "hazard",
e1.breaks = eBLy,
breaks = BLy, pophaz = pmy)
std <- survmean(Surv(FUT, lex.Xst) ~ group, data = xd,
surv.method = "hazard",
e1.breaks = eBLd,
breaks = BLd, pophaz = pmd)
cuy <- data.table(attributes(sty)$survmean.meta$curves)
cud <- data.table(attributes(std)$survmean.meta$curves)
std[, c("est", "exp", "YPLL") := lapply(.SD, function(col) col/yd),
.SDcols = c("est", "exp", "YPLL")]
testthat::expect_equal(sty$est, std$est, scale = 1L, tolerance = 0.0005)
testthat::expect_equal(sty$exp, std$exp, scale = 1L, tolerance = 0.001)
testthat::expect_equal(sty$obs, std$obs)
testthat::expect_equal(cuy$surv, cud$surv, scale = 1L, tolerance = 0.00005)
#### lifetable method
## observed survival & Ederer II
sty <- survmean(Surv(FUT, lex.Xst) ~ group, data = xy,
surv.method = "lifetable",
e1.breaks = eBLy,
breaks = BLy, pophaz = pmy)
std <- survmean(Surv(FUT, lex.Xst) ~ group, data = xd,
surv.method = "lifetable",
e1.breaks = eBLd,
breaks = BLd, pophaz = pmd)
cuy <- data.table(attributes(sty)$survmean.meta$curves)
cud <- data.table(attributes(std)$survmean.meta$curves)
std[, c("est", "exp", "YPLL") := lapply(.SD, function(col) col/yd),
.SDcols = c("est", "exp", "YPLL")]
testthat::expect_equal(sty$est, std$est, scale = 1L, tolerance = 0.0005)
testthat::expect_equal(sty$exp, std$exp, scale = 1L, tolerance = 0.001)
testthat::expect_equal(sty$obs, std$obs)
testthat::expect_equal(cuy$surv, cud$surv, scale = 1L, tolerance = 0.00005)
})
testthat::test_that("updating works for survmean objects", {
popEpi:::skip_normally()
set.seed(1)
sr <- setDT(popEpi::sire[sample(1:.N, 100), ])
x <- Lexis(entry = list(FUT = 0, AGE = dg_age, CAL = get.yrs(dg_date)),
exit = list(CAL = get.yrs(ex_date)),
data = sr,
exit.status = factor(status, levels = 0:2,
labels = c("alive", "canD", "othD")),
entry.status = factor(0, levels = 0:2,
labels = c("alive", "canD", "othD")),
merge = TRUE)
pm <- data.table(popEpi::popmort)
names(pm) <- c("sex", "CAL", "AGE", "haz")
BL <- list(FUT = 1:10)
fo <- Surv(time = FUT, event = lex.Xst != "alive") ~ 1
e <- quote(
survmean(formula = fo, data = x, breaks = BL,pophaz = pm)
)
sm <- eval(e)
fo <- update(fo, ~.+sex)
sm2 <- eval(e)
sm3 <- update(sm, .~ + sex)
testthat::expect_equal(data.frame(sm2), data.frame(sm3))
testthat::expect_equal(formula(sm2), fo)
testthat::expect_equal(formula(sm3), fo)
testthat::expect_equal(getCall(sm2), e)
})
Try the popEpi package in your browser
Any scripts or data that you put into this service are public.
popEpi documentation built on Aug. 23, 2023, 5:08 p.m.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
testthat::context("mean survival testing")
library("data.table")
library("Epi")
testthat::test_that("survmean() agrees with old results", {
popEpi:::skip_normally()
sr <- data.table(popEpi::sire)[dg_date < ex_date, ]
sr$agegr <- cut(sr$dg_age, c(0,45,60,Inf), right=FALSE)
x <- Lexis(entry = list(FUT = 0, AGE = dg_age, CAL = get.yrs(dg_date)),
exit = list(CAL = get.yrs(ex_date)),
data = sr,
exit.status = factor(status, levels = 0:2,
labels = c("alive", "canD", "othD")),
entry.status = factor(0, levels = 0:2,
labels = c("alive", "canD", "othD")),
merge = TRUE)
## observed survival
pm <- data.table(popEpi::popmort)
names(pm) <- c("sex", "CAL", "AGE", "haz")
sm <- survmean(Surv(time = FUT, event = lex.Xst != "alive") ~ agegr,
pophaz = pm, data = x,
breaks = list(FUT = seq(0, 10, 1/12)),
e1.breaks = list(FUT = c(seq(0, 10, 1/12), 11:100)))
## values to test against computed on 2016-03-04;
## git ref: 5077677
testthat::expect_equal(sm$est, c(33.951439, 21.611419, 7.604318), tol = 0.005, scale = 1)
testthat::expect_equal(sm$exp, c(45.25686, 31.22712, 13.06725), tol = 0.005, scale = 1)
})
# test_that("survmean() agrees with results computed using pkg survival", {
# popEpi:::skip_normally()
#
#
# BL <- list(fot= seq(0,15,1/24))
# eBL <- list(fot = unique(c(BL$fot, seq(15, 115,0.5))))
#
# sire2 <- data.table(popEpi::sire)[dg_date<ex_date, ]
# sire2$statusf <- factor(sire2$status, levels = 0:2,
# labels = c("alive", "canD", "othD"))
#
# x <- lexpand(sire2,
# birth = bi_date, entry = dg_date, exit = ex_date,
# status = statusf,
# breaks = NULL)
#
# popmort_sm <- setDT(copy(popEpi::popmort))
# setnames(popmort_sm, c("agegroup", "year"), c("age", "per"))
#
# sm <- survmean(Surv(fot, event = lex.Xst) ~ 1,
# breaks = BL, e1.breaks = eBL, r = 1,
# pophaz = popmort_sm, data = x)
#
# st <- survtab(Surv(fot, event = lex.Xst) ~ 1,
# breaks = BL,
# pophaz = popmort_sm, data = x)
#
# sm_curve <- copy(attr(sm, "survmean.meta")$curve)
# td01 <- readRDS("tests/testthat/survmean_test_data_01.rds")
#
# expect_equal(sm_curve[survmean_type == "est", surv],
# td01$est,
# scale = 1L,
# tolerance = 0.001)
# })
testthat::test_that("survmean expected survival curve corresponds to full Ederer I", {
popEpi:::skip_normally()
sr <- data.table(sire)[dg_date < ex_date, ]
sr$agegr <- cut(sr$dg_age, c(0,45,60,Inf), right=FALSE)
x <- Lexis(entry = list(FUT = 0, AGE = dg_age, CAL = get.yrs(dg_date)),
exit = list(CAL = get.yrs(ex_date)),
data = sr,
exit.status = factor(status, levels = 0:2,
labels = c("alive", "canD", "othD")),
entry.status = factor(0, levels = 0:2,
labels = c("alive", "canD", "othD")),
merge = TRUE)
pm <- copy(popEpi::popmort)
names(pm) <- c("sex", "CAL", "AGE", "haz")
BL <- list(FUT = seq(0, 10, 1/12))
eBL <- list(FUT = c(BL$FUT, seq(11,110,1/2)))
sm <- survmean(Surv(time = FUT, event = lex.Xst != "alive") ~ 1,
pophaz = pm, data = x,
breaks = BL,
e1.breaks = eBL)
## pure Ederer I curve
setDT(x)
x[, lex.dur := 110]
setattr(x, "class", c("Lexis", "data.table", "data.frame"))
e1 <- comp_e1(x = x, breaks = eBL,
pophaz = pm, survScale = "FUT")
setkeyv(e1, "FUT")
e1[, "delta" := diff(eBL$FUT)]
e1[, "l1" := c(1,surv.exp[-.N])]
sm.e1 <- e1[, sum((l1+surv.exp)/2L*delta)]
testthat::expect_equal(sm$exp, sm.e1, tol = 0.0005, scale = 1)
})
testthat::test_that("survmean period method is useful", {
popEpi:::skip_normally()
sr <- data.table(popEpi::sire)[dg_date < ex_date, ]
sr$agegr <- cut(sr$dg_age, c(0,45,60,Inf), right=FALSE)
x <- Lexis(entry = list(FUT = 0, AGE = dg_age, CAL = get.yrs(dg_date)),
exit = list(CAL = get.yrs(ex_date)),
data = sr,
exit.status = factor(status, levels = 0:2,
labels = c("alive", "canD", "othD")),
entry.status = factor(0, levels = 0:2,
labels = c("alive", "canD", "othD")),
merge = TRUE)
pm <- data.table(popEpi::popmort)
names(pm) <- c("sex", "CAL", "AGE", "haz")
BL <- list(FUT = seq(0, 10, 1/12))
eBL <- list(FUT = c(BL$FUT, seq(max(BL$FUT),110,1/2)))
eBL$FUT <- sort(unique(eBL$FUT))
sm <- survmean(Surv(time = FUT, event = lex.Xst != "alive") ~ 1,
subset = dg_date >= "1998-01-01" & dg_date < "2003-01-01",
pophaz = pm, data = x,
r = 1,
breaks = BL,
e1.breaks = eBL)
BL <- list(FUT = seq(0, 5, 1/12), CAL = c(1998,2003))
eBL <- list(FUT = c(BL$FUT, seq(max(BL$FUT), 10, 1/12), seq(10,110,1/2)))
eBL$FUT <- sort(unique(eBL$FUT))
smp <- survmean(Surv(time = FUT, event = lex.Xst != "alive") ~ 1,
pophaz = pm, data = x,
breaks = BL,
r = 1,
e1.breaks = eBL)
testthat::expect_equal(sm$obs, smp$obs)
testthat::expect_equal(sm$exp, smp$exp)
testthat::expect_equal(smp$est, 10.01542, tol = 0.0005, scale = 1)
testthat::expect_equal(sm$est, 10.0216, tol = 0.0005, scale = 1)
})
testthat::test_that("Dates and frac. yrs produce congruent results", {
popEpi:::skip_normally()
x <- data.table(popEpi::sire)
x <- x[dg_date<ex_date]
## phony group variable
set.seed(1L)
x$group <- rbinom(nrow(x), 1, 0.5)
## yrs
xy <- Lexis(entry = list(FUT = 0, AGE = dg_age, CAL = get.yrs(dg_date)),
exit = list(CAL = get.yrs(ex_date)),
data = x,
exit.status = factor(status, levels = 0:2,
labels = c("alive", "canD", "othD")),
entry.status = factor(0, levels = 0:2,
labels = c("alive", "canD", "othD")),
merge = TRUE)
## dates
xd <- Lexis(entry = list(FUT = 0L, AGE = dg_date-bi_date, CAL = dg_date),
exit = list(CAL = ex_date),
data = x,
exit.status = factor(status, levels = 0:2,
labels = c("alive", "canD", "othD")),
entry.status = factor(0, levels = 0:2,
labels = c("alive", "canD", "othD")),
merge = TRUE)
yd <- 365.242199
BLy <- list(FUT = seq(0, 9, 1/4))
BLd <- lapply(BLy, function(el) el * yd)
eBLy <- list(FUT = c(seq(0, 9, 1/4), 9.5, 10:75))
eBLd <- lapply(eBLy, function(el) el * yd)
pmy <- data.table(popEpi::popmort)
setnames(pmy, c("year", "agegroup"), c("CAL", "AGE"))
pmd <- data.table(pmy)
pmd[, CAL := as.Date(paste0(CAL, "-01-01"))]
pmd[, AGE := AGE * yd]
pmd[, haz := haz/yd]
#### hazard method
## observed survival & Ederer II
sty <- survmean(Surv(FUT, lex.Xst) ~ group, data = xy,
surv.method = "hazard",
e1.breaks = eBLy,
breaks = BLy, pophaz = pmy)
std <- survmean(Surv(FUT, lex.Xst) ~ group, data = xd,
surv.method = "hazard",
e1.breaks = eBLd,
breaks = BLd, pophaz = pmd)
cuy <- data.table(attributes(sty)$survmean.meta$curves)
cud <- data.table(attributes(std)$survmean.meta$curves)
std[, c("est", "exp", "YPLL") := lapply(.SD, function(col) col/yd),
.SDcols = c("est", "exp", "YPLL")]
testthat::expect_equal(sty$est, std$est, scale = 1L, tolerance = 0.0005)
testthat::expect_equal(sty$exp, std$exp, scale = 1L, tolerance = 0.001)
testthat::expect_equal(sty$obs, std$obs)
testthat::expect_equal(cuy$surv, cud$surv, scale = 1L, tolerance = 0.00005)
#### lifetable method
## observed survival & Ederer II
sty <- survmean(Surv(FUT, lex.Xst) ~ group, data = xy,
surv.method = "lifetable",
e1.breaks = eBLy,
breaks = BLy, pophaz = pmy)
std <- survmean(Surv(FUT, lex.Xst) ~ group, data = xd,
surv.method = "lifetable",
e1.breaks = eBLd,
breaks = BLd, pophaz = pmd)
cuy <- data.table(attributes(sty)$survmean.meta$curves)
cud <- data.table(attributes(std)$survmean.meta$curves)
std[, c("est", "exp", "YPLL") := lapply(.SD, function(col) col/yd),
.SDcols = c("est", "exp", "YPLL")]
testthat::expect_equal(sty$est, std$est, scale = 1L, tolerance = 0.0005)
testthat::expect_equal(sty$exp, std$exp, scale = 1L, tolerance = 0.001)
testthat::expect_equal(sty$obs, std$obs)
testthat::expect_equal(cuy$surv, cud$surv, scale = 1L, tolerance = 0.00005)
})
testthat::test_that("updating works for survmean objects", {
popEpi:::skip_normally()
set.seed(1)
sr <- setDT(popEpi::sire[sample(1:.N, 100), ])
x <- Lexis(entry = list(FUT = 0, AGE = dg_age, CAL = get.yrs(dg_date)),
exit = list(CAL = get.yrs(ex_date)),
data = sr,
exit.status = factor(status, levels = 0:2,
labels = c("alive", "canD", "othD")),
entry.status = factor(0, levels = 0:2,
labels = c("alive", "canD", "othD")),
merge = TRUE)
pm <- data.table(popEpi::popmort)
names(pm) <- c("sex", "CAL", "AGE", "haz")
BL <- list(FUT = 1:10)
fo <- Surv(time = FUT, event = lex.Xst != "alive") ~ 1
e <- quote(
survmean(formula = fo, data = x, breaks = BL,pophaz = pm)
)
sm <- eval(e)
fo <- update(fo, ~.+sex)
sm2 <- eval(e)
sm3 <- update(sm, .~ + sex)
testthat::expect_equal(data.frame(sm2), data.frame(sm3))
testthat::expect_equal(formula(sm2), fo)
testthat::expect_equal(formula(sm3), fo)
testthat::expect_equal(getCall(sm2), e)
})
Try the popEpi package in your browser
Any scripts or data that you put into this service are public.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Embedding an R snippet on your website
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.