Markov chains, hysteresis and logistic regression.

Logistic regression models the response of a binary variable to a covariate. The binary data produced by a two state Markov process coupled to an external covariate can be fit using logistic regression. Under what circumstances, if any, is this a valid approach?

At each time step i, the non-stationary Markov process is governed by a transition matrix

(1)   \begin{equation*}P_i=\left(\begin{array}{cc} p_i & 1-p_i \\ 1-q_i & q _i\end{array}\right)\end{equation*}

where 0 \leq p_i,q_i \leq 1P can be coupled to the external covariate x using logistic functions, which preserve the stochastic matrix property. For instance, choose:

(2)   \begin{eqnarray*} p_i &=& \frac{1}{1+e^{2\left(1+2x_i\right)}}\\ q_i &=& \frac{1}{1+e^{2\left(1-2x_i\right)}} \end{eqnarray*}

Also, assume that x varies harmonically with timescale (angular period) T i.e.

(3)   \begin{equation*} x_i= \sin \left( \frac{i}{T} \right) \end{equation*}

The graph below the expected response (a number between 1 & 2) as a function of x obtained by simulating 200 Markov chains of length 10^5 using Equations (1)-(3).

hysteresis

The striking feature is the appearance of an hysteresis loop. Hysteresis, where the mean response has separate branches for increasing and decreasing x, becomes pronounced when x varies rapidly (small T). On the other hand, when x varies slowly, the mean response collapses to a single branch.

What’s going on?

Starting from some initial state, a markov chain approaches it’s equilibrium distribution after a mixing time. The mixing time is longer when there are bottlenecks e.g.when q \simeq 1 in (1). When the covariate varies on a time-scale which is long compared to the mixing time, the system remains close to equilibrium at all times. Hysteresis is absent in this case. Conversely, when the covariate varies rapidly the transition matrix changes too quickly for the system to reach equilibrium. This leads to hysteresis.

Logistic regression is valid as long as the mixing time is short compared to the timescale over which the covariate varies. However, if the mixing time is long, multiple branches of the response are present. This aspect of Markov chain dynamics is missing from ordinary logistic regression. This conclusion generalises to the N-state/multinomial regression case.

Markov chain simulation

A good way simulate non-stationary Markov chains in R is to pass a list of transition matrices to a C++ function (MCsim1) using Rcpp. For example:

where “~/markov.cpp” contains:

 

 

Coffee prices, long range forecasts and drought

The tropical highlands of Minas Gerais, Brazil are responsible for 25% of the world’s Arabica coffee production. In 2014, the region experienced drought during the critical austral summer months of January and February. World coffee prices moved sharply higher at the end of January. By early March, prices had nearly doubled.

In an efficient market the price of a commodity reflects all available information. Did the coffee price assimilate long-range weather forecast information available in 2014 ?

rainfall1

The above chart show monthly rainfall (crosses) and average rainfall (black line). Rainfall totals were extracted from ERA-interim reanalysis. Periods of deficit relative to average rainfall are indicated in red. Rainfall was less than 50% of average in both January and February 2014.

prices

Despite dryness, coffee futures actually trended slightly lower during January 2014 (above). However this situation reversed dramatically after January 29 (indicated by the red arrow). It is as though the market abruptly woke up to the fact that drought would continue well into February and that this would impact Arabica coffee fruit development.

In fact, well in advance of January 29, long range weather forecasts were indicating a high probability of continued drought in February. The graph below shows a large ensemble of CFSv2 rainfall forecasts for Minas Gerais for December, January and February. Such forecasts[*] are made every 6 hours up to 9 months prior to the forecast month. A high probability of anomalous rainfall for the months of January and February is evident some 2-3 weeks in advance.

 

forecast

This analysis points to a surprising conclusion. For perhaps two weeks, world coffee market prices did not properly reflect probabilistic information available from long range weather forecasts.

[*] Raw forecasts, not bias corrected.