Title: | Scale Mixture of Normal Distribution in Linear Mixed-Effects Model |
---|---|
Description: | Bayesian analysis of censored linear mixed-effects models that replace Gaussian assumptions with a flexible class of distributions, such as the scale mixture of normal family distributions, considering a damped exponential correlation structure which was employed to account for within-subject autocorrelation among irregularly observed measures. For more details, see Zhong et al. (2025, forthcoming in Statistics in Medicine). |
Authors: | Kelin Zhong [aut, cre], Fernanda L. Schumacher [aut], Luis M. Castro [aut], Victor H. Lachos [aut], Jalmar M.F. Carrasco [aut] |
Maintainer: | Kelin Zhong <[email protected]> |
License: | GPL-3 |
Version: | 0.1.0 |
Built: | 2024-11-27 07:41:01 UTC |
Source: | https://github.com/kelinzhong/smnlmec |
This function fits left, right censored mixed-effects linear model, with scale mixture of normal distribution errors, using the Stan. It returns estimates, standard errors and LPML, AIC, BIC and DIC.
SMNlmec.est( ID, x_set, z_set, tt, y_complete, censor_vector, dist = "Normal", struc = "UNC", direction = "left", thin_num = 1, chains_num = 1, iter_num = 3000, burn_percen = 0.1, seed_set = NULL, adapt_delta_set = 0.8 )
SMNlmec.est( ID, x_set, z_set, tt, y_complete, censor_vector, dist = "Normal", struc = "UNC", direction = "left", thin_num = 1, chains_num = 1, iter_num = 3000, burn_percen = 0.1, seed_set = NULL, adapt_delta_set = 0.8 )
ID |
Vector |
x_set |
Design matrix of the fixed effects of order |
z_set |
Design matrix of the random effects of order |
tt |
Vector |
y_complete |
Vector |
censor_vector |
Vector |
dist |
Distribution of the random effects and random error. Available options are |
struc |
Structure of the correlation structure. Available options are |
direction |
Direction of censoring type. Available options are |
thin_num |
A positive integer specifying the period for saving samples. The default is 5. See more details in rstan::stan(). |
chains_num |
A positive integer specifying the number of chains generating by rstan::stan(). The default is 3. |
iter_num |
A positive integer specifying the number of iterations for each chain (including warmup). The default is 5000. |
burn_percen |
A percentage of the warm-up iterations in each chain the Stan. The default is 0.2. |
seed_set |
A random seed. The default is NULL. |
adapt_delta_set |
A parameter to control the sampler's behavior. The default is 0.8. See rstan::stan() for more details. |
Return a S4 class SMNlmecfit object. Using function SMNlmec.summary()
to obtain the estimation of parameters and model selection criteria. The SMNlmecfit include:
stan_object |
A stanfit object from rstan::stan(). |
model_criteria |
A list includes LPML, DIC, EAIC, EBIC, K-L divergence. |
dist_set |
The setting of distribution of the stan model. |
struc_set |
The setting of correlation structure of the stan model. |
Kelin Zhong, Fernanda L. Schumacher, Luis M. Castro and Victor H. Lachos. Bayesian analysis of censored linear mixed-effects models for heavy-tailed irregularly observed repeated measures. Statistics in Medicine, 2025. doi:10.1002/sim.10295
require(rstan) require(StanHeaders) require(MASS) require(tmvtnorm) require(mvtnorm) require(mnormt) data("UTIdata_sub") data1 <- UTIdata_sub y1 <- c(log10(data1$RNA)) cc <- (data1$RNAcens==1)+0 y_com<-as.numeric(y1) rho_com<-as.numeric(cc) x <- cbind( (data1$Fup==0)+0, (data1$Fup==1)+0, (data1$Fup==3)+0, (data1$Fup==6)+0, (data1$Fup==9)+0, (data1$Fup==12)+0, (data1$Fup==18)+0, (data1$Fup==24)+0 ) z <- matrix(rep(1, length(y1)), ncol=1) UTI_T_DEC <- SMNlmec.est(ID = data1$Patid, x_set = x, z_set = z, tt = data1$Fup, y_complete = y_com, censor_vector = rho_com, dist = "Student", struc = "DEC", direction = "left", thin_num = 1, chains_num = 1, iter_num = 3000, burn_percen = 0.1, seed_set = 9955, adapt_delta_set = 0.8) SMNlmec.summary(UTI_T_DEC)
require(rstan) require(StanHeaders) require(MASS) require(tmvtnorm) require(mvtnorm) require(mnormt) data("UTIdata_sub") data1 <- UTIdata_sub y1 <- c(log10(data1$RNA)) cc <- (data1$RNAcens==1)+0 y_com<-as.numeric(y1) rho_com<-as.numeric(cc) x <- cbind( (data1$Fup==0)+0, (data1$Fup==1)+0, (data1$Fup==3)+0, (data1$Fup==6)+0, (data1$Fup==9)+0, (data1$Fup==12)+0, (data1$Fup==18)+0, (data1$Fup==24)+0 ) z <- matrix(rep(1, length(y1)), ncol=1) UTI_T_DEC <- SMNlmec.est(ID = data1$Patid, x_set = x, z_set = z, tt = data1$Fup, y_complete = y_com, censor_vector = rho_com, dist = "Student", struc = "DEC", direction = "left", thin_num = 1, chains_num = 1, iter_num = 3000, burn_percen = 0.1, seed_set = 9955, adapt_delta_set = 0.8) SMNlmec.summary(UTI_T_DEC)
Generating Censored UNC, DEC, CAR errors with Mixed Effects, for normal, student's-t and slash distribution.
SMNlmec.sim( m, x, z, tt, nj, beta, sigma2, D, phi, struc = "UNC", typeModel = "Normal", p.cens = 0.1, n.cens = NULL, cens_type = "right", nu_set = NULL )
SMNlmec.sim( m, x, z, tt, nj, beta, sigma2, D, phi, struc = "UNC", typeModel = "Normal", p.cens = 0.1, n.cens = NULL, cens_type = "right", nu_set = NULL )
m |
Number of individuals. |
x |
Design matrix of the fixed effects of order |
z |
Design matrix of the random effects of order |
tt |
Vector |
nj |
Vector |
beta |
Vector of values fixed effects. |
sigma2 |
Values of the scalar of the variance matrix. |
D |
Variance matrix of the random effects of order |
phi |
Vector of parameter in the |
struc |
Structure for the simulated data. Available options are |
typeModel |
Distribution of the simulated data. Available options are |
p.cens |
Percentage of censored measurements in the responses. The default value is 0.1. |
n.cens |
Number of censored measurements in the responses. The default value is NULL. |
cens_type |
The direction of cesoring. Available options are |
nu_set |
degrees of freedom for student's-t or slash simulated data. The default value is NULL. |
return list:
cc |
Vector of censoring indicators. |
y_cc |
Vector of responses censoring. |
p.cens <- 0.1 m <- 50 D <- matrix(c(0.049,0.001,0.001,0.002),2,2) sigma2_set <- 0.15 beta <- c(-2.83,-0.18) nu <- 2 phi <- c(0.6,2) nj <- rep(6,m) tt <- rep(1:6,length(nj)) X1 <- rep(1,sum(nj)) X2 <- tt x <- as.matrix(cbind(X1,X2)) Z1 <- rep(1,sum(nj)) Z2 <- tt z <- as.matrix(cbind(Z1,Z2)) ID_sim <- rep(0,length(tt)) ID_log <- 0 for(i in 1:m) { for(j in 1:nj[i]) { ID_sim[ID_log + j] <- i } ID_log <- ID_log + nj[i] } Slash_DEC_sim <- SMNlmec.sim(m = m,x = x,z = z,tt = tt,nj = nj,beta = beta, sigma2 = sigma2_set,D = D,phi= phi,struc ="DEC", typeModel="Slash",p.cens = p.cens,n.cens = NULL, cens_type="right",nu_set=nu) head(Slash_DEC_sim$cc) sum(Slash_DEC_sim$cc)/length(Slash_DEC_sim$cc) head(Slash_DEC_sim$y_cc) y_com <- as.numeric(Slash_DEC_sim$y_cc) rho_com <- as.numeric(Slash_DEC_sim$cc) tem <- tt Slash_DEC_est <- SMNlmec.est(ID = ID_sim, x_set = x, z_set = z, tt = tem, y_complete = y_com, censor_vector = rho_com, dist = "Slash", struc = "DEC", direction = "right", thin_num = 1, chains_num = 1, iter_num = 3000, burn_percen = 0.1, seed_set = 9955, adapt_delta_set = 0.8) SMNlmec.summary(Slash_DEC_est)
p.cens <- 0.1 m <- 50 D <- matrix(c(0.049,0.001,0.001,0.002),2,2) sigma2_set <- 0.15 beta <- c(-2.83,-0.18) nu <- 2 phi <- c(0.6,2) nj <- rep(6,m) tt <- rep(1:6,length(nj)) X1 <- rep(1,sum(nj)) X2 <- tt x <- as.matrix(cbind(X1,X2)) Z1 <- rep(1,sum(nj)) Z2 <- tt z <- as.matrix(cbind(Z1,Z2)) ID_sim <- rep(0,length(tt)) ID_log <- 0 for(i in 1:m) { for(j in 1:nj[i]) { ID_sim[ID_log + j] <- i } ID_log <- ID_log + nj[i] } Slash_DEC_sim <- SMNlmec.sim(m = m,x = x,z = z,tt = tt,nj = nj,beta = beta, sigma2 = sigma2_set,D = D,phi= phi,struc ="DEC", typeModel="Slash",p.cens = p.cens,n.cens = NULL, cens_type="right",nu_set=nu) head(Slash_DEC_sim$cc) sum(Slash_DEC_sim$cc)/length(Slash_DEC_sim$cc) head(Slash_DEC_sim$y_cc) y_com <- as.numeric(Slash_DEC_sim$y_cc) rho_com <- as.numeric(Slash_DEC_sim$cc) tem <- tt Slash_DEC_est <- SMNlmec.est(ID = ID_sim, x_set = x, z_set = z, tt = tem, y_complete = y_com, censor_vector = rho_com, dist = "Slash", struc = "DEC", direction = "right", thin_num = 1, chains_num = 1, iter_num = 3000, burn_percen = 0.1, seed_set = 9955, adapt_delta_set = 0.8) SMNlmec.summary(Slash_DEC_est)
A generic function to provide a summary for objects of class SMNlmecfit
.
SMNlmec.summary(object) ## S4 method for signature 'SMNlmecfit' SMNlmec.summary(object)
SMNlmec.summary(object) ## S4 method for signature 'SMNlmecfit' SMNlmec.summary(object)
object |
An object of class |
A summary of model estimations, R-hats, standard errors, and criteria.
A printed summary of the SMNlmecfit object.
Get a summary of results from SMNlmec.est.
## S4 method for signature 'SMNlmecfit' SMNlmec.summary(object)
## S4 method for signature 'SMNlmecfit' SMNlmec.summary(object)
object |
An object of class |
The summary of estimations, R hats, standard errors and 95
SMNlmecfit Class
stan_object
stanfit object from rstan.
model_criteria
list, model selection criteria.
dist_set
character, the name of distribution.
struc_set
character, the name of correlation structure.
A function to create objects of class SMNlmecfit
.
SMNlmecfit.creator(stan_object, model_criteria, dist_set, struc_set)
SMNlmecfit.creator(stan_object, model_criteria, dist_set, struc_set)
stan_object |
stanfit object from rstan. |
model_criteria |
list, model selection criteria. |
dist_set |
character, the name of distribution. |
struc_set |
character, the name of correlation structure. |
A SMNlmecfit object containing:
stan_object |
A stanfit object from rstan::stan(). |
model_criteria |
A list includes LPML, DIC, EAIC, EBIC, K-L divergence. |
dist_set |
The setting of distribution of the stan model. |
struc_set |
The setting of correlation structure of the stan model. |
Data set from a study of Unstructured Treatment Interruption in HIV-infected adolescents in four institutions in the US. The main outcome is the HIV-1 RNA viral load, which is subject to censoring below the lower limit of detection of the assay (50 copies/mL). The censored observations are indicated by the variable RNAcens.
data(UTIdata)
data(UTIdata)
A data frame with 373 observations on the following 5 variables.
patient ID
days after treatment interruption.
follow-up months
viral load RNA
censoring indicator for viral load
Saitoh, A., Foca, M, et al. (2008), Clinical outcome in perinatally acquired HIV-infected children and adolescents after unstructured treatment interruption, Pediatrics,121, e513-e521.
Data set from a study of Unstructured Treatment Interruption in HIV-infected adolescents in four institutions in the US. The main outcome is the HIV-1 RNA viral load, which is subject to censoring below the lower limit of detection of the assay (50 copies/mL). The censored observations are indicated by the variable RNAcens. Excluding subjects whose observations are less than 2 and with missing RNA (excluding subject ID C6 T16).
data(UTIdata_sub)
data(UTIdata_sub)
A data frame with 360 observations on the following 5 variables.
patient ID
days after treatment interruption.
follow-up months
viral load RNA
censoring indicator for viral load
Saitoh, A., Foca, M, et al. (2008), Clinical outcome in perinatally acquired HIV-infected children and adolescents after unstructured treatment interruption, Pediatrics,121, e513-e521.