[Vegan-commits] r274 - in pkg: . R inst man
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Thu Mar 20 17:33:41 CET 2008
Author: gsimpson
Date: 2008-03-20 17:33:41 +0100 (Thu, 20 Mar 2008)
New Revision: 274
Added:
pkg/R/allPerms.R
pkg/R/print.allPerms.R
pkg/R/print.summary.allPerms.R
pkg/R/summary.allPerms.R
Modified:
pkg/DESCRIPTION
pkg/inst/ChangeLog
pkg/man/permCheck.Rd
Log:
New function allPerms to return a complete enumeration of all possible permutations for a given scheme and number of observations
Modified: pkg/DESCRIPTION
===================================================================
--- pkg/DESCRIPTION 2008-03-20 16:27:59 UTC (rev 273)
+++ pkg/DESCRIPTION 2008-03-20 16:33:41 UTC (rev 274)
@@ -1,7 +1,7 @@
Package: vegan
Title: Community Ecology Package
-Version: 1.12-4
-Date: Mar 9, 2008
+Version: 1.12-5
+Date: Mar 20, 2008
Author: Jari Oksanen, Roeland Kindt, Pierre Legendre, Bob O'Hara, Gavin L. Simpson,
M. Henry H. Stevens
Maintainer: Jari Oksanen <jari.oksanen at oulu.fi>
Added: pkg/R/allPerms.R
===================================================================
--- pkg/R/allPerms.R (rev 0)
+++ pkg/R/allPerms.R 2008-03-20 16:33:41 UTC (rev 274)
@@ -0,0 +1,172 @@
+`allPerms` <- function(n, control = permControl(), max = 9999,
+ observed = FALSE) {
+ ## in-line functions
+ `all.free` <- function(n, v = 1:n) {
+ if( n == 1 ) {
+ matrix(v, 1, 1)
+ } else {
+ X <- NULL
+ for(i in 1:n)
+ X <- rbind(X, cbind(v[i],
+ Recall(n-1, v[-i])))
+ X
+ }
+ }
+ `all.series` <- function(n, control) {
+ v <- seq_len(n)
+ nperms <- numPerms(v, control = control)
+ X <- matrix(nrow = nperms, ncol = n)
+ for(i in v) {
+ X[i,] <- seq(i, length = n)%%n + 1
+ }
+ ## if mirroring, rev the cols of X[v,]
+ if(control$mirror)
+ X[(n+1):(2*n),] <- X[v, rev(v)]
+ X
+ }
+ `all.grid` <- function(n, control) {
+ v <- seq_len(n)
+ nperms <- numPerms(v, control)
+ nr <- control$nrow
+ nc <- control$ncol
+ X <- matrix(nrow = nperms, ncol = n)
+ idx <- 1
+ ## ncol == 2 is special case
+ if(control$ncol == 2) {
+ X <- all.series(n, permControl(type = "series",
+ mirror = control$mirror,
+ constant = control$constant)
+ )
+ } else {
+ for(i in seq_len(nr)) {
+ for(j in seq_len(nc)) {
+ ir <- seq(i, length = nr)%%nr
+ ic <- seq(j, length = nc)%%nc
+ ## block 1 - no reversals
+ X[idx, ] <- rep(ic, each = nr) * nr +
+ rep(ir, len = nr * nc) + 1
+ if(control$mirror) {
+ ## block 2 - rev rows but not columns
+ X[idx + n, ] <- rep(ic, each = nr) * nr +
+ rep(rev(ir), len = nr * nc) + 1
+ ## block 3 - rev columns but not rows
+ X[idx + (2*n), ] <- rep(rev(ic), each = nr) *
+ nr + rep(ir, len = nr * nc) + 1
+ }
+ idx <- idx + 1
+ }
+ }
+ if(control$mirror) {
+ ## rev columns and rows
+ ## no calculations, just rev cols of block 1
+ v <- seq_len(n)
+ X[((3*n)+1):(4*n), ] <- X[v, rev(v)]
+ }
+ }
+ X
+ }
+ `all.strata` <- function(n, control) {#, nperms) {
+ v <- seq_len(n)
+ nperms <- numPerms(v, control)
+ lev <- length(levels(control$strata))
+ X <- matrix(nrow = nperms, ncol = length(control$strata))
+ perms <- all.free(lev)
+ sp <- split(v, control$strata)
+ for(i in seq_len(nrow(perms)))
+ X[i,] <- unname(do.call(c, sp[perms[i,]]))
+ X
+ }
+ ## recursive fun for perms within strata
+ bar <- function(mat, n) {
+ if(n == 1)
+ mat
+ else
+ mat <- rbind(mat, Recall(mat, n-1))
+ mat
+ }
+ ## start
+ v <- seq_len(n)
+ ## check permutation scheme and update control
+ pcheck <- permCheck(v, control = control)
+ control <- pcheck$control
+ ## get max number of permutations
+ ## originally had:
+ ##nperms <- numPerms(v, control = control)
+ ## but pcheck contains 'n', the result of call to numPerms
+ nperms <- pcheck$n
+ ## sanity check - don't let this run away to infinity
+ ## esp with type = "free"
+ if(nperms > max)
+ stop("Number of possible permutations too big (> 'max')")
+ type <- control$type
+ if(type != "strata" && !is.null(control$strata)) {
+ ## permuting within blocks
+ if(control$constant) {
+ ## same permutation in each block
+ v <- seq_len(n)
+ pg <- unique(table(control$strata))
+ control.wi <- permControl(type = control$type,
+ mirror = control$mirror,
+ nrow = control$nrow,
+ ncol = control$ncol)
+ nperm <- numPerms(v, control)
+ ord <- switch(control$type,
+ free = all.free(pg),
+ series = all.series(pg, control = control.wi),
+ grid = all.grid(pg, control = control.wi))
+ perm.wi <- nrow(ord)
+ sp <- split(v, control$strata)
+ res <- matrix(nrow = nperm, ncol = n)
+ for(i in seq_len(perm.wi))
+ res[i,] <- sapply(sp, function(x) x[ord[i,]])
+ } else {
+ ## different permutations within blocks
+ ng <- length(levels(control$strata))
+ pg <- length(control$strata) / ng
+ control.wi <- permControl(type = control$type,
+ mirror = control$mirror,
+ nrow = control$nrow,
+ ncol = control$ncol)
+ ord <- switch(control$type,
+ free = all.free(pg),
+ series = all.series(pg, control = control.wi),
+ grid = all.grid(pg, control = control.wi)
+ )
+ perm.wi <- nrow(ord)
+ add <- seq(from = 0, by = pg, length.out = ng)
+ res <- vector(mode = "list", length = ng)
+ a <- 1
+ b <- nperms / perm.wi
+ for(i in seq_len(ng)) {
+ res[[i]] <- matrix(rep(bar(ord+add[i], a), each = b),
+ ncol = pg)
+ a <- a*perm.wi
+ b <- b/perm.wi
+ }
+ res <- do.call(cbind, res)
+ }
+ } else {
+ ## no blocks
+ res <- switch(type,
+ free = all.free(n),
+ series = all.series(n, control=control),
+ grid = all.grid(n, control=control),
+ strata = all.strata(n, control=control)
+ )
+ }
+ ## some times storage.mode of res is numeric, sometimes
+ ## it is integer, set to "integer" for comparisons using
+ ## identical to match the observed ordering
+ storage.mode(res) <- "integer"
+ if(!observed) {
+ obs.row <- apply(res, 1, function(x, v) {identical(x, v)}, v)
+ res <- res[!obs.row, ]
+ ## reduce the number of permutations to get rid of the
+ ## observed ordering
+ control$nperm <- control$nperm - 1
+ }
+ class(res) <- "allPerms"
+ attr(res, "control") <- control
+ attr(res, "observed") <- observed
+ return(res)
+}
Added: pkg/R/print.allPerms.R
===================================================================
--- pkg/R/print.allPerms.R (rev 0)
+++ pkg/R/print.allPerms.R 2008-03-20 16:33:41 UTC (rev 274)
@@ -0,0 +1,9 @@
+`print.allPerms` <- function(x, ...) {
+ dims <- dim(x)
+ control <- attr(x, "control")
+ observed <- attr(x, "observed")
+ attributes(x) <- NULL
+ dim(x) <- dims
+ print(x)
+ return(invisible(x))
+}
Added: pkg/R/print.summary.allPerms.R
===================================================================
--- pkg/R/print.summary.allPerms.R (rev 0)
+++ pkg/R/print.summary.allPerms.R 2008-03-20 16:33:41 UTC (rev 274)
@@ -0,0 +1,17 @@
+`print.summary.allPerms` <- function(x, ...) {
+ dims <- dim(x)
+ control <- attr(x, "control")
+ observed <- attr(x, "observed")
+ attributes(x) <- NULL
+ dim(x) <- dims
+ cat("\n")
+ writeLines(strwrap("Complete enumeration of permutations\n",
+ prefix = "\t"))
+ cat("\nPermutation Scheme:\n")
+ print(control)
+ cat(paste("Contains observed ordering?:", ifelse(observed, "Yes", "No"),
+ "\n"))
+ cat("\nAll permutations:\n")
+ print(x)
+ return(invisible(x))
+}
Added: pkg/R/summary.allPerms.R
===================================================================
--- pkg/R/summary.allPerms.R (rev 0)
+++ pkg/R/summary.allPerms.R 2008-03-20 16:33:41 UTC (rev 274)
@@ -0,0 +1,4 @@
+`summary.allPerms` <- function(object, ...) {
+ class(object) <- "summary.allPerms"
+ object
+}
Modified: pkg/inst/ChangeLog
===================================================================
--- pkg/inst/ChangeLog 2008-03-20 16:27:59 UTC (rev 273)
+++ pkg/inst/ChangeLog 2008-03-20 16:33:41 UTC (rev 274)
@@ -2,8 +2,14 @@
VEGAN DEVEL VERSIONS at http://r-forge.r-project.org/
-Version 1.12-4 (opened Mar 9, 2008)
+Version 1.12-5 (opened Mar 20, 2008)
+
+ * allPerms: New function to enumerate all possible permutations
+ for a given permutation scheme and number of samples. Has 'print'
+ and 'summary' methods
+Version 1.12-4 (closed Mar 20, 2008)
+
* permDisper: Added pairwise comparisons of group dispersions
via a classical t test and via permutation test, accessed via
new argument 'pairwise = TRUE'. 'permDisper' is now documented
Modified: pkg/man/permCheck.Rd
===================================================================
--- pkg/man/permCheck.Rd 2008-03-20 16:27:59 UTC (rev 273)
+++ pkg/man/permCheck.Rd 2008-03-20 16:33:41 UTC (rev 274)
@@ -8,14 +8,19 @@
\alias{getNumObs.default}
\alias{getNumObs.integer}
\alias{getNumObs.numeric}
+\alias{allPerms}
+\alias{print.allPerms}
+\alias{summary.allPerms}
+\alias{print.summary.allPerms}
-\title{Number of possible permutations and checking of permutation schemes}
+\title{Utility functions for permutation schemes}
\description{
\code{permCheck} provides checking of permutation schemes for
validity. \code{numPerms} calculates the maximum number of
permutations possible under the current permutation
- scheme. \code{getNumObs} is a utility function to return the number of
- observations for a range of R and ordination objects.
+ scheme. \code{allPerms} enumerates all possible permutations for the
+ given scheme. \code{getNumObs} is a utility function to return the
+ number of observations for a range of R and ordination objects.
}
\usage{
permCheck(object, control = permControl())
@@ -24,6 +29,11 @@
numPerms(object, control = permControl())
+allPerms(n, control = permControl(), max = 9999,
+ observed = FALSE)
+
+\method{summary}{allPerms}(object, \dots)
+
getNumObs(object, \dots)
\method{getNumObs}{default}(object, \dots)
@@ -41,11 +51,17 @@
\item{control}{a list of control values describing properties of the
permutation design, as returned by a call to
\code{\link{permControl}}.}
+ \item{n}{the number of observations.}
+ \item{max}{the maximum number of permutations, below which complete
+ enumeration will be attempted. See Details.}
+ \item{observed}{logical, should the observed ordering of samples be
+ returned as part of the complete enumeration? Default is
+ \code{FALSE} to facilitate usage in higher level functions.}
\item{\dots}{arguments to other methods.}
}
\details{
- \code{permCheck} and \code{numPerms} are utility functions for working
- with the new permutation schemes available in
+ \code{permCheck}, \code{allPerms} and \code{numPerms} are utility
+ functions for working with the new permutation schemes available in
\code{\link{permuted.index2}}.
\code{permCheck} is used to check the current permutation schemes
@@ -66,6 +82,23 @@
than \code{control$minperm}, it is better to enumerate all possible
permutations, and as such complete enumeration of all permutations is
turned on (\code{control$complete} is set to \code{TRUE}).
+
+ Function \code{allPerms} enumerates all possible permutations for the
+ number of observations and the selected permutation scheme. It has
+ \code{\link{print}} and \code{\link{summary}} methods. \code{allPerms}
+ returns a matrix containing all possible permutations, possibly
+ containing the observed ordering (if argument \code{observed} is
+ \code{TRUE}). The rows of this matrix are the various permutations and
+ the columns reflect the number of samples.
+
+ With free permutation designs, and restricted permutation schemes with
+ large numbers of observations, there are a potentially huge number of
+ possible permutations of the samples. It would be inefficient, not to
+ mention incredibly time consuming, to enumerate them all. Storing all
+ possible permutations would also become problematic in such cases. To
+ control this and guard against trying to evaluate too large a number
+ of permutations, if the number of possible permutations is larger than
+ \code{max}, \code{allPerms} exits with an error.
}
\value{
For \code{permCheck} a list containing the maximum number of
@@ -79,17 +112,17 @@
}
%\references{
%}
-\note{
- Currently, complete enumeration is not implemented in
- \code{\link{permuted.index2}}. \code{permCheck} sets complete
- enumeration in the \code{control} component of the returned object
- but this is not acted upon in the permutations.
-
- \code{permCheck} does reduce the \code{nperm} and \code{maxperm}
- components of the returned object to the maximum number of possible
- permutations, and \code{\link{permuted.index2}} obeys these
- directives.
-}
+%\note{
+% Currently, complete enumeration is not implemented in
+% \code{\link{permuted.index2}}. \code{permCheck} sets complete
+% enumeration in the \code{control} component of the returned object
+% but this is not acted upon in the permutations.
+%
+% \code{permCheck} does reduce the \code{nperm} and \code{maxperm}
+% components of the returned object to the maximum number of possible
+% permutations, and \code{\link{permuted.index2}} obeys these
+% directives.
+%}
\author{Gavin Simpson}
\seealso{\code{\link{permuted.index2}} and \code{\link{permControl}}.}
@@ -144,6 +177,22 @@
vec4 <- permCheck(1:100, permControl(type = "series", mirror = TRUE))
all.equal(vec4$n, 200)
+## enumerate all possible permutations
+fac <- gl(3,6)
+ctrl <- permControl(type = "grid", mirror = FALSE, strata = fac,
+ constant = TRUE, nrow = 3, ncol = 2)
+numPerms(1:18, control = ctrl)
+(tmp <- allPerms(18, control = ctrl, observed = TRUE))
+(tmp2 <- allPerms(18, control = ctrl))
+## turn on mirroring
+ctrl$mirror <- TRUE
+numPerms(1:18, control = ctrl)
+(tmp3 <- allPerms(18, control = ctrl, observed = TRUE))
+(tmp4 <- allPerms(18, control = ctrl))
+## prints out details of the permutation scheme as
+## well as the matrix of permutations
+summary(tmp)
+summary(tmp2)
}
\keyword{ utilities }
\keyword{ design }
More information about the Vegan-commits
mailing list