Trading Rules, Aphorisms & Books

An interesting link.

The Perfect Oscillator on S&P E-Mini Data

Following on from my previous posts the chart below is a snapshot of the last 150 days of the S&P E-mini continuous futures contract, up to and including the 9th December.

The top pane is obviously a plot of price showing triangles to indicate cyclic highs and lows when the Perfect Oscillator (middle pane) crosses its SMA (see previous post for explanation). As explained earlier the algorithm used in calculating the PO makes use of projected/predicted prices, and since these "forward prices" are available from the .oct function the bars are coloured according to the 5 and 10 day "forward prices." Based on a very simple if/or statement in the code essentially we should be long when the bars are blue and short when red. The few whites bars that are visible are bars where the conditions for being blue or red are not met.

Superimposed over the middle pane plot of the PO are the full period EMA and plus and minus 1 exponential standard deviation lines, a sort of Bollinger Bands over the indicator. The bottom plot is of the bandpass indicator which is used in the "forward prices" algorithm.

Based on the logic and theory of the PO the rules for use are essentially:-
1) go long/short according to the PO SMA crossings
2) confirm the above by the bar colours
3) exit/close up stops when the PO crosses its exponential standard deviation lines

and based on these simple rules a short entry is indicated for 12th December, particularly since the bandpass is also indicating a short trade with its leading function cross. Will be interesting to see how this indicated trade works out.

Update #2 on Perfect Moving Average and Oscillator

Below are screenshots of the PMA and PO based upon the improved price projection algorithm discussed in the previous post, shown on an idealised time series of a sinewave (plus trend) in cyan in all screenshots.

The first 3 are plots of the 10 bar PMA in red with a normal 10 bar EMA in green shown for comparative purposes.

As the theory suggests it can be seen that the PMA has the equivalent smoothing of the 10 bar EMA but without any lag.

The next series of screenshots are of the PO, where again the "price" is cyan, the various red lines are 2 bar, 3 bar etc. POs up to an 11 bar PO, and the green PO is a simple average of all these POs.

As I had anticipated these POs lead price by a quarter period which means the cycle highs and lows are indicated by the green "average PO" crossing its centre line, which in turn is a full cycle period length simple moving average of the average PO. My leading oscillator indicator should give timely warnings of these upcoming turns in price. It is also interesting to note that the individual, longer length PO crossings (greater amplitude) also give advance warning of the upcoming turn and the fact that the shorter length POs NOT crossing the centre line and being entirely above or below said centre line seem to indicate the strength of any trend.

Update on the Perfect Moving Average and Oscillator

It has been some time since my last post and the reason for this was, in large part, due to my changing over from an Ubuntu to a Kubuntu OS and the resultant teething problems etc. However, all that is over so now I return to the Perfect Moving Average (PMA) and the Perfect Oscillator (PO).

In my previous post I stated that I would show these on real world data, but before doing this there was one thing I wanted to investigate - the accuracy of the price projection algorithm upon which the average and oscillator calculations are based. Below are four sample screenshots of these investigations.

Sideways Market

Up Trending Market

Down Trending Market

The charts show a simulated "ideal" time series (cyan) with two projected price series of 10 points each (red and green) plotted such that the red and green projections are plotted concurrently with the actual series points they are projections, or predictions, of. Both projections/predictions start at the point at which the red "emerges" from the actual series. As can be seen the green projection/prediction is much more accurate than the red, and the green is the improved projection algorithm I have recently been working on.

As a result of this improvement and due to the nature of the internal calculations I expect that the PMA will follow prices much closer and that the PO will lead prices by a quarter cycle compared with my initial price projection algorithm. More on this in the next post.

The Theoretically Perfect Moving Average and Oscillator

Readers of this blog will probably have noticed that I am partial to using concepts from digital signal processing in the development of my trading system. Recently I "rediscovered" this PDF, written by Tim Tillson, on Google Docs, which has a great opening introduction:

" 'Digital filtering includes the process of smoothing, predicting, differentiating,

integrating, separation of signals, and removal of noise from a signal. Thus many people
who do such things are actually using digital filters without realizing that they are; being
unacquainted with the theory, they neither understand what they have done nor the
possibilities of what they might have done.'

This quote from R. W. Hamming applies to the vast majority of indicators in technical
analysis."

The purpose of this blog post is to outline my recent work in applying some of the principles discussed in the linked PDF.

Long time readers of this blog may remember that back in April this year I abandoned work I was doing on the AFIRMA trend line because I was dissatisfied with the results. The code for this required projecting prices forward using my leading signal code, and I now find that I can reuse the bulk of this code to create a theoretically perfect moving average and oscillator. I have projected prices forward by 10 bars and then taken the average of the forwards and backwards exponential moving averages, as per the linked PDF, to create a series of averages that theoretically are in phase with the underlying price, and then taken the mean of all these averages to produce my "perfect" moving average. Similarly, I have done the same to create a "perfect" oscillator and the results are shown in the images below.
Sideways Market

Trending up Market

Trending down Market

The upper panes show "price" and its perfect moving average, the middle panes show the perfect oscillator, its one bar leading signal, and exponential standard deviation bands set at 1 standard deviation above and below an exponential moving average of the oscillator. The lower panes show the %B indicator of the oscillator within the bands, offset so that the exponential standard deviation bands are at 1 and -1 levels.

Even cursory examination of many charts such as these shows the efficacy of the signals generated by the crossovers of price with its moving average, the oscillator with its leading signal and the %B indicator crossing the 1 and -1 levels, when applied to idealised data. My next post will discuss the application to real price series.

Post to TradingBlox Forum

I have recently posted a reply to a thread on the TradingBlox forum here which readers of this blog might be interested in. The Octave code below is that used to generate the graphic in the linked forum reply.

clear all

a_returns = [2,-1,2,-1,2,-1,2,-1,2,-1];
b_returns = [3,-1.5,3,-1.5,3,-1.5,3,-1.5,3,-1.5];
c_returns = shift(b_returns,1);

a_equity = cumsum([1,a_returns]);
b_equity = cumsum([1,b_returns]);
c_equity = cumsum([1,c_returns]);

ab_equity = a_equity .+ b_equity;
ac_equity = a_equity .+ c_equity;

ab_equity_correl = corrcoef(a_equity,b_equity)
ac_equity_correl = corrcoef(a_equity,c_equity)

ab_returns_correl = corrcoef(a_returns,b_returns)
ac_returns_correl = corrcoef(a_returns,c_returns)

ab_centre_returns_correl = corrcoef(center(a_returns),center(b_returns))
ac_centre_returns_correl = corrcoef(center(a_returns),center(c_returns))

plot(a_equity,"k",b_equity,"b",c_equity,"c",ab_equity,"r",ac_equity,"g")
legend("a equity","b equity","c equity","ab equity","ac equity")

Additional Screenshots of Classifier in Action

Gold

10 Year US Treasuries

30 Year US Bonds

EurUsd forex

All charts are daily bars.

Classifier Mark 2 Update

I have now incorporated the revised code for my classifier in my production code and a screen shot of the output on recent S & P 500 data, as of last Friday, is shown below.

The most recent yellow candlesticks indicate that the market is classified as a "down with retracement" market, and the oscillator leading signals (downwards pointing triangles) indicate that the retracement is likely to be over. It will be interesting to see how the classifier performs, in real time, over the coming days.

Naive Bayes Classifier, Mark 2.

It has taken some time, but I have finally been able to incorporate the Trend Vigor indicator into my Naive Bayesian classifier, but with a slight twist. Instead of being purely Bayesian, the classifier has evolved to become a hybrid Bayesian/clustering classifier. The reason for this is that the Trend Vigor indicator has no varying distribution of values but tends to return values that are so close to each other that they can be considered a single value, as mentioned in an earlier post of mine. This can be clearly seen in the short 3D visualisation animation below. The x, y and z axis each represent an input to the classifier, and about 7 seconds into the video you can see the Trend Vigor axis in the foreground with almost straight vertical lines for its "distributions" for each market type. However, it can also be seen that there are spaces in 3D where only combined values for one specific market type appear, particularly evident in the "tails" of the no retracement markets ( the outermost blue and magenta distributions in the video. )



( Non embedded view here )

The revised version of the classifier takes advantage of this fact. Through a series conditional statements each 3D datum point is checked to see if it falls in any of these mutually exclusive spaces and if it does, it is classified as belonging to the market type that has "ownership" of the space in which it lies. If the point cannot be classified via this simple form of clustering then it is assigned a market type through Bayesian analysis.

This Bayesian analysis has also been revised to take into account the value of the Trend Vigor indicator. Since these values have no distribution to speak of a simple linear model is used. If a point is equidistant between two Trend Vigor classifications it is assigned a 0.5 probability of belong to each, this probability rising in linear fashion to 1.0 if it falls exactly on one of the vertical lines mentioned above, with a corresponding decrease in probability assigned to the other market type classification. There is also a boundary condition applied where the probability is set to 0.0 for belonging to a particular market type.

The proof of the pudding is in the eating, and this next chart shows the classification error rate when the classifier is applied to my usual "ideal" time series.

The y axis is the percentage of ideal time series runs in which market type was mis-classified, and the x axis is the period of the cyclic component of the time series being tested. In this test I am only concerned with the results for periods greater than 10 as in real data I have never seen extracted periods less than this. As can be seen the sideways market and both the up and down with no retracement markets have zero mis-classification rates, apart from a small blip at period 12, which is within the 5% mis-classification error rate I had set as my target earlier.

