Clojure: Barnsley Fern

In addition to other pleasurable events, the last month I got as a gift yet another book devoted to the Chaos Theory – ‘Simply Chaos’ by Sergey Demenok (in Russian). This one vaguely reminds me of ‘Chaos: Making a New Science’ by James Gleick, since both aim mostly at attracting readers to science, explaining it to them and showing how beautiful it can be. Still they differ a lot, with ‘Simply Chaos’ being more joyful reading balancing on the boundary between historical and scientifical truth on one side and myth and fantasy on the other. This way the last book was a great and funny refresher, making me remember of some things, that fascinated me while I was reading the Gleick’s ‘Chaos’ a year and a half ago. Most notably, I faced the idea of Barnsley Fern once more.

 

The fern is a mathematical system. Its key feature is the same as with all the chaotic systems, which I have explored so far: being very simple in its essence it produces fascinating and unexpected behavior. Like the logistic or tent maps, it is defined in the form of iterated system, so that starting at some point we plug its coordinates into the equations and arrive at the next point. Feeding the latter to the same equations we get to the third point and so on – the resulting sequence of points forms a trajectory of the system. On the other side, unlike the mentioned chaotic maps, the Barnsley’s system is two-dimensional, although that’s not the most important difference between these systems.

While logistic and tent maps are entirely deterministic and do not involve any coin tossing, the fern is a probabilistic system. In its core lie four different affine transforms, which give way to new points. The choice of the transform to be used on a particular step is made by a fair coin – or, to be more precise, by a random number generator. I lied a bit when I said that the coin is fair, because the probabilities assigned to each of the four transforms are not equal, hence the coin is unfair too – leave the fact that it has four sides. This way, each of the affine transforms constituting the system is characterized by seven numbers. Four of them are the quotients for scaling and rotating, another two – for translation and the last one is the probability of selecting this particular transform to produce the next point. This looks like something great to implement in a simple program, ain’t it?

I started with the affine transforms. Since I wanted the ability to create ones with different quotients, I wrote a Clojure function taking the quotients as arguments and returning a new function, which represents the transform itself. There is no need to pass the probability associated with the affinity to this routine, because it is required for the thing that choses among transforms, but not for the transforms themselves.

(defn affinity [a b c d e f]
 (fn [x y]
  [(+ (* a x) (* b y) e) (+ (* c x) (* d y) f)]))

Another crucial component of our program is the function making the choice. In my case it doesn’t even care about the fact that it picks transforms – it would readily select tomatoes as well. Actually, I had to write three functions to accomplish the task, two of which are helpers.

(defn scale [probs]
 (let [s (apply + probs)]
  (map #(/ % s) probs)))

(defn accumulate [probs]
 (conj
  (vec (take (dec (count probs))
   (drop 1 (reduce 
    #(conj %1 (+ (last %1) %2)) 
    [0] 
    probs))))
  1))

The first one does precisely what its name suggests – it scales probabilities to fit them into zero-one interval. The second helper transforms scaled probabilities into a sequence of thresholds. Any threshold is a sum of probabilities assigned to options up to the current one. For instance, if we have four options with probabilities 10%, 20%, 55% and 15%, the thresholds will be 10%, 30%, 85% and 100%. We use the thresholds to pick transforms in the chooser function.

(defn chooser [probs options]
 (let [threshes (accumulate (scale probs))
       chsr (for [[t o] (map vector threshes options)]
     [#(< % t) o])]
  (fn []
   (let [v (clojure.core/rand)]
      (loop [[checker option] (first chsr) other (next chsr)]
       (if (checker v)
        option
        (recur (first other) (next other))))))))

Several things go on here. The function takes two sequences: one with the probabilities for each option and the other with the options themselves. First, it scales probabilities and computes the thresholds by means of the above helpers. Next, it assembles a sequence of options paired with the functions (created via lambda expression), which simply check whether an argument lies below the threshold associated with an option. Once this last sequence is in place, it is easy to create the so much desired function yielding different options with different probabilities – that’s precisely what our chooser returns. This selector-function generates a random number and then successively applies the threshold-functions to it. Once one yielding truth is reached it returns the corresponding option to the caller. Sincerely speaking, I can’t say what makes me fear the task of creating such a simple thing.

By this moment, we have the functions, responsible for creating affinities and offering us a choice mechanism. That’s enough to accomplish our main goal – create a Barnsley’s system implementation. Moreover, with all this utility in place the last step looks straightforward:

(defn fern [transforms]
 (let [selector (chooser
  (map #(nth % 6) transforms)
  (map #(apply affinity (take 6 %)) transforms))]
   (fn [[x y]]
    ((selector) x y))))

For the sake of simplicity, I've decided that the fern function should accept a sequence of vectors, holding seven numbers each – the affinity quotients, followed by the corresponding probability. Most of work in fern is actually done to separate these from each other and properly call the chooser function. Its first argument (probabilities) is assembled by taking the seventh (sixth zero-based) element of each incoming vector. The second one – the collection of transforms – comes as a result of applying affinity function to the first six elements (quotients) of each incoming sequence. The fern returns one more function (yeah, that’s FP), which, when passed two coordinates of a point, calls the selector to get the transform and applies it to the point, thus returning the next one. That’s all – now a trajectory of our system can be easily generated by means of the built-in iterate function:

(def default-fern
 [[0.0 0.0 0.0 0.16 0.0 0.0 0.01]
  [0.85 0.04 -0.04 0.85 0.0 1.6 0.85]
  [0.2 -0.26 0.23 0.22 0.0 1.6 0.07]
  [-0.15 0.28 0.26 0.24 0.0 0.44 0.07]])

(take 100 (iterate (fern default-fern) [0.0 0.0]))

The above code produces the first 100 two-dimensional points of the original version of the Barnsley fern. The final step should be attaching some visualization to our herd of functions. I haven’t done this in Clojure, since the fern fits nicely into a section of my site devoted to simple chaotic systems. Now it has its own page, so feel free to visit it to experiment with parameters and watch the impressive pictures produced by Barnsley’s system.

CodeProject