mleCOP: Maximum Pseudo-Log-Likelihood Estimation for Copula Parameter Estimation

Description

Perform maximum pseudo-log-likelihood estimation (pMLE) for copula parameters by maximizing the function: $$\mathcal{L}(\Theta_p) = \sum_{i=1}^n \log\bigl[ c(F_x(x_i), F_y(y_i); \Theta_p)\bigr]\mbox{,}$$ where $\mathcal{L}(\Theta_p)$ is the log-likelihood for parameter vector $\Theta_p$ of dimension $p$, and $c(u,v; \Theta_p)$ is the bivariate copula density. The $u$ and $v$ are estimated by the respective empirical cumulative distribution functions $u = F_x(\cdots)$ and $v = F_y(\cdots)$ for each of the joint realizations of a sample of size $n$. The $c(u,v)$ is numerically estimated by the copula using the densityCOP function.

Usage

mleCOP(u, v=NULL, cop=NULL, parafn=function(k) return(k),
          interval=NULL, init.para=NULL, verbose=FALSE, control=list(),
          the.zero=.Machine$double.eps^0.25, ...)

Arguments

Nonexceedance probability $u$ in the $X$ direction;

Nonexceedance probability $v$ in the $Y$ direction and if NULL then u is treated as a two column R data.frame;

cop

A copula function;

parafn

A function responsible for generating the parameters. This is often just a simple return of a parameter vector as copBasic uses this style of parameterization, but this function can take over parameter remapping to handle boundary conditions to benefit the search or provide an interface into other copula packages in R (see Examples);

interval

The search interval for root finding if the parameter dimension of the copula is $p = 1$. The interval is not used for $p \ge 2$;

init.para

The initial guesses for the parameters for the $p$-dimensional optimization for $p \ge 2$. The initial guess is used if the parameter dimension of the copula is $p = 1$ and interval is NULL (see Examples);

verbose

A logical to trigger an minimal reporting scheme within the objective function. This is independent from the control$trace of function optim();

control

This argument is the argument of the same name for optim();

the.zero

The value for “the zero” of the copula density function. This argument is the argument of the same name for densityCOP. The default here is intended to suggest that a tiny nonzero value for density will trap the numerical zero densities; and

...

Additional arguments to pass, see source code for the internally used functions that can pick these additional arguments up.

Value

The value(s) for the estimated parameters are returned within an R list where the elements listed below are populated unique to this package. The other elements of the returned list are generated from either the optimise() or optim() functions of R.

para

The parameter(s) in a canonical element after the one-dimensional root finding ($p = 1$) or multi-dimensional optimization ($p \ge 2$) solutions are passed through parafn;

packagetext

A helpful message unique to the copBasic package;

loglik

The maximum of the log-likelihood matching the name for the same quantity by the function fitCopula in package copula though a separate implementation is used in copBasic;

AIC

Akaike information criterion (AIC) (see also aicCOP): $\mathrm{AIC} = 2p - 2\mathcal{L}(\Theta_p)$; and

BIC

Bayesian information criterion (BIC) (see also bicCOP): $\mathrm{BIC} = p\log(n) - 2\mathcal{L}(\Theta_p)$.

References

Joe, H., 2014, Dependence modeling with copulas: Boca Raton, CRC Press, 462 p.

Examples

Run this code

# NOT RUN {
set.seed(162); UV <- simCOP(n=188, cop=PLACKETTcop, para=5.6)
ML <- mleCOP(UV$U, UV$V, cop=PLACKETTcop, interval=c(0.1,40)) # 5.225459 estimated

Thetas <- 10^(seq(log10(0.001),log10(100), by=.005))
MLs <- sapply(Thetas, function(k)
              densityCOP(UV$U, UV$V, cop=PLACKETTcop, para=k, sumlogs=TRUE))
plot(Thetas, MLs, log="x", type="l", # draw the pMLE solution process
     xlab="Plackett Theta", ylab="sum of log densities")
lines(rep(ML$para, 2), c(ML$objective, par()$usr[3]), col=2) # 
# }
# NOT RUN {
# }
# NOT RUN {
set.seed(149); UV <- simCOP(1000, cop=CLcop, para=pi)
# Warning messages about using optim() for 1D solution
mleCOP(UV, cop=CLcop, init.para=2)$para          # 3.082031
# No warning message, optimise() called instead.
mleCOP(UV, cop=CLcop, interval=c(0,1E2))$para    # 3.081699 
# }
# NOT RUN {
# See extended code listings and discussion in the Note section
# }