[Robast-commits] r630 - in branches/robast-0.9/pkg: RobAStRDA/R RobAStRDA/man RobExtremes/R RobExtremes/inst/AddMaterial/interpolation RobExtremes/man

noreply at r-forge.r-project.org noreply at r-forge.r-project.org
Wed Mar 13 14:33:37 CET 2013


Author: ruckdeschel
Date: 2013-03-13 14:33:37 +0100 (Wed, 13 Mar 2013)
New Revision: 630

Modified:
   branches/robast-0.9/pkg/RobAStRDA/R/interpolAux.R
   branches/robast-0.9/pkg/RobAStRDA/man/internal-interpolate.Rd
   branches/robast-0.9/pkg/RobExtremes/R/SnQn.R
   branches/robast-0.9/pkg/RobExtremes/R/getStartIC.R
   branches/robast-0.9/pkg/RobExtremes/R/internal-getpsi.R
   branches/robast-0.9/pkg/RobExtremes/inst/AddMaterial/interpolation/SnTest.Rdata
   branches/robast-0.9/pkg/RobExtremes/inst/AddMaterial/interpolation/Snplot.pdf
   branches/robast-0.9/pkg/RobExtremes/inst/AddMaterial/interpolation/checkSn.R
   branches/robast-0.9/pkg/RobExtremes/man/internal-interpolate.Rd
Log:
+ RobExtremes: revised Sn diagnostics
+ RobAStRDA:  interpolators (.generateInterpolators, .computeInterpolators) gain possibility to restrict domain of extrapolation by argument extrapol (if beyond, the interpolating functions return NA);
+ RobExtremes: Sn and getStartIC learn to deal with return NAs from interpolation (to this end a new function   .is.na.Psi)

Modified: branches/robast-0.9/pkg/RobAStRDA/R/interpolAux.R
===================================================================
--- branches/robast-0.9/pkg/RobAStRDA/R/interpolAux.R	2013-03-13 09:52:00 UTC (rev 629)
+++ branches/robast-0.9/pkg/RobAStRDA/R/interpolAux.R	2013-03-13 13:33:37 UTC (rev 630)
@@ -35,12 +35,15 @@
 # .generateInterpolators generates the interpolators to a given grid
 #     and returns a list of the given grid and the function list
 ############################################################################
