Lévy Distribution for Stock Price Simulation

[Hans Brinker]

Dear Steve,

I found this site searching for a Monte Carlo Simulator to use for my investment modeling. I am glad I found the one you created, which is using the Normal distribution. I already used it in some of my models, but based on the findings of Benoît Mandelbrot, I would like to use the Lévy distribution, since that much better reflects stock market price changes.
My real question is whether anyone knows whether someone already made a program like yours using the Lévy distribution available over the internet?
Many quants / financial engineers made small programs modeling stock market prices using Lévy distribution based random walks for the companies they work (and worked) for.
I want to use this privately.

I knew that stock prices were not quite normally distributed, but did not know the extend of the errors until I happened to read some examples in the Jan 24 2009 Economist (page 14 of "A Special Report on the Future of Finance").
For example, the Dow Jones should have moved by more than 3.4% on 58 days between 1916 and 2003 based on a normal distribution but actually did so 1001 times.
It should have moved by more than 4.5% on 6 days, it did so on 366.
It should have moved by more than 7% only once in 300,000 years but did so 48 times in the 20th century alone.
For the Lévy distribution see: http://en.wikipedia.org/wiki/Modeling_and_analysis_of_financial_markets

Thanks for your help,

Hans

[SteveW]

A very interesting question, which is taking some time to look into. Here are some links that will allow tracing the paths I followed, and you'll probably get ahead of where I am in the process.

I started by reading the link you posted at Wikipedia, and followed to some pages about the Levy distribution, including http://en.wikipedia.org/wiki/L%C3%A9vy_distribution. I'd recently encountered this before in various contexts but hadn't read everything carefully.

I'd not realized previously (correct me if I'm wrong) that the standard normal distribution actually is a Levy distribution, with alpha = 2. The distribution Mandelbrot suggests is a Levy with alpha = 1.7, so it's an adjustment of a parameter rather than a tremendously different paradigm. At http://en.wikipedia.org/wiki/Stable_distributions is a different-looking graph of a Levy ("symmetric centered") where it is easy to see the difference between alpha = 1.5 vs 2.0. True to what I've seen of Mandelbrot's hypothesis, the 1.5 is slightly more pointy and has more of its area in the tails, i.e. greater leptokurtosis. 1.7 would be in between.

Next was a search at Yahoo on levy distribution algorithm, which turned up many more links than I've been able to follow. Early ones are promising... It is clear that people have created random generators for Levy.

Some potentially useful links. I haven't looked at all of these thoroughly.

--General info and method: http://www.sjsu.edu/faculty/watkins/stablerand.htm
--It looks as though the GNU Scientific Library (GSL) has open source C code for generating random variables with a Levy distribution (and many others, as well), based on the info at http://www.gnu.org/software/gsl/manual/gsl-ref.html#The-Levy-alpha_002dStable-Distributions. Even if you don't wind up using their functions as-is, it may greatly help clarify how to do it.

From this page http://www.gnu.org/software/gsl/ I selected a mirror and ended up at the download page (for GSL 1.12) at http://cudlug.cudenver.edu/GNU/gnu/gsl/. There's a brief description of the GSL at http://directory.fsf.org/project/gsl/. I've thought it would be good to have this library for other reasons recently, so I'll eventually download it and take a look.

However, as I mentioned, I barely got into the search results before following that path, so there may be algorithms posted online that you can find without going for the entire GSL.

That's as far as I've gotten, but I think it should certainly be possible to find the way to do this.

It's getting ahead of things, but since the Levy distribution is being substituted for a regular normal distribution, and simulated stock prices are better generated using a lognormal distribution, it makes sense to me that once the Levy method is found, a lognormal variation of it would be more appropriate for the final calculations.

[Hans Brinker]

Dear Steve,

Thank you for looking into this. Prior to my original posting, I already browsed through much of what wikipedia has to offer and also found and looked through Watkins' site, but missed his CMS algorithm page. It looks like this is implementable in Excel (I used some Basic and Fortran 30 years ago and more recently Visual Basic and C#, but am not really capable in using those languages).
I am not sure what the alpha would be best reflecting the stock market prices, but doubt whether the alpha for the best reflection for stock prices and cotton would be the same.
I pretty much follow Robert Shiller's approach in believing that stock (and other) market prices are very much effected by human emotion and their mass manipulation (read "Irrational Exuberance"). Following that, the emotional reasons for buying and selling cotton may be very different for buying stock. In addition, these emotions and their effects on markets may change over time and with technology. Hence it very well may be that the skewness and kurtosis are different for different products (stocks, cotton, art), markets (Japan, where most large companies hold large fractions of shares of other companies, versus the West where this is not the case) and over time (newspapers in late 1800, radio news in early 1900, TV in 1950+ to online trading now!).
One reported example (I think by Prof. Krugman) was that never before world-wide production was cut as fast as during this recession. This probably is caused by all the linked supply chains, measuring sales at one end and almost instantly adjusting production at the other. This on itself effects salary and revenue streams as well as human emotions, both likely effecting stock prices.
Once modeled data using a log-normal Lévy distribution is available, it would be interesting to see what the best parameter settings are for skewness and kurtosis and whether it is prudent to assume that they are constant over time.

