[Distr-commits] r1268 - in branches/distr-2.8/pkg/distrMod: R man tests/Examples

noreply at r-forge.r-project.org noreply at r-forge.r-project.org
Fri Aug 10 19:47:04 CEST 2018


Author: ruckdeschel
Date: 2018-08-10 19:47:04 +0200 (Fri, 10 Aug 2018)
New Revision: 1268

Modified:
   branches/distr-2.8/pkg/distrMod/R/asCvMVarianceQtl.R
   branches/distr-2.8/pkg/distrMod/man/internals.Rd
   branches/distr-2.8/pkg/distrMod/tests/Examples/distrMod-Ex.Rout.save
Log:
[distrMod] branch 2.8
+ tuned .CvMMDCovariance() in asCvMVarianceQtl.R for speed (like with kStepEstimator timings are 
   taken in comment ##-t-##) as the function .CvMMDCovariance was much slower than 
   .oldCvMMDCovariance for Generalized EVD with Mu Unknown... / now they are at equal there 


Modified: branches/distr-2.8/pkg/distrMod/R/asCvMVarianceQtl.R
===================================================================
--- branches/distr-2.8/pkg/distrMod/R/asCvMVarianceQtl.R	2018-08-10 17:37:37 UTC (rev 1267)
+++ branches/distr-2.8/pkg/distrMod/R/asCvMVarianceQtl.R	2018-08-10 17:47:04 UTC (rev 1268)
@@ -1,5 +1,5 @@
 .CvMMDCovarianceWithMux <- function(L2Fam, param, withplot = FALSE, withpreIC = FALSE,
-                            N = 400, rel.tol=.Machine$double.eps^0.3,
+                            N = 1021, rel.tol=.Machine$double.eps^0.3,
                             TruncQuantile = getdistrOption("TruncQuantile"),
                             IQR.fac = 15, ..., x=NULL){
    mu <- distribution(L2Fam)
@@ -21,7 +21,7 @@
 #             via quantile transformation
 .CvMMDCovariance<- function(L2Fam, param, mu = distribution(L2Fam),
                             withplot = FALSE, withpreIC = FALSE,
-                            N = 400, rel.tol=.Machine$double.eps^0.3,
+                            N = 1021, rel.tol=.Machine$double.eps^0.3,
                             TruncQuantile = getdistrOption("TruncQuantile"),
                             IQR.fac = 15,
                             ...){
@@ -89,13 +89,6 @@
            x.seq.2r <- 1-rev(x.seq.2l)
            x.seq.2 <- (x.seq.2l+x.seq.2r)/2
            x.seq.a <- c(x.seq.1[odd.1],x.seq.2[odd.2],x.seq.3[odd.3])
-#           x.seq.b <- c(x.seq.1,x.seq.2,x.seq.3)
-#           iN.1 <- 1:N1.1
-#           iN.2 <- N1.1+(1:N1.2)
-#           iN.3 <- N1.1+N1.2+(1:N1.3)
-#           riN.3 <- 1:N1.3
-#           riN.2 <- N1.3+1:N1.2
-#           riN.1 <- N1.3+N1.2+1:N1.1
           }else{
            x.seq <- seq(low,up, length = N)
           }
@@ -123,13 +116,10 @@
            x.mu.seq.2r <- 1-rev(x.mu.seq.2l)
            x.mu.seq.2 <- (x.mu.seq.2l+x.mu.seq.2r)/2
            x.mu.seq.a <- c(x.mu.seq.1[odd.1],x.mu.seq.2[odd.2],x.mu.seq.3[odd.3])
-#           x.mu.seq.b <- x.seq.b
-#           iN.mu.1 <- iN.1
-#           iN.mu.2 <- iN.2
-#           iN.mu.3 <- iN.3
-#           riN.mu.1 <- riN.1
-#           riN.mu.2 <- riN.2
-#           riN.mu.3 <- riN.3
+           x.mu.seq.b <- c(x.mu.seq.1,x.mu.seq.2,x.mu.seq.3)
+            iN.mu.1 <- seq(N1.1)
+            iN.mu.2 <- N1.1+seq(N1.2)
+            iN.mu.3 <- N1.1+N1.2+seq(N1.3)
           }else{
            x.mu.seq <- seq(low, up, length = N)
           }
@@ -159,10 +149,6 @@
       fqx <- function(x){qx <- q.l(distr)(x)
                          return(sapply(qx,function(y)evalRandVar(L2deriv.0, y)))
                         }
-      #Delta0x  <- sapply(x.seq.b,fqx)
-      #Delta0x.1 <- Delta0x[iN.1]
-      #Delta0x.2 <- Delta0x[iN.2]
-      #Delta0x.3 <- Delta0x[iN.3]
       Delta0x.1 <- sapply(x.seq.1,fqx)
       Delta0x.2 <- sapply(x.seq.2,fqx)
       Delta0x.3 <- sapply(x.seq.3,fqx)
@@ -174,8 +160,6 @@
       J1 <- do.call(myint, c(list(f=Delta1.q), dotsInt))
       Delta.0 <- function(x) Delta1.q(p(distr)(x))-J1
       J <- do.call(myint, c(list(f=function(x) (Delta1.q(x)-J1)^2),dotsInt))
-#      print(J1)
-#      print(J)
   }else{
       if(is(distr,"DiscreteDistribution")){
          L2x <- sapply(x.seq, function(x) evalRandVar(L2deriv.0, x))
@@ -211,10 +195,6 @@
 
       phiqx <- function(x){qx <- q.l(mu)(x)
                           return(phi(qx))}
-      #psi0qx <- sapply(rev(x.mu.seq.b), phiqx)
-      #psi0qx.1 <- psi0qx[riN.mu.1]
-      #psi0qx.2 <- psi0qx[riN.mu.2]
-      #psi0qx.3 <- psi0qx[riN.mu.3]
       psi0qx.1 <- sapply(rev(x.mu.seq.1), phiqx)
       psi0qx.2 <- sapply(rev(x.mu.seq.2), phiqx)
       psi0qx.3 <- sapply(rev(x.mu.seq.3), phiqx)
@@ -224,7 +204,7 @@
 
       psi0q <-  c(psi0q.1,psi0q.2,psi0q.3)
       psi.q1 <- approxfun(x.mu.seq.a, psi0q, yleft = 0, yright = rev(psi0q)[1])
-      psi <- function(x) psi.q1(p(mu)(x))-psi1
+      psi.fct <- function(x) psi.q1(p(mu)(x))-psi1
    }else{
       if(is(mu,"DiscreteDistribution")&&is(distr,"DiscreteDistribution")){
          if(!all(support(mu)==support(distr))){
@@ -238,10 +218,8 @@
          }
          psi1 <- sum(pprob.mu*Delta.mu*prob.mu)
          psi0 <- cumsum(rev(Delta.mu*prob.mu))
-#         psi1 <- psi0[1]
          psi0 <- rev(psi0)-psi1
-         psi <- approxfun(x.mu.seq, psi0, yleft = -psi1, yright = -psi1)
-#   print(sapply(x.mu.seq,psi))
+         psi.fct <- approxfun(x.mu.seq, psi0, yleft = -psi1, yright = -psi1)
       }else{
    ## integrand phi x Ptheta in formula (51) [ibid]
          phi1 <- function(x) phi(x) * p(distr)(x)
@@ -251,56 +229,57 @@
          psi0 <- sapply(x.mu.seq, function(X){ fct <- function(y) phixy(x=X,y=y)
                                         return(E(object=mu, fun = fct))})
          psi.1 <- approxfun(x.mu.seq, psi0, yleft = 0, yright = rev(psi0)[1])
-         psi <- function(x) psi.1(x)-psi1
-#         if(is(distr,"DiscreteDistribution"))
-#            psi <- function(x) (psi.1(x)-psi1) * (x %in% support(mu))
+         psi.fct <- function(x) psi.1(x)-psi1
       }
    }
  #  print(psi0)
    if(is(distr,"AbscontDistribution")){
-      psi.q <- function(x){qx <- q.l(distr)(x); return(psi(qx))}
+      psi.q <- function(x){qx <- q.l(distr)(x); return(psi.fct(qx))}
+   ## E2 = Cov_mu (psi)
 #      E2 <- do.call(myint, c(list(f=function(x)psi.q(x)^2),dotsInt))
       E1 <- do.call(myint, c(list(f=psi.q),dotsInt))
       E3 <- do.call(myint, c(list(f=function(x){
                                      qx <- q.l(distr)(x)
                                      L2qx <- sapply(qx,function(y)
                                                     evalRandVar(L2deriv.0, y))
-                                     return(psi(qx)*L2qx)
+                                     return(psi.fct(qx)*L2qx)
                                     }), dotsInt))
-      psi.01 <- function(x) (psi(x)-E1)/E3
+      psi.01.f <- function(x) (psi.fct(x)-E1)/E3
       E4 <- do.call(myint, c(list(f=function(x) (psi.q(x)-E1)^2/E3^2),dotsInt))
   }else{
       if(is(distr,"DiscreteDistribution")){
+   ## E2 = Cov_mu (psi)
 #         E2 <- sum(psi0^2*prob)
-         psi0 <- sapply(x.seq, psi)
+         psi0 <- sapply(x.seq, psi.fct)
 
          E1 <- sum(psi0*prob)
          E3 <- sum(psi0*L2x*prob)
          psi.01d <- (psi0-E1)/E3
          E4 <- sum(psi.01d^2*prob)
-         psi.01 <- function(x) (psi(x)-E1)/E3*liesInSupport(distr,x)
+         psi.01.f <- function(x) (psi.fct(x)-E1)/E3*liesInSupport(distr,x)
       }else{
+   ## E2 = Cov_mu (psi)
 #         E2 <- E(object=distr, fun = function(x) psi(x)^2)
          L2x  <- function(x,y)  (x<=y)*evalRandVar(L2deriv.0, x)
-         E1 <- E(object=distr, fun = psi )
-         E3 <- E(object=distr, fun = function(x) psi(x)*evalRandVar(L2deriv.0, x))
-         psi.01 <- function(x) (psi(x) - E1)/E3
-         E4 <- E(object=distr, fun = function(x) psi.01(x)^2)
+         E1 <- E(object=distr, fun = psi.fct )
+         E3 <- E(object=distr, fun = function(x) psi.fct(x)*evalRandVar(L2deriv.0, x))
+         psi.01.f <- function(x) (psi.fct(x) - E1)/E3
+         E4 <- E(object=distr, fun = function(x) psi.01.f(x)^2)
       }
    }
-   ## E2 = Cov_mu (psi)
 
 #   ### control: centering & standardization
    if(withplot)
        { dev.new() #windows()
          x0.seq <- x.seq
          if(is(distr,"AbscontDistribution")) x0.seq <- q.l(distr)(x.seq)
-         plot(x0.seq, psi.01(x0.seq),
+         plot(x0.seq, psi.01.f(x0.seq),
                      type = if(is(distr,"DiscreteDistribution")) "p" else "l")
        }
-   psi.01 <- EuclRandVariable(Map = list(psi.01), Domain = Reals())
 
+   psi.01 <- EuclRandVariable(Map = list(psi.01.f), Domain = Reals())
 
+
       }else{
 
    ## multivariate case
@@ -314,32 +293,18 @@
 
    for(i in 1:Dim)
        { if(is(distr,"AbscontDistribution")){
-            #fct0.q <- sapply(x.seq.b, function(x){qx <- q.l(distr)(x); return(L2deriv.0 at Map[[i]](qx))})
-            #fct0.q1 <- fct0.q[iN.1]
-            #fct0.q2 <- fct0.q[iN.2]
-            #fct0.q3 <- fct0.q[iN.3]
             fct0.q1 <- sapply(x.seq.1, function(x){qx <- q.l(distr)(x); return(L2deriv.0 at Map[[i]](qx))})
             fct0.q2 <- sapply(x.seq.2, function(x){qx <- q.l(distr)(x); return(L2deriv.0 at Map[[i]](qx))})
             fct0.q3 <- sapply(x.seq.3, function(x){qx <- q.l(distr)(x); return(L2deriv.0 at Map[[i]](qx))})
-            #print(fct0)
             Delta0.q1 <-  h0/100*.csimpsum(fct0.q1)
             Delta0.q2 <-  rev(Delta0.q1)[1]+h0*.csimpsum(fct0.q2)
             Delta0.q3 <-  rev(Delta0.q2)[1]+h0/100*.csimpsum(fct0.q3)
             Delta0.q <- c(Delta0.q1,Delta0.q2,Delta0.q3)
             Delta1.q <- approxfun(x.seq.a, Delta0.q, yleft = 0, yright = 0)
-            Delta <- function(x) Delta1.q(p(distr)(x))
-            Map.Delta[[i]] <- Delta
+            Map.Delta[[i]] <- function(x) Delta1.q(p(distr)(x))
             env.i <- environment(Map.Delta[[i]]) <- new.env()
             assign("i", i, envir=env.i)
-            assign("fct0.q1", fct0.q1, envir=env.i)
-            assign("fct0.q2", fct0.q2, envir=env.i)
-            assign("fct0.q3", fct0.q3, envir=env.i)
-            assign("Delta0.q1", Delta0.q1, envir=env.i)
-            assign("Delta0.q2", Delta0.q2, envir=env.i)
-            assign("Delta0.q3", Delta0.q3, envir=env.i)
-            assign("Delta0.q", Delta0.q, envir=env.i)
             assign("Delta1.q", Delta1.q, envir=env.i)
-            assign("Delta", Delta, envir=env.i)
          }else{
             if(is(distr,"DiscreteDistribution")){
                L2x <- sapply(x.seq, function(x) evalRandVar(L2deriv.0, x)[i])
@@ -384,9 +349,14 @@
 
 
    ## J = Var_Ptheta Delta
+##-t-## print(system.time({
    J1 <- E(object=distr, fun = Delta)
+##-t-## }))
    Delta.0 <- Delta - J1
+
+##-t-## print(system.time({
    J <- E(object=distr, fun = Delta.0 %*%t(Delta.0))
+##-t-## }))
    ### CvM-IC phi
    phi <- as(solve(J)%*%Delta.0,"EuclRandVariable")
 
@@ -394,25 +364,21 @@
 
    Map.phi1 <- vector("list",Dim)
    for(i in 1:Dim)
-       { Map.phi1[[i]] <- function(x) evalRandVar(phi,x)[i] * p(distr)(x)
+       { Map.phi1[[i]] <- function(x) phi at Map[[i]](x)* p(distr)(x)
          env.i <- environment(Map.phi1[[i]]) <- new.env()
          assign("i", i, envir=env.i)
          }
 
    phi1 <- EuclRandVariable(Map = Map.phi1, Domain = Reals())
+##-t-## print(system.time({
    psi1 <- E(object=mu, fun = phi1)
+##-t-## }))
 
-#   for(i in 1:Dim)
-#       { Map.phi1[[i]] <- function(x) evalRandVar(phi,x)[i]
-#         env.i <- environment(Map.phi1[[i]]) <- new.env()
-#         assign("i", i, envir=env.i)
-#         }
-
-
    ## obtaining IC psi  (formula (51))
    Map.psi <- vector("list",Dim)
 
 
+##-t-##  print(system.time({
    for(i in 1:Dim)
      {
 
@@ -420,29 +386,20 @@
        assign("i", i, envir=env.i)
 
        if(is(mu,"AbscontDistribution")){
-            fct01.q <- function(x){qx <- q.l(mu)(x);
-                                   return(evalRandVar(phi,qx)[i])}
-            #fct0.q <- sapply(rev(x.mu.seq.b),fct01.q)
-            #fct0.q1 <-  fct0.q[riN.mu.1]
-            #fct0.q2 <-  fct0.q[riN.mu.2]
-            #fct0.q3 <-  fct0.q[riN.mu.3]
-            fct0.q1 <-  sapply(rev(x.mu.seq.1),fct01.q)
-            fct0.q2 <-  sapply(rev(x.mu.seq.2),fct01.q)
-            fct0.q3 <-  sapply(rev(x.mu.seq.3),fct01.q)
+            qxm <- q.l(mu)(x.mu.seq.b)
+
+##-t-##  print(system.time({
+            fct0.qq <- sapply(qxm, phi at Map[[i]])
+##-t-##   }))
+            fct0.q1 <-  rev(fct0.qq[iN.mu.1])
+            fct0.q2 <-  rev(fct0.qq[iN.mu.2])
+            fct0.q3 <-  rev(fct0.qq[iN.mu.3])
             phi0.q3 <-  h0.mu/100*rev(.csimpsum(fct0.q3))
             phi0.q2 <-  phi0.q3[1]+h0.mu*rev(.csimpsum(fct0.q2))
             phi0.q1 <-  phi0.q2[1]+h0.mu/100*rev(.csimpsum(fct0.q1))
             phi0.q <- c(phi0.q1,phi0.q2,phi0.q3)
             phi0a.q <- approxfun(x.mu.seq.a, phi0.q, yleft = 0, yright = rev(phi0.q)[1])
             psi0 <- function(x) {phi0a.q(p(mu)(x))-psi1[i]}
-
-            assign("fct01.q", fct01.q, envir=env.i)
-            assign("fct0.q1", fct0.q1, envir=env.i)
-            assign("fct0.q2", fct0.q2, envir=env.i)
-            assign("fct0.q3", fct0.q3, envir=env.i)
-            assign("phi0.q1", phi0.q1, envir=env.i)
-            assign("phi0.q2", phi0.q2, envir=env.i)
-            assign("phi0.q3", phi0.q3, envir=env.i)
             assign("phi0.q", phi0.q, envir=env.i)
             assign("phi0a.q", phi0a.q, envir=env.i)
             assign("psi0", psi0, envir=env.i)
@@ -476,23 +433,25 @@
             }
        }
 
-#       env.i0 <- environment(phi1) <- new.env()
-#       assign("i", i, envir=env.i0)
-
        Map.psi[[i]] <- psi0
        environment(Map.psi[[i]]) <- env.i
 
     }
+##-t-##  }))
 #   print(Map.psi)
    psi <-  EuclRandVariable(Map = Map.psi, Domain = Reals())
 
-   E2 <- E(object=distr, fun = psi %*%t(psi))
+#   E2 <- E(object=distr, fun = psi %*%t(psi))
    ## E2 = Cov_mu (psi)
 
    ### control: centering & standardization
    L2deriv.0 <- L2Fam at L2deriv[[1]]
+##-t-##  print(system.time({
    E1 <- E(object=distr, fun = psi )
+##-t-##  }))
+##-t-##  print(system.time({
    E3 <- E(object=distr, fun = psi %*%t(L2deriv.0))
+##-t-##  }))
    psi.0 <- psi - E1
    psi.01 <- as(solve(E3)%*%psi.0,"EuclRandVariable")
    if(withplot)
@@ -503,15 +462,22 @@
            plot(x0.mu.seq, sapply(x0.mu.seq,psi.01 at Map[[i]]),
                      type = if(is(distr,"DiscreteDistribution")) "p" else "l")
          }}
+##-t-##  print(system.time({
    E4 <- E(object=distr, fun = psi.01 %*%t(psi.01))
+##-t-##  }))
    }
   E4 <- PosSemDefSymmMatrix(E4)
 
   psi <-  EuclRandVarList(psi.01)
-
   if(onedim){
-     fct1 <- psi[[1]]@Map[[1]]
-     psi[[1]]@Map[[1]] <- function(x) fct1(x)*liesInSupport(distr,x,checkFin=TRUE)
+     fctl <- vector("list",1)
+     fct1 <- psi.01 at Map[[1]]
+     fctl[[1]] <- function(x) fct1(x)*liesInSupport(distr,x,checkFin=TRUE)
+     env.1 <- environment(fctl[[1]]) <- new.env()
+     assign("distr", distr, env.1)
+     assign("fct1", fct1, env.1)
+     assign("psi", psi, env.1)
+     psi[[1]]@Map <- fctl
   }else{
      fctl <- vector("list",Dim)
      for(i in 1:Dim){
@@ -573,7 +539,6 @@
    up1.mu  <- m0.mu + IQR.fac * s0.mu
    low.mu <- max(low0.mu,low1.mu); up.mu <- min(up0.mu,up1.mu)
 
-
    if(is(distr,"DiscreteDistribution"))
        x.seq <-support(distr)
    else
@@ -738,9 +703,13 @@
 
 
    ## J = Var_Ptheta Delta
+##-t-##  print(system.time({
    J1 <- E(object=distr, fun = Delta)
+##-t-##  }))
    Delta.0 <- Delta - J1
+##-t-##  print(system.time({
    J <- E(object=distr, fun = Delta.0 %*%t(Delta.0))
+##-t-##  }))
    ### CvM-IC phi
    phi <- as(solve(J)%*%Delta.0,"EuclRandVariable")
 
@@ -754,7 +723,9 @@
          }
 
    phi1 <- EuclRandVariable(Map = Map.phi1, Domain = Reals())
+##-t-##  print(system.time({
    psi1 <- E(object=mu, fun = phi1)
+##-t-##  }))
 
 
    ## obtaining IC psi  (formula (51))
@@ -793,13 +764,19 @@
     }
    psi <-  EuclRandVariable(Map = Map.psi, Domain = Reals())
 
+##-t-##  print(system.time({
    E2 <- E(object=distr, fun = psi %*%t(psi))
+##-t-##  }))
    ## E2 = Cov_mu (psi)
 
    ### control: centering & standardization
    L2deriv <- L2Fam at L2deriv[[1]]
+##-t-##  print(system.time({
    E1 <- E(object=distr, fun = psi )
+##-t-##  }))
+##-t-##  print(system.time({
    E3 <- E(object=distr, fun = psi %*%t(L2deriv))
+##-t-##  }))
    psi.0 <- psi - E1
    psi.01 <- as(solve(E3)%*%psi.0,"EuclRandVariable")
    if(withplot)
@@ -808,7 +785,9 @@
            plot(x.mu.seq, sapply(x.mu.seq,psi.01 at Map[[i]]),
                      type = if(is(distr,"DiscreteDistribution")) "p" else "l")
          }}
