#' @title Estimation of a non-central squared copula model #' #' @description This function estimates the copula parameter and the non-centrality parameters of a #' non-central squared copula #' #' @param y (nx2) data matrix (observations or residuals) that will be transformed to pseudo-observations #' @param family 'Gaussian' , 't' , 'Clayton' , 'Frank' , 'Gumbel' #' @param p number of non-centrality parameters to be estimated (p = 0,1,2) #' @param InitialValues initial values c(a1,a2,tau) to start the estimation; otherwise pre-selected values will be used #' #' @author Bouchra R. Nasri, August 14, 2019 #' @return \item{theta}{Estimated parameter of the copula according to CRAN copula package} #' @return \item{dof}{Estimated degrees of freedom, only for the Student copula} #' @return \item{tau}{Estimated theoretical Kendall tau for the copula family} #' #' @references Section 5.1 of Nasri, RĂ©millard & Bouezmarni (2019). Semi-parametric copula-based models under non-stationarity, #'Journal of Multivariate Analysis, 173, pages 347-365. #' #' @examples #' \donttest{ #' param <- c(0.8, 2.5, 0.7) ; #'U <- SimNCSCop('Clayton', 250, param) #'estimation <- EstNCSCop(U,'Clayton') #'} #' #' #' @importFrom copula iTau pcopula dcopula tau normalCopula tCopula frankCopula gumbelCopula claytonCopula #' @importFrom stats qnorm pnorm dnorm qchisq #' @export #' EstNCSCop<-function(y,family,p=2,InitialValues=NULL){ n = dim(y)[1]; d = dim(y)[2] U = floor( apply(y,2,rank) )/ (n + 1) #InitialValues alpha0 = ParamTau(family,0.5); #unconstrained parameters if (is.null(InitialValues)) { if (p==2){ initialGuess = c(finv( 1.0),finv(1.5), alpha0, 0); #a1 = 1.0, a2 = 1.5, tau = 0.5, nu = 12.5 } else if (p==1){ initialGuess =c(finv(1.5), alpha0, 0); #a1 = 1.0, a2 = 1.5, tau = 0.5, nu = 12.5 } else if (p==0){ initialGuess = c(alpha0, 0); # tau = 0.5, nu = 12.5 } } else{ a10 = InitialValues[1]; a20 = InitialValues[2]; tau0 = InitialValues[3]; nu0 = 12.5; initialGuess = c(finv(a10), finv(a20),ParamTau(family,tau0),-log(25/nu0 -1)); } # Optimisation alpha = stats::optim(initialGuess,function(x){ -sum(LoglikNCSCop(x,U,family, p))}, method = "Nelder-Mead",control=list(maxit=20000))$par # alpha_new = stats::optim(par= alpha,fun, method = "L-BFGS-B")$par dof = NaN; #theta = ParamCop(family,alpha[(p+1):(p+2)]); out = KendallTau(family,alpha[(p+1):(p+2)]); tau = out$tau; theta0 = out$theta; theta = theta0[1]; dof = theta0[2]; if(p>0){ a=f(alpha[1:p])} else{a=NaN} LL = sum(LoglikNCSCop(alpha,U,family, p)); out = list(tau=tau, a=a, theta=theta, dof=dof,LL=LL,alpha=alpha); return(out) } finv = function(x){ return( log(x/(3-x)) )} f = function(x){ return(3./(1+exp(-x)))}