Nothing
tests/testthat/test_loggompertz.R
In drda: Dose-Response Data Analysis
test_that("Constructor", {
x <- lltd$D$x
y <- lltd$D$y
m <- length(unique(x))
n <- length(y)
w <- rep(1, n)
max_iter <- 10000
stats <- lltd$stats_1
start <- c(0, 1, 1, 1)
lower_bound <- c(0, -1, 0.5, 1)
upper_bound <- c(3, 2, 2, 5)
object <- loggompertz_new(x, y, w, NULL, max_iter, NULL, NULL)
expect_true(inherits(object, "loggompertz"))
expect_equal(object$x, x)
expect_equal(object$y, y)
expect_equal(object$w, w)
expect_equal(object$n, n)
expect_equal(object$m, m)
expect_equal(object$stats, stats)
expect_false(object$constrained)
expect_equal(object$max_iter, max_iter)
expect_null(object$start)
expect_null(object$lower_bound)
expect_null(object$upper_bound)
object <- loggompertz_new(x, y, w, start, max_iter, lower_bound, upper_bound)
i <- c(1, 2)
expect_true(inherits(object, "loggompertz"))
expect_equal(object$x, x)
expect_equal(object$y, y)
expect_equal(object$w, w)
expect_equal(object$n, n)
expect_equal(object$m, m)
expect_equal(object$stats, stats)
expect_true(object$constrained)
expect_equal(object$max_iter, max_iter)
expect_equal(object$start, c(start[i], log(start[-i])))
expect_equal(object$lower_bound, c(lower_bound[i], log(lower_bound[-i])))
expect_equal(object$upper_bound, c(upper_bound[i], log(upper_bound[-i])))
w <- lltd$D$w
stats <- lltd$stats_2
object <- loggompertz_new(x, y, w, NULL, max_iter, NULL, NULL)
expect_true(inherits(object, "loggompertz"))
expect_equal(object$x, x)
expect_equal(object$y, y)
expect_equal(object$w, w)
expect_equal(object$n, n)
expect_equal(object$m, m)
expect_equal(object$stats, stats)
expect_false(object$constrained)
expect_equal(object$max_iter, max_iter)
expect_null(object$start)
expect_null(object$lower_bound)
expect_null(object$upper_bound)
object <- loggompertz_new(x, y, w, start, max_iter, lower_bound, upper_bound)
expect_true(inherits(object, "loggompertz"))
expect_equal(object$x, x)
expect_equal(object$y, y)
expect_equal(object$w, w)
expect_equal(object$n, n)
expect_equal(object$m, m)
expect_equal(object$stats, stats)
expect_true(object$constrained)
expect_equal(object$max_iter, max_iter)
expect_equal(object$start, c(start[i], log(start[-i])))
expect_equal(object$lower_bound, c(lower_bound[i], log(lower_bound[-i])))
expect_equal(object$upper_bound, c(upper_bound[i], log(upper_bound[-i])))
})
test_that("Constructor: errors", {
x <- lltd$D$x
y <- lltd$D$y
w <- lltd$D$w
max_iter <- 10000
expect_error(
loggompertz_new(x, y, w, c(0, 1, 1), max_iter, NULL, NULL),
"'start' must be of length 4"
)
expect_error(
loggompertz_new(x, y, w, c(0, 1, 1, 1, 1), max_iter, NULL, NULL),
"'start' must be of length 4"
)
expect_error(
loggompertz_new(x, y, w, c(0, 1, 0, 1), max_iter, NULL, NULL),
"parameter 'eta' cannot be negative nor zero"
)
expect_error(
loggompertz_new(x, y, w, c(0, 1, -1, 1), max_iter, NULL, NULL),
"parameter 'eta' cannot be negative nor zero"
)
expect_error(
loggompertz_new(x, y, w, c(0, 1, 1, 0), max_iter, NULL, NULL),
"parameter 'phi' cannot be negative nor zero"
)
expect_error(
loggompertz_new(x, y, w, c(0, 1, 1, -1), max_iter, NULL, NULL),
"parameter 'phi' cannot be negative nor zero"
)
expect_error(
loggompertz_new(x, y, w, NULL, max_iter, rep(-Inf, 3), rep(Inf, 3)),
"'lower_bound' must be of length 4"
)
expect_error(
loggompertz_new(x, y, w, NULL, max_iter, rep(-Inf, 3), rep(Inf, 4)),
"'lower_bound' must be of length 4"
)
expect_error(
loggompertz_new(x, y, w, NULL, max_iter, rep(-Inf, 4), rep(Inf, 3)),
"'upper_bound' must be of length 4"
)
expect_error(
loggompertz_new(x, y, w, NULL, max_iter, rep(-Inf, 4), c(1, 1, 0, Inf)),
"'upper_bound[3]' cannot be negative nor zero",
fixed = TRUE
)
expect_error(
loggompertz_new(x, y, w, NULL, max_iter, rep(-Inf, 4), c(1, 1, -1, Inf)),
"'upper_bound[3]' cannot be negative nor zero",
fixed = TRUE
)
expect_error(
loggompertz_new(x, y, w, NULL, max_iter, rep(-Inf, 4), c(1, 1, Inf, 0)),
"'upper_bound[4]' cannot be negative nor zero",
fixed = TRUE
)
expect_error(
loggompertz_new(x, y, w, NULL, max_iter, rep(-Inf, 4), c(1, 1, Inf, -1)),
"'upper_bound[4]' cannot be negative nor zero",
fixed = TRUE
)
})
test_that("Function value", {
x <- lltd$stats_1[, 1]
theta <- lltd$theta_g
m <- length(x)
true_value <- c(
1, 0.99835911398420645, 0.82183032092156685, 0.57555097969061526,
0.42486123076253040, 0.33801933438930586, 0.15212234596215889,
0.15002124973437721
)
value <- loggompertz_fn(x, theta)
expect_type(value, "double")
expect_length(value, m)
expect_equal(value, true_value)
object <- structure(list(stats = lltd$stats_1), class = "loggompertz")
value <- fn(object, object$stats[, 1], theta)
expect_type(value, "double")
expect_length(value, m)
expect_equal(value, true_value)
object <- structure(list(stats = lltd$stats_1), class = "loggompertz_fit")
value <- fn(object, object$stats[, 1], theta)
expect_type(value, "double")
expect_length(value, m)
expect_equal(value, true_value)
})
test_that("Gradient (1)", {
x <- lltd$stats_1[, 1]
theta <- lltd$theta_g
m <- length(x)
true_gradient <- matrix(
c(
# alpha
rep(1, m),
# delta
0, 0.0019304541362277092, 0.20961138715109782, 0.49935178859927617,
0.67663384616172894, 0.77880078307140487, 0.99750312239746012,
0.99997500031249740,
# eta
0, 0.0093970540520846307, 0.062120960821991612, -0.053740420946163813,
-0.10559269877403762, -0.11471250798831215, -0.0063500361305660149,
-0.00011258642934322865,
# phi
0, 0.0041022150394838821, 0.11135604942402072, 0.11790250564149576,
0.089865432693354624, 0.066198066561069414, 0.00084787765403784111,
8.4997875026562279e-06
),
nrow = m,
ncol = 4
)
G <- loggompertz_gradient(x, theta)
expect_type(G, "double")
expect_length(G, m * 4)
expect_equal(G, true_gradient)
})
test_that("Hessian (1)", {
x <- lltd$stats_1[, 1]
theta <- lltd$theta_g
m <- length(x)
true_hessian <- array(
c(
# (alpha, alpha)
rep(0, m),
# (alpha, delta)
rep(0, m),
# (alpha, eta)
rep(0, m),
# (alpha, phi)
rep(0, m),
# (delta, alpha)
rep(0, m),
# (delta, delta)
rep(0, m),
# (delta, eta)
0, -0.011055357708334860, -0.073083483319990132, 0.063224024642545663,
0.12422670444004426, 0.13495589175095547, 0.0074706307418423704,
0.00013245462275673958,
# (delta, phi)
0, -0.0048261353405692731, -0.13100711696943614, -0.13870883016646560,
-0.10572403846277015, -0.077880078307140487, -0.00099750312239746012,
-9.9997500031249740e-06,
# (eta, alpha)
rep(0, m),
# (eta, delta)
0, -0.011055357708334860, -0.073083483319990132, 0.063224024642545663,
0.12422670444004426, 0.13495589175095547, 0.0074706307418423704,
0.00013245462275673958,
# (eta, eta)
0, -0.045204776057939509, -0.0077973141424894308, 0.0029938447029538161,
0.030242642409067612, 0.059634488615294076, 0.018975450654134089,
0.00059650372086101613,
# (eta, phi)
0, -0.017682705989635783, 0.041700808527062246, 0.052382979149550082,
0.019194496020505641, -0.0013147191159589372, -0.0021097255890769194,
-0.000040783552121669912,
# (phi, alpha)
rep(0, m),
# (phi, delta)
0, -0.0048261353405692731, -0.13100711696943614, -0.13870883016646560,
-0.10572403846277015, -0.077880078307140487, -0.00099750312239746012,
-9.9997500031249740e-06,
# (phi, eta)
0, -0.017682705989635783, 0.041700808527062246, 0.052382979149550082,
0.019194496020505641, -0.0013147191159589372, -0.0021097255890769194,
-0.000040783552121669912,
# (phi, phi)
0, -0.0094350945908129289, -0.047326321005208805, -0.0091701948832274482,
0.0039316126803342648, 0.0066198066561069414, 0.00016872765315353038,
1.6998725026562190e-06
),
dim = c(m, 4, 4)
)
H <- loggompertz_hessian(x, theta)
expect_type(H, "double")
expect_length(H, m * 4 * 4)
expect_equal(H, true_hessian)
})
test_that("Gradient and Hessian (1)", {
x <- lltd$stats_1[, 1]
theta <- lltd$theta_g
m <- length(x)
true_gradient <- matrix(
c(
# alpha
rep(1, m),
# delta
0, 0.0019304541362277092, 0.20961138715109782, 0.49935178859927617,
0.67663384616172894, 0.77880078307140487, 0.99750312239746012,
0.99997500031249740,
# eta
0, 0.0093970540520846307, 0.062120960821991612, -0.053740420946163813,
-0.10559269877403762, -0.11471250798831215, -0.0063500361305660149,
-0.00011258642934322865,
# phi
0, 0.0041022150394838821, 0.11135604942402072, 0.11790250564149576,
0.089865432693354624, 0.066198066561069414, 0.00084787765403784111,
8.4997875026562279e-06
),
nrow = m,
ncol = 4
)
true_hessian <- array(
c(
# (alpha, alpha)
rep(0, m),
# (alpha, delta)
rep(0, m),
# (alpha, eta)
rep(0, m),
# (alpha, phi)
rep(0, m),
# (delta, alpha)
rep(0, m),
# (delta, delta)
rep(0, m),
# (delta, eta)
0, -0.011055357708334860, -0.073083483319990132, 0.063224024642545663,
0.12422670444004426, 0.13495589175095547, 0.0074706307418423704,
0.00013245462275673958,
# (delta, phi)
0, -0.0048261353405692731, -0.13100711696943614, -0.13870883016646560,
-0.10572403846277015, -0.077880078307140487, -0.00099750312239746012,
-9.9997500031249740e-06,
# (eta, alpha)
rep(0, m),
# (eta, delta)
0, -0.011055357708334860, -0.073083483319990132, 0.063224024642545663,
0.12422670444004426, 0.13495589175095547, 0.0074706307418423704,
0.00013245462275673958,
# (eta, eta)
0, -0.045204776057939509, -0.0077973141424894308, 0.0029938447029538161,
0.030242642409067612, 0.059634488615294076, 0.018975450654134089,
0.00059650372086101613,
# (eta, phi)
0, -0.017682705989635783, 0.041700808527062246, 0.052382979149550082,
0.019194496020505641, -0.0013147191159589372, -0.0021097255890769194,
-0.000040783552121669912,
# (phi, alpha)
rep(0, m),
# (phi, delta)
0, -0.0048261353405692731, -0.13100711696943614, -0.13870883016646560,
-0.10572403846277015, -0.077880078307140487, -0.00099750312239746012,
-9.9997500031249740e-06,
# (phi, eta)
0, -0.017682705989635783, 0.041700808527062246, 0.052382979149550082,
0.019194496020505641, -0.0013147191159589372, -0.0021097255890769194,
-0.000040783552121669912,
# (phi, phi)
0, -0.0094350945908129289, -0.047326321005208805, -0.0091701948832274482,
0.0039316126803342648, 0.0066198066561069414, 0.00016872765315353038,
1.6998725026562190e-06
),
dim = c(m, 4, 4)
)
gh <- loggompertz_gradient_hessian(x, theta)
expect_type(gh, "list")
expect_type(gh$G, "double")
expect_type(gh$H, "double")
expect_length(gh$G, m * 4)
expect_length(gh$H, m * 4 * 4)
expect_equal(gh$G, true_gradient)
expect_equal(gh$H, true_hessian)
})
test_that("Gradient (2)", {
x <- lltd$stats_1[, 1]
theta <- lltd$theta_g
m <- length(x)
true_gradient <- matrix(
c(
# alpha
rep(1, m),
# delta
0, 0.0019304541362277092, 0.20961138715109782, 0.49935178859927617,
0.67663384616172894, 0.77880078307140487, 0.99750312239746012,
0.99997500031249740,
# log_eta
0, 0.018794108104169261, 0.12424192164398322, -0.10748084189232763,
-0.21118539754807525, -0.22942501597662429, -0.012700072261132030,
-0.00022517285868645729,
# log_phi
0, 0.020511075197419411, 0.55678024712010359, 0.58951252820747881,
0.44932716346677312, 0.33099033280534707, 0.0042393882701892055,
0.000042498937513281139
),
nrow = m,
ncol = 4
)
G <- loggompertz_gradient_2(x, theta)
expect_type(G, "double")
expect_length(G, m * 4)
expect_equal(G, true_gradient)
})
test_that("Hessian (2)", {
x <- lltd$stats_1[, 1]
theta <- lltd$theta_g
m <- length(x)
true_hessian <- array(
c(
# (alpha, alpha)
rep(0, m),
# (alpha, delta)
rep(0, m),
# (alpha, log_eta)
rep(0, m),
# (alpha, log_phi)
rep(0, m),
# (delta, alpha)
rep(0, m),
# (delta, delta)
rep(0, m),
# (delta, log_eta)
0, -0.022110715416669719, -0.14616696663998026, 0.12644804928509133,
0.24845340888008853, 0.26991178350191093, 0.014941261483684741,
0.00026490924551347917,
# (delta, log_phi)
0, -0.024130676702846366, -0.65503558484718070, -0.69354415083232801,
-0.52862019231385073, -0.38940039153570243, -0.0049875156119873006,
-0.000049998750015624870,
# (log_eta, alpha)
rep(0, m),
# (log_eta, delta)
0, -0.022110715416669719, -0.14616696663998026, 0.12644804928509133,
0.24845340888008853, 0.26991178350191093, 0.014941261483684741,
0.00026490924551347917,
# (log_eta, log_eta)
0, -0.16202499612758877, 0.093052665074025501, -0.095505463080512362,
-0.090214827911804800, 0.0091129384845520113, 0.063201730355404326,
0.0021608420247576072,
# (log_eta, log_phi)
0, -0.17682705989635783, 0.41700808527062246, 0.52382979149550082,
0.19194496020505641, -0.013147191159589372, -0.021097255890769194,
-0.00040783552121669912,
# (log_phi, alpha)
rep(0, m),
# (log_phi, delta)
0, -0.024130676702846366, -0.65503558484718070, -0.69354415083232801,
-0.52862019231385073, -0.38940039153570243, -0.0049875156119873006,
-0.000049998750015624870,
# (log_phi, log_eta)
0, -0.17682705989635783, 0.41700808527062246, 0.52382979149550082,
0.19194496020505641, -0.013147191159589372, -0.021097255890769194,
-0.00040783552121669912,
# (log_phi, log_phi)
0, -0.21536628957290381, -0.62637777801011654, 0.36025765612679261,
0.54761748047512974, 0.49648549920802060, 0.0084575795990274650,
0.000084995750079686615
),
dim = c(m, 4, 4)
)
H <- loggompertz_hessian_2(x, theta)
expect_type(H, "double")
expect_length(H, m * 4 * 4)
expect_equal(H, true_hessian)
})
test_that("Gradient and Hessian (2)", {
x <- lltd$stats_1[, 1]
theta <- lltd$theta_g
m <- length(x)
true_gradient <- matrix(
c(
# alpha
rep(1, m),
# delta
0, 0.0019304541362277092, 0.20961138715109782, 0.49935178859927617,
0.67663384616172894, 0.77880078307140487, 0.99750312239746012,
0.99997500031249740,
# log_eta
0, 0.018794108104169261, 0.12424192164398322, -0.10748084189232763,
-0.21118539754807525, -0.22942501597662429, -0.012700072261132030,
-0.00022517285868645729,
# log_phi
0, 0.020511075197419411, 0.55678024712010359, 0.58951252820747881,
0.44932716346677312, 0.33099033280534707, 0.0042393882701892055,
0.000042498937513281139
),
nrow = m,
ncol = 4
)
true_hessian <- array(
c(
# (alpha, alpha)
rep(0, m),
# (alpha, delta)
rep(0, m),
# (alpha, log_eta)
rep(0, m),
# (alpha, log_phi)
rep(0, m),
# (delta, alpha)
rep(0, m),
# (delta, delta)
rep(0, m),
# (delta, log_eta)
0, -0.022110715416669719, -0.14616696663998026, 0.12644804928509133,
0.24845340888008853, 0.26991178350191093, 0.014941261483684741,
0.00026490924551347917,
# (delta, log_phi)
0, -0.024130676702846366, -0.65503558484718070, -0.69354415083232801,
-0.52862019231385073, -0.38940039153570243, -0.0049875156119873006,
-0.000049998750015624870,
# (log_eta, alpha)
rep(0, m),
# (log_eta, delta)
0, -0.022110715416669719, -0.14616696663998026, 0.12644804928509133,
0.24845340888008853, 0.26991178350191093, 0.014941261483684741,
0.00026490924551347917,
# (log_eta, log_eta)
0, -0.16202499612758877, 0.093052665074025501, -0.095505463080512362,
-0.090214827911804800, 0.0091129384845520113, 0.063201730355404326,
0.0021608420247576072,
# (log_eta, log_phi)
0, -0.17682705989635783, 0.41700808527062246, 0.52382979149550082,
0.19194496020505641, -0.013147191159589372, -0.021097255890769194,
-0.00040783552121669912,
# (log_phi, alpha)
rep(0, m),
# (log_phi, delta)
0, -0.024130676702846366, -0.65503558484718070, -0.69354415083232801,
-0.52862019231385073, -0.38940039153570243, -0.0049875156119873006,
-0.000049998750015624870,
# (log_phi, log_eta)
0, -0.17682705989635783, 0.41700808527062246, 0.52382979149550082,
0.19194496020505641, -0.013147191159589372, -0.021097255890769194,
-0.00040783552121669912,
# (log_phi, log_phi)
0, -0.21536628957290381, -0.62637777801011654, 0.36025765612679261,
0.54761748047512974, 0.49648549920802060, 0.0084575795990274650,
0.000084995750079686615
),
dim = c(m, 4, 4)
)
gh <- loggompertz_gradient_hessian_2(x, theta)
expect_type(gh, "list")
expect_type(gh$G, "double")
expect_type(gh$H, "double")
expect_length(gh$G, m * 4)
expect_length(gh$H, m * 4 * 4)
expect_equal(gh$G, true_gradient)
expect_equal(gh$H, true_hessian)
object <- structure(list(stats = lltd$stats_1), class = "loggompertz")
gh <- gradient_hessian(object, theta)
expect_type(gh, "list")
expect_type(gh$G, "double")
expect_type(gh$H, "double")
expect_length(gh$G, m * 4)
expect_length(gh$H, m * 4 * 4)
expect_equal(gh$G, true_gradient)
expect_equal(gh$H, true_hessian)
})
test_that("Value of the RSS", {
theta <- lltd$theta_g
theta[3:4] <- log(theta[3:4])
true_value <- 0.32075900098848013
object <- structure(
list(stats = lltd$stats_1, m = nrow(lltd$stats_1)),
class = "loggompertz"
)
rss_fn <- rss(object)
expect_type(rss_fn, "closure")
value <- rss_fn(theta)
expect_type(value, "double")
expect_length(value, 1)
expect_equal(value, true_value)
known_param <- c(theta[1], NA, NA, theta[4])
rss_fn <- rss_fixed(object, known_param)
expect_type(rss_fn, "closure")
value <- rss_fn(theta[2:3])
expect_type(value, "double")
expect_length(value, 1)
expect_equal(value, true_value)
})
test_that("Gradient and Hessian of the RSS", {
theta <- lltd$theta_g
theta[3:4] <- log(theta[3:4])
true_gradient <- c(
2.2653108698938061, 0.73068466091804248, -0.075702268537242359,
0.57850423102806821
)
true_hessian <- matrix(
c(
# alpha
19, 8.2261048577719928, -1.0136637259706368, 5.2206369313459388,
# delta
8.2261048577719928, 6.0552995869783359, -0.70635145977914647,
1.9088026101488051,
# log_eta
-1.0136637259706368, -0.70635145977914647, 0.17854179243952868,
-0.20592244468846133,
# log_phi
5.2206369313459388, 1.9088026101488051, -0.20592244468846133,
2.7235196164555188
),
nrow = 4,
ncol = 4
)
object <- structure(
list(stats = lltd$stats_1, m = nrow(lltd$stats_1)),
class = "loggompertz"
)
rss_gh <- rss_gradient_hessian(object)
expect_type(rss_gh, "closure")
gh <- rss_gh(theta)
expect_type(gh$G, "double")
expect_type(gh$H, "double")
expect_length(gh$G, 4)
expect_length(gh$H, 4 * 4)
expect_equal(gh$G, true_gradient)
expect_equal(gh$H, true_hessian)
known_param <- c(theta[1], NA, NA, theta[4])
rss_gh <- rss_gradient_hessian_fixed(object, known_param)
expect_type(rss_gh, "closure")
gh <- rss_gh(theta[2:3])
expect_type(gh$G, "double")
expect_type(gh$H, "double")
expect_length(gh$G, 2)
expect_length(gh$H, 2 * 2)
expect_equal(gh$G, true_gradient[2:3])
expect_equal(gh$H, true_hessian[2:3, 2:3])
})
test_that("mle_asy", {
x <- lltd$D$x
y <- lltd$D$y
w <- rep(1, length(y))
max_iter <- 10000
theta <- c(0, 1, 0.63524784841453301, 1.6101961281319852)
true_value <- c(
0.83777215113750696, -0.75729507682865963, 0.63524784841453301,
1.6101961281319852
)
object <- loggompertz_new(x, y, w, NULL, max_iter, NULL, NULL)
result <- mle_asy(object, theta)
expect_type(result, "double")
expect_length(result, 4)
expect_equal(result, true_value)
})
test_that("fit", {
x <- lltd$D$x
y <- lltd$D$y
n <- length(y)
w <- rep(1, n)
k <- as.numeric(table(x))
max_iter <- 10000
estimated <- c(alpha = TRUE, delta = TRUE, eta = TRUE, phi = TRUE)
rss_value <- 0.073630385904407914
theta <- c(
alpha = 0.83777215113750696, delta = -0.75729507682865963,
eta = exp(0.63524784841453301), phi = exp(1.6101961281319852)
)
fitted_values <- rep(
c(
0.83777215113750696, 0.83509727314890451, 0.67312406791826807,
0.46539726317063078, 0.33641138140141234, 0.26005233684193214,
0.083128285791155688, 0.080511485993225808
),
k
)
residuals <- c(
0.014927848862493037, -0.080672151137506963, 0.099227848862493037,
-0.0075972731489045085, -0.057297273148904508, 0.049502726851095492,
-0.11702406791826807, 0.031675932081731933, 0.079902736829369223,
0.020602736829369223, 0.069702736829369223, -0.025497263170630777,
0.018888618598587661, -0.023211381401412339, 0.013488618598587661,
-0.11795233684193214, -0.066328285791155688, 0.055471714208844312,
0.042188514006774192
)
object <- loggompertz_new(x, y, w, NULL, max_iter, NULL, NULL)
result <- fit(object)
expect_true(inherits(result, "loggompertz_fit"))
expect_true(inherits(result, "loglogistic"))
expect_true(result$converged)
expect_false(result$constrained)
expect_equal(result$estimated, estimated)
expect_equal(result$coefficients, theta, tolerance = 1.0e-6)
expect_equal(result$rss, rss_value)
expect_equal(result$df.residual, object$n - 4)
expect_equal(result$fitted.values, fitted_values, tolerance = 1.0e-6)
expect_equal(result$residuals, residuals, tolerance = 1.0e-6)
expect_equal(result$weights, w)
object <- loggompertz_new(x, y, w, c(0, 1, 1, 1), max_iter, NULL, NULL)
result <- fit(object)
expect_true(inherits(result, "loggompertz_fit"))
expect_true(inherits(result, "loglogistic"))
expect_true(result$converged)
expect_false(result$constrained)
expect_equal(result$estimated, estimated)
expect_equal(result$coefficients, theta, tolerance = 1.0e-6)
expect_equal(result$rss, rss_value)
expect_equal(result$df.residual, object$n - 4)
expect_equal(result$fitted.values, fitted_values, tolerance = 1.0e-6)
expect_equal(result$residuals, residuals, tolerance = 1.0e-6)
expect_equal(result$weights, w)
})
test_that("fit_constrained: inequalities", {
x <- lltd$D$x
y <- lltd$D$y
n <- length(y)
w <- rep(1, n)
k <- as.numeric(table(x))
max_iter <- 10000
estimated <- c(alpha = TRUE, delta = TRUE, eta = TRUE, phi = TRUE)
rss_value <- 0.073630385904407914
theta <- c(
alpha = 0.83777215113750696, delta = -0.75729507682865963,
eta = exp(0.63524784841453301), phi = exp(1.6101961281319852)
)
fitted_values <- rep(
c(
0.83777215113750696, 0.83509727314890451, 0.67312406791826807,
0.46539726317063078, 0.33641138140141234, 0.26005233684193214,
0.083128285791155688, 0.080511485993225808
),
k
)
residuals <- c(
0.014927848862493037, -0.080672151137506963, 0.099227848862493037,
-0.0075972731489045085, -0.057297273148904508, 0.049502726851095492,
-0.11702406791826807, 0.031675932081731933, 0.079902736829369223,
0.020602736829369223, 0.069702736829369223, -0.025497263170630777,
0.018888618598587661, -0.023211381401412339, 0.013488618598587661,
-0.11795233684193214, -0.066328285791155688, 0.055471714208844312,
0.042188514006774192
)
object <- loggompertz_new(
x, y, w, NULL, max_iter,
c(0.5, -1, 1, 3), c(1, -0.5, 5, 12)
)
result <- fit_constrained(object)
expect_true(inherits(result, "loggompertz_fit"))
expect_true(inherits(result, "loglogistic"))
expect_true(result$converged)
expect_true(result$constrained)
expect_equal(result$estimated, estimated)
expect_equal(result$coefficients, theta, tolerance = 1.0e-6)
expect_equal(result$rss, rss_value)
expect_equal(result$df.residual, object$n - 4)
expect_equal(result$fitted.values, fitted_values, tolerance = 1.0e-6)
expect_equal(result$residuals, residuals, tolerance = 1.0e-6)
expect_equal(result$weights, w)
# initial values within the boundaries
object <- loggompertz_new(
x, y, w, c(0.6, -0.6, 4, 8), max_iter,
c(0.5, -1, 1, 3), c(1, -0.5, 5, 12)
)
result <- fit_constrained(object)
expect_true(inherits(result, "loggompertz_fit"))
expect_true(inherits(result, "loglogistic"))
expect_true(result$converged)
expect_true(result$constrained)
expect_equal(result$estimated, estimated)
expect_equal(result$coefficients, theta, tolerance = 1.0e-6)
expect_equal(result$rss, rss_value)
expect_equal(result$df.residual, object$n - 4)
expect_equal(result$fitted.values, fitted_values, tolerance = 1.0e-6)
expect_equal(result$residuals, residuals, tolerance = 1.0e-6)
expect_equal(result$weights, w)
# initial values outside the boundaries
object <- loggompertz_new(
x, y, w, c(-2, 2, 7, 1), max_iter,
c(0.5, -1, 1, 3), c(1, -0.5, 5, 12)
)
result <- fit_constrained(object)
expect_true(inherits(result, "loggompertz_fit"))
expect_true(inherits(result, "loglogistic"))
expect_true(result$converged)
expect_true(result$constrained)
expect_equal(result$estimated, estimated)
expect_equal(result$coefficients, theta, tolerance = 1.0e-6)
expect_equal(result$rss, rss_value)
expect_equal(result$df.residual, object$n - 4)
expect_equal(result$fitted.values, fitted_values, tolerance = 1.0e-6)
expect_equal(result$residuals, residuals, tolerance = 1.0e-6)
expect_equal(result$weights, w)
})
test_that("fit_constrained: equalities", {
x <- lltd$D$x
y <- lltd$D$y
n <- length(y)
w <- rep(1, n)
k <- as.numeric(table(x))
max_iter <- 10000
estimated <- c(alpha = FALSE, delta = FALSE, eta = TRUE, phi = TRUE)
rss_value <- 0.17654034036023363
theta <- c(
alpha = 0.8, delta = -0.9, eta = exp(0.34256135717772634),
phi = exp(1.8046122081510387)
)
fitted_values <- rep(
c(
0.8, 0.79248376718906222, 0.65161512055560499, 0.47490216690684499,
0.34359567627384339, 0.25186908386433202, -0.082746689754940315,
-0.099320242477951822
),
k
)
residuals <- c(
0.0527, -0.0429, 0.137, 0.035016232810937784, -0.014683767189062216,
0.092116232810937784, -0.095515120555604989, 0.053184879444395011,
0.070397833093155009, 0.011097833093155009, 0.060197833093155009,
-0.035002166906844991, 0.011704323726156605, -0.030395676273843395,
0.0063043237261566051, -0.10976908386433202, 0.099546689754940315,
0.22134668975494032, 0.22202024247795182
)
object <- loggompertz_new(
x, y, w, NULL, max_iter,
c(0.8, -0.9, rep(-Inf, 2)), c(0.8, -0.9, rep(Inf, 2))
)
result <- fit_constrained(object)
expect_true(inherits(result, "loggompertz_fit"))
expect_true(inherits(result, "loglogistic"))
expect_true(result$converged)
expect_false(result$constrained)
expect_equal(result$estimated, estimated)
expect_equal(result$coefficients, theta, tolerance = 1.0e-6)
expect_equal(result$rss, rss_value)
expect_equal(result$df.residual, object$n - 2)
expect_equal(result$fitted.values, fitted_values, tolerance = 1.0e-6)
expect_equal(result$residuals, residuals, tolerance = 1.0e-6)
expect_equal(result$weights, w)
# initial values with same equalities
object <- loggompertz_new(
x, y, w, c(0.8, -0.9, 1, 1), max_iter,
c(0.8, -0.9, rep(-Inf, 2)), c(0.8, -0.9, rep(Inf, 2))
)
result <- fit_constrained(object)
expect_true(inherits(result, "loggompertz_fit"))
expect_true(inherits(result, "loglogistic"))
expect_true(result$converged)
expect_false(result$constrained)
expect_equal(result$estimated, estimated)
expect_equal(result$coefficients, theta, tolerance = 1.0e-6)
expect_equal(result$rss, rss_value)
expect_equal(result$df.residual, object$n - 2)
expect_equal(result$fitted.values, fitted_values, tolerance = 1.0e-6)
expect_equal(result$residuals, residuals, tolerance = 1.0e-6)
expect_equal(result$weights, w)
# initial values with different equalities
object <- loggompertz_new(
x, y, w, c(0, 1, 1, 1), max_iter,
c(0.8, -0.9, rep(-Inf, 2)), c(0.8, -0.9, rep(Inf, 2))
)
result <- fit_constrained(object)
expect_true(inherits(result, "loggompertz_fit"))
expect_true(inherits(result, "loglogistic"))
expect_true(result$converged)
expect_false(result$constrained)
expect_equal(result$estimated, estimated)
expect_equal(result$coefficients, theta, tolerance = 1.0e-6)
expect_equal(result$rss, rss_value)
expect_equal(result$df.residual, object$n - 2)
expect_equal(result$fitted.values, fitted_values, tolerance = 1.0e-6)
expect_equal(result$residuals, residuals, tolerance = 1.0e-6)
expect_equal(result$weights, w)
})
test_that("fit_constrained: equalities and inequalities", {
x <- lltd$D$x
y <- lltd$D$y
n <- length(y)
w <- rep(1, n)
k <- as.numeric(table(x))
max_iter <- 10000
estimated <- c(alpha = FALSE, delta = FALSE, eta = TRUE, phi = TRUE)
rss_value <- 0.17654034036023363
theta <- c(
alpha = 0.8, delta = -0.9, eta = exp(0.34256135717772634),
phi = exp(1.8046122081510387)
)
fitted_values <- rep(
c(
0.8, 0.79248376718906222, 0.65161512055560499, 0.47490216690684499,
0.34359567627384339, 0.25186908386433202, -0.082746689754940315,
-0.099320242477951822
),
k
)
residuals <- c(
0.0527, -0.0429, 0.137, 0.035016232810937784, -0.014683767189062216,
0.092116232810937784, -0.095515120555604989, 0.053184879444395011,
0.070397833093155009, 0.011097833093155009, 0.060197833093155009,
-0.035002166906844991, 0.011704323726156605, -0.030395676273843395,
0.0063043237261566051, -0.10976908386433202, 0.099546689754940315,
0.22134668975494032, 0.22202024247795182
)
object <- loggompertz_new(
x, y, w, NULL, max_iter,
c(0.8, -0.9, 1, 3), c(0.8, -0.9, 5, 12)
)
result <- fit_constrained(object)
expect_true(inherits(result, "loggompertz_fit"))
expect_true(inherits(result, "loglogistic"))
expect_true(result$converged)
expect_true(result$constrained)
expect_equal(result$estimated, estimated)
expect_equal(result$coefficients, theta, tolerance = 1.0e-6)
expect_equal(result$rss, rss_value)
expect_equal(result$df.residual, object$n - 2)
expect_equal(result$fitted.values, fitted_values, tolerance = 1.0e-6)
expect_equal(result$residuals, residuals, tolerance = 1.0e-6)
expect_equal(result$weights, w)
# initial values within the boundaries
object <- loggompertz_new(
x, y, w, c(0.8, -0.9, 3, 7), max_iter,
c(0.8, -0.9, 1, 3), c(0.8, -0.9, 5, 12)
)
result <- fit_constrained(object)
expect_true(inherits(result, "loggompertz_fit"))
expect_true(inherits(result, "loglogistic"))
expect_true(result$converged)
expect_true(result$constrained)
expect_equal(result$estimated, estimated)
expect_equal(result$coefficients, theta, tolerance = 1.0e-6)
expect_equal(result$rss, rss_value)
expect_equal(result$df.residual, object$n - 2)
expect_equal(result$fitted.values, fitted_values, tolerance = 1.0e-6)
expect_equal(result$residuals, residuals, tolerance = 1.0e-6)
expect_equal(result$weights, w)
# initial values outside the boundaries
object <- loggompertz_new(
x, y, w, c(0, 1, 8, 1), max_iter,
c(0.8, -0.9, 1, 3), c(0.8, -0.9, 5, 12)
)
result <- fit_constrained(object)
expect_true(inherits(result, "loggompertz_fit"))
expect_true(inherits(result, "loglogistic"))
expect_true(result$converged)
expect_true(result$constrained)
expect_equal(result$estimated, estimated)
expect_equal(result$coefficients, theta, tolerance = 1.0e-6)
expect_equal(result$rss, rss_value)
expect_equal(result$df.residual, object$n - 2)
expect_equal(result$fitted.values, fitted_values, tolerance = 1.0e-6)
expect_equal(result$residuals, residuals, tolerance = 1.0e-6)
expect_equal(result$weights, w)
})
test_that("fit (weighted)", {
x <- lltd$D$x
y <- lltd$D$y
w <- lltd$D$w
n <- length(y)
k <- as.numeric(table(x))
max_iter <- 10000
estimated <- c(alpha = TRUE, delta = TRUE, eta = TRUE, phi = TRUE)
rss_value <- 0.035138156284360913
theta <- c(
alpha = 0.84575428432937011, delta = -0.76656134360156587,
eta = exp(0.89238310459026989), phi = exp(1.6622042309903930)
)
fitted_values <- rep(
c(
0.84575428432937011, 0.84573608713048233, 0.73789083646643889,
0.47593010582855891, 0.31155933858113756, 0.22406801104864969,
0.079774485527624186, 0.079195048422139125
),
k
)
residuals <- c(
0.0069457156706298950, -0.088654284329370105, 0.091245715670629895,
-0.018236087130482331, -0.067936087130482331, 0.038863912869517669,
-0.18179083646643889, -0.033090836466438890, 0.069369894171441088,
0.010069894171441088, 0.059169894171441088, -0.036030105828558912,
0.043740661418862440, 0.0016406614188624399, 0.038340661418862440,
-0.081968011048649687, -0.062974485527624186, 0.058825514472375814,
0.043504951577860875
)
object <- loggompertz_new(x, y, w, NULL, max_iter, NULL, NULL)
result <- fit(object)
expect_true(inherits(result, "loggompertz_fit"))
expect_true(inherits(result, "loglogistic"))
expect_true(result$converged)
expect_false(result$constrained)
expect_equal(result$estimated, estimated)
expect_equal(result$coefficients, theta, tolerance = 1.0e-6)
expect_equal(result$rss, rss_value)
expect_equal(result$df.residual, object$n - 4)
expect_equal(result$fitted.values, fitted_values, tolerance = 1.0e-6)
expect_equal(result$residuals, residuals, tolerance = 1.0e-6)
expect_equal(result$weights, w)
object <- loggompertz_new(x, y, w, c(0, 1, 1, 1), max_iter, NULL, NULL)
result <- fit(object)
expect_true(inherits(result, "loggompertz_fit"))
expect_true(inherits(result, "loglogistic"))
expect_true(result$converged)
expect_false(result$constrained)
expect_equal(result$estimated, estimated)
expect_equal(result$coefficients, theta, tolerance = 1.0e-6)
expect_equal(result$rss, rss_value)
expect_equal(result$df.residual, object$n - 4)
expect_equal(result$fitted.values, fitted_values, tolerance = 1.0e-6)
expect_equal(result$residuals, residuals, tolerance = 1.0e-6)
expect_equal(result$weights, w)
})
test_that("fit_constrained (weighted): inequalities", {
x <- lltd$D$x
y <- lltd$D$y
w <- lltd$D$w
n <- length(y)
k <- as.numeric(table(x))
max_iter <- 10000
estimated <- c(alpha = TRUE, delta = TRUE, eta = TRUE, phi = TRUE)
rss_value <- 0.035138156284360913
theta <- c(
alpha = 0.84575428432937011, delta = -0.76656134360156587,
eta = exp(0.89238310459026989), phi = exp(1.6622042309903930)
)
fitted_values <- rep(
c(
0.84575428432937011, 0.84573608713048233, 0.73789083646643889,
0.47593010582855891, 0.31155933858113756, 0.22406801104864969,
0.079774485527624186, 0.079195048422139125
),
k
)
residuals <- c(
0.0069457156706298950, -0.088654284329370105, 0.091245715670629895,
-0.018236087130482331, -0.067936087130482331, 0.038863912869517669,
-0.18179083646643889, -0.033090836466438890, 0.069369894171441088,
0.010069894171441088, 0.059169894171441088, -0.036030105828558912,
0.043740661418862440, 0.0016406614188624399, 0.038340661418862440,
-0.081968011048649687, -0.062974485527624186, 0.058825514472375814,
0.043504951577860875
)
object <- loggompertz_new(
x, y, w, NULL, max_iter,
c(0.5, -1, 1, 3), c(1, -0.5, 5, 12)
)
result <- fit_constrained(object)
expect_true(inherits(result, "loggompertz_fit"))
expect_true(inherits(result, "loglogistic"))
expect_true(result$converged)
expect_true(result$constrained)
expect_equal(result$estimated, estimated)
expect_equal(result$coefficients, theta, tolerance = 1.0e-6)
expect_equal(result$rss, rss_value)
expect_equal(result$df.residual, object$n - 4)
expect_equal(result$fitted.values, fitted_values, tolerance = 1.0e-6)
expect_equal(result$residuals, residuals, tolerance = 1.0e-6)
expect_equal(result$weights, w)
# initial values within the boundaries
object <- loggompertz_new(
x, y, w, c(0.6, -0.6, 4, 8), max_iter,
c(0.5, -1, 1, 3), c(1, -0.5, 5, 12)
)
result <- fit_constrained(object)
expect_true(inherits(result, "loggompertz_fit"))
expect_true(inherits(result, "loglogistic"))
expect_true(result$converged)
expect_true(result$constrained)
expect_equal(result$estimated, estimated)
expect_equal(result$coefficients, theta, tolerance = 1.0e-6)
expect_equal(result$rss, rss_value)
expect_equal(result$df.residual, object$n - 4)
expect_equal(result$fitted.values, fitted_values, tolerance = 1.0e-6)
expect_equal(result$residuals, residuals, tolerance = 1.0e-6)
expect_equal(result$weights, w)
# initial values outside the boundaries
object <- loggompertz_new(
x, y, w, c(-2, 2, 7, 1), max_iter,
c(0.5, -1, 1, 3), c(1, -0.5, 5, 12)
)
result <- fit_constrained(object)
expect_true(inherits(result, "loggompertz_fit"))
expect_true(inherits(result, "loglogistic"))
expect_true(result$converged)
expect_true(result$constrained)
expect_equal(result$estimated, estimated)
expect_equal(result$coefficients, theta, tolerance = 1.0e-6)
expect_equal(result$rss, rss_value)
expect_equal(result$df.residual, object$n - 4)
expect_equal(result$fitted.values, fitted_values, tolerance = 1.0e-6)
expect_equal(result$residuals, residuals, tolerance = 1.0e-6)
expect_equal(result$weights, w)
})
test_that("fit_constrained (weighted): equalities", {
x <- lltd$D$x
y <- lltd$D$y
w <- lltd$D$w
n <- length(y)
k <- as.numeric(table(x))
max_iter <- 10000
estimated <- c(alpha = FALSE, delta = FALSE, eta = TRUE, phi = TRUE)
rss_value <- 0.12152790634234742
theta <- c(
alpha = 0.8, delta = -0.9, eta = exp(0.67344700327430744),
phi = exp(1.8267797417585254)
)
fitted_values <- rep(
c(
0.8, 0.79991222146934487, 0.71604576658544182, 0.49162778848529489,
0.31063536742142226, 0.19268868846616760, -0.096135387777275743,
-0.099957631352575174
),
k
)
residuals <- c(
0.0527, -0.0429, 0.137, 0.027587778530655130, -0.022112221469344870,
0.084687778530655130, -0.15994576658544182, -0.011245766585441821,
0.053672211514705107, -0.0056277884852948934, 0.043472211514705107,
-0.051727788485294893, 0.044664632578577739, 0.0025646325785777394,
0.039264632578577739, -0.050588688466167597, 0.11293538777727574,
0.23473538777727574, 0.22265763135257517
)
object <- loggompertz_new(
x, y, w, NULL, max_iter,
c(0.8, -0.9, rep(-Inf, 2)), c(0.8, -0.9, rep(Inf, 2))
)
result <- fit_constrained(object)
expect_true(inherits(result, "loggompertz_fit"))
expect_true(inherits(result, "loglogistic"))
expect_true(result$converged)
expect_false(result$constrained)
expect_equal(result$estimated, estimated)
expect_equal(result$coefficients, theta, tolerance = 1.0e-6)
expect_equal(result$rss, rss_value)
expect_equal(result$df.residual, object$n - 2)
expect_equal(result$fitted.values, fitted_values, tolerance = 1.0e-6)
expect_equal(result$residuals, residuals, tolerance = 1.0e-6)
expect_equal(result$weights, w)
# initial values with same equalities
object <- loggompertz_new(
x, y, w, c(0.8, -0.9, 1, 1), max_iter,
c(0.8, -0.9, rep(-Inf, 2)), c(0.8, -0.9, rep(Inf, 2))
)
result <- fit_constrained(object)
expect_true(inherits(result, "loggompertz_fit"))
expect_true(inherits(result, "loglogistic"))
expect_true(result$converged)
expect_false(result$constrained)
expect_equal(result$estimated, estimated)
expect_equal(result$coefficients, theta, tolerance = 1.0e-6)
expect_equal(result$rss, rss_value)
expect_equal(result$df.residual, object$n - 2)
expect_equal(result$fitted.values, fitted_values, tolerance = 1.0e-6)
expect_equal(result$residuals, residuals, tolerance = 1.0e-6)
expect_equal(result$weights, w)
# initial values with different equalities
object <- loggompertz_new(
x, y, w, c(0, 1, 1, 1), max_iter,
c(0.8, -0.9, rep(-Inf, 2)), c(0.8, -0.9, rep(Inf, 2))
)
result <- fit_constrained(object)
expect_true(inherits(result, "loggompertz_fit"))
expect_true(inherits(result, "loglogistic"))
expect_true(result$converged)
expect_false(result$constrained)
expect_equal(result$estimated, estimated)
expect_equal(result$coefficients, theta, tolerance = 1.0e-6)
expect_equal(result$rss, rss_value)
expect_equal(result$df.residual, object$n - 2)
expect_equal(result$fitted.values, fitted_values, tolerance = 1.0e-6)
expect_equal(result$residuals, residuals, tolerance = 1.0e-6)
expect_equal(result$weights, w)
})
test_that("fit_constrained (weighted): equalities and inequalities", {
x <- lltd$D$x
y <- lltd$D$y
w <- lltd$D$w
n <- length(y)
k <- as.numeric(table(x))
max_iter <- 10000
estimated <- c(alpha = FALSE, delta = FALSE, eta = TRUE, phi = TRUE)
rss_value <- 0.12152790634234742
theta <- c(
alpha = 0.8, delta = -0.9, eta = exp(0.67344700327430744),
phi = exp(1.8267797417585254)
)
fitted_values <- rep(
c(
0.8, 0.79991222146934487, 0.71604576658544182, 0.49162778848529489,
0.31063536742142226, 0.19268868846616760, -0.096135387777275743,
-0.099957631352575174
),
k
)
residuals <- c(
0.0527, -0.0429, 0.137, 0.027587778530655130, -0.022112221469344870,
0.084687778530655130, -0.15994576658544182, -0.011245766585441821,
0.053672211514705107, -0.0056277884852948934, 0.043472211514705107,
-0.051727788485294893, 0.044664632578577739, 0.0025646325785777394,
0.039264632578577739, -0.050588688466167597, 0.11293538777727574,
0.23473538777727574, 0.22265763135257517
)
object <- loggompertz_new(
x, y, w, NULL, max_iter,
c(0.8, -0.9, 1, 3), c(0.8, -0.9, 5, 12)
)
result <- fit_constrained(object)
expect_true(inherits(result, "loggompertz_fit"))
expect_true(inherits(result, "loglogistic"))
expect_true(result$converged)
expect_true(result$constrained)
expect_equal(result$estimated, estimated)
expect_equal(result$coefficients, theta, tolerance = 1.0e-6)
expect_equal(result$rss, rss_value)
expect_equal(result$df.residual, object$n - 2)
expect_equal(result$fitted.values, fitted_values, tolerance = 1.0e-6)
expect_equal(result$residuals, residuals, tolerance = 1.0e-6)
expect_equal(result$weights, w)
# initial values within the boundaries
object <- loggompertz_new(
x, y, w, c(0.8, -0.9, 3, 7), max_iter,
c(0.8, -0.9, 1, 3), c(0.8, -0.9, 5, 12)
)
result <- fit_constrained(object)
expect_true(inherits(result, "loggompertz_fit"))
expect_true(inherits(result, "loglogistic"))
expect_true(result$converged)
expect_true(result$constrained)
expect_equal(result$estimated, estimated)
expect_equal(result$coefficients, theta, tolerance = 1.0e-6)
expect_equal(result$rss, rss_value)
expect_equal(result$df.residual, object$n - 2)
expect_equal(result$fitted.values, fitted_values, tolerance = 1.0e-6)
expect_equal(result$residuals, residuals, tolerance = 1.0e-6)
expect_equal(result$weights, w)
# initial values outside the boundaries
object <- loggompertz_new(
x, y, w, c(0, 1, 8, 1), max_iter,
c(0.8, -0.9, 1, 3), c(0.8, -0.9, 5, 12)
)
result <- fit_constrained(object)
expect_true(inherits(result, "loggompertz_fit"))
expect_true(inherits(result, "loglogistic"))
expect_true(result$converged)
expect_true(result$constrained)
expect_equal(result$estimated, estimated)
expect_equal(result$coefficients, theta, tolerance = 1.0e-6)
expect_equal(result$rss, rss_value)
expect_equal(result$df.residual, object$n - 2)
expect_equal(result$fitted.values, fitted_values, tolerance = 1.0e-6)
expect_equal(result$residuals, residuals, tolerance = 1.0e-6)
expect_equal(result$weights, w)
})
test_that("fisher_info", {
x <- lltd$D$x
y <- lltd$D$y
w <- lltd$D$w
max_iter <- 10000
theta <- lltd$theta_g
names(theta) <- c("alpha", "delta", "eta", "phi")
sigma <- lltd$sigma
true_value <- matrix(c(
# alpha
6206.96, 2748.9449629599373, -189.18362192800753, 344.04436069543435,
-28355.127259362132,
# delta
2748.9449629599373, 2023.9485432497464, -124.69711327516109,
134.06347723771206, -10502.299351663864,
# eta
-189.18362192800753, -124.69711327516109, 23.773250144905426,
-11.783571918074014, 850.51618712697142,
# phi
344.04436069543435, 134.06347723771206, -11.783571918074014,
29.686894742510176, -1441.3142462178932,
# sigma
-28355.127259362132, -10502.299351663864, 850.51618712697142,
-1441.3142462178932, 113909.14953385276
),
nrow = 5,
ncol = 5
)
rownames(true_value) <- colnames(true_value) <- c(
"alpha", "delta", "eta", "phi", "sigma"
)
object <- loggompertz_new(x, y, w, NULL, max_iter, NULL, NULL)
fim <- fisher_info(object, theta, sigma)
expect_type(fim, "double")
expect_length(fim, 5 * 5)
expect_equal(fim, true_value)
})
test_that("drda: 'lower_bound' argument errors", {
x <- lltd$D$x
y <- lltd$D$y
expect_error(
drda(
y ~ x, mean_function = "loggompertz",
lower_bound = c("a", "b", "c", "d")
),
"'lower_bound' must be a numeric vector"
)
expect_error(
drda(
y ~ x, mean_function = "loggompertz",
lower_bound = matrix(-Inf, nrow = 4, ncol = 2),
upper_bound = rep(Inf, 4)
),
"'lower_bound' must be a numeric vector"
)
expect_error(
drda(
y ~ x, mean_function = "loggompertz",
lower_bound = rep(-Inf, 5),
upper_bound = rep(Inf, 4)
),
"'lower_bound' and 'upper_bound' must have the same length"
)
expect_error(
drda(
y ~ x, mean_function = "loggompertz",
lower_bound = c( 0, -Inf, -Inf, -Inf),
upper_bound = c(-1, Inf, Inf, Inf)
),
"'lower_bound' cannot be larger than 'upper_bound'"
)
expect_error(
drda(
y ~ x, mean_function = "loggompertz",
lower_bound = c(Inf, -Inf, -Inf, -Inf),
upper_bound = c(Inf, Inf, Inf, Inf)
),
"'lower_bound' cannot be equal to infinity"
)
expect_error(
drda(
y ~ x, mean_function = "loggompertz",
lower_bound = rep(-Inf, 5),
upper_bound = rep(Inf, 5)
),
"'lower_bound' must be of length 4"
)
})
test_that("drda: 'upper_bound' argument errors", {
x <- lltd$D$x
y <- lltd$D$y
expect_error(
drda(
y ~ x, mean_function = "loggompertz",
upper_bound = c("a", "b", "c", "d")
),
"'upper_bound' must be a numeric vector"
)
expect_error(
drda(
y ~ x, mean_function = "loggompertz",
lower_bound = rep(-Inf, 4),
upper_bound = matrix(Inf, nrow = 4, ncol = 2)
),
"'upper_bound' must be a numeric vector"
)
expect_error(
drda(
y ~ x, mean_function = "loggompertz",
lower_bound = c(-Inf, -Inf, -Inf, -Inf),
upper_bound = c(-Inf, Inf, Inf, Inf)
),
"'upper_bound' cannot be equal to -infinity"
)
expect_error(
drda(
y ~ x, mean_function = "loggompertz",
lower_bound = rep(-Inf, 5),
upper_bound = rep(Inf, 5)
),
"'lower_bound' must be of length 4"
)
})
test_that("drda: 'start' argument errors", {
x <- lltd$D$x
y <- lltd$D$y
expect_error(
drda(
y ~ x, mean_function = "loggompertz",
start = c("a", "b", "c", "d")
),
"'start' must be a numeric vector"
)
expect_error(
drda(
y ~ x, mean_function = "loggompertz",
start = c(0, Inf, 1, 1)
),
"'start' must be finite"
)
expect_error(
drda(
y ~ x, mean_function = "loggompertz",
start = c(-Inf, 1, 1, 1)
),
"'start' must be finite"
)
expect_error(
drda(
y ~ x, mean_function = "loggompertz",
start = rep(1, 5)
),
"'start' must be of length 4"
)
expect_error(
drda(
y ~ x, mean_function = "loggompertz",
start = c(0, 1, -1, 1)
),
"parameter 'eta' cannot be negative nor zero"
)
expect_error(
drda(
y ~ x, mean_function = "loggompertz",
start = c(0, 1, 0, 1)
),
"parameter 'eta' cannot be negative nor zero"
)
expect_error(
drda(
y ~ x, mean_function = "loggompertz",
start = c(0, 1, 1, -1)
),
"parameter 'phi' cannot be negative nor zero"
)
expect_error(
drda(
y ~ x, mean_function = "loggompertz",
start = c(0, 1, 1, 0)
),
"parameter 'phi' cannot be negative nor zero"
)
})
test_that("nauc: decreasing", {
x <- lltd$D$x
y <- lltd$D$y
w <- lltd$D$w
result <- drda(y ~ x, weights = w, mean_function = "loggompertz")
expect_equal(nauc(result), 0.085299967613769906, tolerance = 1.0e-6)
expect_equal(nauc(result, xlim = c(0, 2)), 0.84575366076886644)
expect_equal(nauc(result, ylim = c(0.3, 0.7)), 0.0059283502711447546)
expect_equal(nauc(result, xlim = c(0, 2), ylim = c(0.3, 0.7)), 1.0)
expect_equal(
nauc(result, xlim = c(5, 8), ylim = c(0.3, 0.7)), 0.33844687600885550
)
expect_equal(nauc(result, xlim = c(10, 15), ylim = c(0.3, 0.7)), 0.0)
})
test_that("naac: decreasing", {
x <- lltd$D$x
y <- lltd$D$y
w <- lltd$D$w
result <- drda(y ~ x, weights = w, mean_function = "loggompertz")
expect_equal(naac(result), 1 - 0.085299967613769906)
expect_equal(naac(result, xlim = c(0, 2)), 1 - 0.84575366076886644)
expect_equal(naac(result, ylim = c(0.3, 0.7)), 1 - 0.0059283502711447546)
expect_equal(naac(result, xlim = c(0, 2), ylim = c(0.3, 0.7)), 0.0)
expect_equal(
naac(result, xlim = c(5, 8), ylim = c(0.3, 0.7)), 1 - 0.33844687600885550
)
expect_equal(naac(result, xlim = c(10, 15), ylim = c(0.3, 0.7)), 1.0)
})
test_that("nauc: increasing", {
x <- lltd$D$x
y <- rev(lltd$D$y)
w <- lltd$D$w
result <- drda(y ~ x, weights = w, mean_function = "loggompertz")
expect_equal(nauc(result), 0.88478280795251264)
expect_equal(nauc(result, xlim = c(0, 2)), 0.129049490857619206)
expect_equal(nauc(result, ylim = c(0.3, 0.7)), 0.99551714173529925)
expect_equal(nauc(result, xlim = c(0, 2), ylim = c(0.3, 0.7)), 0.0)
expect_equal(
nauc(result, xlim = c(5, 8), ylim = c(0.3, 0.7)), 0.85655631688941921
)
expect_equal(nauc(result, xlim = c(9, 12), ylim = c(0.3, 0.7)), 1.0)
})
test_that("naac: increasing", {
x <- lltd$D$x
y <- rev(lltd$D$y)
w <- lltd$D$w
result <- drda(y ~ x, weights = w, mean_function = "loggompertz")
expect_equal(naac(result), 1 - 0.88478280795251264)
expect_equal(naac(result, xlim = c(0, 2)), 1 - 0.129049490857619206)
expect_equal(naac(result, ylim = c(0.3, 0.7)), 1 - 0.99551714173529925)
expect_equal(naac(result, xlim = c(0, 2), ylim = c(0.3, 0.7)), 1.0)
expect_equal(
naac(result, xlim = c(5, 8), ylim = c(0.3, 0.7)), 1 - 0.85655631688941921
)
expect_equal(naac(result, xlim = c(9, 12), ylim = c(0.3, 0.7)), 0.0)
})
Try the drda package in your browser
Any scripts or data that you put into this service are public.
drda documentation built on April 3, 2025, 6 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.
test_that("Constructor", {
x <- lltd$D$x
y <- lltd$D$y
m <- length(unique(x))
n <- length(y)
w <- rep(1, n)
max_iter <- 10000
stats <- lltd$stats_1
start <- c(0, 1, 1, 1)
lower_bound <- c(0, -1, 0.5, 1)
upper_bound <- c(3, 2, 2, 5)
object <- loggompertz_new(x, y, w, NULL, max_iter, NULL, NULL)
expect_true(inherits(object, "loggompertz"))
expect_equal(object$x, x)
expect_equal(object$y, y)
expect_equal(object$w, w)
expect_equal(object$n, n)
expect_equal(object$m, m)
expect_equal(object$stats, stats)
expect_false(object$constrained)
expect_equal(object$max_iter, max_iter)
expect_null(object$start)
expect_null(object$lower_bound)
expect_null(object$upper_bound)
object <- loggompertz_new(x, y, w, start, max_iter, lower_bound, upper_bound)
i <- c(1, 2)
expect_true(inherits(object, "loggompertz"))
expect_equal(object$x, x)
expect_equal(object$y, y)
expect_equal(object$w, w)
expect_equal(object$n, n)
expect_equal(object$m, m)
expect_equal(object$stats, stats)
expect_true(object$constrained)
expect_equal(object$max_iter, max_iter)
expect_equal(object$start, c(start[i], log(start[-i])))
expect_equal(object$lower_bound, c(lower_bound[i], log(lower_bound[-i])))
expect_equal(object$upper_bound, c(upper_bound[i], log(upper_bound[-i])))
w <- lltd$D$w
stats <- lltd$stats_2
object <- loggompertz_new(x, y, w, NULL, max_iter, NULL, NULL)
expect_true(inherits(object, "loggompertz"))
expect_equal(object$x, x)
expect_equal(object$y, y)
expect_equal(object$w, w)
expect_equal(object$n, n)
expect_equal(object$m, m)
expect_equal(object$stats, stats)
expect_false(object$constrained)
expect_equal(object$max_iter, max_iter)
expect_null(object$start)
expect_null(object$lower_bound)
expect_null(object$upper_bound)
object <- loggompertz_new(x, y, w, start, max_iter, lower_bound, upper_bound)
expect_true(inherits(object, "loggompertz"))
expect_equal(object$x, x)
expect_equal(object$y, y)
expect_equal(object$w, w)
expect_equal(object$n, n)
expect_equal(object$m, m)
expect_equal(object$stats, stats)
expect_true(object$constrained)
expect_equal(object$max_iter, max_iter)
expect_equal(object$start, c(start[i], log(start[-i])))
expect_equal(object$lower_bound, c(lower_bound[i], log(lower_bound[-i])))
expect_equal(object$upper_bound, c(upper_bound[i], log(upper_bound[-i])))
})
test_that("Constructor: errors", {
x <- lltd$D$x
y <- lltd$D$y
w <- lltd$D$w
max_iter <- 10000
expect_error(
loggompertz_new(x, y, w, c(0, 1, 1), max_iter, NULL, NULL),
"'start' must be of length 4"
)
expect_error(
loggompertz_new(x, y, w, c(0, 1, 1, 1, 1), max_iter, NULL, NULL),
"'start' must be of length 4"
)
expect_error(
loggompertz_new(x, y, w, c(0, 1, 0, 1), max_iter, NULL, NULL),
"parameter 'eta' cannot be negative nor zero"
)
expect_error(
loggompertz_new(x, y, w, c(0, 1, -1, 1), max_iter, NULL, NULL),
"parameter 'eta' cannot be negative nor zero"
)
expect_error(
loggompertz_new(x, y, w, c(0, 1, 1, 0), max_iter, NULL, NULL),
"parameter 'phi' cannot be negative nor zero"
)
expect_error(
loggompertz_new(x, y, w, c(0, 1, 1, -1), max_iter, NULL, NULL),
"parameter 'phi' cannot be negative nor zero"
)
expect_error(
loggompertz_new(x, y, w, NULL, max_iter, rep(-Inf, 3), rep(Inf, 3)),
"'lower_bound' must be of length 4"
)
expect_error(
loggompertz_new(x, y, w, NULL, max_iter, rep(-Inf, 3), rep(Inf, 4)),
"'lower_bound' must be of length 4"
)
expect_error(
loggompertz_new(x, y, w, NULL, max_iter, rep(-Inf, 4), rep(Inf, 3)),
"'upper_bound' must be of length 4"
)
expect_error(
loggompertz_new(x, y, w, NULL, max_iter, rep(-Inf, 4), c(1, 1, 0, Inf)),
"'upper_bound[3]' cannot be negative nor zero",
fixed = TRUE
)
expect_error(
loggompertz_new(x, y, w, NULL, max_iter, rep(-Inf, 4), c(1, 1, -1, Inf)),
"'upper_bound[3]' cannot be negative nor zero",
fixed = TRUE
)
expect_error(
loggompertz_new(x, y, w, NULL, max_iter, rep(-Inf, 4), c(1, 1, Inf, 0)),
"'upper_bound[4]' cannot be negative nor zero",
fixed = TRUE
)
expect_error(
loggompertz_new(x, y, w, NULL, max_iter, rep(-Inf, 4), c(1, 1, Inf, -1)),
"'upper_bound[4]' cannot be negative nor zero",
fixed = TRUE
)
})
test_that("Function value", {
x <- lltd$stats_1[, 1]
theta <- lltd$theta_g
m <- length(x)
true_value <- c(
1, 0.99835911398420645, 0.82183032092156685, 0.57555097969061526,
0.42486123076253040, 0.33801933438930586, 0.15212234596215889,
0.15002124973437721
)
value <- loggompertz_fn(x, theta)
expect_type(value, "double")
expect_length(value, m)
expect_equal(value, true_value)
object <- structure(list(stats = lltd$stats_1), class = "loggompertz")
value <- fn(object, object$stats[, 1], theta)
expect_type(value, "double")
expect_length(value, m)
expect_equal(value, true_value)
object <- structure(list(stats = lltd$stats_1), class = "loggompertz_fit")
value <- fn(object, object$stats[, 1], theta)
expect_type(value, "double")
expect_length(value, m)
expect_equal(value, true_value)
})
test_that("Gradient (1)", {
x <- lltd$stats_1[, 1]
theta <- lltd$theta_g
m <- length(x)
true_gradient <- matrix(
c(
# alpha
rep(1, m),
# delta
0, 0.0019304541362277092, 0.20961138715109782, 0.49935178859927617,
0.67663384616172894, 0.77880078307140487, 0.99750312239746012,
0.99997500031249740,
# eta
0, 0.0093970540520846307, 0.062120960821991612, -0.053740420946163813,
-0.10559269877403762, -0.11471250798831215, -0.0063500361305660149,
-0.00011258642934322865,
# phi
0, 0.0041022150394838821, 0.11135604942402072, 0.11790250564149576,
0.089865432693354624, 0.066198066561069414, 0.00084787765403784111,
8.4997875026562279e-06
),
nrow = m,
ncol = 4
)
G <- loggompertz_gradient(x, theta)
expect_type(G, "double")
expect_length(G, m * 4)
expect_equal(G, true_gradient)
})
test_that("Hessian (1)", {
x <- lltd$stats_1[, 1]
theta <- lltd$theta_g
m <- length(x)
true_hessian <- array(
c(
# (alpha, alpha)
rep(0, m),
# (alpha, delta)
rep(0, m),
# (alpha, eta)
rep(0, m),
# (alpha, phi)
rep(0, m),
# (delta, alpha)
rep(0, m),
# (delta, delta)
rep(0, m),
# (delta, eta)
0, -0.011055357708334860, -0.073083483319990132, 0.063224024642545663,
0.12422670444004426, 0.13495589175095547, 0.0074706307418423704,
0.00013245462275673958,
# (delta, phi)
0, -0.0048261353405692731, -0.13100711696943614, -0.13870883016646560,
-0.10572403846277015, -0.077880078307140487, -0.00099750312239746012,
-9.9997500031249740e-06,
# (eta, alpha)
rep(0, m),
# (eta, delta)
0, -0.011055357708334860, -0.073083483319990132, 0.063224024642545663,
0.12422670444004426, 0.13495589175095547, 0.0074706307418423704,
0.00013245462275673958,
# (eta, eta)
0, -0.045204776057939509, -0.0077973141424894308, 0.0029938447029538161,
0.030242642409067612, 0.059634488615294076, 0.018975450654134089,
0.00059650372086101613,
# (eta, phi)
0, -0.017682705989635783, 0.041700808527062246, 0.052382979149550082,
0.019194496020505641, -0.0013147191159589372, -0.0021097255890769194,
-0.000040783552121669912,
# (phi, alpha)
rep(0, m),
# (phi, delta)
0, -0.0048261353405692731, -0.13100711696943614, -0.13870883016646560,
-0.10572403846277015, -0.077880078307140487, -0.00099750312239746012,
-9.9997500031249740e-06,
# (phi, eta)
0, -0.017682705989635783, 0.041700808527062246, 0.052382979149550082,
0.019194496020505641, -0.0013147191159589372, -0.0021097255890769194,
-0.000040783552121669912,
# (phi, phi)
0, -0.0094350945908129289, -0.047326321005208805, -0.0091701948832274482,
0.0039316126803342648, 0.0066198066561069414, 0.00016872765315353038,
1.6998725026562190e-06
),
dim = c(m, 4, 4)
)
H <- loggompertz_hessian(x, theta)
expect_type(H, "double")
expect_length(H, m * 4 * 4)
expect_equal(H, true_hessian)
})
test_that("Gradient and Hessian (1)", {
x <- lltd$stats_1[, 1]
theta <- lltd$theta_g
m <- length(x)
true_gradient <- matrix(
c(
# alpha
rep(1, m),
# delta
0, 0.0019304541362277092, 0.20961138715109782, 0.49935178859927617,
0.67663384616172894, 0.77880078307140487, 0.99750312239746012,
0.99997500031249740,
# eta
0, 0.0093970540520846307, 0.062120960821991612, -0.053740420946163813,
-0.10559269877403762, -0.11471250798831215, -0.0063500361305660149,
-0.00011258642934322865,
# phi
0, 0.0041022150394838821, 0.11135604942402072, 0.11790250564149576,
0.089865432693354624, 0.066198066561069414, 0.00084787765403784111,
8.4997875026562279e-06
),
nrow = m,
ncol = 4
)
true_hessian <- array(
c(
# (alpha, alpha)
rep(0, m),
# (alpha, delta)
rep(0, m),
# (alpha, eta)
rep(0, m),
# (alpha, phi)
rep(0, m),
# (delta, alpha)
rep(0, m),
# (delta, delta)
rep(0, m),
# (delta, eta)
0, -0.011055357708334860, -0.073083483319990132, 0.063224024642545663,
0.12422670444004426, 0.13495589175095547, 0.0074706307418423704,
0.00013245462275673958,
# (delta, phi)
0, -0.0048261353405692731, -0.13100711696943614, -0.13870883016646560,
-0.10572403846277015, -0.077880078307140487, -0.00099750312239746012,
-9.9997500031249740e-06,
# (eta, alpha)
rep(0, m),
# (eta, delta)
0, -0.011055357708334860, -0.073083483319990132, 0.063224024642545663,
0.12422670444004426, 0.13495589175095547, 0.0074706307418423704,
0.00013245462275673958,
# (eta, eta)
0, -0.045204776057939509, -0.0077973141424894308, 0.0029938447029538161,
0.030242642409067612, 0.059634488615294076, 0.018975450654134089,
0.00059650372086101613,
# (eta, phi)
0, -0.017682705989635783, 0.041700808527062246, 0.052382979149550082,
0.019194496020505641, -0.0013147191159589372, -0.0021097255890769194,
-0.000040783552121669912,
# (phi, alpha)
rep(0, m),
# (phi, delta)
0, -0.0048261353405692731, -0.13100711696943614, -0.13870883016646560,
-0.10572403846277015, -0.077880078307140487, -0.00099750312239746012,
-9.9997500031249740e-06,
# (phi, eta)
0, -0.017682705989635783, 0.041700808527062246, 0.052382979149550082,
0.019194496020505641, -0.0013147191159589372, -0.0021097255890769194,
-0.000040783552121669912,
# (phi, phi)
0, -0.0094350945908129289, -0.047326321005208805, -0.0091701948832274482,
0.0039316126803342648, 0.0066198066561069414, 0.00016872765315353038,
1.6998725026562190e-06
),
dim = c(m, 4, 4)
)
gh <- loggompertz_gradient_hessian(x, theta)
expect_type(gh, "list")
expect_type(gh$G, "double")
expect_type(gh$H, "double")
expect_length(gh$G, m * 4)
expect_length(gh$H, m * 4 * 4)
expect_equal(gh$G, true_gradient)
expect_equal(gh$H, true_hessian)
})
test_that("Gradient (2)", {
x <- lltd$stats_1[, 1]
theta <- lltd$theta_g
m <- length(x)
true_gradient <- matrix(
c(
# alpha
rep(1, m),
# delta
0, 0.0019304541362277092, 0.20961138715109782, 0.49935178859927617,
0.67663384616172894, 0.77880078307140487, 0.99750312239746012,
0.99997500031249740,
# log_eta
0, 0.018794108104169261, 0.12424192164398322, -0.10748084189232763,
-0.21118539754807525, -0.22942501597662429, -0.012700072261132030,
-0.00022517285868645729,
# log_phi
0, 0.020511075197419411, 0.55678024712010359, 0.58951252820747881,
0.44932716346677312, 0.33099033280534707, 0.0042393882701892055,
0.000042498937513281139
),
nrow = m,
ncol = 4
)
G <- loggompertz_gradient_2(x, theta)
expect_type(G, "double")
expect_length(G, m * 4)
expect_equal(G, true_gradient)
})
test_that("Hessian (2)", {
x <- lltd$stats_1[, 1]
theta <- lltd$theta_g
m <- length(x)
true_hessian <- array(
c(
# (alpha, alpha)
rep(0, m),
# (alpha, delta)
rep(0, m),
# (alpha, log_eta)
rep(0, m),
# (alpha, log_phi)
rep(0, m),
# (delta, alpha)
rep(0, m),
# (delta, delta)
rep(0, m),
# (delta, log_eta)
0, -0.022110715416669719, -0.14616696663998026, 0.12644804928509133,
0.24845340888008853, 0.26991178350191093, 0.014941261483684741,
0.00026490924551347917,
# (delta, log_phi)
0, -0.024130676702846366, -0.65503558484718070, -0.69354415083232801,
-0.52862019231385073, -0.38940039153570243, -0.0049875156119873006,
-0.000049998750015624870,
# (log_eta, alpha)
rep(0, m),
# (log_eta, delta)
0, -0.022110715416669719, -0.14616696663998026, 0.12644804928509133,
0.24845340888008853, 0.26991178350191093, 0.014941261483684741,
0.00026490924551347917,
# (log_eta, log_eta)
0, -0.16202499612758877, 0.093052665074025501, -0.095505463080512362,
-0.090214827911804800, 0.0091129384845520113, 0.063201730355404326,
0.0021608420247576072,
# (log_eta, log_phi)
0, -0.17682705989635783, 0.41700808527062246, 0.52382979149550082,
0.19194496020505641, -0.013147191159589372, -0.021097255890769194,
-0.00040783552121669912,
# (log_phi, alpha)
rep(0, m),
# (log_phi, delta)
0, -0.024130676702846366, -0.65503558484718070, -0.69354415083232801,
-0.52862019231385073, -0.38940039153570243, -0.0049875156119873006,
-0.000049998750015624870,
# (log_phi, log_eta)
0, -0.17682705989635783, 0.41700808527062246, 0.52382979149550082,
0.19194496020505641, -0.013147191159589372, -0.021097255890769194,
-0.00040783552121669912,
# (log_phi, log_phi)
0, -0.21536628957290381, -0.62637777801011654, 0.36025765612679261,
0.54761748047512974, 0.49648549920802060, 0.0084575795990274650,
0.000084995750079686615
),
dim = c(m, 4, 4)
)
H <- loggompertz_hessian_2(x, theta)
expect_type(H, "double")
expect_length(H, m * 4 * 4)
expect_equal(H, true_hessian)
})
test_that("Gradient and Hessian (2)", {
x <- lltd$stats_1[, 1]
theta <- lltd$theta_g
m <- length(x)
true_gradient <- matrix(
c(
# alpha
rep(1, m),
# delta
0, 0.0019304541362277092, 0.20961138715109782, 0.49935178859927617,
0.67663384616172894, 0.77880078307140487, 0.99750312239746012,
0.99997500031249740,
# log_eta
0, 0.018794108104169261, 0.12424192164398322, -0.10748084189232763,
-0.21118539754807525, -0.22942501597662429, -0.012700072261132030,
-0.00022517285868645729,
# log_phi
0, 0.020511075197419411, 0.55678024712010359, 0.58951252820747881,
0.44932716346677312, 0.33099033280534707, 0.0042393882701892055,
0.000042498937513281139
),
nrow = m,
ncol = 4
)
true_hessian <- array(
c(
# (alpha, alpha)
rep(0, m),
# (alpha, delta)
rep(0, m),
# (alpha, log_eta)
rep(0, m),
# (alpha, log_phi)
rep(0, m),
# (delta, alpha)
rep(0, m),
# (delta, delta)
rep(0, m),
# (delta, log_eta)
0, -0.022110715416669719, -0.14616696663998026, 0.12644804928509133,
0.24845340888008853, 0.26991178350191093, 0.014941261483684741,
0.00026490924551347917,
# (delta, log_phi)
0, -0.024130676702846366, -0.65503558484718070, -0.69354415083232801,
-0.52862019231385073, -0.38940039153570243, -0.0049875156119873006,
-0.000049998750015624870,
# (log_eta, alpha)
rep(0, m),
# (log_eta, delta)
0, -0.022110715416669719, -0.14616696663998026, 0.12644804928509133,
0.24845340888008853, 0.26991178350191093, 0.014941261483684741,
0.00026490924551347917,
# (log_eta, log_eta)
0, -0.16202499612758877, 0.093052665074025501, -0.095505463080512362,
-0.090214827911804800, 0.0091129384845520113, 0.063201730355404326,
0.0021608420247576072,
# (log_eta, log_phi)
0, -0.17682705989635783, 0.41700808527062246, 0.52382979149550082,
0.19194496020505641, -0.013147191159589372, -0.021097255890769194,
-0.00040783552121669912,
# (log_phi, alpha)
rep(0, m),
# (log_phi, delta)
0, -0.024130676702846366, -0.65503558484718070, -0.69354415083232801,
-0.52862019231385073, -0.38940039153570243, -0.0049875156119873006,
-0.000049998750015624870,
# (log_phi, log_eta)
0, -0.17682705989635783, 0.41700808527062246, 0.52382979149550082,
0.19194496020505641, -0.013147191159589372, -0.021097255890769194,
-0.00040783552121669912,
# (log_phi, log_phi)
0, -0.21536628957290381, -0.62637777801011654, 0.36025765612679261,
0.54761748047512974, 0.49648549920802060, 0.0084575795990274650,
0.000084995750079686615
),
dim = c(m, 4, 4)
)
gh <- loggompertz_gradient_hessian_2(x, theta)
expect_type(gh, "list")
expect_type(gh$G, "double")
expect_type(gh$H, "double")
expect_length(gh$G, m * 4)
expect_length(gh$H, m * 4 * 4)
expect_equal(gh$G, true_gradient)
expect_equal(gh$H, true_hessian)
object <- structure(list(stats = lltd$stats_1), class = "loggompertz")
gh <- gradient_hessian(object, theta)
expect_type(gh, "list")
expect_type(gh$G, "double")
expect_type(gh$H, "double")
expect_length(gh$G, m * 4)
expect_length(gh$H, m * 4 * 4)
expect_equal(gh$G, true_gradient)
expect_equal(gh$H, true_hessian)
})
test_that("Value of the RSS", {
theta <- lltd$theta_g
theta[3:4] <- log(theta[3:4])
true_value <- 0.32075900098848013
object <- structure(
list(stats = lltd$stats_1, m = nrow(lltd$stats_1)),
class = "loggompertz"
)
rss_fn <- rss(object)
expect_type(rss_fn, "closure")
value <- rss_fn(theta)
expect_type(value, "double")
expect_length(value, 1)
expect_equal(value, true_value)
known_param <- c(theta[1], NA, NA, theta[4])
rss_fn <- rss_fixed(object, known_param)
expect_type(rss_fn, "closure")
value <- rss_fn(theta[2:3])
expect_type(value, "double")
expect_length(value, 1)
expect_equal(value, true_value)
})
test_that("Gradient and Hessian of the RSS", {
theta <- lltd$theta_g
theta[3:4] <- log(theta[3:4])
true_gradient <- c(
2.2653108698938061, 0.73068466091804248, -0.075702268537242359,
0.57850423102806821
)
true_hessian <- matrix(
c(
# alpha
19, 8.2261048577719928, -1.0136637259706368, 5.2206369313459388,
# delta
8.2261048577719928, 6.0552995869783359, -0.70635145977914647,
1.9088026101488051,
# log_eta
-1.0136637259706368, -0.70635145977914647, 0.17854179243952868,
-0.20592244468846133,
# log_phi
5.2206369313459388, 1.9088026101488051, -0.20592244468846133,
2.7235196164555188
),
nrow = 4,
ncol = 4
)
object <- structure(
list(stats = lltd$stats_1, m = nrow(lltd$stats_1)),
class = "loggompertz"
)
rss_gh <- rss_gradient_hessian(object)
expect_type(rss_gh, "closure")
gh <- rss_gh(theta)
expect_type(gh$G, "double")
expect_type(gh$H, "double")
expect_length(gh$G, 4)
expect_length(gh$H, 4 * 4)
expect_equal(gh$G, true_gradient)
expect_equal(gh$H, true_hessian)
known_param <- c(theta[1], NA, NA, theta[4])
rss_gh <- rss_gradient_hessian_fixed(object, known_param)
expect_type(rss_gh, "closure")
gh <- rss_gh(theta[2:3])
expect_type(gh$G, "double")
expect_type(gh$H, "double")
expect_length(gh$G, 2)
expect_length(gh$H, 2 * 2)
expect_equal(gh$G, true_gradient[2:3])
expect_equal(gh$H, true_hessian[2:3, 2:3])
})
test_that("mle_asy", {
x <- lltd$D$x
y <- lltd$D$y
w <- rep(1, length(y))
max_iter <- 10000
theta <- c(0, 1, 0.63524784841453301, 1.6101961281319852)
true_value <- c(
0.83777215113750696, -0.75729507682865963, 0.63524784841453301,
1.6101961281319852
)
object <- loggompertz_new(x, y, w, NULL, max_iter, NULL, NULL)
result <- mle_asy(object, theta)
expect_type(result, "double")
expect_length(result, 4)
expect_equal(result, true_value)
})
test_that("fit", {
x <- lltd$D$x
y <- lltd$D$y
n <- length(y)
w <- rep(1, n)
k <- as.numeric(table(x))
max_iter <- 10000
estimated <- c(alpha = TRUE, delta = TRUE, eta = TRUE, phi = TRUE)
rss_value <- 0.073630385904407914
theta <- c(
alpha = 0.83777215113750696, delta = -0.75729507682865963,
eta = exp(0.63524784841453301), phi = exp(1.6101961281319852)
)
fitted_values <- rep(
c(
0.83777215113750696, 0.83509727314890451, 0.67312406791826807,
0.46539726317063078, 0.33641138140141234, 0.26005233684193214,
0.083128285791155688, 0.080511485993225808
),
k
)
residuals <- c(
0.014927848862493037, -0.080672151137506963, 0.099227848862493037,
-0.0075972731489045085, -0.057297273148904508, 0.049502726851095492,
-0.11702406791826807, 0.031675932081731933, 0.079902736829369223,
0.020602736829369223, 0.069702736829369223, -0.025497263170630777,
0.018888618598587661, -0.023211381401412339, 0.013488618598587661,
-0.11795233684193214, -0.066328285791155688, 0.055471714208844312,
0.042188514006774192
)
object <- loggompertz_new(x, y, w, NULL, max_iter, NULL, NULL)
result <- fit(object)
expect_true(inherits(result, "loggompertz_fit"))
expect_true(inherits(result, "loglogistic"))
expect_true(result$converged)
expect_false(result$constrained)
expect_equal(result$estimated, estimated)
expect_equal(result$coefficients, theta, tolerance = 1.0e-6)
expect_equal(result$rss, rss_value)
expect_equal(result$df.residual, object$n - 4)
expect_equal(result$fitted.values, fitted_values, tolerance = 1.0e-6)
expect_equal(result$residuals, residuals, tolerance = 1.0e-6)
expect_equal(result$weights, w)
object <- loggompertz_new(x, y, w, c(0, 1, 1, 1), max_iter, NULL, NULL)
result <- fit(object)
expect_true(inherits(result, "loggompertz_fit"))
expect_true(inherits(result, "loglogistic"))
expect_true(result$converged)
expect_false(result$constrained)
expect_equal(result$estimated, estimated)
expect_equal(result$coefficients, theta, tolerance = 1.0e-6)
expect_equal(result$rss, rss_value)
expect_equal(result$df.residual, object$n - 4)
expect_equal(result$fitted.values, fitted_values, tolerance = 1.0e-6)
expect_equal(result$residuals, residuals, tolerance = 1.0e-6)
expect_equal(result$weights, w)
})
test_that("fit_constrained: inequalities", {
x <- lltd$D$x
y <- lltd$D$y
n <- length(y)
w <- rep(1, n)
k <- as.numeric(table(x))
max_iter <- 10000
estimated <- c(alpha = TRUE, delta = TRUE, eta = TRUE, phi = TRUE)
rss_value <- 0.073630385904407914
theta <- c(
alpha = 0.83777215113750696, delta = -0.75729507682865963,
eta = exp(0.63524784841453301), phi = exp(1.6101961281319852)
)
fitted_values <- rep(
c(
0.83777215113750696, 0.83509727314890451, 0.67312406791826807,
0.46539726317063078, 0.33641138140141234, 0.26005233684193214,
0.083128285791155688, 0.080511485993225808
),
k
)
residuals <- c(
0.014927848862493037, -0.080672151137506963, 0.099227848862493037,
-0.0075972731489045085, -0.057297273148904508, 0.049502726851095492,
-0.11702406791826807, 0.031675932081731933, 0.079902736829369223,
0.020602736829369223, 0.069702736829369223, -0.025497263170630777,
0.018888618598587661, -0.023211381401412339, 0.013488618598587661,
-0.11795233684193214, -0.066328285791155688, 0.055471714208844312,
0.042188514006774192
)
object <- loggompertz_new(
x, y, w, NULL, max_iter,
c(0.5, -1, 1, 3), c(1, -0.5, 5, 12)
)
result <- fit_constrained(object)
expect_true(inherits(result, "loggompertz_fit"))
expect_true(inherits(result, "loglogistic"))
expect_true(result$converged)
expect_true(result$constrained)
expect_equal(result$estimated, estimated)
expect_equal(result$coefficients, theta, tolerance = 1.0e-6)
expect_equal(result$rss, rss_value)
expect_equal(result$df.residual, object$n - 4)
expect_equal(result$fitted.values, fitted_values, tolerance = 1.0e-6)
expect_equal(result$residuals, residuals, tolerance = 1.0e-6)
expect_equal(result$weights, w)
# initial values within the boundaries
object <- loggompertz_new(
x, y, w, c(0.6, -0.6, 4, 8), max_iter,
c(0.5, -1, 1, 3), c(1, -0.5, 5, 12)
)
result <- fit_constrained(object)
expect_true(inherits(result, "loggompertz_fit"))
expect_true(inherits(result, "loglogistic"))
expect_true(result$converged)
expect_true(result$constrained)
expect_equal(result$estimated, estimated)
expect_equal(result$coefficients, theta, tolerance = 1.0e-6)
expect_equal(result$rss, rss_value)
expect_equal(result$df.residual, object$n - 4)
expect_equal(result$fitted.values, fitted_values, tolerance = 1.0e-6)
expect_equal(result$residuals, residuals, tolerance = 1.0e-6)
expect_equal(result$weights, w)
# initial values outside the boundaries
object <- loggompertz_new(
x, y, w, c(-2, 2, 7, 1), max_iter,
c(0.5, -1, 1, 3), c(1, -0.5, 5, 12)
)
result <- fit_constrained(object)
expect_true(inherits(result, "loggompertz_fit"))
expect_true(inherits(result, "loglogistic"))
expect_true(result$converged)
expect_true(result$constrained)
expect_equal(result$estimated, estimated)
expect_equal(result$coefficients, theta, tolerance = 1.0e-6)
expect_equal(result$rss, rss_value)
expect_equal(result$df.residual, object$n - 4)
expect_equal(result$fitted.values, fitted_values, tolerance = 1.0e-6)
expect_equal(result$residuals, residuals, tolerance = 1.0e-6)
expect_equal(result$weights, w)
})
test_that("fit_constrained: equalities", {
x <- lltd$D$x
y <- lltd$D$y
n <- length(y)
w <- rep(1, n)
k <- as.numeric(table(x))
max_iter <- 10000
estimated <- c(alpha = FALSE, delta = FALSE, eta = TRUE, phi = TRUE)
rss_value <- 0.17654034036023363
theta <- c(
alpha = 0.8, delta = -0.9, eta = exp(0.34256135717772634),
phi = exp(1.8046122081510387)
)
fitted_values <- rep(
c(
0.8, 0.79248376718906222, 0.65161512055560499, 0.47490216690684499,
0.34359567627384339, 0.25186908386433202, -0.082746689754940315,
-0.099320242477951822
),
k
)
residuals <- c(
0.0527, -0.0429, 0.137, 0.035016232810937784, -0.014683767189062216,
0.092116232810937784, -0.095515120555604989, 0.053184879444395011,
0.070397833093155009, 0.011097833093155009, 0.060197833093155009,
-0.035002166906844991, 0.011704323726156605, -0.030395676273843395,
0.0063043237261566051, -0.10976908386433202, 0.099546689754940315,
0.22134668975494032, 0.22202024247795182
)
object <- loggompertz_new(
x, y, w, NULL, max_iter,
c(0.8, -0.9, rep(-Inf, 2)), c(0.8, -0.9, rep(Inf, 2))
)
result <- fit_constrained(object)
expect_true(inherits(result, "loggompertz_fit"))
expect_true(inherits(result, "loglogistic"))
expect_true(result$converged)
expect_false(result$constrained)
expect_equal(result$estimated, estimated)
expect_equal(result$coefficients, theta, tolerance = 1.0e-6)
expect_equal(result$rss, rss_value)
expect_equal(result$df.residual, object$n - 2)
expect_equal(result$fitted.values, fitted_values, tolerance = 1.0e-6)
expect_equal(result$residuals, residuals, tolerance = 1.0e-6)
expect_equal(result$weights, w)
# initial values with same equalities
object <- loggompertz_new(
x, y, w, c(0.8, -0.9, 1, 1), max_iter,
c(0.8, -0.9, rep(-Inf, 2)), c(0.8, -0.9, rep(Inf, 2))
)
result <- fit_constrained(object)
expect_true(inherits(result, "loggompertz_fit"))
expect_true(inherits(result, "loglogistic"))
expect_true(result$converged)
expect_false(result$constrained)
expect_equal(result$estimated, estimated)
expect_equal(result$coefficients, theta, tolerance = 1.0e-6)
expect_equal(result$rss, rss_value)
expect_equal(result$df.residual, object$n - 2)
expect_equal(result$fitted.values, fitted_values, tolerance = 1.0e-6)
expect_equal(result$residuals, residuals, tolerance = 1.0e-6)
expect_equal(result$weights, w)
# initial values with different equalities
object <- loggompertz_new(
x, y, w, c(0, 1, 1, 1), max_iter,
c(0.8, -0.9, rep(-Inf, 2)), c(0.8, -0.9, rep(Inf, 2))
)
result <- fit_constrained(object)
expect_true(inherits(result, "loggompertz_fit"))
expect_true(inherits(result, "loglogistic"))
expect_true(result$converged)
expect_false(result$constrained)
expect_equal(result$estimated, estimated)
expect_equal(result$coefficients, theta, tolerance = 1.0e-6)
expect_equal(result$rss, rss_value)
expect_equal(result$df.residual, object$n - 2)
expect_equal(result$fitted.values, fitted_values, tolerance = 1.0e-6)
expect_equal(result$residuals, residuals, tolerance = 1.0e-6)
expect_equal(result$weights, w)
})
test_that("fit_constrained: equalities and inequalities", {
x <- lltd$D$x
y <- lltd$D$y
n <- length(y)
w <- rep(1, n)
k <- as.numeric(table(x))
max_iter <- 10000
estimated <- c(alpha = FALSE, delta = FALSE, eta = TRUE, phi = TRUE)
rss_value <- 0.17654034036023363
theta <- c(
alpha = 0.8, delta = -0.9, eta = exp(0.34256135717772634),
phi = exp(1.8046122081510387)
)
fitted_values <- rep(
c(
0.8, 0.79248376718906222, 0.65161512055560499, 0.47490216690684499,
0.34359567627384339, 0.25186908386433202, -0.082746689754940315,
-0.099320242477951822
),
k
)
residuals <- c(
0.0527, -0.0429, 0.137, 0.035016232810937784, -0.014683767189062216,
0.092116232810937784, -0.095515120555604989, 0.053184879444395011,
0.070397833093155009, 0.011097833093155009, 0.060197833093155009,
-0.035002166906844991, 0.011704323726156605, -0.030395676273843395,
0.0063043237261566051, -0.10976908386433202, 0.099546689754940315,
0.22134668975494032, 0.22202024247795182
)
object <- loggompertz_new(
x, y, w, NULL, max_iter,
c(0.8, -0.9, 1, 3), c(0.8, -0.9, 5, 12)
)
result <- fit_constrained(object)
expect_true(inherits(result, "loggompertz_fit"))
expect_true(inherits(result, "loglogistic"))
expect_true(result$converged)
expect_true(result$constrained)
expect_equal(result$estimated, estimated)
expect_equal(result$coefficients, theta, tolerance = 1.0e-6)
expect_equal(result$rss, rss_value)
expect_equal(result$df.residual, object$n - 2)
expect_equal(result$fitted.values, fitted_values, tolerance = 1.0e-6)
expect_equal(result$residuals, residuals, tolerance = 1.0e-6)
expect_equal(result$weights, w)
# initial values within the boundaries
object <- loggompertz_new(
x, y, w, c(0.8, -0.9, 3, 7), max_iter,
c(0.8, -0.9, 1, 3), c(0.8, -0.9, 5, 12)
)
result <- fit_constrained(object)
expect_true(inherits(result, "loggompertz_fit"))
expect_true(inherits(result, "loglogistic"))
expect_true(result$converged)
expect_true(result$constrained)
expect_equal(result$estimated, estimated)
expect_equal(result$coefficients, theta, tolerance = 1.0e-6)
expect_equal(result$rss, rss_value)
expect_equal(result$df.residual, object$n - 2)
expect_equal(result$fitted.values, fitted_values, tolerance = 1.0e-6)
expect_equal(result$residuals, residuals, tolerance = 1.0e-6)
expect_equal(result$weights, w)
# initial values outside the boundaries
object <- loggompertz_new(
x, y, w, c(0, 1, 8, 1), max_iter,
c(0.8, -0.9, 1, 3), c(0.8, -0.9, 5, 12)
)
result <- fit_constrained(object)
expect_true(inherits(result, "loggompertz_fit"))
expect_true(inherits(result, "loglogistic"))
expect_true(result$converged)
expect_true(result$constrained)
expect_equal(result$estimated, estimated)
expect_equal(result$coefficients, theta, tolerance = 1.0e-6)
expect_equal(result$rss, rss_value)
expect_equal(result$df.residual, object$n - 2)
expect_equal(result$fitted.values, fitted_values, tolerance = 1.0e-6)
expect_equal(result$residuals, residuals, tolerance = 1.0e-6)
expect_equal(result$weights, w)
})
test_that("fit (weighted)", {
x <- lltd$D$x
y <- lltd$D$y
w <- lltd$D$w
n <- length(y)
k <- as.numeric(table(x))
max_iter <- 10000
estimated <- c(alpha = TRUE, delta = TRUE, eta = TRUE, phi = TRUE)
rss_value <- 0.035138156284360913
theta <- c(
alpha = 0.84575428432937011, delta = -0.76656134360156587,
eta = exp(0.89238310459026989), phi = exp(1.6622042309903930)
)
fitted_values <- rep(
c(
0.84575428432937011, 0.84573608713048233, 0.73789083646643889,
0.47593010582855891, 0.31155933858113756, 0.22406801104864969,
0.079774485527624186, 0.079195048422139125
),
k
)
residuals <- c(
0.0069457156706298950, -0.088654284329370105, 0.091245715670629895,
-0.018236087130482331, -0.067936087130482331, 0.038863912869517669,
-0.18179083646643889, -0.033090836466438890, 0.069369894171441088,
0.010069894171441088, 0.059169894171441088, -0.036030105828558912,
0.043740661418862440, 0.0016406614188624399, 0.038340661418862440,
-0.081968011048649687, -0.062974485527624186, 0.058825514472375814,
0.043504951577860875
)
object <- loggompertz_new(x, y, w, NULL, max_iter, NULL, NULL)
result <- fit(object)
expect_true(inherits(result, "loggompertz_fit"))
expect_true(inherits(result, "loglogistic"))
expect_true(result$converged)
expect_false(result$constrained)
expect_equal(result$estimated, estimated)
expect_equal(result$coefficients, theta, tolerance = 1.0e-6)
expect_equal(result$rss, rss_value)
expect_equal(result$df.residual, object$n - 4)
expect_equal(result$fitted.values, fitted_values, tolerance = 1.0e-6)
expect_equal(result$residuals, residuals, tolerance = 1.0e-6)
expect_equal(result$weights, w)
object <- loggompertz_new(x, y, w, c(0, 1, 1, 1), max_iter, NULL, NULL)
result <- fit(object)
expect_true(inherits(result, "loggompertz_fit"))
expect_true(inherits(result, "loglogistic"))
expect_true(result$converged)
expect_false(result$constrained)
expect_equal(result$estimated, estimated)
expect_equal(result$coefficients, theta, tolerance = 1.0e-6)
expect_equal(result$rss, rss_value)
expect_equal(result$df.residual, object$n - 4)
expect_equal(result$fitted.values, fitted_values, tolerance = 1.0e-6)
expect_equal(result$residuals, residuals, tolerance = 1.0e-6)
expect_equal(result$weights, w)
})
test_that("fit_constrained (weighted): inequalities", {
x <- lltd$D$x
y <- lltd$D$y
w <- lltd$D$w
n <- length(y)
k <- as.numeric(table(x))
max_iter <- 10000
estimated <- c(alpha = TRUE, delta = TRUE, eta = TRUE, phi = TRUE)
rss_value <- 0.035138156284360913
theta <- c(
alpha = 0.84575428432937011, delta = -0.76656134360156587,
eta = exp(0.89238310459026989), phi = exp(1.6622042309903930)
)
fitted_values <- rep(
c(
0.84575428432937011, 0.84573608713048233, 0.73789083646643889,
0.47593010582855891, 0.31155933858113756, 0.22406801104864969,
0.079774485527624186, 0.079195048422139125
),
k
)
residuals <- c(
0.0069457156706298950, -0.088654284329370105, 0.091245715670629895,
-0.018236087130482331, -0.067936087130482331, 0.038863912869517669,
-0.18179083646643889, -0.033090836466438890, 0.069369894171441088,
0.010069894171441088, 0.059169894171441088, -0.036030105828558912,
0.043740661418862440, 0.0016406614188624399, 0.038340661418862440,
-0.081968011048649687, -0.062974485527624186, 0.058825514472375814,
0.043504951577860875
)
object <- loggompertz_new(
x, y, w, NULL, max_iter,
c(0.5, -1, 1, 3), c(1, -0.5, 5, 12)
)
result <- fit_constrained(object)
expect_true(inherits(result, "loggompertz_fit"))
expect_true(inherits(result, "loglogistic"))
expect_true(result$converged)
expect_true(result$constrained)
expect_equal(result$estimated, estimated)
expect_equal(result$coefficients, theta, tolerance = 1.0e-6)
expect_equal(result$rss, rss_value)
expect_equal(result$df.residual, object$n - 4)
expect_equal(result$fitted.values, fitted_values, tolerance = 1.0e-6)
expect_equal(result$residuals, residuals, tolerance = 1.0e-6)
expect_equal(result$weights, w)
# initial values within the boundaries
object <- loggompertz_new(
x, y, w, c(0.6, -0.6, 4, 8), max_iter,
c(0.5, -1, 1, 3), c(1, -0.5, 5, 12)
)
result <- fit_constrained(object)
expect_true(inherits(result, "loggompertz_fit"))
expect_true(inherits(result, "loglogistic"))
expect_true(result$converged)
expect_true(result$constrained)
expect_equal(result$estimated, estimated)
expect_equal(result$coefficients, theta, tolerance = 1.0e-6)
expect_equal(result$rss, rss_value)
expect_equal(result$df.residual, object$n - 4)
expect_equal(result$fitted.values, fitted_values, tolerance = 1.0e-6)
expect_equal(result$residuals, residuals, tolerance = 1.0e-6)
expect_equal(result$weights, w)
# initial values outside the boundaries
object <- loggompertz_new(
x, y, w, c(-2, 2, 7, 1), max_iter,
c(0.5, -1, 1, 3), c(1, -0.5, 5, 12)
)
result <- fit_constrained(object)
expect_true(inherits(result, "loggompertz_fit"))
expect_true(inherits(result, "loglogistic"))
expect_true(result$converged)
expect_true(result$constrained)
expect_equal(result$estimated, estimated)
expect_equal(result$coefficients, theta, tolerance = 1.0e-6)
expect_equal(result$rss, rss_value)
expect_equal(result$df.residual, object$n - 4)
expect_equal(result$fitted.values, fitted_values, tolerance = 1.0e-6)
expect_equal(result$residuals, residuals, tolerance = 1.0e-6)
expect_equal(result$weights, w)
})
test_that("fit_constrained (weighted): equalities", {
x <- lltd$D$x
y <- lltd$D$y
w <- lltd$D$w
n <- length(y)
k <- as.numeric(table(x))
max_iter <- 10000
estimated <- c(alpha = FALSE, delta = FALSE, eta = TRUE, phi = TRUE)
rss_value <- 0.12152790634234742
theta <- c(
alpha = 0.8, delta = -0.9, eta = exp(0.67344700327430744),
phi = exp(1.8267797417585254)
)
fitted_values <- rep(
c(
0.8, 0.79991222146934487, 0.71604576658544182, 0.49162778848529489,
0.31063536742142226, 0.19268868846616760, -0.096135387777275743,
-0.099957631352575174
),
k
)
residuals <- c(
0.0527, -0.0429, 0.137, 0.027587778530655130, -0.022112221469344870,
0.084687778530655130, -0.15994576658544182, -0.011245766585441821,
0.053672211514705107, -0.0056277884852948934, 0.043472211514705107,
-0.051727788485294893, 0.044664632578577739, 0.0025646325785777394,
0.039264632578577739, -0.050588688466167597, 0.11293538777727574,
0.23473538777727574, 0.22265763135257517
)
object <- loggompertz_new(
x, y, w, NULL, max_iter,
c(0.8, -0.9, rep(-Inf, 2)), c(0.8, -0.9, rep(Inf, 2))
)
result <- fit_constrained(object)
expect_true(inherits(result, "loggompertz_fit"))
expect_true(inherits(result, "loglogistic"))
expect_true(result$converged)
expect_false(result$constrained)
expect_equal(result$estimated, estimated)
expect_equal(result$coefficients, theta, tolerance = 1.0e-6)
expect_equal(result$rss, rss_value)
expect_equal(result$df.residual, object$n - 2)
expect_equal(result$fitted.values, fitted_values, tolerance = 1.0e-6)
expect_equal(result$residuals, residuals, tolerance = 1.0e-6)
expect_equal(result$weights, w)
# initial values with same equalities
object <- loggompertz_new(
x, y, w, c(0.8, -0.9, 1, 1), max_iter,
c(0.8, -0.9, rep(-Inf, 2)), c(0.8, -0.9, rep(Inf, 2))
)
result <- fit_constrained(object)
expect_true(inherits(result, "loggompertz_fit"))
expect_true(inherits(result, "loglogistic"))
expect_true(result$converged)
expect_false(result$constrained)
expect_equal(result$estimated, estimated)
expect_equal(result$coefficients, theta, tolerance = 1.0e-6)
expect_equal(result$rss, rss_value)
expect_equal(result$df.residual, object$n - 2)
expect_equal(result$fitted.values, fitted_values, tolerance = 1.0e-6)
expect_equal(result$residuals, residuals, tolerance = 1.0e-6)
expect_equal(result$weights, w)
# initial values with different equalities
object <- loggompertz_new(
x, y, w, c(0, 1, 1, 1), max_iter,
c(0.8, -0.9, rep(-Inf, 2)), c(0.8, -0.9, rep(Inf, 2))
)
result <- fit_constrained(object)
expect_true(inherits(result, "loggompertz_fit"))
expect_true(inherits(result, "loglogistic"))
expect_true(result$converged)
expect_false(result$constrained)
expect_equal(result$estimated, estimated)
expect_equal(result$coefficients, theta, tolerance = 1.0e-6)
expect_equal(result$rss, rss_value)
expect_equal(result$df.residual, object$n - 2)
expect_equal(result$fitted.values, fitted_values, tolerance = 1.0e-6)
expect_equal(result$residuals, residuals, tolerance = 1.0e-6)
expect_equal(result$weights, w)
})
test_that("fit_constrained (weighted): equalities and inequalities", {
x <- lltd$D$x
y <- lltd$D$y
w <- lltd$D$w
n <- length(y)
k <- as.numeric(table(x))
max_iter <- 10000
estimated <- c(alpha = FALSE, delta = FALSE, eta = TRUE, phi = TRUE)
rss_value <- 0.12152790634234742
theta <- c(
alpha = 0.8, delta = -0.9, eta = exp(0.67344700327430744),
phi = exp(1.8267797417585254)
)
fitted_values <- rep(
c(
0.8, 0.79991222146934487, 0.71604576658544182, 0.49162778848529489,
0.31063536742142226, 0.19268868846616760, -0.096135387777275743,
-0.099957631352575174
),
k
)
residuals <- c(
0.0527, -0.0429, 0.137, 0.027587778530655130, -0.022112221469344870,
0.084687778530655130, -0.15994576658544182, -0.011245766585441821,
0.053672211514705107, -0.0056277884852948934, 0.043472211514705107,
-0.051727788485294893, 0.044664632578577739, 0.0025646325785777394,
0.039264632578577739, -0.050588688466167597, 0.11293538777727574,
0.23473538777727574, 0.22265763135257517
)
object <- loggompertz_new(
x, y, w, NULL, max_iter,
c(0.8, -0.9, 1, 3), c(0.8, -0.9, 5, 12)
)
result <- fit_constrained(object)
expect_true(inherits(result, "loggompertz_fit"))
expect_true(inherits(result, "loglogistic"))
expect_true(result$converged)
expect_true(result$constrained)
expect_equal(result$estimated, estimated)
expect_equal(result$coefficients, theta, tolerance = 1.0e-6)
expect_equal(result$rss, rss_value)
expect_equal(result$df.residual, object$n - 2)
expect_equal(result$fitted.values, fitted_values, tolerance = 1.0e-6)
expect_equal(result$residuals, residuals, tolerance = 1.0e-6)
expect_equal(result$weights, w)
# initial values within the boundaries
object <- loggompertz_new(
x, y, w, c(0.8, -0.9, 3, 7), max_iter,
c(0.8, -0.9, 1, 3), c(0.8, -0.9, 5, 12)
)
result <- fit_constrained(object)
expect_true(inherits(result, "loggompertz_fit"))
expect_true(inherits(result, "loglogistic"))
expect_true(result$converged)
expect_true(result$constrained)
expect_equal(result$estimated, estimated)
expect_equal(result$coefficients, theta, tolerance = 1.0e-6)
expect_equal(result$rss, rss_value)
expect_equal(result$df.residual, object$n - 2)
expect_equal(result$fitted.values, fitted_values, tolerance = 1.0e-6)
expect_equal(result$residuals, residuals, tolerance = 1.0e-6)
expect_equal(result$weights, w)
# initial values outside the boundaries
object <- loggompertz_new(
x, y, w, c(0, 1, 8, 1), max_iter,
c(0.8, -0.9, 1, 3), c(0.8, -0.9, 5, 12)
)
result <- fit_constrained(object)
expect_true(inherits(result, "loggompertz_fit"))
expect_true(inherits(result, "loglogistic"))
expect_true(result$converged)
expect_true(result$constrained)
expect_equal(result$estimated, estimated)
expect_equal(result$coefficients, theta, tolerance = 1.0e-6)
expect_equal(result$rss, rss_value)
expect_equal(result$df.residual, object$n - 2)
expect_equal(result$fitted.values, fitted_values, tolerance = 1.0e-6)
expect_equal(result$residuals, residuals, tolerance = 1.0e-6)
expect_equal(result$weights, w)
})
test_that("fisher_info", {
x <- lltd$D$x
y <- lltd$D$y
w <- lltd$D$w
max_iter <- 10000
theta <- lltd$theta_g
names(theta) <- c("alpha", "delta", "eta", "phi")
sigma <- lltd$sigma
true_value <- matrix(c(
# alpha
6206.96, 2748.9449629599373, -189.18362192800753, 344.04436069543435,
-28355.127259362132,
# delta
2748.9449629599373, 2023.9485432497464, -124.69711327516109,
134.06347723771206, -10502.299351663864,
# eta
-189.18362192800753, -124.69711327516109, 23.773250144905426,
-11.783571918074014, 850.51618712697142,
# phi
344.04436069543435, 134.06347723771206, -11.783571918074014,
29.686894742510176, -1441.3142462178932,
# sigma
-28355.127259362132, -10502.299351663864, 850.51618712697142,
-1441.3142462178932, 113909.14953385276
),
nrow = 5,
ncol = 5
)
rownames(true_value) <- colnames(true_value) <- c(
"alpha", "delta", "eta", "phi", "sigma"
)
object <- loggompertz_new(x, y, w, NULL, max_iter, NULL, NULL)
fim <- fisher_info(object, theta, sigma)
expect_type(fim, "double")
expect_length(fim, 5 * 5)
expect_equal(fim, true_value)
})
test_that("drda: 'lower_bound' argument errors", {
x <- lltd$D$x
y <- lltd$D$y
expect_error(
drda(
y ~ x, mean_function = "loggompertz",
lower_bound = c("a", "b", "c", "d")
),
"'lower_bound' must be a numeric vector"
)
expect_error(
drda(
y ~ x, mean_function = "loggompertz",
lower_bound = matrix(-Inf, nrow = 4, ncol = 2),
upper_bound = rep(Inf, 4)
),
"'lower_bound' must be a numeric vector"
)
expect_error(
drda(
y ~ x, mean_function = "loggompertz",
lower_bound = rep(-Inf, 5),
upper_bound = rep(Inf, 4)
),
"'lower_bound' and 'upper_bound' must have the same length"
)
expect_error(
drda(
y ~ x, mean_function = "loggompertz",
lower_bound = c( 0, -Inf, -Inf, -Inf),
upper_bound = c(-1, Inf, Inf, Inf)
),
"'lower_bound' cannot be larger than 'upper_bound'"
)
expect_error(
drda(
y ~ x, mean_function = "loggompertz",
lower_bound = c(Inf, -Inf, -Inf, -Inf),
upper_bound = c(Inf, Inf, Inf, Inf)
),
"'lower_bound' cannot be equal to infinity"
)
expect_error(
drda(
y ~ x, mean_function = "loggompertz",
lower_bound = rep(-Inf, 5),
upper_bound = rep(Inf, 5)
),
"'lower_bound' must be of length 4"
)
})
test_that("drda: 'upper_bound' argument errors", {
x <- lltd$D$x
y <- lltd$D$y
expect_error(
drda(
y ~ x, mean_function = "loggompertz",
upper_bound = c("a", "b", "c", "d")
),
"'upper_bound' must be a numeric vector"
)
expect_error(
drda(
y ~ x, mean_function = "loggompertz",
lower_bound = rep(-Inf, 4),
upper_bound = matrix(Inf, nrow = 4, ncol = 2)
),
"'upper_bound' must be a numeric vector"
)
expect_error(
drda(
y ~ x, mean_function = "loggompertz",
lower_bound = c(-Inf, -Inf, -Inf, -Inf),
upper_bound = c(-Inf, Inf, Inf, Inf)
),
"'upper_bound' cannot be equal to -infinity"
)
expect_error(
drda(
y ~ x, mean_function = "loggompertz",
lower_bound = rep(-Inf, 5),
upper_bound = rep(Inf, 5)
),
"'lower_bound' must be of length 4"
)
})
test_that("drda: 'start' argument errors", {
x <- lltd$D$x
y <- lltd$D$y
expect_error(
drda(
y ~ x, mean_function = "loggompertz",
start = c("a", "b", "c", "d")
),
"'start' must be a numeric vector"
)
expect_error(
drda(
y ~ x, mean_function = "loggompertz",
start = c(0, Inf, 1, 1)
),
"'start' must be finite"
)
expect_error(
drda(
y ~ x, mean_function = "loggompertz",
start = c(-Inf, 1, 1, 1)
),
"'start' must be finite"
)
expect_error(
drda(
y ~ x, mean_function = "loggompertz",
start = rep(1, 5)
),
"'start' must be of length 4"
)
expect_error(
drda(
y ~ x, mean_function = "loggompertz",
start = c(0, 1, -1, 1)
),
"parameter 'eta' cannot be negative nor zero"
)
expect_error(
drda(
y ~ x, mean_function = "loggompertz",
start = c(0, 1, 0, 1)
),
"parameter 'eta' cannot be negative nor zero"
)
expect_error(
drda(
y ~ x, mean_function = "loggompertz",
start = c(0, 1, 1, -1)
),
"parameter 'phi' cannot be negative nor zero"
)
expect_error(
drda(
y ~ x, mean_function = "loggompertz",
start = c(0, 1, 1, 0)
),
"parameter 'phi' cannot be negative nor zero"
)
})
test_that("nauc: decreasing", {
x <- lltd$D$x
y <- lltd$D$y
w <- lltd$D$w
result <- drda(y ~ x, weights = w, mean_function = "loggompertz")
expect_equal(nauc(result), 0.085299967613769906, tolerance = 1.0e-6)
expect_equal(nauc(result, xlim = c(0, 2)), 0.84575366076886644)
expect_equal(nauc(result, ylim = c(0.3, 0.7)), 0.0059283502711447546)
expect_equal(nauc(result, xlim = c(0, 2), ylim = c(0.3, 0.7)), 1.0)
expect_equal(
nauc(result, xlim = c(5, 8), ylim = c(0.3, 0.7)), 0.33844687600885550
)
expect_equal(nauc(result, xlim = c(10, 15), ylim = c(0.3, 0.7)), 0.0)
})
test_that("naac: decreasing", {
x <- lltd$D$x
y <- lltd$D$y
w <- lltd$D$w
result <- drda(y ~ x, weights = w, mean_function = "loggompertz")
expect_equal(naac(result), 1 - 0.085299967613769906)
expect_equal(naac(result, xlim = c(0, 2)), 1 - 0.84575366076886644)
expect_equal(naac(result, ylim = c(0.3, 0.7)), 1 - 0.0059283502711447546)
expect_equal(naac(result, xlim = c(0, 2), ylim = c(0.3, 0.7)), 0.0)
expect_equal(
naac(result, xlim = c(5, 8), ylim = c(0.3, 0.7)), 1 - 0.33844687600885550
)
expect_equal(naac(result, xlim = c(10, 15), ylim = c(0.3, 0.7)), 1.0)
})
test_that("nauc: increasing", {
x <- lltd$D$x
y <- rev(lltd$D$y)
w <- lltd$D$w
result <- drda(y ~ x, weights = w, mean_function = "loggompertz")
expect_equal(nauc(result), 0.88478280795251264)
expect_equal(nauc(result, xlim = c(0, 2)), 0.129049490857619206)
expect_equal(nauc(result, ylim = c(0.3, 0.7)), 0.99551714173529925)
expect_equal(nauc(result, xlim = c(0, 2), ylim = c(0.3, 0.7)), 0.0)
expect_equal(
nauc(result, xlim = c(5, 8), ylim = c(0.3, 0.7)), 0.85655631688941921
)
expect_equal(nauc(result, xlim = c(9, 12), ylim = c(0.3, 0.7)), 1.0)
})
test_that("naac: increasing", {
x <- lltd$D$x
y <- rev(lltd$D$y)
w <- lltd$D$w
result <- drda(y ~ x, weights = w, mean_function = "loggompertz")
expect_equal(naac(result), 1 - 0.88478280795251264)
expect_equal(naac(result, xlim = c(0, 2)), 1 - 0.129049490857619206)
expect_equal(naac(result, ylim = c(0.3, 0.7)), 1 - 0.99551714173529925)
expect_equal(naac(result, xlim = c(0, 2), ylim = c(0.3, 0.7)), 1.0)
expect_equal(
naac(result, xlim = c(5, 8), ylim = c(0.3, 0.7)), 1 - 0.85655631688941921
)
expect_equal(naac(result, xlim = c(9, 12), ylim = c(0.3, 0.7)), 0.0)
})
Try the drda 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.