Pensions and tax you think to yourself - he’s managed to find the two most boring subjects possible and now he appears to be planning to write about a combination of them! Is this some sort of attempt to win an award for the most boring and tedious blog post every?

Well what if I told you this blog post, whilst being about pension and tax was also about a burning social justice, about high finance, and about death? Would you find that interesting? 

Lifetime Allowance

​Before diving into these deeps social issues, we need to briefly set the scene. Let’s think about all the types of tax that the government levies on us; we get taxed on our income (income tax), we get taxed on our property wealth (council tax), we get taxed on our wealth when we die (inheritance tax), and we get taxed when we buy things (VAT and other sales taxes). Since the government taxes just about every financial transaction we make, they can tweak all these different types of tax to incentivise and disincentivizes certain behaviours. One thing governments have traditionally been very keen on incentivising is getting its citizens saving enough for retirement. The way this has historically been done is to allow income tax relief on the proportion of an employee’s gross salary which is put into a pension. i.e. If you earn £30,000 and put £5,000 of that into a pension, then you pay no tax on the £5,000 and are then taxed as if you only earn £25,000.
​
You might think to yourself great, so now people are saving more for their retirement and everyone lives happily ever after, the end. Well not quite. The problem with this arrangement is certain unscrupulous individuals started to use this mechanism to avoid paying tax. Let’s say you’re 64, quite well off, and hoping to retire next year, why not chuck most of your money into your pension? You can live off your savings for a year and then all that money you would have taken as income and therefore been taxed on can instead be deferred for a year, taken as pension when you retire, and you would pay no income tax on that year's salary. Or let’s suppose an individual come from a very wealthy family, and doesn't really need to get paid a salary to live on, then they could just put their entire salary into a pension, not take anything until they retire (at 55), and then not pay any tax at all on their earnings from their job across their entire career.

Okay, fine, it sounds like we’re going to need some caps on the maximum pension tax relief to stop people doing this. Let’s say – £225k per annum is the maximum that can be put into a pension pot per year, and £1.5m is the total size of a pension pot you can build up before we start making people pay tax. That’s a large amount of money, who’s going to complain about that?

These were the levels set by Tony Blair’s government between 2004 and 2006. These caps were then slowly increased over time, reaching a height of around £250k pa, and £1.8m as a lifetime allowance in 2011. These are very big figures, and it’s hard to imagine anyone objecting to a tax on an annual pension contribution above £250k.

In 2011, with a rapidly expanding budget deficit to deal with and a promise to balance the books, the Conservative government instigated their infamous austerity programme. The tories had pledged to ring-fence certain areas of government spending (the NHS, and Education for example), they therefore started to look around for other areas to raise additional income which had not traditionally been taped. Tax relief for high earners seemed like an obvious target, and who can really say that the rich should not have to pay additional tax given the situation to help try to protect spending on areas such as policing or social care.

By 2014 the annual allowance (the amount that can be saved per annum) had been reduced to £40k per annum, and the lifetime allowance was down to £1m. These values still seem high, however with the economy in the state it was, they are perhaps not as high as they first appear. This reduction came during a period of all-time low interest rates, I'm going to demonstrate which this caused issued by examining a series of graphs.

The following graph displays the Bank of England base rate for the last twenty years, notice the significant reduction in 2007 following the global financial crisis. Note we have been at this very low level ever since.

As a proxy for the investment return a pension or life insurance company could receive on a long term, almost risk-free investment I’m going to use Base rate + 2%, I've plotted this as well on the graph below:

Okay, so why is this important?

The interesting idea is that the ‘price’ of an annuity can be approximated using quite a simple formula. It is (very) roughly equal to:

$ A = \frac{1}{i} $

Where i is our investment return.

What does this mean in real numbers? If I want to retire and get a pension of £10,000 per annum, Then the price of that is going to be 10,000 * A = 10,000 / i. If i=10%, then I will need to have a pot of £100,000. If i = 1%, then I will need a pot of £1m.
​
Let’s plot the implied annuity price on the same chart as a secondary axis:

So we see we have this reciprocal relationship, as interest rates go down, the cost of an annuity goes up. At the moment, due to the very low interest rates, it costs something like £37 to purchase an annuity of £1 per annum for life.

Interesting you might (or might not) think. But this annuity price becomes useful if we combine it with the lifetime tax-free allowance. If we assume a potential retiree is going to use their pension pot to purchase an annuity, then to calculate the level of pension they will be able to get, we need to take their pension pot and divide by the cost of the annuity.

The graph below shows the approximate tax-free pension one could expect if retiring with a pension pot equal to the full Lifetime allowance, using the prevailing annuity rates which are in turn based on the prevailing interest rates.
​

