Title: | Model Selection in Multivariate Longitudinal Data Analysis |
---|---|
Description: | An efficient Gibbs sampling algorithm is developed for Bayesian multivariate longitudinal data analysis with the focus on selection of important elements in the generalized autoregressive matrix. It provides posterior samples and estimates of parameters. In addition, estimates of several information criteria such as Akaike information criterion (AIC), Bayesian information criterion (BIC), deviance information criterion (DIC) and prediction accuracy such as the marginal predictive likelihood (MPL) and the mean squared prediction error (MSPE) are provided for model selection. |
Authors: | Kuo-Jung Lee |
Maintainer: | Kuo-Jung Lee <[email protected]> |
License: | GPL-2 |
Version: | 1.0 |
Built: | 2024-11-08 06:32:57 UTC |
Source: | CRAN |
Using MCMC procedure to generate posterior samples and provide AIC, BIC, DIC, MPL, MSPE, and predicted values.
MLModelSelectionMCMC(Num.of.iterations, list.Data, list.InitialValues, list.HyperPara, list.UpdatePara, list.TuningPara)
MLModelSelectionMCMC(Num.of.iterations, list.Data, list.InitialValues, list.HyperPara, list.UpdatePara, list.TuningPara)
Num.of.iterations |
Number of iterations. |
list.Data |
List of data set containing response |
list.InitialValues |
List of initial values for parameters. |
list.HyperPara |
List of given hyperparameters in priors. |
list.UpdatePara |
Determine which parameter will be updated. |
list.TuningPara |
Provide turning parameters in proposal distributions. |
We set the subject (
) has
continuous responses at each time point
(
). Assume that the measurement times are common across subjects, but not necessarily equally-spaced. Let
denote the response vector containing
continuous responses for
th subject at time
along with a
vectof of covariates,
. An efficient Gibbs sampling algorithm is developed for model estimation in the multivariate longitudinal model given by
where is a vector of regression coefficients of length
,
is a generalized autoregressive parameter (GARP) to explain the serial dependence of responses across time. Moreover,
The priors for the parameters in the model given by
where ,
, and
are prespecified values. For
and
, we further assume
where is prespecified value and
is the point mass at 0.
Lists of posterior samples, parameters estimates, AIC, BIC, DIC, MPL, MSPE, and predicted values are returned
We'll provide the reference for details of the model and the algorithm for performing model estimation whenever the manuscript is accepted.
Kuo-Jung Lee
Keunbaik Lee et al. (2015) Estimation of covariance matrix of multivariate longitudinal data using modified Choleksky and hypersphere decompositions. Biometrics. 75-86, 2020. doi:10.1111/biom.13113.
library(MASS) library(MLModelSelection) AR.Order = 6 #denote \phi_{itj, kg} = \alpha_{kg} \mathbf{1}{|t-j|=1} ISD.Model = 1 #denote \log(\sigma_{itk}) = \lambda_{k0} + \lambda_{k1} h_{it} data(SimulatedData) N = dim(SimulatedData$Y)[1] # the number of subjects T = dim(SimulatedData$Y)[2] # time points K = dim(SimulatedData$Y)[3] # the number of attributes P = dim(SimulatedData$X)[3] # the number of covariates M = AR.Order # the demension of alpha nlamb = ISD.Model + 1 # the dimension of lambda Data = list(Y = SimulatedData$Y, X = SimulatedData$X, TimePointsAvailable = SimulatedData$TimePointsAvailable, AR.Order = AR.Order, ISD.Model = ISD.Model) beta.ini = matrix(rnorm(P*K), P, K) delta.ini = array(rbinom(K*K*M, 1, 0.1), c(K, K, M)) alpha.ini = array(runif(K*K*M, -1, 1), c(K, K, M)) lambda.ini = matrix(rnorm(nlamb*K), K, nlamb, byrow=T) nu.ini = rnorm(K) InitialValues = list(beta = beta.ini, delta = delta.ini, alpha = alpha.ini, lambda = lambda.ini, nu = nu.ini) # Hyperparameters in priors sigma2.beta = 1 sigma2.alpha = 10 sigma2.lambda = 0.01 sigma2.nu = 0.01 # Whehter the parameter will be updated UpdateBeta = TRUE UpdateDelta = TRUE UpdateAlpha = TRUE UpdateLambda = TRUE UpdateNu = TRUE HyperPara = list(sigma2.beta = sigma2.beta, sigma2.alpha=sigma2.alpha, sigma2.lambda=sigma2.lambda, sigma2.nu=sigma2.nu) UpdatePara = list(UpdateBeta = UpdateBeta, UpdateAlpha = UpdateAlpha, UpdateDelta = UpdateDelta, UpdateLambda = UpdateLambda, UpdateNu = UpdateNu) # Tuning parameters in proposal distribution within MCMC TuningPara = list(TuningAlpha = 0.01, TuningLambda = 0.005, TuningNu = 0.005) num.of.iter = 100 start.time <- Sys.time() PosteriorSamplesEstimation = MLModelSelectionMCMC(num.of.iter, Data, InitialValues, HyperPara, UpdatePara, TuningPara) end.time <- Sys.time() cat("Estimate of beta\n") print(PosteriorSamplesEstimation$PosteriorEstimates$beta.mean)
library(MASS) library(MLModelSelection) AR.Order = 6 #denote \phi_{itj, kg} = \alpha_{kg} \mathbf{1}{|t-j|=1} ISD.Model = 1 #denote \log(\sigma_{itk}) = \lambda_{k0} + \lambda_{k1} h_{it} data(SimulatedData) N = dim(SimulatedData$Y)[1] # the number of subjects T = dim(SimulatedData$Y)[2] # time points K = dim(SimulatedData$Y)[3] # the number of attributes P = dim(SimulatedData$X)[3] # the number of covariates M = AR.Order # the demension of alpha nlamb = ISD.Model + 1 # the dimension of lambda Data = list(Y = SimulatedData$Y, X = SimulatedData$X, TimePointsAvailable = SimulatedData$TimePointsAvailable, AR.Order = AR.Order, ISD.Model = ISD.Model) beta.ini = matrix(rnorm(P*K), P, K) delta.ini = array(rbinom(K*K*M, 1, 0.1), c(K, K, M)) alpha.ini = array(runif(K*K*M, -1, 1), c(K, K, M)) lambda.ini = matrix(rnorm(nlamb*K), K, nlamb, byrow=T) nu.ini = rnorm(K) InitialValues = list(beta = beta.ini, delta = delta.ini, alpha = alpha.ini, lambda = lambda.ini, nu = nu.ini) # Hyperparameters in priors sigma2.beta = 1 sigma2.alpha = 10 sigma2.lambda = 0.01 sigma2.nu = 0.01 # Whehter the parameter will be updated UpdateBeta = TRUE UpdateDelta = TRUE UpdateAlpha = TRUE UpdateLambda = TRUE UpdateNu = TRUE HyperPara = list(sigma2.beta = sigma2.beta, sigma2.alpha=sigma2.alpha, sigma2.lambda=sigma2.lambda, sigma2.nu=sigma2.nu) UpdatePara = list(UpdateBeta = UpdateBeta, UpdateAlpha = UpdateAlpha, UpdateDelta = UpdateDelta, UpdateLambda = UpdateLambda, UpdateNu = UpdateNu) # Tuning parameters in proposal distribution within MCMC TuningPara = list(TuningAlpha = 0.01, TuningLambda = 0.005, TuningNu = 0.005) num.of.iter = 100 start.time <- Sys.time() PosteriorSamplesEstimation = MLModelSelectionMCMC(num.of.iter, Data, InitialValues, HyperPara, UpdatePara, TuningPara) end.time <- Sys.time() cat("Estimate of beta\n") print(PosteriorSamplesEstimation$PosteriorEstimates$beta.mean)
A simulated multivariate longitudinal data for demonstration.
data("SimulatedData")
data("SimulatedData")
A list
consists of Y
the observations 100 subjects in 3 attributes along 10 time points,
X
the design matrix with 6 covariate including the intercept, TimePointsAvailable
the avilable time points for each subject.
Y
The response variables.
X
The design matrix.
TimePointsAvailable
The available time points for each subject.
library(MLModelSelection) data(SimulatedData) SimulatedData = data(SimulatedData)
library(MLModelSelection) data(SimulatedData) SimulatedData = data(SimulatedData)