Description Usage Arguments Details Value Note Author(s) References Examples
Using MCMC procedure to generate posterior samples and provide AIC, BIC, DIC, MPL, MSPE, and predicted values.
1 2 | 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 Y, design matrix X, avialable time points for each subject, GARP model, and ISD model. |
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 i (i=1, …, N) has K continuous responses at each time point t (t=1, …, n_i). Assume that the measurement times are common across subjects, but not necessarily equally-spaced. Let {y}_{it} = (y_{it1}, …, y_{itK}) denote the response vector containing K continuous responses for ith subject at time t along with a p\times 1 vectof of covariates, {x}_{it} = (x_{it1}, …, x_{itp}). An efficient Gibbs sampling algorithm is developed for model estimation in the multivariate longitudinal model given by
y_{i1k} = {x}'_{it}{β}_k + e_{i1k}, t=1;
y_{itk} = {x}'_{it}{β}_k + ∑_{g=1}^K∑_{j=1}^{t-1} φ_{itj, kg} (y_{ijg}-x'_{ij}{β}_g)+ e_{itk}, t≥q 2,
where {β}_k = (β_{k1}, …, β_{kp})' is a vector of regression coefficients of length p, φ_{itj, kg} is a generalized autoregressive parameter (GARP) to explain the serial dependence of responses across time. Moreover,
φ_{itj, kg} = α_{kg} \mathbf{1}\{|t-j|=1\} ,\; \log(σ_{itk}) = λ_{k0} + λ_{k1} h_{it}, \; \log≤ft(\frac{ω_{ilm}}{π-ω_{ilm}}\right) = ν_l + ν_m.
The priors for the parameters in the model given by
{β} \sim \mathcal{N}(0, σ_β^2 I);
{λ}_k \sim \mathcal{N}(0, σ_λ^2 I);
{ν}_k \sim \mathcal{N}(0, σ_ν^2 I), \quad k=1, …, K,
where σ_β^2, σ_λ^2, and σ_ν^2 are prespecified values. For k, g = 1, …, K and m=1, …, a, we further assume
α_{kgm} \sim δ_{kgm} \mathcal{N}(0, σ^2_δ) + (1-δ_{kgm})η_0,
where σ^2_δ is prespecified value and η_0 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.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 | 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)
|
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.