Nothing
context("Distribution functions and utilities")
## note - standard q fns in R return zero for p=0 for positive dists, but -Inf for real dists.
tol <- 1e-06
test_that("Exponential hazards",{
expect_equal(hexp(c(-Inf, NaN, Inf, NA, -1, 0, 1), c(1,2,3)), c(0,NaN,3,NA,0,3,1))
expect_equal(hexp(c(-Inf, NaN, Inf, NA, -1, 0, 1), c(1,2,3), log=TRUE), log(c(0,NaN,3,NA,0,3,1)))
expect_equal(Hexp(c(-Inf, NaN, Inf, NA, -1, 0, 1, 2, 4), c(1,2,3)), c(0,NaN,Inf, NA,0,0, 1, 4, 12))
expect_equal(Hexp(c(-Inf, NaN, Inf, NA, -1, 0, 1, 2, 4), c(1,2,3), log=TRUE), log(c(0,NaN,Inf, NA,0,0, 1, 4, 12)))
expect_equal(dexp(c(1,2,3), c(2,3,4)) / (1 - pexp(c(1,2,3), c(2,3,4))), hexp(c(1,2,3), c(2,3,4)))
expect_equal(-log(1 - pexp(c(1,2,3), c(2,3,4))), Hexp(c(1,2,3), c(2,3,4)))
expect_equal(log(-log(1 - pexp(c(1,2,3), c(2,3,4)))), Hexp(c(1,2,3), c(2,3,4), log=TRUE))
expect_warning(hexp(c(1,1,1), c(-1, 0, 2)), "Negative rate")
expect_warning(dexp(c(1,1,1), c(-1, 0, 2)), "NaNs produced")
})
test_that("Weibull hazards",{
expect_equal(hweibull(c(1,1,1,1), c(1,1,2,2), c(1,3,1,3)),
dweibull(c(1,1,1,1), c(1,1,2,2), c(1,3,1,3)) / (1 - pweibull(c(1,1,1,1), c(1,1,2,2), c(1,3,1,3))))
expect_equal(Hweibull(c(1,1,1,1), c(1,1,2,2), c(1,3,1,3)),
-pweibull(c(1,1,1,1), c(1,1,2,2), c(1,3,1,3), lower.tail=FALSE, log.p=TRUE))
## exponential reduction
expect_equal(hweibull(c(-Inf, NaN, Inf, NA, -1, 0, 1), c(1,1,1), c(1,2,3)), c(0,NaN,1/3,NA,0,1/3,1))
### positive shape, increasing
hweibull(c(-Inf, -1, 0, 1, 2, Inf), 2, 2.5)
### neg shape, decreasing
hweibull(c(-Inf, -1, 0, 1, 2, Inf), 0.5, 1)
expect_equal(Hexp(c(-Inf, NaN, Inf, NA, -1, 0, 1, 2, 4), c(1,2,3)), c(0,NaN,Inf, NA,0,0, 1, 4, 12))
})
test_that("Log-normal hazards",{
expect_equal(hlnorm(1, 2, 3, log=TRUE), log(hlnorm(1,2,3,log=FALSE)))
expect_equal(Hlnorm(1, 2, 3, log=TRUE), log(-log(1 - plnorm(1,2,3,log=FALSE))))
})
test_that("Gamma hazards",{
expect_equal(hgamma(1, 2, 3, log=TRUE), log(hgamma(1,2,3,log=FALSE)))
expect_equal(Hgamma(1, 2, 3, log=TRUE), log(-log(1 - pgamma(1,2,3,log=FALSE))))
})
test_that("WeibullPH",{
a <- 0.1; m <- 2
b <- m^(-1/a)
x <- c(-Inf, NaN, NA, -1, 0, 1, 2, Inf)
expect_equal(dweibullPH(x, a, m), dweibull(x, a, b))
expect_equal(dweibullPH(x, a, m), dweibull(x, a, b), log=TRUE)
expect_equal(pweibullPH(x, a, m), pweibull(x, a, b))
expect_equal(pweibullPH(x, a, m), pweibull(x, a, b), log.p=TRUE)
expect_equal(pweibullPH(x, a, m), pweibull(x, a, b), lower.tail=FALSE)
qq <- c(0, 0.5, 0.7, 1)
expect_equal(qweibullPH(qq, a, m), qweibull(qq, a, b))
expect_equal(hweibullPH(x, a, m), hweibull(x, a, b))
expect_equal(hweibullPH(x, a, m), hweibull(x, a, b), log=TRUE)
expect_equal(HweibullPH(x, a, m), Hweibull(x, a, b))
expect_equal(HweibullPH(x, a, m), Hweibull(x, a, b), log=TRUE)
set.seed(1); x1 <- rweibull(10, a, b)
set.seed(1); x2 <- rweibullPH(10, a, m)
expect_equal(x1, x2)
fitw <- flexsurvreg(formula = Surv(ovarian$futime, ovarian$fustat) ~ rx, data = ovarian, dist="weibull")
fitwp <- flexsurvreg(formula = Surv(ovarian$futime, ovarian$fustat) ~ rx, data = ovarian, dist="weibullPH")
expect_equal(fitw$res["shape","est"], fitwp$res["shape","est"], tolerance=1e-06)
expect_equal(fitw$res["scale","est"], fitwp$res["scale","est"]^(-1/fitwp$res["shape","est"]), tolerance=1e-05)
expect_equal(coef(fitw)["rx"], -coef(fitwp)["rx"] / fitwp$res["shape","est"])
})
test_that("rmst_generic",{
expect_equal(rmst_lnorm(500, start=250, meanlog=7.4225, sdlog = 1.1138),
rmst_generic(plnorm, t=500, start=250, meanlog=7.4225, sdlog = 1.1138))
})
#implemented in C++
test_that("exph",{
exph_r <- function(y)(y + sqrt(y^2 + 1))
dexph_r <- function(y)(1 + y / sqrt(y^2 + 1))
expect_equal(exph(1.1), exph_r(1.1))
expect_equal(dexph(1.1), dexph_r(1.1))
})
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.