+##-t-##  print(system.time({
    E4 <- E(object=distr, fun = psi.01 %*%t(psi.01))
+##-t-##  }))
    }
   E4 <- PosSemDefSymmMatrix(E4)
 
@@ -860,6 +839,21 @@
 .oldCvMMDCovariance(GF,par=ParamFamParameter(main=c(scale=2.3,shape=0.3)), withplot=TRUE, N = 100)
 .CvMMDCovariance(GF,par=ParamFamParameter(main=c(scale=2.3,shape=0.3)), withplot=TRUE, N = 100)
 
+P0 <- PoisFamily()
+B0 <- BinomFamily(size=8, prob=0.3)
+N0 <- NormLocationFamily();
+C0 <- CauchyLocationFamily()
+cls <- CauchyLocationScaleFamily();
+N1 <- NormScaleFamily()
+NS <- NormLocationScaleFamily(); paramP <- ParamFamParameter(name = "locscale", main = c("loc"=0,"scale"=1),trafo = diag(2));
+Els <- L2LocationScaleFamily(loc = 0, scale = 1,
+                  name = "Laplace Location and scale family",
+                  centraldistribution = DExp(),
+                  LogDeriv = function(x)  sign(x),
+                  FisherInfo = diag(2),
+                  trafo = diag(2))
+Nb <- NbinomwithSizeFamily()
+GF <- GammaFamily()
 system.time(print(.oldCvMMDCovariance(P0,par=ParamFamParameter("lambda",1))))
 system.time(print(.CvMMDCovariance(P0,par=ParamFamParameter("lambda",1))))
 system.time(print(.oldCvMMDCovariance(B0,par=ParamFamParameter("",.3))))
