[Distr-commits] r1276 - in branches/distr-2.8/pkg/distrEx: . R inst man

noreply at r-forge.r-project.org noreply at r-forge.r-project.org
Wed Aug 15 15:49:37 CEST 2018


Author: ruckdeschel
Date: 2018-08-15 15:49:37 +0200 (Wed, 15 Aug 2018)
New Revision: 1276

Modified:
   branches/distr-2.8/pkg/distrEx/NAMESPACE
   branches/distr-2.8/pkg/distrEx/R/AsymTotalVarDist.R
   branches/distr-2.8/pkg/distrEx/R/CvMDist.R
   branches/distr-2.8/pkg/distrEx/R/Expectation.R
   branches/distr-2.8/pkg/distrEx/R/GammaWeibullExpectation.R
   branches/distr-2.8/pkg/distrEx/R/HellingerDist.R
   branches/distr-2.8/pkg/distrEx/R/Internalfunctions.R
   branches/distr-2.8/pkg/distrEx/R/OAsymTotalVarDist.R
   branches/distr-2.8/pkg/distrEx/R/TotalVarDist.R
   branches/distr-2.8/pkg/distrEx/R/distrExIntegrate.R
   branches/distr-2.8/pkg/distrEx/inst/NEWS
   branches/distr-2.8/pkg/distrEx/man/0distrEx-package.Rd
   branches/distr-2.8/pkg/distrEx/man/AsymTotalVarDist.Rd
   branches/distr-2.8/pkg/distrEx/man/CvMDist.Rd
   branches/distr-2.8/pkg/distrEx/man/E.Rd
   branches/distr-2.8/pkg/distrEx/man/HellingerDist.Rd
   branches/distr-2.8/pkg/distrEx/man/OAsymTotalVarDist.Rd
   branches/distr-2.8/pkg/distrEx/man/TotalVarDist.Rd
   branches/distr-2.8/pkg/distrEx/man/distrExIntegrate.Rd
   branches/distr-2.8/pkg/distrEx/man/internals.Rd
Log:
[distrEx] branch 2.8: 
+ introduced filter functions .filterEargs and .filterFunargs to warrant some safety that only those args from   "..." become arguments 
  of the integrand, of distrExIntegrate, of E(), of integrate, of GLIntegrate, of quantiles and IQR, which are within the formals of 
  the respective function...
+ .filterEargs in fact is (almost) .filterEargs from RobAStBase version 1.2, but:
    (a) it takes the formals automatically and not by an explicit string
    (b) .filterEargs in RobAStBase also filters in the items of "E.argList" 
+ the return values of distrExIntegrate and all E()-methods gain an 
  optional attribute "diagnostic" which is filled if argument
  diagnostic is TRUE (the E()-methods whenever they use distrExIntegrate
  in (parts of) their computation.



Modified: branches/distr-2.8/pkg/distrEx/NAMESPACE
===================================================================
--- branches/distr-2.8/pkg/distrEx/NAMESPACE	2018-08-14 11:40:31 UTC (rev 1275)
+++ branches/distr-2.8/pkg/distrEx/NAMESPACE	2018-08-15 13:49:37 UTC (rev 1276)
@@ -53,4 +53,4 @@
        "distrExMASK", "distrExoptions", "distrExMOVED")
 export("make01","PrognCondDistribution",
        "PrognCondition")
-export(".getIntbounds", ".qtlIntegrate")
+export(".getIntbounds", ".qtlIntegrate", ".filterEargs", ".filterFunargs")

