[Distr-commits] r1417 - in pkg/distrEx: . R inst
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Mon Nov 27 07:08:27 CET 2023
Author: stamats
Date: 2023-11-27 07:08:26 +0100 (Mon, 27 Nov 2023)
New Revision: 1417
Modified:
pkg/distrEx/DESCRIPTION
pkg/distrEx/R/distrExIntegrate.R
pkg/distrEx/inst/NEWS
Log:
changed encoding from latin1 to UTF-8 as requested by CRAN
Modified: pkg/distrEx/DESCRIPTION
===================================================================
--- pkg/distrEx/DESCRIPTION 2023-09-19 15:57:19 UTC (rev 1416)
+++ pkg/distrEx/DESCRIPTION 2023-11-27 06:08:26 UTC (rev 1417)
@@ -1,6 +1,6 @@
Package: distrEx
-Version: 2.9.0
-Date: 2022-11-12
+Version: 2.9.1
+Date: 2023-11-27
Title: Extensions of Package 'distr'
Description: Extends package 'distr' by functionals, distances, and conditional distributions.
Depends: R(>= 3.4), methods, distr(>= 2.8.0)
@@ -11,7 +11,7 @@
"cph")))
ByteCompile: yes
License: LGPL-3
-Encoding: latin1
+Encoding: UTF-8
URL: http://distr.r-forge.r-project.org/
LastChangedDate: {$LastChangedDate$}
LastChangedRevision: {$LastChangedRevision$}
Modified: pkg/distrEx/R/distrExIntegrate.R
===================================================================
--- pkg/distrEx/R/distrExIntegrate.R 2023-09-19 15:57:19 UTC (rev 1416)
+++ pkg/distrEx/R/distrExIntegrate.R 2023-11-27 06:08:26 UTC (rev 1417)
@@ -1,240 +1,240 @@
-# Gauss-Legendre abscissas and weights
-# cf. for example Numerical Recipies in C (1992), p. 152
-
-#implementation in S:
-.GLawOld <- function(n){
- if(n %% 2 == 1) stop("n has to be an even number")
-
- m <- (n + 1)/2
- A <- numeric(n)
- W <- numeric(n)
- for(i in 1:floor(m)){
- z <- cos(pi*(i - 0.25)/(n + 0.5))
- repeat{
- p1 <- 1
- p2 <- 0
- for(j in 1:n){
- p3 <- p2
- p2 <- p1
- p1 <- ((2*j - 1)*z*p2 - (j-1)*p3)/j
- }
- pp <- n*(z*p1 - p2)/(z^2 - 1)
- z1 <- z
- z <- z - p1/pp
- if(abs(z - z1) < .Machine$double.eps) break
- }
- A[i] <- -z
- A[n + 1 - i] <- z
- W[i] <- 2/((1 - z^2)*pp^2)
- W[n + 1 - i] <- W[i]
- }
- cbind(A, W)
-}
-
-## new: interface to C P.R. 01-04-06
-
-.GLaw <- function(n){
- if(n %% 2 == 1) stop("n has to be an even number")
-
- A <- numeric(n)
- W <- numeric(n)
-
-# mm<-dyn.load("G:/rtest/GLaw.dll")
- erg<-.C(C_gauleg,n = as.integer(n),eps = as.double(.Machine$double.eps),
- A = as.double(A),W = as.double(W)) #, PACKAGE = "distrEx")
- ### PACKAGE ARGUMENT added P.R. 270507
- #### removed again 200417 / used registered symbol instead
-# dyn.unload("G:/rtest/GLaw.dll")
-#
-# P.R. 20140810: .Call interface instead of .C interface
-#
-# erg0 <- .Call("Gauleg", n, eps, PACKAGE="distrEx")
-# erg <- matrix(erg0,n,2); colnames(erg) <- c("A","W")
-#
- cbind(A=erg$A, W=erg$W)
-}
-
-
-if(FALSE){
-# code to produce the AW values stored in the namespace of distrEx
-##
-
-## timing code borrowed from base::system.time
-
- ppt <- function(y) {
- if (!is.na(y[4L]))
- y[1L] <- y[1L] + y[4L]
- if (!is.na(y[5L]))
- y[2L] <- y[2L] + y[5L]
- paste(formatC(y[1L:3L]), collapse = " ")
- }
-
-
-
-todo <- c(50, 100, 400, 500, 800, 1000, 4000, 5000, 8000, 10000, 40000, 50000, 80000, 100000)
-l <- length(todo)
-nE <- new.env()
-svncheckout <- "C:/rtest/distr"
-pkg <- file.path(svncheckout, "branches/distr-2.8/pkg/distrEx")
-sysdataFilename <- file.path(pkg, "R/sysdata.rda")
-load(sysdataFilename,envir=nE)
-
-gc()
-starttime <- proc.time()
-on.exit(message("Timing stopped at: ", ppt(proc.time() - starttime)))
-
-lasttime <- starttime
-for(gridsize.i in seq(todo)){
- cat("Gridpoint i =", gridsize.i, ", order = ", todo[gridsize.i],", time needed: ")
- res <- distrEx:::.GLaw(todo[gridsize.i])
- newtime <- proc.time()
- timN <- structure(newtime - lasttime, class = "proc_time")
- lasttime <- newtime
- cat(paste(round(timN,3)), "\n")
- nam <- paste(".AW",as.character(todo[gridsize.i]), sep = ".")
- assign(x=nam, value=res, envir=nE)
-}
-
- timN <- structure(proc.time() - starttime, class = "proc_time")
- cat("Time altogether:", paste(round(timN,3)), "\n")
-
-rm(".AW.100000", envir=nE)
-what <- ls(all=TRUE, env=nE)
-for(item in what) {cat(item, ":\n");print(object.size(get(item, envir=nE)))}
-on.exit()
-
-save(list=what,file=sysdataFilename,envir=nE)
-rm(nE)
-}
-
-GLIntegrate <- function(f, lower, upper, order = 500, ...){
- if(order %in% c(50, 100, 400, 500, 800, 1000, 4000, 5000, 8000, 10000,
- 40000, 50000, 80000, 100000))
- AW <- getFromNamespace(paste(".AW", as.character(order),
- sep = "."), ns = "distrEx")
- else
- AW <- .GLaw(order)
-
- # transformation to [lower, upper]
- xl <- (upper - lower)/2
- W <- xl*AW[,2]
- A <- xl*AW[,1] + (lower + upper)/2
-
- res <- W*c(f(A, ...))
- sum(res)
-}
-
-distrExIntegrate <- function(f, lower, upper, subdivisions = 100,
- rel.tol = .Machine$double.eps^0.25,
- abs.tol = rel.tol, stop.on.error = TRUE,
- distr, order = .distrExOptions$GLIntegrateOrder,
- ..., diagnostic = FALSE){
- mc <- match.call()
- time <- proc.time()
-
- ## taken from base::system.time
- ppt <- function(y) {
- if (!is.na(y[4L]))
- y[1L] <- y[1L] + y[4L]
- if (!is.na(y[5L]))
- y[2L] <- y[2L] + y[5L]
- paste(formatC(y[1L:3L]), collapse = " ")
- }
-
- dots <- list(...)
- dotsFun <- .filterFunargs(dots,f)
- funwD <- function(x) do.call(f,c(list(x), dotsFun))
-
- on.exit(message("Timing stopped at: ", ppt(proc.time() -
- time)))
- dotsInt <- if(length(names(dots))) dots[names(dots)%in% names(formals(integrate))] else NULL
- res <- try(do.call(integrate, c(list(funwD, lower = lower, upper = upper, rel.tol = rel.tol,
- abs.tol = abs.tol, stop.on.error = stop.on.error), dotsInt)),
- silent = TRUE)
-
- # if integrate fails => Gauß-Legendre integration
- if(!is(res,"try-error")){
- val <- res$value
- if(diagnostic){
- diagn <- list(call = mc, method = "integrate",
- args = c(list(lower=lower, upper = upper, rel.tol = rel.tol,
- abs.tol = abs.tol,
- stop.on.error = stop.on.error),list(...)),
- result = res)
- res <- val
- }else res <- val
- }else{
- errmess <- res
- Zi <- 1
- if(lower >= upper){
- lo <- lower
- lower <- upper
- upper <- lo
- Zi <- -1
- }
- if(!missing(distr)){
- q.lDots <- NULL
- if(length(names(dots))) {
- q.lDots <- dots[names(dots) %in% names(formals(q.l(distr)))]
- q.lDots[["p"]] <- q.lDots[["lower.tail"]] <- NULL
- }
- }
-
- if(!is.finite(lower))
- if(missing(distr)) stop(res)
- else{
- lower <- do.call(q.l(distr),
- c(list(.distrExOptions$GLIntegrateTruncQuantile),q.lDots))
-# value <- c(...)
-# if(is(distr, "UnivariateCondDistribution"))
-# lower <- q.l(distr)(.distrExOptions$GLIntegrateTruncQuantile,
-# cond = value$cond)
-# else
-# lower <- q.l(distr)(.distrExOptions$GLIntegrateTruncQuantile)
- }
- if(!is.finite(upper))
- if(missing(distr)) stop(res)
- else{
- q.lArgs <- if("lower.tail" %in% names(formals(distr at q)))
- list(p=.distrExOptions$GLIntegrateTruncQuantile, lower.tail=FALSE) else
- list(p=1-.distrExOptions$GLIntegrateTruncQuantile)
- q.lArgs <- c(q.lArgs, q.lDots)
- upper <- do.call(q.l(distr),q.lArgs)
-# value <- c(...)
-# if(is(distr, "UnivariateCondDistribution")){
-# if("lower.tail" %in% names(formals(distr at q)))
-# upper <- q.l(distr)(.distrExOptions$GLIntegrateTruncQuantile,
-# cond = value$cond, lower.tail = FALSE)
-# else
-# upper <- q.l(distr)(1 - .distrExOptions$GLIntegrateTruncQuantile,
-# cond = value$cond)
-# }else{
-# if("lower.tail" %in% names(formals(distr at q)))
-# upper <- q.l(distr)(.distrExOptions$GLIntegrateTruncQuantile,
-# lower.tail = FALSE)
-# else
-# upper <- q.l(distr)(1 - .distrExOptions$GLIntegrateTruncQuantile)
-# }
- }
- dotsGLInt <- NULL
- if(length(names(dots))) dotsGLInt <- dots[names(dots)%in% names(formals(GLIntegrate))]
- res <- Zi* do.call(GLIntegrate,c(list(f = funwD, lower = lower, upper = upper,
- order = order),dotsGLInt))
- if(diagnostic){
- diagn <- list(call = mc, method = "GLIntegrate",
- args = c(list(lower=lower, upper=upper, order=order),
- list(...)),
- result = list(GLIresult = res, errorMessage = errmess),
- distrExOptions = .distrExOptions)
- }
- }
-
- new.time <- proc.time()
- on.exit()
- if(diagnostic){
- diagn$time <- structure(new.time - time, class = "proc_time")
- attr(res,"diagnostic") <- diagn
- class(attr(res,"diagnostic"))<- "DiagnosticClass"
- }
- return(res)
-}
+# Gauss-Legendre abscissas and weights
+# cf. for example Numerical Recipies in C (1992), p. 152
+
+#implementation in S:
+.GLawOld <- function(n){
+ if(n %% 2 == 1) stop("n has to be an even number")
+
+ m <- (n + 1)/2
+ A <- numeric(n)
+ W <- numeric(n)
+ for(i in 1:floor(m)){
+ z <- cos(pi*(i - 0.25)/(n + 0.5))
+ repeat{
+ p1 <- 1
+ p2 <- 0
+ for(j in 1:n){
+ p3 <- p2
+ p2 <- p1
+ p1 <- ((2*j - 1)*z*p2 - (j-1)*p3)/j
+ }
+ pp <- n*(z*p1 - p2)/(z^2 - 1)
+ z1 <- z
+ z <- z - p1/pp
+ if(abs(z - z1) < .Machine$double.eps) break
+ }
+ A[i] <- -z
+ A[n + 1 - i] <- z
+ W[i] <- 2/((1 - z^2)*pp^2)
+ W[n + 1 - i] <- W[i]
+ }
+ cbind(A, W)
+}
+
+## new: interface to C P.R. 01-04-06
+
+.GLaw <- function(n){
+ if(n %% 2 == 1) stop("n has to be an even number")
+
+ A <- numeric(n)
+ W <- numeric(n)
+
+# mm<-dyn.load("G:/rtest/GLaw.dll")
+ erg<-.C(C_gauleg,n = as.integer(n),eps = as.double(.Machine$double.eps),
+ A = as.double(A),W = as.double(W)) #, PACKAGE = "distrEx")
+ ### PACKAGE ARGUMENT added P.R. 270507
+ #### removed again 200417 / used registered symbol instead
+# dyn.unload("G:/rtest/GLaw.dll")
+#
+# P.R. 20140810: .Call interface instead of .C interface
+#
+# erg0 <- .Call("Gauleg", n, eps, PACKAGE="distrEx")
+# erg <- matrix(erg0,n,2); colnames(erg) <- c("A","W")
+#
+ cbind(A=erg$A, W=erg$W)
+}
+
+
+if(FALSE){
+# code to produce the AW values stored in the namespace of distrEx
+##
+
+## timing code borrowed from base::system.time
+
+ ppt <- function(y) {
+ if (!is.na(y[4L]))
+ y[1L] <- y[1L] + y[4L]
+ if (!is.na(y[5L]))
+ y[2L] <- y[2L] + y[5L]
+ paste(formatC(y[1L:3L]), collapse = " ")
+ }
+
+
+
+todo <- c(50, 100, 400, 500, 800, 1000, 4000, 5000, 8000, 10000, 40000, 50000, 80000, 100000)
+l <- length(todo)
+nE <- new.env()
+svncheckout <- "C:/rtest/distr"
+pkg <- file.path(svncheckout, "branches/distr-2.8/pkg/distrEx")
+sysdataFilename <- file.path(pkg, "R/sysdata.rda")
+load(sysdataFilename,envir=nE)
+
+gc()
+starttime <- proc.time()
+on.exit(message("Timing stopped at: ", ppt(proc.time() - starttime)))
+
+lasttime <- starttime
+for(gridsize.i in seq(todo)){
+ cat("Gridpoint i =", gridsize.i, ", order = ", todo[gridsize.i],", time needed: ")
+ res <- distrEx:::.GLaw(todo[gridsize.i])
+ newtime <- proc.time()
+ timN <- structure(newtime - lasttime, class = "proc_time")
+ lasttime <- newtime
+ cat(paste(round(timN,3)), "\n")
+ nam <- paste(".AW",as.character(todo[gridsize.i]), sep = ".")
+ assign(x=nam, value=res, envir=nE)
+}
+
+ timN <- structure(proc.time() - starttime, class = "proc_time")
+ cat("Time altogether:", paste(round(timN,3)), "\n")
+
+rm(".AW.100000", envir=nE)
+what <- ls(all=TRUE, env=nE)
+for(item in what) {cat(item, ":\n");print(object.size(get(item, envir=nE)))}
+on.exit()
+
+save(list=what,file=sysdataFilename,envir=nE)
+rm(nE)
+}
+
+GLIntegrate <- function(f, lower, upper, order = 500, ...){
+ if(order %in% c(50, 100, 400, 500, 800, 1000, 4000, 5000, 8000, 10000,
+ 40000, 50000, 80000, 100000))
+ AW <- getFromNamespace(paste(".AW", as.character(order),
+ sep = "."), ns = "distrEx")
+ else
+ AW <- .GLaw(order)
+
+ # transformation to [lower, upper]
+ xl <- (upper - lower)/2
+ W <- xl*AW[,2]
+ A <- xl*AW[,1] + (lower + upper)/2
+
+ res <- W*c(f(A, ...))
+ sum(res)
+}
+
+distrExIntegrate <- function(f, lower, upper, subdivisions = 100,
+ rel.tol = .Machine$double.eps^0.25,
+ abs.tol = rel.tol, stop.on.error = TRUE,
+ distr, order = .distrExOptions$GLIntegrateOrder,
+ ..., diagnostic = FALSE){
+ mc <- match.call()
+ time <- proc.time()
+
+ ## taken from base::system.time
+ ppt <- function(y) {
+ if (!is.na(y[4L]))
+ y[1L] <- y[1L] + y[4L]
+ if (!is.na(y[5L]))
+ y[2L] <- y[2L] + y[5L]
+ paste(formatC(y[1L:3L]), collapse = " ")
+ }
+
+ dots <- list(...)
+ dotsFun <- .filterFunargs(dots,f)
+ funwD <- function(x) do.call(f,c(list(x), dotsFun))
+
+ on.exit(message("Timing stopped at: ", ppt(proc.time() -
+ time)))
+ dotsInt <- if(length(names(dots))) dots[names(dots)%in% names(formals(integrate))] else NULL
+ res <- try(do.call(integrate, c(list(funwD, lower = lower, upper = upper, rel.tol = rel.tol,
+ abs.tol = abs.tol, stop.on.error = stop.on.error), dotsInt)),
+ silent = TRUE)
+
+ # if integrate fails => Gauss-Legendre integration
+ if(!is(res,"try-error")){
+ val <- res$value
+ if(diagnostic){
+ diagn <- list(call = mc, method = "integrate",
+ args = c(list(lower=lower, upper = upper, rel.tol = rel.tol,
+ abs.tol = abs.tol,
+ stop.on.error = stop.on.error),list(...)),
+ result = res)
+ res <- val
+ }else res <- val
+ }else{
+ errmess <- res
+ Zi <- 1
+ if(lower >= upper){
+ lo <- lower
+ lower <- upper
+ upper <- lo
+ Zi <- -1
+ }
+ if(!missing(distr)){
+ q.lDots <- NULL
+ if(length(names(dots))) {
+ q.lDots <- dots[names(dots) %in% names(formals(q.l(distr)))]
+ q.lDots[["p"]] <- q.lDots[["lower.tail"]] <- NULL
+ }
+ }
+
+ if(!is.finite(lower))
+ if(missing(distr)) stop(res)
+ else{
+ lower <- do.call(q.l(distr),
+ c(list(.distrExOptions$GLIntegrateTruncQuantile),q.lDots))
+# value <- c(...)
+# if(is(distr, "UnivariateCondDistribution"))
+# lower <- q.l(distr)(.distrExOptions$GLIntegrateTruncQuantile,
+# cond = value$cond)
+# else
+# lower <- q.l(distr)(.distrExOptions$GLIntegrateTruncQuantile)
+ }
+ if(!is.finite(upper))
+ if(missing(distr)) stop(res)
+ else{
+ q.lArgs <- if("lower.tail" %in% names(formals(distr at q)))
+ list(p=.distrExOptions$GLIntegrateTruncQuantile, lower.tail=FALSE) else
+ list(p=1-.distrExOptions$GLIntegrateTruncQuantile)
+ q.lArgs <- c(q.lArgs, q.lDots)
+ upper <- do.call(q.l(distr),q.lArgs)
+# value <- c(...)
+# if(is(distr, "UnivariateCondDistribution")){
+# if("lower.tail" %in% names(formals(distr at q)))
+# upper <- q.l(distr)(.distrExOptions$GLIntegrateTruncQuantile,
+# cond = value$cond, lower.tail = FALSE)
+# else
+# upper <- q.l(distr)(1 - .distrExOptions$GLIntegrateTruncQuantile,
+# cond = value$cond)
+# }else{
+# if("lower.tail" %in% names(formals(distr at q)))
+# upper <- q.l(distr)(.distrExOptions$GLIntegrateTruncQuantile,
+# lower.tail = FALSE)
+# else
+# upper <- q.l(distr)(1 - .distrExOptions$GLIntegrateTruncQuantile)
+# }
+ }
+ dotsGLInt <- NULL
+ if(length(names(dots))) dotsGLInt <- dots[names(dots)%in% names(formals(GLIntegrate))]
+ res <- Zi* do.call(GLIntegrate,c(list(f = funwD, lower = lower, upper = upper,
+ order = order),dotsGLInt))
+ if(diagnostic){
+ diagn <- list(call = mc, method = "GLIntegrate",
+ args = c(list(lower=lower, upper=upper, order=order),
+ list(...)),
+ result = list(GLIresult = res, errorMessage = errmess),
+ distrExOptions = .distrExOptions)
+ }
+ }
+
+ new.time <- proc.time()
+ on.exit()
+ if(diagnostic){
+ diagn$time <- structure(new.time - time, class = "proc_time")
+ attr(res,"diagnostic") <- diagn
+ class(attr(res,"diagnostic"))<- "DiagnosticClass"
+ }
+ return(res)
+}
Modified: pkg/distrEx/inst/NEWS
===================================================================
--- pkg/distrEx/inst/NEWS 2023-09-19 15:57:19 UTC (rev 1416)
+++ pkg/distrEx/inst/NEWS 2023-11-27 06:08:26 UTC (rev 1417)
@@ -23,6 +23,8 @@
+ KolmogorovDist methods gain argument "..." in order to enable them to
digest argument validity.check in MCEstimator
+ fixed some broken URLs and changed URLs from http to https where possible
++ changed latin1 to UTF-8 encoding as requested by CRAN since latin1
+ encoding will be deprecated.
##############
v 2.8
More information about the Distr-commits
mailing list