Conditional density function of a Generalized Additive model for the copula parameter or Kendall's tau

Source: R/gamBiCopPDF.R

This function returns the density of a bivariate conditional copula, where either the copula parameter or the Kendall's tau is modeled as a function of the covariates.

gamBiCopPDF(object, newdata = NULL)

Arguments

object

gamBiCop-class object.
newdata

(Same as in predict.gam from the mgcv package) A matrix or data frame containing the values of the model covariates at which predictions are required, along with two columns named `"u1"` and `"u2"`. If this is not provided then the density corresponding to the original data are returned. If newdata is provided then it should contain all the variables needed for prediction: a warning is generated if not.

Value

The conditional density.

See also

gamBiCop and gamBiCopPredict.

Examples

require(copula)
set.seed(0)

## Simulation parameters (sample size, correlation between covariates,
## Gaussian copula family)
n <- 2e2
rho <- 0.5
fam <- 1

## A calibration surface depending on three variables
eta0 <- 1
calib.surf <- list(
  calib.quad <- function(t, Ti = 0, Tf = 1, b = 8) {
    Tm <- (Tf - Ti) / 2
    a <- -(b / 3) * (Tf^2 - 3 * Tf * Tm + 3 * Tm^2)
    return(a + b * (t - Tm)^2)
  },
  calib.sin <- function(t, Ti = 0, Tf = 1, b = 1, f = 1) {
    a <- b * (1 - 2 * Tf * pi / (f * Tf * pi +
      cos(2 * f * pi * (Tf - Ti))
      - cos(2 * f * pi * Ti)))
    return((a + b) / 2 + (b - a) * sin(2 * f * pi * (t - Ti)) / 2)
  },
  calib.exp <- function(t, Ti = 0, Tf = 1, b = 2, s = Tf / 8) {
    Tm <- (Tf - Ti) / 2
    a <- (b * s * sqrt(2 * pi) / Tf) * (pnorm(0, Tm, s) - pnorm(Tf, Tm, s))
    return(a + b * exp(-(t - Tm)^2 / (2 * s^2)))
  }
)

## 3-dimensional matrix X of covariates
covariates.distr <- mvdc(normalCopula(rho, dim = 3),
  c("unif"), list(list(min = 0, max = 1)),
  marginsIdentical = TRUE
)
X <- rMvdc(n, covariates.distr)
colnames(X) <- paste("x", 1:3, sep = "")

## U in [0,1]x[0,1] with copula parameter depending on X
U <- condBiCopSim(fam, function(x1, x2, x3) {
  eta0 + sum(mapply(function(f, x)
    f(x), calib.surf, c(x1, x2, x3)))
}, X[, 1:3], par2 = 6, return.par = TRUE)

## Merge U and X
data <- data.frame(U$data, X)
names(data) <- c(paste("u", 1:2, sep = ""), paste("x", 1:3, sep = ""))

## Model fit with penalized cubic splines (via min GCV)
basis <- c(3, 10, 10)
formula <- ~ s(x1, k = basis[1], bs = "cr") +
  s(x2, k = basis[2], bs = "cr") +
  s(x3, k = basis[3], bs = "cr")
system.time(fit <- gamBiCopFit(data, formula, fam))
#>    user  system elapsed 
#>   0.149   0.000   0.149 

## Evaluate the conditional density
gamBiCopPDF(fit$res)
#>   [1]   1.53897409   0.65895255   7.54308919   7.29902794   7.98200466
#>   [6]   0.91280983   0.65051146   1.08273677   0.13773165   8.67048283
#>  [11]   4.95057615   3.20943703   1.00633039   0.84382926   0.40244706
#>  [16]   2.76005256   1.96292286  11.68003966   1.47383017   1.07385649
#>  [21]   0.99519798   1.18493994   1.44425297   0.83340597   2.40090338
#>  [26]   6.05695932   1.41528126   3.00012259   1.74635340   0.89680036
#>  [31]  10.49237658   1.42704608   1.01410397   1.20991114   0.31228006
#>  [36]   1.94781172   1.08920172  10.33179260   0.96340114   1.10652903
#>  [41]   1.19280224   1.53050049   1.15345886   1.49250150   3.65567450
#>  [46]   1.15697132  10.71967518   0.69294323   1.64953687   1.04327879
#>  [51]  88.38436240   3.62295213  10.47617993   0.77327325   2.40824831
#>  [56]   0.08006449   1.37469193   3.00543370   0.85171208   1.12234724
#>  [61]   1.60635642   0.83929746   2.41797654   0.81561166   4.07324027
#>  [66]  32.31728518   0.77131922   1.03527579   2.91140932   2.50619363
#>  [71]   1.14147003   1.34993051   5.42784075   1.65326684   1.27372418
#>  [76]   1.50209872   1.03634247   1.00381308   1.44370479   5.84889100
#>  [81]   5.19282176   0.76370453 100.00992066   2.79801297   1.26603622
#>  [86]   0.99666029   5.95339495   0.77555019   1.29134500   0.85042384
#>  [91]   1.68855016  13.22072485   1.30870031   1.09627773   0.98611410
#>  [96]   1.64269684   1.25826086  18.73796411   1.26322251   1.20958580
#> [101]   1.14939531   0.56989490   1.78403902   1.20540159   9.17700886
#> [106]   1.86256652   2.08354774   1.27445574   0.81208319   1.09719242
#> [111]   0.70119548   1.42652282   1.74619315  31.69259035   1.57919837
#> [116]   0.07155906   3.75765970   1.06833619   1.06277686   2.75424876
#> [121]   1.03289717   0.98409803   1.06658308   4.13468165   1.43183856
#> [126]   1.14460591   1.62000488   1.72714023   1.29943373   2.45416520
#> [131]   1.39534538   1.55503242  10.17514164   1.39641877   0.96021880
#> [136]   6.23700638   0.96829628   2.18752443   1.04595150   1.03916368
#> [141]   0.31373842   1.91360632   1.04706871   1.27716309  12.45134113
#> [146]   0.55571760   1.05842882   1.18931050   1.15449172   2.97966281
#> [151]   0.98788715   2.10243551   1.19538471   1.03388254   0.75967887
#> [156]   0.91491599   0.61518975   5.93053002   1.03114386   2.00393186
#> [161]  22.53256040   0.96323476   7.46971594   1.82261251   2.34808682
#> [166]   2.74248586   1.36414169   3.88879034   1.18696751   1.61579038
#> [171]   1.49307845  88.82655883   1.35422005   0.95393911   1.77950413
#> [176]   0.91590531   2.24707243   3.45363709   3.73011797  11.21433096
#> [181]   1.21093188   2.96206743   5.20979379   1.41870963   4.06028494
#> [186]   0.82568046   1.28401255   1.14938890   6.06026205   7.85370444
#> [191]   1.55187793   2.27710006   1.03309736   1.98248935   2.15382610
#> [196]   1.10177792   1.62369088   2.79315172   2.37038226   0.62678583