So we see that through a combination of the LTA being reduced significantly and interest rates bottoming out, we are now in the position where the maximum tax-free pension (based on approximate annuity rates) is as little as £25k. A far cry from the approximately £100k when the rule was first implemented.

What about that burning social injustice you mentioned?

​But at least we are all treated the same right? Well not exactly… all of the above glossed over one quite important point, it assumes that the individual in question has a defined contribution pension. For individuals with a final salary pension (which includes MPs, Judges, a majority of CEOs and senior public sector employees), the LTA is calculated with reference to a fixed annuity factor of 20. Given the actual average annuity price is more like 40, this has a massive impact. Let’s plot this on a graph as well.

​

Hmmm, so we see that in 2006 when the LTA was first set up, this factor of 20 was quite generous to Defined Contribution pension holders, they could often retire with a pension of around £120k without paying any tax, whereas a defined benefit pension holder could only retire with a maximum of approx. £80k. This has now flipped, due to changes in annuity rates, and the fact that this factor has not been adjusted for so long. Under current conditions, the maximum pension that can be taken tax free is around £25k pa for a definied contribution pension, but around £50k pa for a definied benefit pension. I told you this blog post wasn't just about pensions and tax, but was really about a burning social injustice (albeit one that only affected the highest earners, so maybe I overplayed this aspect slightly)

So what is the point of all this? What am I advocating?

To me it seems clear that we need (approximate) equality of outcome in terms of taxation for defined contribution pensions compared to defined benefit schemes. This could be achieved in one of two ways, either we apply some sort of floor to how low the tax free pension can reduced to for defined contribution pensions, or we adjust the defined benefit annuity factor  periodically to keep it roughly in line with market conditions and therefore apply a lower maximum tax free rate to defined benefit pensioners. 

 Have you ever tried to simulate a negative binomial random variable in a Spreadsheet?

If the answer to that is ‘nope – I’d just use Igloo/Metarisk/Remetrica’ then consider yourself lucky! Unfortunately not every actuary has access to a decent software package, and for those muddling through in Excel, this is not a particularly easy task.

If on the other hand your answer is ‘nope – I’d use Python/R, welcome to the 21st century’. I’d say great, I like using those programs as well, but sometimes for reasons our of your control, things just have to be done in Excel.
​
This is the situation I found myself in recently, and here is my attempt to solve it:

Attempt 0

The first step I took in attempting to solve the problem was of course to Google it, then cross my fingers and hope that someone else has already solved it and this is just going to be a simple copy and paste. Unfortunately when I did search for VBA code to generate a negative binomial random variable, nothing comes up. In fact, nothing comes up when searching for code to simulate a Poisson random variable in VBA. Hopefully if you've found your way here, looking for this exact thing, then you're in luck, just scroll to the bottom and copy and paste my code.

When I Googled it, there were a few solutions that almost solved the problem; there is a really useful Excel add-in called ‘Real statistics’ which I’ve used a few times:
http://www.real-statistics.com/

It's a free excel add-in, and it does have functionality to simulate negative bimonials. If however you need someone to be able to re-run the Spreadsheet, they also will need to have it installed. In that case, you might as well use Python, and then hard code the frequency numbers. Also, I have had issues with it slowing Excel down considerably, so I decided not to use this in this case.

I realised I’d have to come up with something myself, which ideally would meet the following criteria

• Is free (aren’t all the best things)
• Can be bundled into a Spreadsheet
• Runs in a reasonable amount of time

How hard can that be?

Attempt 1

I’d seen a trick before (from Pietro Parodi’s excellent book ‘Pricing in General Insurance’) that a negative binomial can be thought of as a Poisson distribution with a Gamma distribution as the conjugate prior. See the link below for more details:

https://en.wikipedia.org/wiki/Conjugate_prior#Table_of_conjugate_distributions

Since Excel has a built in Gamma inverse, we have simplified the problem to needing to write our own Poisson inverse. We can then easily generate negative binomials using a two step process:

1. Generate the Gamma sample which we will use as our lambda
2. Generate a Poisson sample, where the Poisson dist has mean equal to our gamma sample

Great, so we’ve reduced our problem to just being able to simulate a Poisson in VBA. Unfortunately there’s still no built in Poisson inverse in Excel (or at least the version I have), so we now need a VBA based method to generate this.
​
There is another trick we can use for this - which is also taken from Pietro Parodi - the waiting time for a Poisson dist is an Exponential Dist. And the CDF of an Exponential dist is simple enough that we can just invert it and come up with a formula for generating an Exponential sample. We then set up a loop and add together exponential values, to arrive at Poisson sample. The code for this is give below:

Function Poisson_Inv(Lambda As Double)
    
    s = 0
    N = 0
    Do While s < 1
        u = Rnd()
        s = s - (Application.WorksheetFunction.Ln(u) / Lambda)
        k = k + 1
    Loop
    BH_Poisson_Inv = (k - 1)

End Function

​The VBA code for our negative binomial is therefore:

Function NegBinOld2(b, r)

    Dim Lambda As Double
    Dim N As Long
    
    u = Rnd()

    Lambda = Application.WorksheetFunction.Gamma_Inv(u, r, b)

    N = Poisson_Inv(Lambda)

    NegBinOld2 = N

End Function

Does this do everything we want?

• It’s free (yay!)
• It can be included as a VBA module of a Spreadsheet
• It does run reasonably quickly

 
There are a couple of downside of though:

• Custom VBA functions referenced directly in the Spreadsheet can be a bit buggy when refreshing - it would be good to have more control over this
• Every time we refresh we will get a different result, and there is no way to set a seed to prevent this.

 
This leads us on to Attempt 2

Attempt 2

If we pass the VBA a random uniform sample, then whenever we hit refresh in the Spreadsheet the random sample will refresh, which will force the Negative Binomial to resample. Without this, sometimes the VBA will function will not reload. i.e. we can use the sample to force a refresh whenever we like. Adapting the code gives the following:
​
​

Function NegBinOld(b, r, Rnd1 As Double)

    Dim Lambda As Double
    Dim N As Long
    
    u = Rnd1

    Lambda = Application.WorksheetFunction.Gamma_Inv(u, r, b)

    N = Poisson_Inv(Lambda)

    NegBinOld = N

End Function

So this solves the refresh problem. What about the random seed problem? Even though we now always get the same lambda for a given rand – and personally I quite like to hardcode these in the Spreadsheet once I’m happy with the model, just to speed things up. We still use the VBA rand function to generate the Poisson, this means everytime we refresh, even when passing it the same rand, we will get a different answer and this answer will be non-replicable. This suggests we should somehow use the first random uniform sample to generate all the others in a deterministic (but still pseudo-random) way.

Attempt 3

The way I implemented this was to the set the seed in VBA to be equal to the uniform random we are passing the function, and then using the VBA random number generator (which works deterministically for a given seed) after that. This gives the following code:
​

Function NegBin(b, r, Rnd1 As Double)
    
    Rnd (-1)
    Randomize (Rnd1)
    
    Dim Lambda As Double
    Dim N As Long
       
    u = Rnd()

    Lambda = Application.WorksheetFunction.Gamma_Inv(u, r, b)

    N = Poisson_Inv(Lambda)

    NegBin = N

End Function

So we seem to have everything we want – a free, quick, solution that can be bundled in a Spreadsheet, which allows other people to rerun without installing any software, and we’ve also eliminated the forced refresh issue. What more could we want?

The only slight issue with the last version of the negative binomial is that our parameters are still specified in terms of ‘b’ and ‘r’. Now what exactly are ‘b’ and ‘r’ and how do we relate them to our sample data? I’m not quite sure.... The next trick is shamelessly taken from a conversation I had with Guy Carp’s chief Actuary about their implementation of severity distributions in MetaRisk.

Attempt 4

​Why can't we reparameterise the distribution using parameters that we find useful, instead of feeling bound by using the standard statistics textbook definition (or even more specifically the list given in the appendix to ‘Loss Models – from data to decisions’, which seems to be somewhat of an industry standard), why can't we redefine all the parameters from all common actuarial distributions using a systematic approach for parameters?

Let's imagine a framework where no matter which specific severity distribution you are looking at, the first parameter contains information about the mean (even better if it is literally scaled to the mean in some way), the second contains information about the shape or volatility, the third contains information about the tail weight, and so on. This makes fitting distributions easier, it makes comparing the goodness of fit of different distributions easier, and it make sense checking our fit much easier.

I took this idea, and tied this in neatly to a method of moments parameterisation, whereby the first value is simply the mean of the distribution, and the second is the variance over the mean. This gives us our final version:
​

Function NegBin(Mean, VarOMean, Rnd1 As Double)
    
    Rnd (-1)
    Randomize (Rnd1)
    
    Dim Lambda As Double
    Dim N As Long
    
    b = VarOMean - 1
    r = Mean / b
    
    u = Rnd()

    Lambda = Application.WorksheetFunction.Gamma_Inv(u, r, b)

    N = Poisson_Inv(Lambda)

    NegBin = N

