| Title: | Bayesian Inference for the Multivariate Skew-t Model |
|---|---|
| Description: | Estimates the multivariate skew-t and nested models, as described in the articles Liseo, B., Parisi, A. (2013). Bayesian inference for the multivariate skew-normal model: a population Monte Carlo approach. Comput. Statist. Data Anal. <doi:10.1016/j.csda.2013.02.007> and in Parisi, A., Liseo, B. (2017). Objective Bayesian analysis for the multivariate skew-t model. Statistical Methods & Applications <doi: 10.1007/s10260-017-0404-0>. |
| Authors: | Antonio Parisi [aut, cre], Brunero Liseo [aut], Dirk Eddelbuettel [ctb], Romain Francois [ctb] |
| Maintainer: | Antonio Parisi <[email protected]> |
| License: | GPL-3 |
| Version: | 1.1.1 |
| Built: | 2024-11-27 06:29:06 UTC |
| Source: | CRAN |
Estimates the multivariate skew-t and nested models, as described in the articles Liseo, B., Parisi, A. (2013). Bayesian inference for the multivariate skew-normal model: a population Monte Carlo approach. Comput. Statist. Data Anal. <doi:10.1016/j.csda.2013.02.007> and in Parisi, A., Liseo, B. (2017). Objective Bayesian analysis for the multivariate skew-t model. Statistical Methods & Applications <doi: 10.1007/s10260-017-0404-0>.
| Package: | mvst |
| Type: | Package |
| Title: | Bayesian Inference for the Multivariate Skew-t Model |
| Version: | 1.1.1 |
| Date: | 2023-12-05 |
| Author: | Antonio Parisi [aut, cre], Brunero Liseo [aut], Dirk Eddelbuettel [ctb], Romain Francois [ctb] |
| Maintainer: | Antonio Parisi <[email protected]> |
| Description: | Estimates the multivariate skew-t and nested models, as described in the articles Liseo, B., Parisi, A. (2013). Bayesian inference for the multivariate skew-normal model: a population Monte Carlo approach. Comput. Statist. Data Anal. <doi:10.1016/j.csda.2013.02.007> and in Parisi, A., Liseo, B. (2017). Objective Bayesian analysis for the multivariate skew-t model. Statistical Methods & Applications <doi: 10.1007/s10260-017-0404-0>. |
| SystemRequirements: | GNU Scientific Library |
| License: | GPL-3 |
| Imports: | MCMCpack, mvtnorm, mnormt |
| NeedsCompilation: | yes |
| Packaged: | 2023-12-05 13:29:21 UTC; user |
| Repository: | CRAN |
| Date/Publication: | 2023-12-05 23:10:02 UTC |
| Config/pak/sysreqs: | libgsl0-dev |
Index of help topics:
MNmargLike Marginal Likelihood for the Multivariate Normal
Model.
bivPlot Marginal and joint plots for bivariate data.
cmlSE CML for the parameters of a p-variate
Skew-Elliptical model.
coef.mcSEsummary Extract mcSE Model Coefficients.
dmvSE Density function for the SE distributions.
mcSE MC sampler for a p-variate Skew-Elliptical
model.
mvst-package Bayesian Inference for the Multivariate Skew-t
Model
rmvSE Random generation from a SE distribution.
summary.mcSE Summary function for mcSE objects.
Antonio Parisi [aut, cre], Brunero Liseo [aut], Dirk Eddelbuettel [ctb], Romain Francois [ctb]
Maintainer: Antonio Parisi <[email protected]>
Parisi A, Liseo B (2017). Objective Bayesian Analysis for the Multivariate Skew-t Model. Statistical Methods & Applications. ISSN 1613-981X. doi:10.1007/s10260-017-0404-0 Parisi, A., Liseo, B. (2018). Statistical Inference with Skew-t Distributions: the MVST R Package. Annali del Dipartimento di Metodi e Modelli per l'Economia il Territorio e la Finanza. ISSN 2385-0825. Liseo, B., Parisi, A. (2013). Bayesian inference for the multivariate skew-normal model: a population Monte Carlo approach. Comput. Statist. Data Anal. <doi:10.1016/j.csda.2013.02.007>
Scatterplot and marginal histograms for bivariate data. If theta is provided, the joint and marginal densities will be superimposed.
bivPlot(y, modelType=NULL, theta=NULL)bivPlot(y, modelType=NULL, theta=NULL)
y |
Data matrix or vector. |
modelType |
(optional) model specification. Already implemented modelTypes are 'N' (Normal), 'SN' (skew-normal), 'T' (Student T), and 'ST' (skew-t). |
theta |
(optional) list of parameters. The list should contain elements named 'xi' (a numeric vector) or 'B' (a kxp matrix), 'G' (a pxp matrix), 'psi' (a numeric vetor, optional) and 'nu' (a scalar, optional). See Details. |
The function draws a plot for bivariate data.
# Define the parameters' list pars = list(xi=c(5,2), G=diag(2), psi=rep(1,2), nu=4) # Generate data values = rmvSE(200, 2, NULL, 'ST', theta=pars) y = values$y # Draw the data points. bivPlot(y) # Draw the data points and the density function. bivPlot(y, modelType='ST', theta=pars)# Define the parameters' list pars = list(xi=c(5,2), G=diag(2), psi=rep(1,2), nu=4) # Generate data values = rmvSE(200, 2, NULL, 'ST', theta=pars) y = values$y # Draw the data points. bivPlot(y) # Draw the data points and the density function. bivPlot(y, modelType='ST', theta=pars)
Complete Maximum Likelihood for the parameters of a p-variate Skew-Elliptical model.
cmlSE(modelType, y, z=NULL, v=NULL, X=NULL)cmlSE(modelType, y, z=NULL, v=NULL, X=NULL)
modelType |
character indicating the model. Possible values are "N", "SN", "T", "ST" |
y |
data matrix. |
z |
(optional) numeric vector. |
v |
(optional) numeric vector. |
X |
(optional) design matrix. |
Given the value of the latent variables z and v, the function returns a list containing the estimates for the required model. Where available, a design matrix with the value of the covariates can be provided. In this case, the parameters of a regression model with skewed errors are estimated.
Parisi, A. and Liseo, B. (2017) "Objective Bayesian Analysis for the Multivariate Skew-t Model" Statistical Methods & Applications
## Generate artificial data pars = list(xi=c(3,5), psi=c(2,4), G=diag(2), nu=6) values = rmvSE(n=20, p=2, modelType='ST', theta=pars) ## CML estimates for pars thetaHat = cmlSE(modelType='ST', y=values$y, z=values$z, v=values$v)## Generate artificial data pars = list(xi=c(3,5), psi=c(2,4), G=diag(2), nu=6) values = rmvSE(n=20, p=2, modelType='ST', theta=pars) ## CML estimates for pars thetaHat = cmlSE(modelType='ST', y=values$y, z=values$z, v=values$v)
The point estimates for the model parameters are obtained from mcSE summary objects.
## S3 method for class 'mcSEsummary' coef(object, ...)## S3 method for class 'mcSEsummary' coef(object, ...)
object |
an object of class mcSEsummary. |
... |
further arguments passed to or from other methods. |
A list containing the point estimates for the estimated model.
# Generate ST-distributed data (including the value of the latent variables) pars = list(xi=c(2,2), G=diag(2), psi=c(0.3,0.5), nu=5) values = rmvSE(n=100, p=2, modelType='ST', theta=pars) # Estimate a Skew-t model (not run) # fit = mcSE(y=values$y, X=NULL, N=20000, Ti=3, modelType='ST', warmUp=TRUE) # stats = summary(fit) # coef(stats)# Generate ST-distributed data (including the value of the latent variables) pars = list(xi=c(2,2), G=diag(2), psi=c(0.3,0.5), nu=5) values = rmvSE(n=100, p=2, modelType='ST', theta=pars) # Estimate a Skew-t model (not run) # fit = mcSE(y=values$y, X=NULL, N=20000, Ti=3, modelType='ST', warmUp=TRUE) # stats = summary(fit) # coef(stats)
This function computes the density function for p-variate Skew-Elliptical variables.
dmvSE(y, X=NULL, modelType, theta, LOG=FALSE)dmvSE(y, X=NULL, modelType, theta, LOG=FALSE)
y |
vector of quantiles. If y is a matrix, each row is considered as a quantile. |
X |
(optional) a design matrix. |
modelType |
Model specification. Already implemented modelTypes are 'N' (Normal), 'SN' (skew-normal), 'T' (Student T), and 'ST' (skew-t). |
theta |
list of parameters. The list should contain elements named 'xi' (a numeric vector) or 'B' (a kxp matrix), 'G' (a pxp matrix), 'psi' (a numeric vetor, optional) and 'nu' (a scalar, optional). See Details. |
LOG |
logical; if TRUE, log-densities are returned. |
A numeric vector with n values of the density function, one for each row in y.
Azzalini, A. and Capitanio, A. (2003) "Distributions generated by perturbation of symmetry with emphasis on a multivariate skew t distribution", JRSSB.
# Define the parameters' list pars = list(xi=c(5,2), G=diag(2), psi=rep(1,2), nu=4) # Generate data value = rmvSE(1, 2, NULL, 'ST', theta=pars) # Compute the density function in the point y dmvSE(y=value$y, X=NULL, modelType='ST', theta=pars, LOG=FALSE)# Define the parameters' list pars = list(xi=c(5,2), G=diag(2), psi=rep(1,2), nu=4) # Generate data value = rmvSE(1, 2, NULL, 'ST', theta=pars) # Compute the density function in the point y dmvSE(y=value$y, X=NULL, modelType='ST', theta=pars, LOG=FALSE)
MonteCarlo sampler for a p-variate Skew-Elliptical model.
mcSE(y, X=NULL, N, Ti, modelType='ST', warmUp=FALSE, control=list())mcSE(y, X=NULL, N, Ti, modelType='ST', warmUp=FALSE, control=list())
y |
data matrix. |
X |
design matrix. (default: NULL) |
N |
number of particles per iteration. |
Ti |
number of iterations. |
modelType |
model to be estimated (default: ST). See 'Details'. |
warmUp |
logical flag (default: FALSE). See 'Deatils' |
control |
(optional) a list of control parameters. See 'Details'. |
Already implemented modelTypes are 'N' (Normal), 'SN' (skew-normal), 'T' (Student T), and 'ST' (skew-t, the default). To estimate a regression model, an 'X' should be added: for example, 'STX' stands for a regression model with ST errors. For these models, the argument parTypes in the control list is overridden. It is however possible to implement other models; in this case, parType is required and should contain the names of the parameters of the model. The argument warmUp allows to run preliminary iterations with a smaller number of particles, in order to speed up the algorithm. The number of these iterations, and the number of particles can be controlled using the Nwu argument in the control list.
To estimate regression models with skewed errors, it is sufficient to specify the argument X, which should contain the design matrix.
The (optional) argument control can provide a list with the following elements
seedif different from NULL, sets the random seed for replicability purposes.
parInfodata.frame containing the informations about the model parameters. Each row of the data.frame should contain the names of the parameters, the type (u: 'univariate', m: 'multivariate', M: 'matrix-variate', SM: 'symmetric matrix-variate'), the number nCols of elements, or columns, of the parameter and the number nRows of rows (eventually 1). Required if the modelType is not already implemented.
propFuncsnamed character vector with the names of the functions for the proposal distributions. The names of the elements of propFuncs should be the relevant elements of the set ('z', 'v', 'xi', 'psi', 'G', 'nu'). Custom proposal functions should require at least four arguments: y, the data matrix, X, the covariates (NULL if not relevant), particles, that is the list of current values of the parameters, for each particle, and priorList, a list containing the hyperparameters of the prior distributions. It should return two objects: values (the proposed values for the parameter, for each particle) and log.dq (vector with the N (log-)values of the proposal density).
logPriorFuncname of the file containing the function to compute the value (in logarithms) of the posterior density for all the particles. It should depend on the objects y, particles and priorList, while it should return N values of the posterior density, in logarithms.
Nwunumeric vector with the number of particles for each warm-up iteration. Default value is rep(2000, 3). It is however overridden if warmUp is FALSE.
priorListlist of hyperparameters.
saveParticleslogical flag (default: FALSE) indicating whether the value of the particles proposed in each iteration should be saved. If TRUE, the folder 'Iterations' is created in the current directory.
outFolderthe folder in which the outputs are saved (if saveParticles is TRUE). The default folder is '/Output/Iterations'.
verboselogical flag (default: TRUE). If TRUE, details about the progress of the algorithm are printed.
The function returns
particles |
a list with the resampled particles of the last iteration, |
log.py |
the vector of the estimates of the marginal likelihood (one value per iteration), |
nResampled |
the vector of the number of resampled particles (one value per iteration). |
If saveParticles is TRUE, the lists of the sampled particles, the importance weights, and the indices of the resampled particles are saved in the folder specified in outFolder, or in the default folder '/Output'. If outFolder doesn't already exists, it will be created.
Parisi, A. and Liseo, B. (2017) "Objective Bayesian Analysis for the Multivariate Skew-t Model" Statistical Methods & Applications
Azzalini, A. and Arellano-Valle, R.B. (2013) "Maximum Penalized Likelihood Estimation for Skew-normal and Skew-t Distributions" J. Statist. Plann. Inference, 143 (2), 419–433.
## Generate artificial data pars = list(xi=c(3,5), psi=c(2,4), G=diag(2), nu=6) values = rmvSE(n=60, p=2, modelType='ST', theta=pars) ## Estimate a Skew-t model (not run) # fit = mcSE(y=values$y, N=20000, Ti=3, modelType='ST') # stats = summary(fit) # coef(stats)## Generate artificial data pars = list(xi=c(3,5), psi=c(2,4), G=diag(2), nu=6) values = rmvSE(n=60, p=2, modelType='ST', theta=pars) ## Estimate a Skew-t model (not run) # fit = mcSE(y=values$y, N=20000, Ti=3, modelType='ST') # stats = summary(fit) # coef(stats)
This function computes the exact marginal likelihood for Normally distributed data, under the default priors.
MNmargLike(y, X=NULL, LOG=FALSE)MNmargLike(y, X=NULL, LOG=FALSE)
y |
data matrix. |
X |
(optional) a design matrix. |
LOG |
logical; if TRUE, the log-marginal likelihood is returned. |
A scalar representing the marginal likelihood of a (multivariate) Normal model under the default priors for data y. If the design matrix X is provided, the function returns the marginal likelihood of a (multivariate) regression model with Normally distributed errors.
Liseo B, Parisi A (2013). Bayesian Inference for the Multivariate Skew-Normal Model: A Population Monte Carlo approach. Comput. Statist. Data Anal., 63, 125-138. ISSN 0167-9473. doi:10.1016/j.csda.2013.02.007.
# Generate Normally distributed data require(mvtnorm) y = rmvnorm(100, rep(2,2), diag(2)) # Marginal likelihood (exact value) MNmargLike(y, X=NULL, LOG=TRUE)# Generate Normally distributed data require(mvtnorm) y = rmvnorm(100, rep(2,2), diag(2)) # Marginal likelihood (exact value) MNmargLike(y, X=NULL, LOG=TRUE)
This function generates draws from a p-variate Skew-Elliptical distribution.
rmvSE(n, p, X=NULL, modelType, theta)rmvSE(n, p, X=NULL, modelType, theta)
n |
number of draws. |
p |
dimension of the drawn vectors. |
X |
a design matrix. |
modelType |
generating distribution. Already implemented modelTypes are 'N' (Normal), 'SN' (skew-normal), 'T' (Student T), and 'ST' (skew-t). |
theta |
list of parameters. The list should contain elements named 'xi' (numeric), 'G' (pxp matrix), 'psi' (numeric, optional) and 'nu' (scalar, optional). |
A list with three elements
yn x p matrix of the random draws from a p-variate SE distribution.
zvector of the latent values z (NULL for symmetric models)
vvector of the latent values v (NULL for the N and SN models)
Azzalini, A. and Capitanio, A. (2003) "Distributions generated by perturbation of symmetry with emphasis on a multivariate skew t distribution", JRSSB (see eq. 25).
## Generate artificial data pars = list(xi=c(3,5), psi=c(2,4), G=diag(2), nu=6) values = rmvSE(n=200, p=2, modelType='ST', theta=pars) ## X contains the data matrix and the vectors z and v of latent variables: y = values$y z = values$z v = values$v## Generate artificial data pars = list(xi=c(3,5), psi=c(2,4), G=diag(2), nu=6) values = rmvSE(n=200, p=2, modelType='ST', theta=pars) ## X contains the data matrix and the vectors z and v of latent variables: y = values$y z = values$z v = values$v
summary method for class "mcSE".
## S3 method for class 'mcSE' summary(object, ...)## S3 method for class 'mcSE' summary(object, ...)
object |
an object of class mcSE. |
... |
further arguments passed to or from other methods. |
A numeric vector with n values of the density function, one for each row in y.
Parisi A, Liseo B (2017). Objective Bayesian Analysis for the Multivariate Skew-t Model. Statistical Methods & Applications.
# Generate Normally distributed data pars = list(xi=c(2,2), G=diag(2), psi=c(0.3,0.5), nu=5) values = rmvSE(n=100, p=2, modelType='N', theta=pars) # Estimate a Skew-t model (not run) # fit = mcSE(y=values$y, X=NULL, N=20000, Ti=3, modelType='ST', warmUp=FALSE) # summary(fit)# Generate Normally distributed data pars = list(xi=c(2,2), G=diag(2), psi=c(0.3,0.5), nu=5) values = rmvSE(n=100, p=2, modelType='N', theta=pars) # Estimate a Skew-t model (not run) # fit = mcSE(y=values$y, X=NULL, N=20000, Ti=3, modelType='ST', warmUp=FALSE) # summary(fit)