[Robast-commits] r572 - in branches/robast-0.9/pkg: ROptEst ROptEst/R ROptEst/inst/scripts ROptEst/man RobExtremes/inst/scripts
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Sun Jan 27 20:10:04 CET 2013
Author: ruckdeschel
Date: 2013-01-27 20:10:04 +0100 (Sun, 27 Jan 2013)
New Revision: 572
Added:
branches/robast-0.9/pkg/ROptEst/R/getInfRad.R
branches/robast-0.9/pkg/ROptEst/R/getRadius.R
branches/robast-0.9/pkg/ROptEst/man/getInfRad.Rd
branches/robast-0.9/pkg/ROptEst/man/getRadius.Rd
branches/robast-0.9/pkg/RobExtremes/inst/scripts/GumbelLocationModel.R
Removed:
branches/robast-0.9/pkg/ROptEst/inst/scripts/GumbelLocationModel.R
Modified:
branches/robast-0.9/pkg/ROptEst/NAMESPACE
branches/robast-0.9/pkg/ROptEst/R/AllGeneric.R
branches/robast-0.9/pkg/ROptEst/R/getInfClip.R
branches/robast-0.9/pkg/ROptEst/R/getRiskIC.R
branches/robast-0.9/pkg/ROptEst/man/getInfClip.Rd
branches/robast-0.9/pkg/ROptEst/man/getRiskIC.Rd
Log:
RobExtremes/ROptEst: Moved script GumbelLocationModel.R from the latter to the former pkg.
ROptEst: + new functionality to determine "optimal" radius (getRadius, getInfRad (the latter an S4method));
+"getInfClip", signature(clip = "numeric", L2deriv = "UnivariateDistribution",risk = "asSemivar", neighbor = "ContNeighborhood") regains argument biastype as it is called this way in getInfRobIC_asGRisk.R
+ getRiskIC gains example and for signature(IC = "HampIC", risk = "asCov", neighbor = "missing", L2Fam = "L2ParamFamily") in 1dim case now uses getInfV() .
Modified: branches/robast-0.9/pkg/ROptEst/NAMESPACE
===================================================================
--- branches/robast-0.9/pkg/ROptEst/NAMESPACE 2013-01-27 15:54:29 UTC (rev 571)
+++ branches/robast-0.9/pkg/ROptEst/NAMESPACE 2013-01-27 19:10:04 UTC (rev 572)
@@ -11,6 +11,7 @@
"getFixRobIC",
"getAsRisk",
"getFiRisk",
+ "getInfRad",
"getInfClip",
"getFixClip",
"getInfGamma",
@@ -31,7 +32,7 @@
"comparePlot", "getRiskFctBV")
export("getL2normL2deriv",
"asAnscombe", "asL1", "asL4",
- "getReq", "getMaxIneff")
+ "getReq", "getMaxIneff", "getRadius")
export("roptest","roptest.old", "robest",
"getLagrangeMultByOptim","getLagrangeMultByIter")
export("genkStepCtrl", "genstartCtrl", "gennbCtrl")
Modified: branches/robast-0.9/pkg/ROptEst/R/AllGeneric.R
===================================================================
--- branches/robast-0.9/pkg/ROptEst/R/AllGeneric.R 2013-01-27 15:54:29 UTC (rev 571)
+++ branches/robast-0.9/pkg/ROptEst/R/AllGeneric.R 2013-01-27 19:10:04 UTC (rev 572)
@@ -18,6 +18,10 @@
setGeneric("getInfClip",
function(clip, L2deriv, risk, neighbor, ...) standardGeneric("getInfClip"))
}
+if(!isGeneric("getInfRad")){
+ setGeneric("getInfRad",
+ function(clip, L2deriv, risk, neighbor, ...) standardGeneric("getInfRad"))
+}
if(!isGeneric("getFixClip")){
setGeneric("getFixClip",
function(clip, Distr, risk, neighbor, ...) standardGeneric("getFixClip"))
Modified: branches/robast-0.9/pkg/ROptEst/R/getInfClip.R
===================================================================
--- branches/robast-0.9/pkg/ROptEst/R/getInfClip.R 2013-01-27 15:54:29 UTC (rev 571)
+++ branches/robast-0.9/pkg/ROptEst/R/getInfClip.R 2013-01-27 19:10:04 UTC (rev 572)
@@ -155,7 +155,7 @@
L2deriv = "UnivariateDistribution",
risk = "asSemivar",
neighbor = "ContNeighborhood"),
- function(clip, L2deriv, risk, neighbor, cent, symm, trafo){
+ function(clip, L2deriv, risk, neighbor, biastype, cent, symm, trafo){
biastype <- if(sign(risk)==1) positiveBias() else negativeBias()
z0 <- getInfCent(L2deriv = L2deriv, risk = risk, neighbor = neighbor,
biastype = biastype,
Added: branches/robast-0.9/pkg/ROptEst/R/getInfRad.R
===================================================================
--- branches/robast-0.9/pkg/ROptEst/R/getInfRad.R (rev 0)
+++ branches/robast-0.9/pkg/ROptEst/R/getInfRad.R 2013-01-27 19:10:04 UTC (rev 572)
@@ -0,0 +1,175 @@
+###############################################################################
+## optimal radius for given clipping bound for asymptotic MSE
+###############################################################################
+
+setMethod("getInfRad", signature(clip = "numeric",
+ L2deriv = "UnivariateDistribution",
+ risk = "asMSE",
+ neighbor = "ContNeighborhood"),
+ function(clip, L2deriv, risk, neighbor, biastype,
+ cent, symm, trafo){
+ gamm <- getInfGamma(L2deriv = L2deriv, risk = risk, neighbor = neighbor,
+ biastype = biastype, cent = cent, clip = clip)
+ return((-gamm/clip)^.5)
+ })
+
+setMethod("getInfRad", signature(clip = "numeric",
+ L2deriv = "UnivariateDistribution",
+ risk = "asMSE",
+ neighbor = "TotalVarNeighborhood"),
+ function(clip, L2deriv, risk, neighbor, biastype,
+ cent, symm, trafo){
+ gamm <- getInfGamma(L2deriv = sign(as.vector(trafo))*L2deriv, risk = risk,
+ neighbor = neighbor, biastype = biastype,
+ cent = if(symm) -clip/2 else cent , clip = clip)
+ return((-gamm/clip)^.5)
+ })
+
+setMethod("getInfRad", signature(clip = "numeric",
+ L2deriv = "EuclRandVariable",
+ risk = "asMSE",
+ neighbor = "UncondNeighborhood"),
+ function(clip, L2deriv, risk, neighbor, biastype,
+ Distr, stand, cent, trafo){
+ gamm <- getInfGamma(L2deriv = L2deriv, risk = risk, neighbor = neighbor,
+ biastype = biastype, Distr = Distr, stand = stand,
+ cent = cent, clip = clip)
+ return((-gamm/clip)^.5)
+ })
+
+###############################################################################
+## optimal radius for given clipping bound for asymptotic L1risk
+###############################################################################
+
+setMethod("getInfRad", signature(clip = "numeric",
+ L2deriv = "UnivariateDistribution",
+ risk = "asL1",
+ neighbor = "ContNeighborhood"),
+ function(clip, L2deriv, risk, neighbor, biastype,
+ cent, symm, trafo){
+ s <- getInfV(L2deriv, neighbor, biastype, clip, cent, stand=1)
+ gamm <- getInfGamma(L2deriv = L2deriv, risk = risk, neighbor = neighbor,
+ biastype = biastype, cent = cent, clip = clip)
+ solvfct <- function(r){
+ w <- r * clip / s^.5
+ dp <- 2*dnorm(w)
+ pp <- 2*pnorm(w)-1
+ lhs <- s^.5*r*pp/dp
+ return(lhs + gamm)
+ }
+ r <- try(uniroot(solvfct, lower=1e-5, upper = 10)$root,silent=TRUE)
+ if(is(r, "try-error")) return(NA)
+ return(r)
+ })
+
+setMethod("getInfRad", signature(clip = "numeric",
+ L2deriv = "UnivariateDistribution",
+ risk = "asL1",
+ neighbor = "TotalVarNeighborhood"),
+ function(clip, L2deriv, risk, neighbor, biastype,
+ cent, symm, trafo){
+ s <- getInfV(L2deriv, neighbor, biastype, clip, cent, stand=1)
+ gamm <- getInfGamma(L2deriv = sign(as.vector(trafo))*L2deriv, risk = risk,
+ neighbor = neighbor, biastype = biastype,
+ cent = if(symm) -clip/2 else cent , clip = clip)
+ solvfct <- function(r){
+ w <- r * clip / s^.5
+ dp <- 2*dnorm(w)
+ pp <- 2*pnorm(w)-1
+ lhs <- s^.5*r*pp/dp
+ return(lhs + gamm)
+ }
+ r <- try(uniroot(solvfct, lower=1e-5, upper = 10)$root,silent=TRUE)
+ if(is(r, "try-error")) return(NA)
+ return(r)
+ })
+
+
+###############################################################################
+## optimal radius for given clipping bound for asymptotic L4 risk
+###############################################################################
+
+setMethod("getInfRad", signature(clip = "numeric",
+ L2deriv = "UnivariateDistribution",
+ risk = "asL4",
+ neighbor = "ContNeighborhood"),
+ function(clip, L2deriv, risk, neighbor, biastype,
+ cent, symm, trafo){
+ s <- getInfV(L2deriv, neighbor, biastype, clip, cent, stand=1)
+ gamm <- getInfGamma(L2deriv = L2deriv, risk = risk, neighbor = neighbor,
+ biastype = biastype, cent = cent, clip = clip)
+ solvfct <- function(r){
+ mse <- r^2 *clip^2 + s
+ mse4 <- (r^2 *clip^2/3 + s)/mse
+ r^2*clip*mse4 + gamm }
+ r <- try(uniroot(solvfct, lower=1e-5, upper = 10)$root,silent=TRUE)
+ if(is(r, "try-error")) return(NA)
+ return(r)
+ })
+
+setMethod("getInfRad", signature(clip = "numeric",
+ L2deriv = "UnivariateDistribution",
+ risk = "asL4",
+ neighbor = "TotalVarNeighborhood"),
+ function(clip, L2deriv, risk, neighbor, biastype,
+ cent, symm, trafo){
+ s <- getInfV(L2deriv, neighbor, biastype, clip, cent, stand=1)
+ gamm <- getInfGamma(L2deriv = sign(as.vector(trafo))*L2deriv, risk = risk,
+ neighbor = neighbor, biastype = biastype,
+ cent = if(symm) -clip/2 else cent , clip = clip)
+ solvfct <- function(r){
+ mse <- r^2 *clip^2 + s
+ mse4 <- (r^2 *clip^2/3 + s)/mse
+ r^2*clip*mse4 + gamm }
+ r <- try(uniroot(solvfct, lower=1e-5, upper = 10)$root,silent=TRUE)
+ if(is(r, "try-error")) return(NA)
+ return(r)
+ })
+
+
+
+###############################################################################
+## optimal radius for given clipping bound for asymptotic under-/overshoot risk
+###############################################################################
+setMethod("getInfRad", signature(clip = "numeric",
+ L2deriv = "UnivariateDistribution",
+ risk = "asUnOvShoot",
+ neighbor = "UncondNeighborhood"),
+ function(clip, L2deriv, risk, neighbor, biastype,
+ cent, symm, trafo){
+ gamm <- getInfGamma(L2deriv = sign(as.vector(trafo))*L2deriv, risk = risk,
+ neighbor = neighbor, biastype = biastype,
+ cent = if(symm) -clip/2 else cent , clip = clip)
+ return( -gamm*risk at width)
+ })
+
+###############################################################################
+## optimal radius for given clipping bound for asymptotic semivariance
+###############################################################################
+setMethod("getInfRad", signature(clip = "numeric",
+ L2deriv = "UnivariateDistribution",
+ risk = "asSemivar",
+ neighbor = "ContNeighborhood"),
+ function(clip, L2deriv, risk, neighbor, biastype, cent, symm, trafo){
+ biastype <- if(sign(risk)==1) positiveBias() else negativeBias()
+ z0 <- getInfCent(L2deriv = L2deriv, risk = risk, neighbor = neighbor,
+ biastype = biastype,
+ clip = max(clip, 1e-4), cent = 0, trafo = trafo,
+ symm = symm, tol.z = 1e-6)
+
+ ga <- getInfGamma(L2deriv = L2deriv, risk = risk, neighbor = neighbor,
+ biastype = biastype, cent = cent, clip = clip)
+
+
+ solvfct <- function(r){
+ si <- if (sign(risk)>0) 1 else -1
+ v0 <- E(L2deriv, function(x) pmin( x-z0, si*clip)^2 )
+ s0 <- sqrt(v0)
+ sv <- r * clip / s0
+ er <- r^2 * clip + r * s0 * dnorm(sv) / pnorm(sv) + ga
+ return(er)
+ }
+ r <- try(uniroot(solvfct, lower=1e-5, upper = 10)$root,silent=TRUE)
+ if(is(r, "try-error")) return(NA)
+ return(r)
+ })
Added: branches/robast-0.9/pkg/ROptEst/R/getRadius.R
===================================================================
--- branches/robast-0.9/pkg/ROptEst/R/getRadius.R (rev 0)
+++ branches/robast-0.9/pkg/ROptEst/R/getRadius.R 2013-01-27 19:10:04 UTC (rev 572)
@@ -0,0 +1,33 @@
+getRadius <- function(IC, risk, neighbor, L2Fam){
+ if(!is(IC, "HampIC")) stop("'IC' must be of class 'HampIC'.")
+ if(!is(risk,"asGRisk")) stop("'risk' must be of class 'asGRisk'.")
+ if(!is(neighbor,"UncondNeighborhood"))
+ stop("'neighbor' must be of class 'UncondNeighborhood'.")
+ if (missing(L2Fam)) L2Fam <- force(eval(IC at CallL2Fam))
+ if(!is(L2Fam,"L2ParamFamily")) stop("'L2Fam' must be of class 'L2ParamFamily'.")
+ L2deriv.0 <- L2Fam at L2deriv
+ if(numberOfMaps(L2deriv.0)==1){ ## L2derivDim <- L2Fam at L2deriv
+ z <- cent(IC)/as.vector(stand(IC))
+ c0 <- clip(IC)/abs(as.vector(stand(IC)))
+ symm <- FALSE
+ if(is(L2Fam at L2derivDistrSymm[[1]], "SphericalSymmetry"))
+ symm <- L2Fam at L2derivDistrSymm[[1]]@SymmCenter == 0
+ r <- getInfRad(clip = c0, L2deriv = L2Fam at L2derivDistr[[1]],
+ risk = risk, neighbor = neighbor, biastype = biastype(risk),
+ cent = z, symm = symm, trafo = trafo(L2Fam at param))
+ }else{
+ if(!is(L2Fam at distribution, "UnivariateDistribution"))
+ stop("not yet implemented")
+ if((length(L2Fam at L2deriv) == 1) & is(L2Fam at L2deriv[[1]], "RealRandVariable")){
+ L2deriv <- L2Fam at L2deriv[[1]]
+ }else{
+ L2deriv <- diag(dimension(L2Fam at L2deriv)) %*% L2Fam at L2deriv
+ }
+ z <- solve(stand(IC),cent(IC))
+ r <- getInfRad(clip = clip(IC), L2deriv = L2deriv,
+ risk = risk, neighbor = neighbor, biastype = biastype(risk),
+ Distr = L2Fam at distribution, stand = stand(IC),
+ cent = z , trafo = trafo(L2Fam at param))
+ }
+ return(r)
+}
Modified: branches/robast-0.9/pkg/ROptEst/R/getRiskIC.R
===================================================================
--- branches/robast-0.9/pkg/ROptEst/R/getRiskIC.R 2013-01-27 15:54:29 UTC (rev 571)
+++ branches/robast-0.9/pkg/ROptEst/R/getRiskIC.R 2013-01-27 19:10:04 UTC (rev 572)
@@ -16,9 +16,20 @@
L2Fam = "L2ParamFamily"),
function(IC, risk, L2Fam){
Cov <- IC at Risks[["asCov"]]
- if(is.null(Cov))
+ if(is.null(Cov)){
+ if(numberOfMaps(L2Fam at L2deriv)==1){ ## L2derivDim <- L2Fam at L2deriv
+ L2deriv <- L2Fam at L2derivDistr[[1]]
+ A <- as.vector(IC at stand)
+ c0 <- IC at clip/abs(A)
+ z <- IC at cent/A
+ neighbor <- ContNeighborhood(1)
+ Cov <- getInfV(L2deriv = L2deriv, neighbor = neighbor,
+ biastype = biastype(IC), clip = c0, cent = z, stand = A)
+ return(list(asCov = list(distribution = .getDistr(L2Fam),
+ value = Cov)))
+ }
return(getRiskIC(as(IC, "IC"), risk = risk, L2Fam = L2Fam))
- else
+ }else
return(list(asCov = list(distribution = .getDistr(L2Fam), value = Cov)))
})
@@ -30,9 +41,9 @@
Cov <- IC at Risks[["asCov"]]
if (is.null(Cov)){
L2deriv <- L2Fam at L2derivDistr[[1]]
- A <- IC at stand
- c0 <- (IC at clipUp-IC@clipLo)/A
- z <- IC at clipLo/A
+ A <- as.vector(IC at stand)
+ c0 <- (IC at clipUp-IC@clipLo)/abs(A)
+ z <- IC at clipLo/abs(A)
neighbor <- TotalVarNeighborhood(1)
Cov <- getInfV(L2deriv = L2deriv, neighbor = neighbor,
biastype = biastype(IC), clip = c0, cent = z, stand = A)
Deleted: branches/robast-0.9/pkg/ROptEst/inst/scripts/GumbelLocationModel.R
===================================================================
--- branches/robast-0.9/pkg/ROptEst/inst/scripts/GumbelLocationModel.R 2013-01-27 15:54:29 UTC (rev 571)
+++ branches/robast-0.9/pkg/ROptEst/inst/scripts/GumbelLocationModel.R 2013-01-27 19:10:04 UTC (rev 572)
@@ -1,175 +0,0 @@
-###############################################################################
-## Example: Gumbel Location Family
-## computations numerically less stable than in case of the
-## Exponential Scale Family
-###############################################################################
-require(ROptEst)
-options("newDevice"=TRUE)
-
-## generates Gumbel Location Family with loc = 0
-## (known scale = 1)
-distrExOptions(ElowerTruncQuantile = 1e-7) # non-finite function value in integrate
-G0 <- GumbelLocationFamily(loc = 0, scale = 1)
-G0 # show G0
-plot(G0) # plot of Gumbel(loc = 0, scale = 1) and L_2 derivative
-checkL2deriv(G0)
-
-## classical optimal IC
-G0.IC0 <- optIC(model = G0, risk = asCov())
-G0.IC0 # show IC
-plot(G0.IC0) # plot IC
-checkIC(G0.IC0)
-
-## L_2 family + infinitesimal neighborhood
-G0.Rob1 <- InfRobModel(center = G0, neighbor = ContNeighborhood(radius = 0.5))
-G0.Rob1 # show G0.Rob1
-G0.Rob2 <- InfRobModel(center = G0, neighbor = TotalVarNeighborhood(radius = 0.5))
-
-## MSE solution
-E1.Rob1 <- InfRobModel(center = ExpScaleFamily(), neighbor = ContNeighborhood(radius = 0.5))
-(E1.IC1 <- optIC(model=E1.Rob1, risk=asMSE()))
-G0.IC1 <- optIC(model=G0.Rob1, risk=asMSE())
-checkIC(G0.IC1)
-Risks(G0.IC1)
-clip(E1.IC1)
-cent(E1.IC1)
-stand(E1.IC1)
-clip(G0.IC1)
-cent(G0.IC1)
-stand(G0.IC1)
-
-## alternatively
-G0.IC11 <- E1.IC1 # rate = 1!
-CallL2Fam(G0.IC11) <- call("GumbelLocationFamily")
-cent(G0.IC11) <- -cent(E1.IC1)
-G0.IC11
-checkIC(G0.IC11)
-Risks(G0.IC11)
-
-E1.Rob2 <- InfRobModel(center = ExpScaleFamily(), neighbor = TotalVarNeighborhood(radius = 0.5))
-E1.IC2 <- optIC(model=E1.Rob2, risk=asMSE())
-G0.IC2 <- optIC(model=G0.Rob2, risk=asMSE())
-checkIC(G0.IC2)
-Risks(G0.IC2)
-clipLo(E1.IC2)
-clipUp(E1.IC2)
-stand(E1.IC2)
-clipLo(G0.IC2)
-clipUp(G0.IC2)
-stand(G0.IC2)
-## alternatively
-G0.IC21 <- E1.IC2 # rate = 1!
-CallL2Fam(G0.IC21) <- call("GumbelLocationFamily")
-clipLo(G0.IC21) <- -clipUp(E1.IC2)
-clipUp(G0.IC21) <- -clipLo(E1.IC2)
-G0.IC21
-checkIC(G0.IC21)
-Risks(G0.IC21)
-
-## lower case solutions
-(G0.IC3 <- optIC(model=G0.Rob1, risk=asBias()))
-checkIC(G0.IC3)
-Risks(G0.IC3)
-(G0.IC4 <- optIC(model=G0.Rob2, risk=asBias()))
-checkIC(G0.IC4)
-Risks(G0.IC4)
-
-## Hampel solution
-(G0.IC5 <- optIC(model=G0.Rob1, risk=asHampel(bound=clip(G0.IC1))))
-checkIC(G0.IC5)
-Risks(G0.IC5)
-(G0.IC6 <- optIC(model=G0.Rob2, risk=asHampel(bound=Risks(G0.IC2)$asBias$value), maxiter = 100))
-checkIC(G0.IC6)
-Risks(G0.IC6)
-
-## radius minimax IC
-## numerically instable for small 'loRad'!
-## => use connection to ExpScaleFamily for computations
-(G0.IC7 <- radiusMinimaxIC(L2Fam=G0, neighbor=ContNeighborhood(),
- risk=asMSE(), loRad=0.5, upRad=1.0))
-checkIC(G0.IC7)
-Risks(G0.IC7)
-(G0.IC8 <- radiusMinimaxIC(L2Fam=G0, neighbor=TotalVarNeighborhood(),
- risk=asMSE(), loRad=0.5, upRad=1.5))
-checkIC(G0.IC8)
-Risks(G0.IC8)
-
-## least favorable radius
-## numerically instable!
-## => use connection to ExpScaleFamily for computations
-(G0.r.rho1 <- leastFavorableRadius(L2Fam=G0, neighbor=ContNeighborhood(),
- risk=asMSE(), rho=0.5))
-(G0.r.rho2 <- leastFavorableRadius(L2Fam=G0, neighbor=TotalVarNeighborhood(),
- risk=asMSE(), rho=1/3))
-
-## one-step estimation
-## 1. generate a contaminated sample
-ind <- rbinom(1e2, size=1, prob=0.05)
-G0.x <- rgumbel(1e2, loc=(1-ind)*0.5+ind*1)
-
-## 2. Kolmogorov(-Smirnov) minimum distance estimator
-(G0.est01 <- MDEstimator(x=G0.x, GumbelLocationFamily()))
-
-## 3. Cramer von Mises minimum distance estimator
-(G0.est02 <- MDEstimator(x=G0.x, GumbelLocationFamily(), distance = CvMDist))
-
-## 4. one-step estimation: radius known
-G0.Rob3 <- InfRobModel(center=GumbelLocationFamily(loc=estimate(G0.est02)),
- neighbor=ContNeighborhood(radius=0.5))
-G0.IC9 <- optIC(model=G0.Rob3, risk=asMSE())
-(G0.est1 <- oneStepEstimator(G0.x, IC=G0.IC9, start=G0.est02))
-
-## 5. k-step estimation: radius known
-(G0.est2 <- kStepEstimator(G0.x, IC=G0.IC9, start=G0.est02, steps = 3))
-
-## 6. M estimation: radius known
-G0.Rob31 <- InfRobModel(center=GumbelLocationFamily(loc=0),
- neighbor=ContNeighborhood(radius=0.5))
-G0.IC91 <- optIC(model=G0.Rob31, risk=asMSE())
-(G0.est11 <- locMEstimator(G0.x, IC=G0.IC91))
-
-## 7. one-step estimation: radius interval
-G0.IC10 <- radiusMinimaxIC(L2Fam=GumbelLocationFamily(loc=estimate(G0.est02)),
- neighbor=ContNeighborhood(), risk=asMSE(), loRad=0.5, upRad=1)
-(G0.est3 <- oneStepEstimator(G0.x, IC=G0.IC10, start=G0.est02))
-
-## 8. k-step estimation: radius interval
-(G0.est4 <- kStepEstimator(G0.x, IC=G0.IC10, start=G0.est02, steps = 3))
-
-## 9. M estimation: radius interval
-G0.IC101 <- radiusMinimaxIC(L2Fam=GumbelLocationFamily(),
- neighbor=ContNeighborhood(), risk=asMSE(), loRad=0.5, upRad=1)
-(G0.est21 <- locMEstimator(G0.x, IC=G0.IC101))
-
-## 10. It's easier to use function roptest for k-step (k >= 1) estimation!
-sqrtn <- sqrt(length(G0.x))
-G0.est5 <- roptest(G0.x, eps.lower = 0.5/sqrtn, eps.upper = 1/sqrtn,
- L2Fam = GumbelLocationFamily())
-G0.est6 <- roptest(G0.x, eps.lower = 0.5/sqrtn, eps.upper = 1/sqrtn,
- L2Fam = GumbelLocationFamily(), steps = 3)
-
-## comparison - radius known
-estimate(G0.est1)
-estimate(G0.est2)
-estimate(G0.est11)
-
-## confidence intervals
-confint(G0.est1, symmetricBias())
-confint(G0.est2, symmetricBias())
-confint(G0.est11, symmetricBias())
-
-## comparison - radius interval
-estimate(G0.est3)
-estimate(G0.est5)
-estimate(G0.est4)
-estimate(G0.est6)
-estimate(G0.est21)
-
-## confidence intervals
-confint(G0.est3, symmetricBias())
-confint(G0.est5, symmetricBias())
-confint(G0.est4, symmetricBias())
-confint(G0.est6, symmetricBias())
-confint(G0.est21, symmetricBias())
-
-distrExOptions(ElowerTruncQuantile=0) # default
Modified: branches/robast-0.9/pkg/ROptEst/man/getInfClip.Rd
===================================================================
--- branches/robast-0.9/pkg/ROptEst/man/getInfClip.Rd 2013-01-27 15:54:29 UTC (rev 571)
+++ branches/robast-0.9/pkg/ROptEst/man/getInfClip.Rd 2013-01-27 19:10:04 UTC (rev 572)
@@ -45,7 +45,7 @@
risk, neighbor, biastype, cent, symm, trafo)
\S4method{getInfClip}{numeric,UnivariateDistribution,asSemivar,ContNeighborhood}(clip, L2deriv,
- risk, neighbor, cent, symm, trafo)
+ risk, neighbor, biastype, cent, symm, trafo)
}
\arguments{
\item{clip}{ positive real: clipping bound }
Added: branches/robast-0.9/pkg/ROptEst/man/getInfRad.Rd
===================================================================
--- branches/robast-0.9/pkg/ROptEst/man/getInfRad.Rd (rev 0)
+++ branches/robast-0.9/pkg/ROptEst/man/getInfRad.Rd 2013-01-27 19:10:04 UTC (rev 572)
@@ -0,0 +1,128 @@
+\name{getInfRad}
+\alias{getInfRad}
+\alias{getInfRad-methods}
+\alias{getInfRad,numeric,UnivariateDistribution,asMSE,ContNeighborhood-method}
+\alias{getInfRad,numeric,UnivariateDistribution,asMSE,TotalVarNeighborhood-method}
+\alias{getInfRad,numeric,UnivariateDistribution,asL1,ContNeighborhood-method}
+\alias{getInfRad,numeric,UnivariateDistribution,asL1,TotalVarNeighborhood-method}
+\alias{getInfRad,numeric,UnivariateDistribution,asL4,ContNeighborhood-method}
+\alias{getInfRad,numeric,UnivariateDistribution,asL4,TotalVarNeighborhood-method}
+\alias{getInfRad,numeric,EuclRandVariable,asMSE,UncondNeighborhood-method}
+\alias{getInfRad,numeric,UnivariateDistribution,asUnOvShoot,UncondNeighborhood-method}
+\alias{getInfRad,numeric,UnivariateDistribution,asSemivar,ContNeighborhood-method}
+
+\title{Generic Function for the Computation of the Optimal Radius for Given Clipping Bound}
+\description{
+ The usual robust optimality problem for given asGRisk searches the optimal
+ clipping height b of a Hampel-type IC to given radius of the neighborhood.
+ Instead, again for given asGRisk and for given Hampel-Type IC with
+ given clipping height b we may determine the radius of the neighborhood
+ for which it is optimal in the sense of the first sentence. This
+ radius is determined by \code{getInfRad}. This function is rarely called
+ directly. It is used withing \code{\link{getRadius}}.
+}
+\usage{
+getInfRad(clip, L2deriv, risk, neighbor, ...)
+
+\S4method{getInfRad}{numeric,UnivariateDistribution,asMSE,ContNeighborhood}(clip, L2deriv,
+ risk, neighbor, biastype, cent, symm, trafo)
+
+\S4method{getInfRad}{numeric,UnivariateDistribution,asMSE,TotalVarNeighborhood}(clip, L2deriv,
+ risk, neighbor, biastype, cent, symm, trafo)
+
+\S4method{getInfRad}{numeric,UnivariateDistribution,asL1,ContNeighborhood}(clip, L2deriv,
+ risk, neighbor, biastype, cent, symm, trafo)
+
+\S4method{getInfRad}{numeric,UnivariateDistribution,asL1,TotalVarNeighborhood}(clip, L2deriv,
+ risk, neighbor, biastype, cent, symm, trafo)
+
+\S4method{getInfRad}{numeric,UnivariateDistribution,asL4,ContNeighborhood}(clip, L2deriv,
+ risk, neighbor, biastype, cent, symm, trafo)
+
+\S4method{getInfRad}{numeric,UnivariateDistribution,asL4,TotalVarNeighborhood}(clip, L2deriv,
+ risk, neighbor, biastype, cent, symm, trafo)
+
+\S4method{getInfRad}{numeric,EuclRandVariable,asMSE,UncondNeighborhood}(clip, L2deriv, risk,
+ neighbor, biastype, Distr, stand, cent, trafo)
+
+\S4method{getInfRad}{numeric,UnivariateDistribution,asUnOvShoot,UncondNeighborhood}(clip, L2deriv,
+ risk, neighbor, biastype, cent, symm, trafo)
+
+\S4method{getInfRad}{numeric,UnivariateDistribution,asSemivar,ContNeighborhood}(clip, L2deriv,
+ risk, neighbor, biastype, cent, symm, trafo)
+}
+\arguments{
+ \item{clip}{ positive real: clipping bound }
+ \item{L2deriv}{ L2-derivative of some L2-differentiable family
+ of probability measures. }
+ \item{risk}{ object of class \code{"RiskType"}. }
+ \item{neighbor}{ object of class \code{"Neighborhood"}. }
+ \item{\dots}{ additional parameters. }
+ \item{biastype}{ object of class \code{"BiasType"} }
+ \item{cent}{ optimal centering constant. }
+ \item{stand}{ standardizing matrix. }
+ \item{Distr}{ object of class \code{"Distribution"}. }
+ \item{symm}{ logical: indicating symmetry of \code{L2deriv}. }
+ \item{trafo}{ matrix: transformation of the parameter. }
+}
+%\details{}
+\value{The optimal clipping bound is computed.}
+\section{Methods}{
+\describe{
+ \item{clip = "numeric", L2deriv = "UnivariateDistribution",
+ risk = "asMSE", neighbor = "ContNeighborhood"}{
+ optimal clipping bound for asymtotic mean square error. }
+
+
+ \item{clip = "numeric", L2deriv = "UnivariateDistribution",
+ risk = "asMSE", neighbor = "TotalVarNeighborhood"}{
+ optimal clipping bound for asymtotic mean square error. }
+
+ \item{clip = "numeric", L2deriv = "EuclRandVariable",
+ risk = "asMSE", neighbor = "UncondNeighborhood"}{
+ optimal clipping bound for asymtotic mean square error. }
+
+ \item{clip = "numeric", L2deriv = "UnivariateDistribution",
+ risk = "asL1", neighbor = "ContNeighborhood"}{
+ optimal clipping bound for asymtotic mean absolute error. }
+
+ \item{clip = "numeric", L2deriv = "UnivariateDistribution",
+ risk = "asL1", neighbor = "TotalVarNeighborhood"}{
+ optimal clipping bound for asymtotic mean absolute error. }
+
+ \item{clip = "numeric", L2deriv = "UnivariateDistribution",
+ risk = "asL4", neighbor = "ContNeighborhood"}{
+ optimal clipping bound for asymtotic mean power 4 error. }
+
+ \item{clip = "numeric", L2deriv = "UnivariateDistribution",
+ risk = "asL4", neighbor = "TotalVarNeighborhood"}{
+ optimal clipping bound for asymtotic mean power 4 error. }
+
+ \item{clip = "numeric", L2deriv = "UnivariateDistribution",
+ risk = "asUnOvShoot", neighbor = "UncondNeighborhood"}{
+ optimal clipping bound for asymtotic under-/overshoot risk. }
+
+ \item{clip = "numeric", L2deriv = "UnivariateDistribution",
+ risk = "asSemivar", neighbor = "ContNeighborhood"}{
+ optimal clipping bound for asymtotic semivariance.}
+}}
+\references{
+ Rieder, H. (1980) Estimates derived from robust tests. Ann. Stats. \bold{8}: 106--115.
+
+ Rieder, H. (1994) \emph{Robust Asymptotic Statistics}. New York: Springer.
+
+ Ruckdeschel, P. and Rieder, H. (2004) Optimal Influence Curves for
+ General Loss Functions. Statistics & Decisions \emph{22}, 201-223.
+
+ Ruckdeschel, P. (2005) Optimally One-Sided Bounded Influence Curves.
+ Mathematical Methods in Statistics \emph{14}(1), 105-131.
+
+ Kohl, M. (2005) \emph{Numerical Contributions to the Asymptotic Theory of Robustness}.
+ Bayreuth: Dissertation.
+}
+\author{Peter Ruckdeschel \email{Peter.Ruckdeschel at itwm.fraunhofer.de}}
+%\note{}
+\seealso{\code{\link[RobAStBase]{ContIC-class}}, \code{\link[RobAStBase]{TotalVarIC-class}}}
+%\examples{}
+\concept{influence curve}
+\keyword{robust}
Added: branches/robast-0.9/pkg/ROptEst/man/getRadius.Rd
===================================================================
--- branches/robast-0.9/pkg/ROptEst/man/getRadius.Rd (rev 0)
+++ branches/robast-0.9/pkg/ROptEst/man/getRadius.Rd 2013-01-27 19:10:04 UTC (rev 572)
@@ -0,0 +1,56 @@
+\name{getRadius}
+\alias{getRadius}
+
+\title{Computation of the Optimal Radius for Given Clipping Bound}
+\description{
+ The usual robust optimality problem for given asGRisk searches the optimal
+ clipping height b of a Hampel-type IC to given radius of the neighborhood.
+ Instead, again for given asGRisk and for given Hampel-Type IC with
+ given clipping height b we may determine the radius of the neighborhood
+ for which it is optimal in the sense of the first sentence.
+}
+\usage{
+getRadius(IC, risk, neighbor, L2Fam)
+}
+\arguments{
+ \item{IC}{ an IC of class \code{"HampIC"}.}
+ \item{risk}{ object of class \code{"RiskType"}. }
+ \item{neighbor}{ object of class \code{"Neighborhood"}. }
+ \item{L2Fam}{ object of class \code{"L2FamParameter"}. Can be missing;
+ in this case it is constructed from slot \code{CallL2Fam} of
+ \code{IC}. }
+}
+\value{The optimal radius is computed.}
+\references{
+ Rieder, H. (1980) Estimates derived from robust tests. Ann. Stats. \bold{8}: 106--115.
+
+ Rieder, H. (1994) \emph{Robust Asymptotic Statistics}. New York: Springer.
+
+ Ruckdeschel, P. and Rieder, H. (2004) Optimal Influence Curves for
+ General Loss Functions. Statistics & Decisions \emph{22}, 201-223.
+
+ Ruckdeschel, P. (2005) Optimally One-Sided Bounded Influence Curves.
+ Mathematical Methods in Statistics \emph{14}(1), 105-131.
+
+ Kohl, M. (2005) \emph{Numerical Contributions to the Asymptotic Theory of Robustness}.
+ Bayreuth: Dissertation.
+}
+\author{Peter Ruckdeschel \email{Peter.Ruckdeschel at itwm.fraunhofer.de}}
+%\note{}
+\seealso{\code{\link[RobAStBase]{ContIC-class}}, \code{\link[RobAStBase]{TotalVarIC-class}}}
+\examples{
+N <- NormLocationFamily(mean=0, sd=1)
+nb <- ContNeighborhood(); ri <- asMSE()
+radIC <- radiusMinimaxIC(L2Fam=N, neighbor=nb, risk=ri, loRad=0.1, upRad=0.5)
+getRadius(radIC, L2Fam=N, neighbor=nb, risk=ri)
+
+## taken from script NormalScaleModel.R in folder scripts
+N0 <- NormScaleFamily(mean=0, sd=1)
+(N0.IC7 <- radiusMinimaxIC(L2Fam=N0, neighbor=nb, risk=ri, loRad=0, upRad=Inf))
+##
+getRadius(N0.IC7, risk=asMSE(), neighbor=nb, L2Fam=N0)
+getRadius(N0.IC7, risk=asL4(), neighbor=nb, L2Fam=N0)
+}
+
+\concept{influence curve}
+\keyword{robust}
Modified: branches/robast-0.9/pkg/ROptEst/man/getRiskIC.Rd
===================================================================
--- branches/robast-0.9/pkg/ROptEst/man/getRiskIC.Rd 2013-01-27 15:54:29 UTC (rev 571)
+++ branches/robast-0.9/pkg/ROptEst/man/getRiskIC.Rd 2013-01-27 19:10:04 UTC (rev 572)
@@ -58,6 +58,12 @@
\author{Peter Ruckdeschel \email{Peter.Ruckdeschel at itwm.fraunhofer.de}}
\note{This generic function is still under construction.}
\seealso{\code{\link[ROptEst]{getRiskIC}}, \code{\link[RobAStBase]{InfRobModel-class}}}
-%\examples{}
+\examples{
+B <- BinomFamily(size = 25, prob = 0.25)
+
+## classical optimal IC
+IC0 <- optIC(model = B, risk = asCov())
+getRiskIC(IC0, asCov())
+}
\concept{influence curve}
\keyword{robust}
Copied: branches/robast-0.9/pkg/RobExtremes/inst/scripts/GumbelLocationModel.R (from rev 565, branches/robast-0.9/pkg/ROptEst/inst/scripts/GumbelLocationModel.R)
===================================================================
--- branches/robast-0.9/pkg/RobExtremes/inst/scripts/GumbelLocationModel.R (rev 0)
+++ branches/robast-0.9/pkg/RobExtremes/inst/scripts/GumbelLocationModel.R 2013-01-27 19:10:04 UTC (rev 572)
@@ -0,0 +1,175 @@
+###############################################################################
+## Example: Gumbel Location Family
+## computations numerically less stable than in case of the
+## Exponential Scale Family
+###############################################################################
+require(RobExtremes)
+options("newDevice"=TRUE)
+
+## generates Gumbel Location Family with loc = 0
+## (known scale = 1)
+distrExOptions(ElowerTruncQuantile = 1e-7) # non-finite function value in integrate
+G0 <- GumbelLocationFamily(loc = 0, scale = 1)
+G0 # show G0
+plot(G0) # plot of Gumbel(loc = 0, scale = 1) and L_2 derivative
+checkL2deriv(G0)
+
+## classical optimal IC
+G0.IC0 <- optIC(model = G0, risk = asCov())
+G0.IC0 # show IC
+plot(G0.IC0) # plot IC
[TRUNCATED]
To get the complete diff run:
svnlook diff /svnroot/robast -r 572
More information about the Robast-commits
mailing list