[Robast-commits] r1131 - in branches/robast-1.2/pkg/ROptEst: . R inst inst/scripts man
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Sun Aug 12 10:54:13 CEST 2018
Author: ruckdeschel
Date: 2018-08-12 10:54:12 +0200 (Sun, 12 Aug 2018)
New Revision: 1131
Added:
branches/robast-1.2/pkg/ROptEst/R/CheckMakeContIC.R
branches/robast-1.2/pkg/ROptEst/man/checkmakeIC.Rd
Modified:
branches/robast-1.2/pkg/ROptEst/NAMESPACE
branches/robast-1.2/pkg/ROptEst/R/L1L2normL2deriv.R
branches/robast-1.2/pkg/ROptEst/R/LowerCaseMultivariate.R
branches/robast-1.2/pkg/ROptEst/R/getAsRisk.R
branches/robast-1.2/pkg/ROptEst/R/getComp.R
branches/robast-1.2/pkg/ROptEst/R/getInfCent.R
branches/robast-1.2/pkg/ROptEst/R/getInfClip.R
branches/robast-1.2/pkg/ROptEst/R/getInfGamma.R
branches/robast-1.2/pkg/ROptEst/R/getInfLM.R
branches/robast-1.2/pkg/ROptEst/R/getInfRad.R
branches/robast-1.2/pkg/ROptEst/R/getInfRobIC_asAnscombe.R
branches/robast-1.2/pkg/ROptEst/R/getInfRobIC_asBias.R
branches/robast-1.2/pkg/ROptEst/R/getInfRobIC_asGRisk.R
branches/robast-1.2/pkg/ROptEst/R/getInfRobIC_asHampel.R
branches/robast-1.2/pkg/ROptEst/R/getInfStand.R
branches/robast-1.2/pkg/ROptEst/R/getInfV.R
branches/robast-1.2/pkg/ROptEst/R/getMaxIneff.R
branches/robast-1.2/pkg/ROptEst/R/getModifyIC.R
branches/robast-1.2/pkg/ROptEst/R/getReq.R
branches/robast-1.2/pkg/ROptEst/R/getRiskIC.R
branches/robast-1.2/pkg/ROptEst/R/getStartIClcsc.R
branches/robast-1.2/pkg/ROptEst/R/internal.roptest.R
branches/robast-1.2/pkg/ROptEst/R/leastFavorableRadius.R
branches/robast-1.2/pkg/ROptEst/R/optIC.R
branches/robast-1.2/pkg/ROptEst/R/radiusMinimaxIC.R
branches/robast-1.2/pkg/ROptEst/R/roptest.new.R
branches/robast-1.2/pkg/ROptEst/inst/NEWS
branches/robast-1.2/pkg/ROptEst/inst/scripts/MBRE.R
branches/robast-1.2/pkg/ROptEst/man/getBiasIC.Rd
branches/robast-1.2/pkg/ROptEst/man/getInfCent.Rd
branches/robast-1.2/pkg/ROptEst/man/getInfClip.Rd
branches/robast-1.2/pkg/ROptEst/man/getInfGamma.Rd
branches/robast-1.2/pkg/ROptEst/man/getInfRad.Rd
branches/robast-1.2/pkg/ROptEst/man/getInfStand.Rd
branches/robast-1.2/pkg/ROptEst/man/getInfV.Rd
branches/robast-1.2/pkg/ROptEst/man/getMaxIneff.Rd
branches/robast-1.2/pkg/ROptEst/man/getReq.Rd
branches/robast-1.2/pkg/ROptEst/man/getRiskIC.Rd
branches/robast-1.2/pkg/ROptEst/man/getinfLM.Rd
branches/robast-1.2/pkg/ROptEst/man/inputGenerator.Rd
branches/robast-1.2/pkg/ROptEst/man/internals.Rd
branches/robast-1.2/pkg/ROptEst/man/leastFavorableRadius.Rd
branches/robast-1.2/pkg/ROptEst/man/minmaxBias.Rd
branches/robast-1.2/pkg/ROptEst/man/optIC.Rd
branches/robast-1.2/pkg/ROptEst/man/robest.Rd
branches/robast-1.2/pkg/ROptEst/man/roptest.Rd
Log:
[ROptEst] branch 1.2
+ getBiasIC for signature {HampIC,UncondNeighborhood} no longer has argument withCheck
+ getRiskIC gains argument "..." to pass on arguments to E()
+ roptest gains an argument E.argList for arguments to be passed to E() from (a)
\code{MDEstimator} (here this additional argument is only used if \code{initial.est}
is missing), (b) \code{getStartIC}, and (c) \code{kStepEstimator}. Potential clashes with
arguments of the same name in \code{\dots} are resolved by inserting the items of argument
list \code{E.argList} as named items, so in case of collisions the item of \code{E.argList}
overwrites the existing one from \code{\dots}.
+ the input generators genkStepCtrl, genstartCtrl, genstartICCtrl gain argument
E.argList to pass on arguments to E()
+ the help to robest is more explicit about the usage of argument E.argList at
various places in the generator
+ the help to roptest has an extra paragraph explicating the usage of initial.est
and start.Par and the different steps done in roptest
+ getBiasIC for signature {HampIC,UncondNeighborhood} no longer has argument withCheck
+ several getRiskIC and optIC methods gain argument "..." to pass on arguments to E()
+ revised script MBRE.R
----------------------------------------------------------------------------------
+ an embarrassing error when computing E(X t(X)) (saving the lower half by symmetry)
in getInfStand: a[row(a)>col(a)] <- a[row(a)<col(a)] in general is not symmetric...
did not bite us so far as the highest dim we considered so far was 3, and the
recipe is exact exactly up to dim 3... the exact formula would have been (and is now
in the code) a[row(a)>col(a)] <- t(a)[row(a)>col(a)]
+ fix an argument collision in MBRE-method / in call to getInfRobIC_asGRisk, when called
from OptIC through getInfRobIC_asBias ("upper" got into ...)
----------------------------------------------------------------------------------
+ internal function .getComp, determining by symmetry slots which entries in LMs a and A
have to be computed, now fills the lower triangle of A with FALSE (was not used so far,
but can be used in a faster computation method for checkIC makeIC to determine whether
it is cleverer to integrate in k or in p space)
+ particular checkIC and makeIC methods for ContICs which allow for speed up
if in k space many entries of the LMs can be skipped due to symmetry
+ several methods (getRiskIC, getBiasIC, getBoundedIC, makeIC, checkIC, modifyIC)
gain argument "..." to pass on arguments to E(). This holds in particular for the
functions used to compute the optimally-robust ICs, i.e. getInfRob_asBias, getInfRob_asHampel,
getInfRob_asGRisk, getInfRob_asAnscombe, leastFavorableRadius, radiusMinimaxIC, and
geInfGamma, getInfRad, getInfClip, getInfCent, getInfStand, getInfV, getAsRisk, getReq,
(at least as far as multivariate ICs are concerned), .LowerCaseMultivariate, getMaxIneff,
getInfRad, getLagrangeMultByIter and getLagrangeMultByOptim
+ new internal helper function .filterEarg by means of constant ..IntegrateArgs (exported
from package RobAStBase) is used to filter out arguments from dots which are meant for E()
+ the local .modifyIC0 functions only used to produce the new IC but not for filling
slot modifyIC loose argument withMakeIC (and dots) -- this is now done in the outer
modifyIC function
Modified: branches/robast-1.2/pkg/ROptEst/NAMESPACE
===================================================================
--- branches/robast-1.2/pkg/ROptEst/NAMESPACE 2018-08-12 08:15:22 UTC (rev 1130)
+++ branches/robast-1.2/pkg/ROptEst/NAMESPACE 2018-08-12 08:54:12 UTC (rev 1131)
@@ -39,7 +39,8 @@
"getModifyIC")
exportMethods("updateNorm", "scaleUpdateIC", "eff",
"get.asGRisk.fct", "getStartIC", "plot",
- "comparePlot", "getRiskFctBV", "roptestCall")
+ "comparePlot", "getRiskFctBV", "roptestCall",
+ "checkIC", "makeIC")
export("getL2normL2deriv",
"asAnscombe", "asL1", "asL4",
"getReq", "getMaxIneff", "getRadius")
Added: branches/robast-1.2/pkg/ROptEst/R/CheckMakeContIC.R
===================================================================
--- branches/robast-1.2/pkg/ROptEst/R/CheckMakeContIC.R (rev 0)
+++ branches/robast-1.2/pkg/ROptEst/R/CheckMakeContIC.R 2018-08-12 08:54:12 UTC (rev 1131)
@@ -0,0 +1,209 @@
+#if(FALSE){
+## faster check for ContICs
+
+setMethod("checkIC", signature(IC = "ContIC", L2Fam = "L2ParamFamily"),
+ function(IC, L2Fam, out = TRUE, ...){
+
+ D1 <- L2Fam at distribution
+ if( dimension(Domain(IC at Curve[[1]])) != dimension(img(D1)))
+ stop("dimension of 'Domain' of 'Curve' != dimension of 'img' of 'distribution' of 'L2Fam'")
+
+ res <- .prepareCheckMakeIC(L2Fam, w = IC at weight, ...)
+ ## if it pays off to use symmetry/ to compute integrals in L2deriv space
+ ## we compute the following integrals:
+ ## G1 = E w, G2 = E Lambda w, G3 = E Lambda Lambda' w
+ ## we want to compute:
+ ## Delta1 = E (A Lambda-a) w, Delta2 = E (A Lambda-a) Lambda' w
+ ## where A = stand(IC), a=cent(IC)
+ ## hence Delta1 = A G2 - a G1, Delta2 = A G3 - a G2'
+ ### otherwise the return value is NULL and we use the standard method
+
+ if(is.null(res))
+ return(getMethod("checkIC", signature(IC = "IC",
+ L2Fam = "L2ParamFamily"))(IC,L2Fam, out = out, ...))
+
+
+ A <- stand(IC); a <- cent(IC)
+ G1 <- res$G1; G2 <- res$G2; G3 <- res$G3
+ Delta1 <- A%*%G2- a*G1
+ Delta2 <- A%*%G3 - a%*%t(G2) - trafo(L2Fam at param)
+
+ if(out)
+ cat("precision of centering:\t", Delta1, "\n")
+
+ if(out){
+ cat("precision of Fisher consistency:\n")
+ print(Delta2)
+ cat("precision of Fisher consistency - relative error [%]:\n")
+ print(100*Delta2/trafo)
+ }
+
+ prec <- max(abs(Delta1), abs(Delta2))
+ names(prec) <- "maximum deviation"
+
+ return(prec)
+ })
+
+## make some L2function a pIC at a model
+setMethod("makeIC", signature(IC = "ContIC", L2Fam = "L2ParamFamily"),
+ function(IC, L2Fam, ...){
+
+ D1 <- L2Fam at distribution
+ if( dimension(Domain(IC at Curve[[1]])) != dimension(img(D1)))
+ stop("dimension of 'Domain' of 'Curve' != dimension of 'img' of 'distribution' of 'L2Fam'")
+
+ if(dimension(IC at Curve) != dims)
+ stop("Dimension of IC and parameter must be equal")
+
+ res <- .prepareCheckMakeIC(L2Fam, w = IC at weight, ...)
+
+ ## if it pays off to use symmetry/ to compute integrals in L2deriv space
+ ## we compute the following integrals:
+ ## G1 = E w, G2 = E Lambda w, G3 = E Lambda Lambda' w
+ ## we want to compute:
+ ## Delta1 = E (A Lambda-a) w, Delta2 = E (A Lambda-a) Lambda' w
+ ## where A = stand(IC), a=cent(IC)
+ ## hence Delta1 = A G2 - a G1, Delta2 = A G3 - a G2'
+ ### otherwise the return value is NULL and we use the standard method
+
+ if(is.null(res))
+ return(getMethod("makeIC", signature(IC = "IC",
+ L2Fam = "L2ParamFamily"))(IC,L2Fam))
+
+ G1 <- res$G1; G2 <- res$G2; G3 <- res$G3
+ trafo <- trafo(L2Fam at param)
+ nrvalues <- nrow(trafo)
+ dims <- ncol(trafo)
+
+ cent0 <- G2/G1
+ stand1 <- trafo%*%distr::solve(G3-cent0%*%t(G2))
+ cent1 <- stand1%*%cent0
+
+ L2deriv <- as(diag(dims) %*% L2Fam at L2deriv, "EuclRandVariable")
+ D1 <- L2Fam at distribution
+
+ IC1.0 <- stand1%*%L2deriv
+ IC1.1 <- IC1.0 -cent1
+ IC1.f <- function(x) evalRandVar(IC1.1,x)
+
+ IC1.l <- vector("list",nrvalues)
+ for(i in 1:nrvalues){
+ IC1.l[[i]] <- function(x){}
+ body(IC.l[[i]]) <- substitute({indS <- liesInSupport(D0,x,checkFin=TRUE)
+ indS*((IC1.s(x))[i])
+ }, list(IC1.s=IC1.f, D0=D1, i=i))
+ }
+ IC1.c <- EuclRandVariable(Map = IC1.l, Domain = IC at Curve[[1]],
+ Range = Reals())
+
+ cIC1 <- new("ContIC")
+ cIC1 at name <- name
+ cIC1 at Curve <- IC1.c
+ cIC1 at Risks <- IC at Risks
+ cIC1 at Infos <- IC at Infos
+ cIC1 at CallL2Fam <- L2Fam at fam.call
+ cIC1 at clip <- IC at clip
+ cIC1 at cent <- cent1
+ cIC1 at stand <- stand1
+ cIC1 at lowerCase <- IC at lowerCase
+ cIC1 at neighborRadius <- IC at neighborRadius
+ cIC1 at weight <- IC at weight
+ cIC1 at biastype <- IC at biastype
+ cIC1 at normtype <- IC at normtype
+ cIC1 at modifyIC <- IC at modifyIC
+ addInfo(cIC1) <- c("IC<-",
+ "generated by affine linear trafo to enforce consistency")
+ return(cIC1)
+ })
+
+.prepareCheckMakeIC <- function(L2Fam, w, ...){
+
+ dims <- length(L2Fam at param)
+ trafo <- trafo(L2Fam at param)
+ nrvalues <- nrow(trafo)
+
+ z.comp <- rep(TRUE,dims)
+ A.comp <- matrix(rep(TRUE,dims^2),nrow=dims)
+ to.comp.i <- (dims+1)*(dims+2)/2
+ to.comp.a <- (dims+1)*nrvalues
+
+ L2deriv <- as(diag(dims) %*% L2Fam at L2deriv, "EuclRandVariable")
+
+ z.comp <- rep(TRUE,dims)
+ A.comp <- matrix(TRUE, dims, dims)
+ # otherwise if trafo == unitMatrix may use symmetry info
+ if(.isUnitMatrix(trafo)){
+ comp <- .getComp(L2deriv, L2Fam at distrSymm, L2Fam at L2derivSymm, L2Fam at L2derivDistrSymm)
+ z.comp <- comp$"z.comp"
+ A.comp <- comp$"A.comp"
+ t.comp.i <- sum(z.comp)+sum(A.comp)+1
+ }
+
+ if(to.comp.a < to.comp.i) return(NULL)
+
+
+ res <- .getG1G2G3Stand(L2deriv = L2deriv, Distr = L2Fam at distribution,
+ A.comp = A.comp, z.comp = z.comp, w = w, ...)
+ return(res)
+}
+
+
+
+.getG1G2G3Stand <- function(L2deriv, Distr, A.comp, z.comp, w, ...){
+
+ dotsI <- .filterEargs(list(...))
+ if(is.null(dotsI$useApply)) dotsI$useApply <- FALSE
+
+ w.fct <- function(x){
+ weight(w)(evalRandVar(L2deriv, as.matrix(x)) [,,1])
+ }
+
+
+ integrand2 <- function(x, L2.i){
+ return(L2.i(x)*w.fct(x))
+ }
+
+ Eargs <- c(list(object = Distr, fun = w.fct), dotsI)
+ res1 <- do.call(E,Eargs)
+
+ nrvalues <- length(L2deriv)
+ res2 <- numeric(nrvalues)
+ for(i in 1:nrvalues){
+ if(z.comp[i]){
+ Eargs <- c(list(object = Distr, fun = integrand2,
+ L2.i = L2deriv at Map[[i]]), dotsI)
+ res2[i] <- do.call(E,Eargs)
+ }else{
+ res2[i] <- 0
+ }
+ }
+
+ cent <- res2/res1
+
+ integrandA <- function(x, L2.i, L2.j, i, j){
+ return((L2.i(x) - cent[i])*(L2.j(x) - cent[j])*w.fct(x = x))
+ }
+
+ nrvalues <- length(L2deriv)
+ erg <- matrix(0, ncol = nrvalues, nrow = nrvalues)
+
+ for(i in 1:nrvalues){
+ for(j in i:nrvalues){
+ if(A.comp[i,j]){
+ Eargs <- c(list(object = Distr, fun = integrandA,
+ L2.i = L2deriv at Map[[i]],
+ L2.j = L2deriv at Map[[j]], i = i, j = j), dotsI)
+ erg[i, j] <- do.call(E,Eargs)
+ }
+ }
+ }
+ erg[col(erg) < row(erg)] <- t(erg)[col(erg) < row(erg)]
+
+ return(list(G1=res1,G2=res2, G3=erg))
+ }
+
+
+
+
+
+#}
Modified: branches/robast-1.2/pkg/ROptEst/R/L1L2normL2deriv.R
===================================================================
--- branches/robast-1.2/pkg/ROptEst/R/L1L2normL2deriv.R 2018-08-12 08:15:22 UTC (rev 1130)
+++ branches/robast-1.2/pkg/ROptEst/R/L1L2normL2deriv.R 2018-08-12 08:54:12 UTC (rev 1131)
@@ -8,6 +8,10 @@
setMethod("getL1normL2deriv", signature(L2deriv = "RealRandVariable"),
function(L2deriv, cent, stand, Distr, normtype, ...){
+
+ dotsI <- .filterEargs(list(...))
+ if(is.null(dotsI$useApply)) dotsI$useApply <- FALSE
+
integrandG <- function(x, L2, stand, cent){
X <- evalRandVar(L2, as.matrix(x))[,,1] - cent
Y <- apply(X, 2, "%*%", t(stand))
@@ -15,6 +19,7 @@
return((res > 0)*res)
}
- return(E(object = Distr, fun = integrandG, L2 = L2deriv,
- stand = stand, cent = cent, useApply = FALSE))
+
+ return(do.call(E, c(list(object = Distr, fun = integrandG, L2 = L2deriv,
+ stand = stand, cent = cent),dotsI)))
})
Modified: branches/robast-1.2/pkg/ROptEst/R/LowerCaseMultivariate.R
===================================================================
--- branches/robast-1.2/pkg/ROptEst/R/LowerCaseMultivariate.R 2018-08-12 08:15:22 UTC (rev 1130)
+++ branches/robast-1.2/pkg/ROptEst/R/LowerCaseMultivariate.R 2018-08-12 08:54:12 UTC (rev 1131)
@@ -1,8 +1,11 @@
.LowerCaseMultivariate <- function(L2deriv, neighbor, biastype,
normtype, Distr, Finfo, trafo, z.start = NULL,
A.start = NULL, z.comp = NULL, A.comp = NULL, maxiter, tol,
- verbose = NULL){
+ verbose = NULL, ...){
+ dotsI <- .filterEargs(list(...))
+ if(is.null(dotsI$useApply)) dotsI$useApply <- FALSE
+
if(missing(verbose)|| is.null(verbose))
verbose <- getRobAStBaseOption("all.verbose")
@@ -59,8 +62,8 @@
w <<- w0
}
- E1 <- E(object = Distr, fun = abs.fct, L2 = L2deriv, stand = A,
- cent = z, normtype.0 = normtype, useApply = FALSE)
+ E1 <- do.call(E,c(list(object = Distr, fun = abs.fct, L2 = L2deriv, stand = A,
+ cent = z, normtype.0 = normtype), dotsI))
stA <- if (is(normtype,"QFNorm"))
QuadForm(normtype)%*%A else A
# erg <- E1/sum(diag(stA %*% t(trafo)))
@@ -101,8 +104,11 @@
.LowerCaseMultivariateTV <- function(L2deriv, neighbor, biastype,
normtype, Distr, Finfo, trafo,
A.start = NULL, maxiter, tol,
- verbose = NULL){
+ verbose = NULL, ...){
+ dotsI <- .filterEargs(list(...))
+ if(is.null(dotsI$useApply)) dotsI$useApply <- FALSE
+
if(missing(verbose)|| is.null(verbose))
verbose <- getRobAStBaseOption("all.verbose")
@@ -124,8 +130,8 @@
p <- 1
A <- matrix(param, ncol = k, nrow = 1)
# print(A)
- E1 <- E(object = Distr, fun = pos.fct, L2 = L2deriv, stand = A,
- useApply = FALSE)
+ E1 <- do.call(E, c(list( object = Distr, fun = pos.fct,
+ L2 = L2deriv, stand = A), dotsI))
erg <- E1/sum(diag(A %*% t(trafo)))
return(erg)
}
@@ -144,10 +150,10 @@
Y <- as.numeric(A %*% X)
return(as.numeric(pr.sign*Y>0))
}
- p.p <- E(object = Distr, fun = pr.fct, L2 = L2deriv,
- useApply = FALSE, pr.sign = 1)
- m.p <- E(object = Distr, fun = pr.fct, L2 = L2deriv,
- useApply = FALSE, pr.sign = -1)
+ p.p <- do.call(E, c(list( object = Distr, fun = pr.fct, L2 = L2deriv,
+ pr.sign = 1), dotsI))
+ m.p <- do.call(E, c(list( object = Distr, fun = pr.fct, L2 = L2deriv,
+ pr.sign = -1), dotsI))
a <- -b * p.p/(p.p+m.p)
Modified: branches/robast-1.2/pkg/ROptEst/R/getAsRisk.R
===================================================================
--- branches/robast-1.2/pkg/ROptEst/R/getAsRisk.R 2018-08-12 08:15:22 UTC (rev 1130)
+++ branches/robast-1.2/pkg/ROptEst/R/getAsRisk.R 2018-08-12 08:54:12 UTC (rev 1131)
@@ -86,7 +86,7 @@
eerg <- .LowerCaseMultivariate(L2deriv = L2deriv, neighbor = neighbor,
biastype = biastype, normtype = normtype, Distr = Distr, Finfo = Finfo,
trafo = trafo, z.start = z.start, A.start = A.start, z.comp = z.comp,
- A.comp = DA.comp, maxiter = maxiter, tol = tol, verbose = verbose)
+ A.comp = DA.comp, maxiter = maxiter, tol = tol, verbose = verbose, ...)
erg <- eerg$erg
bias <- 1/erg$value
@@ -112,7 +112,7 @@
neighbor = neighbor, biastype = biastype,
normtype = normtype, Distr = Distr, Finfo = Finfo, trafo = trafo,
A.start = A.start, maxiter = maxiter,
- tol = tol, verbose = verbose)
+ tol = tol, verbose = verbose, ...)
erg <- eerg$b
bias <- 1/erg$value
@@ -167,7 +167,7 @@
Cov <- getInfV(L2deriv = L2deriv, neighbor = neighbor,
biastype = biastype, Distr = Distr,
V.comp = V.comp, cent = cent,
- stand = stand, w = w)
+ stand = stand, w = w, ...)
if(!is.null(trafo)) Cov <- trafo%*%Cov%*%t(trafo)
return(list(asCov = Cov))
})
@@ -219,7 +219,7 @@
Cov <- getAsRisk(risk = asCov(), L2deriv = L2deriv, neighbor = neighbor,
biastype = biastype, Distr = Distr, clip = clip,
cent = cent, stand = stand, trafo = trafo,
- V.comp = V.comp, w = w)$asCov
+ V.comp = V.comp, w = w, ...)$asCov
p <- nrow(stand)
std <- if(is(normtype,"QFNorm")) QuadForm(normtype) else diag(p)
@@ -255,7 +255,7 @@
biastype = biastype, normtype = normtype,
Distr = Distr, clip = clip,
cent = cent, stand = stand, V.comp = V.comp,
- w = w)$trAsCov
+ w = w, ...)$trAsCov
return(list(asAnscombe = FI/trAsCov.0))
})
Modified: branches/robast-1.2/pkg/ROptEst/R/getComp.R
===================================================================
--- branches/robast-1.2/pkg/ROptEst/R/getComp.R 2018-08-12 08:15:22 UTC (rev 1130)
+++ branches/robast-1.2/pkg/ROptEst/R/getComp.R 2018-08-12 08:54:12 UTC (rev 1131)
@@ -8,6 +8,7 @@
if(L2derivDistrSymm[[i]]@SymmCenter == 0)
z.comp[i] <- FALSE
}
+ if(nrvalues>1){
for(i in 1:(nrvalues-1))
for(j in (i+1):nrvalues){
if(is(DistrSymm, "SphericalSymmetry")){
@@ -19,6 +20,7 @@
}
}
A.comp[col(A.comp) < row(A.comp)] <- FALSE
+ }
return(list(A.comp = A.comp, z.comp = z.comp))
}
Modified: branches/robast-1.2/pkg/ROptEst/R/getInfCent.R
===================================================================
--- branches/robast-1.2/pkg/ROptEst/R/getInfCent.R 2018-08-12 08:15:22 UTC (rev 1130)
+++ branches/robast-1.2/pkg/ROptEst/R/getInfCent.R 2018-08-12 08:54:12 UTC (rev 1131)
@@ -39,7 +39,11 @@
neighbor = "TotalVarNeighborhood",
biastype = "BiasType"),
function(L2deriv, neighbor, biastype, Distr, z.comp, w,
- tol.z = .Machine$double.eps^.5){
+ tol.z = .Machine$double.eps^.5, ...){
+
+ dotsI <- .filterEargs(list(...))
+ if(is.null(dotsI$useApply)) dotsI$useApply <- FALSE
+
stand <- stand(w)
clip <- clip(w)
b <- clip[2]-clip[1]
@@ -51,7 +55,7 @@
Y <- as.numeric(stand%*%Lx)
pmin(pmax(g,Y),g+c0)
}
- return(E(object = Distr, fun = fct, useApply = FALSE))
+ return(do.call(E, c(list(object = Distr, fun = fct), dotsI)))
}
lower <- -b
upper <- 0
@@ -64,7 +68,11 @@
neighbor = "ContNeighborhood",
biastype = "BiasType"),
function(L2deriv, neighbor, biastype, Distr, z.comp, w,
- tol.z = .Machine$double.eps^.5){
+ tol.z = .Machine$double.eps^.5, ...){
+
+ dotsI <- .filterEargs(list(...))
+ if(is.null(dotsI$useApply)) dotsI$useApply <- FALSE
+
integrand1 <- function(x){
weight(w)(evalRandVar(L2deriv, as.matrix(x)) [,,1])
}
@@ -72,13 +80,13 @@
return(L2.i(x)*integrand1(x))
}
- res1 <- E(object = Distr, fun = integrand1, useApply = FALSE)
+ res1 <- do.call(E, c(list(object = Distr, fun = integrand1),dotsI))
nrvalues <- length(L2deriv)
res2 <- numeric(nrvalues)
for(i in 1:nrvalues){
if(z.comp[i]){
- res2[i] <- E(object = Distr, fun = integrand2,
- L2.i = L2deriv at Map[[i]], useApply = FALSE)
+ res2[i] <- do.call(E, c(list(object = Distr, fun = integrand2,
+ L2.i = L2deriv at Map[[i]]), dotsI))
}else{
res2[i] <- 0
}
Modified: branches/robast-1.2/pkg/ROptEst/R/getInfClip.R
===================================================================
--- branches/robast-1.2/pkg/ROptEst/R/getInfClip.R 2018-08-12 08:15:22 UTC (rev 1130)
+++ branches/robast-1.2/pkg/ROptEst/R/getInfClip.R 2018-08-12 08:54:12 UTC (rev 1131)
@@ -35,11 +35,11 @@
risk = "asMSE",
neighbor = "UncondNeighborhood"),
function(clip, L2deriv, risk, neighbor, biastype,
- Distr, stand, cent, trafo){
+ Distr, stand, cent, trafo, ...){
return(neighbor at radius^2*clip +
getInfGamma(L2deriv = L2deriv, risk = risk, neighbor = neighbor,
biastype = biastype, Distr = Distr, stand = stand,
- cent = cent, clip = clip))
+ cent = cent, clip = clip, ...))
})
###############################################################################
@@ -155,7 +155,11 @@
L2deriv = "UnivariateDistribution",
risk = "asSemivar",
neighbor = "ContNeighborhood"),
- function(clip, L2deriv, risk, neighbor, biastype, cent, symm, trafo){
+ function(clip, L2deriv, risk, neighbor, biastype, cent, symm, trafo, ...){
+
+ dotsI <- .filterEargs(list(...))
+ if(is.null(dotsI$useApply)) dotsI$useApply <- FALSE
+
biastype <- if(sign(risk)==1) positiveBias() else negativeBias()
z0 <- getInfCent(L2deriv = L2deriv, risk = risk, neighbor = neighbor,
biastype = biastype,
@@ -168,9 +172,9 @@
r <- neighbor at radius
if (sign(risk)>0)
- v0 <- E(L2deriv, function(x) pmin( x-z0, clip)^2 )
+ v0 <- do.call(E,c(list(L2deriv, function(x) pmin( x-z0, clip)^2),dotsI))
else
- v0 <- E(L2deriv, function(x) pmax( x-z0, -clip)^2 )
+ v0 <- do.call(E,c(list(L2deriv, function(x) pmax( x-z0, -clip)^2),dotsI))
s0 <- sqrt(v0)
sv <- r * clip / s0
Modified: branches/robast-1.2/pkg/ROptEst/R/getInfGamma.R
===================================================================
--- branches/robast-1.2/pkg/ROptEst/R/getInfGamma.R 2018-08-12 08:15:22 UTC (rev 1130)
+++ branches/robast-1.2/pkg/ROptEst/R/getInfGamma.R 2018-08-12 08:54:12 UTC (rev 1131)
@@ -29,7 +29,11 @@
neighbor = "ContNeighborhood",
biastype = "BiasType"),
function(L2deriv, risk, neighbor, biastype, Distr,
- stand, cent, clip, power = 1L){
+ stand, cent, clip, power = 1L, ...){
+
+ dotsI <- .filterEargs(list(...))
+ if(is.null(dotsI$useApply)) dotsI$useApply <- FALSE
+
integrandG <- function(x, L2, stand, cent, clip){
X <- evalRandVar(L2, as.matrix(x))[,,1] - cent
Y <- stand %*% X
@@ -38,8 +42,9 @@
return((res > 0)*res^power)
}
- return(-E(object = Distr, fun = integrandG, L2 = L2deriv,
- stand = stand, cent = cent, clip = clip, useApply = FALSE))
+ res <- do.call(E, c(list(object = Distr, fun = integrandG, L2 = L2deriv,
+ stand = stand, cent = cent, clip = clip),dotsI))
+ return(-res)
})
setMethod("getInfGamma", signature(L2deriv = "RealRandVariable",
@@ -47,7 +52,11 @@
neighbor = "TotalVarNeighborhood",
biastype = "BiasType"),
function(L2deriv, risk, neighbor, biastype, Distr,
- stand, cent, clip, power = 1L){
+ stand, cent, clip, power = 1L, ...){
+
+ dotsI <- .filterEargs(list(...))
+ if(is.null(dotsI$useApply)) dotsI$useApply <- FALSE
+
integrandG <- function(x, L2, stand, cent, clip){
X <- evalRandVar(L2, as.matrix(x))[,,1] - cent
Y <- stand %*% X
@@ -56,8 +65,9 @@
return((res > 0)*res^power)
}
- return(-E(object = Distr, fun = integrandG, L2 = L2deriv,
- stand = stand, cent = cent, clip = clip, useApply = FALSE))
+ res <- do.call(E, c(list(object = Distr, fun = integrandG, L2 = L2deriv,
+ stand = stand, cent = cent, clip = clip),dotsI))
+ return(-res)
})
###############################################################################
## gamma in case of asymptotic under-/overshoot risk
Modified: branches/robast-1.2/pkg/ROptEst/R/getInfLM.R
===================================================================
--- branches/robast-1.2/pkg/ROptEst/R/getInfLM.R 2018-08-12 08:15:22 UTC (rev 1130)
+++ branches/robast-1.2/pkg/ROptEst/R/getInfLM.R 2018-08-12 08:54:12 UTC (rev 1131)
@@ -6,7 +6,7 @@
neighbor, biastype, normtype, Distr,
a.start, z.start, A.start, w.start, std,
z.comp, A.comp, maxiter, tol,
- verbose = NULL, warnit = TRUE){
+ verbose = NULL, warnit = TRUE, ...){
if(missing(verbose)|| is.null(verbose))
verbose <- getRobAStBaseOption("all.verbose")
if(missing(warnit)|| is.null(warnit)) warnit <- TRUE
@@ -43,7 +43,7 @@
## update centering
z <- getInfCent(L2deriv = L2deriv, neighbor = neighbor,
biastype = biastype, Distr = Distr, z.comp = z.comp,
- w = w, tol.z = .Machine$double.eps^.5)
+ w = w, tol.z = .Machine$double.eps^.5, ...)
# print(c("z"=z))
if(is(neighbor,"TotalVarNeighborhood")){
a <- z
@@ -56,7 +56,7 @@
# update standardization
A <- getInfStand(L2deriv = L2deriv, neighbor = neighbor,
biastype = biastype, Distr = Distr, A.comp = A.comp,
- cent = zc, trafo = trafo, w = w)
+ cent = zc, trafo = trafo, w = w, ...)
# print(c("A"=A))
## in case of self-standardization: update norm
@@ -106,11 +106,16 @@
a.start, z.start, A.start, w.start, std, z.comp,
A.comp, maxiter, tol, verbose = NULL, ...){
+
if(missing(verbose)|| is.null(verbose))
verbose <- getRobAStBaseOption("all.verbose")
LMcall <- match.call()
### manipulate dots in call -> set control argument for optim
dots <- list(...)
+
+ dotsI <- .filterEargs(dots)
+ if(is.null(dotsI$useApply)) dotsI$useApply <- FALSE
+
if(is.null(dots$method)) dots$method <- "L-BFGS-B"
if(!is.null(dots$control)){
@@ -171,7 +176,7 @@
-getInfGamma(L2deriv = L2deriv, risk = risk0,
neighbor = neighbor, biastype = biastype,
Distr = Distr, stand = A0, cent = z0, clip = b1,
- power = 2)+radius(neighbor)^2*b1^2
+ power = 2,...)+radius(neighbor)^2*b1^2
}
b0 <- optimize(funint.opt, interval=c(1e-8,1e8))$minimum
@@ -223,7 +228,7 @@
getInfGamma(L2deriv = L2deriv, risk = riskA,
neighbor = neighbor, biastype = biastype,
Distr = Distr, stand = A0, cent = z0, clip = b0,
- power = 2)/2 -
+ power = 2,...)/2 -
# ~ - E[|Y_A|_Q^2 (1-w_b(|Y_A|_Q))^2]/2
sum(diag(std0%*%A0%*%t(trafo)) ))
## ~tr_Q AD'
@@ -231,32 +236,30 @@
## in case TotalVarNeighborhood additional correction term:
if(is(neighbor,"TotalVarNeighborhood"))
val <- (val -a0^2/2 -
- E(Distr, fun = function(x){ ## ~ - E Y_-^2/2
+ do.call(E, c(list(Distr, fun = function(x){ ## ~ - E Y_-^2/2
L2 <- evalRandVar(L2deriv, as.matrix(x)) [,,1]- z0
Y <- A0 %*% L2
return(Y^2*(Y<0))
- }, useApply = FALSE)/2)
+ }),dotsI))/2)
}else if(is(riskA,"asMSE")){
- val <- (E(object = Distr, fun = function(x){
+ val <- (do.call(E, c(list(object = Distr, fun = function(x){
X <- evalRandVar(L2deriv, as.matrix(x))[,,1] - z0
Y <- A0 %*% X
nY <- norm(risk0)(Y)
return(nY^2*weight(w0)(X))
- }, # E|Y|^2 w
- useApply=FALSE) /2 -
+ }),dotsI))/2 - # E|Y|^2 w
sum(diag(std0%*%A0%*%t(trafo)) ))
## ~tr_Q AD'
## in case TotalVarNeighborhood additional correction term:
if(is(neighbor,"TotalVarNeighborhood"))
val <- (val -a0^2/2 -
- E(Distr, fun = function(x){
+ do.call(E,c(list(Distr, fun = function(x){
X <- evalRandVar(L2deriv, as.matrix(x))[,,1] - z0
Y <- A0 %*% X
return(Y^2*(Y<0))
- },
- useApply=FALSE)/2)
+ }),dotsI))/2)
}
## if this is the current optimum
Modified: branches/robast-1.2/pkg/ROptEst/R/getInfRad.R
===================================================================
--- branches/robast-1.2/pkg/ROptEst/R/getInfRad.R 2018-08-12 08:15:22 UTC (rev 1130)
+++ branches/robast-1.2/pkg/ROptEst/R/getInfRad.R 2018-08-12 08:54:12 UTC (rev 1131)
@@ -30,10 +30,10 @@
risk = "asMSE",
neighbor = "UncondNeighborhood"),
function(clip, L2deriv, risk, neighbor, biastype,
- Distr, stand, cent, trafo){
+ Distr, stand, cent, trafo, ...){
gamm <- getInfGamma(L2deriv = L2deriv, risk = risk, neighbor = neighbor,
biastype = biastype, Distr = Distr, stand = stand,
- cent = cent, clip = clip)
+ cent = cent, clip = clip, ...)
return((-gamm/clip)^.5)
})
Modified: branches/robast-1.2/pkg/ROptEst/R/getInfRobIC_asAnscombe.R
===================================================================
--- branches/robast-1.2/pkg/ROptEst/R/getInfRobIC_asAnscombe.R 2018-08-12 08:15:22 UTC (rev 1130)
+++ branches/robast-1.2/pkg/ROptEst/R/getInfRobIC_asAnscombe.R 2018-08-12 08:54:12 UTC (rev 1131)
@@ -109,6 +109,9 @@
OptOrIter = "iterate", maxiter, tol, warn,
verbose = NULL, checkBounds = TRUE, ...){
+ dotsI <- .filterEargs(list(...))
+ if(is.null(dotsI$useApply)) dotsI$useApply <- FALSE
+
if(missing(verbose)|| is.null(verbose))
verbose <- getRobAStBaseOption("all.verbose")
@@ -156,7 +159,7 @@
trafo = trafo, maxiter = maxiter,
tol = tol,
warn = FALSE, Finfo = Finfo,
- QuadForm = std, verbose = verbose)
+ QuadForm = std, verbose = verbose,...)
if(is.null(lower)||(lower< lowBerg$b))
{lower <- lowBerg$b
@@ -212,12 +215,12 @@
chkbd <- if(it.erg<25) FALSE else checkBounds
verbL <- if(it.erg<25) FALSE else verbose
- erg <<- getInfRobIC(L2deriv, risk.b, neighbor,
+ erg <<- do.call(getInfRobIC, c(list(L2deriv, risk.b, neighbor,
Distr, DistrSymm, L2derivSymm,
L2derivDistrSymm, Finfo, trafo, onesetLM = onesetLM,
z.start, A.start, upper = upper, lower = lower,
OptOrIter = OptOrIter, maxiter = maxi, tol = toli , warn = warn,
- verbose = verbL, checkBounds = chkbd, ...)
+ verbose = verbL, checkBounds = chkbd), dotsI))
trV <- erg$risk$trAsCov$value
if(verbose) cat("Outer iteration:", it.erg," b_0=", round(b0,3),
" eff=", round(trV.ML/trV,3), "\n")
Modified: branches/robast-1.2/pkg/ROptEst/R/getInfRobIC_asBias.R
===================================================================
--- branches/robast-1.2/pkg/ROptEst/R/getInfRobIC_asBias.R 2018-08-12 08:15:22 UTC (rev 1130)
+++ branches/robast-1.2/pkg/ROptEst/R/getInfRobIC_asBias.R 2018-08-12 08:54:12 UTC (rev 1131)
@@ -32,6 +32,9 @@
A.start, Finfo, trafo, maxiter, tol, warn,
verbose = NULL, ...){
+ dots <- list(...)
+ dotsnames <- names(dots)
+
if(missing(verbose)|| is.null(verbose))
verbose <- getRobAStBaseOption("all.verbose")
@@ -52,14 +55,20 @@
))
if(warn) cat(warntxt)
neighbor at radius <- 15
- res <- getInfRobIC(L2deriv = L2deriv,
- risk = asMSE(normtype = normtype),
- neighbor = neighbor, Distr = Distr,
- DistrSymm = DistrSymm, L2derivSymm = L2derivSymm,
- L2derivDistrSymm = L2derivDistrSymm, Finfo = Finfo,
- trafo = trafo, onesetLM = FALSE, z.start = z.start,
- A.start = A.start, upper = 1e4, maxiter = maxiter,
+ getInfRobICargList <- list(L2deriv = L2deriv,
+ risk = asMSE(normtype = normtype),
+ neighbor = neighbor, Distr = Distr,
+ DistrSymm = DistrSymm, L2derivSymm = L2derivSymm,
+ L2derivDistrSymm = L2derivDistrSymm, Finfo = Finfo,
+ trafo = trafo, onesetLM = FALSE, z.start = z.start,
+ A.start = A.start, upper = 1e4, maxiter = maxiter,
tol = tol, warn = warn, verbose = verbose)
+
+ for(item in dotsnames)
+ if(!item %in% names(getInfRobICargList))
+ getInfRobICargList[[item]] <- dots[[item]]
+
+ res <- do.call(getInfRobIC,getInfRobICargList)
A.max <- max(abs(res$A))
res$A <- res$A/A.max
@@ -189,7 +198,7 @@
biastype = "BiasType"),
function(L2deriv, neighbor, biastype, normtype, Distr,
z.start, A.start, z.comp, A.comp, Finfo, trafo, maxiter, tol,
- verbose = NULL){
+ verbose = NULL, ...){
if(missing(verbose)|| is.null(verbose))
verbose <- getRobAStBaseOption("all.verbose")
@@ -197,7 +206,7 @@
eerg <- .LowerCaseMultivariate(L2deriv = L2deriv, neighbor = neighbor,
biastype = biastype, normtype = normtype, Distr = Distr,
Finfo = Finfo, trafo, z.start, A.start = A.start, z.comp = z.comp,
- A.comp = DA.comp, maxiter = maxiter, tol = tol, verbose = verbose)
+ A.comp = DA.comp, maxiter = maxiter, tol = tol, verbose = verbose, ...)
erg <- eerg$erg
b <- 1/erg$value
@@ -229,7 +238,7 @@
Cov <- getInfV(L2deriv = L2deriv, neighbor = neighbor,
biastype = biastype, Distr = Distr,
[TRUNCATED]
To get the complete diff run:
svnlook diff /svnroot/robast -r 1131
More information about the Robast-commits
mailing list