Package {splitpopsurv}

splitpopsurv: Split-Population (Cure / Mover-Stayer) Survival Models

Description

Maximum-likelihood estimation of split-population survival models, which combine an accelerated failure-time regression for event timing among "movers" with a logistic regression on the probability of belonging to the immune "stayer" population. See Schmidt & Witte (1989) and Yamaguchi (1992, 1998) for the underlying theory, and the package README for a full manual including a likelihood correction relative to the Stata 'ml' programs this package translates.

Main functions

[fit_splitpop_loglogistic()], [fit_splitpop_weibull()], [fit_splitpop_lognormal()], [fit_splitpop_gamma()], [fit_splitpop_ggamma()].

Author(s)

Maintainer: Nobutaka Fukuda nobutaka.fukuda@tohoku.ac.jp

Authors:

• Nobutaka Fukuda nobutaka.fukuda@tohoku.ac.jp

See Also

Useful links:

• https://github.com/nobifukuda/splitpopsurv

• Report bugs at https://github.com/nobifukuda/splitpopsurv/issues

Fit a split-population gamma survival model

Description

Maximum-likelihood estimation of a split-population (cure) survival model with a gamma baseline for the "mover" population. R translation of the Stata program 'SphGam'; see 'docs/manual.html' for the model and a note on a likelihood correction relative to the original Stata code.

Usage

fit_splitpop_gamma(
  hform,
  pform,
  data,
  time,
  event,
  group,
  start = NULL,
  method = "NR"
)

Arguments

Details

Note: unlike the other four models, the H_regression linear index here is used directly as a rate (no 'exp()' transform), so it must stay positive; the default 'start' sets a positive intercept and zero slopes to keep it safe regardless of covariate values.

Value

A 'maxLik' object; use 'summary()', 'coef()', 'logLik()', 'vcov()'.

Examples

set.seed(1)
n <- 150
x1 <- rnorm(n); x2 <- rbinom(n, 1, 0.5)
cured <- rbinom(n, 1, plogis(0.2 + 0.5 * x2))
t_latent <- (-log(runif(n)) / exp(0.5 - 0.4 * x1))^(1 / 1.3)
censor_time <- rexp(n, rate = 0.2)
time  <- pmax(ifelse(cured == 1, censor_time, pmin(t_latent, censor_time)), 1e-3)
event <- ifelse(cured == 1, 0, as.numeric(t_latent <= censor_time))
mydata <- data.frame(time = time, event = event, x1 = x1, x2 = x2, group = 1)

fit <- fit_splitpop_gamma(~x1, ~x2, mydata,
                           "time", "event", "group", method = "BFGS")
summary(fit)

Fit a split-population generalized gamma survival model

Description

Maximum-likelihood estimation of a split-population (cure) survival model with a generalized gamma baseline (Prentice 1974 / Yamaguchi & Ferguson 1995, note 10) for the "mover" population – the family that nests the Weibull (kappa=1), log-normal (kappa=0), and gamma (sigma=1) models above. R translation of the Stata program 'SphGGam' ('d0' method); see 'docs/manual.html' for the model and a note on a likelihood correction relative to the original Stata code.

Usage

fit_splitpop_ggamma(
  hform,
  pform,
  data,
  time,
  event,
  group,
  start = NULL,
  method = "NR"
)

Arguments

Details

Note: '|kappa| < 0.01' switches to a log-normal-equivalent branch, a genuine kink in the likelihood surface (preserved from the Stata guard). If 'summary()' reports a non-finite standard error, the optimizer likely landed near that threshold – try a different starting 'kappa'.

Value

A 'maxLik' object; use 'summary()', 'coef()', 'logLik()', 'vcov()'.

Examples

set.seed(1)
n <- 150
x1 <- rnorm(n); x2 <- rbinom(n, 1, 0.5)
cured <- rbinom(n, 1, plogis(0.2 + 0.5 * x2))
t_latent <- (-log(runif(n)) / exp(0.5 - 0.4 * x1))^(1 / 1.3)
censor_time <- rexp(n, rate = 0.2)
time  <- pmax(ifelse(cured == 1, censor_time, pmin(t_latent, censor_time)), 1e-3)
event <- ifelse(cured == 1, 0, as.numeric(t_latent <= censor_time))
mydata <- data.frame(time = time, event = event, x1 = x1, x2 = x2, group = 1)