@@ -886,6 +880,35 @@
 system.time(print(.CvMMDCovariance(P0,par=ParamFamParameter("lambda",1),mu=Norm())))
 system.time(print(.oldCvMMDCovariance(Nb,par=ParamFamParameter(main=c(size=2.3,prob=0.3)),mu=Norm())))
 system.time(print(.CvMMDCovariance(Nb,par=ParamFamParameter(main=c(size=2.3,prob=0.3)),mu=Norm())))
+
+" ## Code capsulated as string to avoid interpretation during R CMD Check
+
+require(RobExtremes)
+Pf <- ParetoFamily()
+.oldCvMMDCovariance(Pf,par=param(Pf), withplot=TRUE)
+.CvMMDCovariance(Pf,par=param(Pf), withplot=TRUE) ## better trust this one
+GPDf <- GParetoFamily()
+.oldCvMMDCovariance(GPDf,par=param(GPDf), withplot=TRUE)
+.CvMMDCovariance(GPDf,par=param(GPDf), withplot=TRUE) ## better trust this one
+GEVf <- GEVFamily()
+.oldCvMMDCovariance(GEVf,par=param(GEVf), withplot=TRUE)
+.CvMMDCovariance(GEVf,par=param(GEVf), withplot=TRUE) ## better trust this one
+GEVuf <- GEVFamilyMuUnknown()
+.oldCvMMDCovariance(GEVuf,par=param(GEVuf), withplot=TRUE)
+.CvMMDCovariance(GEVuf,par=param(GEVuf), withplot=TRUE) ## better trust this one
+
+
+system.time(print(.oldCvMMDCovariance(Pf,par=param(Pf))))
+system.time(print(.CvMMDCovariance(Pf,par=param(Pf)))) ## better trust this one
+system.time(print(.oldCvMMDCovariance(GPDf,par=param(GPDf))))
+system.time(print(.CvMMDCovariance(GPDf,par=param(GPDf)))) ## better trust this one
+system.time(print(.oldCvMMDCovariance(GEVf,par=param(GEVf))))
+system.time(print(.CvMMDCovariance(GEVf,par=param(GEVf)))) ## better trust this one
+system.time(print(.oldCvMMDCovariance(GEVuf,par=param(GEVuf))))
+system.time(print(.CvMMDCovariance(GEVuf,par=param(GEVuf)))) ## better trust this one
+
+"
+
 set.seed(123)
 x <- rnorm(100)
 diF <- DiscreteDistribution(x)

