[Distr-commits] r809 - in branches/distr-2.4/pkg/distrMod: . R inst inst/scripts man vignettes
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Sat May 19 18:16:12 CEST 2012
Author: ruckdeschel
Date: 2012-05-19 18:16:12 +0200 (Sat, 19 May 2012)
New Revision: 809
Removed:
branches/distr-2.4/pkg/distrMod/R/GParetoFamily.R
branches/distr-2.4/pkg/distrMod/R/LDEstimator.R
branches/distr-2.4/pkg/distrMod/man/GParetoFamily.Rd
branches/distr-2.4/pkg/distrMod/man/GumbelLocationFamily.Rd
branches/distr-2.4/pkg/distrMod/man/LDEstimator.Rd
branches/distr-2.4/pkg/distrMod/man/internalldeHelpers.Rd
Modified:
branches/distr-2.4/pkg/distrMod/NAMESPACE
branches/distr-2.4/pkg/distrMod/R/0distrModOptions.R
branches/distr-2.4/pkg/distrMod/R/0distrModUtils.R
branches/distr-2.4/pkg/distrMod/R/AllClass.R
branches/distr-2.4/pkg/distrMod/R/AllGeneric.R
branches/distr-2.4/pkg/distrMod/R/AllReturnClasses.R
branches/distr-2.4/pkg/distrMod/R/AllShow.R
branches/distr-2.4/pkg/distrMod/R/L2GroupFamilies-methods.R
branches/distr-2.4/pkg/distrMod/R/MLEstimator.R
branches/distr-2.4/pkg/distrMod/R/SimpleL2ParamFamilies.R
branches/distr-2.4/pkg/distrMod/R/confint.R
branches/distr-2.4/pkg/distrMod/inst/NEWS
branches/distr-2.4/pkg/distrMod/inst/scripts/example_CvMMDE.R
branches/distr-2.4/pkg/distrMod/man/0distrMod-package.Rd
branches/distr-2.4/pkg/distrMod/man/InternalReturnClasses-class.Rd
branches/distr-2.4/pkg/distrMod/man/L2LocationFamily-class.Rd
branches/distr-2.4/pkg/distrMod/man/L2LocationFamily.Rd
branches/distr-2.4/pkg/distrMod/man/MDEstimator.Rd
branches/distr-2.4/pkg/distrMod/man/ParamFamily.Rd
branches/distr-2.4/pkg/distrMod/man/internalClassUnions-class.Rd
branches/distr-2.4/pkg/distrMod/man/validParameter-methods.Rd
branches/distr-2.4/pkg/distrMod/vignettes/distrMod.Rnw
Log:
distrMod: started moving functionality for extreme value distributions (in this case GumbelLocationFamily and temporarily stored LDEstimators) from package distrMod to new package RobExtremes developed in robast family on r-forge
Modified: branches/distr-2.4/pkg/distrMod/NAMESPACE
===================================================================
--- branches/distr-2.4/pkg/distrMod/NAMESPACE 2012-05-19 16:07:43 UTC (rev 808)
+++ branches/distr-2.4/pkg/distrMod/NAMESPACE 2012-05-19 16:16:12 UTC (rev 809)
@@ -22,10 +22,11 @@
exportClasses("L2GroupParamFamily", "L2LocationFamily",
"L2ScaleFamily", "L2LocationScaleFamily")
exportClasses("L2LocationScaleUnion")
+exportClasses("L2ScaleShapeUnion")
exportClasses("BinomFamily","PoisFamily", "NormLocationFamily",
- "GumbelLocationFamily", "NormScaleFamily", "ExpScaleFamily",
+ "NormScaleFamily", "ExpScaleFamily",
"LnormScaleFamily", "GammaFamily", "BetaFamily", "NormLocationScaleFamily",
- "CauchyLocationScaleFamily", "GParetoFamily")
+ "CauchyLocationScaleFamily")
exportClasses("NormType", "QFNorm", "InfoNorm", "SelfNorm")
exportClasses("Estimate", "MCEstimate")
exportClasses("Confint")
@@ -65,10 +66,9 @@
export("NonSymmetric", "EvenSymmetric", "OddSymmetric", "FunSymmList")
export("ParamFamParameter", "ParamFamily", "L2ParamFamily",
"BinomFamily", "PoisFamily", "NbinomFamily", "NormLocationFamily",
- "GumbelLocationFamily", "NormScaleFamily", "ExpScaleFamily",
+ "NormScaleFamily", "ExpScaleFamily",
"LnormScaleFamily", "GammaFamily", "BetaFamily", "NormLocationScaleFamily",
- "CauchyLocationScaleFamily", "NbinomwithSizeFamily", "NbinomMeanSizeFamily",
- "GParetoFamily")
+ "CauchyLocationScaleFamily", "NbinomwithSizeFamily", "NbinomMeanSizeFamily")
export("asCov", "trAsCov", "asHampel", "asBias", "asMSE", "asUnOvShoot",
"fiCov", "trFiCov", "fiHampel", "fiMSE", "fiBias", "fiUnOvShoot")
export("positiveBias", "negativeBias", "symmetricBias",
@@ -80,4 +80,3 @@
export("NormScaleUnknownLocationFamily", "NormLocationUnknownScaleFamily")
export("L2LocationUnknownScaleFamily", "L2ScaleUnknownLocationFamily")
export("meRes", "get.criterion.fct")
-export("LDEstimator", "medkMAD", "medSn", "medQn", "medkMADhybr")
Modified: branches/distr-2.4/pkg/distrMod/R/0distrModOptions.R
===================================================================
--- branches/distr-2.4/pkg/distrMod/R/0distrModOptions.R 2012-05-19 16:07:43 UTC (rev 808)
+++ branches/distr-2.4/pkg/distrMod/R/0distrModOptions.R 2012-05-19 16:16:12 UTC (rev 809)
@@ -25,4 +25,4 @@
}
getdistrModOption <- function(x) distrModOptions(x)[[1]]
-distrModoptions <- distrModOptions
\ No newline at end of file
+distrModoptions <- distrModOptions
Modified: branches/distr-2.4/pkg/distrMod/R/0distrModUtils.R
===================================================================
--- branches/distr-2.4/pkg/distrMod/R/0distrModUtils.R 2012-05-19 16:07:43 UTC (rev 808)
+++ branches/distr-2.4/pkg/distrMod/R/0distrModUtils.R 2012-05-19 16:16:12 UTC (rev 809)
@@ -18,7 +18,7 @@
upp1 <- me + IQR.fac * s1
low <- max(low0,low1); upp <- min(upp0, upp1)
xs <- seq(low, upp, length = getdistrOption("DefaultNrGridPoints"))
- m <- getdistrOption("DefaultNrGridPoints")%/%100
+ m <- getdistrOption("DefaultNrGridPoints")%/%100+1
dxs<- -d(distr)(xs, log = TRUE)
# plot(xs, dxs,type="l")
x1 <- xs[1]; xn <- (rev(xs)[1])
Modified: branches/distr-2.4/pkg/distrMod/R/AllClass.R
===================================================================
--- branches/distr-2.4/pkg/distrMod/R/AllClass.R 2012-05-19 16:07:43 UTC (rev 808)
+++ branches/distr-2.4/pkg/distrMod/R/AllClass.R 2012-05-19 16:16:12 UTC (rev 809)
@@ -218,6 +218,11 @@
contains = c("L2GroupParamFamily","VIRTUAL")
)
+## virtual in-between class for common parts in modifyModel - method
+setClass("L2ScaleShapeUnion",
+ representation(withPos = "logical"),
+ contains = c("L2GroupParamFamily","VIRTUAL")
+ )
## L2-differentiable (univariate) location family
setClass("L2LocationFamily",
Modified: branches/distr-2.4/pkg/distrMod/R/AllGeneric.R
===================================================================
--- branches/distr-2.4/pkg/distrMod/R/AllGeneric.R 2012-05-19 16:07:43 UTC (rev 808)
+++ branches/distr-2.4/pkg/distrMod/R/AllGeneric.R 2012-05-19 16:16:12 UTC (rev 809)
@@ -278,3 +278,9 @@
if(!isGeneric("optimwarn")){
setGeneric("optimwarn", function(object) standardGeneric("optimwarn"))
}
+if(!isGeneric("withPos")){
+ setGeneric("withPos", function(object) standardGeneric("withPos"))
+}
+if(!isGeneric("withPos<-")){
+ setGeneric("withPos<-", function(object,value) standardGeneric("withPos<-"))
+}
Modified: branches/distr-2.4/pkg/distrMod/R/AllReturnClasses.R
===================================================================
--- branches/distr-2.4/pkg/distrMod/R/AllReturnClasses.R 2012-05-19 16:07:43 UTC (rev 808)
+++ branches/distr-2.4/pkg/distrMod/R/AllReturnClasses.R 2012-05-19 16:16:12 UTC (rev 809)
@@ -28,8 +28,8 @@
contains = "L2ParamFamily")
## Gamma family
-setClass("GammaFamily",
- contains = "L2ParamFamily")
+setClass("GammaFamily", prototype=prototype(withPos=TRUE),
+ contains = "L2ScaleShapeUnion")
## Beta family
setClass("BetaFamily",
@@ -39,10 +39,6 @@
setClass("NormLocationFamily",
contains = "L2LocationFamily")
-## Gumbel location family
-setClass("GumbelLocationFamily",
- contains = "L2LocationFamily")
-
## Normal scale family
setClass("NormScaleFamily",
contains = "L2ScaleFamily")
@@ -62,7 +58,5 @@
## Cauchy location scale family
setClass("CauchyLocationScaleFamily",
contains = "L2LocationScaleFamily")
-## class
-setClass("GParetoFamily", contains="L2ParamFamily")
Modified: branches/distr-2.4/pkg/distrMod/R/AllShow.R
===================================================================
--- branches/distr-2.4/pkg/distrMod/R/AllShow.R 2012-05-19 16:07:43 UTC (rev 808)
+++ branches/distr-2.4/pkg/distrMod/R/AllShow.R 2012-05-19 16:16:12 UTC (rev 809)
@@ -338,4 +338,4 @@
# options("digits" = digits)
# show(object = x)
# options("digits" = oldDigits)
-# })
\ No newline at end of file
+# })
Deleted: branches/distr-2.4/pkg/distrMod/R/GParetoFamily.R
===================================================================
--- branches/distr-2.4/pkg/distrMod/R/GParetoFamily.R 2012-05-19 16:07:43 UTC (rev 808)
+++ branches/distr-2.4/pkg/distrMod/R/GParetoFamily.R 2012-05-19 16:16:12 UTC (rev 809)
@@ -1,362 +0,0 @@
-#################################
-##
-## Class: GParetoFamily
-##
-################################
-
-
-## methods
-setMethod("validParameter",signature(object="GParetoFamily"),
- function(object, param, tol =.Machine$double.eps){
- if (is(param, "ParamFamParameter"))
- param <- main(param)
- if (!all(is.finite(param)))
- return(FALSE)
- if (any(param[1] <= tol))
- return(FALSE)
- if (any(param[2] <= tol))
- return(FALSE)
- return(TRUE)
- })
-
-
-## generating function
-## loc: known/fixed threshold/location parameter
-## scale: scale parameter
-## shape: shape parameter
-## of.interest: which parameters, transformations are of interest
-## posibilites are: scale, shape, quantile, expected loss, expected shortfall
-## a maximum number of two of these may be selected
-## p: probability needed for quantile and expected shortfall
-## N: expected frequency for expected loss
-## trafo: optional parameter transformation
-## start0Est: startEstimator for MLE and MDE --- if NULL HybridEstimator is used;
-
-GParetoFamily <- function(loc = 0, scale = 1, shape = 0.5,
- of.interest = c("scale", "shape"),
- p = NULL, N = NULL, trafo = NULL,
- start0Est = NULL){
- if(is.null(trafo)){
- of.interest <- unique(of.interest)
- if(length(of.interest) > 2)
- stop("A maximum number of two parameters resp. parameter transformations may be selected.")
- if(!all(of.interest %in% c("scale", "shape", "quantile", "expected loss", "expected shortfall")))
- stop("Parameters resp. transformations of interest have to be selected from: ",
- "'scale', 'shape', 'quantile', 'expected loss', 'expected shortfall'.")
-
- ## reordering of of.interest
- if(("scale" %in% of.interest) && ("scale" != of.interest[1])){
- of.interest[2] <- of.interest[1]
- of.interest[1] <- "scale"
- }
- if(!("scale" %in% of.interest) && ("shape" %in% of.interest) && ("shape" != of.interest[1])){
- of.interest[2] <- of.interest[1]
- of.interest[1] <- "shape"
- }
- if(!any(c("scale", "shape") %in% of.interest) && ("quantile" %in% of.interest)
- && ("quantile" != of.interest[1])){
- of.interest[2] <- of.interest[1]
- of.interest[1] <- "quantile"
- }
- if(!any(c("scale", "shape", "quantile") %in% of.interest)
- && ("expected shortfall" %in% of.interest)
- && ("expected shortfall" != of.interest[1])){
- of.interest[2] <- of.interest[1]
- of.interest[1] <- "expected shortfall"
- }
- }
- theta <- c(loc, scale, shape)
-
- ##symmetry
- distrSymm <- NoSymmetry()
-
- ## parameters
- names(theta) <- c("loc", "scale", "shape")
-
- if(is.null(trafo)){
- tau <- NULL
- if("scale" %in% of.interest){
- tau <- function(theta){
- th <- theta[1]
- names(th) <- "scale"
- th
- }
- Dtau <- function(theta){
- D <- t(c(1, 0))
- rownames(D) <- "scale"
- D
- }
- }
- if("shape" %in% of.interest){
- if(is.null(tau)){
- tau <- function(theta){
- th <- theta[2]
- names(th) <- "shape"
- th
- }
- Dtau <- function(theta){
- D <- t(c(0, 1))
- rownames(D) <- "shape"
- D
- }
- }else{
- tau <- function(theta){
- th <- theta
- names(th) <- c("scale", "shape")
- th
- }
- Dtau <- function(theta){
- D <- diag(2)
- rownames(D) <- c("scale", "shape")
- D
- }
- }
- }
- if("quantile" %in% of.interest){
- if(is.null(p)) stop("Probability 'p' has to be specified.")
- if(is.null(tau)){
- tau <- function(theta){ }
- body(tau) <- substitute({ q <- loc0 + theta[1]*((1-p0)^(-theta[2])-1)/theta[2]
- names(q) <- "quantile"
- q },
- list(loc0 = loc, p0 = p))
- Dtau <- function(theta){ }
- body(Dtau) <- substitute({ scale <- theta[1]
- shape <- theta[2]
- D1 <- ((1-p0)^(-shape)-1)/shape
- D2 <- -scale/shape*(D1 + log(1-p0)*(1-p0)^(-shape))
- D <- t(c(D1, D2))
- rownames(D) <- "quantile"
- colnames(D) <- NULL
- D },
- list(p0 = p))
- }else{
- tau1 <- tau
- tau <- function(theta){ }
- body(tau) <- substitute({ q <- loc0 + theta[1]*((1-p0)^(-theta[2])-1)/theta[2]
- names(q) <- "quantile"
- c(tau0(theta), q) },
- list(tau0 = tau1, loc0 = loc, p0 = p))
- Dtau1 <- Dtau
- Dtau <- function(theta){}
- body(Dtau) <- substitute({ scale <- theta[1]
- shape <- theta[2]
- D1 <- ((1-p0)^(-shape)-1)/shape
- D2 <- -scale/shape*(D1 + log(1-p0)*(1-p0)^(-shape))
- D <- t(c(D1, D2))
- rownames(D) <- "quantile"
- colnames(D) <- NULL
- rbind(Dtau0(theta), D) },
- list(Dtau0 = Dtau1, p0 = p))
- }
- }
- if("expected shortfall" %in% of.interest){
- if(is.null(p)) stop("Probability 'p' has to be specified.")
- if(is.null(tau)){
- tau <- function(theta){ }
- body(tau) <- substitute({ q <- loc0 + theta[1]*((1-p0)^(-theta[2])-1)/theta[2]
- es <- (q + theta[1] - theta[2]*loc0)/(1-theta[2])
- names(es) <- "expected shortfall"
- es },
- list(loc0 = loc, p0 = p))
- Dtau <- function(theta){ }
- body(Dtau) <- substitute({ scale <- theta[1]
- shape <- theta[2]
- q <- loc0 + theta[1]*((1-p0)^(-theta[2])-1)/theta[2]
- dq1 <- ((1-p0)^(-shape)-1)/shape
- dq2 <- -scale/shape*(dq1 + log(1-p0)*(1-p0)^(-shape))
- D1 <- (dq1 + 1)/(1-shape)
- D2 <- (dq2 - loc0)/(1-shape) + (q + scale - loc0*shape)/(1-shape)^2
- D <- t(c(D1, D2))
- rownames(D) <- "expected shortfall"
- colnames(D) <- NULL
- D },
- list(loc0 = loc, p0 = p))
- }else{
- tau1 <- tau
- tau <- function(theta){ }
- body(tau) <- substitute({ q <- loc0 + theta[1]*((1-p0)^(-theta[2])-1)/theta[2]
- es <- (q + theta[1] - theta[2]*loc0)/(1-theta[2])
- names(es) <- "expected shortfall"
- c(tau0(theta), es) },
- list(tau0 = tau1, loc0 = loc, p0 = p))
- Dtau1 <- Dtau
- Dtau <- function(theta){}
- body(Dtau) <- substitute({ scale <- theta[1]
- shape <- theta[2]
- q <- loc0 + theta[1]*((1-p0)^(-theta[2])-1)/theta[2]
- dq1 <- ((1-p0)^(-shape)-1)/shape
- dq2 <- -scale/shape*(dq1 + log(1-p0)*(1-p0)^(-shape))
- D1 <- (dq1 + 1)/(1-shape)
- D2 <- (dq2 - loc0)/(1-shape) + (q + scale - loc0*shape)/(1-shape)^2
- D <- t(c(D1, D2))
- rownames(D) <- "expected shortfall"
- colnames(D) <- NULL
- rbind(Dtau0(theta), D) },
- list(Dtau0 = Dtau1, loc0 = loc, p0 = p))
- }
- }
- if("expected loss" %in% of.interest){
- if(is.null(N)) stop("Expected frequency 'N' has to be specified.")
- if(is.null(tau)){
- tau <- function(theta){ }
- body(tau) <- substitute({ el <- N0*(loc0 + theta[1]*gamma(1/theta[2]-1)/(theta[2]^2*gamma(1/theta[2]+1)))
- names(el) <- "expected loss"
- el },
- list(loc0 = loc,N0 = N))
- Dtau <- function(theta){ }
- body(Dtau) <- substitute({ scale <- theta[1]
- shape <- theta[2]
- Gneg <- gamma(1/shape-1)
- Gpos <- gamma(1/shape+1)
- D1 <- N0*Gneg/(shape^2*Gpos)
- D2 <- N0*scale*Gneg*(digamma(1/shape+1) - 2*shape - digamma(1/shape-1))/(shape^4*Gpos)
- D <- t(c(D1, D2))
- rownames(D) <- "expected loss"
- colnames(D) <- NULL
- D },
- list(loc0 = loc, N0 = N))
- }else{
- tau1 <- tau
- tau <- function(theta){ }
- body(tau) <- substitute({ el <- N0*(loc0 + theta[1]*gamma(1/theta[2]-1)/(theta[2]^2*gamma(1/theta[2]+1)))
- names(el) <- "expected loss"
- c(tau0(theta), el) },
- list(tau0 = tau1, loc0 = loc,N0 = N))
- Dtau1 <- Dtau
- Dtau <- function(theta){}
- body(Dtau) <- substitute({ scale <- theta[1]
- shape <- theta[2]
- Gneg <- gamma(1/shape-1)
- Gpos <- gamma(1/shape+1)
- D1 <- N0*Gneg/(shape^2*Gpos)
- D2 <- N0*scale*Gneg*(digamma(1/shape+1) - 2*shape - digamma(1/shape-1))/(shape^4*Gpos)
- D <- t(c(D1, D2))
- rownames(D) <- "expected loss"
- colnames(D) <- NULL
- rbind(Dtau0(theta), D) },
- list(Dtau0 = Dtau1, loc0 = loc, N0 = N))
- }
- }
- trafo <- function(x){ list(fval = tau(x), mat = Dtau(x)) }
- }else{
- if(is.matrix(trafo) & nrow(trafo) > 2) stop("number of rows of 'trafo' > 2")
- }
- param <- ParamFamParameter(name = "theta", main = theta[2:3],
- fixed = theta[1],
- trafo = trafo)
-
- ## distribution
- distribution <- GPareto(loc = loc, scale = scale, shape = shape)
-
- ## starting parameters
- startPar <- function(x,...){
- tr <- theta[1]
-
- if(any(x < tr)) stop("some data smaller than 'loc' parameter")
-
- ## Pickand estimator
- if(is.null(start0Est)){
- e0 <- estimate(medkMADhybr(x, k=10, ParamFamily=GParetoFamily(loc = theta[1],
- scale = theta[2], shape = theta[3]),
- q.lo = 1e-3, q.up = 15))
- }else{
- if(is(start0Est,"function")){
- e1 <- start0Est(x, ...)
- e0 <- if(is(e1,"Estimate")) estimate(e1) else e1
- }
- if(!is.null(names(e0)))
- e0 <- e0[c("scale", "shape")]
- }
- names(e0) <- NULL
- return(e0)
- }
-
- modifyPar <- function(theta){
- if(!is.null(names)){
- sc <- theta["scale"]
- sh <- theta["shape"]
- }else{
- theta <- abs(theta)
- sc <- theta[1]
- sh <- theta[2]
- }
- GPareto(loc = loc, scale = sc, shape = sh)
- }
-
- ## what to do in case of leaving the parameter domain
- makeOKPar <- function(theta) {
- if(!is.null(names)){
- sc <- theta["scale"]
- sh <- theta["shape"]
- }else{
- theta <- abs(theta)
- sc <- theta[1]
- sh <- theta[2]
- }
- theta[2] <- pmin(sh,10)
- return(theta)
- }
-
- ## L2-derivative of the distribution
- L2deriv.fct <- function(param) {
- sc <- force(main(param)[1])
- k <- force(main(param)[2])
- tr <- fixed(param)[1]
-
- Lambda1 <- function(x) {
- y <- x*0
- x0 <- (x-tr)/sc
- x1 <- x0[x0>0]
- y[x0>0] <- -1/sc + (1+k)/(1+k*x1)*x1/sc
- return(y)
- }
- Lambda2 <- function(x) {
- y <- x*0
- x0 <- (x-tr)/sc
- x1 <- x0[x0>0]
- y[x0>0] <- log(1+k*x1)/k^2 - (1/k+1)*x1/(1+k*x1)
- return(y)
- }
- ## additional centering of scores to increase numerical precision!
- z1 <- E(distribution, fun=Lambda1)
- z2 <- E(distribution, fun=Lambda2)
- return(list(function(x){ Lambda1(x)-z1 },function(x){ Lambda2(x)-z2 }))
- }
-
- ## Fisher Information matrix as a function of parameters
- FisherInfo.fct <- function(param) {
- sc <- force(main(param)[1])
- k <- force(main(param)[2])
-# tr <- force(fixed(param)[1])
-# fct <- L2deriv.fct(param)
-# P2 <- GPareto(loc = tr, scale = sc, shape = k)
- E11 <- sc^-2
- E12 <- (sc*(1+k))^-1
- E22 <- 2/(1+k)
- return(PosSemDefSymmMatrix(matrix(c(E11,E12,E12,E22)/(1+2*k),2,2)))
- }
-
- FisherInfo <- FisherInfo.fct(param)
- name <- "Generalized Pareto Family"
-
- ## initializing the GPareto family with components of L2-family
- res <- L2ParamFamily(name = name, param = param,
- distribution = distribution,
- L2deriv.fct = L2deriv.fct,
- FisherInfo.fct = FisherInfo.fct,
- FisherInfo = FisherInfo,
- startPar = startPar,
- makeOKPar = makeOKPar,
- modifyParam = modifyPar,
- .returnClsName = "GParetoFamily")
- f.call <- substitute(GParetoFamily(loc = loc0, scale = scale0, shape = shape0,
- of.interest = of.interest0, p = p0,
- N = N0, trafo = trafo0),
- list(loc0 = loc, scale0 = scale, shape0 = shape,
- of.interest0 = of.interest, p0 = p, N0 = N,
- trafo0 = trafo))
- res at fam.call <- f.call
- return(res)
-}
-
Modified: branches/distr-2.4/pkg/distrMod/R/L2GroupFamilies-methods.R
===================================================================
--- branches/distr-2.4/pkg/distrMod/R/L2GroupFamilies-methods.R 2012-05-19 16:07:43 UTC (rev 808)
+++ branches/distr-2.4/pkg/distrMod/R/L2GroupFamilies-methods.R 2012-05-19 16:16:12 UTC (rev 809)
@@ -7,6 +7,9 @@
setMethod("locscalename", signature(object = "L2LocationScaleUnion"),
function(object) object at locscalename)
+setMethod("withPos", signature(object = "L2ScaleShapeUnion"),
+ function(object) object at withPos)
+
setReplaceMethod("LogDeriv", "L2GroupParamFamily",
function(object, value){
object at LogDeriv <- value
@@ -20,3 +23,11 @@
object at locscalename <- value
object
})
+
+setReplaceMethod("withPos", "L2ScaleShapeUnion",
+ function(object, value){
+ if(length(value)!=1)
+ stop("value of slot 'withPos' must be of length one")
+ object at withPos <- value
+ object
+ })
Deleted: branches/distr-2.4/pkg/distrMod/R/LDEstimator.R
===================================================================
--- branches/distr-2.4/pkg/distrMod/R/LDEstimator.R 2012-05-19 16:07:43 UTC (rev 808)
+++ branches/distr-2.4/pkg/distrMod/R/LDEstimator.R 2012-05-19 16:16:12 UTC (rev 809)
@@ -1,206 +0,0 @@
-.prepend <- function(prep, list0, dots = NULL){
- if(length(list0)+length(dots)==0) return(list(prep))
- n <- length(list0) + 1
- list1 <- vector("list",n)
- list1[[1]] <- prep
- names(list1)[1] <- "x"
- if(n>1) for(i in 2:n) {
- list1[[i]] <- list0[[i-1]]
- names(list1)[i] <- names(list0)[i-1]}
- ldots <- length(dots)
- l1 <- length(list1)
- if(ldots) {
- for( i in 1:ldots){
- list1[[l1+i]] <- dots[[i]]
- names(list1)[l1+i] <- names(dots)[i]
- }
- }
- return(list1)
-}
-
-.LDMatch <- function(x.0, loc.est.0,disp.est.0,
- loc.fctal.0, disp.fctal.0, ParamFamily.0,
- loc.est.ctrl.0 = NULL, loc.fctal.ctrl.0=NULL,
- disp.est.ctrl.0 = NULL, disp.fctal.ctrl.0=NULL,
- q.lo.0 =0, q.up.0=Inf, log.q.0 =TRUE, ...
- ){
- dots <- list(...)
- loc.emp <- do.call(loc.est.0, args = .prepend(x.0,loc.est.ctrl.0, dots))
- disp.emp <- do.call(disp.est.0, args = .prepend(x.0,disp.est.ctrl.0, dots))
- q.emp <- if(log.q.0) log(loc.emp)-log(disp.emp) else loc.emp/disp.emp
- q.f <- function(xi){
- distr.new <- ParamFamily.0 at modifyParam(theta=c("scale"=1,"shape"=xi))
- loc.th <- do.call(loc.fctal.0, args = .prepend(distr.new,loc.fctal.ctrl.0, dots))
- sc.th <- do.call(disp.fctal.0, args = .prepend(distr.new,disp.fctal.ctrl.0, dots))
- val <- if(log.q.0) log(loc.th)-log(sc.th) - q.emp else
- loc.th/sc.th-q.emp
- return(val)
- }
- xi.0 <- uniroot(q.f,lower=q.lo.0,upper=q.up.0)$root
- distr.new.0 <- ParamFamily.0 at modifyParam(theta=c("scale"=1,"shape"=xi.0))
- m1xi <- do.call(loc.fctal.0, args = .prepend(distr.new.0,loc.fctal.ctrl.0, dots))
- val <- c("shape"=xi.0,"scale"=loc.emp/m1xi)
- return(val)
-}
-
-LDEstimator <- function(x, loc.est, disp.est,
- loc.fctal, disp.fctal, ParamFamily,
- loc.est.ctrl = NULL, loc.fctal.ctrl=NULL,
- disp.est.ctrl = NULL, disp.fctal.ctrl=NULL,
- q.lo =1e-3, q.up=15, log.q =TRUE,
- name, Infos, asvar = NULL, nuis.idx = NULL,
- trafo = NULL, fixed = NULL, asvar.fct = NULL, na.rm = TRUE,
- ...){
- param0 <- main(param(ParamFamily))
- if(!all(c("shape","scale") %in% names(param0)))
- stop("LDEstimators expect shape-scale models.")
- name.est <- "LDEstimator"
- es.call <- match.call()
- if(missing(name))
- name <- "Some estimator"
- LDnames <- paste("Location:",
- paste(deparse(substitute(loc.fctal))),
- " ","Dispersion:",
- paste(deparse(substitute(disp.fctal))))
- estimator <- function(x,...){
- .LDMatch(x.0= x,
- loc.est.0 = loc.est, disp.est.0 = disp.est,
- loc.fctal.0 = loc.fctal, disp.fctal.0 = disp.fctal,
- ParamFamily.0 = ParamFamily,
- loc.est.ctrl.0 = loc.est.ctrl,
- loc.fctal.ctrl.0 = loc.fctal.ctrl,
- disp.est.ctrl.0 = disp.est.ctrl,
- disp.fctal.ctrl.0 = disp.fctal.ctrl,
- q.lo.0 = q.lo, q.up.0 = q.up, log.q.0 = log.q)
- }
-
-
- asvar.0 <- asvar
- nuis.idx.0 <- nuis.idx
- trafo.0 <- trafo
- fixed.0 <- fixed
- na.rm.0 <- na.rm
-
- estimate <- Estimator(x, estimator, name, Infos,
- asvar = asvar.0, nuis.idx = nuis.idx.0,
- trafo = trafo.0, fixed = fixed.0,
- na.rm = na.rm.0, ...)
- if(missing(asvar)) asvar <- NULL
- if(is.null(asvar))
- if(!missing(asvar.fct))
- if(!is.null(asvar.fct))
- asvar <- asvar.fct(ParamFamily, estimate, ...)
-
- estimate at untransformed.asvar <- asvar
-
- l.e <- length(estimate at untransformed.estimate)
- idx <- NULL
- idm <- 1:l.e
- if(!is.null(nuis.idx))
- {idx <- nuis.idx
- idm <- idm[-idx]
- mat <- diag(length(idm))}
-
- if(!.isUnitMatrix(estimate at trafo$mat)){
- estimate at estimate <- estimate at trafo$fct(estimate)
- if(!is.null(asvar))
- estimate at asvar <- estimate at trafo$mat%*%asvar[idm,idm]%*%t(estimate at trafo$mat)
- }
-
- estimate at estimate.call <- es.call
-
- if(missing(Infos))
- Infos <- matrix(c("LDEstimator", LDnames),
- ncol=2, dimnames=list(character(0), c("method", "message")))
- else{
- Infos <- matrix(c(rep("LDEstimator", length(Infos)+1), c(LDnames,Infos)),
- ncol = 2)
- colnames(Infos) <- c("method", "message")
- }
- estimate at Infos <- Infos
- return(estimate)
-}
-
-
-medkMAD <- function(x, k=1, ParamFamily, q.lo =1e-3, q.up=15, nuis.idx = NULL,
- trafo = NULL, fixed = NULL, asvar.fct = NULL, na.rm = TRUE,
- ...){
- es.call <- match.call()
- if(missing(k)) k <- 1
- es <- LDEstimator(x, loc.est = median, disp.est = kMAD,
- loc.fctal = median, disp.fctal = kMAD,
- ParamFamily = ParamFamily,
- loc.est.ctrl = NULL, loc.fctal.ctrl = NULL,
- disp.est.ctrl = list(k=k, na.rm = na.rm),
- disp.fctal.ctrl=list(k=k),
- q.lo =q.lo, q.up=q.up, log.q=TRUE,
- name = "medkMAD", Infos="medkMAD",
- asvar = NULL, nuis.idx = nuis.idx, trafo = trafo, fixed = fixed,
- asvar.fct = asvar.fct, na.rm = na.rm, ...)
- es at estimate.call <- es.call
- return(es)
- }
-
-medQn <- function(x, ParamFamily, q.lo =1e-3, q.up=15, nuis.idx = NULL,
- trafo = NULL, fixed = NULL, asvar.fct = NULL, na.rm = TRUE,
- ...){
- es.call <- match.call()
- es <- LDEstimator(x, loc.est = median, disp.est = Qn,
- loc.fctal = median, disp.fctal = Qn,
- ParamFamily = ParamFamily,
- loc.est.ctrl = NULL, loc.fctal.ctrl = NULL,
- disp.est.ctrl = list(constant=1,na.rm = na.rm),
- disp.fctal.ctrl = NULL,
- q.lo =q.lo, q.up=q.up, log.q=TRUE,
- name = "medQn", Infos="medQn",
- asvar = NULL, nuis.idx = nuis.idx, trafo = trafo, fixed = fixed,
- asvar.fct = asvar.fct, na.rm = na.rm, ...)
- es at estimate.call <- es.call
- return(es)
- }
-
-medSn <- function(x, ParamFamily, q.lo =1e-3, q.up=15, nuis.idx = NULL,
- trafo = NULL, fixed = NULL, asvar.fct = NULL, na.rm = TRUE,
- accuracy = 100, ...){
- es.call <- match.call()
- es <- LDEstimator(x, loc.est = median, disp.est = Sn,
- loc.fctal = median, disp.fctal = Sn,
- ParamFamily = ParamFamily,
- loc.est.ctrl = NULL, loc.fctal.ctrl = NULL,
- disp.est.ctrl = list(constant=1,na.rm = na.rm),
- disp.fctal.ctrl = list(accuracy=accuracy),
- q.lo =q.lo, q.up=q.up, log.q=TRUE,
- name = "medSn", Infos="medSn",
- asvar = NULL, nuis.idx = nuis.idx, trafo = trafo, fixed = fixed,
- asvar.fct = asvar.fct, na.rm = na.rm, ...)
- es at estimate.call <- es.call
- return(es)
- }
-
-medkMADhybr <- function(x, k=1, ParamFamily, q.lo =1e-3, q.up=15,
- KK=20, nuis.idx = NULL,
- trafo = NULL, fixed = NULL, asvar.fct = NULL, na.rm = TRUE,
- ...){
- i <- 1
- es <- try(medkMAD(x, k = k, ParamFamily = ParamFamily,
- q.lo = q.lo, q.up = q.up,
- nuis.idx = nuis.idx, trafo = trafo,
- fixed = fixed, asvar.fct = asvar.fct, na.rm = na.rm,
- ...), silent=TRUE)
- if(! any(is.na(es)) && !is(es,"try-error"))
- {return(es)}
-
- k1 <- 3.23
- while(i<KK){
- i <- i + 1
- es <- try(medkMAD(x, k = k1, ParamFamily = ParamFamily,
- q.lo = q.lo, q.up = q.up,
- nuis.idx = nuis.idx, trafo = trafo,
- fixed = fixed, asvar.fct = asvar.fct, na.rm = na.rm,
- ...), silent=TRUE)
- k1 <- k1 * 3
- if(! any(is.na(es)) && !is(es,"try-error"))
- {return(es)}
- }
- return(c("scale"=NA,"shape"=NA))
-}
Modified: branches/distr-2.4/pkg/distrMod/R/MLEstimator.R
===================================================================
--- branches/distr-2.4/pkg/distrMod/R/MLEstimator.R 2012-05-19 16:07:43 UTC (rev 808)
+++ branches/distr-2.4/pkg/distrMod/R/MLEstimator.R 2012-05-19 16:16:12 UTC (rev 809)
@@ -52,5 +52,3 @@
return(res)
}
-
-
\ No newline at end of file
[TRUNCATED]
To get the complete diff run:
svnlook diff /svnroot/distr -r 809
More information about the Distr-commits
mailing list