Modified: branches/distr-2.8/pkg/distrEx/R/AsymTotalVarDist.R
===================================================================
--- branches/distr-2.8/pkg/distrEx/R/AsymTotalVarDist.R	2018-08-14 11:40:31 UTC (rev 1275)
+++ branches/distr-2.8/pkg/distrEx/R/AsymTotalVarDist.R	2018-08-15 13:49:37 UTC (rev 1276)
@@ -11,11 +11,14 @@
     function(e1, e2, rho = 1,  
              rel.tol = .Machine$double.eps^0.3, maxiter=1000, Ngrid = 10000,
              TruncQuantile = getdistrOption("TruncQuantile"),
-             IQR.fac = 15){ 
+             IQR.fac = 15, ..., diagnostic = FALSE){
 
         ## if we have to recall this method with a smaller TruncQuantile arg:
         mc <-  as.list(match.call(call = sys.call(sys.parent(1)))[-1])
         mc$TruncQuantile <- TruncQuantile * 1.8
+        dots <- list(...)
+        dotsn <- names(dots)
+        dotsI <- dots[dotsn %in% c("order","subdivisions", "stop.on.error")]
 
         #block warnings:
         o.warn <- getOption("warn"); options(warn = -1)
@@ -36,9 +39,11 @@
         d1 <- d(e1); d2 <- d(e2)
         #
         ### integration as a function of c:
-        Eip <- function(f, c00)
-                distrExIntegrate(f, lower = low, upper = up, 
-                                    rel.tol = rel.tol, c00 = c00)
+        Eip <- function(f, c00, diagnostic0 = FALSE)
+                do.call(distrExIntegrate,c(list(f, lower = low, upper = up,
+                                    rel.tol = rel.tol, c00 = c00,
+                                    diagnostic = diagnostic0),dotsI))
+
        # positive part
        integ.p <- function(x,c00)  pmax(d2(x)-c00*d1(x),0)
        # negative part
@@ -46,8 +51,8 @@
 
        ## function zero c(rho) of which is to be found 
        fct <- function(c0, rho = rho){
-           e.p <- Eip(f=integ.p, c00=c0)
-           e.m <- Eip(f=integ.m, c00=c0)
+           e.p <- Eip(f=integ.p, c00=c0, diagnostic0 = FALSE)
+           e.m <- Eip(f=integ.m, c00=c0, diagnostic0 = FALSE)
            e.p*rho - e.m
            } 
        
@@ -73,8 +78,13 @@
        tef <- fct(1, rho = rho)
        ### if c=1 is already a zero:
        if(tef == 0){
-          res <- Eip(f=integ.p,c00=1)
+          res <- Eip(f=integ.p,c00=1, diagnostic0 = diagnostic)
           names(res) <- "asym. total variation distance"
+          if(diagnostic){
+             diagn <- attr(res,"diagnostic")
+             diagn[["call"]] <- match.call()
+             attr(res,"diagnostic") <- diagn
+          }
           return(res)
        }   
        # else: only have to search in c in [low1;1] resp [1;up1]
@@ -82,7 +92,7 @@
 
        c.rho <- try(uniroot(fct, lower = low1, upper = up1, 
                             rho = rho, tol = rel.tol,
-                            maxiter = maxiter)$root, 
+                            maxiter = maxiter)$root,
                     silent = TRUE)
 
        ## if does not give reasonable solution recall function with 
@@ -93,7 +103,12 @@
                                         e2 = "AbscontDistribution")), 
                               args = mc))             
        ## else:
-       res <- Eip(f=integ.p, c00=c.rho)
+       res <- Eip(f=integ.p, c00=c.rho, diagnostic0=diagnostic)
+       if(diagnostic){
+          diagn <- attr(res,"diagnostic")
+          diagn[["call"]] <- match.call()
+          attr(res,"diagnostic") <- diagn
+       }
        names(res) <- "asym. total variation distance"
        return(res)
     })
@@ -186,11 +201,19 @@
              up.discr = getUp(e2), h.smooth = getdistrExOption("hSmooth"),
              rel.tol = .Machine$double.eps^0.3, maxiter=1000, Ngrid = 10000,
              TruncQuantile = getdistrOption("TruncQuantile"),