Modified: branches/distr-2.8/pkg/distrMod/man/internals.Rd
===================================================================
--- branches/distr-2.8/pkg/distrMod/man/internals.Rd	2018-08-10 17:37:37 UTC (rev 1267)
+++ branches/distr-2.8/pkg/distrMod/man/internals.Rd	2018-08-10 17:47:04 UTC (rev 1268)
@@ -23,7 +23,7 @@
 
 .CvMMDCovariance(L2Fam, param, mu = distribution(L2Fam),
                  withplot = FALSE, withpreIC = FALSE,
-                 N = 400, rel.tol=.Machine$double.eps^0.3,
+                 N = 1021, rel.tol=.Machine$double.eps^0.3,
                  TruncQuantile = getdistrOption("TruncQuantile"), 
                  IQR.fac = 15, ...)
 .oldCvMMDCovariance(L2Fam, param, mu = distribution(L2Fam),
@@ -33,7 +33,7 @@
                  TruncQuantile = getdistrOption("TruncQuantile"),
                  IQR.fac = 15, ...)
 .CvMMDCovarianceWithMux(L2Fam, param, withplot = FALSE, withpreIC = FALSE,
-                        N = 400, rel.tol=.Machine$double.eps^0.3,
+                        N = 1021, rel.tol=.Machine$double.eps^0.3,
                         TruncQuantile = getdistrOption("TruncQuantile"),
                         IQR.fac = 15, ..., x=NULL)
 