-.generateInterpolators <- function(Grid, approxOrspline = "spline"){
+.generateInterpolators <- function(Grid, approxOrspline = "spline",
+                                   extrapol = c(NA,NA)){
   thGrid <- Grid[,1]
   LMGrid <- Grid[,-1,drop=FALSE]
   fctL <- vector("list",ncol(LMGrid))
   xm <- thGrid[1]
   xM <- (rev(thGrid))[1]
+  xm0 <- if(is.na(extrapol[1])) -Inf else extrapol[1]
+  xM0 <- if(is.na(extrapol[2]))  Inf else extrapol[2]
   for(i in 1:ncol(LMGrid)){
        LMG <- LMGrid[,i]
        fct <- if(approxOrspline=="spline")
@@ -52,8 +55,10 @@
        fctX <- function(x){
             y0 <- fct(x)
             y1 <- y0
-            y1[x<xm] <- ym+dym*(x[x<xm]-xm)
-            y1[x>xM] <- yM+dyM*(x[x>xM]-xM)
+            y1[x < xm & x >= xm0] <- ym+dym*(x[x < xm & x >= xm0]-xm)
+            y1[x > xM & x <= xM0] <- yM+dyM*(x[x > xM & x <= xM0]-xM)
+            y1[x < xm0] <- NA
+            y1[x > xM0] <- NA
             if(any(is.na(y0)))
                warning(paste("There have been xi-values out of range ",
                              "of the interpolation grid.", sep = ""))
@@ -191,7 +196,8 @@
                                    includeGrids = NULL, includeNams = NULL,
                                    excludeGrids = NULL, excludeNams = NULL,
                                    withPrint = TRUE, withSmoothFct = FALSE,
-                                   approxOrspline = "spline"){
+                                   approxOrspline = "spline",
+                                   extrapol = c(NA,NA)){
 
   wprint <- function(...){ if (withPrint) print(...)}
 

Modified: branches/robast-0.9/pkg/RobAStRDA/man/internal-interpolate.Rd
===================================================================
--- branches/robast-0.9/pkg/RobAStRDA/man/internal-interpolate.Rd	2013-03-13 09:52:00 UTC (rev 629)
+++ branches/robast-0.9/pkg/RobAStRDA/man/internal-interpolate.Rd	2013-03-13 13:33:37 UTC (rev 630)
@@ -27,7 +27,7 @@
 
 .readGridFromCSV(fromFileCSV)
 
-.generateInterpolators(Grid, approxOrspline = "spline")
+.generateInterpolators(Grid, approxOrspline = "spline", extrapol = c(NA,NA))
 
 .saveGridToRda(fromFileCSV, toFileRDA = "sysdata.rda", withMerge = FALSE,
                withPrint = TRUE, withSmooth = TRUE, df = NULL)
@@ -37,7 +37,8 @@
 .computeInterpolators(sysdataFiles, toFileRDA = "sysdata.rda",
       includeGrids = NULL, includeNams = NULL,
       excludeGrids = NULL, excludeNams = NULL,
-      withPrint = TRUE, withSmoothFct = FALSE, approxOrspline = "spline")
+      withPrint = TRUE, withSmoothFct = FALSE,
+      approxOrspline = "spline", extrapol = c(NA,NA))
 
 .mergeF(file,envir, includeGrids = NULL, includeNams = NULL,
         excludeGrids = NULL, excludeNams = NULL)
@@ -65,6 +66,10 @@
   \item{approxOrspline}{character; if \code{approxOrspline=="spline"} (default),
     \code{\link{splinefun}} is used for generating the interpolators, otherwise
     we use \code{\link{approxfun}}. }
+  \item{extrapol}{numeric of length 2; lower and upper bound, upto which
+       extrapolation is done; beyond, the interpolator returns \code{NA};
+       if one (or both) entries of \code{extrapol} are \code{NA}, we extrapolate
+       beyond limit. }
   \item{toFileRDA}{character; the \file{.rda}-file to which the interpolators
         are saved. }
   \item{withMerge}{logical of length 1: in case a respective grid already

Modified: branches/robast-0.9/pkg/RobExtremes/R/SnQn.R
===================================================================
--- branches/robast-0.9/pkg/RobExtremes/R/SnQn.R	2013-03-13 09:52:00 UTC (rev 629)
+++ branches/robast-0.9/pkg/RobExtremes/R/SnQn.R	2013-03-13 13:33:37 UTC (rev 630)
@@ -111,6 +111,7 @@
        if(!nam %in% names(sng)) return(Sn(as(x,"AbscontDistribution")))
        snf <- sng[[nam]][[.versionSuff("fun")]]
        ret <- snf(shape(x))
+       if(is.na(ret)) return(Sn(as(x,"AbscontDistribution")))
     }else ret <- scale(x)*Sn(x=x/scale(x))
     return(ret)
 }

Modified: branches/robast-0.9/pkg/RobExtremes/R/getStartIC.R
===================================================================
--- branches/robast-0.9/pkg/RobExtremes/R/getStartIC.R	2013-03-13 09:52:00 UTC (rev 629)
+++ branches/robast-0.9/pkg/RobExtremes/R/getStartIC.R	2013-03-13 13:33:37 UTC (rev 630)
@@ -2,10 +2,16 @@
            function(model, risk, ...){
 
     mc <- match.call(expand.dots=TRUE)
+    mc$risk <- if(type(risk)==".MBRE") asMSE() else asBias()
+    mc$neighbor <- ContNeighborhood(radius=0.5)
 
     gridn <- type(risk)
     nam <- gsub(" ","",name(model))
     param1 <- param(model)
+
+    scshnm <- scaleshapename(model)
+    shnam <- scshnm["shape"]
+
     nsng <- character(0)
     sng <- try(getFromNamespace(gridn, ns = "RobAStRDA"), silent=TRUE)
     if(!is(sng,"try-error")) nsng <- names(sng)
@@ -14,15 +20,20 @@
           fctN <- .versionSuff("fun")
           interpolfct <- sng[[nam]][[fctN]]
           .modifyIC <- function(L2Fam, IC){
-                   para <- param(L2Fam)
-                   .getPsi(para, interpolfct, L2Fam, type(risk))}
-          IC0 <- .getPsi(param1, interpolfct, model, type(risk))
-          IC0 at modifyIC <- .modifyIC
-          return(IC0)
+                    para <- param(L2Fam)
+                    if(!.is.na.Psi(para, interpolfct, shnam))
+                       return(.getPsi(para, interpolfct, L2Fam, type(risk)))
+                    else
+                       return(do.call(getStartIC, as.list(mc[-1]),
+                              envir=parent.frame(2)))
+          }
+          if(!.is.na.Psi(param1, interpolfct, shnam)){
+             IC0 <- .getPsi(param1, interpolfct, model, type(risk))
+             IC0 at modifyIC <- .modifyIC
+             return(IC0)
+          }
        }
     }
-    mc$risk <- if(type(risk)==".MBRE") asMSE() else asBias()
-    mc$neighbor <- ContNeighborhood(radius=0.5)
     return(do.call(getStartIC, as.list(mc[-1]), envir=parent.frame(2)))
     })
 

