Monte Carlo consistency Part II a continuation of Part I

So far we have:

  • MC(D, W, i, k) is the Monte Carlo probability that a $1.00 portfolio with initial withdrawal W, increasing annually with inflation at the rate i, will survive for k years if the annual returns are selected at random from a distribution D.
  • Starting at t = 0 with P portfolios:
    [a]   N = MC(D, Wo, i, n)P is the number of portfolios that survive n years (each starting with initial withdrawal Wo), and
    [b]   M = MC(D, Wo, i, m)P is the number of portfolios that survive m years (each starting with initial withdrawal Wo).
  • The Monte Carlo function MC should satisfy:
      [MC(D, W1, i, n-m)+MC(D, W2, i, n-m)+MC(D, W3, i, n-m)+ ... +MC(D, WM, i, n-m)] / M = N / M.
      where W1, W2, ... WM are the withdrawal rates for the M portfolios that have survived m years.

Since all portfolios, every year, have their returns selected from the same distribution and subject to the same inflation rate i, we'll delete these parameters and let
[0]       MC(D, W, i, k) = f(W,k).

The above condition then becomes:
[1]       [f(W1, n-m)+f(W2, n-m)+f(W3, n-m)+ ... +f(WM, n-m)] / M = N / M.

Now let's order the withdrawal rates (hence the terms in the sum) in order of increasing withdrawal rate and write:
[2]       [f(x1, n-m)+f(x2, n-m)+f(x3, n-m)+ ... +f(xM, n-m)] / M = N / M.
where x1 < x2 < ... < xM.

>No two rates are the same?
Uh, good point. It may be that x1 = x2 = x3 (for example), so let's rewrite [2] assuming that
  N(x1,m) is the number of terms with withdrawal rate x1
  N(x2,m) is the number of terms with withdrawal rate x2
  N(x3,m) is the number of terms with withdrawal rate x3
  etc..

Then we can replace [2] by:
[3]       [N(x1,m)f(x1, n-m)+N(x2,m)f(x2, n-m)+N(x3,m)f(x3, n-m)+ ... ] / M = N / M.

Note that N(x1,m)/M, N(x2,m)/M, ... are the fractions of the M portfolios (the ones that have survived m years) that have withdrawal rate x1, x2, etc.

Indeed, let's divide the entire spectrum of possible withdrawal rates into intervals of width Δx and replace N(x1,m)/M, N(x2,m)/M etc. by the fraction of portfolios that lies in each interval. For example, N(x1,m)/M is the fraction that lies in (x1,x1+Δx).

That is, we introduce a distribution of withdrawal rates (for those portfolios that have survived m years), namely p(x,m), where
[4]       p(x,m)Δx is the fraction of m-year survivors with withdrawal rates in the interval (x, x+Δx).

Using this withdrawal rate distribution function, we (again!) rewrite [3] as:
[5]       [p(x1,m)f(x1, n-m)+p(x2,m)f(x2, n-m)+p(x3,m)f(x3, n-m)+ ... ]Δx = N / M.

However, using [a], [b] and [0] above, we note that N / M = f(Wo, n) / f(Wo, m) = (n-year survival rate) / (m-year survival rate).

Then [5] becomes:
[6]       [p(x1,m)f(x1, n-m)+p(x2,m)f(x2, n-m)+p(x3,m)f(x3, n-m)+ ... ]Δx = f(Wo, n) / f(Wo, m)

>zzzZZZ
Wait! We're almost there!

We rewrite (for the last time!) [6] by letting the number of terms become infinite and Δx0 and let m = n-1.
We get:
[B]       p(x,n-1)f(x,1)dx = f(Wo, n) / f(Wo, n-1)

>Huh? The number of terms become infinite?
Remember that we initially started with P portfolios. To do a proper Monte Carlo we should consider an infinite number of starting portfolios. The number that have survived m years would also be infinite as well as the number surviving n years and ...

>What if none survive 40 years?
Are you saying that there's a zero probability that there's a 40-year survivor? Remember, the range of annual returns is infinite. There's a finite possibility that the returns are all HUGE. The lucky investor that gets these will most certainly survive 40 years, even if her initial withdrawal rate is huge.

>Not if the intial withdrawal rate is 100%!
Good point, but I'm assuming that the first withdrawal occurs at the end of year 1 and the withdrawal is a fraction of the starting portfolio. After one year, there's a possibility that that portfolio might have doubled or tripled or ...

>I doubt if there are many of those.
Precisely! The withdrawal probability distribution we introduced, namely p(x,m), would be small for x very small ... and x very small means a very small withdrawal rate and that means a HUGE portfolio. Similarly, at year m, there are some tiny portfolios (soon to become negative!) so the withdrawal rate, namely (withdrawal)/(tiny portfolio), would be large. We would expect a few portfolios would be at both ends of this withdrawal spectrum and ...
>And that picture is accurate?
No. It's a guess. To get the real picture we'd start a jillion portfolios with an initial withdrawal rate Wo and follow each for m years ... and see what fraction have (at year m) a withdrawal rate of Δx, getting p(Δx,m), then what fraction have a withdrawal rate of 2Δx then 3Δx etc. etc. and plot p(x,m) versus x.
>And have you done that?
Uh ... no. Life is far too short for that.
>So now what?
Our magic equation [B] gives the n-year survival rate f(Wo, n) in terms of the (n-1)-year survival rate f(Wo, n-1).
We just multiply by some curious factor: p(x,n-1)f(x,1)dx which is essentially the sum of terms like:
    (probability that portfolio has lasted n-1 years)*(probability that it lasts one more year).

Figure 1