Modified: branches/distr-2.8/pkg/distrMod/tests/Examples/distrMod-Ex.Rout.save
===================================================================
--- branches/distr-2.8/pkg/distrMod/tests/Examples/distrMod-Ex.Rout.save	2018-08-10 17:37:37 UTC (rev 1267)
+++ branches/distr-2.8/pkg/distrMod/tests/Examples/distrMod-Ex.Rout.save	2018-08-10 17:47:04 UTC (rev 1268)
@@ -280,6 +280,63 @@
 > base::assign(".dptime", (proc.time() - get(".ptime", pos = "CheckExEnv")), pos = "CheckExEnv")
 > base::cat("BinomFamily", base::get(".format_ptime", pos = 'CheckExEnv')(get(".dptime", pos = "CheckExEnv")), "\n", file=base::get(".ExTimings", pos = 'CheckExEnv'), append=TRUE, sep="\t")
 > cleanEx()
+> nameEx("CauchyLocationFamily")
+> ### * CauchyLocationFamily
+> 
+> flush(stderr()); flush(stdout())
+> 
+> base::assign(".ptime", proc.time(), pos = "CheckExEnv")
+> ### Name: CauchyLocationFamily
+> ### Title: Generating function for Cauchy location families
+> ### Aliases: CauchyLocationFamily
+> ### Keywords: models
+> 
+> ### ** Examples
+> 
+> (C1 <- CauchyLocationFamily())
+An object of class "CauchyLocationFamily"
+### name:	Cauchy location family
+
+### distribution:	Distribution Object of Class: Cauchy
+ location: 0
+ scale: 1
+
+### param:	An object of class "ParamFamParameter"
+name:	loc
+loc:	0
+trafo:
+    loc
+loc   1
+
+### props:
+[1] "The Cauchy location family is invariant under"
+[2] "the group of transformations 'g(x) = x + loc'"
+[3] "with location parameter 'loc'"                
+> plot(C1)
+> FisherInfo(C1)
+An object of class "PosDefSymmMatrix"
+    loc
+loc 0.5
+> ### need smaller integration range:
+> checkL2deriv(C1)
+precision of centering:	 9.75782e-18 
+precision of Fisher information:
+              loc
+loc -1.061373e-13
+precision of Fisher information - relativ error [%]:
+              loc
+loc -2.122746e-11
+condition of Fisher information:
+[1] 1
+$maximum.deviation
+[1] 1.061373e-13
+
+> 
+> 
+> 
+> base::assign(".dptime", (proc.time() - get(".ptime", pos = "CheckExEnv")), pos = "CheckExEnv")
+> base::cat("CauchyLocationFamily", base::get(".format_ptime", pos = 'CheckExEnv')(get(".dptime", pos = "CheckExEnv")), "\n", file=base::get(".ExTimings", pos = 'CheckExEnv'), append=TRUE, sep="\t")
+> cleanEx()
 > nameEx("CauchyLocationScaleFamily")
 > ### * CauchyLocationScaleFamily
 > 