-             IQR.fac = 15){
-        .asis.smooth.discretize.distance(e1, e2, asis.smooth.discretize, n.discr,
-                 low.discr, up.discr, h.smooth, AsymTotalVarDist, rho = rho,
+             IQR.fac = 15, ..., diagnostic = FALSE){
+        res <- .asis.smooth.discretize.distance(x = e1, Distribution = e2,
+                asis.smooth.discretize=asis.smooth.discretize, n.discr = n.discr,
+                 low.discr= low.discr, up.discr = up.discr, h.smooth = h.smooth,
+                 distance = AsymTotalVarDist, rho = rho,
                  rel.tol = rel.tol, maxiter = maxiter, Ngrid = Ngrid,
-                 TruncQuantile = TruncQuantile, IQR.fac = IQR.fac)
+                 TruncQuantile = TruncQuantile, IQR.fac = IQR.fac, ..., diagnostic = diagnostic)
+       if(diagnostic){
+          diagn <- attr(res,"diagnostic")
+          diagn[["call"]] <- match.call()
+          attr(res,"diagnostic") <- diagn
+       }
+        return(res)
      })
 setMethod("AsymTotalVarDist", signature(e1 = "AbscontDistribution",
                                      e2 = "numeric"),
@@ -199,12 +222,18 @@
              up.discr = getUp(e1), h.smooth = getdistrExOption("hSmooth"),
              rel.tol = .Machine$double.eps^0.3, maxiter=1000, Ngrid = 10000,
              TruncQuantile = getdistrOption("TruncQuantile"),
-             IQR.fac = 15){
-        return(AsymTotalVarDist(e2, e1, rho= 1/rho, 
+             IQR.fac = 15, ..., diagnostic = FALSE){
+        res <- AsymTotalVarDist(e1=e2, e2=e1, rho= 1/rho,
                   asis.smooth.discretize = asis.smooth.discretize, 
                   low.discr = low.discr, up.discr = up.discr, h.smooth = h.smooth,
                   rel.tol = rel.tol, maxiter = maxiter, Ngrid = Ngrid,
-                  TruncQuantile = TruncQuantile, IQR.fac = IQR.fac))
+                  TruncQuantile = TruncQuantile, IQR.fac = IQR.fac, ..., diagnostic = diagnostic)
+       if(diagnostic){
+          diagn <- attr(res,"diagnostic")
+          diagn[["call"]] <- match.call()
+          attr(res,"diagnostic") <- diagn
+       }
+        return(res)
     })
 
 setMethod("AsymTotalVarDist",  signature(e1 = "AcDcLcDistribution",
@@ -212,12 +241,15 @@
            function(e1, e2, rho = 1,  
              rel.tol = .Machine$double.eps^0.3, maxiter=1000, Ngrid = 10000,
              TruncQuantile = getdistrOption("TruncQuantile"),
-             IQR.fac = 15){
+             IQR.fac = 15, ..., diagnostic = FALSE){
         ## if we have to recall this method with a smaller TruncQuantile arg:
         mc <-  as.list(match.call(call = sys.call(sys.parent(1)))[-1])
         mc$TruncQuantile <- TruncQuantile * 1.8
+        dots <- list(...)
+        dotsn <- names(dots)
+        dotsI <- dots[dotsn %in% c("order","subdivisions", "stop.on.error")]
+        #block warnings:
 
-        #block warnings:
         o.warn <- getOption("warn"); options(warn = -1)
         on.exit(options(warn=o.warn))
 
@@ -254,9 +286,10 @@
         low <- max(low,low0); up <- min(up,up0)
         #
         ### integration as a function of c:
-        Eip <- function(f, c00)
-                distrExIntegrate(f, lower = low, upper = up, 
-                                    rel.tol = rel.tol, c00 = c00)
+        Eip <- function(f, c00, diagnostic0 = FALSE)
+                do.call(distrExIntegrate, c(list(f, lower = low, upper = up,
+                                    rel.tol = rel.tol, c00 = c00,
+                                    diagnostic = diagnostic0),dotsI))
        # positive part
        integ.p.c <- function(x,c00)  pmax(ac2.d(x)-c00*ac1.d(x),0)
        # negative part
@@ -275,8 +308,8 @@
 
        ## function zero c(rho) of which is to be found 
        fct <- function(c0, rho = rho){
-           e.p.c <- Eip(f=integ.p.c, c00=c0)
-           e.m.c <- Eip(f=integ.m.c, c00=c0)
+           e.p.c <- Eip(f=integ.p.c, c00=c0, diagnostic0 = FALSE)
+           e.m.c <- Eip(f=integ.m.c, c00=c0, diagnostic0 = FALSE)
            e.p.d <- sum(integ.p.d(c0))
            e.m.d <- sum(integ.m.d(c0))
            e.p <- e.p.c + e.p.d
@@ -310,11 +343,17 @@
        ## gives range for c:
        low1 <- min(d.range,dx.range); up1 <- max(d.range,dx.range)
        ## in any case compare with c=1
-       tef <- fct(1, rho = rho)
+       tef <- fct(1, rho = rho, ...)
        ### if c=1 is already a zero:
        if(tef == 0){
-          res <- Eip(f=integ.p.c,c00=1)+sum(integ.p.d(1))
+          res <- Eip(f=integ.p.c,c00=1, diagnostic0 = diagnostic)
+          if(diagnostic){
+             diagn <- attr(res, "diagnostic")
+             diagn[["call"]] <-  match.call()
+          }
+          res <- res +sum(integ.p.d(1))
           names(res) <- "asym. total variation distance"
+          if(diagnostic) attr(res, "diagnostic") <- diagn
           return(res)
        }   
        # else: only have to search in c in [low1;1] resp [1;up1]
@@ -322,7 +361,7 @@
 
        c.rho <- try(uniroot(fct, lower = low1, upper = up1, 
                             rho = rho, tol = rel.tol,
-                            maxiter = maxiter)$root, 
+                            maxiter = maxiter)$root,
                     silent = TRUE)
 
        ## if does not give reasonable solution recall function with 
@@ -331,12 +370,22 @@
            return(do.call(getMethod("AsymTotalVarDist", 
                            signature(e1 = "AcDcLcDistribution",
                                      e2 = "AcDcLcDistribution")), 
-                              args = mc))             
+                           args = mc))
        }
-       res <- Eip(f=integ.p.c, c00=c.rho)+sum(integ.p.d(c.rho))
+       res <- Eip(f=integ.p.c, c00=c.rho, diagnostic0 = diagnostic)
+       if(diagnostic){
+          diagn <- attr(res, "diagnostic")
+          diagn[["call"]] <-  match.call()
+       }
+       res <- res +sum(integ.p.d(c.rho))
        names(res) <- "asym. total variation distance"
+       if(diagnostic){
+          diagn <- attr(res,"diagnostic")
+          diagn[["call"]] <- match.call()
+          attr(res,"diagnostic") <- diagn
+       }
        return(res)
