Fractal Landscapes in R

The past few days I’ve gone on a bit of a fractals kick. It’s partly been motivated by the recent passing of Benoit Mandelbrot, the iconoclastic mathematical genius who coined the term “fractal” and profoundly shook up a number of fields by introducing fractal concepts into them. Many of you have no doubt seen the Mandelbrot Set before—it’s probably the most famous mathematical fractal shape out there. If you haven’t, check out this video (via

Mandelbrot Fractal Set Trip To e214 HD from teamfresh on Vimeo.

The day he died, I received a link to his obituary in an email from a certain humorist I know, with the following in the subject line: “founder of fractal geometry dead, stoners everywhere in mourning.”

But aside from psychedelic entertainment value, fractal geometry turns out to be tremendously useful in describing patterns found in nature. It is especially useful compared with Euclidean geometry, the familiar system of straight lines, squares, triangles, circles, etc. As anyone who has ever looked at nature can tell you, nothing is straight. Crookedness, roughness, and convoluted patterns are the norm. Using euclidean geometry, there is no simple way to describe a coastline, or mountain range, or the fluctuations of animal populations. Fractal geometry, on the other hand, allows us to approximate many of these shapes surprisingly well using very simple equations.

One neat fractal trick is generating realistic-looking mountainscapes using only a few simple operations repeated over and over on themselves (generally speaking, repeating an operation over and over on itself is the way all fractal patterns arise). Consider a two-dimensional case first, for simplicity: we’ll generate the profile of a mountain range, as seen on the horizon. Code to generate this pattern is below. I’ve written my examples here in R, since I’m familiar with it’s graphics system and since my mostly ecological audience is more likely to be familiar with this language than anything else (I presume mostly ecological…who are you people, anyway?).

mountainize <- function(a, roughness, sdev) {
    # Recursive function to create a one-dimensional "mountain range."
    # Arguments:
    #    a: numeric vector
    #    roughness: number between 0 and 1, controlling how 
    #    quickly the random perturbations tail off as the scale decreases.
    #    sdev: standard deviation of the random perturbations.
    # Value:
    #    A numeric vector the same length as a.
    size <- length(a)
    if (size < 2) {
    } else {
        middle <- (size + 1) / 2
        a[middle] <- mean(c(a[1], a[size])) + rnorm(1, 0, sdev)
        a[1:middle] <- mountainize(a[1:middle], 
                                    roughness, sdev * roughness)
        a[middle:size] <- mountainize(a[middle:size], 
                                    roughness, sdev * roughness)

iterations <- 8  # set the number of times to repeat the algorithm
mountains <- rep(0, 2 ^ iterations + 1)  # make an array of the correct size
set.seed(1)   # delete or change the seed to get other mountain ranges
mountains <- mountainize(mountains, 0.5, 1)
plot(mountains, ty='l')

The function mountainize takes a vector of numbers (in this case, just a bunch of zeros to start) and makes the middle number equal to the mean of the two ends, plus a random number from the normal distribution. It then divides the vector into two parts—from the beginning to the middle, and the middle to the end—and repeats the procedure on each half, but with a smaller random deviation. This recursion continues until the sub-section of the vector is only one element long, at which point it unwinds and returns the whole (now jagged) array. The picture below shows what this looks like step-by-step. With each additional iteration, the number of points doubles. After the seventh or eighth time, it looks pretty good. Voila! A mountain range!

Some technical notes. First, the array has to be an appropriate length for the divide-and-conquer algorithm to work—namely, 2i + 1, where i is the desired number of iterations (i.e. the depth of recursion). Second, there is an argument to the mountainize function called roughness. This parameter controls how much the standard deviation of the random numbers decreases from iteration to iteration, and hence how rough or smooth the mountain range ends up. Run the code yourself and play around with it to see how it affects the finished product.

I need to go to sleep, but I will leave off for now with a teaser: With just a little more effort, it is possible to use the same principle to create awesome looking three-dimensional terrain maps.

Don’t you want to know how to make something like this? Tune in tomorrow to find out how…

Tweet about this on TwitterShare on FacebookShare on RedditShare on Google+Pin on PinterestShare on TumblrEmail this to someone
This entry was posted in Uncategorized and tagged , , , , , . Bookmark the permalink.

4 Responses to Fractal Landscapes in R

  1. Pingback: Fractal Landscapes in R: Part Two « Oceanographer’s Choice

  2. Kellen Proctor says:

    Hi Sam,

    I trust you are doing well. I came across your blog on a random google search for “r fractals”. Fascinating. Completely fascinating.

    Anyway, I was trying to run your code, and ran into trouble with the divide function in mountainize. For some reason, R doesn’t recognize the function. Is there a specific package that you have installed and running that includes divide?

    Best Regards,

    Kellen Proctor

    • Sam says:

      That’s a mistake on my part–I changed the name of the function from “divide” to “mountainize” without changing the two recursive calls inside it. If you replace each “divide” with “mountainize” it will work. I’ve corrected the code block in the post as well.

  3. Gordon says:

    Nicely implemented! Though I have a silly question: Lets say I have 2 points on a plot p1 and p2, obviously they form a line segment. How would it be possible to ‘connect’ these 2 points with mountainize() by making the line that connects them ‘jagged’.

    I think I’m missing something in the code… I have already tried to use as input vector ‘a’ the 4 coordinates but it doesn’t work.

Leave a Reply

Your email address will not be published. Required fields are marked *