@@ -323,19 +380,19 @@
 > ### need smaller integration range:
 > distrExoptions("ElowerTruncQuantile"=1e-4,"EupperTruncQuantile"=1e-4)
 > checkL2deriv(C1)
-precision of centering:	 0 -0.02119711 
+precision of centering:	 9.97466e-18 -2e-04 
 precision of Fisher information:
+                loc         scale
+loc   -2.631895e-11 -3.577867e-17
+scale -3.577867e-17 -2.000000e-04
+precision of Fisher information - relativ error [%]:
                 loc       scale
-loc   -3.137524e-05  0.00000000
-scale  0.000000e+00 -0.02118143
-precision of Fisher information - relativ error [%]:
-               loc     scale
-loc   -0.006275047       NaN
-scale          NaN -4.236286
+loc   -5.263789e-09        -Inf
+scale          -Inf -0.03999999
 condition of Fisher information:
 [1] 1
 $maximum.deviation
-[1] 0.02119711
+[1] 2e-04
 
 > distrExoptions("ElowerTruncQuantile"=1e-7,"EupperTruncQuantile"=1e-7)
 > 
@@ -402,7 +459,7 @@
         dimnames = list(nms, nms0))
     list(fval = fval0, mat = mat0)
 }
-<bytecode: 0x10138e10>
+<bytecode: 0x0e596b80>
 Trafo / derivative matrix at which estimate was produced:
        scale shape
 shape  0.000     1
@@ -613,12 +670,12 @@
 function (x) 
 {
     y <- 0 * x
-    inS <- liesInSupport(distr.0, x)
+    inS <- liesInSupport(distr.0, x, checkFin = TRUE)
     y[inS] <- ((x[inS] - 0)/scale * LogDeriv((x[inS] - 0)/c(scale = 1)) - 
         1)/c(scale = 1)
     return(y)
 }
-<environment: 0x0b361f70>
+<environment: 0x0ed32a30>
 
 > checkL2deriv(E1)
 precision of centering:	 -2.04266e-06 
@@ -741,7 +798,7 @@
  shape: 1
  scale: 1
 
-### param:	An object of class "ParamFamParameter"
+### param:	An object of class "ParamWithScaleAndShapeFamParameter"
 name:	scale and shape
 scale:	1
 shape:	1
@@ -806,8 +863,8 @@
 Slot "fct":
 function (x) 
 QuadFormNorm(x, A = A)
-<bytecode: 0x0b47c878>
-<environment: 0x0b47cab8>
+<bytecode: 0x0f550478>
+<environment: 0x0f5502b8>
 
 > 
 > ## The function is currently defined as
@@ -1120,12 +1177,12 @@
 function (x) 
 {
     y <- 0 * x
-    inS <- liesInSupport(distr.0, x)
+    inS <- liesInSupport(distr.0, x, checkFin = TRUE)
     y[inS] <- ((x[inS] - 0)/scale * LogDeriv((x[inS] - 0)/c(meanlog = 1)) - 
         1)/c(meanlog = 1)
     return(y)
 }
