1.listMat <- function(M, f) {
2  l <- vector("list", ncol(M))
3  names(l) <- dimnames(M)[[2]]
4  for(i in 1:ncol(M))
5      l[[i]] <- f(M[, i])
6  l
7}
8
9.paste0 <- function(...) paste(..., sep = "")
10
11.zf <- function(x) {
12  x <- as.character(x)
13  m <- max(n <- nchar(x))
14  z <- paste(rep(0, m), collapse="")
15  .paste0(substring(z, 0, m-n), x)
16}
17
18conf.design <- function(G, p, block.name = "Blocks",
19                        treatment.names = NULL) {
20  if(!is.matrix(G))
21        G <- rbind(G)
22  if(is.null(treatment.names))
23      treatment.names <- if(is.null(nam <- dimnames(G)[[2]]))
24          .paste0("T", .zf(1:ncol(G))) else nam
25
26  stopifnot(is.character(treatment.names), length(treatment.names) == ncol(G),
27            is.numeric(G), all(G >= 0), all(G %% 1 == 0),
28            is.numeric(p), length(p) == 1, p > 0, p %% 1 ==0,
29            p %in% primes(p),
30            is.character(block.name), length(block.name) == 1)
31
32  D <- as.matrix(expand.grid(rep(list(0:(p-1)), ncol(G))))
33  B <- .listMat((D %*% t(G)) %% p, format)
34  o <- do.call(order, B)
35  B <- do.call(".paste0", B)
36  D <- cbind(B, format(D))[o,  ]
37  dimnames(D) <- list(B[o], c(block.name, treatment.names))
38  data.frame(.listMat(D, factor))
39}
40
41.space <- function(G, p) {
42###
43### Generate all distinct linear combinations of the rows of G over GF(p)
44###
45  x <- 0:(p - 1)
46  M <- as.matrix(x)
47  k <- nrow(G)
48  if(k > 1) {
49    for(i in 2:k) {
50      N <- NULL
51      for(j in x)
52          N <- rbind(N, cbind(M, j))
53      M <- N
54    }
55  }
56  M <- (M %*% G) %% p
57###
58### if the rows of G can be assumed linearly independent the rest can
59### be omitted.
60###
61  m <- 0
62  for(j in 1:ncol(M))
63      m <- p * m + M[, j]
64  M[!duplicated(m),  , drop = FALSE]
65}
66
67conf.set <- function(G, p) {
68  S <- .space(G, p)
69  dimnames(S)[[2]] <- dimnames(G)[[2]]
70  S[apply(S, 1, function(x)
71          any(t <- x > 0) && x[t][1] == 1),  ]
72}
73
74direct.sum <- function(D1, ..., tiebreak = letters) {
75  l <- list(...)
76  if(length(l) == 0) return(D1)
77  D2 <- if(length(l) == 1) l[[1]] else {
78    Recall(..., tiebreak = c(tiebreak[-1], tiebreak[1]))
79  }
80  stopifnot(is.data.frame(D1), is.data.frame(D2))
81  n1 <- nrow(D1)
82  n2 <- nrow(D2)
83  D1 <- D1[rep(seq_len(n1), each  = n2), ]
84  D2 <- D2[rep(seq_len(n2), times = n1), ]
85  if(any(dups <- names(D2) %in% names(D1)))
86      names(D2)[dups] <- .paste0(names(D2)[dups], tiebreak[1])
87  D <- cbind(D1, D2)
88  row.names(D) <- format(seq_len(nrow(D)))
89  D
90}
91
92factorize <- function(x, ...)
93  UseMethod("factorize")
94
95factorize.default <- function(x, divisors = primes(max(x)), ...) {
96  stopifnot(is.numeric(x))
97  if (length(x) > 1) {
98    l <- vector("list", length(x))
99    names(l) <- as.character(x)
100    for (i in seq_along(x))
101      l[[i]] <- Recall(x[i], divisors = divisors, ...)
102    return(l)
103  }
104  if (x %% 1 > 0 || x < 2)
105    return(x)
106  tab <- divisors
107  fac <- numeric(0)
108  while(length(tab <- tab[x %% tab == 0]) > 0) {
109    x <- x/prod(tab)
110    fac <- c(fac, tab)
111  }
112  sort(fac)
113}
114
115factorize.factor <- function(x, name = deparse(substitute(x)),
116                             extension = letters, ...) {
117  llev <- factorize.default(length(levels(x)))
118  if(length(llev) == 1)
119    return(structure(data.frame(x), names = name))
120
121  D <- expand.grid(lapply(llev, seq_len))[x, ] - 1
122
123  D <- lapply(D, factor)
124  names(D) <- .paste0(name, rep(extension, length.out = length(D)))
125  D <- data.frame(D)
126  row.names(D) <- format(seq_len(nrow(D)))
127  D
128}
129
130.makeList <- function(x) if(is.list(x)) x else list(x)
131
132join <- function(...) {
133  l <- list(...)
134  l <- do.call(c, lapply(l, .makeList))
135  l <- lapply(l, function(f)
136              format(as.character(f)))
137  factor(do.call(paste, c(l, list(sep = ":"))))
138}
139
140rjoin <- function(..., part.name = "Part") {
141  l <- lapply(list(...), as.data.frame)	### for some safety...
142  checkNames <- all(sapply(l, function(bit)
143                           all(names(bit) == names(l[[1]]))))
144  if(!checkNames) stop("The column names differ between components!")
145  bf <- factor(.paste0(part.name, .zf(rep(1:length(l), sapply(l, nrow)))))
146  D <- data.frame(Part = bf, do.call("rbind", l))
147  names(D)[1] <- part.name
148  D
149}
150
151primes <- function(n) {
152### Find all primes less than n (or max(n) if length(n) > 1).
153### Uses an obvious sieve method.  Nothing flash.
154###
155### 2013: This now uses a slightly improved coding of the version of
156###       the algorithm used in the pracma package primes() function
157###
158  stopifnot(is.numeric(n), all(n %% 1 == 0))
159  if ((M2 <- max(n)) <= 1)
160    return(numeric(0))
161  x <- seq(1, M2, 2)
162  np <- length(x)
163  x[1] <- 2
164  if(M2 > 8) {
165    top <- floor(sqrt(M2))
166    p <- 1
167    while((p <- p+2) <= top)
168        if(x[(p + 1)/2] > 0)
169            x[seq((p*p + 1)/2, np, p)] <- 0
170  }
171  x[x > 0]
172}
173
174
175
176