Cellular Automaton Models for Market Prices

Dr Michael Kelly

It has long been established by Brealey [1] and Cootner [2] amongst others that stock prices exhibit random behaviour which satisfy the Markov property of depending only upon their last values and not their previous history. This property is also consistent with the weak version of the Efficient Market Hypothesis, which asserts that the present price of a stock impounds all the information contained in its record of past prices. These two properties are consistent with models of stocks described by one dimensional cellular automata since they also depend upon rules which describe the next state of a process in terms of the previous state without relying on past information contained in earlier times. Furthermore the Efficient Market Hypothesis implies that all the information contained in the transactions at any time can be used to construct a decision rule that will encapsulate the local behaviour of all the particpants in the market at that time in such a way that it completely determines the market state at the next moment. This pattern of behaviour has already been observed by Wolfram [3] in his book A New Kind of Science where he posits rule 90 as a possible model for stock prices. 

In this paper we shall examine a number of one dimensional CA and compare them with popular stock models which are contained in the financial literature. We shall attempt to discern which properties of stock behaviour are best described by simple functions of CA output.

Initialization Cells

Needs["NonlinearRegression`"];Needs["MultivariateStatistics`"];
Needs["BarCharts`"];Needs["Histograms`"];

NormalDistributionFit[data_List] := Module[{mu, sigma,q,quantiles, cdfs, paramrules},quantiles = Table[Quantile[data,q],{q,.01,.99,.01}];
cdfs = Transpose[{Range[.01,.99,.01], quantiles}];
paramrules = BestFitParameters/.NonlinearRegress[cdfs,mu+(Sqrt[2]) sigma InverseErf[0,-1+2 q],{mu, sigma},{q},RegressionReport->BestFitParameters];
NormalDistribution[mu,sigma]/.paramrules]

LogNormalDistributionFit[data_List] := Module[{mu, sigma,q,quantiles, cdfs, paramrules},quantiles = Table[Quantile[Log[data],q],{q,.01,.99,.01}];
cdfs = Transpose[{Range[.01,.99,.01], quantiles}];
paramrules = BestFitParameters/.NonlinearRegress[cdfs,mu+(Sqrt[2]) sigma InverseErf[0,-1+2 q],{mu, sigma},{q},RegressionReport->BestFitParameters];
LogNormalDistribution[mu,sigma]/.paramrules]

NoncentralStudentTDistributionFit[data_List] := Module[{n,lambda, paramrules},
paramrules = FindRoot[{((lambda (Sqrt[n]) Gamma[(1/2) (-1+n)])/((Sqrt[2]) Gamma[(n/2)])) == Mean[data],(((1+(lambda^2)) n)/(-2+n))-(((lambda^2) n ((Gamma[(1/2) (-1+n)])^2))/(2 ((Gamma[(n/2)])^2))) == Variance[data]},{{n,Length[data]},{lambda,Mean[data]}}];
Re /@ NoncentralStudentTDistribution[n,lambda]/.paramrules]

Cellular Automata models for equities

We will first develop a basic model for stock prices along the lines first introduced by Wolfram in his book NKS. Here the emphasis will be on finding the correspondence between code and the intuitive idea of local market behaviour as expressed by simple transformation rules. Because CA are inherently stochastic and are even used to evaluate random values of probability distributions, then it follows that the traditional stochastic models of finance must be capable of being characterised in terms of CA. 
The Wolfram model for stock behaviour is based upon cellular automaton 90, which is shown below.

ArrayPlot[CellularAutomaton[90,RandomInteger[1,101],100]]

When expressed as an explicit rule for r=1 neighbourhoods this is expressed as:

ArrayPlot[CellularAutomaton[{{1,_,1}->0,{1,_,0}->1,{0,_,1}->1,{0,_,0}->0},RandomInteger[1,101],100]]

The Wolfram model for stock behaviour is illustrated below. It has the fractal pattern of stocks but exhibits none of the other properties such as unimodality, fat tails and negligible skewness.