End Function

Function Poisson_Inv(Lambda As Double)
    
    s = 0
    N = 0
    Do While s < 1
        u = Rnd()
        s = s - (Application.WorksheetFunction.Ln(u) / Lambda)
        k = k + 1
    Loop
    BH_Poisson_Inv = (k - 1)

End Function

I was asked an interesting question a couple of weeks ago when talking through some modelling with a client.
​
We were modelling an airline account, and for various reasons we had decided to base our large loss modelling on a very basic top-down allocation method. We would take a view of the market losses at a few different return periods, and then using a scenario approach,  would allocate losses to our client proportionately. Using this method, the frequency of losses is then scaled down by the % of major policies written, and the severity of losses is scaled down by the average line size.

To give some concrete numbers (which I’ve made up as I probably shouldn’t go into exactly what the client’s numbers were), let's say the company was planning on taking a line on around 10% of the Major Airline Risks, and their average line was around 1%. We came up with a table of return periods for market level losses. The table looked something like following (the actual one was also different to the table below, but not miles off):

​Then applying the 10% hit factor if there is a loss, and the 1% line written, we get the following table of return periods for our client:

Hopefully all quite straightforward so far. As an aside, it is quite interesting to sometimes pare back all the assumptions to come up with something transparent and simple like the above. For airline risks, the largest single policy limit is around USD 2.5bn, so we are saying our worst case scenario is a single full limit loss, and that each year this has around a 1 in 50 chance of occurring. We can then directly translate that into an expected loss, in this case it equates to 50m (i.e. 2.5bn *0.02) of pure loss cost. If we don't think the market is paying this level of premium for this type of risk, then we better have a good reason for why we are writing the policy!

So all of this is interesting (I hope), but what was the original question the client asked me?

We can see from the chart that for the market level the highest return period we have listed is 1 in 50. Clearly this does translate to a much longer return period at the client level, but in the meeting where I was asked the original question, we were just talking about the market level. The client was interested in what the 1 in 200 at the market level was and what was driving this in the modelling.

The way I had structured the model was to use four separate risk sources, each with a Poisson frequency (lambda set to be equal to the relevant return period), and a fixed severity. So what this question translates to is, for small Lambdas $(<<1)$, what is the probability that $n=2$, $n=3$, etc.? And at what return period is the $n=2$ driving the $1$ in $200$?

Let’s start with the definition of the Poisson distribution:

Let $N \sim Poi(\lambda)$, then:

$$P(N=n) = e^{-\lambda} \frac{ \lambda ^ n}{ n !} $$

We are interested in small $\lambda$ – note that for large $\lambda$ we can use a different approach and apply sterling’s approximation instead. Which if you are interested, I’ve written about here:
www.lewiswalsh.net/blog/poisson-distribution-what-is-the-probability-the-distribution-is-equal-to-the-mean

$$ e^{-\lambda} = \sum_{i=0}^{\infty} \frac{\lambda^i}{ i!} = 1  - \lambda + \frac{\lambda^2}{2} + o(\lambda^2) $$
​
We can then examine the pdf of the Poisson distribution using this approximation:

$$P(N=1) =\lambda e^{-\lambda} = \lambda ( 1 – \lambda + \frac{\lambda^2}{2} + o(\lambda^2) ) = \lambda - \lambda^2 +o(\lambda^2)$$

​Every year we have an office sweepstake for the Cheltenham horse racing festival. Like most sweepstakes, this one attempts to remove the skill, allowing everyone to take part without necessarily needing to know much about horse racing.

In case you’re not familiar with a Sweepstake, here’s a simple example of how one based on the World Cup might work:

• Entry costs £5
• Each player randomly selects a team
• The player who picked the team that goes on to win, gets all the money

This allows anyone to enter and have an equal chance of winning (the chance is equal prior to drawing, but not equal once the draw has been made), otherwise everyone would just want Spain or Germany, or whichever team is currently considered the favourite.
 
Note that in order for this to work properly, you need to ensure that the team draw is not carried out until everyone who wants to play has put their money in – otherwise you introduce some optionality and people can then decide whether to enter based on the teams that are still left. i.e. if you know that Germany and Spain have already been picked then there is less value in entering the competition.
 
The rules for our Cheltenham sweepstake were as follows:
 
The Rules
 
The festival comprises 7 races a day for 4 days, for a total of 28 races.
 
The sweepstake costs £20 to enter, and the winnings are calculated as follows:
 
3rd place in competition = 10% of pot
2nd place in competition = 20% of pot
1st place in competition = 70% of pot
 
