[Robast-commits] r747 - in branches/robast-1.0/pkg/RobExtremes: R man
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Wed Apr 9 23:20:48 CEST 2014
Author: ruckdeschel
Date: 2014-04-09 23:20:47 +0200 (Wed, 09 Apr 2014)
New Revision: 747
Modified:
branches/robast-1.0/pkg/RobExtremes/R/GEVFamily.R
branches/robast-1.0/pkg/RobExtremes/R/GEVFamilyMuUnknown.R
branches/robast-1.0/pkg/RobExtremes/R/startEstGEV.R
branches/robast-1.0/pkg/RobExtremes/man/GEVFamily.Rd
branches/robast-1.0/pkg/RobExtremes/man/GEVFamilyMuUnknown.Rd
Log:
RobExtremes: startEstGEV.R now works with soft bound 1+ xi (x-mu)/sigma > 0 (only to hold for lower quantile...) controlled by argument secLevel
Modified: branches/robast-1.0/pkg/RobExtremes/R/GEVFamily.R
===================================================================
--- branches/robast-1.0/pkg/RobExtremes/R/GEVFamily.R 2014-04-07 13:25:11 UTC (rev 746)
+++ branches/robast-1.0/pkg/RobExtremes/R/GEVFamily.R 2014-04-09 21:20:47 UTC (rev 747)
@@ -146,6 +146,7 @@
of.interest = c("scale", "shape"),
p = NULL, N = NULL, trafo = NULL,
start0Est = NULL, withPos = TRUE,
+ secLevel = 0.7,
withCentL2 = FALSE,
withL2derivDistr = FALSE,
..ignoreTrafo = FALSE,
@@ -237,6 +238,8 @@
## starting parameters
startPar <- function(x,...){
mu <- theta[1]
+ n <- length(x)
+ epsn <- min(floor(secLevel*sqrt(n))+1,n)
## Pickand estimator
if(is.null(start0Est)){
@@ -256,12 +259,12 @@
e0 <- e0[c("scale", "shape")]
}
# print(e0); print(str(x)); print(head(summary(x))); print(mu)
- if(e0[2]>0) if(any(x < mu-e0[1]/e0[2]))
+ if(quantile(e0[2]/e0[1]*(x-mu), epsn/n)< (-1)){
+ if(e0[2]>0)
stop("shape is positive and some data smaller than 'loc-scale/shape' ")
-
- if(e0[2]<0) if(any(x > mu-e0[1]/e0[2]))
+ else
stop("shape is negative and some data larger than 'loc-scale/shape' ")
-
+ }
names(e0) <- NULL
return(e0)
}
Modified: branches/robast-1.0/pkg/RobExtremes/R/GEVFamilyMuUnknown.R
===================================================================
--- branches/robast-1.0/pkg/RobExtremes/R/GEVFamilyMuUnknown.R 2014-04-07 13:25:11 UTC (rev 746)
+++ branches/robast-1.0/pkg/RobExtremes/R/GEVFamilyMuUnknown.R 2014-04-09 21:20:47 UTC (rev 747)
@@ -32,6 +32,7 @@
of.interest = c("scale", "shape"),
p = NULL, N = NULL, trafo = NULL,
start0Est = NULL, withPos = TRUE,
+ secLevel = 0.7,
withCentL2 = FALSE,
withL2derivDistr = FALSE,
..ignoreTrafo = FALSE,
@@ -125,6 +126,9 @@
## starting parameters
startPar <- function(x,...){
+ n <- length(x)
+ epsn <- min(floor(secLevel*sqrt(n))+1,n)
+
## Pickand estimator
if(is.null(start0Est)){
### replaced 20140402: CvMMDE-with xi on Grid
@@ -142,12 +146,12 @@
e0 <- e0[c("loc","scale", "shape")]
}
# print(e0); print(str(x)); print(head(summary(x))); print(mu)
- if(e0[3]>0) if(any(x < e0[1]-e0[2]/e0[3]))
+ if(quantile(e0[3]/e0[2]*(x-e0[1]), epsn/n)< (-1)){
+ if(e0[3]>0)
stop("shape is positive and some data smaller than 'loc-scale/shape' ")
-
- if(e0[3]<0) if(any(x > e0[1]-e0[2]/e0[3]))
+ else if(e0[3]<0)
stop("shape is negative and some data larger than 'loc-scale/shape' ")
-
+ }
names(e0) <- NULL
return(e0)
}
Modified: branches/robast-1.0/pkg/RobExtremes/R/startEstGEV.R
===================================================================
--- branches/robast-1.0/pkg/RobExtremes/R/startEstGEV.R 2014-04-07 13:25:11 UTC (rev 746)
+++ branches/robast-1.0/pkg/RobExtremes/R/startEstGEV.R 2014-04-09 21:20:47 UTC (rev 747)
@@ -1,7 +1,11 @@
-.getXiGrid <- function(){seq(-0.48,5,by=0.5)}
+.getXiGrid <- function(){c(0.1,seq(-0.48,5,by=0.5))}
-.getBetaXiGEV <- function(x, mu, xiGrid = .getXiGrid(), withPos=TRUE){
+.getBetaXiGEV <- function(x, mu, xiGrid = .getXiGrid(), withPos=TRUE, secLevel = 0.7){
+
+ n <- length(x)
+ epsn <- min(floor(secLevel*sqrt(n))+1,n)
+
x0 <- x-mu
s0 <- max(x0)-min(x0)
crit0 <- Inf
@@ -11,7 +15,7 @@
mygev1 <- GEV(loc=0,scale=sig,shape=xi)
CvMDist(x0,mygev1)
}
- intlo <- max(-xi*(x-mu))
+ intlo <- quantile(-xi*(x-mu),1-epsn/n)
intv <- c(max(1e-5,intlo), s0)
sigCvMMD1 <- optimize(funl, interval=intv)$minimum
mygev <- GEVFamily(loc=0,scale=sigCvMMD1,shape=xi, withPos=withPos,
@@ -21,7 +25,7 @@
if(!is(mde0,"try-error")){
es <- estimate(mde0)
crit1 <- criterion(mde0)
- if(1+min(es[2]*x0/es[1])>0){
+ if(quantile(1+es[2]*x0/es[1], epsn/n)>0){
if(crit1<crit0){
mdeb <- mde0
crit0 <- crit1
@@ -42,9 +46,9 @@
return(es)
}
-.getMuBetaXiGEV <- function(x, xiGrid = .getXiGrid(), withPos=TRUE){
+.getMuBetaXiGEV <- function(x, xiGrid = .getXiGrid(), withPos=TRUE, secLevel = 0.7){
mu <- quantile(x,exp(-1))
- es <- .getBetaXiGEV(x=x, mu=mu, xiGrid=xiGrid, withPos=withPos)
+ es <- .getBetaXiGEV(x=x, mu=mu, xiGrid=xiGrid, withPos=withPos, secLevel = secLevel)
es0 <- c(mu,es)
names(es0) <- c("loc","scale","shape")
return(es0)
Modified: branches/robast-1.0/pkg/RobExtremes/man/GEVFamily.Rd
===================================================================
--- branches/robast-1.0/pkg/RobExtremes/man/GEVFamily.Rd 2014-04-07 13:25:11 UTC (rev 746)
+++ branches/robast-1.0/pkg/RobExtremes/man/GEVFamily.Rd 2014-04-09 21:20:47 UTC (rev 747)
@@ -9,8 +9,8 @@
\usage{
GEVFamily(loc = 0, scale = 1, shape = 0.5, of.interest = c("scale", "shape"),
p = NULL, N = NULL, trafo = NULL, start0Est = NULL, withPos = TRUE,
- withCentL2 = FALSE, withL2derivDistr = FALSE, ..ignoreTrafo = FALSE,
- ..withWarningGEV = TRUE)
+ secLevel = 0.7, withCentL2 = FALSE, withL2derivDistr = FALSE,
+ ..ignoreTrafo = FALSE, ..withWarningGEV = TRUE)
}
\arguments{
\item{loc}{ real: known/fixed threshold/location parameter }
@@ -24,6 +24,13 @@
\item{trafo}{ matrix or NULL: transformation of the parameter }
\item{start0Est}{ startEstimator --- if \code{NULL} \code{\link{PickandsEstimator}} is used }
\item{withPos}{ logical of length 1: Is shape restricted to positive values? }
+ \item{secLevel}{ a numeric of length 1:
+ In the ideal GEV model, for each observastion \eqn{X_i}{Xi}, the expression
+ \eqn{1+\frac{{\rm shape}(X_i-{\rm loc})}{{\rm scale}}}{1+shape(Xi-loc)/scale}
+ must be positive, which in principle could be attacked by a single outlier.
+ Hence for sample size \eqn{n} we allow for \eqn{\varepsilon n}{eps n}
+ violations, interpreting the violations as outliers. Here
+ \eqn{\varepsilon = {\tt secLevel}/\sqrt{n}}{eps = secLevel/sqrt(n)}. }
\item{withCentL2}{logical: shall L2 derivative be centered by substracting
the E()? Defaults to \code{FALSE}, but higher accuracy can be achieved
when set to \code{TRUE}.}
Modified: branches/robast-1.0/pkg/RobExtremes/man/GEVFamilyMuUnknown.Rd
===================================================================
--- branches/robast-1.0/pkg/RobExtremes/man/GEVFamilyMuUnknown.Rd 2014-04-07 13:25:11 UTC (rev 746)
+++ branches/robast-1.0/pkg/RobExtremes/man/GEVFamilyMuUnknown.Rd 2014-04-09 21:20:47 UTC (rev 747)
@@ -9,8 +9,8 @@
\usage{
GEVFamilyMuUnknown(loc = 0, scale = 1, shape = 0.5, of.interest = c("scale", "shape"),
p = NULL, N = NULL, trafo = NULL, start0Est = NULL, withPos = TRUE,
- withCentL2 = FALSE, withL2derivDistr = FALSE, ..ignoreTrafo = FALSE,
- ..withWarningGEV = TRUE)
+ secLevel = 0.7, withCentL2 = FALSE, withL2derivDistr = FALSE,
+ ..ignoreTrafo = FALSE, ..withWarningGEV = TRUE)
}
\arguments{
\item{loc}{ real: known/fixed threshold/location parameter }
@@ -24,6 +24,13 @@
\item{trafo}{ matrix or NULL: transformation of the parameter }
\item{start0Est}{ startEstimator --- if \code{NULL} \code{\link{PickandsEstimator}} is used }
\item{withPos}{ logical of length 1: Is shape restricted to positive values? }
+ \item{secLevel}{ a numeric of length 1:
+ In the ideal GEV model, for each observastion \eqn{X_i}{Xi}, the expression
+ \eqn{1+\frac{{\rm shape}(X_i-{\rm loc})}{{\rm scale}}}{1+shape(Xi-loc)/scale}
+ must be positive, which in principle could be attacked by a single outlier.
+ Hence for sample size \eqn{n} we allow for \eqn{\varepsilon n}{eps n}
+ violations, interpreting the violations as outliers. Here
+ \eqn{\varepsilon = {\tt secLevel}/\sqrt{n}}{eps = secLevel/sqrt(n)}. }
\item{withCentL2}{logical: shall L2 derivative be centered by substracting
the E()? Defaults to \code{FALSE}, but higher accuracy can be achieved
when set to \code{TRUE}.}
More information about the Robast-commits
mailing list