I will follow up here with relevant new info or results and invite anyone else to join in.

Regards,

Hans

[SteveW]

Hans,

I downloaded the GNU Scientific Library, converted two relevant functions to C++, and added them to mylib.cpp. They are at the bottom of the page at http://25yearsofprogramming.com/msvcpp/mylibcpp.htm. If you have some experience with C#, you'll probably find the code comprehensible. It might be possible to create an equivalent calculation method in Excel. I don't find the algorithm comprehensible, but I'm not worrying about that for now.

The Levy function returns a random double from the Levy distribution. (The exponential random function is needed by the Levy function.)

The GSL "Bugs" file has a report from one user who believes that the Levy functions have some glitches or inaccuracies. I'm in no position to judge, and, if it's true, in even less of a position to attempt troubleshooting, but thought I should mention it. I generated some sets of random numbers, transferred them to Excel, and graphed them with Tools > Data Analysis > Histogram. With alpha=2, it does indeed look like a normal curve, and other alpha values seem to produce graphs having the shapes they should, and with long tails.

This is only a step in the right direction. The next step is to use the random numbers to calculate consecutive daily price changes, hopefully using a value of "c" that correctly produces the Levy equivalent of "standard deviation" of price fluctuations. I'm experimenting with how to do that. That is, just having any old Levy distribution isn't enough. You need to be able to simulate the price pattern of a stock that has, say, a 40% volatility. The resulting price series must have the desired 40% volatility when calculated the usual way, except that the distribution is Levy and not lognormal.

After that, I plan to translate the C++ code to JavaScript in hopes of being able to create 4 calculation options on the Monte Carlo page: arithmetic normal, lognormal, arithmetic Levy, and log-Levy, all with user-selectable values of alpha and c (standard deviation). That code will also be GPL open source, so, if successful, there should eventually be an open source library of these methods for people to use for experimentation.

----