Each participant picks one horse per race each morning. Points are then calculated using the following scoring system:

• 1st place – 5 x odds
• 2nd place – 3 x odds
• 3rd place – 1 x odds

 
A running total is kept throughout the competition, and the winner is determined after the final race.
 
The odds the scoring are based on are set using the odds printed in the Metro Newspaper on the morning of the races.
 
(as an example, for a horse which has odds of 11/2 in the Metro - if the horse then places 1st, if we selected this horse, we would win (1+11/2)*5 = 32.5 points)
 
Initial thoughts
 
Any set of betting odds can be converted to an implied probability of winning, these would be the odds which over the long run would cause you to breakeven if the race were repeated multiple times with each horse winning a proportion equal to its probability of winning.

Because the scoring in our sweepstake is based on the betting odds, using implied probabilities derived from the odds to help select our horses ends up cancelling itself out (which was the intention when designing the rules). The implied probability can be calculated as one over the odds As an aside, the bookie is indifferent to whether this is the correct probability of winning, they structures the odds purely on the ratio of how their customers are betting. They then structure the odds so that they make money on each race, regardless of which horse wins, for an explanation of this, see my post on creating a Dutchbook:
​
​www.lewiswalsh.net/blog/archives/12-2017
 
Here is an example showing how we would calculate the implied probabilities using some made up odds:
 

​

We can then use the implied probabilities we just calculated to see what would happen if each horse finished in the first three positions. Once we have done this, we can then calculate the Expected Value of betting on each horse:

We see that the payout varies for each horse, but the EV is the same. This is by design, the intention is that it should not matter which horse you bet on, the sweepstake rules equalise everyone.

So what can we - an actuary who knows nothing much about horse racing - do? I don’t have any special knowledge that would allow me to select a horse which would beat the odds listed in the Metro, we appear to be at an impasse.

I could attempt to build a model of which horse will win, and then select my horse based on that, but unless it proves to be more accurate than the Metro odds, I might as well just pick at random. Furthermore, if I could build such a model, then I could just start betting actual money. This probably shows you what a difficult problem this is! There's no such thing as free money. It would be a cool project to try, and it’s something I’ve been meaning to attempt for a while, but that’s best saved for another day.

Attempt 1 - Metro vs pre-race odds

My first thought was that we can exploit the difference in odds between those published in the Metro in the morning, and the latest odds published closer to the start of the race. It seems reasonable to assume that the odds just before the race should be more accurate than the metro odds. There will have been time for additional information to be included in the more up to date odds, e.g. the weather is worse than expected therefore horse x is expected to be slower than usual. Since the payout will be based on the Metro, we will then be able to maximise our EV by exploiting this differential. 

Our table will end up looking something like this:
​

 
We see that we have very small variations in the EVs for some of the horses. It looks like according to this analysis Horse 3 would be the best selection as it has the highest EVs for 1st, 2nd, and 3rd. Based on this strategy, we would then go through each race and select the horse with the highest EV.

Is this what I did? No, for a couple of reasons. The biggest issue was that the Metro did not publish odds in the morning for all races, meaning we couldn’t use the Metro, and therefore the rules of the sweepstake were amended to use the official pre-race odds to calculate the payout instead. This meant there was only one set of odds used, and our edge disappeared!
​
Even if we had used this method, there was a more fundamental issue - the margins we ended up with were tiny anyway. The Metro vs pre-race odds did not swing wildly, meaning that even selecting the horse with the highest EV were only marginally better than picking at random.

So, was there an alternative?

Attempt 2 - 2nd and 3rd place odds
​
My next attempt at an exploitative strategy was based on the insight that the payout multiplier for 2nd and 3rd place was based on the odds of the horse coming 1st, rather than the odds of the horse coming 2nd or 3rd.

The expected value of a horse was not quite as I calculated above, it was actually:
$$EV = P(1)*p_1 + P(2)*p_2 + P(3)*p_3$$
Above, we were using the implied probability of the horse coming first as a proxy for the probability it would come second and third. This is not the same, and some betting websites do allow you to bet on whether a horse will come 2nd or 3rd.

For websites that do not allow you to bet directly on this, then we may still be able to calculate it from the odds for whether a horse finishes in the top 2 or 3 places. We just need to subtract out the implied probability of coming 1st from the probability of coming in the top 2, and then subtracting this out from coming top 3 etc.

​I therefore added some more columns to my table above, corresponding to the probability of the horses coming 2nd and 3rd, and then used this to calculate the EV instead.