Modified: branches/robast-0.9/pkg/RobExtremes/R/internal-getpsi.R
===================================================================
--- branches/robast-0.9/pkg/RobExtremes/R/internal-getpsi.R	2013-03-13 09:52:00 UTC (rev 629)
+++ branches/robast-0.9/pkg/RobExtremes/R/internal-getpsi.R	2013-03-13 13:33:37 UTC (rev 630)
@@ -1,3 +1,7 @@
+.is.na.Psi <- function(param, fct, nam){
+   xi <- main(param)[nam]
+   return(is.na(fct[[1]](xi)))
+}
 .getPsi <- function(param, fct, L2Fam , type){
 
    scshnm <- scaleshapename(L2Fam)

Modified: branches/robast-0.9/pkg/RobExtremes/inst/AddMaterial/interpolation/SnTest.Rdata
===================================================================
(Binary files differ)

Modified: branches/robast-0.9/pkg/RobExtremes/inst/AddMaterial/interpolation/Snplot.pdf
===================================================================
(Binary files differ)

Modified: branches/robast-0.9/pkg/RobExtremes/inst/AddMaterial/interpolation/checkSn.R
===================================================================
--- branches/robast-0.9/pkg/RobExtremes/inst/AddMaterial/interpolation/checkSn.R	2013-03-13 09:52:00 UTC (rev 629)
+++ branches/robast-0.9/pkg/RobExtremes/inst/AddMaterial/interpolation/checkSn.R	2013-03-13 13:33:37 UTC (rev 630)
@@ -29,10 +29,24 @@
   S3g <- sapply(xig, gSn, Gammad)
   S4g <- sapply(xig, gSn, Weibull)
 })
+##
+#   user  system elapsed
+#   2.31    0.00    2.32
+#
 system.time({S1ga <- sapply(xig, gSna, GPareto)
 S2ga <- sapply(xig, gSna, GEV)
 S3ga <- sapply(xig, gSna, Gammad)
 S4ga <- sapply(xig, gSna, Weibull)})
+##
+#    user   system  elapsed
+#  966.03     0.83   979.77
+#
+setwd("C:/rtest/RobASt/branches/robast-0.9/pkg/RobExtremes/inst/AddMaterial/interpolation")
+#
+save(S1g, S1ga, S2g, S2ga, S3g, S3ga, S4g, S4ga, file="SnTest.Rdata")
+#
+#
+pdf("Snplot.pdf")
 par(mfrow=c(2,2))
 plot(xig, S1g, type="l")
 lines(xig, S1ga, col="red")
@@ -43,3 +57,4 @@
 plot(xig, S4g, type="l")
 lines(xig, S4ga, col="red")
 par(mfrow=c(1,1))
+dev.off()

Modified: branches/robast-0.9/pkg/RobExtremes/man/internal-interpolate.Rd
===================================================================
--- branches/robast-0.9/pkg/RobExtremes/man/internal-interpolate.Rd	2013-03-13 09:52:00 UTC (rev 629)
+++ branches/robast-0.9/pkg/RobExtremes/man/internal-interpolate.Rd	2013-03-13 13:33:37 UTC (rev 630)
@@ -6,6 +6,7 @@
 \alias{.MBRE.xi}
 \alias{.getLMGrid}
 \alias{.getPsi}
+\alias{.is.na.Psi}
 \alias{.svInt}
 \alias{.generateInterpGridSn}
 
@@ -20,6 +21,8 @@
 \usage{
 .getPsi(param, fct, L2Fam , type)
 
+.is.na.Psi(param, fct, nam = "shape")
+
 .modify.xi.PFam.call(xi, PFam)
 
 .RMXE.xi(xi, PFam)
@@ -43,6 +46,7 @@
               at which to evaluate the Lagrange multipliers or LDEstimators;
               in our use case, it is a shape-scale model, hence the respective
               (main) parameter must contain \code{"scale"} and \code{"shape"}. }
+  \item{nam}{character; name of the shape parameter. }
   \item{type}{type of the optimality: one of ".OMSE" for maxMSE,
        ".RMXE" for rmx, and  ".MBRE" for MBRE. }
   \item{xi}{numeric of length 1; shape value. }
@@ -62,6 +66,11 @@
     from an object from \file{sysdata.rda} and generates a respective
    \code{HampelIC} object by a call to  \code{generateIC}.
 
+   \code{.is.na.Psi} checks whether the shape parameter already lies
+   beyond the range for which inter-/extrapolation is admitted
+   (and, correspondingly, returns \code{TRUE} if one has to compute the
+   IC completely anew.).
+
   \code{.MBRE.xi} computes the Lagrange multipliers for the MBRE estimator,
   \code{.OMSE.xi} for the OMSE estimator at radius \code{r=0.5},
   and \code{.RMXE.xi} the RMXE estimator.
@@ -83,7 +92,8 @@
       version.
       }
 \value{
-  \item{.getpsi}{an IC.}
+  \item{.getpsi}{an IC. }
+  \item{.is.na.Psi}{logical of length 1. }
   \item{.modify.xi.PFam.call}{A call to evaluate the parametric
             family at the new parameter value. }
   \item{.MBRE.xi}{A list with items \code{b} (a number; clipping height),
@@ -101,12 +111,5 @@
   \item{.svInt}{  \code{invisible(NULL)}}
 }
 \seealso{\code{\link{interpolateSn}}}
-\examples{
-\dontrun{
-### code to produce grid for GPareto:
-  .saveInterpGrid(sysRdaFolder = "C:/rtest/RobASt/branches/RobASt-0.9/pkg/RobExtremes/R",
-              accuracy = 800)
-}
-}
 \keyword{internal}
 \concept{utilities}



More information about the Robast-commits mailing list