lmrm {BNSP} | R Documentation |
Bayesian semiparametric modelling of covariance matrices for multivariate longitudinal data
Description
Implements an MCMC algorithm for posterior sampling based on a semiparametric model for continuous longitudinal multivariate responses. The overall model consists of 5 regression submodels and it utilizes spike-slab priors for variable selection and function regularization. See ‘Details’ section for a full description of the model.
Usage
lmrm(formula, data = list(), centre=TRUE, id, time,
sweeps, burn = 0, thin = 1, seed, StorageDir,
c.betaPrior = "IG(0.5,0.5*n*p)", pi.muPrior = "Beta(1,1)",
c.alphaPrior = "IG(1.1,1.1)", pi.phiPrior = "Beta(1,1)", c.psiPrior = "HN(2)",
sigmaPrior = "HN(2)", pi.sigmaPrior = "Beta(1,1)",
corr.Model = c("common", nClust = 1), DP.concPrior = "Gamma(5,2)",
c.etaPrior = "IG(0.5,0.5*samp)", pi.nuPrior = "Beta(1,1)",
pi.fiPrior = "Beta(1,1)", c.omegaPrior = "IG(1.1,1.1)", sigmaCorPrior = "HN(2)",
tuneCalpha, tuneSigma2, tuneCbeta, tuneAlpha, tuneSigma2R, tuneR, tuneCpsi,
tuneCbCor, tuneOmega, tuneComega, tau, FT = 1,...)
Arguments
formula |
a formula defining the responses and the covariates in the 5 regression models e.g. |
data |
a data frame. |
centre |
Binary indicator. If set equal to TRUE, the design matrices are centred, to have column mean equl to zero, otherwise, if set to FALSE, the columns are not centred. |
id |
identifiers of the individuals or other sampling units that are observed over time. |
time |
a vector input that specifies the time of observation |
sweeps |
total number of posterior samples, including those discarded in burn-in period
(see argument |
burn |
length of burn-in period. |
thin |
thinning parameter. |
seed |
optional seed for the random generator. |
StorageDir |
a required directory to store files with the posterior samples of models parameters. |
c.betaPrior |
The inverse Gamma prior of |
pi.muPrior |
The Beta prior of |
c.alphaPrior |
The inverse Gamma prior of |
pi.phiPrior |
The Beta prior of |
c.psiPrior |
The prior of |
sigmaPrior |
The prior of |
pi.sigmaPrior |
The Beta prior of |
corr.Model |
Specifies the model for the correlation matrices |
DP.concPrior |
The Gamma prior for the Dirichlet process concentration parameter. |
c.etaPrior |
The inverse Gamma prior of |
pi.nuPrior |
The Beta prior of |
pi.fiPrior |
The Beta prior of |
c.omegaPrior |
The prior of |
sigmaCorPrior |
The prior of |
tuneCalpha |
Starting value of the tuning parameter for sampling |
tuneSigma2 |
Starting value of the tuning parameter for sampling |
tuneCbeta |
Starting value of the tuning parameter for sampling |
tuneAlpha |
Starting value of the tuning parameter for sampling regression coefficients of the variance models
|
tuneSigma2R |
Starting value of the tuning parameter for sampling |
tuneR |
Starting value of the tuning parameter for sampling correlation matrices.
Defaults at |
tuneCpsi |
Starting value of the tuning parameter for sampling variances |
.
tuneCbCor |
Starting value of the tuning parameter for sampling |
tuneOmega |
Starting value of the tuning parameter for sampling regression coefficients of the variance models
|
tuneComega |
Starting value of the tuning parameter for sampling |
tau |
The |
FT |
Binary indicator. If set equal to 1, the Fisher's z transform of the correlations is modelled, otherwise if set equal to 0, the untransformed correlations are modelled. |
... |
Other options that will be ignored. |
Details
Function lmrm
returns samples from the posterior distributions of the parameters of
a regression model with normally distributed multivariate longitudinal responses. To describe the model,
let denote the vector of
responses
observed on individual
at time point
.
The observational time points
are allowed to be unequally spaced.
Further, let
denote the covariate vector that is observed along with
and that may include time,
other time-dependent covariates and time-independent ones.
In addition, let
denote the
th response vector.
With
and
cov
, the model assumes multivariate normality,
.
The means
and covariance matrices
are modelled semiparametrically in terms of covariates.
For the means one can specify semiparametric models,
This is the first of the 5 regression submodels.
To model the covariance matrix, first consider the modified Cholesky decomposition,
where
is a unit block lower triangular matrix
and
is a block diagonal matrix,
For modelling in terms of covariates,
first we separate the variances and the correlations
.
It is easy to model matrix
in terms of covariates as the only requirement on its diagonal
elements is that they are nonnegative,
This is the second of the 5 regression submodels.
For , the
element of
, one can specify semiparametric models
This is the third of the 5 regression submodels.
The elements of the correlations matrices are modelled in terms of covariate time only, hence they are denoted by
.
Subject to the positive definiteness constraint, the elements of
are modelled using a normal distribution
with location and scale parameters,
and
, modelled as
and these are the last 2 of the 5 submodels.
Value
Function lmrm
returns samples from the posteriors of the model parameters.
Author(s)
Georgios Papageorgiou gpapageo@gmail.com
References
Papageorgiou, G. (2020). Bayesian semiparametric modelling of covariance matrices for multivariate longitudinal data. arXiv:2012.09833.
Examples
# Fit a joint mean-covariance model on the simulated dataset simD2
require(ggplot2)
data(simD2)
model <- Y1 | Y2 ~ time | sm(time) | sm(lag) | sm(time) | 1
# the above defines the responses and the regression models on the left and
# right of "~", respectively
# the first model, for the mean, is a linear function of time, this is sufficient as
# the 2 responses have constant mean.
# the second model, for the variances, is a smooth function of time
# the third model, for the dependence structure, is a smooth function of lag,
# that lmrm figures out and it does not need to be computed by the user
# the fourth model, for location of the correlations, is a smooth function of time
# the fifth model, for scale of the correlations, is just an intercept model
## Not run:
m1 <- lmrm(formula = model, corr.Model = c("common", nClust = 1), data = simD2,
id = id, time = time, sweeps = 2500, burn = 500, thin = 2,
StorageDir = getwd(), seed = 1)
plot(m1)
## End(Not run)