-              })
+    })
 
 setMethod("AsymTotalVarDist", signature(e1 = "LatticeDistribution", 
                                          e2 = "LatticeDistribution"),

Modified: branches/distr-2.8/pkg/distrEx/R/CvMDist.R
===================================================================
--- branches/distr-2.8/pkg/distrEx/R/CvMDist.R	2018-08-14 11:40:31 UTC (rev 1275)
+++ branches/distr-2.8/pkg/distrEx/R/CvMDist.R	2018-08-15 13:49:37 UTC (rev 1276)
@@ -4,7 +4,7 @@
 ###############################################################################
 setMethod("CvMDist", signature(e1 = "UnivariateDistribution",
                                     e2 = "UnivariateDistribution"),
-    function(e1, e2, mu = e1, useApply = FALSE, ... ){
+    function(e1, e2, mu = e1, useApply = FALSE, ..., diagnostic = FALSE){
         o.warn <- getOption("warn"); options(warn = -1)
         on.exit(options(warn=o.warn))
         if(is.null(e1 at p)){
@@ -17,7 +17,12 @@
         e2 <- new("UnivariateDistribution", r=e2 at r, 
                    p = e2.erg$pfun, d = e2.erg$dfun, q = e2.erg$qfun,  
                    .withSim = TRUE, .withArith = FALSE)}
-        res <- E(mu, fun = function(t) {(p(e1)(t)-p(e2)(t))^2}, useApply = useApply, ...)^.5
+        res <- E(mu, fun = function(t) {(p(e1)(t)-p(e2)(t))^2}, useApply = useApply, ..., diagnostic = diagnostic)^.5
+        if(diagnostic){
+           diagn <- attr(res,"diagnostic")
+           diagn[["call"]] <- match.call()
+           attr(res,"diagnostic") <- diagn
+        }
         names(res) <- "CvM distance"
         return(res)
     })
@@ -25,15 +30,21 @@
 ## CvM distance
 setMethod("CvMDist", signature(e1 = "numeric",
                                     e2 = "UnivariateDistribution"),
-    function(e1, e2, mu = e1, ...)
+    function(e1, e2, mu = e1, ..., diagnostic = FALSE)
         { o.warn <- getOption("warn"); options(warn = -1)
           on.exit(options(warn=o.warn))
           if(identical(mu,e2)){
              if(is(e2, "AbscontDistribution"))
-             return(.newCvMDist(e1,e2)) }   
+             return(.newCvMDist(e1,e2)) }
           e10 <- DiscreteDistribution(e1)       
           if(identical(mu,e1)) mu <- e10
-          CvMDist(e1 = e10, e2 = e2, mu = mu, ...)
+          res <- CvMDist(e1 = e10, e2 = e2, mu = mu, ..., diagnostic = diagnostic)
+          if(diagnostic){
+             diagn <- attr(res,"diagnostic")
+             diagn[["call"]] <- match.call()
+             attr(res,"diagnostic") <- diagn
+          }
+          return(res)
          }
     )
 

Modified: branches/distr-2.8/pkg/distrEx/R/Expectation.R
===================================================================
--- branches/distr-2.8/pkg/distrEx/R/Expectation.R	2018-08-14 11:40:31 UTC (rev 1275)
+++ branches/distr-2.8/pkg/distrEx/R/Expectation.R	2018-08-15 13:49:37 UTC (rev 1276)
@@ -1,12 +1,41 @@
 ## Helper function:
+.filterEargs <- function(dots, neg=FALSE){
+        if(length(dots)==0) return(NULL)
+        formIntNames <- setdiff(unique(c(names(formals(distrExIntegrate)),
+                            names(formals(integrate)),
+                            names(formals(GLIntegrate)))),c("f","...", "distr"))
 
+        if(neg) return(dots[! names(dots) %in% formIntNames])
+        else return(dots[names(dots) %in% formIntNames])
+}
+.filterFunargs <- function(dots, fun, neg=FALSE){
+        if(length(dots)==0) return(NULL)
+        formFunNames <- names(formals(fun))
+
+        if(neg) return(dots[! names(dots) %in% formFunNames])
+        else return(dots[names(dots) %in% formFunNames])
+}
+
 .getIntbounds <- function(object, low, upp, lowTQ, uppTQ, IQR.fac, ...){
         qx <- q.l(object)
-        low0 <- qx(lowTQ, lower.tail = TRUE, ...) 
-        upp0 <- ifelse( "lower.tail" %in% names(formals(qx)),
-                       qx(uppTQ, lower.tail = FALSE, ...), 
-                       qx(1-uppTQ, ...))        
-        m <- median(object, ...); s <- IQR(object, ...)
+        fqxn <- names(formals(qx))
+        dots <- list(...)
+        dotsqx <- NULL
+        if(length(dots)){
+           dotsqx <- dots[names(dots) %in% fqxn]
+           dotsqx[["lower.tail"]] <- NULL
+           dotsqx[["p"]] <- NULL
+        }
+#        print(c(list(lowTQ, lower.tail = TRUE), dotsqx))
+        low0 <- do.call(qx,c(list(lowTQ, lower.tail = TRUE), dotsqx))
+        upp0 <- if ("lower.tail" %in% fqxn)
+                   do.call(qx,c(list(uppTQ, lower.tail = FALSE), dotsqx)) else
+                   do.call(qx,c(list(1-uppTQ), dotsqx))
+        if("cond" %in% names(dotsqx)){
+           m <- median(object,cond=dots$cond); s <- IQR(object,cond=dots$cond)
+        }else{
+           m <- median(object); s <- IQR(object)
+        }
         low1 <- m - IQR.fac * s 
         upp1 <- m + IQR.fac * s
         low <- max(low0,low1,low) 
@@ -34,8 +63,9 @@
              rel.tol= getdistrExOption("ErelativeTolerance"), 
              lowerTruncQuantile = getdistrExOption("ElowerTruncQuantile"), 
              upperTruncQuantile = getdistrExOption("EupperTruncQuantile"), 
-             IQR.fac = getdistrExOption("IQR.fac"), ...
+             IQR.fac = getdistrExOption("IQR.fac"), ..., diagnostic = FALSE
              ){
+        mc <- match.call()
         if(is(Symmetry(object),"SphericalSymmetry"))
            return(SymmCenter(Symmetry(object)))
         integrand <- function(x, dfun){ x * dfun(x) }
@@ -49,12 +79,17 @@
         upp <- Ib["upp"]
         #print(Ib)
         if(upp<low) return(0)
-
-        return(distrExIntegrate(f = integrand, 
+        res <- distrExIntegrate(f = integrand,
                     lower = low,
                     upper = upp, 
                     rel.tol = rel.tol, 
-                    distr = object, dfun = d(object)))
+                    distr = object, dfun = d(object), diagnostic = diagnostic)
+        if(diagnostic){
+           diagn <- attr(res, "diagnostic")
+           diagn[["call"]] <- mc
+           attr(res, "diagnostic") <- diagn
+        }
+        return(res)
     })
 setMethod("E", signature(object = "DiscreteDistribution", 
                          fun = "missing", 
@@ -81,21 +116,29 @@
 setMethod("E", signature(object = "AffLinDistribution", 
                          fun = "missing", 
                          cond = "missing"),
-    function(object, low = NULL, upp = NULL, ...){
+    function(object, low = NULL, upp = NULL, ..., diagnostic = FALSE){
+             mc <- match.call()
              if(is(Symmetry(object),"SphericalSymmetry"))
                 return(SymmCenter(Symmetry(object)))
              if(is.null(low)) low <- -Inf
              if(is.null(upp)) upp <- Inf
              if(upp<low) return(0)
              if(object at a >= 0){
-                return(object at a * E(object at X0, low = (low-object at b)/object at a,
-                             upp = (upp-object at b)/object at a, ...) +
-                object at b * (p(object)(upp)-p.l(object)(low)))
+                res0 <- E(object at X0, low = (low-object at b)/object at a,
+                             upp = (upp-object at b)/object at a, ...,
+                             diagnostic = diagnostic)
              }else{
-                return(object at a * E(object at X0, low = (upp-object at b)/object at a,
-                             upp = (low-object at b)/object at a, ...) +
-                object at b *    (p(object)(upp)-p.l(object)(low)))
+                res0 <- E(object at X0, low = (upp-object at b)/object at a,
+                             upp = (low-object at b)/object at a, ...,
+                             diagnostic = diagnostic)
              }
+             res1 <- object at a * res0 + object at b * (p(object)(upp)-p.l(object)(low))
+             if(diagnostic){
+                diagn <- attr(res0, "diagnostic")
+                diagn[["call"]] <- mc
+                attr(res1, "diagnostic") <- diagn
+             }
+             return(res1)
     })
 
 setMethod("E", signature(object = "AffLinAbscontDistribution", 
@@ -107,19 +150,23 @@
 setMethod("E", signature(object = "AffLinDiscreteDistribution", 
                          fun = "missing", 
                          cond = "missing"),
-           getMethod("E", signature(object = "AffLinDistribution", 
+    function(object, low = NULL, upp = NULL, ...){
+         getMethod("E", signature(object = "AffLinDistribution",
+                         fun = "missing",
+                         cond = "missing"))(object = object, low = low, upp = upp,
+                                              ..., diagnostic = FALSE)
+    })
+setMethod("E", signature(object = "AffLinLatticeDistribution",
                          fun = "missing", 
-                         cond = "missing")))    
-setMethod("E", signature(object = "AffLinLatticeDistribution", 
-                         fun = "missing", 
                          cond = "missing"),
-           getMethod("E", signature(object = "AffLinDistribution", 
-                         fun = "missing", 
-                         cond = "missing")))    
+    function(object, low = NULL, upp = NULL, ...){
+         getMethod("E", signature(object = "AffLinDistribution",
+                         fun = "missing",
+                         cond = "missing"))(object = object, low = low, upp = upp,
+                                              ..., diagnostic = FALSE)
+    })
 
 
-
-
 setMethod("E", signature(object = "MultivariateDistribution", 
                          fun = "missing", 
                          cond = "missing"),
@@ -151,10 +198,12 @@
         if(is.null(low)) low <- -Inf
         if(is.null(upp)) upp <- Inf
         xsim <- xsim[xsim >= low & xsim <= upp]
-        if(useApply)        
-            return(mean(sapply(xsim, fun, ...)))
+        dotsFun <- .filterFunargs(list(...), fun, neg=FALSE)
+        funwD <- function(x) do.call(fun,c(list(x),dotsFun))
+        if(useApply)
+            return(mean(sapply(X=xsim, FUN=funwD)))
         else
-            return(mean(fun(xsim, ...)))
+            return(mean(fun(xsim)))
     })
 
 setMethod("E", signature(object = "AbscontDistribution", 
@@ -164,17 +213,23 @@
              rel.tol= getdistrExOption("ErelativeTolerance"), 
              lowerTruncQuantile = getdistrExOption("ElowerTruncQuantile"), 
              upperTruncQuantile = getdistrExOption("EupperTruncQuantile"), 
-             IQR.fac = getdistrExOption("IQR.fac"), ...){
+             IQR.fac = getdistrExOption("IQR.fac"), ..., diagnostic = FALSE){
 
+        mc <- match.call()
+
+        dotsFun <- .filterFunargs(list(...), fun, neg=FALSE)
+        funwD <- function(x) do.call(fun,c(list(x),dotsFun))
+
         if(useApply){
-            integrand <- function(x, dfun, fun, ...){ 
-                sapply(x, fun, ...) * dfun(x) 
+            integrand <- function(x){
+                sapply(x, funwD) * d(object)(x)
             }
         }else{
-            integrand <- function(x, dfun, fun, ...){ 
-                fun(x, ...) * dfun(x) 
+            integrand <- function(x){
+                funwD(x) * d(object)(x)
             }
         }
+
         if(is.null(low)) low <- -Inf
         if(is.null(upp)) upp <- Inf
 
@@ -182,12 +237,17 @@
               upperTruncQuantile, IQR.fac)
         low <- Ib["low"]
         upp <- Ib["upp"]
-        
-        return(distrExIntegrate(f = integrand,
-                    lower = low,
-                    upper = upp, 
-                    rel.tol = rel.tol, 
-                    distr = object, fun = fun, dfun = d(object), ...))
+        res <- distrExIntegrate(f = integrand, lower = low,  upper = upp,
+                    rel.tol = rel.tol, distr = object, ...,
+                    diagnostic = diagnostic)
+
+        if(diagnostic){
+           diagn <- attr(res, "diagnostic")
+           diagn[["call"]] <- mc
+           attr(res, "diagnostic") <- diagn
+        }
+
+        return(res)
     })
 
 setMethod("E", signature(object = "DiscreteDistribution", 
@@ -198,16 +258,18 @@
         if(is.null(upp)) upp <- Inf
         supp <- support(object)
         supp <- supp[supp>=low & supp<=upp]
+        dotsFun <- .filterFunargs(list(...), fun, neg=FALSE)
+        funwD <- function(x) do.call(fun,c(list(x),dotsFun))
         if(useApply){
-            integrand <- function(x, dfun, fun, ...){
-                sapply(x, fun, ...) * dfun(x)
+            integrand <- function(x, dfun, fun){
+                sapply(X=x, FUN=fun) * dfun(x)
             }
         }else{
-            integrand <- function(x, dfun, fun, ...){
-                fun(x, ...) * dfun(x)
+            integrand <- function(x, dfun, fun){
+                fun(x) * dfun(x)
             }
         }
-        return(sum(integrand(x = supp, dfun = d(object), fun = fun, ...)))
+        return(sum(integrand(x = supp, dfun = d(object), fun = funwD)))
     })
 setMethod("E", signature(object = "LatticeDistribution", 
                          fun = "function", 
@@ -222,14 +284,16 @@
     function(object, fun, useApply = TRUE, Nsim = getdistrExOption("MCIterations"), 
              ...){
         x <- r(object)(Nsim)
+        dotsFun <- .filterFunargs(list(...), fun, neg=FALSE)
+        funwD <- function(x) do.call(fun,c(list(x),dotsFun))
         if(useApply)
-            erg <- apply(x, 1, fun, ...)
+            erg <- apply(X=x, MARGIN=1, FUN=funwD)
         else
-            erg <- t(fun(x, ...))
+            erg <- t(funwD(x))
         if(is.vector(erg))
             return(mean(erg))
         else{
-            res <- fun(x[1,], ...)
+            res <- funwD(x[1,])
             res[] <- rowMeans(erg)
             return(res)
         }
@@ -239,16 +303,18 @@
                          cond = "missing"),
     function(object, fun, useApply = TRUE, ...){
         supp <- support(object)
+        dotsFun <- .filterFunargs(list(...), fun, neg=FALSE)
+        funwD <- function(x) do.call(fun,c(list(x),dotsFun))
         if(useApply){
-            integrand <- function(x, fun, dfun, ...){ fun(x, ...) * dfun(t(x)) }
-            erg <- apply(supp, 1, integrand, fun = fun, dfun = d(object), ...)
+            integrand <- function(x, fun, dfun){ fun(x) * dfun(t(x)) }
+            erg <- apply(supp, 1, integrand, fun = funwD, dfun = d(object))
         }else{
-            erg <- t(fun(supp, ...) * d(object)(supp))
+            erg <- t(funwD(supp) * d(object)(supp))
         }
         if(is.vector(erg))
             return(sum(erg))
         else{
-            res <- fun(supp[1,], ...)
+            res <- funwD(supp[1,])
             res[] <- rowSums(erg)
             return(res)
         }
@@ -271,11 +337,12 @@
              rel.tol= getdistrExOption("ErelativeTolerance"), 
              lowerTruncQuantile = getdistrExOption("ElowerTruncQuantile"), 
              upperTruncQuantile = getdistrExOption("EupperTruncQuantile"), 
-             IQR.fac = getdistrExOption("IQR.fac"), ...){
-        fct <- function(x, dfun, cond){ x * dfun(x, cond) }
+             IQR.fac = getdistrExOption("IQR.fac"), ..., diagnostic = FALSE){
+        mc <- match.call()
+        fct <- function(x){ x * d(object)(x, cond) }
         if(useApply){
-            integrand <- function(x, dfun, cond){ 
-                return(sapply(x, fct, dfun = dfun, cond = cond))
+            integrand <- function(x){
+                return(sapply(x, fct))
             }
         }else{
             integrand <- fct
@@ -289,10 +356,19 @@
         low <- Ib["low"]
         upp <- Ib["upp"]
 
-        return(distrExIntegrate(integrand, 
-              lower = low, upper = upp, rel.tol = rel.tol, distr = object, 
-              dfun = d(object), cond = cond))
+        res <- distrExIntegrate(f=integrand,
+              lower = low, upper = upp, rel.tol = rel.tol, distr = object,
+              dfun = d(object), cond = cond, diagnostic = FALSE)
+
+        if(diagnostic){
+           diagn <- attr(res, "diagnostic")
+           diagn[["call"]] <- mc
+           attr(res, "diagnostic") <- diagn
+        }
+
+        return(res)
     })
+
 setMethod("E", signature(object = "DiscreteCondDistribution", 
                          fun = "missing",
                          cond = "numeric"),
@@ -313,20 +389,22 @@
                          cond = "numeric"),
     function(object, fun, cond, withCond = FALSE, useApply = TRUE, 
              low = NULL, upp = NULL, Nsim = getdistrExOption("MCIterations"), ...){
+        dotsFun <- .filterFunargs(list(...), fun, neg=FALSE)
+        funwD <- function(x) do.call(fun,c(list(x),dotsFun))
         xsim <- r(object)(Nsim, cond)
         if(is.null(low)) low <- -Inf
         if(is.null(upp)) upp <- Inf
         xsim <- xsim[xsim >= low & xsim <= upp]
         if(withCond){
             if(useApply)
-                res <- mean(sapply(xsim, fun, cond, ...))
+                res <- mean(sapply(xsim, funwD, cond))
             else
-                res <- mean(fun(xsim, cond, ...))
+                res <- mean(funwD(xsim, cond))
         }else{
             if(useApply)
-                res <- mean(sapply(xsim, fun, ...))
+                res <- mean(sapply(xsim, funwD))
             else
-                res <- mean(fun(xsim, ...))                
+                res <- mean(funwD(xsim))
         }
 
         return(res)
@@ -338,28 +416,20 @@
              rel.tol= getdistrExOption("ErelativeTolerance"), 
              lowerTruncQuantile = getdistrExOption("ElowerTruncQuantile"), 
              upperTruncQuantile = getdistrExOption("EupperTruncQuantile"), 
-             IQR.fac = getdistrExOption("IQR.fac"), ...){
-        if(withCond)
-            if(useApply){
-                integrand <- function(x, dfun, fun, cond, ...){ 
-                    sapply(x, fun, cond, ...) * dfun(x, cond) 
-                }        
-            }else{
-                integrand <- function(x, dfun, fun, cond, ...){ 
-                    fun(x, cond, ...) * dfun(x, cond) 
-                }        
-            }
-        else
-            if(useApply){
-                integrand <- function(x, dfun, fun, cond, ...){ 
-                    sapply(x, fun, ...) * dfun(x, cond) 
-                }
-            }else{
-                integrand <- function(x, dfun, fun, cond, ...){ 
-                    fun(x, ...) * dfun(x, cond) 
-                }
-            }
+             IQR.fac = getdistrExOption("IQR.fac"), ..., diagnostic = FALSE){
 
+        mc <- match.call()
+        dotsFun <- .filterFunargs(list(...), fun, neg=FALSE)
+
+        CondArg <- if(withCond) list(cond=cond) else NULL
+        funwD <- function(x) do.call(fun,c(list(x), CondArg,dotsFun))
+
+        if(useApply){
+           integrand <- function(x) sapply(x, funwD) * d(object)(x, cond)
+        }else{
+           integrand <- function(x) funwD(x) * d(object)(x, cond)
+        }
+
         if(is.null(low)) low <- -Inf
         if(is.null(upp)) upp <- Inf
 
@@ -368,41 +438,39 @@
         low <- Ib["low"]
         upp <- Ib["upp"]
         
-        return(distrExIntegrate(integrand, 
-                lower = low, upper = upp, rel.tol = rel.tol, distr = object, 
-                dfun = d(object), fun = fun, cond = cond, ...))
+        res <- distrExIntegrate(f=integrand,
+                lower = low, upper = upp, rel.tol = rel.tol, distr = object,
+                fun = funwD, ..., diagnostic = diagnostic)
+
+        if(diagnostic){
+           diagn <- attr(res, "diagnostic")
+           diagn[["call"]] <- mc
+           attr(res, "diagnostic") <- diagn
+        }
+
+        return(res)
     })
 setMethod("E", signature(object = "DiscreteCondDistribution", 
                          fun = "function",
                          cond = "numeric"),
     function(object, fun, cond, withCond = FALSE, useApply = TRUE, low = NULL, upp = NULL, ...){
+
+        dotsFun <- .filterFunargs(list(...), fun, neg=FALSE)
+
+        CondArg <- if(withCond) list(cond=cond) else NULL
+        funwD <- function(x) do.call(fun,c(list(x), CondArg,dotsFun))
+
         if(is.null(low)) low <- -Inf
         if(is.null(upp)) upp <- Inf
         supp <- support(object)(cond)
         supp <- supp[supp>=low & supp<=upp]
-        if(withCond){
-            if(useApply){
-                fct <- function(x, dfun, fun, cond, ...){ 
-                    sapply(x, fun, cond, ...) * dfun(x, cond) 
-                }
-            }else{
-                fct <- function(x, dfun, fun, cond, ...){ 
-                    fun(x, cond, ...) * dfun(x, cond) 
-                }
-            }
+
+        if(useApply){
+           integrand <- function(x) sapply(x, funwD) * d(object)(x, cond)
         }else{
-            if(useApply){
-                fct <- function(x, dfun, fun, cond, ...){ 
-                    sapply(x, fun, ...) * dfun(x, cond) 
-                }
-            }else{
-                fct <- function(x, dfun, fun, cond, ...){ 
-                    fun(x, ...) * dfun(x, cond) 
-                }
-            }
+           integrand <- function(x) funwD(x) * d(object)(x, cond)
         }
-        return(sum(fct(x = supp, dfun = d(object), fun = fun, 
-                       cond = cond, ...)))
+        return(sum(integrand(supp)))
     })
 
 
@@ -432,13 +500,23 @@
 setMethod("E", signature(object = "Beta", 
                          fun = "missing", 
                          cond = "missing"),
-    function(object, low = NULL, upp = NULL, ...){
+    function(object, low = NULL, upp = NULL, ..., diagnostic = FALSE){
+        mc <- match.call()
+
         if(!is.null(low)) if(low <= 0) low <- NULL
         if(!is.null(upp)) if(upp >= 1) upp <- NULL
-        if((!isTRUE(all.equal(ncp(object),0)))|| !is.null(low) || !is.null(upp))
-          return(E(as(object,"AbscontDistribution"), low=low, upp=upp, ...))
-        else
-          return(shape1(object)/(shape1(object)+shape2(object)))
+        if((!isTRUE(all.equal(ncp(object),0)))|| !is.null(low) || !is.null(upp)){
+
+          res <- E(as(object,"AbscontDistribution"), low=low, upp=upp, ..., diagnostic = diagnostic)
+
+          if(diagnostic){
+             diagn <- attr(res, "diagnostic")
+             diagn[["call"]] <- mc
+             attr(res, "diagnostic") <- diagn
+          }
+
+          return(res)
+        }else return(shape1(object)/(shape1(object)+shape2(object)))
     })
 ## source: https://mathworld.wolfram.com/BetaDistribution.html
 
@@ -470,7 +548,8 @@
 setMethod("E", signature(object = "Cauchy",
                          fun = "missing",
                          cond = "missing"),
-    function(object, low = NULL, upp = NULL, ...){
+    function(object, low = NULL, upp = NULL, ..., diagnostic = FALSE){
+    mc <- match.call()
     if(is.null(low) && is.null(upp))
         return(NA)
     else{
@@ -480,10 +559,19 @@
            if(upp == Inf) return(NA)
            else return(-Inf)
         }else{
-           return(if(upp == Inf) Inf else{
-                    getMethod("E", signature(object = "Cauchy",
-                                     fun = "function", cond = "missing"))(object,
-                                          fun=function(x)(x<upp & x>low)*1.0,...)})
+           if(upp == Inf) return(Inf) else{
+              res <- E(object, fun = function(x)1, low=low, upp= upp,
[TRUNCATED]

To get the complete diff run:
    svnlook diff /svnroot/distr -r 1276


More information about the Distr-commits mailing list