[Distr-commits] r729 - branches/distr-2.4/pkg/distr/R branches/distr-2.4/pkg/distr/man branches/distr-2.4/pkg/distrEx/src pkg/distr pkg/distr/R pkg/distr/man
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Thu Sep 1 03:35:32 CEST 2011
Author: ruckdeschel
Date: 2011-09-01 03:35:31 +0200 (Thu, 01 Sep 2011)
New Revision: 729
Removed:
branches/distr-2.4/pkg/distrEx/src/distrEx.dll
Modified:
branches/distr-2.4/pkg/distr/R/AllClasses.R
branches/distr-2.4/pkg/distr/R/AllInitialize.R
branches/distr-2.4/pkg/distr/R/DiscreteDistribution.R
branches/distr-2.4/pkg/distr/R/LatticeDistribution.R
branches/distr-2.4/pkg/distr/R/internalUtils.R
branches/distr-2.4/pkg/distr/man/internals.Rd
pkg/distr/DESCRIPTION
pkg/distr/R/AllClasses.R
pkg/distr/R/AllInitialize.R
pkg/distr/R/DiscreteDistribution.R
pkg/distr/R/LatticeDistribution.R
pkg/distr/R/internalUtils.R
pkg/distr/man/internals.Rd
Log:
[pkg distr]
+ revised initialize method convolution method and generating function for LatticeDistribution -> new routine to determine the smallest common grid for convolution
+ commented out require command in .onAttach
+ corrected error in definition of Beta() --- d slot was wrong
+ corrected type (inifinite...)
Modified: branches/distr-2.4/pkg/distr/R/AllClasses.R
===================================================================
--- branches/distr-2.4/pkg/distr/R/AllClasses.R 2011-07-08 13:13:36 UTC (rev 728)
+++ branches/distr-2.4/pkg/distr/R/AllClasses.R 2011-09-01 01:35:31 UTC (rev 729)
@@ -1,8 +1,9 @@
############# preparations ################
# (.onload, .onattach ...)
-.onLoad <- function(lib, pkg) { # extended 03-28-06: P.R.
- require("methods", character = TRUE, quietly = TRUE)
-}
+### left out P.R. 01-09-11
+##.onLoad <- function(lib, pkg) { # extended 03-28-06: P.R.
+## require("methods", character = TRUE, quietly = TRUE)
+##}
@@ -722,7 +723,7 @@
},
d = function(x, log = FALSE){
dbeta(x, shape1 = 1, shape2 = 1, ncp = 0,
- lower.tail = lower.tail, log.p = log.p)
+ log = log)
},
p = function(q, lower.tail = TRUE, log.p = FALSE ){
pbeta(q, shape1 = 1, shape2 = 1, ncp = 0,
Modified: branches/distr-2.4/pkg/distr/R/AllInitialize.R
===================================================================
--- branches/distr-2.4/pkg/distr/R/AllInitialize.R 2011-07-08 13:13:36 UTC (rev 728)
+++ branches/distr-2.4/pkg/distr/R/AllInitialize.R 2011-09-01 01:35:31 UTC (rev 729)
@@ -279,13 +279,14 @@
OS <- D at support
- if(is.null(lattice))
- { if(! .is.vector.lattice(OS))
- stop("Support as given/generated is not a lattice.")
- .Object at lattice <- .make.lattice.es.vector(OS)
- }else{
- .Object at lattice <- lattice
- }
+ #if(is.null(lattice))
+ # { if(! .is.vector.lattice(OS))
+ # stop("Support as given/generated is not a lattice.")
+ # .Object at lattice <- .make.lattice.es.vector(OS)
+ #}else{
+ .Object at lattice <- if(is.null(lattice ))
+ new("Lattice") else lattice
+ #}
.Object at support <- OS
Modified: branches/distr-2.4/pkg/distr/R/DiscreteDistribution.R
===================================================================
--- branches/distr-2.4/pkg/distr/R/DiscreteDistribution.R 2011-07-08 13:13:36 UTC (rev 728)
+++ branches/distr-2.4/pkg/distr/R/DiscreteDistribution.R 2011-09-01 01:35:31 UTC (rev 729)
@@ -15,7 +15,7 @@
if(!is.numeric(supp))
stop("'supp' is no numeric vector")
if(any(!is.finite(supp))) # admit +/- Inf?
- stop("inifinite or missing values in supp")
+ stop("infinite or missing values in supp")
len <- length(supp)
if(missing(prob)){
prob <- rep(1/len, len)
@@ -23,7 +23,7 @@
if(len != length(prob))
stop("'supp' and 'prob' must have equal lengths")
if(any(!is.finite(prob)))
- stop("inifinite or missing values in prob")
+ stop("infinite or missing values in prob")
if(!identical(all.equal(sum(prob), 1,
tolerance = 2*getdistrOption("TruncQuantile")), TRUE))
stop("sum of 'prob' has to be (approximately) 1")
@@ -226,78 +226,15 @@
e1.L <- as(e1, "LatticeDistribution")
e2.L <- as(e2, "LatticeDistribution")
if(is(e1.L, "LatticeDistribution") & is(e2.L, "LatticeDistribution"))
- {w1 <- width(lattice(e1.L))
- w2 <- width(lattice(e2.L))
- W <- sort(abs(c(w1,w2)))
- if (abs(abs(w1)-abs(w2))<getdistrOption("DistrResolution") ||
- W[2] %% W[1] < getdistrOption("DistrResolution") )
- return(e1.L + e2.L)
- }
- convolutedsupport <- rep(support(e1), each = length(support(e2))) +
- support(e2)
+ {w1 <- width(lattice(e1.L))
+ w2 <- width(lattice(e2.L))
+ W <- sort(abs(c(w1,w2)))
+ if (abs(abs(w1)-abs(w2))<getdistrOption("DistrResolution") ||
+ W[2] %% W[1] < getdistrOption("DistrResolution") )
+ return(e1.L + e2.L)
+ }
+ .convDiscrDiscr(e1,e2)})
- gridvalues1 <- d(e1)(support(e1)); gridvalues2 <- d(e2)(support(e2))
- convolutedvalues <- rep(gridvalues1, each = length(support(e2))) *
- gridvalues2
- rm(gridvalues1,gridvalues2)
-
- tmptable <- data.frame(x = convolutedsupport, dx = convolutedvalues)
- rm(convolutedsupport,convolutedvalues)
- tmp <- tapply(tmptable$dx, tmptable$x, sum)
- rm(tmptable)
-
- supp.u <- as.numeric(names(tmp))
- prob.u <- as.numeric(tmp)
-
- o <- order(supp.u)
- supp <- supp.u[o]
- prob <- prob.u[o]
-
- #supp.u <- unique(supp)
-
- len <- length(supp)
-
- if(len > 1){
- if (min(diff(supp))< getdistrOption("DistrResolution")){
- if (getdistrOption("DistrCollapse")){
- erg <- .DistrCollapse(supp, prob,
- getdistrOption("DistrResolution"))
- if ( len > length(erg$prob) &&
- getdistrOption("DistrCollapse.Unique.Warn") )
- warning("collapsing to unique support values")
- prob <- erg$prob
- supp <- erg$supp
- }else
- stop("grid too narrow --> change DistrResolution")
- }
- }
-
- rm(tmp, len)
-
- .withSim <- e1 at .withSim || e2 at .withSim
-
- rfun <- function(n) {}
- body(rfun) <- substitute({ f(n) + g(n) },
- list(f = e1 at r, g = e2 at r))
-
- dfun <- .makeDNew(supp, prob, Cont = FALSE)
- pfun <- .makePNew(supp, prob, .withSim, Cont = FALSE)
- qfun <- .makeQNew(supp, cumsum(prob), rev(cumsum(rev(prob))),
- .withSim, min(supp), max(supp), Cont = FALSE)
-
- object <- new("DiscreteDistribution", r = rfun, d = dfun, p = pfun,
- q = qfun, support = supp,
- .withSim = .withSim, .withArith = TRUE)
- rm(rfun, dfun, qfun, pfun)
-
- if(is(e1 at Symmetry,"SphericalSymmetry")&&
- is(e2 at Symmetry,"SphericalSymmetry"))
- object at Symmetry <- SphericalSymmetry(SymmCenter(e1 at Symmetry)+
- SymmCenter(e2 at Symmetry))
-
- object
- })
-
setMethod("+", c("Dirac","DiscreteDistribution"),
function(e1,e2){e2+location(e1)})
@@ -569,3 +506,6 @@
.lowerExact = .lowerExact(object),
.logExact = .logExact(object)))}
)
+
+
+
\ No newline at end of file
Modified: branches/distr-2.4/pkg/distr/R/LatticeDistribution.R
===================================================================
--- branches/distr-2.4/pkg/distr/R/LatticeDistribution.R 2011-07-08 13:13:36 UTC (rev 728)
+++ branches/distr-2.4/pkg/distr/R/LatticeDistribution.R 2011-09-01 01:35:31 UTC (rev 729)
@@ -12,12 +12,12 @@
{ D <- DiscreteDistribution
if (is(lattice, "Lattice"))
{ ### check consistency with support of DiscreteDistribution}
- if( !.is.consistent(lattice, support(D), eq.space = FALSE)){
- if (check)
+ if (check){
+ if( !.is.consistent(lattice, support(D), eq.space = FALSE))
stop(paste("Argument 'lattice' is inconsistent to",
" the support of argument 'DiscreteDistribution'." ,
sep = ""))
- else return(D)}
+ }
return(new("AffLinLatticeDistribution", r = D at r, d = D at d,
q = D at q, p = D at p, support = D at support,
a = D at a, b = D at b, X0 = D at X0,
@@ -25,11 +25,11 @@
.withSim = .withSim, img = D at img,
param = D at param, Symmetry = Symmetry))
}else{
- if( !.is.vector.lattice(support(D))){
- if (check)
+ if (check){
+ if( !.is.vector.lattice(support(D)))
stop(paste("Support of argument 'DiscreteDistribution' ",
"is not a lattice.", sep = ""))
- else return(D)}
+ }
return(new("AffLinLatticeDistribution", r = D at r, d = D at d,
q = D at q, p = D at p, support = D at support,
lattice = .make.lattice.es.vector(D at support),
@@ -44,23 +44,23 @@
{ D <- DiscreteDistribution
if (is(lattice, "Lattice"))
{ ### check consistency with support of DiscreteDistribution}
- if( !.is.consistent(lattice, support(D), eq.space = FALSE)){
- if (check)
+ if (check){
+ if( !.is.consistent(lattice, support(D), eq.space = FALSE))
stop(paste("Argument 'lattice' is inconsistent to the",
" support of argument 'DiscreteDistribution'." ,
sep = ""))
- else return(D)}
+ }
return(new("LatticeDistribution", r = D at r, d = D at d,
q = D at q, p = D at p, support = D at support,
lattice = lattice, .withArith = .withArith,
.withSim = .withSim, img = D at img,
param = D at param, Symmetry = Symmetry))
}else{
- if( !.is.vector.lattice(support(D))){
- if (check)
+ if (check){
+ if( !.is.vector.lattice(support(D)))
stop(paste("Support of argument 'DiscreteDistribution' is",
"not a lattice.", sep = " "))
- else return(D)}
+ }
return(new("LatticeDistribution", r = D at r, d = D at d,
q = D at q, p = D at p, support = D at support,
@@ -76,10 +76,10 @@
.withArith = .withArith,
.withSim = .withSim, Symmetry = Symmetry )
- if( !.is.consistent(lattice, supp, eq.space = FALSE)){
- if (check)
+ if (check){
+ if( !.is.consistent(lattice, supp, eq.space = FALSE))
stop("Argument 'lattice' is inconsistent to argument 'supp'.")
- else return(D)}
+ }
return(new("LatticeDistribution", r = r(D), d = d(D),
q = q(D), p = p(D), support = supp,
@@ -103,9 +103,10 @@
lattice = lattice, .withArith = .withArith,
.withSim = .withSim, Symmetry = Symmetry))
}else{
- if (check)
+ #if (check)
stop("Lengths of lattice and probabilities differ.")
- else return(D)}
+ #else return(D)
+ }
}else {if (is.null(prob))
stop(paste("Insufficient information given to ",
"determine distribution.", sep = ""))
@@ -126,18 +127,16 @@
{if (is.null(prob)) prob <- supp*0+1/length(supp)
D <- DiscreteDistribution(supp, prob, .withArith = .withArith,
.withSim = .withSim, Symmetry = Symmetry )
- if (!.is.vector.lattice (supp)){
- if (check)
+ if (check){
+ if (!.is.vector.lattice (supp))
stop("Argument 'supp' given is not a lattice.")
- else return (D)
- }else{
- return(new("LatticeDistribution", r = D at r, d = D at d,
+ }
+ return(new("LatticeDistribution", r = D at r, d = D at d,
q = D at q, p = D at p, support = D at support,
lattice = .make.lattice.es.vector(D at support),
.withArith = D at .withArith,
.withSim = D at .withSim, img = D at img,
param = D at param, Symmetry = Symmetry))
- }
}else
stop("Insufficient information given to determine distribution.")
}
@@ -183,79 +182,100 @@
setMethod("+", c("LatticeDistribution", "LatticeDistribution"),
function(e1,e2){
- ### Step 1
+ if(length(support(e1))==1) return(e2+support(e1))
+ if(length(support(e2))==1) return(e1+support(e2))
-# e1 <- as(e1, "LatticeDistribution")
-# e2 <- as(e2, "LatticeDistribution")
-# ### casting necessary due to setIs
+### Lattice calculations:
- ### Lattice Calculations:
- w1 <- width(lattice(e1))
- w2 <- width(lattice(e2))
sup1 <- support(e1)
sup2 <- support(e2)
- maxl <- length(sup1)*length(sup2)
- ### length of product grid
- csup <- unique(sort(c(sup1,sup2)))
- ### grid width of convolution grid
+ # left and right endpoint of convolution support
+ su12.l <- sup1[1]+sup2[1]
+ su12.r <- (rev(sup1))[1]+(rev(sup2))[1]
- w <- min(diff(csup))
- commonsup <- unique(sort(c(outer(sup1,sup2,"+"))))
- ### grid width of convolution grid
- dcs <- abs(diff(commonsup))
- mw <- min(dcs[dcs>getdistrOption("DistrResolution")])
- if (abs(abs(w1)-abs(w2)) > getdistrOption("DistrResolution")){
- W <- sort(abs(c(w1,w2)))
- if (W[2] %% W[1] > getdistrOption("DistrResolution")){
+ l1 <- length(sup1)
+ l2 <- length(sup2)
- ## check whether arrangement on common grid really
- ## saves something
+ lat1 <- lattice(e1)
+ lat2 <- lattice(e2)
+ L1 <- Length(lat1)
+ L2 <- Length(lat2)
+ w1 <- width(lat1)
+ w2 <- width(lat2)
- prob1 <- d(e1)(sup1)
- prob2 <- d(e2)(sup2)
- ### convolutional grid
- comsup <- seq(min(commonsup),max(commonsup), by = mw)
- fct <- function(sup0, prob0, bw){
- ### expand original grid,prob onto new width:
- sup00 <- seq(min(sup0), max(sup0), by = mw)
- prb0 <- 0 * sup00
- ind0 <- sup00 %in% sup0
- prb0[ind0] <- prob0
- return(LatticeDistribution(supp = sup00,
- prob = prb0))
- }
- if(length(comsup) < maxl)
- return( fct(sup1,prob1,mw) + fct(sup2,prob2,mw))
- else
- return(as(e1, "DiscreteDistribution") +
- as(e2, "DiscreteDistribution"))
- }
- else
- w <- mw #generate common lattice / support
+ ### take care if lattice is infinite
+ L.inf <- !(is.finite(L1)&&is.finite(L2))
+ if(L.inf){
+ if(is.finite(L2)){
+ if(w1>0)
+ L.lr <- +1
+ else
+ L.lr <- -1
+ }else{
+ if(is.finite(L1)){
+ if(w2>0)
+ L.lr <- +1
+ else
+ L.lr <- -1
+ }else{
+ if(w1*w2>0) L.lr <- if(w1>0) +1 else -1
+ if(w1*w2<0) L.lr <- if(abs(su12.l)<abs(su12.r)) +1 else -1
}
+ }
+ }
- newlat <- NULL
- ### if any lattice is infinite: see if we can keep this in mind:
- if( ! (is.finite(Length(lattice(e1))) &&
- is.finite(Length(lattice(e2))) ) &&
- width(lattice(e1)) * width(lattice(e2)) > 0 )
- {p1 <- pivot(lattice(e1))
- p2 <- pivot(lattice(e2))
- p <- p1 + p2
- w0 <- if (width(lattice(e1))>0) w else -w
- newlat <- Lattice(pivot = p, width = w0, Length = Inf)
- }
- ### end Lattice Calculations
+ e0 <- NULL
+ tol0 <- .distroptions$DistrResolution/1000
+
+ ## treat case separately when Discr + Discr is "faster"
+ if(l1*l2 < 100){
+ d0 <- .convDiscrDiscr(e1,e2)
+ sup0 <- support(d0)
+ md <- min(diff(sup0))
+ sup00 <- seq(from=min(sup0),to=max(sup0),by=md)
+ sup0s <- intersect(sup00,sup0)
+ sup01 <- .inWithTol(sup00, sup0s, tol=tol0)
+ sup10 <- .inWithTol(sup0, sup0s, tol=tol0)
+ if(!all(sup10)) return(d0)
+ pr0 <- sup00*0
+ pr0[sup01] <- (prob(d0))[sup10]
+ pr0 <- pr0/sum(pr0)
+ lat <- Lattice(pivot = sup00[1], width = md,
+ Length = length(sup00))
+ e0 <- LatticeDistribution(supp = sup00, prob = pr0,
+ lattice = lat, check = FALSE)
+ if(L.inf){
+ e0 at lattice <- if(L.lr>0){
+ Lattice(pivot = su12.l, width = wa, Length = Inf)
+ }else{
+ Lattice(pivot = su12.r, width = -wa, Length = Inf)}
+ }
+ }
+
+ ## step 1 common width
+ wa <- .getCommonWidth(abs(w1),abs(w2),
+ tol=tol0)
+
+
+ ## treat case separately when no common support, i.e. when
+ ## w1/w2 is not "rational" enough
+
+ if(is.null(wa)) return(.convDiscrDiscr(e1,e2))
+
+
+ w0 <- ifelse(w1<0,-1,1) * wa
+ pi1 <- pivot(lat1)
+ pi2 <- pivot(lat2)
### Step 2
- supp0 <- seq(by = abs(w),
- from = min(support(e1), support(e2)),
- to = max(support(e1), support(e2)))
+ supp0 <- seq(by = wa, from = min(sup1-pi1, sup2-pi2), to = max(sup1-pi1, sup2-pi2))
+ s1 <- .inWithTol(supp0,sup1-pi1,tol0)
+ s2 <- .inWithTol(supp0,sup2-pi2,tol0)
+ d1 <- d2 <- 0*supp0
+ d1[s1] <- prob(as(e1,"DiscreteDistribution"))
+ d2[s2] <- prob(as(e2,"DiscreteDistribution"))
- d1 <- d(e1)(supp0)
- d2 <- d(e2)(supp0)
-
L <- length(supp0)
Ln <- 2^(ceiling(log(L)/log(2))+1)
@@ -270,13 +290,14 @@
## convolution theorem for DFTs
newd <- (Re(fft(ftde1*ftde2, inverse = TRUE)) / Ln)[1:(2*L+1)]
- newd <- (newd >= .Machine$double.eps)*newd
+ newd <- (newd >= .Machine$double.eps^1.5)*newd
## reduction to relevant support
- supp1 <- seq(by = abs(w),
- from = 2 * min(support(e1), support(e2)),
- to = 2 * max(support(e1), support(e2)))
+ supp1 <- seq(by = wa,
+ from = 2 * min(sup1-pi1, sup2-pi2),
+ to = 2 * max(sup1-pi1, sup2-pi2))+pi1+pi2
+
L1 <- length(supp1)
newd <- newd[1:L1]
@@ -287,9 +308,6 @@
)
rsum.l <- 1 + sum( cumsum(newd) <
getdistrOption("TruncQuantile")/2)
- newd <- newd[rsum.l:rsum.u]
- newd <- newd/sum(newd)
- supp1 <- supp1[rsum.l:rsum.u]
}else{
rsum.u <- min( sum( rev(cumsum(rev(newd))) >=
.Machine$double.eps),
@@ -297,12 +315,15 @@
)
rsum.l <- 1 + sum( cumsum(newd) < .Machine$double.eps)
- newd <- newd[rsum.l:rsum.u]
- newd <- newd/sum(newd)
- supp1 <- supp1[rsum.l:rsum.u]
}
- supp2 <- supp1[newd > getdistrOption("TruncQuantile")]
- newd2 <- newd[newd > getdistrOption("TruncQuantile")]
+ wi1 <- rsum.l:rsum.u
+ newd <- newd[wi1]
+ newd <- newd/sum(newd)
+ supp1 <- supp1[wi1]
+
+ wi2 <- newd > getdistrOption("TruncQuantile")
+ supp2 <- supp1[wi2]
+ newd2 <- newd[wi2]
newd2 <- newd2/sum(newd2)
Symmetry <- NoSymmetry()
@@ -311,13 +332,24 @@
Symmetry <- SphericalSymmetry(SymmCenter(e1 at Symmetry)+
SymmCenter(e2 at Symmetry))
-
- if( length(supp1) >= 2 * length(supp2))
+ if( length(supp1) >= 2 * length(supp2)){
return(DiscreteDistribution(supp = supp2, prob = newd2,
.withArith = TRUE, Symmetry = Symmetry))
- else
- return(LatticeDistribution(supp = supp1, prob = newd,
- .withArith = TRUE, Symmetry = Symmetry))
+ }else{
+ lat <- Lattice(pivot=supp1[1],width=wa, Length=length(supp1))
+
+ e0 <- LatticeDistribution(supp = supp1, prob = newd,
+ lattice = lat,
+ .withArith = TRUE, Symmetry = Symmetry,
+ check = FALSE)
+ if(L.inf){
+ e0 at lattice <- if(L.lr>0){
+ Lattice(pivot = su12.l, width = wa, Length = Inf)
+ }else{
+ Lattice(pivot = su12.r, width = -wa, Length = Inf)}
+ }
+ return(e0)
+ }
})
## extra methods
@@ -336,7 +368,8 @@
Symmetry <- SphericalSymmetry(SymmCenter(e1 at Symmetry)+e2)
LatticeDistribution(lattice = L,
- DiscreteDistribution = Distr, Symmetry = Symmetry)
+ DiscreteDistribution = Distr, Symmetry = Symmetry,
+ check = FALSE)
})
setMethod("*", c("LatticeDistribution", "numeric"),
@@ -356,7 +389,8 @@
Symmetry <- SphericalSymmetry(SymmCenter(e1 at Symmetry) * e2)
return(LatticeDistribution(lattice = L,
- DiscreteDistribution = Distr, Symmetry = Symmetry))
+ DiscreteDistribution = Distr, Symmetry = Symmetry,
+ check = FALSE))
}
}
)
@@ -371,7 +405,8 @@
LatticeDistribution(lattice = L,
DiscreteDistribution =
as(e1, "AffLinDiscreteDistribution") + e2,
- Symmetry = Symmetry)
+ Symmetry = Symmetry,
+ check = FALSE)
})
setMethod("*", c("AffLinLatticeDistribution", "numeric"),
@@ -388,7 +423,8 @@
return(LatticeDistribution(lattice = L,
DiscreteDistribution =
as(e1, "AffLinDiscreteDistribution") *
- e2, Symmetry = Symmetry))
+ e2, Symmetry = Symmetry,
+ check = FALSE))
}
}
)
Modified: branches/distr-2.4/pkg/distr/R/internalUtils.R
===================================================================
--- branches/distr-2.4/pkg/distr/R/internalUtils.R 2011-07-08 13:13:36 UTC (rev 728)
+++ branches/distr-2.4/pkg/distr/R/internalUtils.R 2011-09-01 01:35:31 UTC (rev 729)
@@ -48,6 +48,114 @@
}
#------------------------------------------------------------------------------
+### internal help function to determin common lattice width
+#------------------------------------------------------------------------------
+
+.EuclidAlgo <- function(n1,n2){
+ r<- 2
+ m <- round(max(n1,n2))
+ n <- round(min(n1,n2))
+ if (n==1) return(1)
+ while (r>1){
+ r <- m %% n
+ m <- n
+ if(r==1) return(1)
+ if(r==0) break
+ n <- r
+ }
+ return(n)
+}
+.getCommonWidth <- function(x1,x2, tol=.Machine$double.eps){
+ rI <- function(x) .isInteger(x, tol = x * tol)
+ if(rI(x1)&&rI(x2)) return(.EuclidAlgo(x1,x2))
+ if(rI(10000*x1) && rI(10000*x2)) return(.EuclidAlgo(10000*x1,10000*x2)/10000)
+ n <- max(x1,x2); m <- min(x1,x2)
+ if(rI(n/m)) return(m)
+ vc <- 1:10000
+ vecT <- sapply(vc,function(x) rI(x*(n/m)))
+ vec <- m/vc
+# print(min(vc[vecT]))
+ if(any(vecT)) return((vec/min(vc[vecT]))[1])
+ return(NULL)
+}
+
+#------------------------------------------------------------------------------
+### %in% for numerics with tolerance
+#------------------------------------------------------------------------------
+.inWithTol <- function(x,y,tol=.Machine$double.eps){
+ sapply(x, function(u) any(abs(u-y)<tol))
+}
+
+#------------------------------------------------------------------------------
+### brute force convolution
+#------------------------------------------------------------------------------
+.convDiscrDiscr <- function(e1,e2){
+ convolutedsupport <- rep(support(e1), each = length(support(e2))) +
+ support(e2)
+
+ gridvalues1 <- d(e1)(support(e1)); gridvalues2 <- d(e2)(support(e2))
+ convolutedvalues <- rep(gridvalues1, each = length(support(e2))) *
+ gridvalues2
+ rm(gridvalues1,gridvalues2)
+
+ tmptable <- data.frame(x = convolutedsupport, dx = convolutedvalues)
+ rm(convolutedsupport,convolutedvalues)
+ tmp <- tapply(tmptable$dx, tmptable$x, sum)
+ rm(tmptable)
+
+ supp.u <- as.numeric(names(tmp))
+ prob.u <- as.numeric(tmp)
+
+ o <- order(supp.u)
+ supp <- supp.u[o]
+ prob <- prob.u[o]
+
+ #supp.u <- unique(supp)
+
+ len <- length(supp)
+
+ if(len > 1){
+ if (min(diff(supp))< getdistrOption("DistrResolution")){
+ if (getdistrOption("DistrCollapse")){
+ erg <- .DistrCollapse(supp, prob,
+ getdistrOption("DistrResolution"))
+ if ( len > length(erg$prob) &&
+ getdistrOption("DistrCollapse.Unique.Warn") )
+ warning("collapsing to unique support values")
+ prob <- erg$prob
+ supp <- erg$supp
+ }else
+ stop("grid too narrow --> change DistrResolution")
+ }
+ }
+
+ rm(tmp, len)
+
+ .withSim <- e1 at .withSim || e2 at .withSim
+
+ rfun <- function(n) {}
+ body(rfun) <- substitute({ f(n) + g(n) },
+ list(f = e1 at r, g = e2 at r))
+
+ dfun <- .makeDNew(supp, prob, Cont = FALSE)
+ pfun <- .makePNew(supp, prob, .withSim, Cont = FALSE)
+ qfun <- .makeQNew(supp, cumsum(prob), rev(cumsum(rev(prob))),
+ .withSim, min(supp), max(supp), Cont = FALSE)
+
+ object <- new("DiscreteDistribution", r = rfun, d = dfun, p = pfun,
+ q = qfun, support = supp,
+ .withSim = .withSim, .withArith = TRUE)
+ rm(rfun, dfun, qfun, pfun)
+
+ if(is(e1 at Symmetry,"SphericalSymmetry")&&
+ is(e2 at Symmetry,"SphericalSymmetry"))
+ object at Symmetry <- SphericalSymmetry(SymmCenter(e1 at Symmetry)+
+ SymmCenter(e2 at Symmetry))
+
+ object
+ }
+
+#------------------------------------------------------------------------------
### .fm, .fM, .fM2 functions
#------------------------------------------------------------------------------
Modified: branches/distr-2.4/pkg/distr/man/internals.Rd
===================================================================
--- branches/distr-2.4/pkg/distr/man/internals.Rd 2011-07-08 13:13:36 UTC (rev 728)
+++ branches/distr-2.4/pkg/distr/man/internals.Rd 2011-09-01 01:35:31 UTC (rev 729)
@@ -17,6 +17,8 @@
\alias{.make.lattice.es.vector}
\alias{.presubs}
\alias{.inArgs}
+\alias{.EuclidAlgo}
+\alias{.getCommonWidth}
\alias{.isEqual}
\alias{.isIn}
\alias{.setEqual}
@@ -55,6 +57,8 @@
\alias{.trunc.low}
\alias{.modifyqgaps}
\alias{.DistrCollapse}
+\alias{.convDiscrDiscr}
+\alias{.inWithTol}
\alias{devNew}
\title{Internal functions of package distr}
@@ -125,6 +129,10 @@
.trunc.low(object, lower)
.modifyqgaps(pfun, qfun, gaps, leftright = "left")
.DistrCollapse(support, prob, eps = getdistrOption("DistrResolution"))
+.EuclidAlgo(n1,n2)
+.getCommonWidth(x1,x2, tol=.Machine$double.eps)
+.convDiscrDiscr(e1,e2)
+.inWithTol(x,y,tol=.Machine$double.eps)
devNew(...)
}
@@ -228,6 +236,10 @@
the left continuous version; for slot \code{p}: if partially
matched to \code{"left"} the left continuous version, else
the right continuous version;}
+ \item{n1}{integer argument for \code{.EuclidAlgo}}
+ \item{n2}{integer argument for \code{.EuclidAlgo}}
+ \item{x1}{width argument for \code{.getCommonWidth}}
+ \item{x2}{width argument for \code{.getCommonWidth}}
\item{...}{arguments passed through to other functions}
}
@@ -366,7 +378,14 @@
\code{gaps} is not \code{NULL}. If argument \code{leftright} does not
partially match \code{"right"} (default) returns the left continuous
version of the quantile function, else the right continuous one.
-
+
+\code{.EuclidAlgo} computes the greatest common divisor of two integers by
+ means of the Euclidean algorithm.
+\code{.getCommonWidth} for two lattices with widths \code{x1} and \code{x2}
+ computes the smallest common lattice width for convolution.
+\code{.convDiscrDiscr} computes the convolution of two discrete distributions by
+ brute force.
+\code{.inWithTol} works like \code{\%in\%} but with a given tolerance.
\code{devNew} opens a new device. This function is for back compatibility
with R versions < 2.8.0.
}
@@ -433,6 +452,11 @@
Criterium for collapsing: a distance smaller than argument
\code{eps}.
}
+\item{.EuclidAlgo}{returns the greatest common divisor (an integer).}
+\item{.getCommonWidth}{returns the smallest common lattice width (a numeric).}
+\item{.convDiscrDiscr}{returns the convolution of two discrete distributions.}
+\item{.inWithTol}{returns a logical vector of same lenght as \code{x} for the
+ matches (up to tolerance) with vector \code{y}.}
\item{devNew}{returns the return value of the device opened,
usually invisible \code{NULL}.}
}
Deleted: branches/distr-2.4/pkg/distrEx/src/distrEx.dll
===================================================================
(Binary files differ)
Modified: pkg/distr/DESCRIPTION
===================================================================
--- pkg/distr/DESCRIPTION 2011-07-08 13:13:36 UTC (rev 728)
+++ pkg/distr/DESCRIPTION 2011-09-01 01:35:31 UTC (rev 729)
@@ -1,6 +1,6 @@
Package: distr
-Version: 2.3.2
-Date: 2011-03-22
+Version: 2.3.3
+Date: 2011-08-31
Title: Object oriented implementation of distributions
Description: S4 Classes and Methods for distributions
Author: Florian Camphausen, Matthias Kohl, Peter Ruckdeschel, Thomas Stabla
Modified: pkg/distr/R/AllClasses.R
[TRUNCATED]
To get the complete diff run:
svnlook diff /svnroot/distr -r 729
More information about the Distr-commits
mailing list