Michael
Michael

Reputation: 111

How to create matrix of all 2^n binary sequences of length n using recursion in R?

I know I can use expand.grid for this, but I am trying to learn actual programming. My goal is to take what I have below and use a recursion to get all 2^n binary sequences of length n.

I can do this for n = 1, but I don't understand how I would use the same function in a recursive way to get the answer for higher dimensions.

Here is for n = 1:

binseq <- function(n){
  binmat <- matrix(nrow = 2^n, ncol = n)
  r <- 0 #row counter
  for (i in 0:1) {
        r <- r + 1
        binmat[r,] <- i
    }
  return(binmat)
  }

I know I have to use probably a cbind in the return statement. My intuition says the return statement should be something like cbind(binseq(n-1), binseq(n)). But, honestly, I'm completely lost at this point.

The desired output should produce something like what expand.grid gives:

n = 5
expand.grid(replicate(n, 0:1, simplify = FALSE))

It should just be a matrix as binmat is being filled recursively.

Upvotes: 2

Views: 354

Answers (2)

Michael
Michael

Reputation: 111

Apparently this was the desired recursive code:

binseq <- function(n){
  if(n == 1){
    binmat <- matrix(c(0,1), nrow = 2, ncol = 1)
  }else if(n > 1){
    A <- binseq(n-1)
    B <- cbind(rep(0, nrow(A)), A)
    C <- cbind(rep(1, nrow(A)), A)
    binmat <- rbind(B,C)
  }
  return(binmat)
  }

Basically for n = 1 we create a [0, 1] matrix. For every n there after we add a column of 0's to the original matrix, and, separately, a column of 1's. Then we rbind the two matrices to get the final product. So I get what the algorithm is doing, but I don't really understand what the recursion is doing. For example, I don't understand the step from n = 2 to n = 3 based on the algorithm.

Upvotes: 0

whuber
whuber

Reputation: 2494

As requested in a comment (below), here is a limited implementation for binary sequences only:

eg.binary <- function(n, digits=0:1) {
  if (n <= 0) return(matrix(0,0,0))
  if (n == 1) return(matrix(digits, 2))
  x <- eg.binary(n-1)
  rbind(cbind(digits[1], x), cbind(digits[2], x))
}

After taking care of an initial case that R cannot handle correctly, it treats the "base case" of n=1 and then recursively obtains all n-1-digit binary strings and prepends each digit to each of them. The digits are prepended so that the binary strings end up in their usual lexicographic order (the same as expand.grid).

Example:

eg.binary(3)
     [,1] [,2] [,3]
[1,]    0    0    0
[2,]    0    0    1
[3,]    0    1    0
[4,]    0    1    1
[5,]    1    0    0
[6,]    1    0    1
[7,]    1    1    0
[8,]    1    1    1

A general explanation (with a more flexible solution) follows.


Distill the problem down to the basic operation of tacking the values of an array y onto the rows of a dataframe X, associating a whole copy of X with each value (via cbind) and appending the whole lot (via rbind):

cross <- function(X, y) {
  do.call("rbind", lapply(y, function(z) cbind(X, z)))
}

For example,

cross(data.frame(A=1:2, b=letters[1:2]), c("X","Y"))
  A b z
1 1 a X
2 2 b X
3 1 a Y
4 2 b Y

(Let's worry about the column names later.)

The recursive solution for a list of such arrays y assumes you have already carried out these operations for all but the last element of the list. It has to start somewhere, which evidently consists of converting an array into a one-column data frame. Thus:

  eg_ <- function(y) {
    n <- length(y)
    if (n <= 1) {
      as.data.frame(y) 
    } else {
      cross(eg_(y[-n]), y[[n]])
    }
  }

Why the funny name? Because we might want to do some post-processing, such as giving the result nice names. Here's a fuller implementation:

eg <- function(y) {
  # (Define `eg_` here to keep it local to `eg` if you like)
  X <- eg_(y)
  names.default <- paste0("Var", seq.int(length(y)))
  if (is.null(names(y))) {
    colnames(X) <- names.default
  } else {
    colnames(X) <- ifelse(names(y)=="", names.default, names(y))
  }
  X
}

For example:

eg(replicate(3, 0:1, simplify=FALSE))
  Var1 Var2 Var3
1    0    0    0
2    1    0    0
3    0    1    0
4    1    1    0
5    0    0    1
6    1    0    1
7    0    1    1
8    1    1    1
eg(list(0:1, B=2:3))
  Var1 B
1    0 2
2    1 2
3    0 3
4    1 3

Upvotes: 3

Related Questions