-<environment: 0x0ecb2318>
+<environment: 0x109f06e8>
 
 > checkL2deriv(L1)
 precision of centering:	 -0.003003394 
@@ -1192,42 +1249,42 @@
 > ### need smaller integration range:
 > distrExoptions("ElowerTruncQuantile"=1e-4,"EupperTruncQuantile"=1e-4)
 > checkL2deriv(L1)
-precision of centering:	 -2.264121e-17 -1.41228 
+precision of centering:	 2.264121e-17 -0.5873198 
 precision of Fisher information:
-              location        scale
-location -1.600022e-01 4.454916e-17
-scale     4.454916e-17 8.149930e-01
+              location         scale
+location -1.600022e-01 -6.600106e-18
+scale    -6.600106e-18 -8.349279e-01
 precision of Fisher information - relativ error [%]:
-          location    scale
-location -48.00065      Inf
-scale          Inf 56.99427
+          location     scale
+location -48.00065      -Inf
+scale         -Inf -58.38836
 condition of Fisher information:
 [1] 3.667949
 $maximum.deviation
-[1] 1.41228
+[1] 0.8349279
 
 > distrExoptions("ElowerTruncQuantile"=1e-7,"EupperTruncQuantile"=1e-7)
 > ##
 > set.seed(123)
 > x <- rlogis(100,location=1,scale=2)
 > CvMMDEstimator(x, L1)
-Evaluations of Minimum CvM distance estimate ( mu = emp. cdf ) :
-----------------------------------------------------------------
-An object of class "Estimate" 
+Evaluations of Minimum CvM distance estimate ( mu = model distr. ) :
+--------------------------------------------------------------------
+An object of class "CvMMDEstimate" 
 generated by call
   CvMMDEstimator(x = x, ParamFamily = L1)
 samplesize:   100
 estimate:
    location      scale  
-  0.9199612   1.9512792 
- (2.2045470) (0.6064406)
+  0.9691129   1.9522622 
+ (0.4035840) (0.1868609)
 asymptotic (co)variance (multiplied with samplesize):
-         location     scale
-location 486.0028 130.60897
-scale    130.6090  36.77702
+              location         scale
+location  1.628800e+01 -3.750942e-06
+scale    -3.750942e-06  3.491701e+00
 Criterion:
 CvM distance 
-  0.01414319 
+  0.01469024 
 > 
 > 
 > 
@@ -1242,7 +1299,8 @@
 > base::assign(".ptime", proc.time(), pos = "CheckExEnv")
 > ### Name: MCEstimate-class
 > ### Title: MCEstimate-class.
-> ### Aliases: MCEstimate-class criterion criterion,MCEstimate-method
+> ### Aliases: MCEstimate-class MDEstimate-class MLEstimate-class
+> ###   CvMMDEstimate-class criterion criterion,MCEstimate-method
 > ###   criterion.fct criterion.fct,MCEstimate-method
 > ###   startPar,MCEstimate-method method method,MCEstimate-method optimwarn
 > ###   optimwarn,MCEstimate-method optimReturn optimReturn,MCEstimate-method
@@ -1262,7 +1320,7 @@
 > MDEstimator(x, G)
 Evaluations of Minimum Kolmogorov distance estimate  :
 ------------------------------------------------------
-An object of class "Estimate" 
+An object of class "MDEstimate" 
 generated by call
   MDEstimator(x = x, ParamFamily = G)
 samplesize:   50
@@ -1275,7 +1333,7 @@
 > (m <- MLEstimator(x, G))
 Evaluations of Maximum likelihood estimate:
 -------------------------------------------
-An object of class "Estimate" 
+An object of class "MLEstimate" 
 generated by call
   MLEstimator(x = x, ParamFamily = G)
 samplesize:   50
@@ -1330,7 +1388,7 @@
 > MCEstimator(x = x, ParamFamily = G, criterion = negLoglikelihood)
 Evaluations of Minimum criterion estimate:
 ------------------------------------------
-An object of class "Estimate" 
+An object of class "MCEstimate" 
 generated by call
   MCEstimator(x = x, ParamFamily = G, criterion = negLoglikelihood)
 samplesize:   50
@@ -1359,7 +1417,7 @@
 
 Evaluations of Minimum Kolmogorov distance estimate:
 ----------------------------------------------------
-An object of class "Estimate" 
+An object of class "MCEstimate" 
 generated by call
   MCEstimator(x = x, ParamFamily = G, criterion = KolmogorovDist, 
     crit.name = "Kolmogorov distance")
@@ -1378,7 +1436,7 @@
 +             crit.name = "Total variation distance")
 Evaluations of Minimum Total variation distance estimate:
 ---------------------------------------------------------
-An object of class "Estimate" 
+An object of class "MCEstimate" 
 generated by call
   MCEstimator(x = x, ParamFamily = G, criterion = TotalVarDist, 
     crit.name = "Total variation distance")
@@ -1402,7 +1460,7 @@
 +             crit.name = "Hellinger Distance", startPar = c(1,2))
 Evaluations of Minimum Hellinger Distance estimate:
 ---------------------------------------------------
-An object of class "Estimate" 
+An object of class "MCEstimate" 
 generated by call
   MCEstimator(x = x, ParamFamily = G, criterion = HellingerDist, 
     crit.name = "Hellinger Distance", startPar = c(1, 2))
@@ -1449,7 +1507,7 @@
 > MDEstimator(x = x, ParamFamily = G, distance = KolmogorovDist)
 Evaluations of Minimum Kolmogorov distance estimate  :
 ------------------------------------------------------
