[Vegan-commits] r1094 - in pkg: . permute permute/R permute/man
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Mon Jan 4 14:44:44 CET 2010
Author: gsimpson
Date: 2010-01-04 14:44:44 +0100 (Mon, 04 Jan 2010)
New Revision: 1094
Added:
pkg/permute/
pkg/permute/DESCRIPTION
pkg/permute/R/
pkg/permute/R/Blocks.R
pkg/permute/R/Within.R
pkg/permute/R/allPerms.R
pkg/permute/R/getNumObs.R
pkg/permute/R/numPerms.R
pkg/permute/R/permCheck.R
pkg/permute/R/permControl.R
pkg/permute/R/permuplot.R
pkg/permute/R/permute.R
pkg/permute/R/permuted.index2.R
pkg/permute/R/print.allPerms.R
pkg/permute/R/print.permCheck.R
pkg/permute/R/print.permControl.R
pkg/permute/R/print.summary.allPerms.R
pkg/permute/R/print.summary.permCheck.R
pkg/permute/R/summary.allPerms.R
pkg/permute/R/summary.permCheck.R
pkg/permute/man/
pkg/permute/man/allPerms.Rd
pkg/permute/man/getNumObs.Rd
pkg/permute/man/permCheck.Rd
pkg/permute/man/permuted.index2.Rd
Log:
New permute package
Added: pkg/permute/DESCRIPTION
===================================================================
--- pkg/permute/DESCRIPTION (rev 0)
+++ pkg/permute/DESCRIPTION 2010-01-04 13:44:44 UTC (rev 1094)
@@ -0,0 +1,19 @@
+Package: permute
+Title: Functions for generating restricted permutations of data
+Version: 0.0-1
+Date: $Date$
+Author: Gavin L. Simpson
+Maintainer: Gavin L. Simpson <gavin.simpson at ucl.ac.uk>
+Suggests: vegan
+Description: Restricted experimental designs are commonly encountered
+ in ecological data. Ordination analyses of such data are
+ generally tested using permutation tests. The 'permute'
+ package implements a set of restricted permutation designs
+ for freely exchangeable, line transects (time series), and
+ spatial grid designs plus permutation of blocks (groups of
+ samples). 'permute' also allows split-plot designs, in which
+ the whole-plots or split-plots or both can be freely-exchangeble
+ or one of the restricted designs.
+License: GPL-2
+URL: http://vegan.r-forge.r-project.org/
+
Added: pkg/permute/R/Blocks.R
===================================================================
--- pkg/permute/R/Blocks.R (rev 0)
+++ pkg/permute/R/Blocks.R 2010-01-04 13:44:44 UTC (rev 1094)
@@ -0,0 +1,13 @@
+`Blocks` <- function(type = c("free","series","grid","none"),
+ mirror = FALSE, ncol = NULL, nrow = NULL)
+{
+ if(missing(type))
+ type <- "none"
+ else
+ type <- match.arg(type)
+ out <- list(type = type, mirror = mirror,
+ ncol = ncol, nrow = nrow)
+ ## keep as list for now
+ ##class(out) <- "Blocks"
+ return(out)
+}
Added: pkg/permute/R/Within.R
===================================================================
--- pkg/permute/R/Within.R (rev 0)
+++ pkg/permute/R/Within.R 2010-01-04 13:44:44 UTC (rev 1094)
@@ -0,0 +1,14 @@
+`Within` <- function(type = c("free","series","grid","none"),
+ constant = FALSE, mirror = FALSE,
+ ncol = NULL, nrow = NULL)
+{
+ if(missing(type))
+ type <- "free"
+ else
+ type <- match.arg(type)
+ out <- list(type = type, constant = constant, mirror = mirror,
+ ncol = ncol, nrow = nrow)
+ ## keep as default list for now
+ ##class(out) <- "Within"
+ return(out)
+}
Added: pkg/permute/R/allPerms.R
===================================================================
--- pkg/permute/R/allPerms.R (rev 0)
+++ pkg/permute/R/allPerms.R 2010-01-04 13:44:44 UTC (rev 1094)
@@ -0,0 +1,219 @@
+`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,]
+ ## but only if n > 2
+ if(control$mirror && (nperms > 2))
+ 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) {
+ v <- seq_len(n)
+ nperms <- numPerms(v, control)
+ lev <- length(levels(control$strata))
+ X <- matrix(nrow = nperms, ncol = length(control$strata))
+ perms <- if(control$type == "free") {
+ all.free(lev)
+ } else if(control$type == "series") {
+ all.series(lev, control = control)
+ } else {
+ all.grid(lev, control = control)
+ }
+ sp <- split(v, control$strata)
+ for(i in seq_len(nrow(perms)))
+ X[i,] <- unname(do.call(c, sp[perms[i,]]))
+ X
+ }
+ ## replacement for recursive function above
+ bar <- function(mat, n) {
+ res <- vector(mode = "list", length = n)
+ for(i in seq_len(n))
+ res[[i]] <- mat
+ do.call(rbind, res)
+ }
+ ## start
+ v <- n
+ ## expand n if a numeric or integer vector of length 1
+ if((is.numeric(n) || is.integer(n)) && (length(n) == 1))
+ v <- seq_len(n)
+ ## number of observations in data
+ n <- getNumObs(v)
+ ## check permutation scheme and update control
+ pcheck <- permCheck(v, control = control, make.all = FALSE)
+ control <- pcheck$control
+ ## get max number of permutations
+ 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)) {
+ if(!control$permute.strata && !is.null(control$strata)) {
+ ## permuting within blocks
+ ## FIXME: allperms expects samples to be arranged
+ ## in order of fac, i.e. all level 1, followed by
+ ## all level 2 - fix to allow them to be in any order:
+ ## see permuted.index2 for how to do this
+ 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)
+ nperms <- 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 = nperms, ncol = n)
+ for(i in seq_len(perm.wi))
+ res[i,] <- sapply(sp, function(x) x[ord[i,]])
+ } else {
+ ## different permutations within blocks
+ tab <- table(control$strata)
+ ng <- length(tab)
+ pg <- unique(tab)
+ if(length(pg) > 1) {
+ ## different number of observations per level of strata
+ if(control$type == "grid")
+ ## FIXME: this should not be needed once all checks are
+ ## in place in permCheck()
+ stop("Unbalanced grid designs are not supported")
+ control.wi <- permControl(type = control$type,
+ mirror = control$mirror)
+ sp <- split(v, control$strata)
+ res <- vector(mode = "list", length = ng)
+ add <- c(0, cumsum(tab)[1:(ng-1)])
+ for(j in seq(along = tab)) {
+ ord <- switch(control.wi$type,
+ free = all.free(tab[j]),
+ series = all.series(tab[j],
+ control=control.wi))
+ perm.wi <- nrow(ord)
+ if(j == 1) {
+ a <- 1
+ b <- nperms / perm.wi
+ } else {
+ b <- b/perm.wi
+ a <- nperms / (b*perm.wi)
+ }
+ res[[j]] <- matrix(rep(bar(ord+add[j], a),
+ each = b),
+ ncol = tab[j])
+ }
+ res <- do.call(cbind, res)
+ } else {
+ ## same number of observations per level of strata
+ 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 {
+ ## not permuting within blocks or are permuting strata
+ 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/permute/R/getNumObs.R
===================================================================
--- pkg/permute/R/getNumObs.R (rev 0)
+++ pkg/permute/R/getNumObs.R 2010-01-04 13:44:44 UTC (rev 1094)
@@ -0,0 +1,16 @@
+`getNumObs` <- function(object, ...) UseMethod("getNumObs")
+
+`getNumObs.default` <- function(object, ...)
+{
+ NROW(scores(object, display = "sites"))
+}
+
+`getNumObs.numeric` <- function(object, ...)
+{
+ length(object)
+}
+
+`getNumObs.integer` <- function(object, ...)
+{
+ getNumObs.numeric(object)
+}
Added: pkg/permute/R/numPerms.R
===================================================================
--- pkg/permute/R/numPerms.R (rev 0)
+++ pkg/permute/R/numPerms.R 2010-01-04 13:44:44 UTC (rev 1094)
@@ -0,0 +1,100 @@
+`numPerms` <- function(object, control = permControl())
+{
+ ## constant holding types where something is permuted
+ PTYPES <- c("free","grid","series")
+ ## expand object if a numeric or integer vector of length 1
+ if((is.numeric(object) || is.integer(object)) && (length(object) == 1))
+ object <- seq_len(object)
+ ## number of observations in data
+ nobs <- getNumObs(object)
+ ## within perms object
+ WITHIN <- control$within
+ ## strata perms object
+ BLOCKS <- control$blocks
+ ## are strata present?
+ STRATA <- !is.null(control$strata)
+ ## check that when permuting strata or constant within strata,
+ ## strata have same number of samples
+ if(STRATA) {
+ tab.strata <- table(control$strata)
+ same.n <- length(unique(tab.strata))
+ if((BLOCKS$type %in% PTYPES || isTRUE(WITHIN$constant)) && same.n > 1)
+ stop("All levels of strata must have same number of samples for chosen scheme")
+ if(BLOCKS$type == "grid" && same.n > 1)
+ stop("Unbalanced grid designs are not supported")
+ }
+ ## generate multiplier for restricted permutations
+ if(WITHIN$type %in% c("series","grid")) {
+ within.multi <- 2
+ if(WITHIN$type == "grid" && WITHIN$ncol > 2) {
+ within.multi <- 4
+ } else {
+ if(nobs == 2)
+ within.multi <- 1
+ }
+ }
+ if(BLOCKS$type %in% c("series","grid")) {
+ blocks.multi <- 2
+ if(BLOCKS$type == "grid" && BLOCKS$ncol > 2) {
+ blocks.multi <- 4
+ } else {
+ if(nobs == 2)
+ blocks.multi <- 1
+ }
+ }
+ ## calculate number of possible permutations
+ ## blocks
+ num.blocks <- 1
+ if(BLOCKS$type %in% PTYPES) {
+ if(BLOCKS$type == "free")
+ num.blocks <- exp(lfactorial(length(levels(control$strata))))
+ else if(BLOCKS$type %in% c("series","grid")) {
+ if(BLOCKS$mirror)
+ num.blocks <- blocks.multi * nobs
+ else
+ num.blocks <- nobs
+ }
+ }
+ ## within
+ num.within <- 1
+ if(WITHIN$type %in% PTYPES) {
+ if(WITHIN$type == "free") {
+ if(STRATA)
+ num.within <- prod(factorial(tab.strata))
+ else
+ num.within <- exp(lfactorial(nobs))
+ } else if(WITHIN$type %in% c("series","grid")) {
+ if(STRATA) {
+ if(same.n > 1) {
+ multi <- rep(2, length = length(tab.strata))
+ multi[which(tab.strata == 2)] <- 1
+ if(WITHIN$mirror) {
+ num.within <- prod(multi * tab.strata)
+ } else {
+ num.within <- prod(tab.strata)
+ }
+ } else {
+ if(WITHIN$mirror) {
+ if(WITHIN$constant)
+ num.within <- within.multi * unique(tab.strata)
+ else
+ num.within <- prod(within.multi * tab.strata)
+ } else {
+ if(WITHIN$constant)
+ num.within <- unique(tab.strata)
+ else
+ num.within <- prod(tab.strata)
+ }
+ }
+ } else {
+ if(WITHIN$mirror)
+ num.within <- within.multi * nobs
+ else
+ num.within <- nobs
+ }
+ } else {
+ stop("Ambiguous permutation type in 'control$type'")
+ }
+ }
+ return(num.blocks * num.within)
+}
Added: pkg/permute/R/permCheck.R
===================================================================
--- pkg/permute/R/permCheck.R (rev 0)
+++ pkg/permute/R/permCheck.R 2010-01-04 13:44:44 UTC (rev 1094)
@@ -0,0 +1,54 @@
+`permCheck` <- function(object, control = permControl(),
+ make.all = TRUE)
+{
+ ## if object is numeric or integer and of length 1,
+ ## extend the object
+ if(length(object) == 1 &&
+ (is.integer(object) || is.numeric(object)))
+ object <- seq_len(object)
+ ## check the number of observations in object
+ nobs <- getNumObs(object)
+ ## sample permutation type
+ type <- control$within$type
+ ## if strata, check nobs == length of strata
+ ## but beware empty levels
+ if(!is.null(control$strata)) {
+ tab <- table(control$strata)
+ if(!identical(as.integer(nobs), as.integer(sum(tab))))
+ stop("Number of observations and length of 'strata' do not match.")
+ ## if "grid", check design balanced?
+ if((bal <- length(unique(tab))) > 1 && type == "grid")
+ stop("Unbalanced 'grid' designs are not supported.")
+ ## if grid design, check nrow*ncol is multiple of nobs
+ if(type == "grid" &&
+ !identical(nobs %% (control$within$ncol *
+ control$within$nrow), 0))
+ stop("'nrow' * 'ncol' not a multiple of number of observations.")
+ ## if constant, check design balanced?
+ if(control$within$constant && bal > 1)
+ stop("Unbalanced designs not allowed with 'constant = TRUE'.")
+ ## if permuting strata, must be balanced
+ if(control$blocks$type != "none" && bal > 1)
+ stop("Design must be balanced if permuting 'strata'.")
+ }
+ ## check allPerms is of correct form
+ if(!is.null(control$all.perms) &&
+ !identical(class(control$all.perms), "allPerms"))
+ stop("'control$all.perms' must be of class 'allPerms'.")
+ ## get number of possible permutations
+ num.pos <- numPerms(object, control)
+ ## if number of possible perms < minperm turn on complete enumeration
+ if(num.pos < control$minperm) {
+ control$nperm <- control$maxperm <- num.pos
+ control$complete <- TRUE
+ }
+ ## if complete enumeration, generate all permutations
+ if(control$complete && make.all) {
+ control$all.perms <- allPerms(nobs, control = control,
+ max = control$maxperm,
+ observed = FALSE)
+ }
+ retval <- list(n = num.pos, control = control)
+ class(retval) <- "permCheck"
+ retval
+}
Added: pkg/permute/R/permControl.R
===================================================================
--- pkg/permute/R/permControl.R (rev 0)
+++ pkg/permute/R/permControl.R 2010-01-04 13:44:44 UTC (rev 1094)
@@ -0,0 +1,14 @@
+`permControl` <- function(strata = NULL, nperm = 199, complete = FALSE,
+ within = Within(),
+ blocks = Blocks(),
+ maxperm = 9999, minperm = 99,
+ all.perms = NULL)
+{
+ out <- list(strata = strata, nperm = nperm, complete = complete,
+ within = within, blocks = blocks,
+ maxperm = maxperm, minperm = minperm,
+ all.perms = all.perms,
+ name.strata = deparse(substitute(strata)))
+ class(out) <- "permControl"
+ return(out)
+}
Added: pkg/permute/R/permuplot.R
===================================================================
--- pkg/permute/R/permuplot.R (rev 0)
+++ pkg/permute/R/permuplot.R 2010-01-04 13:44:44 UTC (rev 1094)
@@ -0,0 +1,192 @@
+`permuplot` <- function(n, control = permControl(),
+ col = par("col"),
+ hcol = "red",
+ shade = "lightgrey",
+ xlim=NULL, ylim=NULL,
+ inset = 0.1,
+ main=NULL, sub=NULL,
+ ann = par("ann"),
+ cex = par("cex"),
+ ...) {
+ xy.series <- function(n) {
+ angle <- seq(0, 2*pi, length = n+1)[-(n+1)]
+ x <- rev(cos(angle))
+ y <- rev(sin(angle))
+ xy <- xy.coords(x, y)
+ return(xy)
+ }
+ xy.free <- function(n) {
+ x <- runif(n)
+ y <- runif(n)
+ xy <- xy.coords(x, y)
+ return(xy)
+ }
+ xy.grid <- function(ncol, nrow) {
+ x <- rep(seq_len(ncol), each = nrow)
+ y <- rev(rep(seq_len(nrow), times = ncol))
+ xy <- xy.coords(x, y)
+ return(xy)
+ }
+ axis.limits <- function(vals, inset) {
+ lim <- range(vals[is.finite(vals)])
+ lim.range <- lim[2] - lim[1]
+ res <- c(lim[1] - (lim.range * inset),
+ lim[2] + (lim.range * inset))
+ return(res)
+ }
+ ## currently doesn't support restricted permutations of strata themselves
+ if(control$permute.strata && control$type != "free")
+ stop("Restricted permutations of strata currently not supported")
+ ## check that n and length of strata are equal
+ if( use.strata <- !is.null(control$strata) ) {
+ tab <- table(control$strata)
+ if(!identical(as.integer(sum(tab)), as.integer(n)))
+ stop("'n' and length of 'strata' don't match.")
+ }
+ ## check the control design
+ control <- permCheck(n, control = control)$control
+ if(use.strata) {
+ n.grp <- length(tab)
+ opar <- par(no.readonly=TRUE, mar=c(2,2,2,1)+0.1,
+ mfrow = n2mfrow(n.grp),
+ oma=c(2.1,0,3.1,0))
+ on.exit(par(opar))
+ ## if permuting strata, only need to draw the sub-plots
+ ## in a different order
+ if(control$permute.strata) {
+ ## expand shade, col
+ if(identical(length(col), 1))
+ col <- rep(col, n.grp)
+ if(identical(length(shade), 1))
+ shade <- rep(shade, n.grp)
+ ord <- sample(names(tab))
+ if(is.null(xlim))
+ xlim <- c(0,1)
+ if(is.null(ylim))
+ ylim <- c(0,1)
+ xy <- xy.coords(0.5, 0.5)
+ string <- paste("Stratum:\n", ord)
+ names(string) <- ord
+ strh <- max(strheight(string, cex = cex))
+ strw <- max(strwidth(string, cex = cex))
+ box.coords <- xy.coords(rep(c(0.5-strw, 0.5+strw), each = 2),
+ c(0.5-strh, 0.5+strh,
+ 0.5+strh, 0.5-strh))
+ for(i in ord) {
+ plot.new()
+ plot.window(xlim, ylim, asp = 1, ...)
+ polygon(box.coords, col = shade, border = hcol, ...)
+ text(xy$x, xy$y, labels = string[i],
+ col = col, cex = cex, ...)
+ box()
+ #if(ann) {
+ # title(main = paste("Original order:",
+ # which(ord == i)))
+ #}
+ }
+ } else {
+ ## if free and constant, only need one set of random coords
+ xy <- if(control$constant && control$type == "free") {
+ ## needs to be a list for the main loop below
+ xy <- xy.free(unique(tab))
+ res <- vector("list", length = length(tab))
+ for(i in seq_along(res)) {
+ res[[i]] <- xy
+ }
+ res
+ } else {
+ switch(control$type,
+ free = lapply(tab, xy.free),
+ series = lapply(tab, xy.series),
+ grid = lapply(tab, function(x) {
+ xy.grid(control$ncol, control$nrow)
+ }),
+ stop("Unsupport permutation 'type'"))
+ }
+ perms <- permuted.index2(n, control = control)
+ perms <- tapply(perms, control$strata, function(x) x)
+ if(is.null(main))
+ main <- paste("Stratum:", names(tab))
+ for(i in seq_along(xy)) {
+ if(is.null(xlim))
+ xlim <- axis.limits(xy[[i]]$x, inset)
+ if(is.null(ylim))
+ ylim <- axis.limits(xy[[i]]$y, inset)
+ plot.new()
+ plot.window(xlim, ylim, asp = 1, ...)
+ cols <- switch(control$type,
+ free = rep(col, tab[i]),
+ series = c(hcol, rep(col, tab[i]-1)),
+ grid = {cols <- rep(col, tab[i])
+ cols[which.min(perms[[i]])] <-
+ hcol
+ cols})
+ text(xy[[i]]$x, xy[[i]]$y, labels = perms[[i]],
+ col = cols, ...)
+ if(ann) {
+ title(main = main[i], ...)
+ title(sub = paste("n in stratum:", tab[i]),
+ line = 0.5, ...)
+ }
+ box()
+ }
+ }
+ if(ann) {
+ sub <- paste(paste("n: ", n, ";", sep = ""),
+ paste("mirror: ", control$mirror, ";",
+ sep = ""),
+ paste("constant: ", control$constant, ";",
+ sep = ""),
+ sep = " ")
+ if(control$type == "grid")
+ sub <- paste(sub, paste("ncol: ",
+ control$ncol, ";",
+ sep = ""),
+ paste("nrow: ", control$nrow, ";",
+ sep = ""),
+ sep = " ")
+ title(main = paste("Permutation type:", control$type),
+ outer = TRUE, cex.main = 1.75, ...)
+ title(sub = sub, outer = TRUE, line = 0.5,
+ cex.sub = 1, ...)
+ }
+ } else {
+ xy <- switch(control$type,
+ free = xy.free(n),
+ series = xy.series(n),
+ grid = xy.grid(control$ncol, control$nrow),
+ stop("Unsupport permutation 'type'"))
+ if(is.null(xlim)) {
+ xlim <- axis.limits(xy$x, inset)
+ }
+ if(is.null(ylim)) {
+ ylim <- axis.limits(xy$y, inset)
+ }
+ opar <- par(no.readonly=TRUE, mar=c(2,1,3,1)+0.1)
+ on.exit(par(opar))
+ if(is.null(main))
+ main <- paste("Permutation type:", control$type)
+ if(is.null(sub))
+ sub <- paste(paste("n: ", n, ";", sep = ""),
+ paste("mirror: ", control$mirror, ";",
+ sep = ""),
+ sep = " ")
+ plot.new()
+ plot.window(xlim, ylim, asp = 1, ...)
+ labs <- permuted.index2(n, control=control)
+ cols <- switch(control$type,
+ free = rep(col, n),
+ series = c(hcol, rep(col, n-1)),
+ grid = {cols <- rep(col, n)
+ cols[which.min(labs)] <- hcol
+ cols})
+ text(xy$x, xy$y, labels = labs,
+ col = cols, ...)
+ if(ann) {
+ title(main = main, ...)
+ title(sub = sub, line = 0.5, ...)
+ }
+ box()
+ }
+ invisible()
+}
Added: pkg/permute/R/permute.R
===================================================================
--- pkg/permute/R/permute.R (rev 0)
+++ pkg/permute/R/permute.R 2010-01-04 13:44:44 UTC (rev 1094)
@@ -0,0 +1,10 @@
+permute <- function(i, n, control) {
+ if(control$complete && !is.null(control$all.perms))
+ perm <- control$all.perms[i,]
+ else {
+ if(control$complete)
+ warning("'$all.perms' is NULL, yet '$complete = TRUE'.\nReturning a random permutation.")
+ perm <- permuted.index(n, control)
+ }
+ return(perm)
+}
Added: pkg/permute/R/permuted.index2.R
===================================================================
--- pkg/permute/R/permuted.index2.R (rev 0)
+++ pkg/permute/R/permuted.index2.R 2010-01-04 13:44:44 UTC (rev 1094)
@@ -0,0 +1,155 @@
+`permuted.index` <-
+ function (n, control = permControl())
+{
+ `pStrata` <- function(strata, type, mirror = FALSE, start = NULL,
+ flip = NULL, nrow, ncol, start.row = NULL,
+ start.col = NULL) {
+ lev <- length(levels(strata))
+ ngr <- length(strata) / lev
+ sp <- split(seq(along = strata), strata)
+ if(type == "free") {
+ unname(do.call(c, sp[pFree(lev, lev)]))
+ } else if(type == "series") {
+ unname(do.call(c,
+ sp[pSeries(seq_len(lev),
+ mirror = mirror,
+ start = start,
+ flip = flip)]))
+ } else if(type == "grid") {
+ unname(do.call(c,
+ sp[pGrid(nrow = nrow, ncol = ncol,
+ mirror = mirror,
+ start.row = start.row,
+ start.col = start.col,
+ flip = flip)]))
+ } else {
+ stop("Invalid permutation type.")
+ }
+ }
+ `pGrid` <- function(nrow, ncol, mirror = FALSE, start.row = NULL,
+ start.col = NULL, flip = NULL) {
+ if(is.null(start.row))
+ start.row <- pFree(nrow, 1L)
+ if(is.null(start.col))
+ start.col <- pFree(ncol, 1L)
+ ir <- seq(start.row, length=nrow) %% nrow
+ ic <- seq(start.col, length=ncol) %% ncol
+ if(!is.null(flip)) {
+ if(any(flip)) {
+ if(flip[1L])
+ ir <- rev(ir)
+ if(flip[2L])
+ ic <- rev(ic)
+ }
+ } else {
+ if (mirror) {
+ if (runif(1L) < 0.5)
+ ir <- rev(ir)
+ if (runif(1L) < 0.5)
+ ic <- rev(ic)
+ }
+ }
+ rep(ic, each=nrow) * nrow + rep(ir, len=nrow*ncol) + 1L
+ }
+ `pSeries` <- function(inds, mirror = FALSE, start = NULL,
+ flip = NULL) {
+ n <- length(inds)
+ if(is.null(start))
+ start <- pFree(n, 1L)
+ out <- seq(start, length = n) %% n + 1L
+ if(!is.null(flip)) {
+ if(flip)
+ out <- rev(out)
+ } else {
+ if(mirror && runif(1L) < 0.5)
+ out <- rev(out)
+ }
+ inds[out]
+ }
+ `pFree` <- function(x, size)
+ .Internal(sample(x, size, FALSE, NULL))
+ ## END in-line Functions ##
+
+ ## If no strata then permute all samples using stated scheme
+ if(is.null(control$strata)) {
+ out <-
+ switch(control$within$type,
+ "free" = pFree(n, n),
+ "series" =
+ pSeries(1:n, mirror = control$within$mirror),
+ "grid" =
+ pGrid(nrow = control$within$nrow,
+ ncol = control$within$ncol,
+ mirror = control$within$mirror)
+ )
+ } else {
+ ## If strata present, either permute samples, strata or both
+
+ ## permute strata?
+ if(control$blocks$type == "none") {
+ out <- 1:n
+ } else {
+ flip <- runif(1L) < 0.5
+ out <- pStrata(control$strata,
+ type = control$blocks$type,
+ mirror = control$blocks$mirror,
+ flip = flip,
+ nrow = control$blocks$nrow,
+ ncol = control$blocks$ncol)
+ }
+ ## permute the samples within strata?
+ if(control$within$type != "none") {
+ tab <- table(control$strata[out])
+ ## the levels of the strata
+ inds <- names(tab)
+ ## same permutation within each level of strata?
+ if(control$within$constant) {
+ if(control$within$type == "free") {
+ n <- unique(tab)[1L]
+ same.rand <- pFree(n, n)
+ } else if(control$within$type == "series") {
+ start <- pFree(n / length(inds), 1L)
+ flip <- runif(1L) < 0.5
+ } else if(control$within$type == "grid") {
+ start.row <- pFree(control$within$nrow, 1L)
[TRUNCATED]
To get the complete diff run:
svnlook diff /svnroot/vegan -r 1094
More information about the Vegan-commits
mailing list