MM algorithm based on the AD method for multivariate compound zero-inflated generalized poisson distribution

Description

Let Z_0 \sim Bernoulli(1- \phi_0), \bf{x} = (X_1,\cdots, X_m)^T, X_i \sim ZIGP(\phi_i, \lambda_i, \theta_i), for i = 1,\cdots,m , and (Z_0,X_1,\cdots, X_m) be mutually independent. A random vector \bf{y}=(Y_1,\cdots, Y_m)^T follows a multivariate compound zero-inflated generalized poisson distribution if

\bf{y} \overset{\rm{d}}= Z_{0}\bf{x}= \left\{ \begin{array}{ll} \bf{0} & \text{with probability} \ \lambda_0 \\ \bf{x} & \text{with probability} \ 1-\lambda_0 \\ \end{array} \right.

where \lambda_0 \in [0,1), \bf{\phi}=(\phi_1,\cdots, \phi_m)^{T} \in [0, 1)^m, \bf{\lambda}=(\lambda_1,\cdots, \lambda_m)^{T} \in \mathbb{R}_{+}^{m} , \bf{\theta}=(\theta_1,\cdots, \theta_m)^{T} \in [0, 1)^m. The CZIGPMM function is used to calculate the multivariate compound ZIGP model.

Usage

CZIGPMM(data, phi0, phi, la, th, Maxiter = 2000, convergence = 1e-06, ...)

Arguments

Details

The CZIGPMM function is used to calculate multivariate compound zero-inflated generalized poisson distribution model using MM algorithms based on AD technology. data is provided by user by default, it can be a data frame or a matrix. In addition, unknown parameters require users to give appropriate initial values, where 0<=phi0<1, each phi should 0<=phi<1, th should 0<=th<1, and each la should be greater than 0.

Value

An object of class CZIGPMM that contains the following fields: total amount of observations, the number of iterations, convergence rate, the log likelihood value, estimated results for the unknown parameters, the standard deviation of estimate for the unknown parameters, the likelihood-based 95% confidence interval for the unknown parameters, information criterion: AIC value and BIC value.

References

Tian G.L., Huang X.F. and Xu, J.(2019). 'An assembly and decomposition approach for constructing separable minorizing functions in a class of MM algorithms.' Statistica Sinica 29(2), 961-982.

Huang X.F., Tian G.L., Zhang, C. and Jiang, X.(2017). 'Type I multivariate zero-inflated generalized Poisson distribution with applications.' Statistics and its Interface 10(2), 291-311.

Examples

x1 <- c(0,35,23,34,8,19,0,0,0,0)
x2 <- c(38,15,0,25,34,0,0,0,0,0)
y <- cbind(x1, x2)
phi0 = 0.5; phi = rep(0.5,2); la = rep(1,2); th = rep(0.1,2)
CZIGPMM(y, phi0, phi, la, th)

MM algorithm based on AD technology for Cox model

Description

Let T_i, C_i and X_i = (x_{i1},\cdots, x_{iq})^T denote the, survival time, the censoring time and a q dimension vector of coefficients for the i-th individual, respectively. And assume the censoring time C_i is independent of the survival time T_i are mutually independent, and I_i = I(T_{i} \leqslant C_{i}) is the censoring indicator. Then the instantaneous hazard rate function of T_i is

\lambda(t|X_i)=\lambda_{0}(t) \exp(X_{i}^{T} \beta)

where \lambda_{0}(.) is a baseline hazard rate and \beta = (\beta_1, \cdots, \beta_q)^{T} is a vector of regression parameters. We denote \Lambda as the accumulative hazard rate. Then the observed data likelihood function is

L(\alpha | Y_{obs}) = \prod_{i=1}^n (\lambda_{0}(t_i) \exp(X_{i}^{T} \beta))^{I_i} \exp(-\Lambda(t_i) \exp(X_{i}^{T} \beta))

where \alpha = (\beta, \Lambda). The CoxMM function is used to calculate the Cox model.

Usage

CoxMM(formula, data, beta = NULL, Maxiter = 2000, convergence = 1e-06, ...)

Arguments

Details

The CoxMM function is used to calculate the Cox model using MM algorithms based on AD technology. EM algorithms rely on the fact that, after profiling out the nonparametric component \Lambda, the resulting function is concave. However, when this assumption does not hold, maximizing the resulting function using Newton’s method becomes difficult, especially when there are a large number of covariates. MM algorithms can avoid the concavity requirement and bypass the need for Newton method and matrix inversion.

Value

An object of class CoxMM that contains the following fields: the Time, total amount of observations, total number of failure events, the variable name, the \beta, the \lambda, the \Lambda, convergence result, the log likelihood value, the standard deviation of the estimated \beta, the likelihood-based 95% confidence interval for the \beta.

References

D.R. Cox.(1972). 'Regression models and life tables.' Journal of the Royal Statistical Society(Series B) 34(2), 187-220.

Zhang L.L. and Huang X.F.(2022). 'On MM algorithms for Cox model with right-censored data.' In International Conference on Cloud Computing, Internet of Things, and Computer Applications (CICA 2022) 12303, 29-38.

Examples

library(survival)
CoxMM(Surv(time, status) ~ age + sex, lung)

MM algorithm based on AD technology for gamma frailty model

Description

Let T_{ij}, C_{ij} and X_{ij} = (X_{ij1},\cdots, X_{ijq})^T denote the survival time, the censoring time and a vector of covariates, respectively. For the j-th individual in the i-th cluster, for j=1, \cdots, M_i and i=1, \cdots, B. And assume the censoring time C_{ij} is independent of the survival time T_{ij} given X_{ij}, and I_{ij} = I(T_{ij} \leqslant C_{ij}) is the censoring indicator. Conditional on a cluster-specific frailty \omega_i, then the frailty model postulates that the instantaneous hazard rate function of T_{ij} is

\lambda(t|X_{ij} , \omega_i )=\lambda_{0}(t) \exp(X_{ij}^{T} \beta) \omega_i

where \lambda_{0}(t) is a baseline hazard rate and \beta is a vector of regression parameters. We assume that the frailty \omega has a gamma distribution with mean 1, variance \theta and density

g(\omega) = \frac{\omega^{\frac{1}{1-\theta}} \exp(-\frac{\omega}{\theta})}{\Gamma( \frac{1}{\theta}) \theta^{\frac{1}{\theta}}}

and we denote \Lambda as the accumulative hazard rate. The GaFrailtyMM function is used to calculate the gamma frailty model.

Usage

GaFrailtyMM(
  formula,
  data,
  beta = NULL,
  theta = NULL,
  lambda = NULL,
  Maxiter = 2000,
  convergence = 1e-06,
  ...
)

Arguments

Details

The GaFrailtyMM function is used to calculate gamma frailty survival model using MM algorithms based on AD technology. EM algorithms relies on the fact that, after profiling out the nonparametric component \Lambda_0, the resulting function is concave. When it does not hold, using Newton method to maximize the resulting function is difficult especially when there exist a large number of covariates. MM algorithms that can avoid the concavity requirement and bypass Newton method and matrix inversion.

Value

An object of class GaFrailtyMM that contains the following fields: total amount of observations, the Time, the \Lambda, the \lambda, total number of failure events, total number of iterations, convergence result, the log likelihood value, the \theta, the standard deviation of the estimated \theta, the likelihood-based 95% confidence interval for the \theta, \beta, the standard deviation of the estimated \beta, the likelihood-based 95% confidence interval for the \beta, the variable name.

References

Huang X.F., Xu J.F. and Tian G.L.(2019). 'On profile MM algorithms for gamma frailty survival models.' Statistica Sinica 29(2), 895-916.

Examples

library(survival)
GaFrailtyMM(Surv(time, status) ~ age + sex + cluster(id), data=kidney)

Control IC2Pro object

Description

Control IC2Pro object

Usage

IC2Control(Maxiter = 2000, convergence = 1e-06, Idigits = 4, Pdigits = 4)

Arguments

Value

list of Maxiter, convergence, Idigits, Pdigits.

Examples

IC2Control()

MM algorithm based on the AD method for case II interval-censored data

Description

The IC2MM function is used to calculate the case II interval-censored data model. A failure time study that consists of n independent subjects from a homogeneous population with survival function S_{(t)}. Let T_{i} denote the survival time, and i=1, \ldots, n. Suppose that interval-censored data on the T_i are observed and given by

Y_{obs} = \{ (L_{i}, R_{i}]; i=1, \ldots, n \}

where T_i \in (L_{i}, R_{i}] . Let \{s_i \}_{j=0}^{m} denote the unique ordered elements of {0, L_{i}, R_{i}, i=1, \ldots, n } . Take \alpha_{ij} = I(s_{j} \in (L_{i}, R_{i}] ) and p_{j}= S(s_{j-1}) - S(s_{j}), j= 1, \ldots, m . The log-likelihood function is

\ell( {p} | Y_{obs}) = \sum_{i=1}^{n} \log (S(L_{i}) - S(R_{i}) ) = \sum_{i=1}^{n} \log \left( \sum_{j=1}^{m} \alpha_{ij} p_{j} \right)

where {p} = (p_1, \ldots, p_m)^{T} and \sum_{j=1}^{m} p_{j} = 1 , p_{j} \geqslant 0.

Usage

IC2MM(formula, data, ...)

Arguments

Details

The IC2MM function allows the distributions for multiple strata of dataset to be stored as one IC2 object, e.g. data=bcos.

Value

An object of class IC2MM that contains the following fields: error: convergence result; strata: dimensions of each df_tab; s: unique ordered elements of {0, L_{i}, R_{i}, Inf} , if more than one strata, elements are concatenated; S: the survival function, if more than one strata, values are concatenated; df_tab: the dataframe of survival intervals and survival probabilities for each interval, if more than one strata, dataframes are concatenated.

References

Tian G.L., Huang X.F. and Xu, J.(2019). 'An assembly and decomposition approach for constructing separable minorizing functions in a class of MM algorithms.' Statistica Sinica 29(2), 961-982.

See Also

IC2Pro

Examples

library(survival)
L <- c(1.5, 0.1, 1.5, 0.5, 0.4, 0.2, 0.9, 0.2, 0.08, 1.9)
R <- c(2.1, 2.9, 2.7, 1.9, 1.3, 1.4, 2.3, 0.5, 1.5, 4.6 )
data <- data.frame(L, R)
IC2MM(Surv(L,R, type = 'interval2') ~ 1, data )

IC2MM(Surv(L,R, type = 'interval2') ~ 1, data, control=IC2Control(Pdigits=2) )

Calculate non-parametric estimate for case II interval censored survival function

Description

Calculate non-parametric estimate for case II interval censored survival function

Usage

IC2Pro(L, R, control = IC2Control(), ...)

Arguments

Value

An object of class IC2Pro that contains the following fields: error: convergence result; strata: dimensions of df_tab; s: unique ordered elements of {0, L_{i}, R_{i}, Inf} ; S: the survival function; df_tab: the data frame of survival intervals and survival probabilities for each interval.

References

Tian G.L., Huang X.F. and Xu, J.(2019). 'An assembly and decomposition approach for constructing separable minorizing functions in a class of MM algorithms.' Statistica Sinica 29(2), 961-982.

See Also

IC2Control

Examples

L <- c(1.4, 1.5, 1.3, 0.9, 0.4, 0.2, 0.5, 0.03, 1.7, 0.2)
R <- c(2.2, 3, 2.4, 1.2, 2.8, 0.3, 1.6, 2.5, 2.6, 3.4)
IC2Pro(L, R, control=IC2Control())

MM algorithm based on the AD method for left-truncated normal distribution

Description

The LTNMM function is used to calculate a left-truncated normal distribution model. A LTN(\mu, \sigma^2; a) has the density function

f(y; \mu, \sigma^2; a) = \frac{1}{c \sqrt{2 \pi \sigma^{2}}} \exp{( -\frac{(y-\mu)^{2}}{2 \sigma^{2}} )} \centerdot I(y \geq a)

where (\mu, \sigma^2) are two unknown parameters, a is a known constant, c = 1- \Phi(\frac{a-u}{ \sigma}) , and \Phi(\centerdot) is the cdf of the standard normal distribution.

Usage

LTNMM(
  formula,
  a,
  mu = NULL,
  sigma = NULL,
  data = sys.frame(sys.parent()),
  Maxiter = 2000,
  convergence = 1e-06,
  ...
)

Arguments

Details

The LTNMM function is used to calculate a left-truncated normal distribution model using MM algorithms based on AD technology. The formula parameter can be used to provide the data that needs to be calculated, such as formula=y~1. By default, the data is provided by the user’s environment. The initial values of the mean and variance of the normal distribution are estimated using OLS.

Value

An object of class LTNMM that contains the following fields: total amount of observations, the number of iterations, convergence rate, the log likelihood value, estimated results for the unknown parameters, the standard deviation of estimate for the unknown parameters, the likelihood-based 95% confidence interval for the unknown parameters, information criterion: AIC value and BIC value.

References

Tian G.L., Huang X.F., and Xu, J.(2019). 'An assembly and decomposition approach for constructing separable minorizing functions in a class of MM algorithms.' Statistica Sinica 29(2), 961-982.

Examples

y=c(8.7, 5.4, 8.9, 5.8, 6.2, 9.9, 7.5, 9.5, 6.5, 6.3); a=5
LTNMM(y~1, a=5)

MM algorithm based on the AD method for type I multivariate zero-inflated generalized poisson distribution

Description

Let Z \sim Bernoulli(1- \phi), {x} = (X_1,\cdots, X_m)^T, X_i \sim GP(\lambda_i, \theta_i), for i = 1,\cdots,m , and (Z,X_1,\cdots, X_m) are mutually independent. An m dimensional discrete random vector {y}=(Y_1,\cdots, Y_m)^T is said to have a Type I multivariate zero-inflated generalized Poisson distribution(ZIGP) distribution if

{y} \overset{\rm{d}}= Z{x}= \left\{ \begin{array}{ll} {0} & \text{with probability} \ \lambda \\ {x} & \text{with probability} \ 1-\lambda \\ \end{array} \right.

where \lambda \in [0,1), {\lambda}=(\lambda_1,\cdots, \lambda_m)^T , {\theta}=(\theta_1,\cdots, \theta_m)^T , max(-1, -\lambda_i/q_i) < \theta_i \leqslant 1 and q_i \geqslant 4 is the largest positive integer for each \lambda_i + \theta_{i}q_{i}>0 when \theta_{i}<0. The ZIGPMM function is used to calculate the Type I multivariate ZIGP model.

Usage

ZIGPMM(data, phi0, la, th, Maxiter = 2000, convergence = 1e-06, ...)

Arguments

Details

The ZIGPMM function is used to calculate Type I multivariate zero-inflated generalized Poisson distribution model using MM algorithms based on AD technology. data is provided by user by default, it can be a data frame or a matrix. In addition, the unknown parameters require users to give appropriate initial values, where 0<=phi0<1, each th should satisfy 0<=th<1, and each la should be greater than 0.

Value

An object of class ZIGPMM that contains the following fields: total amount of observations, the number of iterations, convergence rate, the log likelihood value, estimated results for the unknown parameters, the standard deviation of estimate for the unknown parameters, the likelihood-based 95% confidence interval for the unknown parameters, information criterion: AIC value and BIC value.

References

Tian G.L., Huang X.F. and Xu, J.(2019). 'An assembly and decomposition approach for constructing separable minorizing functions in a class of MM algorithms.' Statistica Sinica 29(2), 961-982.

Huang X.F., Tian G.L., Zhang, C. and Jiang, X.(2017). 'Type I multivariate zero-inflated generalized Poisson distribution with applications.' Statistics and its Interface 10(2), 291-311.

Examples

x1 <- c(0, 0, 0,38, 0,19,25, 0,25, 0)
x2 <- c(0, 0, 0,23, 0,51,24, 0,10, 0)
y <- cbind(x1, x2)
phi0 = 0.5; la = rep(1,2); th = rep(0.1,2)
ZIGPMM(y, phi0, la, th)

Breast Cosmesis Data

Description

The often used data set for interval censored data, described and given in full in Finkelstein and Wolfe (1985).

Usage

bcos

Format

An object of class data.frame with 94 rows and 3 columns.

Value

Breast cosmesis data contains the following fields:

References

Finkelstein D.M. and Wolfe R.A.(1985). "A semiparametric model for regression analysis of interval-censored failure time data." Biometrics 41, 933-945.

Examples

data = data(bcos)

The children’s absenteeism data in Indonesia

Description

In a survey of Indonesian family life conducted by Strauss et al. the participants included 7000 households sampled from 321 communities randomly selected from 13 of the nation’s 26 Provinces, in which 83% of the Indonesian population lived. Among those households with one child per household, 437 household heads were asked questions about the health of their children.

Usage

cadi

Format

An object of class data.frame with 437 rows and 2 columns.

Value

The children’s absenteeism data in Indonesia contains the following fields:

References

Huang X.F., Tian G.L., Zhang, C. and Jiang, X.(2017). "Type I multivariate zero-inflated generalized Poisson distribution with applications." Statistics and its Interface 10(2), 291-311.

Strauss J., Beegle K., Sikoki B., Dawiyanto A., Herawati Y. and Witoelar Y.(2004). "The Third Wave of the Indonesia Family Life Survey (IFLS): Overview and Field Report, WR-144/1-NIA/NICHD, RAND Corporation, Santa Monica, CA."

Examples

data = data(cadi)

Kidney Infection Data

Description

The data consisted of the time to first and second infection relapse in 38 kidney disease patients using a portable dialysis machine. Infection may occur where the catheter was inserted. Catheters are subsequently removed if infection develops and may be removed for other reasons, in which case observations are censored.

Usage

kidney

Format

An object of class data.frame with 76 rows and 7 columns.

Value

Kidney infection data contains the following fields:

References

McGilchrist C.A. and Aisbett C.W.(1991). "Regression with frailty in survival analysis." Biometrics 47, 461-466.

Examples

data = data(Kidney)

NCCTG Lung Cancer Data

Description

Survival in patients with advanced lung cancer from the North Central Cancer Treatment Group. Performance scores rate how well the patient can perform usual daily activities.

Usage

lung

Format

An object of class data.frame with 228 rows and 10 columns.

Value

Kidney infection data contains the following fields:

References

Finkelstein D.M. and Wolfe R.A.(1985). "A semiparametric model for regression analysis of interval-censored failure time data." Biometrics 41, 933-945.

Examples

data = data(lung)

Plot the Cox object

Description

Plot the Cox object

Usage

## S3 method for class 'Cox'
plot(
  x,
  xlab = "Time",
  ylab = "Cumulative hazard",
  type = "s",
  lty = 1,
  lwd = 1,
  col = gray(0),
  digits = 4,
  ...
)

Arguments

Value

the dataframe of 'Time' and accumulative hazard \Lambda.

Examples

library(survival)
result <- CoxMM(Surv(time, status) ~ age + sex, lung)

plot(result)

Plot the GaF object

Description

Plot the GaF object

Usage

## S3 method for class 'GaF'
plot(
  x,
  xlab = "Time",
  ylab = "Cumulative hazard",
  type = "s",
  lty = 1,
  lwd = 1,
  col = gray(0),
  digits = 4,
  ...
)

Arguments

Value

the dataframe of 'Time' and accumulative hazard \Lambda.

Examples

library(survival)
result <- GaFrailtyMM(Surv(time, status) ~ age + sex + cluster(id), data=kidney)
plot(result)

Plot the IC2 object

Description

Plot the IC2 object

Usage

## S3 method for class 'IC2'
plot(
  x,
  xlab = "Time",
  ylab = "Survival",
  legend = NULL,
  main = "Survival Function",
  lty = 1:9,
  lwd = 1,
  xleg = 0,
  yleg = 0.15,
  col = gray(0),
  ...
)

Arguments

Value

A list of arguments for the legend. Values are x, y, legend, fill, lty, bty, col.

Examples

library(survival)
result = IC2MM(Surv(left, right, type = 'interval2') ~ treatment, bcos)

plot(result, col=c('red', 'blue'))

Summary of parameter estimates of a multivariate compound ZIGP model

Description

This function returns the result of the CZIGPMM function

Usage

## S3 method for class 'CZIGP'
summary(object, digits = 4, ...)

Arguments

Value

Summary for CZIGPMM objects.

See Also

CZIGPMM

Examples

x1 <- c(0,35,23,34,8,19,0,0,0,0)
x2 <- c(38,15,0,25,34,0,0,0,0,0)
y <- cbind(x1, x2)
phi0 = 0.5; phi = rep(0.5,2); la = rep(1,2); th = rep(0.1,2)
result <- CZIGPMM(y, phi0, phi, la, th)

summary(result,digits=4)

Summary of parameter estimates of a Cox model

Description

This function returns the result of the CoxMM function

Usage

## S3 method for class 'Cox'
summary(object, digits = 4, ...)

Arguments

Value

Summary for CoxMM objects.

See Also

CoxMM

Examples

library(survival)
result <- CoxMM(Surv(time, status) ~ age + sex, lung)

summary(result,digits=4)

Summary of parameter estimates of a gamma frailty model

Description

This function returns the result of the GaFrailtyMM function

Usage

## S3 method for class 'GaF'
summary(object, digits = 4, ...)

Arguments

Value

Summary for GaFrailtyMM objects.

See Also

GaFrailtyMM

Examples

library(survival)
result <- GaFrailtyMM(Surv(time, status) ~ age + sex + cluster(id), data=kidney)

summary(result,digits=4)

Summary of parameter estimates of a IC2 model

Description

This function returns the result of the IC2MM function

Usage

## S3 method for class 'IC2'
summary(object, ...)

Arguments

Value

Summary for IC2MM objects.

See Also

IC2MM

Examples

library(survival)
result <- IC2MM(Surv(left, right, type = 'interval2') ~ treatment, bcos)

summary(result)

Summary of parameter estimates of a LTN model

Description

This function returns the result of the LTNMM function

Usage

## S3 method for class 'LTN'
summary(object, digits = 4, ...)

Arguments

Value

Summary for LTNMM objects.

See Also

LTNMM

Examples

y=c(8.7, 5.4, 8.9, 5.8, 6.2, 9.9, 7.5, 9.5, 6.5, 6.3); a=5
result <- LTNMM(y~1, a=5)

summary(result,digits=4)

Summary of parameter estimates of a Type I multivariate ZIGP model

Description

This function returns the result of the ZIGPMM function

Usage

## S3 method for class 'ZIGP'
summary(object, digits = 4, ...)

Arguments

Value

Summary for ZIGPMM objects.

See Also

ZIGPMM

Examples

x1 <- c(0, 0, 0,38, 0,19,25, 0,25, 0)
x2 <- c(0, 0, 0,23, 0,51,24, 0,10, 0)
y <- cbind(x1, x2)
phi0 = 0.5; la = rep(1,2); th = rep(0.1,2)
result <- ZIGPMM(y, phi0, la, th)

summary(result,digits=4)

Voluntary and involuntary job changes data

Description

Jung and Winkelmann(1993) provided data on both the numbers of voluntary and involuntary job changes of males during ten period 1974–1984. The samples contain 2124 males who started their working career before or in 1974 and did not retire before 1984.

Usage

vijc

Format

An object of class data.frame with 2124 rows and 2 columns.

Value

Voluntary and involuntary job changes data contains the following fields:

References

Huang X.F., Tian G.L., Zhang, C. and Jiang, X.(2017). "Type I multivariate zero-inflated generalized Poisson distribution with applications." Statistics and its Interface 10(2), 291-311.

Jung R.C. and Winkelmann R.(1993). "Two aspects of labor mobility: A bivariate Poisson regression approach." Empirical Economics 18(3), 543–556.

Examples

data = data(vijc)