Voronoi Tessellation

Was messing around with these last weekend. It basically generates some random points, builds the tessellation around them, then fills each cell with the average color.



Voronoi Evolution from dizzy pete on Vimeo.


Also made some static images

I signed up to twitter

Because I just wasn't spending enough time on the net as it is. Say hello if you are there too dizzy_pete.

Follow me if you're interested in stuff like cool pictures of skateboards and my thoughts about fellow passengers on the bus.

A volatility filter using historical vol

We have been looking at a way to improve risk adjusted returns by using a volatility filter. Although we could use VIX or equivalent, it turns out that historical volatility will work just as well, if not a little better.

You can see part 1 here Digging into the VIX, and part 2 here What can we use VIX for?

Although the mean return of how we slice things is zero, the distribution of returns is wider for higher readings of our relative measure of volatility. High volatility begets high volatility, at least for our purposes.

By staying out during periods of higher relative volatility, we aim to reduce drawdowns and the volatility of our returns, leading to better risk adjusted results.

A plus of using HV over some external measure like VIX is that it is readily available for any underlying. This means such a filtering technique can be applied to whatever it is we are trading.

Performance

Below is a table with two comparisons, the first compares the HV filter to buy and hold.  Although performance is generally better, we still get some pretty big drawdowns.

The second adds in a 200 period moving average, which is a reasonably strong way of protecting against downside. Again we can see lower volatility and smaller drawdowns with the addition of a vol filter.

I used a 3 month/63 day look back for our relative volatility measure. I haven’t really dug in to what happens when volatility remains elevated for extended periods of time.

I also did not experiment much with the threshold for where we draw the line on ‘high’ relative volatility. I use 0.6 as the cutoff because I originally split things into quintiles when making the first charts.

I also ran this for RUT and NDX over the same period

These results are all frictionless, don’t factor in dividends, return on cash, etc, etc. I don’t consider this viable as a standalone indicator, but something that can be used along side other factors like rotational strategies, or as a potential tool if you are looking for lower volatility.

The source is up here, feel free to have a play around with it and see how you go.

Thanks for reading, 'till next time.

What can we use the VIX for?

In part 1, we took a look at VIX and the relationship it had between historical volatility and realized volatility.

Continuing on, I thought I would take a look at next day returns and the VIX. There is a relationship between SPX and VIX in that when SPX drops, VIX typically rises. This leads the question: does a high VIX serve as a useful indicator of what might happen to SPX the next day?

A problem in answering this is quantifying what we mean by a "high" VIX. A VIX of 25 might be high if it has been under 20 for the last 3 months, but might be considered low if it has been around 30.

To deal with this I use a proportional measure based on its previous trading range. If it is at the top of its range, it gets a 1, and at the bottom it gets a 0. Anywhere in between it gets a number between 0 and 1.

I used data from 2000-2013 and a 252 period lookback, which is roughly 1 year, recording the next day return of SPX. I split this data into 5 groups, the first with relative VIX readings up to 0.2, the second between 0.2 and 0.4, and so on up to 1.0. Then I plotted the next day returns as a boxplot, giving the distribution of returns for a VIX reading in its relative quintile.

The box on the far left corresponds to the lowest relative VIX levels, while the box on the far right corresponds to the highest.

One thing that jumps out to me is that the mean of all these returns is more or less zero. Using this relative measure of VIX has no value as a predictor of next day returns, as whatever measure we use, a mean return of zero will see us end up flat (or in practice, probably down a bit). If there was a real edge, we would expect the mean return to be somewhere other than zero.

The second thing that I notice is the distribution of returns is a lot wider for the right most two boxes, which corresponds to relative VIX measure above 0.60, i.e. relatively higher volatility. Now this makes sense, as we would reasonably expect periods of higher volatility to have higher volatility.

It seems that high readings of VIX can serve as a useful predictor of higher periods of next day volatility, and for someone trying to implement a volatility filter, that could be quite useful.

Conclusion

To recap, we have seen that:

1) VIX is not very good at forecasting 30 day realized volatility.
2) Historical volatility can be a reasonable proxy for IV/VIX.
3) Higher periods of volatility do portend continued higher volatility, at least for the short term (1 day).

Note that 3 does contradict the mean reversion of volatility, because it is looking at the next day returns only. Over a longer period we would certainly expect it mean revert, otherwise we'd be experiencing some very, very bumpy rides ...

Also, I did similar analysis using RUT & RVX and NDX & VXN, both of which had very similar results. We will take a look at them in the next post where we use the findings from this post and the previous to implement a volatility filter.

Code for this post and the previous is up: here

Digging into the VIX