fit <- fit_splitpop_ggamma(~x1, ~x2, mydata,
                            "time", "event", "group", method = "BFGS")
summary(fit)

Fit a split-population log-logistic survival model

Description

Maximum-likelihood estimation of a split-population (cure) survival model with a log-logistic baseline for the "mover" (susceptible) population and a logistic regression on the probability of being a "stayer" (immune to the event). R translation of the Stata program 'SphLog'; see 'docs/manual.html' for the model, formulas, and a note on a likelihood correction relative to the original Stata code.

Usage

fit_splitpop_loglogistic(
  hform,
  pform,
  data,
  time,
  event,
  group,
  start = NULL,
  method = "NR"
)

Arguments

Value

A 'maxLik' object; use 'summary()', 'coef()', 'logLik()', 'vcov()'.

Examples

set.seed(1)
n <- 150
x1 <- rnorm(n); x2 <- rbinom(n, 1, 0.5)
cured <- rbinom(n, 1, plogis(0.2 + 0.5 * x2))
t_latent <- (-log(runif(n)) / exp(0.5 - 0.4 * x1))^(1 / 1.3)
censor_time <- rexp(n, rate = 0.2)
time  <- pmax(ifelse(cured == 1, censor_time, pmin(t_latent, censor_time)), 1e-3)
event <- ifelse(cured == 1, 0, as.numeric(t_latent <= censor_time))
mydata <- data.frame(time = time, event = event, x1 = x1, x2 = x2, group = 1)

fit <- fit_splitpop_loglogistic(~x1, ~x2, mydata,
                                 "time", "event", "group", method = "BFGS")
summary(fit)

Fit a split-population log-normal survival model

Description

Maximum-likelihood estimation of a split-population (cure) survival model with a log-normal baseline for the "mover" population. R translation of the Stata program 'SphNom'; see 'docs/manual.html' for the model and a note on a likelihood correction relative to the original Stata code.

Usage

fit_splitpop_lognormal(
  hform,
  pform,
  data,
  time,
  event,
  group,
  start = NULL,
  method = "NR"
)

Arguments

Value

A 'maxLik' object; use 'summary()', 'coef()', 'logLik()', 'vcov()'.

Examples

set.seed(1)
n <- 150
x1 <- rnorm(n); x2 <- rbinom(n, 1, 0.5)
cured <- rbinom(n, 1, plogis(0.2 + 0.5 * x2))
t_latent <- (-log(runif(n)) / exp(0.5 - 0.4 * x1))^(1 / 1.3)
censor_time <- rexp(n, rate = 0.2)
time  <- pmax(ifelse(cured == 1, censor_time, pmin(t_latent, censor_time)), 1e-3)
event <- ifelse(cured == 1, 0, as.numeric(t_latent <= censor_time))
mydata <- data.frame(time = time, event = event, x1 = x1, x2 = x2, group = 1)

fit <- fit_splitpop_lognormal(~x1, ~x2, mydata,
                               "time", "event", "group", method = "BFGS")
summary(fit)

Fit a split-population Weibull survival model

Description

Maximum-likelihood estimation of a split-population (cure) survival model with a Weibull baseline for the "mover" population. R translation of the Stata program 'SphWieb'; see 'docs/manual.html' for the model and a note on a likelihood correction relative to the original Stata code.

Usage

fit_splitpop_weibull(
  hform,
  pform,
  data,
  time,
  event,
  group,
  start = NULL,
  method = "NR"
)

Arguments

Value

A 'maxLik' object; use 'summary()', 'coef()', 'logLik()', 'vcov()'.

Examples