Module[{hist=CellularAutomaton[{{1,_,1}->0,{1,_,0}->1,{0,_,1}->1,{0,_,0}->0},RandomInteger[1,101],1000],diffs,cumpoints,min1},diffs=Apply[Subtract,Outer[Count,hist,{1,0},1],{1}];cumpoints=Accumulate[diffs];ListPlot[cumpoints,Joined->True]]

Plot of price differences

The following is a first modification to at least keep the stock path positive, because S=0 is an absorbing state. Essentially we keep the shape but shift the prices vertically.

Module[{hist=CellularAutomaton[{{1,_,1}->0,{1,_,0}->1,{0,_,1}->1,{0,_,0}->0},RandomInteger[1,101],1000],diffs,cumpoints,min1},diffs=Apply[Subtract,Outer[Count,hist,{1,0},1],{1}];cumpoints=Accumulate[diffs];min1=Min[cumpoints];ListPlot[If[min1≤0,1-min1,0]+cumpoints,Joined->True]]

What is the problem with this plot? It displays an instrument that experiences enormous changes that is inconsistent with real market behaviour. Some equities may experience such changes over large periods of time, but most do not. So this basic model needs to be modified and adapted.

Module[{hist=CellularAutomaton[{{1,_,1}->0,{1,_,0}->1,{0,_,1}->1,{0,_,0}->0},RandomInteger[1,101],1000],diffs,cumpoints,min1},diffs=Apply[Subtract,Outer[Count,hist,{1,0},1],{1}];cumpoints=Accumulate[diffs];min1=Min[cumpoints];Histogram[If[min1≤0,1-min1,0]+cumpoints]]

Histogram of price differences

The difficulty with this simple model is further manifest in its multi-modal behaviour. A phenomena that is never observed. So this basic model has to be adapted. The problem appears to be that while this simple model has the apparent random structure of some financial instrument it needs to be standardised so that it's range of values are realistic.

shiftedData=Module[{hist=CellularAutomaton[{{1,_,1}->0,{1,_,0}->1,{0,_,1}->1,{0,_,0}->0},RandomInteger[1,101],100],diffs,cumpoints,min1},diffs=Apply[Subtract,Outer[Count,hist,{1,0},1],{1}];cumpoints=Accumulate[diffs];min1=Min[cumpoints];If[min1≤0,1-min1,0]+cumpoints]

LogNormalDistributionFit[shiftedData]

Histogram[RandomReal[%,1000]]

Sometimes the T distribution fit fails because there is no solution to the number of degrees of freedom consistent with the data.

NoncentralStudentTDistributionFit[shiftedData]

Histogram[RandomReal[%,1000]]

The problem could be with the rule itself. If we change just one of the rules, the last one, then the behaviour of the model changes radically so that we recover a unimodular distribution consistent with observed data.

Module[{hist=CellularAutomaton[{{1,_,1}->1,{1,_,0}->0,{0,_,1}->0,{0,_,0}->1},RandomInteger[1,101],1000],diffs,cumpoints,min1},diffs=Apply[Subtract,Outer[Count,hist,{1,0},1],{1}];cumpoints=Accumulate[diffs];min1=Min[cumpoints];Histogram[If[min1≤0,1-min1,0]+cumpoints]]

Histogram of alternate prices, obtained from a different rule, 95.

The problem with the model could be with the limitation of the size of the neighbourhood. If we increase the neigbourhood to allow for one more agent then it again exhibits the right behaviour and this is consistent with actual experience as traders and other agents do react in groups larger than 3. But we again encounter the same problem of very large deviations.

Module[{hist=CellularAutomaton[{{1,1,_,1}->0,{0,1,_,1}->1,{1,1,_,0}->1,{0,1,_,0}->0,{1,0,_,1}->1,{1,0,_,0}->1,{0,0,_,1}->0,{0,0,_,0}->0},RandomInteger[1,101],1000],diffs,cumpoints,min1},diffs=Apply[Subtract,Outer[Count,hist,{1,0},1],{1}];cumpoints=Accumulate[diffs];min1=Min[cumpoints];ListPlot[If[min1≤0,1-min1,0]+cumpoints,Joined->True]]

Plot of prices, with an alternate 4-neighbour rule.

The difficulty even with this model is apparent in that if we analyse its underlying probability distribution then it depicts only markets with enormous skew.

