[Vinecopula-commits] r111 - in pkg: R src src/include
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Do Aug 6 17:23:55 CEST 2015
Author: tnagler
Date: 2015-08-06 17:23:55 +0200 (Thu, 06 Aug 2015)
New Revision: 111
Modified:
pkg/R/BiCopCDF.r
pkg/R/BiCopCheck.R
pkg/R/BiCopDeriv.r
pkg/R/BiCopDeriv2.r
pkg/R/BiCopHfunc.r
pkg/R/BiCopHfuncDeriv.r
pkg/R/BiCopHfuncDeriv2.r
pkg/R/BiCopHinv.R
pkg/R/BiCopPDF.r
pkg/R/BiCopPar2Beta.r
pkg/R/BiCopPar2TailDep.r
pkg/R/BiCopPar2Tau.r
pkg/R/BiCopSim.R
pkg/R/BiCopTau2Par.r
pkg/src/deriv.c
pkg/src/deriv2.c
pkg/src/hfunc.c
pkg/src/hfuncderiv.c
pkg/src/hfuncderiv2.c
pkg/src/include/deriv.h
pkg/src/include/deriv2.h
pkg/src/include/hfunc.h
pkg/src/include/likelihood.h
pkg/src/include/logderiv.h
pkg/src/include/tcopuladeriv.h
pkg/src/include/tcopuladeriv_new.h
pkg/src/likelihood.c
pkg/src/logderiv.c
pkg/src/tcopuladeriv.c
pkg/src/tcopuladeriv_new.c
Log:
- correct call for non-t families in BiCopHfuncDeriv2.r(.., deriv ="par1u2"))
- use C function "Hinv1" instead of "pcc" for BiCopSim
- vectorize BiCop-functions w.r.t. family, par, par2:
* in C: BiCopPDF, BiCopHfunc, BiCopHinv, BiCopDeriv, BiCopDeriv2, BiCopHfuncDeriv, BiCopHfuncDeriv2
* in R: BiCopCDF, BiCopPar2Tau, BiCopPar2Beta, BiCopPar2TailDep
- add check.pars parameter to the above functions for the option to omit family/parameter consistency checks (for internal usage)
Modified: pkg/R/BiCopCDF.r
===================================================================
--- pkg/R/BiCopCDF.r 2015-08-04 14:09:18 UTC (rev 110)
+++ pkg/R/BiCopCDF.r 2015-08-06 15:23:55 UTC (rev 111)
@@ -1,13 +1,14 @@
-BiCopCDF <- function(u1, u2, family, par, par2 = 0, obj = NULL) {
+BiCopCDF <- function(u1, u2, family, par, par2 = 0, obj = NULL, check.pars = TRUE) {
## sanity checks for u1, u2
- if (is.null(u1) == TRUE || is.null(u2) == TRUE)
- stop("u1 and/or u2 are not set or have length zero.")
- if (any(u1 > 1) || any(u1 < 0))
- stop("Data has be in the interval [0,1].")
- if (any(u2 > 1) || any(u2 < 0))
- stop("Data has be in the interval [0,1].")
- if (length(u1) != length(u2))
- stop("Lengths of 'u1' and 'u2' do not match.")
+ if (is.null(u1) == TRUE || is.null(u2) == TRUE)
+ stop("u1 and/or u2 are not set or have length zero.")
+ if (any(u1 > 1) || any(u1 < 0))
+ stop("Data has be in the interval [0,1].")
+ if (any(u2 > 1) || any(u2 < 0))
+ stop("Data has be in the interval [0,1].")
+ if (length(u1) != length(u2))
+ stop("Lengths of 'u1' and 'u2' do not match.")
+ n <- length(u1)
## extract family and parameters if BiCop object is provided
if (missing(family))
@@ -24,26 +25,52 @@
par2 <- obj$par2
}
- ## sanity checks for family and parameters
- if (is.na(family) | is.na(par))
+ ## check for reasonable input
+ if (any(is.na(family)) | any(is.na(par)))
stop("Provide either 'family' and 'par' or 'obj'")
- if (family == 2)
+ if (any(family == 2))
stop("The CDF of the t-copula is not implemented.")
- if (!(family %in% c(0, 1, 3, 4, 5, 6, 7, 8, 9, 10, 13, 14, 16, 17, 18, 19, 20,
- 23, 24, 26, 27, 28, 29, 30, 33, 34, 36, 37, 38, 39, 40, 41,
- 51, 61, 71, 104, 114, 124, 134, 204, 214, 224, 234)))
- stop("Copula family not implemented.")
- if (family %in% c(7, 8, 9, 10, 17, 18, 19, 20, 27, 28, 29, 30, 37, 38, 39, 40,
- 104, 114, 124, 134, 204, 214, 224, 234) && par2 == 0)
- stop("For BB1, BB6, BB7, BB8 and Tawn copulas, 'par2' must be set.")
- if (family %in% c(1, 3, 4, 5, 6, 13, 14, 16, 23, 24, 26, 33, 34, 36, 41, 51,
- 61, 71) && length(par) < 1)
- stop("'par' not set.")
- BiCopCheck(family, par, par2)
- res <- rep(NA, length(u1))
+ ## adjust length for parameter vectors; stop if not matching
+ if (any(c(length(family), length(par), length(par2)) == n)) {
+ if (length(family) == 1)
+ family <- rep(family, n)
+ if (length(par) == 1)
+ par <- rep(par, n)
+ if (length(par2) == 1)
+ par2 <- rep(par2, n)
+ }
+ if (!(length(family) %in% c(1, n)))
+ stop("'family' has to be a single number or a size n vector")
+ if (!(length(par) %in% c(1, n)))
+ stop("'par' has to be a single number or a size n vector")
+ if (!(length(par2) %in% c(1, n)))
+ stop("'par2' has to be a single number or a size n vector")
- ## CDFs for the different families
+ ## check for family/parameter consistency
+ if (check.pars)
+ BiCopCheck(family, par, par2)
+
+ ## calculate CDF
+ if (length(par) == 1) {
+ # call for single parameters
+ out <- calcCDF(u1, u2, family, par, par2)
+ } else {
+ # vectorized call
+ out <- vapply(1:length(par),
+ function(i) calcCDF(u1[i],
+ u2[i],
+ family[i],
+ par[i],
+ par2[i]),
+ numeric(1))
+ }
+
+ ## return result
+ out
+}
+
+calcCDF <- function(u1, u2, family, par, par2) {
if (family == 0) {
res <- u1 * u2
} else if (family == 1) {
@@ -134,8 +161,11 @@
par3 <- par2
par2 <- 1
res <- u1 - C(u1, 1 - u2, -par, par2, par3)
+ } else {
+ ## CDFs for the different families
+ res <- rep(NA, length(u1))
}
- }
+ }
## return results
res
Modified: pkg/R/BiCopCheck.R
===================================================================
--- pkg/R/BiCopCheck.R 2015-08-04 14:09:18 UTC (rev 110)
+++ pkg/R/BiCopCheck.R 2015-08-06 15:23:55 UTC (rev 111)
@@ -1,106 +1,97 @@
-BiCopCheck <- function(family, par, par2 = 0, obj = NULL) {
- ## extract family and parameters if BiCop object is provided
- if (missing(family))
- family <- NA
- if (missing(par))
- par <- NA
- # for short hand usage extract obj from family
- if (class(family) == "BiCop")
- obj <- family
- if (!is.null(obj)) {
- stopifnot(class(obj) == "BiCop")
- family <- obj$family
- par <- obj$par
- par2 <- obj$par2
- }
-
- ## sanity checks for family and parameters
- if (is.na(family) | is.na(par))
- stop("Provide either 'family' and 'par' or 'obj'")
- if (!(family %in% c(0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 13, 14, 16, 17, 18, 19,
- 20, 23, 24, 26, 27, 28, 29, 30, 33, 34, 36, 37, 38, 39,
- 40, 41, 42, 51, 52, 61, 62, 71, 72,
- 104, 114, 124, 134, 204, 214, 224, 234)))
+BiCopCheck <- function(family, par, par2) {
+ ## check if all required parameters are set
+ if (!(all(family %in% c(0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 13, 14, 16, 17, 18, 19,
+ 20, 23, 24, 26, 27, 28, 29, 30, 33, 34, 36, 37, 38, 39,
+ 40, 41, 42, 51, 52, 61, 62, 71, 72,
+ 104, 114, 124, 134, 204, 214, 224, 234))))
stop("Copula family not implemented.")
- if (c(2, 7, 8, 9, 10, 17, 18, 19, 20, 27, 28, 29, 30, 37, 38, 39, 40, 42, 52,
- 62, 72, 104, 114, 124, 134, 204, 214, 224, 234) %in% family && par2 == 0)
+ if (any((family %in% c(2, 7, 8, 9, 10, 17, 18, 19, 20, 27, 28, 29, 30, 37, 38, 39, 40, 42, 52,
+ 62, 72, 104, 114, 124, 134, 204, 214, 224, 234)) & (par2 == 0)))
stop("For t-, BB1, BB6, BB7, BB8 and Tawn copulas, 'par2' must be set.")
- if (c(1, 3, 4, 5, 6, 11, 13, 14, 16, 23, 24, 26, 33, 34, 36, 41, 51, 61, 71) %in%
- family && length(par) < 1)
+ if (length(par) < 1)
stop("'par' not set.")
+ stopifnot(length(par) == length(par2))
- if ((family == 1 || family == 2) && abs(par[1]) >= 1)
- stop("The parameter of the Gaussian and t-copula has to be in the interval (-1,1).")
- if (family == 2 && par2 <= 2)
- stop("The degrees of freedom parameter of the t-copula has to be larger than 2.")
- if ((family == 3 || family == 13) && par <= 0)
- stop("The parameter of the Clayton copula has to be positive.")
- if ((family == 4 || family == 14) && par < 1)
- stop("The parameter of the Gumbel copula has to be in the interval [1,oo).")
- if ((family == 6 || family == 16) && par <= 1)
- stop("The parameter of the Joe copula has to be in the interval (1,oo).")
- if (family == 5 && par == 0)
- stop("The parameter of the Frank copula has to be unequal to 0.")
- if ((family == 7 || family == 17) && par <= 0)
- stop("The first parameter of the BB1 copula has to be positive.")
- if ((family == 7 || family == 17) && par2 < 1)
- stop("The second parameter of the BB1 copula has to be in the interval [1,oo).")
- if ((family == 8 || family == 18) && par <= 0)
- stop("The first parameter of the BB6 copula has to be in the interval [1,oo).")
- if ((family == 8 || family == 18) && par2 < 1)
- stop("The second parameter of the BB6 copula has to be in the interval [1,oo).")
- if ((family == 9 || family == 19) && par < 1)
- stop("The first parameter of the BB7 copula has to be in the interval [1,oo).")
- if ((family == 9 || family == 19) && par2 <= 0)
- stop("The second parameter of the BB7 copula has to be positive.")
- if ((family == 10 || family == 20) && par < 1)
- stop("The first parameter of the BB8 copula has to be in the interval [1,oo).")
- if ((family == 10 || family == 20) && (par2 <= 0 || par2 > 1))
- stop("The second parameter of the BB8 copula has to be in the interval (0,1].")
- if ((family == 23 || family == 33) && par >= 0)
- stop("The parameter of the rotated Clayton copula has to be negative.")
- if ((family == 24 || family == 34) && par > -1)
- stop("The parameter of the rotated Gumbel copula has to be in the interval (-oo,-1].")
- if ((family == 26 || family == 36) && par >= -1)
- stop("The parameter of the rotated Joe copula has to be in the interval (-oo,-1).")
- if ((family == 27 || family == 37) && par >= 0)
- stop("The first parameter of the rotated BB1 copula has to be negative.")
- if ((family == 27 || family == 37) && par2 > -1)
- stop("The second parameter of the rotated BB1 copula has to be in the interval (-oo,-1].")
- if ((family == 28 || family == 38) && par >= 0)
- stop("The first parameter of the rotated BB6 copula has to be in the interval (-oo,-1].")
- if ((family == 28 || family == 38) && par2 > -1)
- stop("The second parameter of the rotated BB6 copula has to be in the interval (-oo,-1].")
- if ((family == 29 || family == 39) && par > -1)
- stop("The first parameter of the rotated BB7 copula has to be in the interval (-oo,-1].")
- if ((family == 29 || family == 39) && par2 >= 0)
- stop("The second parameter of the rotated BB7 copula has to be negative.")
- if ((family == 30 || family == 40) && par > -1)
- stop("The first parameter of the rotated BB8 copula has to be in the interval (-oo,-1].")
- if ((family == 30 || family == 40) && (par2 >= 0 || par2 < (-1)))
- stop("The second parameter of the rotated BB8 copula has to be in the interval [-1,0).")
- if ((family == 41 || family == 51) && par <= 0)
- stop("The parameter of the reflection asymmetric copula has to be positive.")
- if ((family == 61 || family == 71) && par >= 0)
- stop("The parameter of the rotated reflection asymmetric copula has to be negative.")
- if (family == 42) {
- a <- par
- b <- par2
- limA <- (b - 3 - sqrt(9 + 6 * b - 3 * b^2))/2
- if (abs(b) > 1)
- stop("The second parameter of the two-parametric asymmetric copulas has to be in the interval [-1,1]")
- if (a > 1 || a < limA)
- stop("The first parameter of the two-parametric asymmetric copula has to be in the interval [limA(par2),1]")
+ ## check for family/parameter consistency
+ checkPars <- function(x) {
+ family <- x[1]
+ par <- x[2]
+ par2 <- x[3]
+ if ((family == 1 || family == 2) && abs(par) >= 1)
+ stop("The parameter of the Gaussian and t-copula has to be in the interval (-1,1).")
+ if (family == 2 && par2 <= 2)
+ stop("The degrees of freedom parameter of the t-copula has to be larger than 2.")
+ if ((family == 3 || family == 13) && par <= 0)
+ stop("The parameter of the Clayton copula has to be positive.")
+ if ((family == 4 || family == 14) && par < 1)
+ stop("The parameter of the Gumbel copula has to be in the interval [1,oo).")
+ if ((family == 6 || family == 16) && par <= 1)
+ stop("The parameter of the Joe copula has to be in the interval (1,oo).")
+ if (family == 5 && par == 0)
+ stop("The parameter of the Frank copula has to be unequal to 0.")
+ if ((family == 7 || family == 17) && par <= 0)
+ stop("The first parameter of the BB1 copula has to be positive.")
+ if ((family == 7 || family == 17) && par2 < 1)
+ stop("The second parameter of the BB1 copula has to be in the interval [1,oo).")
+ if ((family == 8 || family == 18) && par <= 0)
+ stop("The first parameter of the BB6 copula has to be in the interval [1,oo).")
+ if ((family == 8 || family == 18) && par2 < 1)
+ stop("The second parameter of the BB6 copula has to be in the interval [1,oo).")
+ if ((family == 9 || family == 19) && par < 1)
+ stop("The first parameter of the BB7 copula has to be in the interval [1,oo).")
+ if ((family == 9 || family == 19) && par2 <= 0)
+ stop("The second parameter of the BB7 copula has to be positive.")
+ if ((family == 10 || family == 20) && par < 1)
+ stop("The first parameter of the BB8 copula has to be in the interval [1,oo).")
+ if ((family == 10 || family == 20) && (par2 <= 0 || par2 > 1))
+ stop("The second parameter of the BB8 copula has to be in the interval (0,1].")
+ if ((family == 23 || family == 33) && par >= 0)
+ stop("The parameter of the rotated Clayton copula has to be negative.")
+ if ((family == 24 || family == 34) && par > -1)
+ stop("The parameter of the rotated Gumbel copula has to be in the interval (-oo,-1].")
+ if ((family == 26 || family == 36) && par >= -1)
+ stop("The parameter of the rotated Joe copula has to be in the interval (-oo,-1).")
+ if ((family == 27 || family == 37) && par >= 0)
+ stop("The first parameter of the rotated BB1 copula has to be negative.")
+ if ((family == 27 || family == 37) && par2 > -1)
+ stop("The second parameter of the rotated BB1 copula has to be in the interval (-oo,-1].")
+ if ((family == 28 || family == 38) && par >= 0)
+ stop("The first parameter of the rotated BB6 copula has to be in the interval (-oo,-1].")
+ if ((family == 28 || family == 38) && par2 > -1)
+ stop("The second parameter of the rotated BB6 copula has to be in the interval (-oo,-1].")
+ if ((family == 29 || family == 39) && par > -1)
+ stop("The first parameter of the rotated BB7 copula has to be in the interval (-oo,-1].")
+ if ((family == 29 || family == 39) && par2 >= 0)
+ stop("The second parameter of the rotated BB7 copula has to be negative.")
+ if ((family == 30 || family == 40) && par > -1)
+ stop("The first parameter of the rotated BB8 copula has to be in the interval (-oo,-1].")
+ if ((family == 30 || family == 40) && (par2 >= 0 || par2 < (-1)))
+ stop("The second parameter of the rotated BB8 copula has to be in the interval [-1,0).")
+ if ((family == 41 || family == 51) && par <= 0)
+ stop("The parameter of the reflection asymmetric copula has to be positive.")
+ if ((family == 61 || family == 71) && par >= 0)
+ stop("The parameter of the rotated reflection asymmetric copula has to be negative.")
+ if (family == 42) {
+ a <- par
+ b <- par2
+ limA <- (b - 3 - sqrt(9 + 6 * b - 3 * b^2))/2
+ if (abs(b) > 1)
+ stop("The second parameter of the two-parametric asymmetric copulas has to be in the interval [-1,1]")
+ if (a > 1 || a < limA)
+ stop("The first parameter of the two-parametric asymmetric copula has to be in the interval [limA(par2),1]")
+ }
+ if ((family == 104 || family == 114 || family == 204 || family == 214) && par < 1)
+ stop("Please choose 'par' of the Tawn copula in [1,oo).")
+ if ((family == 104 || family == 114 || family == 204 || family == 214) && (par2 < 0 || par2 > 1))
+ stop("Please choose 'par2' of the Tawn copula in [0,1].")
+ if ((family == 124 || family == 134 || family == 224 || family == 234) && par > -1)
+ stop("Please choose 'par' of the Tawn copula in (-oo,-1].")
+ if ((family == 124 || family == 134 || family == 224 || family == 234) && (par2 < 0 || par2 > 1))
+ stop("Please choose 'par2' of the Tawn copula in [0,1].")
}
- if ((family == 104 || family == 114 || family == 204 || family == 214) && par < 1)
- stop("Please choose 'par' of the Tawn copula in [1,oo).")
- if ((family == 104 || family == 114 || family == 204 || family == 214) && (par2 < 0 || par2 > 1))
- stop("Please choose 'par2' of the Tawn copula in [0,1].")
- if ((family == 124 || family == 134 || family == 224 || family == 234) && par > -1)
- stop("Please choose 'par' of the Tawn copula in (-oo,-1].")
- if ((family == 124 || family == 134 || family == 224 || family == 234) && (par2 < 0 || par2 > 1))
- stop("Please choose 'par2' of the Tawn copula in [0,1].")
+ apply(cbind(family, par, par2), 1, checkPars)
+
## return TRUE if all checks pass
TRUE
}
\ No newline at end of file
Modified: pkg/R/BiCopDeriv.r
===================================================================
--- pkg/R/BiCopDeriv.r 2015-08-04 14:09:18 UTC (rev 110)
+++ pkg/R/BiCopDeriv.r 2015-08-06 15:23:55 UTC (rev 111)
@@ -1,4 +1,4 @@
-BiCopDeriv <- function(u1, u2, family, par, par2 = 0, deriv = "par", log = FALSE, obj = NULL) {
+BiCopDeriv <- function(u1, u2, family, par, par2 = 0, deriv = "par", log = FALSE, obj = NULL, check.pars = TRUE) {
## sanity checks for u1, u2
if (is.null(u1) == TRUE || is.null(u2) == TRUE)
stop("u1 and/or u2 are not set or have length zero.")
@@ -8,6 +8,7 @@
stop("Data has be in the interval [0,1].")
if (any(u2 > 1) || any(u2 < 0))
stop("Data has be in the interval [0,1].")
+ n <- length(u1)
## extract family and parameters if BiCop object is provided
if (missing(family))
@@ -24,103 +25,191 @@
par2 <- obj$par2
}
- ## sanity checks for family and parameters
- if (is.na(family) | is.na(par))
+ ## check for reasonable input
+ if (any(is.na(family)) | any(is.na(par)))
stop("Provide either 'family' and 'par' or 'obj'")
- if (!(family %in% c(0, 1, 2, 3, 4, 5, 6, 13, 14, 16, 23, 24, 26, 33, 34, 36)))
+ if (!all(family %in% c(0, 1, 2, 3, 4, 5, 6, 13, 14, 16, 23, 24, 26, 33, 34, 36)))
stop("Copula family not implemented.")
- if (family == 2 && par2 == 0)
+ if (any((family == 2) & (par2 == 0)))
stop("For t-copulas, 'par2' must be set.")
- if (family %in% c(1, 3, 4, 5, 6, 13, 14, 16, 23, 24, 26, 33, 34, 36) && length(par) < 1)
- stop("'par' not set.")
- if (deriv == "par2" && family != 2)
+ if ((deriv == "par2") && any(family != 2))
stop("The derivative with respect to the second parameter can only be derived for the t-copula.")
- if (log == TRUE && (deriv %in% c("u1", "u2")))
- stop("The derivative with respect to one of the arguments are not available in the log case.")
- BiCopCheck(family, par, par2)
+ if ((log == TRUE) && (deriv %in% c("u1", "u2")))
+ stop("The derivative with respect to one of the arguments is not available in the log case.")
+ ## adjust length for parameter vectors; stop if not matching
+ if (any(c(length(family), length(par), length(par2)) == n)) {
+ if (length(family) == 1)
+ family <- rep(family, n)
+ if (length(par) == 1)
+ par <- rep(par, n)
+ if (length(par2) == 1)
+ par2 <- rep(par2, n)
+ }
+ if (!(length(family) %in% c(1, n)))
+ stop("'family' has to be a single number or a size n vector")
+ if (!(length(par) %in% c(1, n)))
+ stop("'par' has to be a single number or a size n vector")
+ if (!(length(par2) %in% c(1, n)))
+ stop("'par2' has to be a single number or a size n vector")
+
+ ## check for family/parameter consistency
+ if (check.pars)
+ BiCopCheck(family, par, par2)
+
## call C routines for specified 'deriv' case
- n <- length(u1)
- if (log == TRUE) {
- if (deriv == "par") {
- if (family == 2) {
- out <- .C("difflPDF_rho_tCopula",
+ if (length(par) == 1) {
+ ## call for single parameters
+ if (log == TRUE) {
+ if (deriv == "par") {
+ if (family == 2) {
+ out <- .C("difflPDF_rho_tCopula",
+ as.double(u1),
+ as.double(u2),
+ as.integer(n),
+ as.double(c(par, par2)),
+ as.integer(2),
+ as.double(rep(0, n)),
+ PACKAGE = "VineCopula")[[6]]
+ } else {
+ out <- .C("difflPDF_mod",
+ as.double(u1),
+ as.double(u2),
+ as.integer(n),
+ as.double(par),
+ as.integer(family),
+ as.double(rep(0, n)),
+ PACKAGE = "VineCopula")[[6]]
+ }
+ } else if (deriv == "par2") {
+ out <- .C("difflPDF_nu_tCopula_new",
+ as.double(u1),
+ as.double(u2),
+ as.integer(n),
+ as.double(c(par, par2)),
+ as.integer(2),
+ as.double(rep(0, n)), PACKAGE = "VineCopula")[[6]]
+ }
+ } else {
+ if (deriv == "par") {
+ if (family == 2) {
+ out <- .C("diffPDF_rho_tCopula",
+ as.double(u1),
+ as.double(u2),
+ as.integer(n),
+ as.double(c(par, par2)),
+ as.integer(2),
+ as.double(rep(0, n)),
+ PACKAGE = "VineCopula")[[6]]
+ } else {
+ out <- .C("diffPDF_mod",
+ as.double(u1),
+ as.double(u2),
+ as.integer(n),
+ as.double(par),
+ as.integer(family),
+ as.double(rep(0, n)),
+ PACKAGE = "VineCopula")[[6]]
+ }
+ } else if (deriv == "par2") {
+ out <- .C("diffPDF_nu_tCopula_new",
as.double(u1),
as.double(u2),
as.integer(n),
as.double(c(par, par2)),
- as.integer(2),
- as.double(rep(0, n)),
+ as.integer(2),
+ as.double(rep(0, n)),
PACKAGE = "VineCopula")[[6]]
- } else {
- out <- .C("difflPDF_mod",
+ } else if (deriv == "u1") {
+ out <- .C("diffPDF_u_mod",
as.double(u1),
as.double(u2),
as.integer(n),
- as.double(par),
+ as.double(c(par, par2)),
as.integer(family),
+ as.double(rep(0, n)),
+ PACKAGE = "VineCopula")[[6]]
+ } else if (deriv == "u2") {
+ out <- .C("diffPDF_v_mod",
+ as.double(u1),
+ as.double(u2),
+ as.integer(n),
+ as.double(c(par, par2)),
+ as.integer(family),
as.double(rep(0, n)),
PACKAGE = "VineCopula")[[6]]
+ } else {
+ stop("This kind of derivative is not implemented")
}
- } else if (deriv == "par2") {
- out <- .C("difflPDF_nu_tCopula_new",
- as.double(u1),
- as.double(u2),
- as.integer(n),
- as.double(c(par, par2)),
- as.integer(2),
- as.double(rep(0, n)), PACKAGE = "VineCopula")[[6]]
}
} else {
- if (deriv == "par") {
- if (family == 2) {
- out <- .C("diffPDF_rho_tCopula",
+ ## vectorized call
+ if (log == TRUE) {
+ if (deriv == "par") {
+ out <- .C("difflPDF_mod_vec",
as.double(u1),
as.double(u2),
- as.integer(n),
- as.double(c(par, par2)),
- as.integer(2),
+ as.integer(n),
+ as.double(par),
+ as.double(par2),
+ as.integer(family),
as.double(rep(0, n)),
- PACKAGE = "VineCopula")[[6]]
- } else {
- out <- .C("diffPDF_mod",
+ PACKAGE = "VineCopula")[[7]]
+ } else if (deriv == "par2") {
+ out <- .C("difflPDF_nu_tCopula_new_vec",
as.double(u1),
+ as.double(u2),
+ as.integer(n),
+ as.double(par),
+ as.double(par2),
+ as.integer(2),
+ as.double(rep(0, n)),
+ PACKAGE = "VineCopula")[[7]]
+ }
+ } else {
+ if (deriv == "par") {
+ out <- .C("diffPDF_mod_vec",
+ as.double(u1),
as.double(u2),
as.integer(n),
as.double(par),
+ as.double(par2),
as.integer(family),
as.double(rep(0, n)),
- PACKAGE = "VineCopula")[[6]]
+ PACKAGE = "VineCopula")[[7]]
+ } else if (deriv == "par2") {
+ out <- .C("diffPDF_nu_tCopula_new_vec",
+ as.double(u1),
+ as.double(u2),
+ as.integer(n),
+ as.double(par),
+ as.double(par2),
+ as.integer(family),
+ as.double(rep(0, n)),
+ PACKAGE = "VineCopula")[[7]]
+ } else if (deriv == "u1") {
+ out <- .C("diffPDF_u_mod_vec",
+ as.double(u1),
+ as.double(u2),
+ as.integer(n),
+ as.double(par),
+ as.double(par2),
+ as.integer(family),
+ as.double(rep(0, n)),
+ PACKAGE = "VineCopula")[[7]]
+ } else if (deriv == "u2") {
+ out <- .C("diffPDF_v_mod_vec",
+ as.double(u1),
+ as.double(u2),
+ as.integer(n),
+ as.double(par),
+ as.double(par2),
+ as.integer(family),
+ as.double(rep(0, n)),
+ PACKAGE = "VineCopula")[[7]]
+ } else {
+ stop("This kind of derivative is not implemented")
}
- } else if (deriv == "par2") {
- out <- .C("diffPDF_nu_tCopula_new",
- as.double(u1),
- as.double(u2),
- as.integer(n),
- as.double(c(par, par2)),
- as.integer(2),
- as.double(rep(0, n)),
- PACKAGE = "VineCopula")[[6]]
- } else if (deriv == "u1") {
- out <- .C("diffPDF_u_mod",
- as.double(u1),
- as.double(u2),
- as.integer(n),
- as.double(c(par, par2)),
- as.integer(family),
- as.double(rep(0, n)),
- PACKAGE = "VineCopula")[[6]]
- } else if (deriv == "u2") {
- out <- .C("diffPDF_v_mod",
- as.double(u1),
- as.double(u2),
- as.integer(n),
- as.double(c(par, par2)),
- as.integer(family),
- as.double(rep(0, n)),
- PACKAGE = "VineCopula")[[6]]
- } else {
- stop("This kind of derivative is not implemented")
}
}
Modified: pkg/R/BiCopDeriv2.r
===================================================================
--- pkg/R/BiCopDeriv2.r 2015-08-04 14:09:18 UTC (rev 110)
+++ pkg/R/BiCopDeriv2.r 2015-08-06 15:23:55 UTC (rev 111)
@@ -1,4 +1,4 @@
-BiCopDeriv2 <- function(u1, u2, family, par, par2 = 0, deriv = "par", obj = NULL) {
+BiCopDeriv2 <- function(u1, u2, family, par, par2 = 0, deriv = "par", obj = NULL, check.pars = TRUE) {
## sanity checks for u1, u2
if (is.null(u1) == TRUE || is.null(u2) == TRUE)
stop("u1 and/or u2 are not set or have length zero.")
@@ -8,6 +8,7 @@
stop("Data has be in the interval [0,1].")
if (any(u2 > 1) || any(u2 < 0))
stop("Data has be in the interval [0,1].")
+ n <- length(u1)
## extract family and parameters if BiCop object is provided
if (missing(family))
@@ -24,136 +25,253 @@
par2 <- obj$par2
}
- ## sanity checks for family and parameters
- if (is.na(family) | is.na(par))
+ ## check for reasonable input
+ if (any(is.na(family)) | any(is.na(par)))
stop("Provide either 'family' and 'par' or 'obj'")
- if (!(family %in% c(0, 1, 2, 3, 4, 5, 6, 13, 14, 16, 23, 24, 26, 33, 34, 36)))
+ if (!all(family %in% c(0, 1, 2, 3, 4, 5, 6, 13, 14, 16, 23, 24, 26, 33, 34, 36)))
stop("Copula family not implemented.")
- if (family == 2 && par2 == 0)
+ if (any((family == 2) & (par2 == 0)))
stop("For t-copulas, 'par2' must be set.")
- if (deriv == "par2" && family != 2)
+ if ((deriv %in% c("par2", "par1par2", "par2u1", "par2u2")) && any(family != 2))
stop("The derivative with respect to the second parameter can only be derived for the t-copula.")
- BiCopCheck(family, par, par2)
- ## calculate derivatives
- n <- length(u1)
- if (deriv == "par") {
- if (family == 2) {
- out <- .C("diff2PDF_rho_tCopula",
+ ## adjust length for parameter vectors; stop if not matching
+ if (any(c(length(family), length(par), length(par2)) == n)) {
+ if (length(family) == 1)
+ family <- rep(family, n)
+ if (length(par) == 1)
+ par <- rep(par, n)
+ if (length(par2) == 1)
+ par2 <- rep(par2, n)
+ }
+ if (!(length(family) %in% c(1, n)))
+ stop("'family' has to be a single number or a size n vector")
+ if (!(length(par) %in% c(1, n)))
+ stop("'par' has to be a single number or a size n vector")
+ if (!(length(par2) %in% c(1, n)))
+ stop("'par2' has to be a single number or a size n vector")
+
+ ## check for family/parameter consistency
+ if (check.pars)
+ BiCopCheck(family, par, par2)
+
+ ## call C routines for specified 'deriv' case
+ if (length(par) == 1) {
+ ## call for single parameters
+ if (deriv == "par") {
+ if (family == 2) {
+ out <- .C("diff2PDF_rho_tCopula",
+ as.double(u1),
+ as.double(u2),
+ as.integer(n),
+ as.double(c(par, par2)),
+ as.integer(2),
+ as.double(rep(0, n)),
+ PACKAGE = "VineCopula")[[6]]
+ } else {
+ out <- .C("diff2PDF_mod",
+ as.double(u1),
+ as.double(u2),
+ as.integer(n),
+ as.double(par),
+ as.integer(family),
+ as.double(rep(0, n)),
+ PACKAGE = "VineCopula")[[6]]
+ }
+ } else if (deriv == "par2") {
+ out <- .C("diff2PDF_nu_tCopula_new",
as.double(u1),
- as.double(u2),
+ as.double(u2),
as.integer(n),
- as.double(c(par, par2)),
+ as.double(c(par, par2)),
as.integer(2),
as.double(rep(0, n)),
PACKAGE = "VineCopula")[[6]]
- } else {
- out <- .C("diff2PDF_mod",
+ } else if (deriv == "u1") {
[TRUNCATED]
To get the complete diff run:
svnlook diff /svnroot/vinecopula -r 111
Mehr Informationen über die Mailingliste Vinecopula-commits