I wanted to revisit using some sort of volatility filter for systematic trading. In particular, if we are trading SPX, can we somehow use the VIX to produce better risk adjust returns? This is not about trading volatility, but more about using additional "out of band" data in our systems.

I thought I would take a look at the VIX, what is it, and can it help us?

In theory, the VIX is the market consensus of what future 30 day realized volatility of SPX will actually be, as derived from option prices of SPX. Anyone who watches or trades SPX has seen that a drop in SPX usually results a rise in VIX, or the IV in SPX option chains.

Occasionally this relationship breaks down, which can be a useful as an indicator in itself, but I had a sneaking suspicion that the VIX would have more in common with historical vol than what realized vol turns out to be over the next 30 calendar days. That is, the VIX is not particularly useful as a longer term forward looking indicator of volatility, at least no more than current historical volatility.

I did a bit of data gathering and calculated realized and historical volatility for SPX. I used 21 periods when calculating vols as the VIX covers 30 calendar days, which is roughly a month, and there are roughly 21 trading days in an average month.

First of all, we can take a look in hindsight at how well the VIX forecast of volatility matched the actual realized volatility over the next 21 days.

In this chart, the red line is what realized volatility turned out to be, going forward from a given day. It is impossible to know this in advance, so we are engaging in some serious data snooping, for educational purposes only. The teal line is the VIX. Below them is a plot of VIX - SPX realized vol, and flagrant abuse of gradients.

You can see that the VIX gives an estimate that is typically higher than what actually eventuated, a well-known phenomenon, which for convenience I will put in "The Volatility Risk Premium" basket.

There is one big divergence around the end of July, 2011, when the VIX was significantly underestimating what realized vol would actually turn out to be.  In August 2011 there was a large selloff and you can see the VIX jumped to a rather ebullient level.

If we moved the red line forward 21 periods, we would get historical volatility, that is, the actually volatility that occurred over the last 21 trading periods.

We can see the relationship is much closer. The VIX still generally reflects a premium over historical volatility, but big moves in the underlying SPX correlate to big moves in VIX/IV.

Relations

The correlations between these three measures, IV, RV and HV are worth taking a look at. This is using data from January 2000 up to about the end of January 2013:

There is a 91% correlation between historical volatility and VIX. This is higher than the correlation between HV and RV, and HV is really just RV shifted back 21 periods. It seems the relation between VIX and HV is something quite strong.

Now because these are all measures of very similar things, we might expect high correlation. Another metric we could use too look at the relationship is R-Squared.

We can see a relatively strong relationship between historical vol and IV/VIX, at 0.60. The relationship between IV and realized vol is looking a lot weaker.

Conclusion

From all of this, I concluded that the VIX is not a particularly useful indicator of 30 day future realized volatility. However, it is a reasonable estimate of what historical volatility was.

As it turns out, this is not particularly useful, as we can easily calculate historical volatility. It does mean that we can flip things around, and use historical volatility as a reasonable proxy for VIX or market expectations of future realized volatility.

This will turn out to be quite useful, as we may find ourselves trading instruments that do not have nice, readily available VIX equivalents, but that is a topic for a subsequent post.

So ends part 1. These posts turned out very long so I decided to split them up, part 2 should be available here.

The R code for this and the next post is up here.

OpenGL & openFrameworks Lighting example

I have been wanting to do more with lighting, but did not really feel like I had a reasonable grasp on How It Worked.

I made this app which has the basic lighting types:

• Ambient Lighting
• Directional Lighting
• Point Lighting
• Spot Lighting

It lets you set the various colours, like the diffuse colour and the specular highlights, and for the material used to colour the sphere you can also set the emissive light.

It is basic but for the basics it is fairly comprehensive, if you aren't quite sure what a diffuse or specular colour is, it will probably be useful. The only thing it doesn't have controls for is light attenuation and I commented out the alpha channel controls as they take up UI space and alpha blending is an app in itself.

There are a couple of oF lighting demo apps but they don't let you easily control the parameters. Mine uses the ofxUI addon so has external dependencies, but should make it a lot easier to see what happens when you shine a blue diffuse light with a red specular highlight on a yellow object.

A screen shot is below. I've made a couple of changes, the biggest is if you are holding down the space bar as you change a color it will change the RGB components together, so you can adjust in terms of black to white.

Source is up here: OpenGL openFrameworks Lighting Demo

These are the instructions from the top of the source:

Spot light defaults to a red light with a blue specular highlight
Point light is similar to a spot light, defaults to green light with white spec highlight
Directional light defaults to blue light with red specular light.
Ambient light defaults to weak grey.
Hold down space when changing a color to adjust black <-> white
Use the 'X Source' toggle to draw the light source
You can also adjust the material properties, shine and the various colours.
Toggles at the bottom of the control panel let you turn lights on/off, and also face culling
It uses ofEasyCam, left click + drag to move around right click + drag to zoom
This app depends on ofxUI and in turn ofXmlSettings