We see that the final column, Total EV, now has quite different values for each horse. In this case Horse 15, Seddon has an EV of 11.72. The favourite horse on the other hand - number 7  - only has an EV of 6.2. The intuitive explanation of this is that the probability of Seddon coming first is very low – this is the reflected in the long odds of 67, this then gives us a large multiplier, but the odds of the horse coming second or third are actually relatively less far out – the fact that it is not the favourite actually increases the odds of it coming in a position which is not 1st. But then if does come in 2nd or 3rd, we would still apply the same large multiplier for the odds of it coming 1st. This then gives us our 2nd edge – we can gain a positive EV by focusing.

As a thought experiment, imagine we have a race with three horses – horse A is a clear favourite, horse B is an average horse, and horse C is clearly the weakest. By betting on horse C – the odds of it winning should be very low, so the multiple should be very high, but then this multiple will be applied even if it comes 2nd or 3rd, which is exactly where it is expected to finish.

This therefore suggests our next potential strategy – select horses which maximise our EV using the implied probabilities of the horses coming 2nd, 3rd etc.

So is this what I did? Well kind of....

The issue with this approach is that typically the horses that provide the best EV also have very long odds. In the race analysed above, our horse has an EV of 11.7, but it only has a 7% chance overall of coming in the top 3. In race two for example, the horse with the best EV actually only had a 2.7% chance of coming in the top 3. Since there are only 28 races in total, if each horse we selected only had a 2.7% chance of coming in, then the probability of us getting 0 points overall in the entire competition would then be:
​
$(1-2.7 \%)^{28} = 48 \%$​

So there is roughly a 50% chance we will get 0 points!

Alternatively, if we selected the favourite every time, we could expect it to come top 3 almost every time, and thereby guarantee ourselves points most races, but it also has the lowest EV.

So we have what appears to be a risk vs reward trade off. Pick outsiders and give ourselves the highest EV overall, or pick the favourites thereby reducing our overall EV but also reducing our volatility.

This leads us  neatly to attempt 3 - trying to think about how to maximise our probability of winning the competition rather than simply maximising EV for each race.
​
Attempt 3 - EP curves
​
From the work above, we now have our model of the position each horse will finish in each race – using the 1st, 2nd, and 3rd implied probabilities - and we have the payout for each horse – using the odds of the horse coming 1st. We can then bring our selection of horse and these probabilities together in a summary tab and simulate our daily score stochastically using a Monte Carlo method. To do this we just need to turn the implied probabilities into a CDF, and lookup the value of each position and repeat 10,000 times. The output for this then ends up looking like the following, where the value is the number of points we will for a given race.

So we see that across this sample of 20 simulations, most days we to end up with 0 points overall, but a few days have very high scores. So far so good!

The next step is to set up an EP table of the following, which looks like the following:

The EP table gives us the probability of exceeding various scores in the competition based on our horse selections. In this case, we see that there is a 1 in 20 chance of getting 453 points or greater in the day.

This is useful even on its own – when I was deciding which horses to bet on, I simply played around with the selections until I got an EP table I was comfortable with. My reasoning was quite basic – I decided I wanted to maximise the 1 in 20 value. I wanted to give myself something like a 1/4 chance of winning the whole competition and a 3/4 chance of getting very few points. Since there were four days of races, dividing the 1/4 by another 4 suggests we should be looking at maximising the 1 in 20 level (I admit this reasoning was a bit rough, but it seemed to serve its purpose)

The insight here is that the payout structure of the sweepstake is such that coming in the top 3 is all that matters, and in particular coming 1st is disproportionately rewarded. To see this, we can think of the daily selection of horses as attempting to maximise the EV of the overall prize rather than the EV of our overall score - maximising the EV of each race is only a means to this end. So we are actually interested in maximising the following:

$0.7 * P(1st) + 0.2 + P(2nd) + 0.1 * P(3rd)$

Which will largely be dominated by P(1st), given the $0.7$ factor.
​
This is largely the strategy I went for in the end.

Attempt 4 - Game theory?

I’ve brushed over one difficulty above; in order to maximise our prize EV we need to consider not just which strategy we should take, but how this strategy will fare against the strategies that other people will take. If everyone is maximising their 1 in 20 return period then there’s little point us doing exactly the same. Luckily for me, most people were doing little more than picking horses randomly. We could then formalise this assumption, and come up with a numerical solution to the problem above.

To do this, we would simulate our returns for each day across 10,000 simulations as above, but this time we would compare ourselves against a ‘base strategy’ of random selection of horses and we would simulate this base strategy for the approximately 40 people who entered.

Each simulation would then give us a ranking we would finish in the competition, Here is an example of what that would look like:

​And we could then convert this into an EP table which would look like the following:

So we see that if we select these horses, we end up having somewhere between a 1 in 10 and a 1 in 5 chance of winning the competition.

Now that we have all of this set up, we can then optimise our horse selection to target a particular return period.

I didn’t actually end up setting up the above for the sweepstake, but I suspect it would have been an improvement on my approach

Attempt 5 - Multiple day?

There is a further refinement we can make to the above. We have so far only really been focusing on maximising our chance of winning by using a fixed strategy throughout the competition. But there is no reason we have to do this. After the first day, we should really be including the current scores of each competitor as part of our calculation of our ranking. i.e. if person 1 had a great day and now has 200 points but we had a bad day and still have 0 points, by accounting for this, the model should automatically increase our volatility i.e. start picking horses with longer odds so as to increase the chance of us catching up. If on the other hand, we had a really good first day and are now in the lead, the model should then automatically reduce our volatility and start selecting the favourites more often to help safely maintain our lead.
​
How did it work in practice?

I ended up placing 2nd, and taking 20% of the prize pot which was great! I was behind for most of the competition but then pulled back on the last day when a 66/1 came in 1st place, and I picked up 330 points off a single race. This may look miraculous, but is exactly how the model is supposed to work. 

Does that have any applicability to gambling generally? Unfortunately not, basically all of the work above is based on exploiting the specific scoring system of the sweepstake. There's no real way of apply this to gambling generally.

I came across a useful way of looking at Rate on Lines last week, I was talking to a broker about what return periods to use in a model for various levels of airline market loss (USD250m, USD500m, etc.). The model was intended to be just a very high level, transparent market level model which we could use as a framework to discuss with an underwriter. We were talking through the reasonableness of the assumptions when the broker came out with the following:
 
'Well, you’d pay about 12.5 on line in the retro market at that attachment level, so that’s a 1 in 7 break-even right?'
 
My response was:

'ummmm, come again?'
 
His reasoning was as follows:
 
Suppose the ILW pays $1$ @ $100$% reinstatements, and that it costs $12.5$​% on line.
Then if the layer suffers a loss, the insured will have a net position on the contract of $75$%. This is the $100$% limit which they receive due to the loss, minus the original $12.5$% Premium, minus an additional $12.5$% reinstatement Premium. The reinsurer will now need another $6$ years at $12.5$% RoL $(0.0125 * 6 = 0.75)$ to recover the limit and be at break-even.
 
Here is a breakdown of the cashflow over the seven years for a $10m$ stretch at $12.5$% RoL:
 

So the loss year plus the six clean years, tells us that if a loss occurs once every 7 years, then the contract is at break-even for this level of RoL.

So this is kind of cool - any time we have a RoL for a retro layer, we can immediately convert it to a Return Period for a loss which would trigger the layer.
 
Generalisation 1 – various rates on line
 
We can then generalise this reasoning to apply to a layer with an arbitrary RoL. Using the same reasoning as above, the break-even return period ends up being:
 
$RP= 1 + \frac{(1-2*RoL)}{RoL}$
 
Inverting this gives:
 
$RoL = \frac{1}{(1 + RP)}$
 
So let's say we have an ILW costing $7.5$% on line, the break-even return period is:
 
$1 + \frac{(1-0.15)}{0.075} = 11.3$
 
Or let’s suppose we have a $1$ in $19$ return period, the RoL will be:
 
$0.05 = \frac{1}{(1 + 19)}$
 
Generalisation 2 – other non-proportional layers
 
The formula we derived above was originally intended to apply to ILWs, but it also holds any time we think the loss to the layer, if it occurs, will be a total loss. This might be the case for a cat layer, or a clash layers (layers which have an attachment above the underwriting limit for a single risk), or any layer with a relatively high attachment point compared to the underwriting limit.
 
Adjustments to the formulas
 
There are a few of adjustments we might need to make to these formulas before using them in practice.
 
Firstly, the RoL above has no allowance for profit or expense loading, we can account for this by converting the market RoL to a technical RoL, this is done by simply dividing the RoL by $120-130$% (or any other appropriate profit/expense loading). This has the effect of increasing the number of years before the loss is expected to occur.
 
Alternately, if layer does not have a paid reinstatement, or has a different factor than $100$%, then we would need to amend the multiple we are multiplying the RoL by in the formula above. For example, with nil paid reinstatements, the formula would be:
 
$RP = 1 + \frac{(1-RoL)}{RoL}$
 
Another refinement we might wish to make would be to weaken the total loss assumption.  We would then need to reduce the RoL by an appropriate amount to account for the possibility of partial losses. It’s going to be quite hard to say how much this should be adjusted for – the lower the layer the more it would need to be.