Of more concern was the apparent large mis-classification error rate of the retracement markets ( the green and black lines in the chart. ) However, further investigation of these errors revealed them not to be "errors" as such but more a quirk of the classifier, which lends itself to exploitation. Almost all of the "errors" occur consecutively at the same phase of the cyclic component, at all periods, and the "error" appears in the same direction. By this I mean that if the true market type is up with retracement, the "error" indicates an up with no retracement market; if the true market type is down with retracement, the "error" indicates a down with no retracement market. The two charts below show this visually for both the up and down with retracement markets and are typical representations of the "error" being discussed.

The first pane in each chart shows one complete cycle in which the whole cycle, including the most recent datum point, are correctly classified as being an up with retracement market ( upper chart ) and a down with retracement market ( lower chart. ) The second pane shows a snapshot of the cycle after it has progressed in time through its phase with the last point being the last point that is mis-classified. The "difference" between each chart's respective two panes at the right hand edge shows the portion of the time series that is mis-classified.

It can be seen that the mis-classification occurs at the end of the retracement, immediately prior to the actual turn. This behaviour could easily be exploited via a trading rule. For example, assume that the market has been classified as an up with retracement market and a retracement short trade has been taken. As the retracement proceeds our trade moves into profit but then the market classification changes to up with no retracement. Remember that the classifier (never?) mis-classifies such no retracement markets. What would one want to do in such a situation? Obviously one would want to exit the current short trade and go long, and in so doing would be exiting the short and initiating the possible long at precisely the right time; just before the market turn upwards! This mis-classification "error" could, on real data, turn out to be very serendipitous.

All in all, I think this revised, Mark 2 version of my market classifier is a marked improvement on its predecessor.

Creation of Synthetic Data

Some time ago (the file was last edited in July 2010) I wrote an Octave .oct function to create synthetic data for testing and optimisation purposes. I was inspired to do so by the December 2005 issue of The Breakout Bulletin and it has recently come to mind again due to a posting on the StackExchange Quantitative Finance Forum here. I have posted the code for my .oct function in the code box below.

In writing this function I wanted to extend the ideas presented in the Breakout Bulletin and make them more applicable for the purposes I had/have in mind. By randomly scrambling the data any bar to bar dependency is destroyed (by design of course), but what if you want to preserve some bar to bar dependencies? This .oct function is my solution to preserving this dependency and a brief discussion of the theory behind it follows.

Firstly there is an assumption that any single bar and the market forces that caused the bar to be formed the way it did (up bar, down bar, doji etc.) are dependent on the immediately preceding market activity and the "current mode" of the market. Implicit in this assumption is that certain "types" of bars are more likely to be seen depending on market "mode," i.e. the types of bar in an up trend are likely to be distinctly different from those in a down trending or sideways trending market, so what is needed is some way to bin the bars which reflects this.

My solution is to apply a 21 bar moving median of the close and median absolute deviations from this median as bands above and below it, similar to Bollinger Bands. There are 3 levels; 1 x MAD, 2 x MAD and 3 x MAD above, and 3 below; to give a total of 8 "zones" as they are called in the code. Furthermore, a 21 bar moving median of the True Range and a 4 bar WMA of the True Range are also calculated. The first part of the code ("Code Block A Loop"), after all the required declarations, loops over the input time series data calculating all the above and assigning each bar to a specific bin based upon the "zone" in which the previous bar resides, and further assignation depends on whether the previous bar is a high or low volatility bar decided by the True Range 4 bar WMA being above or below the True Range 21 bar moving median. This gives a total of 16 different bins to which a bar can be assigned. On assignation to a bin, the open, high, low and close are recorded in that bin by their relation to the previous close thus: log10(close/previous_close), log10(open/previous_open)... etc.

The next part of the code ("Code Block B Loop") actually creates the synthetic data by randomly drawing a bar's relationships to its previous close from the "relevant bin" and calculating a "new" bar based upon these relationships. This "relevant bin" is determined by the "zone" position and volatility of the most recently calculated synthetic "new" bar. After a new, "new bar" has been created, the median, MADs and True Range calculations are updated to include this new, "new bar," which becomes the previous bar on the next iteration of the loop for Code Block B Loop.

Finally, a small part of the code adjusts the input data in the case of negative values due to the possible use of continuous back-adjusted futures contracts as the input data. This is necessary to avoid errors in trying to calculate the log10 of a negative number.

The above method of binning the input data and subsequent randomisation is my attempt to ensure that dependencies/characteristics of the original data are preserved - for example - assume a bar is above the upper 3 x MAD level and is determined to be a high volatility bar, then the next synthetically created bar will be drawn only from the binned distribution of bars that in the real data also follow a bar above the upper 3 x MAD level and is determined to be a high volatility bar.

This code is offered as is and comes with no warranty whatsoever. However, if you like it and use it I would be interested to hear from you. In particular, if you have any suggestions for the code's improvement, extension, optimisation etc. or see any errors in the code, I would really appreciate your feedback. 

#include octave/oct.h
#include octave/dColVector.h
#include algorithm
#include math.h
#include "MersenneTwister.h"