set.seed(1)
n <- 150
x1 <- rnorm(n); x2 <- rbinom(n, 1, 0.5)
cured <- rbinom(n, 1, plogis(0.2 + 0.5 * x2))
t_latent <- (-log(runif(n)) / exp(0.5 - 0.4 * x1))^(1 / 1.3)
censor_time <- rexp(n, rate = 0.2)
time  <- pmax(ifelse(cured == 1, censor_time, pmin(t_latent, censor_time)), 1e-3)
event <- ifelse(cured == 1, 0, as.numeric(t_latent <= censor_time))
mydata <- data.frame(time = time, event = event, x1 = x1, x2 = x2, group = 1)

fit <- fit_splitpop_weibull(~x1, ~x2, mydata,
                             "time", "event", "group", method = "BFGS")
summary(fit)

Log-likelihood: split-population gamma model

Description

Log-likelihood: split-population gamma model

Usage

loglik_splitpop_gamma(par, time, event, group, Xh, Xp)

Arguments

Value

A numeric vector of per-observation log-likelihood contributions.

Examples

Xh <- cbind(1, rnorm(20)); Xp <- cbind(1, rbinom(20, 1, 0.5))
par0 <- c(1, 0, 1, 0, 0)  # positive H_regression intercept: theta1 must stay > 0
loglik_splitpop_gamma(par0, time = runif(20, 0.1, 5),
                       event = rbinom(20, 1, 0.7),
                       group = rep(1, 20), Xh = Xh, Xp = Xp)

Log-likelihood: split-population generalized gamma model

Description

Log-likelihood: split-population generalized gamma model

Usage

loglik_splitpop_ggamma(par, time, event, group, Xh, Xp)

Arguments

Value

A numeric vector of per-observation log-likelihood contributions.

Examples

Xh <- cbind(1, rnorm(20)); Xp <- cbind(1, rbinom(20, 1, 0.5))
par0 <- c(0, 0, 0, 1, 0, 0)
loglik_splitpop_ggamma(par0, time = runif(20, 0.1, 5),
                        event = rbinom(20, 1, 0.7),
                        group = rep(1, 20), Xh = Xh, Xp = Xp)

Log-likelihood: split-population log-logistic model

Description

Internal log-likelihood used by [fit_splitpop_loglogistic()]. Exposed so it can be evaluated by hand at candidate starting values (the R analogue of Stata's 'ml check').

Usage

loglik_splitpop_loglogistic(par, time, event, group, Xh, Xp)

Arguments

Value

A numeric vector of per-observation log-likelihood contributions.

Examples

Xh <- cbind(1, rnorm(20)); Xp <- cbind(1, rbinom(20, 1, 0.5))
par0 <- c(0, 0, 1, 0, 0)
loglik_splitpop_loglogistic(par0, time = runif(20, 0.1, 5),
                             event = rbinom(20, 1, 0.7),
                             group = rep(1, 20), Xh = Xh, Xp = Xp)

Log-likelihood: split-population log-normal model

Description

Log-likelihood: split-population log-normal model

Usage

loglik_splitpop_lognormal(par, time, event, group, Xh, Xp)

Arguments

Value

A numeric vector of per-observation log-likelihood contributions.

Examples

Xh <- cbind(1, rnorm(20)); Xp <- cbind(1, rbinom(20, 1, 0.5))
par0 <- c(0, 0, 0, 0, 0)
loglik_splitpop_lognormal(par0, time = runif(20, 0.1, 5),
                           event = rbinom(20, 1, 0.7),
                           group = rep(1, 20), Xh = Xh, Xp = Xp)

Log-likelihood: split-population Weibull model

Description

Log-likelihood: split-population Weibull model

Usage

loglik_splitpop_weibull(par, time, event, group, Xh, Xp)

Arguments

Value

A numeric vector of per-observation log-likelihood contributions.

Examples

Xh <- cbind(1, rnorm(20)); Xp <- cbind(1, rbinom(20, 1, 0.5))
par0 <- c(0, 0, 0, 0, 0)
loglik_splitpop_weibull(par0, time = runif(20, 0.1, 5),
                         event = rbinom(20, 1, 0.7),
                         group = rep(1, 20), Xh = Xh, Xp = Xp)