Like forward Euler, and like backward Euler, the wolff algorithm is a useful thing to talk about in a class because it helps people understand the basic ideas, but no one really uses it anymore. Although, you will see an occasional exception to that statement, usually in the literature of a different field. The primary drawback of the wolff algorithm is its single point nature - the fact that you are picking a starting point and looking for the nearest neighbor at every step. As you can well imagine, any noise in the data might really trash that. The Kantz algorithm and the family of related methods address that by picking a whole bunch of points, tracking them in parallel, averaging their distance to some central point, and watching how that quantity grows. As in the wolff algorithm, you pick a starting point in the trajectory, then you draw an epsilon ball around it, find all the points in the trajectory that are inside that epsilon ball and measure their distances to that central point. Finally, you average all those distances. Now remember, these are all points in a single trajectory. For each one of them, you know where it goes next. If the center point, for example, is the 17th point in the trajectory, you know where the 18th point is, so what you do is you figure out where each of those points goes next...and I'm going to redo this drawing so it's a little less confusing...and this drawing, I'm going to draw in green where each point goes next. This is where the middle one goes next. This one goes here next. This one will go here next. And this one here next, and so on and so forth. This drawing doesn't exactly match the one on the left by the way, it's just a schematic. Then what I do is I compute the distances between the forward image of that center point - which is where the first one went - and the forward images of all those other points. The ratio of the average of all the red distances in this plot to all the average of all the red distances in this plot is a measure of how much the dynamics is stretching out state space around that central point. It's called the stretching factor and that stretching factor is exactly what the Lyapunov exponent is trying to capture. The Kantz algorithm repeats this calculation for a range of different time lengths and a bunch of different initial points. If you plot the results on a log of delta s vs time curve, you'll see a curve like this. Now what's going on here? When time is small, that means that you're not giving things much time to stretch out, so they won't stretch and the ratio will be 1 and the log will be 0. If you let them stretch for a really really long time however, the points in the epsilon ball will spread out all over the attractor and that's as far apart as they can get and that's as far as the volume can grow. That's this upper flat area. In between those asymptotes, there's a scaling region... I apologize for the weakly writing by the way, I downloaded a new copy of the sketchpad software that I use on the tablet and it's making my writing look awful... Let's think about what a diagonal line on a loglinear plot means. What that means is that the natural log of delta s is proportional to time and the constant of proportionality is a. If we take e to the both sides of this, we'll get... where a is this slope. And state space stretches as e to the lambda t, so the slope of the scaling region, a, is the Lyapunov exponent. You saw a similar thing in the previous segment about calculating the fractal dimension by the way, although the axes were different. There, we were after the relationship between the log of the number of boxes of size epsilon to cover something and the log of one over the size of those boxes, so we plotted log of one over epsilon on the horizontal axis, looked for a scaling region on that curve and took its slope as the value of the capacity dimension. In any calculation like these two things that depends on the existence of a scaling region, it's critically important not to impude the existence of one - that is, to see one there just because you want it there when it's really not. And that's a really subtle problem because scaling regions are very hard to define. You could, for example, fit a line to a chunk of your curve and insist that the r-squared value of that line be above 0.9 or something, but that's an arbitrary threshold. All of this brings up an interesting and important point about nonlinear time series analysis. One that's come up before. This is a very very powerful tool, but it has to be used wisely. And what I mean by wisely, is by a person in the loop. Someone needs to look at those plots and see if they have the right shape before deciding to pick off a fractal dimension or a Lyapunov exponent. If it doesn't have the right shape from your perspective, not mine. Or if that shape doesn't persist as you turn the knobs on the algorithm, like the size of that epsilon ball in the Kantz algorithm, you shouldn't use that curve to pick off those values. And that means, among other things, that these kinds of calculations are all but impossible to automate. By that I mean that if you have 1000 data sets and you want to compute the Lyapunov exponent of all 1000 of them, it's not going to work to simply throw those in the hopper of either of these algorithms and go off and have a beer. The answers won't be right. You need to be involved with the process. That's because the problem is hard and subtle, and you now know how to do that. The reason for this is not that these tools are lousy, but rather that nonlinear time series is a very hard problem. Linear tools are easy to use and they always give you an answer, but if you apply them to a nonlinear system, that answer may well be very wrong. Remember the lamp post. Knowing how these algorithms work should also help you understand and respect the requirements for what kind of data are needed for them to work properly. The data set needs to be long enough to cover the behavior that you are trying to sample. It needs to be quickly enough sampled to catch all the details in that, and it needs to not have so much noise as to obscure those details. And those are tough requirements sometimes. The TISEAN toolset, by the way, includes the Kantz algorithm for computing Lyapunov exponents from data. All of that was about calculating lambda 1, the largest positive exponent. There are n-1 other Lyapunov exponents in an n-dimensional system. How to calculate them? That's much harder. Instead of following the spreading a 2-D volume, as sampled by the data, you have to follow higher- dimensional volumes. The algorithms, and the numerics, and the data requirements were hard enough for following the 2-D volumes. The 3-D version occasionally works, giving you lambda 2, but I've never seen results beyond that that I'm really happy with.