[Distr-commits] r807 - in branches/distr-2.4/pkg: distrEx/R distrEx/man distrMod/R distrMod/man
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Thu May 17 17:05:05 CEST 2012
Author: ruckdeschel
Date: 2012-05-17 17:05:05 +0200 (Thu, 17 May 2012)
New Revision: 807
Modified:
branches/distr-2.4/pkg/distrEx/R/SnQn.R
branches/distr-2.4/pkg/distrEx/man/Var.Rd
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/LDEstimator.Rd
Log:
a little speed up for Sn and medSn
Modified: branches/distr-2.4/pkg/distrEx/R/SnQn.R
===================================================================
--- branches/distr-2.4/pkg/distrEx/R/SnQn.R 2012-05-16 11:28:28 UTC (rev 806)
+++ branches/distr-2.4/pkg/distrEx/R/SnQn.R 2012-05-17 15:05:05 UTC (rev 807)
@@ -42,11 +42,16 @@
setMethod("Sn", signature(x = "UnivariateDistribution"),
- function(x, low=0,upp=20, accuracy = 1000, ...){
-
+ function(x, low=0,upp=1.01, accuracy = 1000, ...){
+ m0 <- median(x)
+ M0 <- mad(x)
g <- function(xx){
fct <- function(m) p(x)(m+xx)-p(x)(-m+xx)-0.5
- m <- try(uniroot(fct, lower = low, upper = upp*q(x)(0.75))$root,
+ up0 <- upp*(M0+abs(m0-xx))
+ m <- try(uniroot(fct, lower = low,
+ upper = up0,
+ f.lower=if(low<1e-12) -0.5 else fct(low),
+ f.upper=max(fct(up0),1e-8))$root,
silent = TRUE)
if(is(m,"try-error")) {
# print("error")
@@ -54,9 +59,8 @@
}else{ return(m) }
}
- x0 <- q(x)(seq(1e-2/accuracy,1-1e-2/accuracy,length=accuracy))
+ x0 <- q(x)(seq(.5/accuracy,1-.5/accuracy,length=accuracy))
y <- sapply(x0,g)
-
c0 <- median(y,na.rm=TRUE)
return(c0)
})
@@ -82,7 +86,7 @@
setMethod("Sn", signature(x = "Norm"),
function(x, ...){
- return(sd(x)* 0.8385098038)
+ return(sd(x)* 0.838504603)
})
setMethod("Qn", signature(x = "Norm"),
Modified: branches/distr-2.4/pkg/distrEx/man/Var.Rd
===================================================================
--- branches/distr-2.4/pkg/distrEx/man/Var.Rd 2012-05-16 11:28:28 UTC (rev 806)
+++ branches/distr-2.4/pkg/distrEx/man/Var.Rd 2012-05-17 15:05:05 UTC (rev 807)
@@ -329,7 +329,7 @@
Sn(x, ...)
\S4method{Sn}{ANY}(x, ...)
-\S4method{Sn}{UnivariateDistribution}(x, low = 0, upp = 20, accuracy = 1000, ...)
+\S4method{Sn}{UnivariateDistribution}(x, low = 0, upp = 1.01, accuracy = 1000, ...)
\S4method{Sn}{DiscreteDistribution}(x, ...)
\S4method{Sn}{AffLinDistribution}(x, ...)
\S4method{Sn}{Norm}(x, ...)
@@ -355,10 +355,12 @@
\item{q00}{numeric or NULL: determines search interval (from \code{-q00} to \code{q00})
for \code{Qn}; if \code{NULL} (default) \code{q00}
is set to \code{10*q(x)(3/4)} internally. }
- \item{low}{numeric; lower bound for search interval; defaults to \code{0}. }
- \item{upp}{numeric; upper bound for search interval; defaults to \code{20}.
- Is used internally as \code{upp*q(x)(0.75)}. }
- \item{accuracy}{numeric; number of grid points; defaults to \code{1000}. }
+ \item{low}{numeric; lower bound for search interval for median(abs(x-Y)) where
+ Y (a real constant) runs over the range of x; defaults to \code{0}. }
+ \item{upp}{numeric; upper bound for search interval for median(abs(x-Y)) where
+ Y (a real constant) runs over the range of x; defaults to \code{1.01}.
+ Is used internally as \code{upp*(mad(x)+abs(median(x)-Y))}. }
+ \item{accuracy}{numeric; number of grid points for \code{Sn}; defaults to \code{1000}. }
}
\value{
The value of the corresponding functional at the distribution in the argument is computed.
Modified: branches/distr-2.4/pkg/distrMod/R/GParetoFamily.R
===================================================================
--- branches/distr-2.4/pkg/distrMod/R/GParetoFamily.R 2012-05-16 11:28:28 UTC (rev 806)
+++ branches/distr-2.4/pkg/distrMod/R/GParetoFamily.R 2012-05-17 15:05:05 UTC (rev 807)
@@ -257,9 +257,9 @@
## Pickand estimator
if(is.null(start0Est)){
- e0 <- estimate(medkMADhybr(x, ParamFamily=GParetoFamily(loc = theta[1],
+ e0 <- estimate(medkMADhybr(x, k=10, ParamFamily=GParetoFamily(loc = theta[1],
scale = theta[2], shape = theta[3]),
- q.lo = 1e-6, q.up = 20))
+ q.lo = 1e-3, q.up = 15))
}else{
if(is(start0Est,"function")){
e1 <- start0Est(x, ...)
Modified: branches/distr-2.4/pkg/distrMod/R/LDEstimator.R
===================================================================
--- branches/distr-2.4/pkg/distrMod/R/LDEstimator.R 2012-05-16 11:28:28 UTC (rev 806)
+++ branches/distr-2.4/pkg/distrMod/R/LDEstimator.R 2012-05-17 15:05:05 UTC (rev 807)
@@ -47,7 +47,7 @@
loc.fctal, disp.fctal, ParamFamily,
loc.est.ctrl = NULL, loc.fctal.ctrl=NULL,
disp.est.ctrl = NULL, disp.fctal.ctrl=NULL,
- q.lo =0, q.up=Inf, log.q =TRUE,
+ 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,
...){
@@ -122,7 +122,7 @@
}
-medkMAD <- function(x, k=1, ParamFamily, q.lo =1e-6, q.up=10, nuis.idx = NULL,
+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()
@@ -141,7 +141,7 @@
return(es)
}
-medQn <- function(x, ParamFamily, q.lo =1e-6, q.up=10, nuis.idx = NULL,
+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()
@@ -159,9 +159,9 @@
return(es)
}
-medSn <- function(x, ParamFamily, q.lo =1e-6, q.up=10, nuis.idx = NULL,
+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 = 1000, ...){
+ accuracy = 100, ...){
es.call <- match.call()
es <- LDEstimator(x, loc.est = median, disp.est = Sn,
loc.fctal = median, disp.fctal = Sn,
@@ -177,7 +177,7 @@
return(es)
}
-medkMADhybr <- function(x, k=1, ParamFamily, q.lo =1e-6, q.up=10,
+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,
...){
Modified: branches/distr-2.4/pkg/distrMod/man/GParetoFamily.Rd
===================================================================
--- branches/distr-2.4/pkg/distrMod/man/GParetoFamily.Rd 2012-05-16 11:28:28 UTC (rev 806)
+++ branches/distr-2.4/pkg/distrMod/man/GParetoFamily.Rd 2012-05-17 15:05:05 UTC (rev 807)
@@ -15,12 +15,12 @@
\item{scale}{ positive real: scale parameter }
\item{shape}{ positive real: shape parameter }
\item{of.interest}{ character: which parameters, transformations are of interest.\cr
- posibilites are: "scale", "shape", "quantile", "expected loss",
+ possibilites are: "scale", "shape", "quantile", "expected loss",
"expected shortfall"; a maximum number of two of these may be selected }
\item{p}{real or NULL: probability needed for quantile and expected shortfall }
\item{N}{real or NULL: expected frequency for expected loss }
\item{trafo}{ matrix or NULL: transformation of the parameter }
- \item{start0Est}{ startEstimator for MLE and MDE --- if NULL HybridEstimator is used }
+ \item{start0Est}{ startEstimator --- if \code{NULL} \code{\link{medkMADhybr}} is used }
}
\details{
The slots of the corresponding L2 differentiable
@@ -29,7 +29,18 @@
\value{Object of class \code{"L2ParamFamily"}}
\references{
Kohl, M. (2005) \emph{Numerical Contributions to
- the Asymptotic Theory of Robustness}. Bayreuth: Dissertation.}
+ the Asymptotic Theory of Robustness}. Bayreuth: Dissertation.\cr
+M.~Kohl, P. Ruckdeschel, H.~Rieder (2010):
+Infinitesimally Robust Estimation in General Smoothly Parametrized Models.
+\emph{Stat. Methods Appl.}, \bold{19}, 333--354.\cr
+ Ruckdeschel, P. and Horbenko, N. (2011): Optimally-Robust Estimators in Generalized
+ Pareto Models. ArXiv 1005.1476. To appear at \emph{Statistics}.
+ DOI: 10.1080/02331888.2011.628022. \cr
+ Ruckdeschel, P. and Horbenko, N. (2011): Yet another breakdown point notion:
+EFSBP --illustrated at scale-shape models. ArXiv 1005.1480. To appear at \emph{Metrika}.
+DOI: 10.1007/s00184-011-0366-4
+ }
+
\author{Matthias Kohl \email{Matthias.Kohl at stamats.de}\cr
Peter Ruckdeschel \email{peter.ruckdeschel at itwm.fraunhofer.de}\cr
Nataliya Horbenko \email{nataliya.horbenko at itwm.fraunhofer.de}}
Modified: branches/distr-2.4/pkg/distrMod/man/LDEstimator.Rd
===================================================================
--- branches/distr-2.4/pkg/distrMod/man/LDEstimator.Rd 2012-05-16 11:28:28 UTC (rev 806)
+++ branches/distr-2.4/pkg/distrMod/man/LDEstimator.Rd 2012-05-17 15:05:05 UTC (rev 807)
@@ -17,20 +17,20 @@
LDEstimator(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 =0, q.up=Inf, log.q =TRUE,
+ 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,
...)
-medkMAD(x, k=1, ParamFamily, q.lo =1e-6, q.up=10, nuis.idx = NULL,
+medkMAD(x, k=1, ParamFamily, q.lo =1e-3, q.up=15, nuis.idx = NULL,
trafo = NULL, fixed = NULL, asvar.fct = NULL, na.rm = TRUE,
...)
-medkMADhybr(x, k=1, ParamFamily, q.lo =1e-6, q.up=10, KK = 20, nuis.idx = NULL,
+medkMADhybr(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,
...)
-medSn(x, ParamFamily, q.lo =1e-6, q.up=10, nuis.idx = NULL,
+medSn(x, ParamFamily, q.lo =1e-3, q.up=15, nuis.idx = NULL,
trafo = NULL, fixed = NULL, asvar.fct = NULL, na.rm = TRUE,
- accuracy = 1000, ...)
-medQn(x, ParamFamily, q.lo =1e-6, q.up=10, nuis.idx = NULL,
+ accuracy = 100, ...)
+medQn(x, ParamFamily, q.lo =1e-3, q.up=15, nuis.idx = NULL,
trafo = NULL, fixed = NULL, asvar.fct = NULL, na.rm = TRUE,
...)
}
@@ -98,7 +98,13 @@
Uses \code{\link{.LDMatch}} internally.
}
\note{The values for \code{q.lo} and \code{q.up} are a bit delicate and
- have to be found, model by model, by try and error.}
+ have to be found, model by model, by try and error.
+ As a rule, \code{medSn} is rather slow, as the evaluation of the \code{Sn}
+ functional is quite expensive. So if \code{medSn} is the estimator of choice,
+ it pays off, for a given shape-scale family, to evaluate \code{medSn} on a
+ grid of shape-values (with scale 1) and then to use an interpolation techniques
+ in a particular method to replace the default one for this shape-scale family.
+ As an example, we have done so for the GPD family.}
\value{
An object of S4-class \code{"Estimate"}.
}
@@ -108,7 +114,7 @@
DOI: 10.1007/s00184-011-0366-4.
A. Marazzi and C. Ruffieux (1999): The truncated mean of asymmetric distribution.
-Computational Statistics and Data Analysis 32: 79–100 .
+Computational Statistics and Data Analysis 32: 79-100 .
}
%\references{ }
@@ -124,15 +130,15 @@
## parametric family of probability measures
G <- GammaFamily(scale = 1, shape = 2)
-medQn(x = x, ParamFamily = G, q.lo = 1e-3, q.up = 15)
-medSn(x = x, ParamFamily = G, q.lo = 1e-3, q.up = 15)
-medkMAD(x = x, ParamFamily = G, q.lo = 1e-3, q.up = 15)
-medkMADhybr(x = x, ParamFamily = G, q.lo = 1e-3, q.up = 15)
-medkMAD(x = x, k=10, ParamFamily = G, q.lo = 1e-3, q.up = 15)
+medQn(x = x, ParamFamily = G)
+medSn(x = x, ParamFamily = G)
+medkMAD(x = x, ParamFamily = G)
+medkMADhybr(x = x, ParamFamily = G)
+medkMAD(x = x, k=10, ParamFamily = G)
##not at all robust:
LDEstimator(x, loc.est = mean, disp.est = sd,
loc.fctal = E, disp.fctal = sd,
- ParamFamily = G, q.lo = 1e-3, q.up = 15)
+ ParamFamily = G)
}
\keyword{univar}
More information about the Distr-commits
mailing list