R/gamBiCopCDF.R
gamBiCopCDF.RdThis function returns the distribution function of a bivariate conditional copula, where either the copula parameter or the Kendall's tau is modeled as a function of the covariates.
gamBiCopCDF(object, newdata = NULL)gamBiCop-class object.
(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. If this is not provided then the distribution 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.
The conditional density.
gamBiCop and gamBiCopPredict.
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.182 0.001 0.182
## Evaluate the conditional density
gamBiCopCDF(fit$res)
#> [1] 0.0198474375 0.6717537658 0.6685440525 0.8345873257 0.8288132359
#> [6] 0.1235910766 0.3441990142 0.0931108177 0.5371458487 0.0326774029
#> [11] 0.0259161934 0.0377845745 0.4509528719 0.0138858221 0.3480052303
#> [16] 0.4402536121 0.2161366785 0.0016261750 0.1521913659 0.2440282674
#> [21] 0.3294484333 0.3719879929 0.2957705670 0.5313398868 0.7566311632
#> [26] 0.0243851130 0.0955629741 0.8173055113 0.1199735852 0.4009586262
#> [31] 0.0622281033 0.5543639528 0.1511671716 0.2644151328 0.5553504298
#> [36] 0.1661333398 0.3235422624 0.8569511730 0.1808727492 0.6168688631
#> [41] 0.0099209913 0.4246454533 0.6223504187 0.1060099016 0.6193154126
#> [46] 0.1098135079 0.8964376697 0.0536940127 0.3790518202 0.1528715314
#> [51] 0.0194698895 0.0880012963 0.9738834597 0.1024866326 0.4712303592
#> [56] 0.0048217453 0.3219375695 0.6458725253 0.3022159033 0.1280493367
#> [61] 0.3938457411 0.0965470296 0.2201418364 0.3321879357 0.0119953113
#> [66] 0.8382381448 0.5141623637 0.3232156953 0.5513177846 0.9094837475
#> [71] 0.0417468334 0.1491592997 0.0848606389 0.0184211738 0.1475870185
#> [76] 0.4827166412 0.0277150043 0.0442711944 0.1607491179 0.0528071527
#> [81] 0.9123747087 0.2834527256 0.0346594904 0.2685566180 0.1052068757
#> [86] 0.2012448481 0.0802050066 0.3464284716 0.7555537014 0.4262663030
#> [91] 0.5681282335 0.6173471660 0.0305412671 0.0856821681 0.5049142300
#> [96] 0.4368720322 0.0040903705 0.9445065259 0.1772638256 0.6066537144
#> [101] 0.1132383913 0.6463693253 0.2849581116 0.6005626829 0.0220945385
#> [106] 0.0257241932 0.2874825659 0.3583115960 0.5341512750 0.3833252263
#> [111] 0.3749509560 0.6787883055 0.1557347164 0.0092346646 0.3573093619
#> [116] 0.0512652071 0.7721101236 0.3568639660 0.3086761572 0.9001894702
#> [121] 0.1824546365 0.0687378171 0.6644245674 0.2206306100 0.1114148343
#> [126] 0.2450142303 0.1195683798 0.6348976301 0.0736029194 0.0089034132
#> [131] 0.5231052513 0.4191334309 0.3087568699 0.0944031781 0.5283637504
#> [136] 0.1040881634 0.5251563778 0.0154027157 0.2684306439 0.2635202325
#> [141] 0.3683145107 0.0030515943 0.3481829420 0.4151551369 0.6928313906
#> [146] 0.2197741788 0.1571349383 0.3791930871 0.6612060001 0.8779142171
#> [151] 0.0298212350 0.7401627441 0.3025652182 0.3399796660 0.1900857820
#> [156] 0.0231475565 0.6479278942 0.8367285201 0.2972243388 0.2672207409
#> [161] 0.4121781158 0.3797378066 0.0112251727 0.1262510428 0.7545334331
#> [166] 0.8872254648 0.2747487605 0.0522422264 0.0779799630 0.4700939776
#> [171] 0.0706914601 0.0006675355 0.7895479317 0.3743260487 0.3768860792
#> [176] 0.1698229332 0.0897216653 0.2157991970 0.7839786382 0.9695606129
#> [181] 0.5005304386 0.7630667689 0.0212698039 0.7220052649 0.8567693316
#> [186] 0.0233138132 0.2111615534 0.5491288271 0.9636177902 0.1002917085
#> [191] 0.6855137451 0.8972648574 0.5182638278 0.6045079364 0.3383737628
#> [196] 0.2850117547 0.6727221045 0.0420619042 0.0013328707 0.2277666265