shiftedData=Module[{hist=CellularAutomaton[{{1,1,_,1}->0,{0,1,_,1}->1,{1,1,_,0}->1,{0,1,_,0}->0,{1,0,_,1}->1,{1,0,_,0}->1,{0,0,_,1}->0,{0,0,_,0}->0},RandomInteger[1,101],100],diffs,cumpoints,min1},diffs=Apply[Subtract,Outer[Count,hist,{1,0},1],{1}];cumpoints=Accumulate[diffs];min1=Min[cumpoints];If[min1≤0,1-min1,0]+cumpoints]

LogNormalDistributionFit[shiftedData]

Histogram[RandomReal[%,1000]]

NoncentralStudentTDistributionFit[shiftedData]

Histogram[RandomReal[%,1000]]

What about totalistic CAs?. It might be expected that not only should there be more agents but that also some agents contribute to the market prices greater than others. Here we consider a totalistic k = 3 colours, r = 2 neighbourhood radius.  Again it is shown that it suffers from large volatility and skewness.

Totalistic cellular automata rule for interacting market agents.

CAStep[TotalisticCARule[rule_List,1],a_List]:=rule[[-1-(RotateLeft[a]+a+RotateRight[a])]]
CAStep[TotalisticCARule[rule_List,r_Integer],a_List]:=rule[[-1-Sum[RotateLeft[a,i],{i,-r,r}]]]
ToTotalisticCARule[num_Integer,k_Integer,r_Integer]:=TotalisticCARule[IntegerDigits[num,k,1+(k-1) (2 r+1)],r]
CATotalisticEvolveList[rule_List,init_List,t_Integer]:=NestList[CAStep[TotalisticCARule[rule,1],#]&,init,t]
CATotalisticEvolveList[rule_List,r_Integer,init_List,t_Integer]:=NestList[CAStep[TotalisticCARule[rule,r],#]&,init,t]

FromDigits[{2,0,1,0,2,0,2,0,1,2,0},3]

ToTotalisticCARule[126294,3,2]

Module[{hist=CATotalisticEvolveList[{2,0,1,0,2,0,2,0,1,2,0},2,RandomInteger[{0,2},101],1000],diffs,cumpoints,min1},diffs=Apply[Subtract,Outer[Count,hist,{2,0},1],{1}];cumpoints=Accumulate[diffs];min1=Min[cumpoints];ListPlot[If[min1≤0,1-min1,0]+cumpoints,Joined->True]]

shiftedData=Module[{hist=CATotalisticEvolveList[{2,0,1,0,2,0,2,0,1,2,0},2,RandomInteger[{0,2},101],100],diffs,cumpoints,min1},diffs=Apply[Subtract,Outer[Count,hist,{2,0},1],{1}];cumpoints=Accumulate[diffs];min1=Min[cumpoints];If[min1≤0,1-min1,0]+cumpoints]

LogNormalDistributionFit[shiftedData]

Histogram[RandomReal[%,1000]]

NoncentralStudentTDistributionFit[shiftedData]

Histogram[RandomReal[%,1000]]

Stochastic models for equities

We need to revisit the simplest model and ask how it can be modified, in this case standardised so that it behaves more like real market instruments. We shall start by looking at the traditional lognormal model for the stocks which is based upon the notion that the smallest unit of noise or the underlying measure that enters the market is described by the Wiener martingale (Subscript[W, t]) ~ 𝒩(0, (Sqrt[t])). The use of this stochastic measure is linked to the completeness of the market, but a lot of empirical evidence suggest that the market is incomplete and this implies that (Subscript[W, t]) may no longer be an appropriate representation for noise in equities.

The random walk characterization of equities results in their observed lognormal distribution. This is described by the following stochastic difference equation where (Subscript[S, t]) is the stock price at time t, r is the risk-neutral interest rate, q the continuous dividend yield of the stock, Δ t is the small time increment and Z~ 𝒩(0,1) is the standard normal random variable with mean 0 and standard deviation 1.

((d (Subscript[S, t]))/(Subscript[S, t]))=(r-q) Δ t+Z σ (Sqrt[Δ t])

The above equation leads to a simple updating rule where the next value for the stock is obtained from the previous value

(Subscript[S, t+Δt])=(Subscript[S, t])+d (Subscript[S, t]) = (Subscript[S, t])(1+(r-q) Δ t+Z σ (Sqrt[Δ t]))

This is easily transferred into a Mathematica function called path[s, r, q, σ, t, navg]where the inputs have the same meanings as described above and n is the number of sample points along the stock path process as it evolves over a period of time t.

path[s_,r_,q_,σ_,t_,navg_Integer]:=Module[{mpath=1+(((r-q) t)/(navg-1)),spath=σ (Sqrt[(t/(navg-1))]),randnorms},randnorms=RandomReal[NormalDistribution[mpath,spath],navg-1];FoldList[#1 #2&,s,randnorms]];

path[100, 0.1, 0.02, 0.2, 1, 10]

ListLinePlot[path[100,0.1,0.02,0.2,1,100]]

Stock prices generated by CAs

CA values are now modified and standardised by introducing the CA differences as the variation that occurs for stocks so that the stock path mean is 1+(r-q) Δ T where Δ T is the time increment and the stock path standard deviation is σ (Sqrt[Δ T]), where σ is the stock volatility. But instead of using the standard normal random variate Z ~ 𝒩(0, 1), we use a standardised CA.

Clear[CAStep]
CenterList[n_Integer]:=ReplacePart[Table[0,{n}],1,Ceiling[n/2]]
CAStep[rule_List,a_List]:=Partition[a,3,1,2]/.rule
CAEvolveList[rule_List,init_List,t_Integer]:=NestList[CAStep[rule,#]&,init,t]
CAGraphics[history_List]:=Graphics[Raster[1-Reverse[history]],AspectRatio->Automatic]
MaxRuleNumber[k_Integer,r_Integer]:=k^(k^(2 r+1))

CApath[s_,r_,q_,σ_,t_,navg_Integer,n_Integer]:=Module[{mpath,spath,hist,diffs,newdiffs},
mpath=1+(r-q)*t/(navg-1);spath=σ Sqrt[t/(navg-1)];
hist=CAEvolveList[{{1,_,1}->0,{1,_,0}->1,{0,_,1}->1,{0,_,0}->0},RandomInteger[{0,1},n],navg-2];diffs=Subtract@@@Outer[Count,hist,{1,0},1]/Sqrt[n];newdiffs=mpath+spath*diffs;FoldList[#1*#2&,s,newdiffs]]

CApath2[s_,r_,q_,σ_,t_,navg_Integer,n_Integer]:=Module[{mpath,spath,hist,diffs,newdiffs},
mpath=1+(r-q)*t/(navg-1);spath=σ Sqrt[t/(navg-1)];
hist=CellularAutomaton[{{1,_,1}->0,{1,_,0}->1,{0,_,1}->1,{0,_,0}->0},RandomInteger[{0,1},n],navg-2];diffs=Subtract@@@Outer[Count,hist,{1,0},1]/Sqrt[n];newdiffs=mpath+spath*diffs;FoldList[#1*#2&,s,newdiffs] ]

By introducing a path mean and standard deviation with noise generated by a 'standardised' CA then we now get a model that corresponds to real market behaviour with appropriate magnitudes, pattern and skew.

ListLinePlot[path[100, 0.1, 0.02, 0.2, 1, 200]]

ListLinePlot[CApath[100, 0.1, 0.02, 0.2, 1, 200,1001]]

ListLinePlot[CApath2[100, 0.1, 0.02, 0.2, 1, 200,101]]

shiftedData=CApath[100, 0.1, 0.02, 0.2, 1, 100,101]

LogNormalDistributionFit[shiftedData]

Histogram[RandomReal[%,1000]]

NoncentralStudentTDistributionFit[shiftedData]

Histogram[RandomReal[%,1000]]

Option prices generated by CAs ?

The introduction of this new model suggests a possible new pricing formula for options, again in which the standard normal random variate is replaced by a 'standardised' CA. First we introduce the traditional Black-Scholes model but as before, replace the CDF of the standard normal random variable by the frequency CDF for the corresponding CA.

CDF[NormalDistribution[0,1],z]

{BlackScholesCall[100,100, 0.1, 0.02, 0.2, 1], BlackScholesPut[100,100, 0.1, 0.02, 0.2, 1]}

The frequency CDF for the modified CA used to describe the stock is calculated below and example numbers generated. This approach is consistent with the traditional statistical tables for the standard normal random variate.

CANormProby[navg_Integer,n_Integer]:=Module[{hist=CellularAutomaton[{{1,_,1}->0,{1,_,0}->1,{0,_,1}->1,{0,_,0}->0},RandomInteger[{0,1},n],navg-2]},
Interpolation[Transpose[{Range[-4.9,5.0,0.1],Accumulate[BinCounts[Subtract@@@Outer[Count,hist,{1,0},1]/Sqrt[n],{-5,5,0.1}]/navg]//N }]] ]

Here we compare the frequency CDF for the modified CA to the CDF of the standard normal random variate for the events which correspond to the stock being in and out of the money.

η /@ {done[100,100, 0.1, 0.02, 0.2, 1],-dtwo[100,100, 0.1, 0.02, 0.2, 1]}

CANormProby[100,101] /@ {done[100,100, 0.1, 0.02, 0.2, 1],-dtwo[100,100, 0.1, 0.02, 0.2, 1]}

Here we compare the frequency CDF for the modified CA to the CDF of the standard normal random variate over the range of possible values for the standardized differences for a very large number of sample points and initial array. Note the close approximation between the CA CDF and the Normal CDF.

CAcdf=CANormProby[1000,1001];
Plot[{η[x],CAcdf[x]},{x,-4.9,5},AxesLabel->{"Standardized Diffs", "CDFs"}]

Whereas here we compare the frequency CDF for the modified CA to the CDF of the standard normal random variate over the range of possible values for the standardized differences for a smaller number of sample points and initial array. Note the small oscillating deviations of the CA CDF from the Normal CDF.

CAcdf=CANormProby[100,101];
Plot[{η[x],CAcdf[x]},{x,-4.9,5},AxesLabel->{"Standardized Diffs", "CDFs"}]

So if we can use CAs to represent equities, can we also use the same approach to determine options on equities?

{CACall[100,100, 0.1, 0.02, 0.2, 1,101], CAPut[100,100, 0.1, 0.02, 0.2, 1,101]}

{BlackScholesCall[100,100, 0.1, 0.02, 0.2, 1], BlackScholesPut[100,100, 0.1, 0.02, 0.2, 1]}

Comparison of the option values for the Black-Scholes and the corresponding CA model we see that the CA based probability distributions do not transfer to option pricing. The reason is undoubtedly due to the fact that CAs are not martingale processes and therefore will not satisfy the same assumptions that underlie the standard options pricing models.

(Subscript[X, T])=𝔼((Subscript[X, t]))

References

Black,F. and Scholes,M., The Pricing of Options and Corporate Liabilities, Journal of Political Economy, 81, 659-673, 1973.
Brealey, R.A., An Introduction to Risk and Return from Common Stock, Cambridge, MA. MIT Press, 1983.
Cootner, P.H. (ed.), The Random Character of Stock market Prices, Cambridge, MA. MIT Press, 1964.
Cox,J., Ross,S. and Rubinstein,M., Option Pricing: A Simplified Approach, Journal of Financial Economics, 7, p.229, 1979.
Harrison, J.M., Brownian Motion and Stochastic Flow Systems, Wiley New York, 1985.
Hull,J., Options, Futures and Other Derivatives, Prentice Hall, 6th Edn 2005.
Kelly, M., Evaluation Of Financial Options using Radial Basis Functions in Mathematica, The Mathematica Journal, Vol. 11, No 1, Summer 2007.
MathConsult. UnRisk Pricing Engine. http://www.unriskderivatives.com/
Shaw,W., Modelling Financial Derivatives with Mathematica, Cambridge University Press, 1998.
Trott,M., The Mathematica Guidebook for Numerics, Springer, 2006.
Wolfram, S., A New Kind of Science, Wolfram Media, 2002.