2D lattice random walk plots in functional style
Looks good.  The code is functional with the exception of using
return, which is not necessary as the last value will be returned
anyway, so I'd simplify to the following and also in the interested of
efficiency only compute the coords once and without seq, as : is
shorter:

step &lt;- function(x) {
        unlist(sample(possibleSteps, 1))
}

takeRandomWalk &lt;- function(nSteps) {
        coordPairs &lt;- Reduce(`+`, lapply(1:nSteps, step), accumulate = T)
        x &lt;- sapply(coordPairs, `[`, 1)
        y &lt;- sapply(coordPairs, `[`, 2)
        list(x, y)
}

plotRandomWalk &lt;- function(nSteps, margins) {
        walkObj &lt;- takeRandomWalk(nSteps)
        coords &lt;- -margins:margins
        plot(coords, coords, type = 'n', xlab = "", ylab = "")
        lines(walkObj[[1]], walkObj[[2]])
}


Reduce with accumulate is pretty cool, I think this is the first
time I've seen an accumulate parameter on a reduce-like function
anywhere.

Now I think you know that, but usually R code like this would be written
with less lists and more data frames or matrixes.  That doesn't mean the
code isn't functional:  As long as it's using functions, not mutating
things and using functional abstractions instead of imperative
constructs it's fine.

So the following is probably more idiomatic and at the same time uses
more efficient representations for the data:

possibleSteps &lt;- matrix(c(-1, 0, 0, -1, 1, 0, 0, 1), ncol = 2, byrow = TRUE,
                        dimnames = list(NULL, c("X", "Y")))


Gives a matrix with named columns like so, which will be the representation
throughout the rest of the code (which is good for consistency):

&gt; possibleSteps
      X  Y
[1,] -1  0
[2,]  0 -1
[3,]  1  0
[4,]  0  1


The random walk will still be created by sample, except it's sampling
indexes.  cbind will create a new matrix from the two accumulated sums
(cumsum):

takeRandomWalk &lt;- function(nSteps) {
        indexes &lt;- sample(1:dim(possibleSteps)[1], nSteps, TRUE)
        walk &lt;- possibleSteps[indexes,]
        cbind(X = cumsum(walk[,1]), Y = cumsum(walk[,2]))
}


And finally plotRandomWalk is a bit simpler as we can just give
lines the constructed coordinates matrix:

plotRandomWalk &lt;- function(nSteps, margins) {
        coords &lt;- -margins:margins
        plot(coords, coords, type = 'n', xlab = "", ylab = "")
        lines(takeRandomWalk(nSteps))
}