Tracking down errors in R

It's that moment we all know and love, somewhere in our code something has gone wrong. We think we have done everything right, but instead of expected glory we find only terse red text lain below our lintel. 

This can be very frustrating, and trouble shooting these issues can often be very time consuming. 

All is not lost. There are a few bits of R that can greatly help finding out what exactly has gone wrong and where, which in turn should suggest a reasonable course of action.

First we will look at some simple methods we can use to track down issues, namely the warn option and traceback, and then we will look at stepping through functions with debug. 

First Steps

I'm going to use an example that has stuck with me from when I first started using R, using neural nets for classification with the iris data.

Let's take a look at the error:

library(nnet)

X <- iris[,1:4]

Y <- iris[,5]

mod <- nnet(X, Y, size=2)

# weights: 19

Error in nnet.default(X, Y, size = 2) : 

  NA/NaN/Inf in foreign function call (arg 2)

In addition: Warning message:

In nnet.default(X, Y, size = 2) : NAs introduced by coercion

Urgh. Really? What the hell does that mean?

We can use the built in traceback() function to see where this error occured

> traceback()

2: nnet.default(X, Y, size = 2)

1: nnet(X, Y, size = 2)

> 

We can see our call to nnet(), which in turn has called nnet.default() and which is where our error has come from. 

In the error output, we can see there was also a warning "NAs introduced by coercion". As we weren't expecting any warnings, let's track down that, as errors tend to compound. 

Warnings

To find out where that message was coming from, we will use options(warn = 2) which will turn warning messages into errors. We can do this by setting the warn option to a specific level, in this case 2. 

The default is warn = 0, which means warnings will be stored until the top level function returns. We could use warn = 1 which wil print the warning as it is encountered, but in this case we want to stop straight away, so will set it to 2.

options(warn = 2)

Let's try again:

> mod <- nnet(X, Y, size=2)

# weights: 19

Error in nnet.default(X, Y, size = 2) : 

  (converted from warning) NAs introduced by coercion

>

Hmm, still coming from nnet.default, let's see if traceback() is offering any new information.

> traceback()

6: doWithOneRestart(return(expr), restart)

5: withOneRestart(expr, restarts[[1L]])

4: withRestarts({

       .Internal(.signalCondition(simpleWarning(msg, call), msg, 

           call))

       .Internal(.dfltWarn(msg, call))

   }, muffleWarning = function() NULL)

3: .signalSimpleWarning("NAs introduced by coercion", quote(nnet.default(X, 

       Y, size = 2)))

2: nnet.default(X, Y, size = 2)

1: nnet(X, Y, size = 2)

> 

We see a whole bunch of extra stuff in the traceback(), but from 3 onwards it appears to be sideffects of having set warn = 2. We do however see our warning has come from nnet.default again, so we will dig into that to see if we can find out what is going on.

Getting dirty with debug()

To do this, we can use the debug function. We will turn on debugging for nnet, which will let us step through the code line by line as it is executed.

> debug(nnet)

> mod <- nnet(X, Y, size=2)

debugging in: nnet(X, Y, size = 2)

debug: UseMethod("nnet")

Browse[2]> 

The Browse> prompt tells us we are in the debugger. The debug: UseMethod("nnet") tells us the next line of code to be executed is UseMethod("nnet"). We could enter 'n' here to continue to the next line, however a convenient default is just hitting enter (i.e. an empty line). 

Browse[2]> 

debugging in: nnet.default(X, Y, size = 2)

debug: {

    net <- NULL

...

    class(net) <- "nnet"

    net

}

Browse[3]> 

Here R has printed out the R source for the function we just entered. We can see we are still at our Browse> prompt, so will continue on by hitting enter again and again:

Browse[3]> 

debug: net <- NULL

Browse[3]> 

debug: x <- as.matrix(x)

Browse[3]> 

debug: y <- as.matrix(y)

Browse[3]> 

debug: if (any(is.na(x))) stop("missing values in 'x'")

Browse[3]> 

debug: NULL

Browse[3]> 

Well this is neat. Hitting enter at the prompt, R shows us each line that is about to be executed. We will continue on hitting enter until we see our error message.

Browse[3]> 

… 

Browse[3]>

debug: if (length(weights) != ntr || any(weights < 0)) stop("invalid weights vector")

Browse[3]> 

debug: NULL

Browse[3]>

One thing of interest, we can see a conditional if statement is about to be run. When the conditional is evaluated as false, meaning the conditioning code won't be executed, we will see the debug: NULL printed out. 

Browse[3]> 

