Fitting functions for selm models


A call to selm activates a call to and from here to some other function which actually performs the parameter search, among those listed below. These lower-level functions can be called directly for increased efficiency, at the expense of some more programming effort and lack of methods for the returned object.


22, y, family = "SN", start = NULL, w, fixed.param = list(), 
   offset = NULL, selm.control)

sn.mple(x, y, cp = NULL, w, penalty = NULL, trace = FALSE, opt.method =
   c("nlminb",  "Nelder-Mead", "BFGS", "CG", "SANN"), control = list()) 

st.mple(x, y, dp = NULL, w, = NULL, symmetr = FALSE, penalty = NULL, 
   trace = FALSE, opt.method = c("nlminb", "Nelder-Mead", "BFGS", "CG", "SANN"), 
   control = list()) 

msn.mle(x, y, start = NULL, w, trace = FALSE, opt.method = c("nlminb", 
   "Nelder-Mead", "BFGS", "CG", "SANN"), control = list())

msn.mple(x, y, start = NULL, w, trace = FALSE, penalty = NULL, 
   opt.method = c("nlminb", "Nelder-Mead", "BFGS", "CG", "SANN"), 
   control = list()) 

mst.mple(x, y, start = NULL, w, = NULL, symmetr=FALSE, 
   penalty = NULL, trace = FALSE, 
   opt.method = c("nlminb", "Nelder-Mead", "BFGS", "CG", "SANN"), 
   control = list()) 



a full-rank design matrix with the first column of all 1's.


a vector or a matrix of response values such that NROW(y)=nrow(x).


a character string which selects the parametric family of distributions assumed for the error term of the regression model. It must one of "SN" (default), "ST" or "SC", which correspond to the skew-normal, the skew-t and the skew-Cauchy family, respectively. See makeSECdistr for more information on these families and the skew-elliptically contoured (SEC) distributions; notice that family "ESN" is not allowed here.

start, dp, cp

a vector or a list of initial parameter values, depeding whether y is a vector or a matrix. It is assumed that cp is given in the CP parameterization, dp and start in the DP parameterization.


a vector of non-negative integer weights of length equal to NROW(y); if missing, a vector of all 1's is generated.


a list of assignments of parameter values to be kept fixed during the optimization process. Currently, there is only one such option, namely fixed.param=list(nu='value'), to fix the degrees of freedom at the named 'value' when family="ST", for instance list(nu=3). Setting fixed.param=list(nu=1) is equivalent to select family="SC".


an optional character string with the name of the penalty function of the log-likelihood; default value NULL corresponds to no penalty.


this can be used to specify an a priori known component to be included in the linear predictor during fitting. This should be NULL or a numeric vector of length equal to the number of cases. One or more offset terms can be included in the formula instead or as well, and if more than one are specified their sum is used.


a logical value which regulates printing of successive calls to the target function; default value is FALSE which suppresses printing.

a positive value to keep fixed the parameter nu of the ST distribution in the optimization process; with default value NULL, nu is estimated like the other parameters.


a logical flag indicating whether a contraint of symmetry is imposed on the slant parameter; default is symmetr=FALSE.


a character string which selects the optimization method within the set c("nlminb", "Nelder-Mead", "BFGS", "CG", "SANN"); the last four of these are "methods" of function optim.


a list whose components regulate the working of; see ‘Details’ for their description;


a list of control items passed to the optimization function.


A call to selm produces a call to which selects the appropriate function among sn.mple, st.mple, msn.mle, msn.mple, mst.mple, depending on the arguments of the calling statement. In the adopted scheme for function names, msn refers to a multivariate skew-normal distribution and mst refers to a multivariate skew-t distribution, while mle and mple refers to maximum likelihood and maximum penalized likelihood estimation, respectively. Of these functions, sn.mple works in CP space; the others in the DP space. In all cases, a correspondig mapping to the alternative parameter space is performed before exiting, in addition to the selected parameter set.

The components of selm.control are as follows:

  • method: the estimation method, "MLE" or "MPLE".

  • penalty: a string with the name of the penalty function.

  • info.type: a string with the name of the information matrix, "observed" or "expected"; currently fixed at "observed".

  • opt.method: a character string which selects the optimization method.

  • opt.control: a list of control parameters of opt.method.

Function msn.mle, for MLE estimation of linear models with SN errors, is unchanged from version 0.4-x of the package. Function msn.mple is similar to msn.mle but allows to introduce a penalization of the log-likelihood; when penalty=NULL, a call to msn.mle is more efficient. Functions sn.mple and mst.mple work like sn.mle and mst.mle in version 0.4-x if the argument penalty is not set or it is set to NULL, except that mst.mple does not handle a univariate response (use st.mple for that).


A list whose specific components depend on the named function. Typical components are:


the calling statement


vector or list of estimated DP parameters


vector or list of estimated CP parameters


the maximized (penalized) log-likelihood


a list with auxiliary output values, depending on the function


a list produced by the numerical opt.method


Computational aspects of maximum likelihood estimation for univariate SN distributions are discussed in Section 3.1.7 of Azzalini and Capitanio (2014). The working of sn.mple follows these lines; maximization is performed in the CP space. All other functions operate on the DP space.

The technique underlying msn.mle is based on a partial analytical maximization, leading implicitly to a form of profile log-likelihood. This scheme is formulated in detail in Section 6.1 of Azzalini and Capitanio (1999) and summarized in Section 5.2.1 of Azzalini and Capitanio (2014). The same procedure is not feasible when one adopts MPLE; hence function msn.mple has to maximize over a larger parameter space.

Maximization of the univariate ST log-likelihood is speeded-up by using the expressions of the gradient given by DiCiccio and Monti (2011), reproduced with inessential variants in Section 4.3.3 of Azzalini and Capitanio (2014).

The working of mst.mple is based on a re-parameterization described in Section 5.1 of Azzalini and Capitanio (2003). The expressions of the corresponding log-likelihood derivatives are given in Appendix B of the full version of the paper.


Adelchi Azzalini


Azzalini, A. and Capitanio, A. (1999). Statistical applications of the multivariate skew normal distribution. J.Roy.Statist.Soc. B 61, 579–602. Full-length version available at

Azzalini, A. and Capitanio, A. (2003). Distributions generated by perturbation of symmetry with emphasis on a multivariate skew t distribution. J.Roy. Statist. Soc. B 65, 367–389. Full-length version available at

Azzalini, A. with the collaboration of Capitanio, A. (2014). The Skew-Normal and Related Families. Cambridge University Press, IMS Monographs series.

DiCiccio, T. J. and Monti, A. C. (2011). Inferential aspects of the skew t-distribution. Quaderni di Statistica 13, 1–21.

See Also

selm for a comprehensive higher level fitting function,

Qpenalty for specification of a penalty function


data(wines, package="sn")
X <- model.matrix(~ phenols + wine, data=wines)
fit <- msn.mle(x=X, y=cbind(wines$acidity, wines$alcohol), opt.method="BFGS")
fit <- st.mple(x=X, y = wines$acidity,, penalty="Qpenalty")

Want to suggest features or report bugs for Use the GitHub issue tracker.

comments powered by Disqus