If you have ever generated Random Variables stochastically using a Gaussian Copula, you may have noticed that the correlation of the generated sample ends up being lower than the value of the Covariance matrix of the underlying multivariate Gaussian Distribution. For an explanation of why this happens you can check out a previous post of mine:
www.lewiswalsh.net/blog/correlations-friedrich-gauss-and-copula.

It would be nice if we could amend our method to compensate for this drop. As a quick fix, we can simply run the model a few times and fudge the Covariance input until we get the desired Correlation value. If the model runs quickly, this is quite easy to do, but as soon as the model starts to get bigger and slower, it quickly becomes impractical to run it three of four times just to get the output Correlation we desire.

We can do better than this. The insight we rely on is that for a Gaussian Copula, the Pearson Correlation in the generated sample just depends on the Covariance Value. We can therefore create a precomputed table of Input and Output values, and use this to select the correct input value for the desired output. I wrote some R code to do just that, we compute a table of Pearson's Correlations obtained for various Input Covariance values when using the Gaussian Copula.

a <- library(MASS)
library(psych)

set.seed(100)

m <- 2
n <- 10^6

OutputCor <- 0
InputCor <- 0

for (i in 1:100) {
  
  sigma <- matrix(c(1, i/100,
                    i/100, 1), 
                  nrow=2)
  z <- mvrnorm(n,mu=rep(0, m),Sigma=sigma,empirical=T)
  
  u <- pnorm(z)
  
  OutputCor[i] <- cor(u,method='pearson')[1,2]
  InputCor[i] <- i/10
  
}

OutputCor
InputCor

Here is a sample from the table of results. You can see that the drop is relatively modest, but it does apply consistent across the whole table.

Here is a graph showing the drop in values:
​

Updated Algorithm

We can then use the pre-computed table, interpolating where necessary, to give us a Covariance value for our Multivariate Gaussian Distribution which will generate the desired Pearson Product Moment Correlation Value. So for example, if we would like to generate a sample with a Pearson Product Moment value of $0.5$, according to our table, we would need to use $0.517602$ as an input Covariance. We can test these values using the following code:

a <- library(MASS)
library(psych)

set.seed(100)

m <- 2
n <- 5000000
sigma <- matrix(c(1, 0.517602,
                  0.517602, 1), 
                nrow=2)
z <- mvrnorm(n,mu=rep(0, m),Sigma=sigma,empirical=T)
u <- pnorm(z)

cor(u,method='pearson')

​Analytic Formulas

I tried to find an analytic formula for the Product Moment values obtained in this manner, but I couldn't find anything online, and I also wasn't able to derive one myself. If we could find one, then instead of using the precompued table, we would be able to simply calculate the correct value. While searching, I did come across a number of interesting analytic formulas linking the values of Kendall's Tau, Spearman's Rank, and the input Covariance..

All the formulas below are from Fang, Fang, Kotz (2002) Link to paper: 
www.sciencedirect.com/science/article/pii/S0047259X01920172

The paper gives the following two results, where $\rho$ is the Pearson's Product Moment

We can then use these formulas to extend our method above further to calculate an input Covariance to give any desired Kendall Tau, or Spearman's Rank. I initially thought that they would link the Pearson Product Moment value with Kendall or Spearman's measure, in which case we would still have to use the precomputed table. After testing it I realised that it is actually linking the Covariance to Kendall and Spearman's measures.

Thinking about it, Kendall's Tau, and Spearman's Rank are both invariant to the reverse Guassian transformation when moving from $z$ to $u$ in the algorithm. Therefore the problem of deriving an analytic formula for them is much simpler as one only has to link their values for a multivariate Guassian Distribution. Pearson's however does change, therefore it is a completely different problem and may not even have a closed form solution.

As an example of how to use the above formula, suppose we'd like our generated data to have a Kendall's Tau of $0.4$. First we need to invert the Kendall's Tau formula:

$$ \rho = sin ( \frac{ \tau \pi }{2} ) $$

We then plug in $\rho = 0.4 $ giving:

​We can then test this value with the following R code:

a <- library(MASS)
library(psych)

set.seed(100)

m <- 2
n <- 50000
sigma <- matrix(c(1, 0.587785,
                  0.587785, 1), 
                nrow=2)
z <- mvrnorm(n,mu=rep(0, m),Sigma=sigma,empirical=T)
u <- pnorm(z)
cor(z,method='kendall')

Which we see gives us the value of $\tau$ we want. In this case the difference between the input Covariance $0.587785$, and the value of Kendall's Tau $0.4$ is actually quite significant.