If you're not familiar with the GNU General Public License (what the GSL is licensed under, and what I also use), a brief summary is: as long as what you're doing is for personal use, you can copy the code and do whatever you want with it. If you redistribute it (as I am doing by publishing a modified version of the code), you must retain all copyright notices (and add yours, with notations of what changes you've made) and also license it under GPL3, meaning that the entire package it is part of must be open source.

----

With regard to alpha, what little reading I've done on it from time to time has left me with the (vague) impression that Mandelbrot has been fairly consistent about the 1.7 figure, and not necessarily just for cotton. I could speculate that he might say it results more from collective human behavioral characteristics (aggregate behavior of many people) than from any characteristics of the asset being traded. That's one guess, anyway.

It would be interesting to see what the best parameter settings are for skewness and kurtosis and whether it is prudent to assume that they are constant over time.

The kurtosis is alpha. Skewness is a measure called beta, and would be the trend. This function doesn't allow for beta. There is a related GSL function that does, but I figure it is safest to assume beta is always 0 (no trend).

[SteveW]

All done. I added new calculation options to the page at http://25yearsofprogramming.com/blog/20070412c-montecarlostockprices.htm.

[Hans Brinker]

Thanks Steve,

You did a lot of work!
It will take me a while to go through it and do some checking. Following that I will do extensive testing.
I will come back as soon as I have something worthwhile to say.

Regards,

Hans
 

[Hans Brinker]

Steve,

Finally I got so far to look at the data generated by the Lévy distribution based stock price data generator which you so generously created.
I will look into the distribution aspects later, but so far I see a discrepancy between the volatility % calculated using the same date set in Excel and the number displayed at the bottom of the "Input" sheet.
The Excel calculated numbers are much smaller than the numbers displayed at the bottom of your sheet.
I did the calcs for 2500, 250 and 10 days (latter results displayed below) and the differences appear to be a function of the number days for which data are requested.

I calculate the % volatility as STDEV() / Average().
Am I interpreting any of your data or the definition of volatility incorrectly?


Regards,

Hans

*****************
Inputs Instructions
To calculate by the most common method, simply enter Starting Price and Volatility, and click Calculate.
 Output Prices

Starting Price:
 This will be the starting point from which successive values are calculated, but the start value will not be in the output list. Thus, if you want to extend an existing time series, you can enter its last value as the start value here. 
102.02
102.16
103.41
102.47
101.21
101.57
102.99
103.05
102.54
103.19
 
Volatility (%):
% The same volatility figure required for Black-Scholes option valuation. Enter as a percent, not a decimal. Ex: 62, not 0.62. If you have a set of daily prices for a real stock but don't know its volatility, you can use the volatility calculator. Historical price data is available at Yahoo! Finance. 
Price calculation method:

 Arithmetic
 Logarithmic
 None
 Arithmetic: Up/down percentage point moves are equally likely even though up and down are not symmetrical. For example, a stock that loses 50% must gain 100% to regain the lost ground.

Logarithmic: Up and down moves are symmetrical; the probabilities of the stock doubling or halving are equal. Although it is nonintuitive, the most common model used in economics and finance is Logarithmic + Normal (below) = "lognormal".

None: Uses the calculator as a random number generator. The raw random numbers are not turned into simulated daily prices. 
 
Probability distribution:

 Normal
 Levy
 Uniform
 Normal (aka gaussian): Mean=0, StdDev=Volatility / 100. 

Levy: Mean=0, Beta (skewness)=0, Alpha (kurtosis) is used as entered below. c (variability, "StdDev") = (Volatility / 100) / Vol.adj.factor.

Uniform: Mean=0. For raw random numbers, > -Volatility / 100 and < Volatility / 100. For simulated prices, equally distributed between ±1.733421 * Volatility / 100, which forces calculated volatility into desired range9.
 
For Levy only:
Alpha:


Vol. adj. factor:

 Alpha must be > 0 and <= 2. Default value of 1.7 is probably too low. Alpha = 1 is a Cauchy distribution. Alpha = 2 is a Gaussian normal distribution. 

Volatility adjustment (aka fudge) factor is any number that brings the calculated volatility (in Verify:, below) closer to the requested value. With Alpha = 1.7, an adjustment factor of 2 seems to work ok.
 
# of days:
 The number of values to calculate. Default = 252 = 1 trading year. 
 Click Calculate. 
(Reverse) The output default sort order is from oldest date at the top to most recent at the bottom. To reverse, click the Reverse button before copying the data. 
Copy
and
Paste Click anywhere in the output box. Type Ctrl+A (Select All), then Ctrl+C (Copy to clipboard). Paste the numbers in the application where you want to use them (Ctrl+V).   
14.21824689% Verify: at left is the actual calculated annualized volatility of the generated numbers. It won't exactly match the requested value because the output prices are random numbers.   
 

[SteveW]

Hans,

The most concise definition of financial volatility (as required for calculations such as Black-Scholes option valuation) is: "the standard deviation of the stock's daily logarithmic returns, annualized".

The standard deviation figure is actually the standard deviation of an exponent, so it isn't in any sense a percentage or fraction of the stock price, which makes it a nonintuitive statistic. The "annualization" calculation that is applied to the standard deviation is also a nonintuitive one.

The methods for the two calculations are described on my Volatility Calculator page at http://25yearsofprogramming.com/blog/20070412a-volatility.htm, which has links to Wikipedia articles where I found the methods and equations.

The way you'd get the volatility calculation in Excel would be like this:

Make a column containing all the daily closing prices. Let's assume it's Column A. The prices must be daily.

In Column B, starting in B2, enter the formula: =LN(A2/A1). (You can't do this in B1 because it would need a reference to A0, which doesn't exist.) Be sure to use LN() and not LOG().

Drag B2 down the column so its formula is used in all the cells of Column B. It's the natural logarithm of (today's price / yesterday's price), for each day.

Now in C1 (or any cell), enter the formula =STDEV(B2:Bxx) where xx is the last cell of column B that has a value. This is the standard deviation of all the values in Column B.

Now the annualization:

In D1 (or any cell), enter the formula =C1 * SQRT(252). This is the annualized standard deviation, which is the Volatility. Except for rounding errors, this should match the Black-Scholes value on my Volatility Calculator web page, which matches the "Verify" figure on the Monte Carlo simulator page.

[Hans Brinker]

Steve,

Yes, doing that I get the same number, at least for the first four significant numbers.
Thanks,

Hans