-An object of class "Estimate" 
+An object of class "MDEstimate" 
 generated by call
   MDEstimator(x = x, ParamFamily = G, distance = KolmogorovDist)
 samplesize:   50
@@ -1463,7 +1521,7 @@
 > KolmogorovMDEstimator(x = x, ParamFamily = G)
 Evaluations of Minimum Kolmogorov distance estimate  :
 ------------------------------------------------------
-An object of class "Estimate" 
+An object of class "MDEstimate" 
 generated by call
   KolmogorovMDEstimator(x = x, ParamFamily = G)
 samplesize:   50
@@ -1474,11 +1532,11 @@
 Kolmogorov distance 
          0.07111522 
 > 
-> ## von Mises minimum distance estimator with default mu
+> ## von Mises minimum distance estimator with default mu = Mod
 > MDEstimator(x = x, ParamFamily = G, distance = CvMDist)
 Evaluations of Minimum CvM distance estimate ( mu = emp. cdf )  :
 -----------------------------------------------------------------
-An object of class "Estimate" 
+An object of class "CvMMDEstimate" 
 generated by call
   MDEstimator(x = x, ParamFamily = G, distance = CvMDist)
 samplesize:   50
@@ -1518,7 +1576,7 @@
 > MLEstimator(x, BinomFamily(size = 25))
 Evaluations of Maximum likelihood estimate:
 -------------------------------------------
-An object of class "Estimate" 
+An object of class "MLEstimate" 
 generated by call
   MLEstimator(x = x, ParamFamily = BinomFamily(size = 25))
 samplesize:   100
@@ -1548,7 +1606,7 @@
 > MLEstimator(x, PoisFamily())
 Evaluations of Maximum likelihood estimate:
 -------------------------------------------
-An object of class "Estimate" 
+An object of class "MLEstimate" 
 generated by call
   MLEstimator(x = x, ParamFamily = PoisFamily())
 samplesize:   2608
@@ -1573,7 +1631,7 @@
 > MLEstimator(x, NormLocationScaleFamily())
 Evaluations of Maximum likelihood estimate:
 -------------------------------------------
-An object of class "Estimate" 
+An object of class "MLEstimate" 
 generated by call
   MLEstimator(x = x, ParamFamily = NormLocationScaleFamily())
 samplesize:   100
@@ -1606,7 +1664,7 @@
 > (res <- MLEstimator(x = x, ParamFamily = G))
 Evaluations of Maximum likelihood estimate:
 -------------------------------------------
-An object of class "Estimate" 
+An object of class "MLEstimate" 
 generated by call
   MLEstimator(x = x, ParamFamily = G)
 samplesize:   50
@@ -1813,20 +1871,20 @@
 > FisherInfo(N1.w)
 An object of class "PosSemDefSymmMatrix"
             size     prob
-size  0.03044946  -4.0000
+size  0.03044974  -4.0000
 prob -4.00000000 533.3333
 > checkL2deriv(N1.w)
 precision of centering:	 -6.245978e-06 0.001177892 
 precision of Fisher information:
               size          prob
-size -4.182531e-06  0.0008481424
+size -4.462854e-06  0.0008481424
 prob  8.481424e-04 -0.1601189384
 precision of Fisher information - relativ error [%]:
             size        prob
-size -0.01373598 -0.02120356
+size -0.01465646 -0.02120356
 prob -0.02120356 -0.03002230
 condition of Fisher information:
-[1] 1195572
+[1] 1194827
 $maximum.deviation
 [1] 0.1601189
 
@@ -1849,20 +1907,20 @@
 > FisherInfo(N2.w)
 An object of class "PosSemDefSymmMatrix"
              size         mean
-size 3.044946e-02     1600.091
+size 3.044974e-02     1600.091
 mean 1.600091e+03 85342933.607
 > checkL2deriv(N2.w)
 precision of centering:	 -6.245978e-06 -0.4711755 
 precision of Fisher information:
               size          mean
-size -4.182531e-06 -3.392695e-01
-mean -3.392695e-01 -2.562107e+04
+size -4.462854e-06 -3.392704e-01
+mean -3.392704e-01 -2.562107e+04
 precision of Fisher information - relativ error [%]:
             size        mean
-size -0.01373598 -0.02120313
-mean -0.02120313 -0.03002131
+size -0.01465646 -0.02120319
+mean -0.02120319 -0.03002131
 condition of Fisher information:
-[1] 1.89903e+11
+[1] 189784664099
 $maximum.deviation
 [1] 25621.07
 
@@ -2232,7 +2290,7 @@
         return(abs(x))
     else return(sqrt(colSums(x^2)))
 }
-<bytecode: 0x0d2cc2c0>
+<bytecode: 0x0c5c3cf0>
 <environment: namespace:distrMod>
 > name(EuclNorm)
 [1] "EuclideanNorm"
@@ -2267,7 +2325,7 @@
         return(abs(x))
     else return(sqrt(colSums(x^2)))
 }
-<bytecode: 0x0d2cc2c0>
+<bytecode: 0x0c5c3cf0>
 <environment: namespace:distrMod>
 
 > 
@@ -2663,7 +2721,7 @@
  shape: 2
  scale: 1
 
-### param:	An object of class "ParamFamParameter"
+### param:	An object of class "ParamWithScaleAndShapeFamParameter"
 name:	scale and shape
 scale:	1
 shape:	2
@@ -2750,8 +2808,8 @@
 Slot "fct":
 function (x) 
 QuadFormNorm(x, A = A0)
-<bytecode: 0x0ef13a38>
[TRUNCATED]

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


More information about the Distr-commits mailing list