CondiCopLocFun {LocalCop}

Create a TMB local likelihood function.

Description

Wraps a call to TMB::MakeADFun().

Usage

CondiCopLocFun(u1, u2, family, x, x0, wgt, degree = 1, eta, nu)

Arguments

u1

Vector of first uniform response.

u2

Vector of second uniform response.

family

An integer defining the bivariate copula family to use. See ConvertPar().

x

Vector of observed covariate values.

x0

Scalar covariate value at which to evaluate the local likelihood. Does not have to be a subset of x.

wgt

Vector of positive kernel weights.

degree

Integer specifying the polynomial order of the local likelihood function. Currently only 0 and 1 are supported.

eta

Value of the copula dependence parameter. Scalar or vector of length two, depending on whether degree is 0 or 1.

nu

Value of the other copula parameter. Scalar or vector of same length as u1. Ignored if family != 2.

Value

A list as returned by a call to TMB::MakeADFun(). In particular, this contains elements fun and gr for the negative local likelihood and its gradient with respect to eta.

Examples

# the following example shows how to create
# an unconditional copula likelihood function

# simulate data
n <- 1000 # sample size
family <- 2 # Student-t copula
rho <- runif(1, -1, 1) # unconditional dependence parameter
nu <- runif(1, 4, 20)# degrees of freedom parameter
udata <- VineCopula::BiCopSim(n, family = family, par = rho, par2 = nu)

# create likelihood function

# parameter conversion: equivalent to BiCopPar2Eta(family = 2, ...)
rho2eta <- function(rho) .5 * log((1+rho)/(1-rho))
nll_obj <- CondiCopLocFun(u1 = udata[,1], u2 = udata[,2], family = family,
                          x = rep(0, n), x0 = 0, # centered covariate x - x0 == 0
                          wgt = rep(1, n), # unweighted
                          degree = 0, # zero-order fit
                          eta = c(rho2eta(rho), 0),
                          nu = nu)

# likelihood function: recall that TMB requires a _negative_ ll
stucop_lik <- function(rho) {
  -nll_obj$fn(c(rho2eta(rho), 0))
}

# compare to VineCopula.
rhovec <- runif(50, -1, 1)
system.time({
  ll1 <- sapply(rhovec, stucop_lik) # LocalCop
})
system.time({
  ll2 <- sapply(rhovec, function(rho) {
    # VineCopula
    sum(log(VineCopula::BiCopPDF(u1 = udata[,1], u2 = udata[,2],
                                 family = family,
                                 par = rho, par2 = nu)))
  })
})

# difference between the two
range(ll1 - ll2)