>Do you know what you're doing?
That remains to be seen, but it's fun, eh?
>For whom?

Okay, so we do a bit of Monte Carlo:
Starting with 10,000 portfolios withdrawing (initially) at 4%, increasing with 3% inflation, we follow them for 40 years and look at the 10,000 withdrawal rates at year 10, 20, 30 and 40 and plot the distribution ... and get Figure 2 (which is much like Figure 1, eh?)

>A lucky guess. Besides, what returns did you use?
Oh, I forgot. A lognormal distribution with Mean = 10%, SD = 20%.

>And those red lines at x = 0?
They're portfolios that didn't survive.

>It looks like the average withdrawal rate keeps getting smaller, eh?
Sure. The withdrawal rate is (withdrawal amount) / (portfolio value) and although the withdrawal amount is the same for all portfolios, the average portfolio value is going up.

>So Figure 2 shows the distribution of (withdrawal amount) / (portfolio value)?
Yes. It shows the fraction of portfolios (at year m) which have withdrawal rate
x = (withdrawal amount) / (portfolio value)   for 0 < x < infinity.


Figure 2

So, at year 30 for example, (withdrawal amount) = 0.04(1.03)30 = 0.097 or 9.7 cents (assuming an initial withdrawal rate of 4% and 3% inflation and an initial portfolio of $1.00). So the m = 30 distribution is the distribution of the quantity 0.097 / (portfolio values at year 30). As it turns out (for the example we're considering), the median portfolio has a withdrawal rate of 1.2% (as is suggested from the m = 30 chart in Figure 2). That means that half the withdrawal rates are less than 0.012 (that's 1.2%) so half the portfolios are greater than 0.097/0.012 = 8.08 or $8.08.

>And half are less than $8.08 ... including those that have dropped dead?
Yes, including those that haven't survived 30 years - 1912 of them (or 19.12%) in this 10,000 Monte Carlo run, as indicated by the red line at x = 0.
Since 1912 haven't survived to 30 years, then 8088 portfolios have survived and, of those, 7208 survive the addtional 10 years.
>That's a 72.08% 40-year survival rate?
Yes, 72.08% of the original, starting portfolios ... for this particular Monte Carlo run.
That suggests that, of the 8088 30-year survivors, 7208 survive to year 40 and 7208/8088 = 0.89 so 89% of these survive to year 40.
>Is that what your formulas give?
In order to answer that question I'd have to run a jillion MC simulations for each of the 8088 30-year survivors.
>I'll wait ...

Let's look at [3], with m = 30 and n = 40:
      [N(x1,30)f(x1,10)+N(x2,30)f(x2,10)+N(x3,30)f(x3,10)+ ... ] / M = N / M.
where

  • N(x,30) gives the number of 30-year survivors which (at year 30) have a withdrawal rate x (as in Figure 3a below - derived from 2c above)
  • f(x,10) gives the fraction of portfolios that survive10 years if the starting withdrawal rate is x (as given in Figure 3b below)
  • M = 8088 = number of 30-year survivors
  • N = 7208 = number of 40-year survivors.


Figure 3
The product N(x,30)f(x,10) looks like Figure 3c (where only the 30-year survivors are considered):

The area under the graph in Figure 3c is represented by the sum:
N(x1,30)f(x1,10)+N(x2,30)f(x2,10)+N(x3,30)f(x3,10)+ ...

... or, for an infinite number of Monte Carlo iterations (and a smooooth curve), it's represented by p(x,n-1)f(x,1)dx.
>And the vertical axis measures the number of ...

Sorry. We'd have to change to the fraction of rather than the number of because nobuddy wants infinite numbers.
>So ... what now?
I'm thinking, however let's look at the term f(x,1) which gives the survival rate for portfolios that start with a withdrawal rate x and run for one year.

After one year (with a return of R), a $1.00 portfolio is worth (1+R) after which we withdraw $x, increased by one year's inflation.
That gives (1+R) - x(1+i) which will be negative (meaning the portfolio does NOT survive) if (1+R) < x(1+i).
But we're assuming a lognormal distribution of returns, R, so 1+R = ez for some normally distributed variable z.
Hence we can write ez < x(1+i) or z < log[x(1+i)] and since x(1+i) is something like 0.04 (for a 4% withdrawal),
then log[x(1+i)] will be something like log[0.04] = -3.2 ... so z is going to be negative in order to satisfy this inequality!

However, if R is distributed lognormally with Mean M and Standard Deviation S, then z will have a normal distribution with Mean
    M = (1/2)log[(1+M)2/(1 + S2/(1+M)2)]
and Standard Deviation
    S = SQRT[log{1+S2 / (1+M)2]     as noted here.

If we let N(z) represent the "standard" Normal cumulative distribution (with Mean = 0 and SD = 1) as in Figure 4, then our z distribution
(with Mean = M and SD = S) will be N[(z-M)/S] and we're looking at the probability that our portfolio will NOT survive
... meaning that z < log[x(1+i)] ... and that probability is:   N[{log[x(1+i)]-M}/S].
The probability that our portfolio WILL survive is then 1 - N[{log[x(1+i)]-M}/S] ... and that's f(x,1) ... and that might look like Figure 5, below (similar Figure 3b, above):

Figure 5

Figure 4

>What! A 155% withdrawal rate?
Remember that we're considering withdrawing after one year and it's possible (though not probable) that we get a HUGE return that first year and we're withdrawing a percentage of the initial portfolio so ...
>Why not withdraw immediately? Why wait a year?
We could do that. It'd just mean (with a 4% withdrawal) we start with a reduced portfolio - for example $0.96, instead of $1.00.
>What now?
I'm thinking ...