DEFUN_DLD (syn_3, args, , "Inputs are OHLC column vectors and tick size. Output is MC generated synthetic OHLC")
{
octave_value_list retval_list;

    if (args(0).length () < 5)
    {
        error ("Invalid arguments");
        return retval_list;
    }

    ColumnVector open = args(0).column_vector_value (); // open vector
    ColumnVector high = args(1).column_vector_value (); // high vector
    ColumnVector low = args(2).column_vector_value (); // low vector
    ColumnVector close = args(3).column_vector_value (); // close vector
    double tick = args(4).double_value(); // Tick size

    if (error_state)
    {
        error ("Invalid arguments");
        return retval_list;
    }

// Check for negative or zero price values due to continuous back- adjusting of price series & compensate if necessary
// Note: the "&" below acts as Address-of operator: p = &x; Read: Assign to p (a pointer) the address of x.
    double lowest_low = *std::min_element( &low(0), &low(low.length()) );
    double correction_factor = 0.0;
    if ( lowest_low <= 0.0 )
    {
    correction_factor = fabs(lowest_low) + tick;
	for (octave_idx_type ii (0); ii < args(0).length (); ii++)
	{
	open (ii) = open (ii) + correction_factor;
	high (ii) = high (ii) + correction_factor;
	low (ii) = low (ii) + correction_factor;
	close (ii) = close (ii) + correction_factor;
	}
    }

    ColumnVector moving_median_window (21);
    ColumnVector moving_MAD_window (21);
    ColumnVector moving_true_range_window (21);
    ColumnVector moving_median = args(0).column_vector_value ();
    ColumnVector moving_MAD = args(0).column_vector_value ();
    ColumnVector true_range = args(0).column_vector_value ();
//  declare and "pre-reserve" enough space for the various categorised bins for the MC proceedure
//  first, bins where curerent_4_bar_wma of true_range >= current_median_true_range
    ColumnVector zone_1_open ( args(0).length () ); // zone_1 >= median & < median + MAD 
    ColumnVector zone_1_high ( args(0).length () );
    ColumnVector zone_1_low ( args(0).length () );
    ColumnVector zone_1_close ( args(0).length () ); 
    ColumnVector zone_2_open ( args(0).length () ); // >= median + MAD & < median + 2*MAD
    ColumnVector zone_2_high ( args(0).length () );
    ColumnVector zone_2_low ( args(0).length () );
    ColumnVector zone_2_close ( args(0).length () );
    ColumnVector zone_3_open ( args(0).length () ); // >= median + 2*MAD & < median + 3*MAD
    ColumnVector zone_3_high ( args(0).length () );
    ColumnVector zone_3_low ( args(0).length () );
    ColumnVector zone_3_close ( args(0).length () );
    ColumnVector zone_4_open ( args(0).length () ); // >= median + 3*MAD
    ColumnVector zone_4_high ( args(0).length () );
    ColumnVector zone_4_low ( args(0).length () );
    ColumnVector zone_4_close ( args(0).length () );
    ColumnVector zone_5_open ( args(0).length () ); // < median & >= median - MAD 
    ColumnVector zone_5_high ( args(0).length () );
    ColumnVector zone_5_low ( args(0).length () );
    ColumnVector zone_5_close ( args(0).length () );
    ColumnVector zone_6_open ( args(0).length () ); // < median - MAD & >= median - 2*MAD
    ColumnVector zone_6_high ( args(0).length () );
    ColumnVector zone_6_low ( args(0).length () );
    ColumnVector zone_6_close ( args(0).length () );
    ColumnVector zone_7_open ( args(0).length () ); // < median - 2*MAD & >= median - 3*MAD
    ColumnVector zone_7_high ( args(0).length () );
    ColumnVector zone_7_low ( args(0).length () );
    ColumnVector zone_7_close ( args(0).length () );
    ColumnVector zone_8_open ( args(0).length () ); // < median - 3*MAD
    ColumnVector zone_8_high ( args(0).length () );
    ColumnVector zone_8_low ( args(0).length () );
    ColumnVector zone_8_close ( args(0).length () );
    int zone_1_access_int = 0;
    int zone_2_access_int = 0;
    int zone_3_access_int = 0;
    int zone_4_access_int = 0;
    int zone_5_access_int = 0;
    int zone_6_access_int = 0;
    int zone_7_access_int = 0;
    int zone_8_access_int = 0;
//  second, bins where curerent_4_bar_wma of true_range < current_median_true_range
    ColumnVector zone_1_lv_open ( args(0).length () ); // zone_1 >= median & < median + MAD 
    ColumnVector zone_1_lv_high ( args(0).length () );
    ColumnVector zone_1_lv_low ( args(0).length () );
    ColumnVector zone_1_lv_close ( args(0).length () ); 
    ColumnVector zone_2_lv_open ( args(0).length () ); // >= median + MAD & < median + 2*MAD
    ColumnVector zone_2_lv_high ( args(0).length () );
    ColumnVector zone_2_lv_low ( args(0).length () );
    ColumnVector zone_2_lv_close ( args(0).length () );
    ColumnVector zone_3_lv_open ( args(0).length () ); // >= median + 2*MAD & < median + 3*MAD
    ColumnVector zone_3_lv_high ( args(0).length () );
    ColumnVector zone_3_lv_low ( args(0).length () );
    ColumnVector zone_3_lv_close ( args(0).length () );
    ColumnVector zone_4_lv_open ( args(0).length () ); // >= median + 3*MAD
    ColumnVector zone_4_lv_high ( args(0).length () );
    ColumnVector zone_4_lv_low ( args(0).length () );
    ColumnVector zone_4_lv_close ( args(0).length () );
    ColumnVector zone_5_lv_open ( args(0).length () ); // < median & >= median - MAD 
    ColumnVector zone_5_lv_high ( args(0).length () );
    ColumnVector zone_5_lv_low ( args(0).length () );
    ColumnVector zone_5_lv_close ( args(0).length () );
    ColumnVector zone_6_lv_open ( args(0).length () ); // < median - MAD & >= median - 2*MAD
    ColumnVector zone_6_lv_high ( args(0).length () );
    ColumnVector zone_6_lv_low ( args(0).length () );
    ColumnVector zone_6_lv_close ( args(0).length () );
    ColumnVector zone_7_lv_open ( args(0).length () ); // < median - 2*MAD & >= median - 3*MAD
    ColumnVector zone_7_lv_high ( args(0).length () );
    ColumnVector zone_7_lv_low ( args(0).length () );
    ColumnVector zone_7_lv_close ( args(0).length () );
    ColumnVector zone_8_lv_open ( args(0).length () ); // < median - 3*MAD
    ColumnVector zone_8_lv_high ( args(0).length () );
    ColumnVector zone_8_lv_low ( args(0).length () );
    ColumnVector zone_8_lv_close ( args(0).length () );
    int zone_1_lv_access_int = 0;
    int zone_2_lv_access_int = 0;
    int zone_3_lv_access_int = 0;
    int zone_4_lv_access_int = 0;
    int zone_5_lv_access_int = 0;
    int zone_6_lv_access_int = 0;
    int zone_7_lv_access_int = 0;
    int zone_8_lv_access_int = 0;
    double current_median_true_range;
    double current_4_wma_true_range;

// loop to fill the first 20 spaces of true_range vector
    true_range (0) = high (0) - low (0);
    for (octave_idx_type ii (1); ii < 20; ii++)
    {
    true_range (ii) = fmax ( fmax(high(ii)-low(ii),fabs(high(ii)-close(ii-1))) , fabs(low(ii)-close(ii-1)) );
    }

// following code loops to create a 21 period moving median and a 21 period moving median MAD (Median Absolute Deviation), based upon closing price
// At the same time as these are calculated, each previous bar's close is inspected to place that close in relation to the previous values of the moving 
// median and moving MAD. Also calculated is the 21 period median true range and 4 bar WMA of the true range. Based upon this, the current bar's stats 
// are binned into the appropriate zone vector column for later MC use.

    for (octave_idx_type ii (20); ii < args(0).length (); ii++) // Code Block A loop
    {
	for (octave_idx_type jj (0); jj < 21; jj++) // loop to fill the moving_median_window
	{
        moving_median_window (jj) = close (ii-jj);
   	} // end of loop to fill the moving_median_window
        // Note: the "&" below acts as Address-of operator: p = &x; Read: Assign to p (a pointer) the address of x.
        std::nth_element( &moving_median_window(0), &moving_median_window(10), &moving_median_window(21) );
        moving_median (ii) = moving_median_window(10);

	for (octave_idx_type jj (0); jj < 21; jj++) // loop to fill the moving_MAD_window
	{
        moving_MAD_window (jj) = fabs( close (ii-jj) - moving_median (ii) );
   	} // end of loop to fill the moving_MAD_window
        // Note: the "&" below acts as Address-of operator: p = &x; Read: Assign to p (a pointer) the address of x.
        std::nth_element( &moving_MAD_window(0), &moving_MAD_window(10), &moving_MAD_window(21) );
        moving_MAD (ii) = moving_MAD_window(10);

        // true range calculations
        true_range (ii) = fmax ( fmax(high(ii)-low(ii),fabs(high(ii)-close(ii-1))) , fabs(low(ii)-close(ii-1)) );
	for (octave_idx_type jj (0); jj < 21; jj++) // loop to fill the moving_true_range_window
	{
        moving_true_range_window (jj) = true_range (ii-jj);
   	} // end of loop to fill the moving_true_range_window
        // Note: the "&" below acts as Address-of operator: p = &x; Read: Assign to p (a pointer) the address of x.
        std::nth_element( &moving_true_range_window(0), &moving_true_range_window(10), &moving_true_range_window(21) );
        current_median_true_range = moving_true_range_window(10);
        current_4_wma_true_range = ( 4*true_range(ii) + 3*true_range(ii-1) + 2*true_range(ii-2) + true_range(ii-3) ) / 10 ;

// now analyise the positions of the bar closes in relation to the moving median and moving MAD, fill the relevant MC bins and adjust bin counts

	if ( ii >= 21 & current_4_wma_true_range >= current_median_true_range ) // bin selection based on volatility loop
	{
		if ( close (ii-1) >= moving_median (ii-1) & close (ii-1) < (moving_median (ii-1) + moving_MAD (ii-1)) ) // zone_1
		{
		zone_1_open (zone_1_access_int) = log10 ( open (ii) / close (ii-1) );
		zone_1_high (zone_1_access_int) = log10 ( high (ii) / close (ii-1) );
		zone_1_low (zone_1_access_int) = log10 ( low (ii) / close (ii-1) );
		zone_1_close (zone_1_access_int) = log10 ( close (ii) / close (ii-1) );
		zone_1_access_int = zone_1_access_int + 1;
		}

		else if ( close (ii-1) >= (moving_median (ii-1) + moving_MAD (ii-1)) & close (ii-1) < (moving_median (ii-1) + 2*moving_MAD (ii-1)) ) // zone_2
		{
		zone_2_open (zone_2_access_int) = log10 ( open (ii) / close (ii-1) );
		zone_2_high (zone_2_access_int) = log10 ( high (ii) / close (ii-1) );
		zone_2_low (zone_2_access_int) = log10 ( low (ii) / close (ii-1) );
		zone_2_close (zone_2_access_int) = log10 ( close (ii) / close (ii-1) );
		zone_2_access_int = zone_2_access_int + 1;	
		}

		else if ( close (ii-1) >= (moving_median (ii-1) + 2*moving_MAD (ii-1)) & close (ii-1) < (moving_median (ii-1) + 3*moving_MAD (ii-1)) ) // zone_3
		{
		zone_3_open (zone_3_access_int) = log10 ( open (ii) / close (ii-1) );
		zone_3_high (zone_3_access_int) = log10 ( high (ii) / close (ii-1) );
		zone_3_low (zone_3_access_int) = log10 ( low (ii) / close (ii-1) );
		zone_3_close (zone_3_access_int) = log10 ( close (ii) / close (ii-1) );
		zone_3_access_int = zone_3_access_int + 1;	
		}

		else if ( close (ii-1) >= (moving_median (ii-1) + 3*moving_MAD (ii-1)) ) // zone_4
		{
		zone_4_open (zone_4_access_int) = log10 ( open (ii) / close (ii-1) );
		zone_4_high (zone_4_access_int) = log10 ( high (ii) / close (ii-1) );
		zone_4_low (zone_4_access_int) = log10 ( low (ii) / close (ii-1) );
		zone_4_close (zone_4_access_int) = log10 ( close (ii) / close (ii-1) );
		zone_4_access_int = zone_4_access_int + 1;	
		}

		else if ( close (ii-1) < moving_median (ii-1) & close (ii-1) >= (moving_median (ii-1) - moving_MAD (ii-1)) ) // zone_5
		{
		zone_5_open (zone_5_access_int) = log10 ( open (ii) / close (ii-1) );
		zone_5_high (zone_5_access_int) = log10 ( high (ii) / close (ii-1) );
		zone_5_low (zone_5_access_int) = log10 ( low (ii) / close (ii-1) );
		zone_5_close (zone_5_access_int) = log10 ( close (ii) / close (ii-1) );
		zone_5_access_int = zone_5_access_int + 1;	
		}

		else if ( close (ii-1) < (moving_median (ii-1) - moving_MAD (ii-1)) & close (ii-1) >= (moving_median (ii-1) - 2*moving_MAD (ii-1)) ) // zone_6
		{
		zone_6_open (zone_6_access_int) = log10 ( open (ii) / close (ii-1) );
		zone_6_high (zone_6_access_int) = log10 ( high (ii) / close (ii-1) );
		zone_6_low (zone_6_access_int) = log10 ( low (ii) / close (ii-1) );
		zone_6_close (zone_6_access_int) = log10 ( close (ii) / close (ii-1) );
		zone_6_access_int = zone_6_access_int + 1;	
		}

		else if ( close (ii-1) < (moving_median (ii-1) - 2*moving_MAD (ii-1)) & close (ii-1) >= (moving_median (ii-1) - 3*moving_MAD (ii-1)) ) // zone_7
		{
		zone_7_open (zone_7_access_int) = log10 ( open (ii) / close (ii-1) );
		zone_7_high (zone_7_access_int) = log10 ( high (ii) / close (ii-1) );
		zone_7_low (zone_7_access_int) = log10 ( low (ii) / close (ii-1) );
		zone_7_close (zone_7_access_int) = log10 ( close (ii) / close (ii-1) );
		zone_7_access_int = zone_7_access_int + 1;	
		}

		else // zone_8
		{
		zone_8_open (zone_8_access_int) = log10 ( open (ii) / close (ii-1) );
		zone_8_high (zone_8_access_int) = log10 ( high (ii) / close (ii-1) );
		zone_8_low (zone_8_access_int) = log10 ( low (ii) / close (ii-1) );
		zone_8_close (zone_8_access_int) = log10 ( close (ii) / close (ii-1) );
		zone_8_access_int = zone_8_access_int + 1;	
		}
	}
	else
	{
		if ( close (ii-1) >= moving_median (ii-1) & close (ii-1) < (moving_median (ii-1) + moving_MAD (ii-1)) ) // zone_1
		{
		zone_1_lv_open (zone_1_lv_access_int) = log10 ( open (ii) / close (ii-1) );
		zone_1_lv_high (zone_1_lv_access_int) = log10 ( high (ii) / close (ii-1) );
		zone_1_lv_low (zone_1_lv_access_int) = log10 ( low (ii) / close (ii-1) );
		zone_1_lv_close (zone_1_lv_access_int) = log10 ( close (ii) / close (ii-1) );
		zone_1_lv_access_int = zone_1_lv_access_int + 1;
		}

		else if ( close (ii-1) >= (moving_median (ii-1) + moving_MAD (ii-1)) & close (ii-1) < (moving_median (ii-1) + 2*moving_MAD (ii-1)) ) // zone_2
		{
		zone_2_lv_open (zone_2_lv_access_int) = log10 ( open (ii) / close (ii-1) );
		zone_2_lv_high (zone_2_lv_access_int) = log10 ( high (ii) / close (ii-1) );
		zone_2_lv_low (zone_2_lv_access_int) = log10 ( low (ii) / close (ii-1) );
		zone_2_lv_close (zone_2_lv_access_int) = log10 ( close (ii) / close (ii-1) );
		zone_2_lv_access_int = zone_2_lv_access_int + 1;	
		}

		else if ( close (ii-1) >= (moving_median (ii-1) + 2*moving_MAD (ii-1)) & close (ii-1) < (moving_median (ii-1) + 3*moving_MAD (ii-1)) ) // zone_3
		{
		zone_3_lv_open (zone_3_lv_access_int) = log10 ( open (ii) / close (ii-1) );
		zone_3_lv_high (zone_3_lv_access_int) = log10 ( high (ii) / close (ii-1) );
		zone_3_lv_low (zone_3_lv_access_int) = log10 ( low (ii) / close (ii-1) );
		zone_3_lv_close (zone_3_lv_access_int) = log10 ( close (ii) / close (ii-1) );
		zone_3_lv_access_int = zone_3_lv_access_int + 1;	
		}

		else if ( close (ii-1) >= (moving_median (ii-1) + 3*moving_MAD (ii-1)) ) // zone_4
		{
		zone_4_lv_open (zone_4_lv_access_int) = log10 ( open (ii) / close (ii-1) );
		zone_4_lv_high (zone_4_lv_access_int) = log10 ( high (ii) / close (ii-1) );
		zone_4_lv_low (zone_4_lv_access_int) = log10 ( low (ii) / close (ii-1) );
		zone_4_lv_close (zone_4_lv_access_int) = log10 ( close (ii) / close (ii-1) );
		zone_4_lv_access_int = zone_4_lv_access_int + 1;	
		}

		else if ( close (ii-1) < moving_median (ii-1) & close (ii-1) >= (moving_median (ii-1) - moving_MAD (ii-1)) ) // zone_5
		{
		zone_5_lv_open (zone_5_lv_access_int) = log10 ( open (ii) / close (ii-1) );
		zone_5_lv_high (zone_5_lv_access_int) = log10 ( high (ii) / close (ii-1) );
		zone_5_lv_low (zone_5_lv_access_int) = log10 ( low (ii) / close (ii-1) );
		zone_5_lv_close (zone_5_lv_access_int) = log10 ( close (ii) / close (ii-1) );
		zone_5_lv_access_int = zone_5_lv_access_int + 1;	
		}

		else if ( close (ii-1) < (moving_median (ii-1) - moving_MAD (ii-1)) & close (ii-1) >= (moving_median (ii-1) - 2*moving_MAD (ii-1)) ) // zone_6
		{
		zone_6_lv_open (zone_6_lv_access_int) = log10 ( open (ii) / close (ii-1) );
		zone_6_lv_high (zone_6_lv_access_int) = log10 ( high (ii) / close (ii-1) );
		zone_6_lv_low (zone_6_lv_access_int) = log10 ( low (ii) / close (ii-1) );
		zone_6_lv_close (zone_6_lv_access_int) = log10 ( close (ii) / close (ii-1) );
		zone_6_lv_access_int = zone_6_lv_access_int + 1;	
		}

		else if ( close (ii-1) < (moving_median (ii-1) - 2*moving_MAD (ii-1)) & close (ii-1) >= (moving_median (ii-1) - 3*moving_MAD (ii-1)) ) // zone_7
		{
		zone_7_lv_open (zone_7_lv_access_int) = log10 ( open (ii) / close (ii-1) );
		zone_7_lv_high (zone_7_lv_access_int) = log10 ( high (ii) / close (ii-1) );
		zone_7_lv_low (zone_7_lv_access_int) = log10 ( low (ii) / close (ii-1) );
		zone_7_lv_close (zone_7_lv_access_int) = log10 ( close (ii) / close (ii-1) );
		zone_7_lv_access_int = zone_7_lv_access_int + 1;	
		}

		else // zone_8
		{
		zone_8_lv_open (zone_8_lv_access_int) = log10 ( open (ii) / close (ii-1) );
		zone_8_lv_high (zone_8_lv_access_int) = log10 ( high (ii) / close (ii-1) );
		zone_8_lv_low (zone_8_lv_access_int) = log10 ( low (ii) / close (ii-1) );
		zone_8_lv_close (zone_8_lv_access_int) = log10 ( close (ii) / close (ii-1) );
		zone_8_lv_access_int = zone_8_lv_access_int + 1;	
		}

	} // end of ( ii >= 21 & current_4_wma_true_range >= current_median_true_rang ) if conditional for bin selection based on volatility

    } // end of Code Block A loop

// the next Code Block B performs the MC randomisation routine
    MTRand mtrand1; // Declare the Mersenne Twister Class - will seed from system time
    int random_zone_access_int;

// first, reset the volatility calculations to those at the beginning of the original series
    for (octave_idx_type ii (0); ii < 21; ii++) // loop to fill the moving_true_range_window
    {
    moving_true_range_window (ii) = true_range (ii);
    } // end of loop to fill the moving_true_range_window
    // Note: the "&" below acts as Address-of operator: p = &x; Read: Assign to p (a pointer) the address of x.
    std::nth_element( &moving_true_range_window(0), &moving_true_range_window(10), &moving_true_range_window(21) );
    current_median_true_range = moving_true_range_window(10);
    current_4_wma_true_range = ( 4*true_range(20) + 3*true_range(19) + 2*true_range(18) + true_range(17) ) / 10 ;

// now the MC synthetic routine
    for (octave_idx_type ii (21); ii < args(0).length (); ii++) // Code Block B loop
    {

// identify where previous close is vis a vis previous median and MAD and create a new "current bar" by MC selection from relevant zone bin
	if ( current_4_wma_true_range >= current_median_true_range )
	{
		if ( close (ii-1) >= moving_median (ii-1) & close (ii-1) < (moving_median (ii-1) + moving_MAD (ii-1)) ) // zone_1
		{
		random_zone_access_int = mtrand1.randInt( zone_1_access_int - 1 ); // generate random access int

			if ( random_zone_access_int < 0 ) // check random access int doesn't exceed lower boundary for zone vector column  
			{
			random_zone_access_int = 0; 
			}
			if ( random_zone_access_int > zone_1_access_int - 1 ) // check random access int doesn't exceed upper boundary for zone vector column  
			{
			random_zone_access_int = zone_1_access_int - 1; 
			}

		open (ii) = ( floor( ( close (ii-1) * pow(10, zone_1_open (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		high (ii) = ( floor( ( close (ii-1) * pow(10, zone_1_high (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		low (ii) = ( floor( ( close (ii-1) * pow(10, zone_1_low (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		close (ii) = ( floor( ( close (ii-1) * pow(10, zone_1_close (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		}

		else if ( close (ii-1) >= (moving_median (ii-1) + moving_MAD (ii-1)) & close (ii-1) < (moving_median (ii-1) + 2*moving_MAD (ii-1)) ) // zone_2
		{
		random_zone_access_int = mtrand1.randInt( zone_2_access_int - 1 ); // generate random access int

			if ( random_zone_access_int < 0 ) // check random access int doesn't exceed lower boundary for zone vector column  
			{
			random_zone_access_int = 0; 
			}
			if ( random_zone_access_int > zone_2_access_int - 1 ) // check random access int doesn't exceed upper boundary for zone vector column  
			{
			random_zone_access_int = zone_2_access_int - 1; 
			}

		open (ii) = ( floor( ( close (ii-1) * pow(10, zone_2_open (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		high (ii) = ( floor( ( close (ii-1) * pow(10, zone_2_high (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		low (ii) = ( floor( ( close (ii-1) * pow(10, zone_2_low (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		close (ii) = ( floor( ( close (ii-1) * pow(10, zone_2_close (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		}

		else if ( close (ii-1) >= (moving_median (ii-1) + 2*moving_MAD (ii-1)) & close (ii-1) < (moving_median (ii-1) + 3*moving_MAD (ii-1)) ) // zone_3
		{
		random_zone_access_int = mtrand1.randInt( zone_3_access_int - 1 ); // generate random access int

			if ( random_zone_access_int < 0 ) // check random access int doesn't exceed lower boundary for zone vector column  
			{
			random_zone_access_int = 0; 
			}
			if ( random_zone_access_int > zone_3_access_int - 1 ) // check random access int doesn't exceed upper boundary for zone vector column  
			{
			random_zone_access_int = zone_3_access_int - 1; 
			}

		open (ii) = ( floor( ( close (ii-1) * pow(10, zone_3_open (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		high (ii) = ( floor( ( close (ii-1) * pow(10, zone_3_high (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		low (ii) = ( floor( ( close (ii-1) * pow(10, zone_3_low (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		close (ii) = ( floor( ( close (ii-1) * pow(10, zone_3_close (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		}

		else if ( close (ii-1) >= (moving_median (ii-1) + 3*moving_MAD (ii-1)) ) // zone_4
		{
		random_zone_access_int = mtrand1.randInt( zone_4_access_int - 1 ); // generate random access int

			if ( random_zone_access_int < 0 ) // check random access int doesn't exceed lower boundary for zone vector column  
			{
			random_zone_access_int = 0; 
			}
			if ( random_zone_access_int > zone_4_access_int - 1 ) // check random access int doesn't exceed upper boundary for zone vector column  
			{
			random_zone_access_int = zone_4_access_int - 1; 
			}

		open (ii) = ( floor( ( close (ii-1) * pow(10, zone_4_open (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		high (ii) = ( floor( ( close (ii-1) * pow(10, zone_4_high (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		low (ii) = ( floor( ( close (ii-1) * pow(10, zone_4_low (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		close (ii) = ( floor( ( close (ii-1) * pow(10, zone_4_close (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;	
		}

		else if ( close (ii-1) < moving_median (ii-1) & close (ii-1) >= (moving_median (ii-1) - moving_MAD (ii-1)) ) // zone_5
		{
		random_zone_access_int = mtrand1.randInt( zone_5_access_int - 1 ); // generate random access int

			if ( random_zone_access_int < 0 ) // check random access int doesn't exceed lower boundary for zone vector column  
			{
			random_zone_access_int = 0; 
			}
			if ( random_zone_access_int > zone_5_access_int - 1 ) // check random access int doesn't exceed upper boundary for zone vector column  
			{
			random_zone_access_int = zone_5_access_int - 1; 
			}

		open (ii) = ( floor( ( close (ii-1) * pow(10, zone_5_open (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		high (ii) = ( floor( ( close (ii-1) * pow(10, zone_5_high (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		low (ii) = ( floor( ( close (ii-1) * pow(10, zone_5_low (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		close (ii) = ( floor( ( close (ii-1) * pow(10, zone_5_close (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;	
		}

		else if ( close (ii-1) < (moving_median (ii-1) - moving_MAD (ii-1)) & close (ii-1) >= (moving_median (ii-1) - 2*moving_MAD (ii-1)) ) // zone_6
		{
		random_zone_access_int = mtrand1.randInt( zone_6_access_int - 1 ); // generate random access int

			if ( random_zone_access_int < 0 ) // check random access int doesn't exceed lower boundary for zone vector column  
			{
			random_zone_access_int = 0; 
			}
			if ( random_zone_access_int > zone_6_access_int - 1 ) // check random access int doesn't exceed upper boundary for zone vector column  
			{
			random_zone_access_int = zone_6_access_int - 1; 
			}

		open (ii) = ( floor( ( close (ii-1) * pow(10, zone_6_open (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		high (ii) = ( floor( ( close (ii-1) * pow(10, zone_6_high (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		low (ii) = ( floor( ( close (ii-1) * pow(10, zone_6_low (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		close (ii) = ( floor( ( close (ii-1) * pow(10, zone_6_close (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;	
		}

		else if ( close (ii-1) < (moving_median (ii-1) - 2*moving_MAD (ii-1)) & close (ii-1) >= (moving_median (ii-1) - 3*moving_MAD (ii-1)) ) // zone_7
		{
		random_zone_access_int = mtrand1.randInt( zone_7_access_int - 1 ); // generate random access int

			if ( random_zone_access_int < 0 ) // check random access int doesn't exceed lower boundary for zone vector column  
			{
			random_zone_access_int = 0; 
			}
			if ( random_zone_access_int > zone_7_access_int - 1 ) // check random access int doesn't exceed upper boundary for zone vector column  
			{
			random_zone_access_int = zone_7_access_int - 1; 
			}

		open (ii) = ( floor( ( close (ii-1) * pow(10, zone_7_open (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		high (ii) = ( floor( ( close (ii-1) * pow(10, zone_7_high (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		low (ii) = ( floor( ( close (ii-1) * pow(10, zone_7_low (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		close (ii) = ( floor( ( close (ii-1) * pow(10, zone_7_close (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;		
		}

		else // zone_8
		{
		random_zone_access_int = mtrand1.randInt( zone_8_access_int - 1 ); // generate random access int

			if ( random_zone_access_int < 0 ) // check random access int doesn't exceed lower boundary for zone vector column  
			{
			random_zone_access_int = 0; 
			}
			if ( random_zone_access_int > zone_8_access_int - 1 ) // check random access int doesn't exceed upper boundary for zone vector column  
			{
			random_zone_access_int = zone_8_access_int - 1; 
			}

		open (ii) = ( floor( ( close (ii-1) * pow(10, zone_8_open (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		high (ii) = ( floor( ( close (ii-1) * pow(10, zone_8_high (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		low (ii) = ( floor( ( close (ii-1) * pow(10, zone_8_low (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		close (ii) = ( floor( ( close (ii-1) * pow(10, zone_8_close (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;		
		}
	}
	else
	{
		if ( close (ii-1) >= moving_median (ii-1) & close (ii-1) < (moving_median (ii-1) + moving_MAD (ii-1)) ) // zone_1
		{
		random_zone_access_int = mtrand1.randInt( zone_1_lv_access_int - 1 ); // generate random access int

			if ( random_zone_access_int < 0 ) // check random access int doesn't exceed lower boundary for zone vector column  
			{
			random_zone_access_int = 0; 
			}
			if ( random_zone_access_int > zone_1_lv_access_int - 1 ) // check random access int doesn't exceed upper boundary for zone vector column  
			{
			random_zone_access_int = zone_1_lv_access_int - 1; 
			}

		open (ii) = ( floor( ( close (ii-1) * pow(10, zone_1_lv_open (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		high (ii) = ( floor( ( close (ii-1) * pow(10, zone_1_lv_high (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		low (ii) = ( floor( ( close (ii-1) * pow(10, zone_1_lv_low (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		close (ii) = ( floor( ( close (ii-1) * pow(10, zone_1_lv_close (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		}

		else if ( close (ii-1) >= (moving_median (ii-1) + moving_MAD (ii-1)) & close (ii-1) < (moving_median (ii-1) + 2*moving_MAD (ii-1)) ) // zone_2
		{
		random_zone_access_int = mtrand1.randInt( zone_2_lv_access_int - 1 ); // generate random access int

			if ( random_zone_access_int < 0 ) // check random access int doesn't exceed lower boundary for zone vector column  
			{
			random_zone_access_int = 0; 
			}
			if ( random_zone_access_int > zone_2_lv_access_int - 1 ) // check random access int doesn't exceed upper boundary for zone vector column  
			{
			random_zone_access_int = zone_2_lv_access_int - 1; 
			}

		open (ii) = ( floor( ( close (ii-1) * pow(10, zone_2_lv_open (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		high (ii) = ( floor( ( close (ii-1) * pow(10, zone_2_lv_high (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		low (ii) = ( floor( ( close (ii-1) * pow(10, zone_2_lv_low (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		close (ii) = ( floor( ( close (ii-1) * pow(10, zone_2_lv_close (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		}

		else if ( close (ii-1) >= (moving_median (ii-1) + 2*moving_MAD (ii-1)) & close (ii-1) < (moving_median (ii-1) + 3*moving_MAD (ii-1)) ) // zone_3
		{
		random_zone_access_int = mtrand1.randInt( zone_3_lv_access_int - 1 ); // generate random access int

			if ( random_zone_access_int < 0 ) // check random access int doesn't exceed lower boundary for zone vector column  
			{
			random_zone_access_int = 0; 
			}
			if ( random_zone_access_int > zone_3_lv_access_int - 1 ) // check random access int doesn't exceed upper boundary for zone vector column  
			{
			random_zone_access_int = zone_3_lv_access_int - 1; 
			}

		open (ii) = ( floor( ( close (ii-1) * pow(10, zone_3_lv_open (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		high (ii) = ( floor( ( close (ii-1) * pow(10, zone_3_lv_high (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		low (ii) = ( floor( ( close (ii-1) * pow(10, zone_3_lv_low (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		close (ii) = ( floor( ( close (ii-1) * pow(10, zone_3_lv_close (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		}

		else if ( close (ii-1) >= (moving_median (ii-1) + 3*moving_MAD (ii-1)) ) // zone_4
		{
		random_zone_access_int = mtrand1.randInt( zone_4_lv_access_int - 1 ); // generate random access int

			if ( random_zone_access_int < 0 ) // check random access int doesn't exceed lower boundary for zone vector column  
			{
			random_zone_access_int = 0; 
			}
			if ( random_zone_access_int > zone_4_lv_access_int - 1 ) // check random access int doesn't exceed upper boundary for zone vector column  
			{
			random_zone_access_int = zone_4_lv_access_int - 1; 
			}

		open (ii) = ( floor( ( close (ii-1) * pow(10, zone_4_lv_open (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		high (ii) = ( floor( ( close (ii-1) * pow(10, zone_4_lv_high (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		low (ii) = ( floor( ( close (ii-1) * pow(10, zone_4_lv_low (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		close (ii) = ( floor( ( close (ii-1) * pow(10, zone_4_lv_close (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;	
		}

		else if ( close (ii-1) < moving_median (ii-1) & close (ii-1) >= (moving_median (ii-1) - moving_MAD (ii-1)) ) // zone_5
		{
		random_zone_access_int = mtrand1.randInt( zone_5_lv_access_int - 1 ); // generate random access int

			if ( random_zone_access_int < 0 ) // check random access int doesn't exceed lower boundary for zone vector column  
			{
			random_zone_access_int = 0; 
			}
			if ( random_zone_access_int > zone_5_lv_access_int - 1 ) // check random access int doesn't exceed upper boundary for zone vector column  
			{
			random_zone_access_int = zone_5_lv_access_int - 1; 
			}

		open (ii) = ( floor( ( close (ii-1) * pow(10, zone_5_lv_open (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		high (ii) = ( floor( ( close (ii-1) * pow(10, zone_5_lv_high (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		low (ii) = ( floor( ( close (ii-1) * pow(10, zone_5_lv_low (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		close (ii) = ( floor( ( close (ii-1) * pow(10, zone_5_lv_close (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;	
		}

		else if ( close (ii-1) < (moving_median (ii-1) - moving_MAD (ii-1)) & close (ii-1) >= (moving_median (ii-1) - 2*moving_MAD (ii-1)) ) // zone_6
		{
		random_zone_access_int = mtrand1.randInt( zone_6_lv_access_int - 1 ); // generate random access int

			if ( random_zone_access_int < 0 ) // check random access int doesn't exceed lower boundary for zone vector column  
			{
			random_zone_access_int = 0; 
			}
			if ( random_zone_access_int > zone_6_lv_access_int - 1 ) // check random access int doesn't exceed upper boundary for zone vector column  
			{
			random_zone_access_int = zone_6_lv_access_int - 1; 
			}

		open (ii) = ( floor( ( close (ii-1) * pow(10, zone_6_lv_open (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		high (ii) = ( floor( ( close (ii-1) * pow(10, zone_6_lv_high (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		low (ii) = ( floor( ( close (ii-1) * pow(10, zone_6_lv_low (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		close (ii) = ( floor( ( close (ii-1) * pow(10, zone_6_lv_close (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;	
		}

		else if ( close (ii-1) < (moving_median (ii-1) - 2*moving_MAD (ii-1)) & close (ii-1) >= (moving_median (ii-1) - 3*moving_MAD (ii-1)) ) // zone_7
		{
		random_zone_access_int = mtrand1.randInt( zone_7_lv_access_int - 1 ); // generate random access int

			if ( random_zone_access_int < 0 ) // check random access int doesn't exceed lower boundary for zone vector column  
			{
			random_zone_access_int = 0; 
			}
			if ( random_zone_access_int > zone_7_lv_access_int - 1 ) // check random access int doesn't exceed upper boundary for zone vector column  
			{
			random_zone_access_int = zone_7_lv_access_int - 1; 
			}

		open (ii) = ( floor( ( close (ii-1) * pow(10, zone_7_lv_open (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		high (ii) = ( floor( ( close (ii-1) * pow(10, zone_7_lv_high (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		low (ii) = ( floor( ( close (ii-1) * pow(10, zone_7_lv_low (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		close (ii) = ( floor( ( close (ii-1) * pow(10, zone_7_lv_close (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;		
		}

		else // zone_8
		{
		random_zone_access_int = mtrand1.randInt( zone_8_lv_access_int - 1 ); // generate random access int

			if ( random_zone_access_int < 0 ) // check random access int doesn't exceed lower boundary for zone vector column  
			{
			random_zone_access_int = 0; 
			}
			if ( random_zone_access_int > zone_8_lv_access_int - 1 ) // check random access int doesn't exceed upper boundary for zone vector column  
			{
			random_zone_access_int = zone_8_lv_access_int - 1; 
			}

		open (ii) = ( floor( ( close (ii-1) * pow(10, zone_8_lv_open (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		high (ii) = ( floor( ( close (ii-1) * pow(10, zone_8_lv_high (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		low (ii) = ( floor( ( close (ii-1) * pow(10, zone_8_lv_low (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;
		close (ii) = ( floor( ( close (ii-1) * pow(10, zone_8_lv_close (random_zone_access_int)) ) / tick + 0.5 ) ) * tick;		
		}

	} // end of ( current_4_wma_true_range >= current_median_true_rang ) if conditional 

// create new "current bar" moving_median, moving MAD values and volatility values, overloading previous code 

	for (octave_idx_type jj (0); jj < 21; jj++) // loop to fill the moving_median_window
	{
        moving_median_window (jj) = close (ii-jj);
   	} // end of loop to fill the moving_median_window
        // Note: the "&" below acts as Address-of operator: p = &x; Read: Assign to p (a pointer) the address of x.
        std::nth_element( &moving_median_window(0), &moving_median_window(10), &moving_median_window(21) );
        moving_median (ii) = moving_median_window(10);

	for (octave_idx_type jj (0); jj < 21; jj++) // loop to fill the moving_MAD_window
	{
        moving_MAD_window (jj) = fabs( close (ii-jj) - moving_median (ii) );
   	} // end of loop to fill the moving_MAD_window
        // Note: the "&" below acts as Address-of operator: p = &x; Read: Assign to p (a pointer) the address of x.
        std::nth_element( &moving_MAD_window(0), &moving_MAD_window(10), &moving_MAD_window(21) );
        moving_MAD (ii) = moving_MAD_window(10);

        // true range calculations
        true_range (ii) = fmax ( fmax(high(ii)-low(ii),fabs(high(ii)-close(ii-1))) , fabs(low(ii)-close(ii-1)) );
	for (octave_idx_type jj (0); jj < 21; jj++) // loop to fill the moving_true_range_window
	{
        moving_true_range_window (jj) = true_range (ii-jj);
   	} // end of loop to fill the moving_true_range_window
        // Note: the "&" below acts as Address-of operator: p = &x; Read: Assign to p (a pointer) the address of x.
        std::nth_element( &moving_true_range_window(0), &moving_true_range_window(10), &moving_true_range_window(21) );
        current_median_true_range = moving_true_range_window(10);
        current_4_wma_true_range = ( 4*true_range(ii) + 3*true_range(ii-1) + 2*true_range(ii-2) + true_range(ii-3) ) / 10 ;

    } // end of Code Block B loop

// if compensation due to negative or zero original price values due to continuous back- adjusting of price series took place, re-adjust
    if ( correction_factor != 0.0 )
    {
	for (octave_idx_type ii (0); ii < args(0).length (); ii++)
	{
	open (ii) = open (ii) - correction_factor;
	high (ii) = high (ii) - correction_factor;
	low (ii) = low (ii) - correction_factor;
	close (ii) = close (ii) - correction_factor;
	moving_median (ii) = moving_median (ii) - correction_factor;
	moving_MAD (ii) = moving_MAD (ii) - correction_factor;
	}
    }

retval_list(5) = moving_MAD;
retval_list(4) = moving_median;
retval_list(3) = close;
retval_list(2) = low;
retval_list(1) = high;
retval_list(0) = open;
return retval_list;
} 

A final thought: although not implemented in the above code it would be possible to apply some form of "quality control" to the output. Statistical measures of the input time series could be taken and thresholds established and only those synthetic outputs that fall within these threshold conditions could be accepted as a valid synthetic time series output.

Below is a screenshot of a time series and synthetic data generated from it using the above function code. For the moment I won't say which is the original and which is the synthetic data - perhaps readers would like to post their guesses as comments?

Update on the Trend Vigor indicator

I have now completed some preliminary Monte Carlo testing of the Trend Vigor indicator and I thought this would be a good opportunity to share with readers the general approach I have been taking when talking about my Monte Carlo testing on predefined, "ideal" market types.

Firstly, I create my "ideal" market types by using an Octave .oct function, the code for which is given below. In this code can also be seen my implementation of the code for the Trend Vigor indicator, which I think is slightly different from Elher's. The code is commented, so no further description is required here.

.oct function code

// This function takes as arguments a single value input for a sinewave period and a single value input for a degrees increment value. A sinewave of the 
// given period is created and then trends are added to the sinewave such that 5 hypothesised "ideal" market types are created: these are

// 1) a perfectly cyclic sideways market with no trend i.e. just the sinewave component
// 2) an uptrending market with cyclic retracements (uwr) such that retracements are down to the 50% Fibonacci ratio level
// 3) an uptrending market with no retracemnets (unr) i.e. the uptrend completely swamps the downward cyclic component counter to the trend
// 4) a downtrending market with cyclic retracements (dwr) such that retracements are up to the 50% Fibonacci ratio level
// 5) a downtrending market with no retracements (dnr) i.e. the downtrend completely swamps the upward cyclic component counter to the trend

// The cyclic component of these markets is then extracted using the bandpass indicator. These vector lengths are 500 in order to allow time for the 
// bandpass calculations to settle down. The peak to peak amplitude is then recovered from the bandpass at the end of each market vector. The trend slope 
// is also calculated at the end of each market vector. The trend vigor indicator value is then calculated thus
// trend_slope/p_to_p. The bandpass delta is set to 0.2.

// The idea is that 
// for a sideways market the value should be about 0
// for a trending with retracement market the value should be > 0 && < 1 or < 0 && > -1, depending on whether an up trend or down trend
// for a trending with no retracement market the ratio should be > 1  or > -1, depending on whether an up trend or down trend 

// The original sinewave is then repeatedly phase shifted by the degrees increment value, the above process repeated and new values are calculated. 
// All the calculated values for each market type are the vector outputs of the function. 

#include 
#include 
#include 
#define PI 3.14159265

DEFUN_DLD (trend_vigor_dist, args, , "Inputs are period & degrees increment value, outputs are vectors of repeated median slope values")
{
    octave_value_list retval_list;

    if (args(0).length () < 1 | args(0).length () > 1)
    {
        error ("Invalid arguments. Inputs are single value period length and single value degrees increment value");
        return retval_list;
    }

    if (args(1).length () < 1 | args(1).length () > 1)
    {
        error ("Invalid arguments. Inputs are single value period length and single value degrees increment value");
        return retval_list;
    }

    if (error_state)
    {
        error ("Invalid arguments. Inputs are single value period length and single value degrees increment value");
        return retval_list;
    }
    // read inputs
    int period = args(0).int_value ();
    double degrees_inc = args(1).double_value ();
    double period_inc = 360.0 / double(period);
    int length = period + 1; // the length of the "lookback". Is equal to period + 1 to encompass a full period

    // vectors to hold created market values
    ColumnVector sideways_vec(500); // vector to hold sideways market values
    ColumnVector sideways_bandpass_vec(500); // vector to hold bandpass of sideways market values

    ColumnVector uwr_vec(500); // vector to hold uwr market values
    ColumnVector uwr_bandpass_vec(500); // vector to hold bandpass of uwr market values
 
    ColumnVector unr_vec(500); // vector to hold unr market values
    ColumnVector unr_bandpass_vec(500); // vector to hold bandpass of unr market values

    ColumnVector dwr_vec(500); // vector to hold dwr market values
    ColumnVector dwr_bandpass_vec(500); // vector to hold bandpass of dwr market values
 
    ColumnVector dnr_vec(500); // vector to hold dnr market values
    ColumnVector dnr_bandpass_vec(500); // vector to hold bandpass of dnr market values

    // calculate the market trend_incs
    double uwr_trend_inc = 12 / ( 5 * double(period) );
    double unr_trend_inc = 4 / double(period);
    double dwr_trend_inc = -( 12 / ( 5 * double(period) ) );
    double dnr_trend_inc = -( 4 / double(period) );

    // declare variables for bandpass and trend vigor calculations
    double delta = 0.2;
    double beta = cos( (360.0/period)*PI/180.0 );
    double gamma = 1.0 / cos( (720.0*delta/period)*PI/180.0 );
    double alpha = gamma - sqrt(gamma*gamma - 1.0); 
    double power_side;
    double power_uwr;
    double power_unr;
    double power_dwr;
    double power_dnr;
    double rms;
    double p_to_p;

    // create output vectors
    int output_vec_length = int ( 360 / degrees_inc );
    ColumnVector sideways_dist ( output_vec_length ); // create output column of correct length for sideways market
    ColumnVector uwr_dist  ( output_vec_length ); // create output column of correct length for uwr market
    ColumnVector unr_dist  ( output_vec_length ); // create output column of correct length for unr market
    ColumnVector dwr_dist  ( output_vec_length ); // create output column of correct length for dwr market
    ColumnVector dnr_dist  ( output_vec_length ); // create output column of correct length for dnr market

    for (octave_idx_type ii (0); ii < output_vec_length; ii++) 
    {

        // Create the market types and their bandpasses for this ii iteration
 for (octave_idx_type jj (0); jj < 500; jj++) 
        {

        // First create the sideways market type
        sideways_vec(jj) = sin( (degrees_inc*ii + period_inc*jj) * PI / 180 );
                
                if ( jj < 2 )
                {
                sideways_bandpass_vec(jj) = 0;
                }
                else
                {
                sideways_bandpass_vec(jj) = 0.5*(1.0 - alpha)*(sideways_vec(jj) - sideways_vec(jj-2)) + beta*(1.0 + alpha)*sideways_bandpass_vec(jj-1) - alpha*sideways_bandpass_vec(jj-2);
                }

        // next, the uwr retracement market (uwr)
        uwr_vec(jj) = sideways_vec(jj) + jj*uwr_trend_inc;

                if ( jj < 2 )
                {
                uwr_bandpass_vec(jj) = 0;
                }
                else
                {
                uwr_bandpass_vec(jj) = 0.5*(1.0 - alpha)*(uwr_vec(jj) - uwr_vec(jj-2)) + beta*(1.0 + alpha)*uwr_bandpass_vec(jj-1) - alpha*uwr_bandpass_vec(jj-2);
                }

        // next, the unr retracement market (unr)
        unr_vec(jj) = sideways_vec(jj) + jj*unr_trend_inc;

                if ( jj < 2 )
                {
                unr_bandpass_vec(jj) = 0;
                }
                else
                {
                unr_bandpass_vec(jj) = 0.5*(1.0 - alpha)*(unr_vec(jj) - unr_vec(jj-2)) + beta*(1.0 + alpha)*unr_bandpass_vec(jj-1) - alpha*unr_bandpass_vec(jj-2);
                }

        // next, the dwr retracement market (dwr)
        dwr_vec(jj) = sideways_vec(jj) + jj*dwr_trend_inc;

                if ( jj < 2 )
                {
                dwr_bandpass_vec(jj) = 0;
                }
                else
                {
                dwr_bandpass_vec(jj) = 0.5*(1.0 - alpha)*(dwr_vec(jj) - dwr_vec(jj-2)) + beta*(1.0 + alpha)*dwr_bandpass_vec(jj-1) - alpha*dwr_bandpass_vec(jj-2);
                }

        // next, the dnr retracement market (dnr)
        dnr_vec(jj) = sideways_vec(jj) + jj*dnr_trend_inc;

                if ( jj < 2 )
                {
                dnr_bandpass_vec(jj) = 0;
                }
                else
                {
                dnr_bandpass_vec(jj) = 0.5*(1.0 - alpha)*(dnr_vec(jj) - dnr_vec(jj-2)) + beta*(1.0 + alpha)*dnr_bandpass_vec(jj-1) - alpha*dnr_bandpass_vec(jj-2);
                }

        } // end of jj loop to create the different markets and their bandpasses

    // now loop over end of each market vector to create the distributions

     power_side = 0.0;
     power_uwr = 0.0;
     power_unr = 0.0;
     power_dwr = 0.0;
     power_dnr = 0.0;

     for (octave_idx_type jj (0); jj < length; jj++)
         {
         power_side = power_side + sideways_bandpass_vec(499-jj)*sideways_bandpass_vec(499-jj) + sideways_bandpass_vec(499-jj-int(period/4.0))*sideways_bandpass_vec(499-jj-int(period/4.0)) ;
         power_uwr = power_uwr + uwr_bandpass_vec(499-jj)*uwr_bandpass_vec(499-jj) + uwr_bandpass_vec(499-jj-int(period/4.0))*uwr_bandpass_vec(499-jj-int(period/4.0)) ;
         power_unr = power_unr + unr_bandpass_vec(499-jj)*unr_bandpass_vec(499-jj) + unr_bandpass_vec(499-jj-int(period/4.0))*unr_bandpass_vec(499-jj-int(period/4.0)) ;
         power_dwr = power_dwr + dwr_bandpass_vec(499-jj)*dwr_bandpass_vec(499-jj) + dwr_bandpass_vec(499-jj-int(period/4.0))*dwr_bandpass_vec(499-jj-int(period/4.0)) ;
         power_dnr = power_dnr + dnr_bandpass_vec(499-jj)*dnr_bandpass_vec(499-jj) + dnr_bandpass_vec(499-jj-int(period/4.0))*dnr_bandpass_vec(499-jj-int(period/4.0)) ;
         }

     // fill the distribution vectors
     rms = sqrt( power_side / (period+1) ) ;
     p_to_p = 2.0 * 1.414 * rms ;
     sideways_dist(ii) = ( sideways_vec(499) - sideways_vec(499-period) ) / p_to_p ;

     rms = sqrt( power_uwr / (period+1) ) ;
     p_to_p = 2.0 * 1.414 * rms ;
     uwr_dist(ii) = ( uwr_vec(499) - uwr_vec(499-period) ) / p_to_p ;

     rms = sqrt( power_unr / (period+1) ) ;
     p_to_p = 2.0 * 1.414 * rms ;
     unr_dist(ii) = ( unr_vec(499) - unr_vec(499-period) ) / p_to_p ;

     rms = sqrt( power_dwr / (period+1) ) ;
     p_to_p = 2.0 * 1.414 * rms ;
     dwr_dist(ii) = ( dwr_vec(499) - dwr_vec(499-period) ) / p_to_p ;

     rms = sqrt( power_dnr / (period+1) ) ;
     p_to_p = 2.0 * 1.414 * rms ;
     dnr_dist(ii) = ( dnr_vec(499) - dnr_vec(499-period) ) / p_to_p ;

    } // end of main ii loop

    retval_list(4) = dnr_dist;
    retval_list(3) = dwr_dist;
    retval_list(2) = unr_dist;
    retval_list(1) = uwr_dist;
    retval_list(0) = sideways_dist;

    return retval_list;
}
 

This function is called by this simple Octave script

clear all

inc = input( "Enter phase increment: ");

for ii = 6:50

[sideways_dist,uwr_dist,unr_dist,dwr_dist,dnr_dist] = trend_vigor_dist(ii,inc);
A=[sideways_dist,uwr_dist,unr_dist,dwr_dist,dnr_dist];
file = strcat( int2str(ii),"_period_dist" );
dlmwrite(file,A)

endfor 

which writes the output of the tests to named files which are to be used for further analysis in R.

Firstly, using R, I wanted to see what the distribution of the results looks like, so this R script

rm(list=ls())
data <- as.matrix(read.csv(file="20_period_dist",head=FALSE,sep=,))
side <- density(data[,1])
uwr <- density(data[,2])
max_uwr_y <- max(uwr$y)
max_uwr_x <- max(uwr$x)
min_uwr_x <- min(uwr$x)
unr <- density(data[,3])
max_unr_y <- max(unr$y)
max_unr_x <- max(unr$x)
min_unr_x <- min(unr$x)
dwr <- density(data[,4])
max_dwr_y <- max(dwr$y)
max_dwr_x <- max(dwr$x)
min_dwr_x <- min(dwr$x)
dnr <- density(data[,5])
max_dnr_y <- max(dnr$y)
max_dnr_x <- max(dnr$x)
min_dnr_x <- min(dnr$x)
plot_max_y <- max(max_uwr_y,max_unr_y,max_dwr_y,max_dnr_y)
plot_max_x <- max(max_uwr_x,max_unr_x,max_dwr_x,max_dnr_x)
plot_min_x <- min(min_uwr_x,min_unr_x,min_dwr_x,min_dnr_x)
par(mfrow=c(2,1))
plot(uwr,xlim=c(plot_min_x,plot_max_x),ylim=c(0,plot_max_y),col="red")
par(new=TRUE)
plot(unr,xlim=c(plot_min_x,plot_max_x),ylim=c(0,plot_max_y),col="blue")
par(new=TRUE)
plot(dwr,xlim=c(plot_min_x,plot_max_x),ylim=c(0,plot_max_y),col="green")
par(new=TRUE)
plot(dnr,xlim=c(plot_min_x,plot_max_x),ylim=c(0,plot_max_y),col="black")
plot(side)
 

produces plots such as this output

where the top plot shows the distributions of the uwr, unr, dwr and dnr markets, and the lower plot the sideways market. In this particular case, the spread of each distribution is so narrow ( measured differences of the order of thousandths of a decimal place ) that I consider that for practical purposes the distributions can be treated as single values. This simple R boot strap script gets the average value of the distributions to be used as this single point value.

rm(list=ls())

data <- as.matrix(read.csv(file="trend_vigor_dist_results",head=FALSE,sep=,))
side <- data[,1]
uwr <- data[,2]
unr <- data[,3]
dwr <- data[,4]
dnr <- data[,5]
side_samples <- matrix(0,50000)
uwr_samples <- matrix(0,50000)
unr_samples <- matrix(0,50000)
dwr_samples <- matrix(0,50000)
dnr_samples <- matrix(0,50000)

for(ii in 1:50000) {
side_samples[ii] <- mean(sample(side, replace=T))
uwr_samples[ii] <- mean(sample(uwr, replace=T))
unr_samples[ii] <- mean(sample(unr, replace=T))
dwr_samples[ii] <- mean(sample(dwr, replace=T))
dnr_samples[ii] <- mean(sample(dnr, replace=T))
}

side_mean <- mean(side_samples)
uwr_mean <- mean(uwr_samples)
unr_mean <- mean(unr_samples)
dwr_mean <- mean(dwr_samples)
dnr_mean <- mean(dnr_samples) 

For readers' interest, the actual values are 0.829, 1.329, -0.829 and -1.329 with 0 for the sideways market.

This final plot is the same Natural Gas plot as in my previous post, but with the above values substituted for Ehler's default values of 1 and -1.

What I intend to do now is use these values as the means of normal distributions with varying standard deviations as inputs for my Naive Bayes classifier. Further Monte Carlo testing will be done such that values for the standard deviations are obtained that result in the classifier giving false classifications, when tested using the "ideal" markets code above, within acceptable limits, most probably a 5% classification error rate.

A "new" indicator for possible inclusion in my Naive Bayesian Classifier

I have recently come across another interesting indicator on Ehler's website, information about which is available for download in the seminars section. It is the Trend Vigor indicator and below is my Octave .oct function implementation of it, first shown on the continuous, back-adjusted Natural Gas contract

and next the EURUSD spot forex market, both daily bars.

My implementation is slightly different from Ehler's in that there is no smoothing of the trend slope prior to calculation of the indicator. My reason for this is that smoothing will change the probability distribution of the indicator values, and since my Naive Bayes Classifier uses Kernel Density Estimation of a Mixture model, I prefer not to "corrupt" this density calculation with the subjective choice of a smoothing algorithm. As before, all parameters for this part of the classifier will be determined by using Monte Carlo techniques on "ideal," synthetically modelled market prices.

In the Natural Gas chart above it can be seen that the latter half of the chart is mostly sideways action, as determined by the Trend Vigor indicator, whilst my current version of the classifier ( the colour coded candlesticks ) does not give a corresponding similar determination for the whole of this later half of the chart. This is also the case in the second EURUSD chart. It is my hope that including the Trend Vigor indicator as an input to the classifier will improve the classifier's classification ability.

More on this in due course.

Update on Metatrader 4 C++ DLL coding

Coding of my suite of indicators for the Metatrader 4 platform using C++ DLLs is on-going and below is a screen shot of my installation as it stands at the moment.

The indicators are the Cybercycle, with its 1, 2, 3 and 4 day leads in the sub plot, and the instantaneous trendline, in red, and the Tukey contol lines, blue, in the main chart.

One thing that I have come up against in this coding are some of the apparent limitations of the Metatrader 4 platform; for instance there is a limit to how many lines can be drawn, per indicator, in any chart window. Also I am not yet able to access the values of one indicator, in the sub plot for example, for plotting in the main chart, i.e. plotting arrows above/below the candlesticks to indicate where the leading Cybercycle functions cross each other. Either this requires some substantial "hack" or much more advanced Metaquote language programming skills than I possess at the moment.

Metatrader 4 C++ DLL Coding

I have finally managed to teach myself, with help from various online sources including this very helpful forum thread, how to code a C++ DLL for the Metatrader 4 trading platform. As practise I have coded a simple moving average indicator using a recursive algorithm, the code for which I make freely available as a PDF download on the Dekalog website here.

Now that I have learned this I can begin, as I had hoped I would be able to, to "drag and drop" my C++ .oct functions into the Metatrader platform for testing and trading using live intraday forex data.

Completion of Naive Bayesian Classifier coding

I am pleased to say that my .oct coding of the Naive Bayesian classifier is now complete. The purpose of the function is to take as inputs various measurements of the current state of the time series, and then for the classifier to classify the time series as being in one of the five following states:-

• trending sideways with a cyclic action
• trending upwards with 50% retracements of previous up-legs
• trending upwards with no retracements
• trending downwards with 50% retracemenst of previous down-legs
• trending downwards with no retracements

These classifications will determine the most appropriate trading approach to take given the current state of the "tradable," and two screen shots of the classifier are given below, rendered as a "paint bar" study in the upper candlestick chart.

This first image just shows three classifications; trending sideways with a cycle in cyan, trending either up or down with 50% retracement in green, and trending up or down without retracement when the bars are blue or red (blue when close > open, red when close < open). Additionally, the coloured triangles show when the Cybercycle leading functions in the first subgraph cross, giving a count down to the predicted cyclic highs and lows, the blue "0" triangle being the predicted actual high or low. The green "fp high" and red "fp low" lines are the respective full period channel highs and lows of the price series, with the dotted yellow line being the midpoint of the two.

The second subgraph predicts turning points based on zero line crosses of the Cybercycle in the first subgraph. Read more about this indicator here.

The third subgraph is a plot of the sine of the phase of the Cybercycle, with leading signals, superimposed over a stochastic of the price bars. I may discuss this particular indicator in more depth in a future post.

This second screen shot is similar to the first, except that the classifications for retracement and no retracement are separated out into upward and downwards (chart key uwr = up with retracement, unr = up no retracement etc.). In this graph the coloured triangles represent the leading function crosses of the sine in the third subgraph.


Personally I am pleased that the effort to produce this indicator (almost six weeks of Monte Carlo simulation and statistical analysis in R, plus C++ coding) has resulted in a useful addition to my stable of indicators. I think there could still be some tweaks/improvements/additions that could be added, but for the moment I will leave the indicator as it is. The next thing I have set myself to do is implement my collection of C++ .oct functions as DLLs for Metatrader 4 (maybe an EA?) on a demo account so that I can utilise/test them intraday on forex data. Hopefully this will turn out to be a relatively simple case of "dragging and dropping" the already written C++ code.