debug: Z <- as.double(cbind(x, y))

Browse[3]> 

Error in nnet.default(X, Y, size = 2) : 

  (converted from warning) NAs introduced by coercion

> 

Well, there is our warning message. Unfortunately we have lost the Browse> prompt, meaning we are no longer inside the function being debugged, but back at the main prompt. 

This is a side effect of our rather aggressive options("warn") setting. Let's tone it down a bit and set it to 1 so the warnings will be printed as they occur, then jump back into debugging.

> options(warn=1)

> mod <- nnet(X, Y, size=2)

debugging in: nnet(X, Y, size = 2)

debug: UseMethod("nnet")

Browse[2]> 

debugging in: nnet.default(X, Y, size = 2)

debug: {

    net <- NULL

…

Browse[3]> 

…

Browse[3]> 

debug: Z <- as.double(cbind(x, y))

Browse[3]> 

Warning in nnet.default(X, Y, size = 2) : NAs introduced by coercion

debug: storage.mode(weights) <- "double"

Browse[3]> 

After some time, we get our warning message, and we are still in our debugged function. The warning is coming from the line

Z <- as.double(cbind(x, y))

What could be the problem here? Something is going wrong when nnet is converting x and y to doubles. Let's take a look at them to see if there is anything going on.

We can do this by using get("variable"), where variable is the quoted name of the variable. First let's take a look at other two variables in the function to see how it works:

Browse[3]> get("nout")

[1] 1

Browse[3]> get("ntr")

[1] 150

Browse[3]> 

Looking back at the debug output of each line, we can see they were set as some of the dimension info for X and Y, nout has a value of 1 and ntr has a value of 150. 

Let's take a look at x and y now:

Browse[3]> get("x")

       Sepal.Length Sepal.Width Petal.Length Petal.Width

  [1,] 5.1 3.5 1.4 0.2

  [2,] 4.9 3.0 1.4 0.2

  [3,] 4.7 3.2 1.3 0.2

  [4,] 4.6 3.1 1.5 0.2

…

Browse[3]>

This is our input X data.

Browse[3]> get("y")

       [,1]        

  [1,] "setosa"    

  [2,] "setosa"    

  [3,] "setosa"    

  [4,] "setosa"    

...

Browse[3]>

And here is our Y data, just as we passed in.

Recall the warning was triggered by this line:

Z <- as.double(cbind(x, y))

Which is converting x and y to the numeric double type. 

We can see that x is numeric data, while y is character data. What happens if you convert character data to numeric? It doesn't seem to make sense, but let's try:

Browse[3]> as.numeric(get("y"))

Warning: NAs introduced by coercion

  [1] NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA

 [37] NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA

 [73] NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA

[109] NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA NA

[145] NA NA NA NA NA NA

Browse[3]> 

A-ha!

This is just the warning we see, and thinking about it, we can understand why converting strings to numeric data is probably not going to be particularly meaningful.

Let's get out of the debugger and think about what might be going on. Hit 'c' to continue execution 

Browse[3]> c

Error in nnet.default(X, Y, size = 2) : 

  NA/NaN/Inf in foreign function call (arg 2)

>

We are now back at the main prompt. 

Now what??

What is going on? It seems reasonable that we should be able to pass a factor for classification, in fact we are pretty sure that's what we saw being used in the examples in the nnet package documentation.

A careful reading of help(nnet) reveals some details. In particular:

If the response in formula is a factor, an appropriate classification network is constructed ... If the response is not a factor, it is passed on unchanged to nnet.default.

It is possible to pass a factor in Y, but we must use the formula syntax. Looking at the examples, we see the matrix syntax in use as well, however it is transforming the y values using class.ind(). 

Right about now I contemplate calling the police, as I'm pretty sure someone has snuck in and changed the docs while I wasn't looking. There is no way I would make such a simple mistake ...  

Anyway, let's turn off debugging, and see if we can get this working:

> undebug(nnet)

> mod <- nnet(X, class.ind(Y), size=2)

> mod <- nnet(Y~., data=cbind(X, Y), size=2)

Both methods run without error (or warning). Success at last.

Outro

We've seen a few ways we can dig into R and track down where things are going wrong. First is using options(warn=2) to make the R convert warnings to errors, and using traceback() to find out in which function the issue is arising. 

Often, this may be enough to get things back on track, especially if the function causing trouble is small. For more complex issues, we can use debug() which will let us step through the function line by line, inspecting variables and internal state as needed.

It should be said this will only be of use when the source is available to R, i.e. the function in question has been implemented in R and not compiled C/C++/Fortran or whatever that is imported as a shared lib. Deubgging that is more involved and there is some reasonable documentation on this available here.

That's all for now, happy hunting and thanks for stopping by.