Catastrophe theory has been developed as a deterministic theory for systems that may respond to continuous changes in control variables by a discontinuous change from one equilibrium state to another. A key idea is that system under study is driven towards an equilibrium state. The behavior of the dynamical systems under study is completely determined by a so-called potential function, which depends on behavioral and control variables. The behavioral, or state variable describes the state of the system, while control variables determine the behavior of the system. The dynamics under catastrophe models can become extremely complex, and according to the classification theory of Thom, there are seven different families based on the number of control and dependent variables.
Let us suppose that the process yt evolves over t = 1,…, T as
dyt = -dV(yt; α, β)dt/dyt —– (1)
where V (yt; α, β) is the potential function describing the dynamics of the state variable yt controlled by parameters α and β determining the system. When the right-hand side of (1) equals zero, −dV (yt; α, β)/dyt = 0, the system is in equilibrium. If the system is at a non-equilibrium point, it will move back to its equilibrium where the potential function takes the minimum values with respect to yt. While the concept of potential function is very general, i.e. it can be quadratic yielding equilibrium of a simple flat response surface, one of the most applied potential functions in behavioral sciences, a cusp potential function is defined as
−V(yt; α, β) = −1/4yt4 + 1/2βyt2 + αyt —– (2)
with equilibria at
-dV(yt; α, β)dt/dyt = -yt3 + βyt + α —– (3)
being equal to zero. The two dimensions of the control space, α and β, further depend on realizations from i = 1 . . . , n independent variables xi,t. Thus it is convenient to think about them as functions
αx = α0 +α1x1,t +…+ αnxn,t —– (4)
βx = β0 + β1x1,t +…+ βnxn,t —– (5)
The control functions αx and βx are called normal and splitting factors, or asymmetry and bifurcation factors, respectively and they determine the predicted values of yt given xi,t. This means that for each combination of values of independent variables there might be up to three predicted values of the state variable given by roots of
-dV(yt; αx, βx)dt/dyt = -yt3 + βyt + α = 0 —– (6)
This equation has one solution if
δx = 1/4αx2 − 1/27βx3 —– (7)
is greater than zero, δx > 0 and three solutions if δx < 0. This construction can serve as a statistic for bimodality, one of the catastrophe flags. The set of values for which the discriminant is equal to zero, δx = 0 is the bifurcation set which determines the set of singularity points in the system. In the case of three roots, the central root is called an “anti-prediction” and is least probable state of the system. Inside the bifurcation, when δx < 0, the surface predicts two possible values of the state variable which means that the state variable is bimodal in this case.
Most of the systems in behavioral sciences are subject to noise stemming from measurement errors or inherent stochastic nature of the system under study. Thus for a real-world applications, it is necessary to add non-deterministic behavior into the system. As catastrophe theory has primarily been developed to describe deterministic systems, it may not be obvious how to extend the theory to stochastic systems. An important bridge has been provided by the Itô stochastic differential equations to establish a link between the potential function of a deterministic catastrophe system and the stationary probability density function of the corresponding stochastic process. Adding a stochastic Gaussian white noise term to the system
dyt = -dV(yt; αx, βx)dt/dyt + σytdWt —– (8)
where -dV(yt; αx, βx)dt/dyt is the deterministic term, or drift function representing the equilibrium state of the cusp catastrophe, σyt is the diffusion function and Wt is a Wiener process. When the diffusion function is constant, σyt = σ, and the current measurement scale is not to be nonlinearly transformed, the stochastic potential function is proportional to deterministic potential function and probability distribution function corresponding to the solution from (8) converges to a probability distribution function of a limiting stationary stochastic process as dynamics of yt are assumed to be much faster than changes in xi,t. The probability density that describes the distribution of the system’s states at any t is then
fs(y|x) = ψ exp((−1/4)y4 + (βx/2)y2 + αxy)/σ —– (9)
The constant ψ normalizes the probability distribution function so its integral over the entire range equals to one. As bifurcation factor βx changes from negative to positive, the fs(y|x) changes its shape from unimodal to bimodal. On the other hand, αx causes asymmetry in fs(y|x).