[spcopula-commits] r65 - / pkg pkg/R pkg/data pkg/demo
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Fri Oct 19 17:22:31 CEST 2012
Author: ben_graeler
Date: 2012-10-19 17:22:31 +0200 (Fri, 19 Oct 2012)
New Revision: 65
Added:
pkg/data/simulatedTriples.RData
pkg/demo/
pkg/demo/00Index
pkg/demo/MRP.R
Removed:
pkg/data/my.data.rda
spcopula_1.0.59.tar.gz
spcopula_1.0.59.zip
Modified:
pkg/DESCRIPTION
pkg/NAMESPACE
pkg/R/BB7copula.R
pkg/R/linkingCDVine.R
pkg/R/partialDerivatives.R
Log:
- Further bug fixing due to changes in copula
- demo to multivariate return periods MRP.R
Modified: pkg/DESCRIPTION
===================================================================
--- pkg/DESCRIPTION 2012-09-28 14:08:05 UTC (rev 64)
+++ pkg/DESCRIPTION 2012-10-19 15:22:31 UTC (rev 65)
@@ -1,8 +1,8 @@
Package: spcopula
Type: Package
Title: copula driven spatial analysis
-Version: 1.0.64
-Date: 2012-09-28
+Version: 1.0.65
+Date: 2012-10-19
Author: Benedikt Graeler
Maintainer: Benedikt Graeler <ben.graeler at uni-muenster.de>
Description: This package provides a framework to analyse spatial data provided in the format of the spacetime package with copulas. Additionally, support for calculating multivariate return periods is implemented.
Modified: pkg/NAMESPACE
===================================================================
--- pkg/NAMESPACE 2012-09-28 14:08:05 UTC (rev 64)
+++ pkg/NAMESPACE 2012-10-19 15:22:31 UTC (rev 65)
@@ -35,7 +35,7 @@
# MRP functions
export(genEmpKenFun, genInvKenFun)
-export(kendallRP, criticalLevel, kendallDistribution)
+export(kendallRP, criticalLevel, kendallDistribution, getKendallDistr)
## classes
exportClasses(asCopula, cqsCopula, neighbourhood)
Modified: pkg/R/BB7copula.R
===================================================================
--- pkg/R/BB7copula.R 2012-09-28 14:08:05 UTC (rev 64)
+++ pkg/R/BB7copula.R 2012-10-19 15:22:31 UTC (rev 65)
@@ -1,182 +1,245 @@
-#####################
-## ##
-## the BB7 copulas ##
-## ##
-#####################
-# Joe, H., (1997). Multivariate Models and Dependence Concepts. Monogra. Stat. Appl. Probab. 73, London: Chapman and Hall.
-
-validBB7Copula = function(object) {
- if (object at dimension != 2)
- return("Only BB7 copulas of dimension 2 are supported.")
- param <- object at parameters
- upper <- object at param.upbnd
- lower <- object at param.lowbnd
- if (length(param) != length(upper))
- return("Parameter and upper bound have non-equal length")
- if (length(param) != length(lower))
- return("Parameter and lower bound have non-equal length")
- if (any(is.na(param) | param >= upper | param[1] < lower[1] | param[2] <= lower[2]))
- return("Parameter value out of bound.")
- else return (TRUE)
-}
-
-setClass("BB7Copula",
- representation = representation("copula", family="numeric"),
- validity = validBB7Copula,
- contains = list("copula")
-)
-
-# constructor
-BB7Copula <- function (param) {
- val <- new("BB7Copula", dimension = as.integer(2), parameters = param,
- param.names = c("theta", "delta"), param.lowbnd = c(1, 0), param.upbnd = c(Inf, Inf), family=9, fullname = "BB7 copula family. Number 9 in CDVine.")
- val
-}
-
-## density ##
-setMethod("dCopula", signature("numeric","BB7Copula"), linkCDVine.PDF)
-
-## jcdf ##
-setMethod("pCopula", signature("numeric","BB7Copula"), linkCDVine.CDF)
-
-## partial derivatives ##
-# ddu
-setMethod("dduCopula", signature("numeric","BB7Copula"), linkCDVine.ddu)
-
-# ddv
-setMethod("ddvCopula", signature("numeric","BB7Copula"), linkCDVine.ddv)
-
-## random number generator
-setMethod("rCopula", signature("numeric","BB7Copula"), linkCDVine.r)
-
-## kendall distribution/measure, taken from CDVine:::obs.stat
-kendall.BB7 <- function(copula, t){
- theta = copula at parameters[1]
- delta = copula at parameters[2]
-
- kt <- rep(NA,length(t))
- kt <- t + 1/(theta * delta) * ((1 - (1 - t)^theta)^(-delta) - 1)/
- ((1 - t)^(theta - 1) * (1 - (1 - t)^theta)^(-delta - 1))
- kt[t==1] <- 1
- kt[t==0] <- 0
- return(kt)
-}
-
-setMethod("kendallDistribution", signature("BB7Copula"), kendall.BB7)
-
-setMethod("getKendallDistr", signature("BB7Copula"),
- function(copula) return(function(t) kendall.BB7(copula, t)))
-
-#########################
-## BB7 survival copula ##
-#########################
-
-setClass("surBB7Copula",
- representation = representation("copula", family="numeric"),
- validity = validBB7Copula,
- contains = list("copula")
-)
-
-# constructor
-surBB7Copula <- function (param) {
- val <- new("surBB7Copula", dimension = as.integer(2), parameters = param, param.names = c("theta", "delta"), param.lowbnd = c(1, 0), param.upbnd = c(Inf, Inf), family= 19, fullname = "Survival BB7 copula family. Number 19 in CDVine.")
- return(val)
-}
-
-## density ##
-setMethod("dCopula", signature("numeric","surBB7Copula"), linkCDVine.PDF)
-
-## jcdf ##
-setMethod("pCopula", signature("numeric","surBB7Copula"), linkCDVine.surCDF)
-# persp(surBB7Copula(c(5.329995,2.1201476)),pcopula,zlim=c(0,1))
-
-## partial derivatives ##
-# ddu
-setMethod("dduCopula", signature("numeric","surBB7Copula"), linkCDVine.ddu)
-
-# ddv
-setMethod("ddvCopula", signature("numeric","surBB7Copula"), linkCDVine.ddv)
-
-## random number generator
-setMethod("rCopula", signature("numeric","surBB7Copula"), linkCDVine.r)
-
-####################
-## BB7 copula 90� ##
-####################
-
-validRotBB7Copula = function(object) {
- if (object at dimension != 2)
- return("Only BB7 copulas of dimension 2 are supported.")
- param <- object at parameters
- upper <- object at param.upbnd
- lower <- object at param.lowbnd
- if (length(param) != length(upper))
- return("Parameter and upper bound have non-equal length")
- if (length(param) != length(lower))
- return("Parameter and lower bound have non-equal length")
- if (any(is.na(param) | param[1] > upper[1] | param[2] >= upper[2] | param <= lower))
- return("Parameter value out of bound")
- else return (TRUE)
-}
-
-setClass("r90BB7Copula",
- representation = representation("copula", family="numeric"),
- validity = validRotBB7Copula,
- contains = list("copula")
-)
-
-# constructor
-r90BB7Copula <- function (param) {
- val <- new("r90BB7Copula", dimension = as.integer(2), parameters = param, param.names = c("theta", "delta"), param.lowbnd = c(-Inf, -Inf), param.upbnd = c(-1, 0), family=29, fullname = "90 deg rotated BB7 copula family. Number 29 in CDVine.")
- val
-}
-
-## density ##
-setMethod("dCopula", signature("numeric","r90BB7Copula"), linkCDVine.PDF)
-
-## jcdf ##
-setMethod("pCopula", signature("numeric","r90BB7Copula"), linkCDVine.r90CDF)
-
-## partial derivatives ##
-# ddu
-setMethod("dduCopula", signature("numeric","r90BB7Copula"), linkCDVine.ddu)
-
-## ddv
-setMethod("ddvCopula", signature("numeric","r90BB7Copula"), linkCDVine.ddv)
-
-## random number generator
-setMethod("rCopula", signature("numeric","r90BB7Copula"), linkCDVine.r)
-# rcopula(r90BB7Copula(c(-5.329995,-1.1201476)),500)
-
-########################
-## BB7 copula 270 deg ##
-########################
-
-setClass("r270BB7Copula",
- representation = representation("copula", family="numeric"),
- validity = validRotBB7Copula,
- contains = list("copula")
-)
-
-# constructor
-r270BB7Copula <- function (param) {
- val <- new("r270BB7Copula", dimension = as.integer(2), parameters = param, param.names = c("theta", "delta"), param.lowbnd = c(-Inf, -Inf), param.upbnd = c(-1, -1), family=39, fullname = "270 deg rotated BB7 copula family. Number 39 in CDVine.")
- val
-}
-
-## density ##
-setMethod("dCopula", signature("numeric","r270BB7Copula"), linkCDVine.PDF)
-
-## jcdf ##
-setMethod("pCopula", signature("numeric","r270BB7Copula"), linkCDVine.r270CDF)
-
-## partial derivatives ##
-# ddu
-setMethod("dduCopula", signature("numeric","r270BB7Copula"), linkCDVine.ddu)
-
-# ddv
-setMethod("ddvCopula", signature("numeric","r270BB7Copula"), linkCDVine.ddv)
-
-## random number generator
-setMethod("rCopula", signature("numeric","r270BB7Copula"), linkCDVine.r)
-
+#####################
+## ##
+## the BB7 copulas ##
+## ##
+#####################
+# Joe, H., (1997). Multivariate Models and Dependence Concepts. Monogra. Stat. Appl. Probab. 73, London: Chapman and Hall.
+
+validBB7Copula = function(object) {
+ if (object at dimension != 2)
+ return("Only BB7 copulas of dimension 2 are supported.")
+ param <- object at parameters
+ upper <- object at param.upbnd
+ lower <- object at param.lowbnd
+ if (length(param) != length(upper))
+ return("Parameter and upper bound have non-equal length")
+ if (length(param) != length(lower))
+ return("Parameter and lower bound have non-equal length")
+ if (any(is.na(param) | param >= upper | param[1] < lower[1] | param[2] <= lower[2]))
+ return("Parameter value out of bound.")
+ else return (TRUE)
+}
+
+setClass("BB7Copula",
+ representation = representation("copula", family="numeric"),
+ validity = validBB7Copula,
+ contains = list("copula")
+)
+
+# constructor
+BB7Copula <- function (param) {
+ val <- new("BB7Copula", dimension = as.integer(2), parameters = param,
+ param.names = c("theta", "delta"), param.lowbnd = c(1, 0), param.upbnd = c(Inf, Inf), family=9, fullname = "BB7 copula family. Number 9 in CDVine.")
+ val
+}
+
+## density ##
+setMethod("dCopula", signature("numeric","BB7Copula"),
+ function(u, copula, log) {
+ linkCDVine.PDF(matrix(u,ncol=copula at dimension),copula, log)
+ })
+setMethod("dCopula", signature("matrix","BB7Copula"), function(u, copula, log) linkCDVine.PDF(u, copula, log))
+
+## jcdf ##
+setMethod("pCopula", signature("numeric","BB7Copula"),
+ function(u, copula, ...) {
+ linkCDVine.CDF(matrix(u,ncol=copula at dimension),copula)
+ })
+setMethod("pCopula", signature("matrix","BB7Copula"), linkCDVine.CDF)
+
+## partial derivatives ##
+# ddu
+setMethod("dduCopula", signature("numeric","BB7Copula"),
+ function(u, copula, ...) {
+ linkCDVine.ddu(matrix(u,ncol=copula at dimension),copula)
+ })
+setMethod("dduCopula", signature("matrix","BB7Copula"), linkCDVine.ddu)
+
+# ddv
+setMethod("ddvCopula", signature("numeric","BB7Copula"),
+ function(u, copula, ...) {
+ linkCDVine.ddv(matrix(u,ncol=copula at dimension),copula)
+ })
+setMethod("ddvCopula", signature("matrix","BB7Copula"), linkCDVine.ddv)
+
+## random number generator
+setMethod("rCopula", signature("numeric","BB7Copula"), linkCDVine.r)
+
+## kendall distribution/measure, taken from CDVine:::obs.stat
+kendall.BB7 <- function(copula, t){
+ theta = copula at parameters[1]
+ delta = copula at parameters[2]
+
+ kt <- rep(NA,length(t))
+ kt <- t + 1/(theta * delta) * ((1 - (1 - t)^theta)^(-delta) - 1)/
+ ((1 - t)^(theta - 1) * (1 - (1 - t)^theta)^(-delta - 1))
+ kt[t==1] <- 1
+ kt[t==0] <- 0
+ return(kt)
+}
+
+setMethod("kendallDistribution", signature("BB7Copula"), kendall.BB7)
+
+setMethod("getKendallDistr", signature("BB7Copula"),
+ function(copula) return(function(t) kendall.BB7(copula, t)))
+
+#########################
+## BB7 survival copula ##
+#########################
+
+setClass("surBB7Copula",
+ representation = representation("copula", family="numeric"),
+ validity = validBB7Copula,
+ contains = list("copula")
+)
+
+# constructor
+surBB7Copula <- function (param) {
+ val <- new("surBB7Copula", dimension = as.integer(2), parameters = param, param.names = c("theta", "delta"), param.lowbnd = c(1, 0), param.upbnd = c(Inf, Inf), family= 19, fullname = "Survival BB7 copula family. Number 19 in CDVine.")
+ return(val)
+}
+
+## density ##
+setMethod("dCopula", signature("numeric","surBB7Copula"),
+ function(u, copula, ...) {
+ linkCDVine.PDF(matrix(u,ncol=copula at dimension),copula)
+ })
+setMethod("dCopula", signature("matrix","surBB7Copula"), linkCDVine.PDF)
+
+## jcdf ##
+setMethod("pCopula", signature("numeric","surBB7Copula"),
+ function(u, copula, ...) {
+ linkCDVine.CDF(matrix(u,ncol=copula at dimension),copula)
+ })
+setMethod("pCopula", signature("matrix","surBB7Copula"), linkCDVine.surCDF)
+
+## partial derivatives ##
+# ddu
+setMethod("dduCopula", signature("numeric","surBB7Copula"),
+ function(u, copula, ...) {
+ linkCDVine.ddu(matrix(u,ncol=copula at dimension),copula)
+ })
+setMethod("dduCopula", signature("matrix","surBB7Copula"), linkCDVine.ddu)
+
+# ddv
+setMethod("ddvCopula", signature("numeric","surBB7Copula"),
+ function(u, copula, ...) {
+ linkCDVine.ddv(matrix(u,ncol=copula at dimension),copula)
+ })
+setMethod("ddvCopula", signature("matrix","surBB7Copula"), linkCDVine.ddv)
+
+## random number generator
+setMethod("rCopula", signature("numeric","surBB7Copula"), linkCDVine.r)
+
+###################
+## BB7 copula 90 ##
+###################
+
+validRotBB7Copula = function(object) {
+ if (object at dimension != 2)
+ return("Only BB7 copulas of dimension 2 are supported.")
+ param <- object at parameters
+ upper <- object at param.upbnd
+ lower <- object at param.lowbnd
+ if (length(param) != length(upper))
+ return("Parameter and upper bound have non-equal length")
+ if (length(param) != length(lower))
+ return("Parameter and lower bound have non-equal length")
+ if (any(is.na(param) | param[1] > upper[1] | param[2] >= upper[2] | param <= lower))
+ return("Parameter value out of bound")
+ else return (TRUE)
+}
+
+setClass("r90BB7Copula",
+ representation = representation("copula", family="numeric"),
+ validity = validRotBB7Copula,
+ contains = list("copula")
+)
+
+# constructor
+r90BB7Copula <- function (param) {
+ val <- new("r90BB7Copula", dimension = as.integer(2), parameters = param, param.names = c("theta", "delta"), param.lowbnd = c(-Inf, -Inf), param.upbnd = c(-1, 0), family=29, fullname = "90 deg rotated BB7 copula family. Number 29 in CDVine.")
+ val
+}
+
+## density ##
+setMethod("dCopula", signature("numeric","r90BB7Copula"),
+ function(u, copula, ...) {
+ linkCDVine.PDF(matrix(u,ncol=copula at dimension),copula)
+ })
+setMethod("dCopula", signature("matrix","r90BB7Copula"), linkCDVine.PDF)
+
+## jcdf ##
+setMethod("pCopula", signature("numeric","r90BB7Copula"),
+ function(u, copula, ...) {
+ linkCDVine.CDF(matrix(u,ncol=copula at dimension),copula)
+ })
+setMethod("pCopula", signature("matrix","r90BB7Copula"), linkCDVine.r90CDF)
+
+## partial derivatives ##
+# ddu
+setMethod("dduCopula", signature("numeric","r90BB7Copula"),
+ function(u, copula, ...) {
+ linkCDVine.ddu(matrix(u,ncol=copula at dimension),copula)
+ })
+setMethod("dduCopula", signature("matrix","r90BB7Copula"), linkCDVine.ddu)
+
+## ddv
+setMethod("ddvCopula", signature("numeric","r90BB7Copula"),
+ function(u, copula, ...) {
+ linkCDVine.ddv(matrix(u,ncol=copula at dimension),copula)
+ })
+setMethod("ddvCopula", signature("matrix","r90BB7Copula"), linkCDVine.ddv)
+
+## random number generator
+setMethod("rCopula", signature("numeric","r90BB7Copula"), linkCDVine.r)
+# rcopula(r90BB7Copula(c(-5.329995,-1.1201476)),500)
+
+########################
+## BB7 copula 270 deg ##
+########################
+
+setClass("r270BB7Copula",
+ representation = representation("copula", family="numeric"),
+ validity = validRotBB7Copula,
+ contains = list("copula")
+)
+
+# constructor
+r270BB7Copula <- function (param) {
+ val <- new("r270BB7Copula", dimension = as.integer(2), parameters = param, param.names = c("theta", "delta"), param.lowbnd = c(-Inf, -Inf), param.upbnd = c(-1, -1), family=39, fullname = "270 deg rotated BB7 copula family. Number 39 in CDVine.")
+ val
+}
+
+## density ##
+setMethod("dCopula", signature("numeric","r270BB7Copula"),
+ function(u, copula, ...) {
+ linkCDVine.PDF(matrix(u,ncol=copula at dimension),copula)
+ })
+setMethod("dCopula", signature("matrix","r270BB7Copula"), linkCDVine.PDF)
+
+## jcdf ##
+setMethod("pCopula", signature("numeric","r270BB7Copula"),
+ function(u, copula, ...) {
+ linkCDVine.CDF(matrix(u,ncol=copula at dimension),copula)
+ })
+setMethod("pCopula", signature("matrix","r270BB7Copula"), linkCDVine.r270CDF)
+
+## partial derivatives ##
+# ddu
+setMethod("dduCopula", signature("numeric","r270BB7Copula"),
+ function(u, copula, ...) {
+ linkCDVine.ddu(matrix(u,ncol=copula at dimension),copula)
+ })
+setMethod("dduCopula", signature("matrix","r270BB7Copula"), linkCDVine.ddu)
+
+# ddv
+setMethod("ddvCopula", signature("numeric","r270BB7Copula"),
+ function(u, copula, ...) {
+ linkCDVine.ddv(matrix(u,ncol=copula at dimension),copula)
+ })
+setMethod("ddvCopula", signature("matrix","r270BB7Copula"), linkCDVine.ddv)
+
+## random number generator
+setMethod("rCopula", signature("numeric","r270BB7Copula"), linkCDVine.r)
+
Modified: pkg/R/linkingCDVine.R
===================================================================
--- pkg/R/linkingCDVine.R 2012-09-28 14:08:05 UTC (rev 64)
+++ pkg/R/linkingCDVine.R 2012-10-19 15:22:31 UTC (rev 65)
@@ -4,17 +4,17 @@
# density from BiCopPDF
-linkCDVine.PDF <- function (u, copula) {
+linkCDVine.PDF <- function (u, copula, log) {
param <- copula at parameters
if(length(param)==1) param <- c(param,0)
- if (!is.matrix(u)) u <- matrix(u, ncol = 2)
n <- nrow(u)
fam <- copula at family
coplik = .C("LL_mod_seperate", as.integer(fam), as.integer(n), as.double(u[,1]),
as.double(u[,2]), as.double(param[1]), as.double(param[2]),
as.double(rep(0, n)), PACKAGE = "CDVine")[[7]]
- return(exp(coplik))
+ if(log) return(coplik)
+ else return(exp(coplik))
}
# cdf from BiCopCDF
Modified: pkg/R/partialDerivatives.R
===================================================================
--- pkg/R/partialDerivatives.R 2012-09-28 14:08:05 UTC (rev 64)
+++ pkg/R/partialDerivatives.R 2012-10-19 15:22:31 UTC (rev 65)
@@ -33,8 +33,8 @@
warning("Numerical evaluation of invddu takes place.")
res <- NULL
for (i in 1:length(u)) {
- res <- rbind(res, optimize(function(x) (dduCopula(copula,
- cbind(rep(u[i], length(x)), x)) - y[i])^2,
+ res <- rbind(res, optimize( function(x)
+ (dduCopula(cbind(rep(u[i], length(x)), x),copula) - y[i])^2,
interval = c(0, 1))$minimum)
}
return(res)
@@ -49,8 +49,8 @@
warning("Numerical evaluation of invddv takes place.")
res <- NULL
for (i in 1:length(v)) {
- res <- rbind(res, optimize(function(x) (ddvCopula(copula,
- cbind(x, rep(v[i], length(x)))) - y[i])^2,
+ res <- rbind(res, optimize(function(x)
+ (ddvCopula(cbind(x, rep(v[i], length(x))),copula) - y[i])^2,
interval = c(0, 1))$minimum)
}
return(res)
Deleted: pkg/data/my.data.rda
===================================================================
(Binary files differ)
Added: pkg/data/simulatedTriples.RData
===================================================================
(Binary files differ)
Property changes on: pkg/data/simulatedTriples.RData
___________________________________________________________________
Added: svn:mime-type
+ application/octet-stream
Added: pkg/demo/00Index
===================================================================
--- pkg/demo/00Index (rev 0)
+++ pkg/demo/00Index 2012-10-19 15:22:31 UTC (rev 65)
@@ -0,0 +1 @@
+MRP The MRP demo gives insight in the code used in the paper multivariate return periods, a practical ...
Added: pkg/demo/MRP.R
===================================================================
--- pkg/demo/MRP.R (rev 0)
+++ pkg/demo/MRP.R 2012-10-19 15:22:31 UTC (rev 65)
@@ -0,0 +1,46 @@
+## get the data
+data(simulatedTriples)
+
+## rank order transformation
+peakVol <- rankTransform(triples[,1],triples[,3])
+colnames(peakVol) <- c("Qp","Vp")
+plot(peakVol, asp=1)
+
+# Kendall's tau correlation
+cor(triples,method="kendall")
+
+sum(log(dCopula(peakVol,BB7Copula(c(2,14)))))
+
+sum(dCopula(peakVol,BB7Copula(c(2,14)),log=T))
+sum(dCopula(peakVol,BB7Copula(c(2,14)),log=F))
+
+class(peakVol)
+
+loglikCopula(c(2,14),x=peakVol, BB7Copula(c(2,14)))
+
+# estiamte the BB7 copula by means of maximum likelihood
+copQV <- fitCopula(BB7Copula(param=c(2,14)), peakVol, method="ml",start=c(2,14), estimate.variance=F)@copula
+copQV
+
+# we use a design return period of 100 years
+# the MAR-case
+v_MAR <- c(0.99,invdduCopula(0.99,copQV,0.99))
+v_MAR
+
+## the anlytical kendall distribution
+kendallFunQV <- getKendallDistr(copQV)
+
+# the Kendall distribution value (fraction of pairs having a smaller copula value than "t")
+kendallFunQV(t=0.99)
+
+curve(kendallFunQV, from=0, to=1, asp=1)
+
+# the critical level of the KEN2-RP for 100 years
+t_KEN2 <- criticalLevel(kendallFunQV,100,mu=1)
+t_KEN2
+
+# the corresponmding KEN2-RP for the OR-RP
+kendallRP(kendallFunQV,cl=0.99,mu=1,copula=copQV)
+
+# illustrating the critical lines (empirically)
+contour(copQV,pCopula,levels=c(0.99,t_KEN2),xlim=c(0.98,1),ylim=c(0.98,1),n=1000, asp=1, col="blue")
\ No newline at end of file
Deleted: spcopula_1.0.59.tar.gz
===================================================================
(Binary files differ)
Deleted: spcopula_1.0.59.zip
===================================================================
(Binary files differ)
More information about the spcopula-commits
mailing list