Temporal Clustering Times on Forex Majors Pairs

In the following code box there are the results from the temporal clustering routine of my last few posts on the four forex majors pairs of EUR_USD, GBP_USD, USD_CHF and USD_JPY.

###### EUR_USD 10 minute bars #######
## In the following order
## Both Delta turning point filter and "normal" TPF combined ##
## Delta turning point filter only ##
## "Normal" turning point filter only

###################### Monday ##############################################
K_opt == 8, ix values == 13  38    63    89     112    135    162    186    ## averaged over all 15 n_bars 1 to 15 inclusive
                         00  4:10  8:20  12:40  16:30  20:20  00:50  4:50
                         
K_opt == 8, ix values == 13 39 64 89 112 135 161 186                        ## averaged over n_bars 1 to 6 inclusive ( upto and include 1 hour )

K_opt == 5, ix_values == 21 60 97 134 175                                   ## averaged over n_bars 7 to 15 inclusive ( over 1 hour )
K == 6,     ix values == 21 59 94 125 158 184

K_opt == 11, ix values == 9 26 43 60 78 95 113 132 151 169 185              ## averaged over all 15 n_bars 1 to 15 inclusive

K_opt == 8, ix values == 13  36  61  86 111 136 161 186                     ## averaged over n_bars 1 to 6 inclusive ( upto and include 1 hour )

K_opt == 8, ix values == 13  34  61  87 110 137 164 187                     ## averaged over n_bars 7 to 15 inclusive ( over 1 hour )

K_opt == 8, ix values == 13  38  63  88 112 137 162 186                     ## averaged over all 15 n_bars 1 to 15 inclusive

K_opt == 10, ix values == 10  31  52  72  91 112 131 150 169 188            ## averaged over n_bars 1 to 6 inclusive ( upto and include 1 hour )

K_opt == 8, ix values == 12  35  62  88 112 137 164 187                     ## averaged over n_bars 7 to 15 inclusive ( over 1 hour )

###################### Tuesday #############################################
K_opt == 6, ix values == 131   169    206   244    283    322               ## averaged over all 15 n_bars 1 to 15 inclusive
                         19:40 02:00  8:10  14:30  21:00  03:30
                         
K_opt == 6, ix values == 131 170 207 245 284 323                            ## averaged over n_bars 1 to 6 inclusive ( upto and include 1 hour )

K_opt == 7, ix values == 131 168 206 243 274 305 330                        ## averaged over n_bars 7 to 15 inclusive ( over 1 hour )

K_opt == 11, ix values == 124 143 164 184 205 226 247 268 289 310 331       ## averaged over all 15 n_bars 1 to 15 inclusive

K_opt == 11, ix values == 124 144 164 185 204 225 246 267 288 309 332       ## averaged over n_bars 1 to 6 inclusive ( upto and include 1 hour ) 

K_opt == 7, ix values = 133 169 206 241 273 304 329                         ## averaged over n_bars 7 to 15 inclusive ( over 1 hour )

K_opt == 9, ix values == 127 152 175 202 228 253 278 305 330                ## averaged over all 15 n_bars 1 to 15 inclusive

K_opt == 9, ix values == 127 152 177 202 228 253 278 304 329                ## averaged over n_bars 1 to 6 inclusive ( upto and include 1 hour )

K_opt == 7, ix values == 132 168 205 242 273 304 329                        ## averaged over n_bars 7 to 15 inclusive ( over 1 hour )

###################### Wednesday ###########################################
K_opt == 6, ix values == 275    312    351    389    426    465             ## averaged over all 15 n_bars 1 to 15 inclusive
                         19:40  01:50  08:20  14:40  20:50  03:20
                         
K_opt == 6, ix values == 275 313 352 391 428 466                            ## averaged over n_bars 1 to 6 inclusive ( upto and include 1 hour )

K_opt == 6, ix values == 274 312 350 389 424 463                            ## averaged over n_bars 7 to 15 inclusive ( over 1 hour )

K_opt == 9, ix values == 272 299 322 347 372 397 422 449 474                ## averaged over all 15 n_bars 1 to 15 inclusive

K_opt == 11, ix values == 268 288 308 329 348 369 390 411 432 453 476       ## averaged over n_bars 1 to 6 inclusive ( upto and include 1 hour )

K_opt == 6, ix values == 275 312 351 388 424 463                            ## averaged over n_bars 7 to 15 inclusive ( over 1 hour )

K_opt == 9, ix values == 272 297 322 348 373 398 423 449 474                ## averaged over all 15 n_bars 1 to 15 inclusive 

K_opt == 9, ix values == 271 297 322 348 373 398 423 448 473                ## averaged over n_bars 1 to 6 inclusive ( upto and include 1 hour )

K_opt == 6, ix values == 276 311 350 389 426 465                            ## averaged over n_bars 7 to 15 inclusive ( over 1 hour )

####################### Thursday ###########################################
K_opt == 6, ix values == 420    457    495    532    570    609             ## averaged over all 15 n_bars 1 to 15 inclusive
                         19:50  02:00  08:20  14:30  20:50  03:20
                         
K_opt == 6, ix values == 420 457 494 531 570 610                            ## averaged over n_bars 1 to 6 inclusive ( upto and include 1 hour )

K_opt == 6, ix values == 420 457 495 532 568 607                            ## averaged over n_bars 7 to 15 inclusive ( over 1 hour )

K_opt == 9, ix values == 416 443 466 492 518 543 568 593 618                ## averaged over all 15 n_bars 1 to 15 inclusive 

K_opt == 10, ix values == 414 437 460 483 506 527 550 573 596 619           ## averaged over n_bars 1 to 6 inclusive ( upto and include 1 hour )

K_opt == 9, ix values == 416 443 466 493 520 543 568 595 618                ## averaged over n_bars 7 to 15 inclusive ( over 1 hour )

K_opt == 9, ix values == 415 440 465 492 518 543 568 593 618                ## averaged over all 15 n_bars 1 to 15 inclusive 

K_opt == 9, ix values ==  415 440 465 492 518 543 568 593 618               ## averaged over n_bars 1 to 6 inclusive ( upto and include 1 hour )

K_opt == 7, ix values == 420 457 494 529 561 592 617                        ## averaged over n_bars 7 to 15 inclusive ( over 1 hour )

####################### Friday #############################################
K_opt == 5, ix values == 564    599    635    670     703                   ## averaged over all 15 n_bars 1 to 15 inclusive
                         19:50  01:40  07:40  13:30   19:00
                         
K_opt == 6, ix values == 563 596 627 654 680 707                            ## averaged over n_bars 1 to 6 inclusive ( upto and include 1 hour ) 
K == 5,     ix values == 564 599 635 668 703

K_opt == 5, ix values == 564 601 639 674 705                                ## averaged over n_bars 7 to 15 inclusive ( over 1 hour )

K_opt == 9, ix values == 556 575 595 614 633 652 672 691 711                ## averaged over all 15 n_bars 1 to 15 inclusive

K_opt == 11, ix values == 554 570 587 602 619 634 651 667 682 698 713       ## averaged over n_bars 1 to 6 inclusive ( upto and include 1 hour )

K_opt == 9, ix values == 556 575 595 614 633 652 671 691 711                ## averaged over n_bars 1 to 6 inclusive ( upto and include 1 hour )

K_opt == 9, ix values == 556 575 596 613 634 652 672 691 711                ## averaged over n_bars 7 to 15 inclusive ( over 1 hour )

K_opt == 9, ix values == 556 575 594 613 633 652 672 691 710                ## averaged over all 15 n_bars 1 to 15 inclusive 

K_opt == 9, ix values == 556 575 594 613 634 653 672 691 710                ## averaged over n_bars 1 to 6 inclusive ( upto and include 1 hour )

K_opt == 5, ix values == 564 600 637 674 705                                ## averaged over n_bars 7 to 15 inclusive ( over 1 hour )

############################################################################

###### GBP_USD 10 minute bars #######
## In the following order
## Both Delta turning point filter and "normal" TPF combined ##

###################### Monday ##############################################
K_opt = 8, ix_values = 13    36    61    86     111    136    162    186    ## averaged over all 15 n_bars 1 to 15 inclusive
                       0:00  3:50  8:00  12:10  16:20  20:30  0:50   4:50

K_opt = 9, ix_values = 12  34  56  78  99 120 141 164 187                   ## averaged over n_bars 1 to 6 inclusive ( upto and include 1 hour ) 

K_opt = 8, ix_values = 12  35  61  86 110 136 163 186                       ## averaged over n_bars 7 to 15 inclusive ( over 1 hour )

###################### Tuesday #############################################
K_opt = 12, ix_values = 124    143    162   180   199   216   235   254   274   293   312   332     ## averaged over all 15 n_bars 1 to 15 inclusive
                        18:30  21:40  0:50  3:50  7:00  9:50  13:00 16:10 19:30 22:40 1:50  5:10

K_opt = 11, ix_values = 124 143 164 185 206 227 248 269 290 311 332         ## averaged over n_bars 1 to 6 inclusive ( upto and include 1 hour )  

K_opt = 9, ix_values = 128 154 177 205 230 254 279 307 330                  ## averaged over n_bars 7 to 15 inclusive ( over 1 hour )

###################### Wednesday ###########################################
K_opt = 11, ix_values = 269   290   311  331  352  373   394   415   434   455   476   ## averaged over all 15 n_bars 1 to 15 inclusive
                        18:40 22:10 1:40 5:00 8:30 12:00 15:30 19:00 22:10 1:40  5:10

K_opt = 11, ix_values = 269 289 310 330 351 372 393 413 434 455 476         ## averaged over n_bars 1 to 6 inclusive ( upto and include 1 hour ) 

K_opt = 8, ix_values = 275 310 341 367 394 422 451 475                      ## averaged over n_bars 7 to 15 inclusive ( over 1 hour )

###################### Thursday ############################################
K_opt = 9, ix_values = 415   440   465  492  517   542   568   594   618    ## averaged over all 15 n_bars 1 to 15 inclusive
                       19:00 23:10 3:20 7:50 12:00 16:10 20:30 0:50  4:50

K_opt = 9, ix_values = 415 440 465 491 517 542 568 593 618                  ## averaged over n_bars 1 to 6 inclusive ( upto and include 1 hour )  

K_opt = 9, ix_values = 416 441 464 492 519 542 569 596 619                  ## averaged over n_bars 7 to 15 inclusive ( over 1 hour )

###################### Friday ##############################################
K_opt = 9, ix_values = 557   576   595  614  633  652   671   690   711     ## averaged over all 15 n_bars 1 to 15 inclusive
                       18:40 21:50 1:00 4:10 7:20 10:30 13:40 16:50 20:20

K_opt = 9, ix_values = 557 576 595 614 633 652 671 691 711                  ## averaged over n_bars 1 to 6 inclusive ( upto and include 1 hour ) 

K_opt = 8, ix_values = 557 576 599 621 642 665 686 709                      ## averaged over n_bars 7 to 15 inclusive ( over 1 hour )

############################################################################

###### USD_CHF 10 minute bars #######
## In the following order
## Both Delta turning point filter and "normal" TPF combined ##

###################### Monday ##############################################
K_opt = 11, ix_values = 8      25    42   61   79    96    113   131   150   169  188   ## averaged over all 15 n_bars 1 to 15 inclusive
                        23:10  2:00  4:50 8:00 11:00 13:50 16:40 19:40 22:50 2:00 5:10

K_opt = 11, ix_values =  9  26  43  60  79  96 114 133 151 170 189          ## averaged over n_bars 1 to 6 inclusive ( upto and include 1 hour ) 

K_opt = 7, ix_values =  13  38  66  99 127 157 184                          ## averaged over n_bars 7 to 15 inclusive ( over 1 hour )

###################### Tuesday #############################################
K_opt = 9, ix_values = 127   152   177  202  228   253   279   306  330     ## averaged over all 15 n_bars 1 to 15 inclusive
                       19:00 23:10 3:20 7:30 11:50 16:00 20:20 0:50 4:50

K_opt = 11, ix_values = 124 144 165 185 204 225 246 267 288 309 331         ## averaged over n_bars 1 to 6 inclusive ( upto and include 1 hour )  

K_opt = 7, ix_values = 133 170 205 240 270 301 328                          ## averaged over n_bars 7 to 15 inclusive ( over 1 hour )

###################### Wednesday ###########################################
K_opt = 10, ix_values = 270   293   316  342  365   388   411   432   454  475  ## averaged over all 15 n_bars 1 to 15 inclusive
                        18:50 22:40 2:30 6:50 10:40 14:30 18:20 21:50 1:30 5:00

K_opt = 12, ix_values = 268 287 308 327 346 365 384 401 420 439 458 477     ## averaged over n_bars 1 to 6 inclusive ( upto and include 1 hour )  

K_opt = 7, ix_values = 276 313 349 383 414 444 471                          ## averaged over n_bars 7 to 15 inclusive ( over 1 hour )

###################### Thursday ############################################
K_opt = 11, ix_values = 413   432   452  471  491  512   533   554   575   598  619  ## averaged over all 15 n_bars 1 to 15 inclusive
                        18:40 21:50 1:10 4:20 7:40 11:10 14:40 18:10 21:40 1:30 5:00

K_opt = 12, ix_values = 412 431 450 469 488 507 526 545 563 582 601 621     ## averaged over n_bars 1 to 6 inclusive ( upto and include 1 hour ) 

K_opt = 9, ix_values = 415 440 463 491 518 543 570 597 619                  ## averaged over n_bars 7 to 15 inclusive ( over 1 hour )

###################### Friday ##############################################
K_opt = 9, ix_values = 557   576   596  615  634  653   672   691   710     ## averaged over all 15 n_bars 1 to 15 inclusive
                       18:40 21:50 1:10 4:20 7:30 10:40 13:50 17:00 20:10

K_opt = 9, ix_values = 556 575 595 614 633 652 671 690 710                  ## averaged over n_bars 1 to 6 inclusive ( upto and include 1 hour ) 

K_opt = 7, ix_values = 558 579 602 629 652 677 705                          ## averaged over n_bars 7 to 15 inclusive ( over 1 hour )                

############################################################################

###### USD_JPY 10 minute bars #######
## In the following order
## Both Delta turning point filter and "normal" TPF combined ##

###################### Monday ##############################################
K_opt = 12, ix_values = 8     24   41   58   73    90    107   124   141   158  173  190  ## averaged over all 15 n_bars 1 to 15 inclusive
                        23:10 1:50 4:40 7:30 10:00 12:50 15:40 18:30 21:20 0:10 2:40 5:30

K_opt = 12, ix_values = 8  24  41  56  73  90 107 124 141 158 173 190       ## averaged over n_bars 1 to 6 inclusive ( upto and include 1 hour )

K_opt = 5, ix_values = 20  60  99 136 175                                   ## averaged over n_bars 7 to 15 inclusive ( over 1 hour )

###################### Tuesday #############################################
K_opt = 9, ix_values = 128   154   179  204  229   254   279   306  331     ## averaged over all 15 n_bars 1 to 15 inclusive
                       19:10 23:30 3:40 7:50 12:00 16:10 20:20 0:50 5:00

K_opt = 9, ix_values = 128 153 178 203 228 254 279 305 330                  ## averaged over n_bars 1 to 6 inclusive ( upto and include 1 hour )

K_opt = 7, ix_values = 133 168 205 240 271 302 329                          ## averaged over n_bars 7 to 15 inclusive ( over 1 hour )

###################### Wednesday ###########################################
K_opt = 11, ix_values = 269   289   310  331  352  373   394   414   433   454  476  ## averaged over all 15 n_bars 1 to 15 inclusive
                        18:40 22:00 1:30 5:00 8:30 12:00 15:30 18:50 22:00 1:30 5:10

K_opt = 9, ix_values = 272 297 322 348 374 399 424 449 474                  ## averaged over n_bars 1 to 6 inclusive ( upto and include 1 hour )

K_opt = 10, ix_values = 269 288 309 331 352 376 398 423 450 475             ## averaged over n_bars 7 to 15 inclusive ( over 1 hour )

###################### Thursday ############################################
K_opt = 9, ix_values = 416   442   467  492  518   543   568   593  618     ## averaged over all 15 n_bars 1 to 15 inclusive
                       19:10 23:30 3:40 7:50 12:10 16:20 20:30 0:40 4:50

K_opt = 12, ix_values = 412 431 450 469 488 507 526 545 564 583 602 621     ## averaged over n_bars 1 to 6 inclusive ( upto and include 1 hour )

K_opt = 7, ix_values = 420 455 492 527 560 591 618                          ## averaged over n_bars 7 to 15 inclusive ( over 1 hour )

###################### Friday ##############################################
K_opt = 7, 8 or 9
ix_values 7 = 561   588   613       638        663   686    709             ## averaged over all 15 n_bars 1 to 15 inclusive
ix_values 8 = 557   578   599  622  643        666   687    710
ix_values 9 = 557   576   596  616  635  653   672   691    711             ## timings are for this bottom row
              18:40 21:50 1:10 4:30 7:40 10:40 13:50 17:00  20:20

K_opt = 8, ix_values = 558 579 600 621 644 665 687 709                      ## averaged over n_bars 1 to 6 inclusive ( upto and include 1 hour )

K_opt = 6, ix_values = 563 594 621 646 676 705                              ## averaged over n_bars 7 to 15 inclusive ( over 1 hour )

############################################################################

This is based on 10 minute bars over the last year or so. Readers should read my last few previous posts for background.

The first set of results, EUR_USD, are what the charts of my previous posts were based on and include combined results of my "Delta Turning Point Filter" and "Normal Turning Point Filter" and the results for each filter separately. Since there doesn't appear to be significant differences between these, the other three pairs' results are the combined filter results only.

The K_opt variable is the optimal number of clusters (see my temporal-clustering-part-3 post for how "optimal" is decided) and the ix_values are also described in this post. For convenience the first set of ix_values per day have the relevant times anotated underneath and therefore it is a simple matter to count forwards/backwards in 10 minute increments to place times to the other ix_values. The variable n_bars is an input to the turning point filter functions and essentially indicates the lookback/lookforward period (n_bar == 2 would mean 2 x 10 minute periods) used for determining a local high/low according to each function's logic.

As to how to interpret this, a typical sequence of times per day might look like this:

18:40 22:00 1:30 5:00 8:30 12:00 15:30 18:50 22:00 1:30 5:10

where the highlighted times represent the BST times for the period covering the London session open to the New York session close for one day. The preceding and following times are the two "book-ending" Asian sessions. 

Close inspection of these results reveals some surprising regularities. In even just the above single example (an actual copy and paste of a code box example) there appear to be definite times per day at which a local high/low occurs. I hopefully will be able to incorporate this into some type of chart for a nice visual presentation of the data. 

More in due course. Enjoy.

Temporal Clustering on Real Prices, Part 2

Below are some more out of sample plots for the Temporal Clustering solutions of the EUR_USD forex pair for the week just gone. The details of how these solutions are derived is explained in my previous post, Temporal Clustering on Real Prices. First is Tuesday's solution

where the major (blue vertical lines) turns are a combination of optimal K values of 6 and 7 (5 sets of data in total) plus 2 sets of data each for K = 9 and 11 (red and green vertical lines). The price plot is

Next up is Wednesday's solution

where the blue vertical lines represent 5 sets of data with K = 6 and the red and green vertical lines 3 sets and 1 set with K = 9 and K= 11 respectively. The price plot is

Thursday's solution is

where black/blue vertical lines are K values of 9 and 6 respectively, whilst green/red are K values 10 and 7. Thursday's price plot is

Finally, Friday's solution is

where the major blue vertical lines are K = 9 over 5 sets of data, with the remainder being K = 5, 6 and 11 over the last 4 sets of data. Friday's price plot is

The above seems to tie in nicely with my previous post about Forex Intraday Seasonality whereby the above identified turning points signify the end points of said intraday tendencies to trend. Readers might also be interested in another paper I have come across, Segmentation and Time-of-Day Patterns in Foreign Exchange Markets, which gives a possible, theoretical explanation as to why such patterns manifest themselves. In particular, for the EUR_USD pair, the paper states  

• "the US dollar appreciates significantly from 8:00 to 12:00 GMT
  and the euro appreciates significantly from 16:00 to 22:00 GMT"
  

Readers can judge for themselves whether this appears to be true, out of sample, by inspecting the above plots. Enjoy!

Temporal Clustering on Real Prices

Having now had time to run the code shown in my previous post, Temporal Clustering, part 3, in this post I want to show the results on real prices.

Firstly, I have written two functions in Octave to identify market turning points and each function takes as input an n_bar argument which determines the lookback/lookforward length along price series to determine local relative highs and lows. I ran both these for n_bar values of 1 to 15 inclusive on EUR_USD forex 10 minute bars from July 2012 upto and including last week's set of 10 minute bars. I created 3 sets of turning point data per function by averaging the function outputs over n_bar 1 - 15, 1 - 6 and 7 - 15, and also averaged the outputs over the average of the 2 functions over the same ranges. In total this gives 9 slightly different sets of turning point data.

I then ran the optimal K clustering code, shown in previous posts, over each set of data to get the "solutions" per set of data. Six of the sets had an optimal K value of 8 and a combined plot of these is shown below.

For each "solution" turning point ix (ix ranges from 1 to 198) a turning point value of 1 is added to get a sort of spike train plot through time. The ix = 1 value is 22:00 BST on Sunday and ix = 198 is 06:50 BST on Tuesday. I chose this range so that there would be a buffer at each end of the time range I am really interested in: 7:00 BST to 22:00 BST, which covers the time from the London open to the New York close. The vertical blue lines are plotted for clarity to help identify the the turns and are plotted as 3 consecutive lines 10 minutes apart. The added text shows the time of occurence of the first bar of each triplet of lines, the time being London BST. The following second plot is the same as above but with the other 3 "solutions" of K = 5, 10 and 11 added.

For those readers who are familiar with the Delta Phenomenon the main vertical blue lines could conceptually be thought of as MTD lines with the other lines being lower timeframe ITD lines, but on an intraday scale. However, it is important to bear in mind that this is NOT a Delta solution and therefore rules about numbering, alternating highs and lows and inversions etc. do not apply. It is more helpful to think in terms of probability and see the various spikes/lines as indicating times of the day at which there is a higher probability of price making a local high or low. The size of a move after such a high or low is not indicated, and the timings are only approximate or alternatively represent the centre of a window in which the high or low might occur.

The proof of the pudding is in the eating, however, and the following plots are yesterday's (23 November 2020) out of sample EUR_USD forex pair price action with the lines of the above "solution" overlaid. The first plot is just the K = 8 solution plot

whilst this second plot has all lines shown.

Given the above caveats about caution with regards to the lines only being probabilities, it seems uncanny how accurately the major highs and lows of the day are picked out. I only wish I had done this analysis sooner as then yesterday could have been one of my best trading days ever!

More soon.

Temporal Clustering, Part 3

Continuing on with the subject matter of my last post, in the code box below there is R code which is a straight forward refactoring of the Octave code contained in the second code box of my last post. This code is my implementation of the cross validation routine described in the paper Cluster Validation by Prediction Strength, but adapted for use in the one dimensional case. I have refactored this into R code so that I can use the Ckmeans.1d.dp package for optimal, one dimensional clustering.

library( Ckmeans.1d.dp )

## load the training data from Octave output (comment out as necessary )
data = read.csv( "~/path/to//all_data_matrix" , header = FALSE )

## comment out as necessary
adjust = 0 ## default adjust value
sum_seq = seq( from = 1 , to = 198 , by = 1 ) ; adjust = 1 ; sum_seq_l = as.numeric( length( sum_seq ) )## Monday
##sum_seq = seq( from = 115 , to = 342 , by = 1 ) ; sum_seq_l = as.numeric( length( sum_seq ) ) ## Tuesday
##sum_seq = seq( from = 115 , to = 342 , by = 1 ) ; sum_seq_l = as.numeric( length( sum_seq ) ) ## Wednesday
##sum_seq = seq( from = 115 , to = 342 , by = 1 ) ; sum_seq_l = as.numeric( length( sum_seq ) ) ## Thursday
##sum_seq = seq( from = 547 , to = 720 , by = 1 ) ; adjust = 2 ; sum_seq_l = as.numeric( length( sum_seq ) ) ## Friday

## intraday --- commnet out or adjust as necessary
##sum_seq = seq( from = 25 , to = 100 , by = 1 ) ; sum_seq_l = as.numeric( length( sum_seq ) )

upper_tri_mask = 1 * upper.tri( matrix( 0L , nrow = sum_seq_l , ncol = sum_seq_l ) , diag = FALSE )
no_sample_iters = 1000
max_K = 20
all_k_ps = matrix( 0L , nrow = 1 , ncol = max_K )

for ( iters in 1 : no_sample_iters ) {

## sample the data in data by rows
train_ix = sample( nrow( data ) , size = round( nrow( data ) / 2 ) , replace = FALSE )
train_data = data[ train_ix , sum_seq ] ## extract training data using train_ix rows of data
train_data_sum = colSums( train_data )  ## sum down the columns of train_data
test_data = data[ -train_ix , sum_seq ] ## extract test data using NOT train_ix rows of data
test_data_sum = colSums( test_data )    ## sum down the columns of test_data
## adjust for weekend if necessary
if ( adjust == 1 ) { ## Monday, so correct artifacts of weekend gap
  train_data_sum[ 1 : 5 ] = mean( train_data_sum[ 1 : 48 ] )
  test_data_sum[ 1 : 5 ] = mean( test_data_sum[ 1 : 48 ] )   
} else if ( adjust == 2 ) { ## Friday, so correct artifacts of weekend gap
  train_data_sum[ ( sum_seq_l - 4 ) : sum_seq_l ] = mean( train_data_sum[ ( sum_seq_l - 47 ) : sum_seq_l ] )
  test_data_sum[  ( sum_seq_l - 4 ) : sum_seq_l ] = mean( test_data_sum[ ( sum_seq_l - 47 ) : sum_seq_l ] ) 
}

for ( k in 1 : max_K ) {
  
## K segment train_data_sum
train_res = Ckmeans.1d.dp( sum_seq , k , train_data_sum )
train_out_pairs_mat = matrix( 0L , nrow = sum_seq_l , ncol = sum_seq_l )

## K segment test_data_sum
test_res = Ckmeans.1d.dp( sum_seq , k , test_data_sum )
test_out_pairs_mat = matrix( 0L , nrow = sum_seq_l , ncol = sum_seq_l )

  for ( ii in 1 : length( train_res$centers ) ) {
    ix = which( train_res$cluster == ii )
    train_out_pairs_mat[ ix , ix ] = 1 
    ix = which( test_res$cluster == ii )
    test_out_pairs_mat[ ix , ix ] = 1
    }
  ## coerce to upper triangular matrix
  train_out_pairs_mat = train_out_pairs_mat * upper_tri_mask
  test_out_pairs_mat = test_out_pairs_mat * upper_tri_mask
  
  ## get minimum co-membership cluster proportion
  sample_min_vec = matrix( 0L , nrow = 1 , ncol = length( test_res$centers ) )
  for ( ii in 1 : length( test_res$centers ) ) {
    ix = which( test_res$cluster == ii )
    test_cluster_sum = sum( test_out_pairs_mat[ ix , ix ] )
    train_cluster_sum = sum( test_out_pairs_mat[ ix , ix ] * train_out_pairs_mat[ ix , ix ] )
    sample_min_vec[ , ii ] = train_cluster_sum / test_cluster_sum
  }
  
  ## get min of sample_min_vec
  min_val = min( sample_min_vec[ !is.nan( sample_min_vec ) ] ) ## removing any NaN
  all_k_ps[ , k ] = all_k_ps[ , k ] + min_val

} ## end of K for loop

} ## end of sample loop

all_k_ps = all_k_ps / no_sample_iters ## average values
plot( 1 : length( all_k_ps ) , all_k_ps , "b" , xlab = "Number of Clusters K" , ylab = "Prediction Strength Value" )
abline( h = 0.8 , col = "red" )

The purpose of the cross validation routine is to select the number of clusters K, in the model selection sense, that is best supported by the available data. The above linked paper suggests that the optimal number of clusters K is the highest number K that has a prediction strength value over some given threshold (e.g. 0.8 or 0.9). The last part of the code plots the values of prediction strength for K (x-axis) vs. prediction strength (y-axis), along with the threshold value of 0.8 in red. For the particular set of data in question, it can be seen that the optimal K value for the number of clusters is 8.

This second code box shows code, re-using some of the above code, to visualise the clusters for a given K,

library( Ckmeans.1d.dp )

## load the training data from Octave output (comment out as necessary )
data = read.csv( "~/path/to/all_data_matrix" , header = FALSE )
data_sum = colSums( data ) ## sum down the columns of data
data_sum[ 1 : 5 ] = mean( data_sum[ 1 : 48 ] ) ## correct artifacts of weekend gap
data_sum[ 716 : 720 ] = mean( data_sum[ 1 : 48 ] ) ## correct artifacts of weekend gap

## comment out as necessary
adjust = 0 ## default adjust value
sum_seq = seq( from = 1 , to = 198 , by = 1 ) ; sum_seq_l = as.numeric( length( sum_seq ) ) ## Monday
##sum_seq = seq( from = 115 , to = 342 , by = 1 ) ; sum_seq_l = as.numeric( length( sum_seq ) ) ## Tuesday
# sum_seq = seq( from = 115 , to = 342 , by = 1 ) ; sum_seq_l = as.numeric( length( sum_seq ) ) ## Wednesday
# sum_seq = seq( from = 115 , to = 342 , by = 1 ) ; sum_seq_l = as.numeric( length( sum_seq ) ) ## Thursday
##sum_seq = seq( from = 547 , to = 720 , by = 1 ) ; sum_seq_l = as.numeric( length( sum_seq ) ) ## Friday

## intraday --- commnet out or adjust as necessary
##sum_seq = seq( from = 25 , to = 100 , by = 1 ) ; sum_seq_l = as.numeric( length( sum_seq ) )

k = 8
res = Ckmeans.1d.dp( sum_seq , k , data_sum[ sum_seq ] )

plot( sum_seq , data_sum[ sum_seq ], main = "Cluster centres. Cluster centre ix is a predicted turning point",
     col = res$cluster,
     pch = res$cluster, type = "h", xlab = "Count from beginning ix at ix = 1",
     ylab = "Total Counts per ix" )

abline( v = res$centers, col = "chocolate" , lty = "dashed" )

text( res$centers, max(data_sum[sum_seq]) * 0.95, cex = 0.75, font = 2,
      paste( round(res$centers) ) )

a typical plot for which is shown below.

The above plot can be thought of as a clustering at a particular scale, and one can go down in scale by selecting smaller ranges of the data. For example, taking all the datum clustered in the 3 clusters centred at x-axis ix values 38, 63 and 89 and re-running the code in the first code box on just this data gives this prediction strength plot, which suggests a K value of 6.

Re-running the code in the second code box plots these 6 clusters thus.

Looking at this last plot, it can be seen that there is a cluster at x-axis ix value 58, which corresponds to 7.30 a.m. London time, and within this green cluster there are 2 distinct peaks which correspond to 7.00 a.m. and 8.00 a.m. A similar, visual analysis of the far right cluster, centre ix = 94, shows a peak at the time of the New York open.

My hypothesis is that by clustering in the above manner it will be possible to identify distinct, intraday times at which the probability of a market turn is greater than at other times. More in due course.

A Temporal Clustering Function, Part 2

Further to my previous post, below is an extended version of the "blurred_maxshift_1d_linear" function. This updated version has two extra outputs: a vector of the cluster centre index ix values and a vector the same length as the input data with the cluster centres to which each datum has been assigned. These changes have necessitated some extensive re-writing of the function to include various checks contained in nested conditional statements.

## Copyright (C) 2020 dekalog
## 
## This program is free software: you can redistribute it and/or modify it
## under the terms of the GNU General Public License as published by
## the Free Software Foundation, either version 3 of the License, or
## (at your option) any later version.
## 
## This program is distributed in the hope that it will be useful, but
## WITHOUT ANY WARRANTY; without even the implied warranty of
## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
## GNU General Public License for more details.
## 
## You should have received a copy of the GNU General Public License
## along with this program.  If not, see
## .

## -*- texinfo -*- 
## @deftypefn {} {@var{train_vec}, @var{cluster_centre_ix}, @var{assigned_cluster_centre_ix} =} blurred_maxshift_1d_linear_V2 (@var{train_vec}, @var{bandwidth})
##
## @seealso{}
## @end deftypefn

## Author: dekalog 
## Created: 2020-10-21

function [ new_train_vec , cluster_centre_ix , assigned_cluster_centre_ix ] = blurred_maxshift_1d_linear_V2 ( train_vec , bandwidth )

if ( nargin < 2 )
 bandwidth = 1 ;
endif

if ( numel( train_vec ) < 2 * bandwidth + 1 )
 error( 'Bandwidth too wide for length of train_vec.' ) ;
endif

length_train_vec = numel( train_vec ) ;
new_train_vec = zeros( size( train_vec ) ) ;
assigned_cluster_centre_ix = ( 1 : 1 : length_train_vec ) ;

## initialising loop
## do the beginning 
[ ~ , ix ] = max( train_vec( 1 : 2 * bandwidth + 1 ) ) ;
new_train_vec( ix ) = sum( train_vec( 1 : bandwidth + 1 ) ) ;
assigned_cluster_centre_ix( 1 : bandwidth + 1 ) = ix ;

## and end of train_vec first
[ ~ , ix ] = max( train_vec( end - 2 * bandwidth : end ) ) ;
new_train_vec( end - 2 * bandwidth - 1 + ix ) = sum( train_vec( end - bandwidth : end ) ) ;
assigned_cluster_centre_ix( end - bandwidth : end ) = length_train_vec - 2 * bandwidth - 1 + ix ;

for ii = ( bandwidth + 2 ) : ( length_train_vec - bandwidth - 1 )
 [ ~ , ix ] = max( train_vec( ii - bandwidth : ii + bandwidth ) ) ;
 new_train_vec( ii - bandwidth - 1 + ix ) += train_vec( ii ) ;
 assigned_cluster_centre_ix( ii ) = ii - bandwidth - 1 + ix ; 
endfor
## end of initialising loop

train_vec = new_train_vec ;

## initialise the while condition variable
has_converged = 0 ;

while ( has_converged < 1 )

new_train_vec = zeros( size( train_vec ) ) ;

## do the beginning 
[ ~ , ix ] = max( train_vec( 1 : 2 * bandwidth + 1 ) ) ;
new_train_vec( ix ) += sum( train_vec( 1 : bandwidth + 1 ) ) ;
assigned_cluster_centre_ix( 1 : bandwidth + 1 ) = ix ;

## and end of train_vec first
[ ~ , ix ] = max( train_vec( end - 2 * bandwidth : end ) ) ;
new_train_vec( end - 2 * bandwidth - 1 + ix ) += sum( train_vec( end - bandwidth : end ) ) ;
assigned_cluster_centre_ix( end - bandwidth : end ) = length_train_vec - 2 * bandwidth - 1 + ix ;

for ii = ( bandwidth + 2 ) : ( length_train_vec - bandwidth - 1 )

 [ max_val , ix ] = max( train_vec( ii - bandwidth : ii + bandwidth ) ) ;
 ## check for ties in max_val value in window
 no_ties = sum( train_vec( ii - bandwidth : ii + bandwidth ) == max_val ) ;

  if ( no_ties == 1 && max_val == train_vec( ii ) && ix == bandwidth + 1 ) ## main if
   ## value in train_vec(ii) is max val of window, with no ties
   new_train_vec( ii ) += train_vec( ii ) ;
   assigned_cluster_centre_ix( ii ) = ii ;

  elseif ( no_ties == 1 && max_val != train_vec( ii ) && ix != bandwidth + 1 ) ## main if
   ## no ties for max_val, but need to move data at ii and change ix
   ## get assigned_cluster_centre_ix that point to ii, which needs to be updated
   assigned_ix = find( assigned_cluster_centre_ix == ii ) ;
    if ( !isempty( assigned_ix ) ) ## should always be true because at least the one original ii == ii
     assigned_cluster_centre_ix( assigned_ix ) = ii - ( bandwidth + 1 ) + ix ;
    elseif ( isempty( assigned_ix ) ) ## but cheap insurance
     assigned_cluster_centre_ix( ii ) = ii - ( bandwidth + 1 ) + ix ;
    endif
   new_train_vec( ii - ( bandwidth + 1 ) + ix ) += train_vec( ii ) ;

  elseif ( no_ties > 1 && max_val > train_vec( ii ) ) ## main if
   ## 2 ties for max_val, which is > val at ii, need to move data at ii 
   ## to the closer max_val ix and change ix in assigned_cluster_centre_ix
   match_max_val_ix = find( train_vec( ii - bandwidth : ii + bandwidth ) == max_val ) ;

    if ( numel( match_max_val_ix ) == 2 ) ## only 2 matching max vals
     centre_window_dist = ( bandwidth + 1 ) .- match_max_val_ix ;

           if ( abs( centre_window_dist( 1 ) ) == abs( centre_window_dist( 2 ) ) ) 
            ## equally distant from centre ii of moving window

                  assigned_ix = find( assigned_cluster_centre_ix == ii ) ;
                   if ( !isempty( assigned_ix ) ) ## should always be true because at least the one original ii == ii
                    ix_before = find( assigned_ix < ii ) ;
                    ix_after = find( assigned_ix > ii ) ;
                    new_train_vec( ii - ( bandwidth + 1 ) + match_max_val_ix( 1 ) ) += train_vec( ii ) / 2 ;
                    new_train_vec( ii - ( bandwidth + 1 ) + match_max_val_ix( 2 ) ) += train_vec( ii ) / 2 ;
                    assigned_cluster_centre_ix( assigned_ix( ix_before ) ) = ii - ( bandwidth + 1 ) + match_max_val_ix( 1 ) ;
                    assigned_cluster_centre_ix( assigned_ix( ix_after ) ) = ii - ( bandwidth + 1 ) + match_max_val_ix( 2 ) ;
                    assigned_cluster_centre_ix( ii ) = ii ; ## bit of a kluge
                   elseif ( isempty( assigned_ix ) ) ## but cheap insurance
                    ## no other assigned_cluster_centre_ix values to account for, so just split equally
                    new_train_vec( ii - ( bandwidth + 1 ) + match_max_val_ix( 1 ) ) += train_vec( ii ) / 2 ;
                    new_train_vec( ii - ( bandwidth + 1 ) + match_max_val_ix( 2 ) ) += train_vec( ii ) / 2 ;
                    assigned_cluster_centre_ix( ii ) = ii ; ## bit of a kluge
                   else
                    error( 'There is an unknown error in instance ==2 matching max_vals with equal distances to centre of moving window with assigned_ix. Write code to deal with this edge case.' ) ;
                   endif

           else ## not equally distant from centre ii of moving window

                  assigned_ix = find( assigned_cluster_centre_ix == ii ) ;
                   if ( !isempty( assigned_ix ) ) ## should always be true because at least the one original ii == ii
                    ## There is an instance == 2 matching max_vals with non equal distances to centre of moving window with previously assigned_ix to ii ix
                    ## Assign all assigned_ix to the nearest max value ix
                    [ ~ , min_val_ix ] = min( [ abs( centre_window_dist( 1 ) ) abs( centre_window_dist( 2 ) ) ] ) ;
                    new_train_vec( ii - ( bandwidth + 1 ) + match_max_val_ix( min_val_ix ) ) += train_vec( ii ) ;
                    assigned_cluster_centre_ix( ii ) = ii - ( bandwidth + 1 ) + match_max_val_ix( min_val_ix ) ;
                    assigned_cluster_centre_ix( assigned_ix ) = ii - ( bandwidth + 1 ) + match_max_val_ix( min_val_ix ) ;
                   elseif ( isempty( assigned_ix ) ) ## but cheap insurance
                    [ ~ , min_val_ix ] = min( abs( centre_window_dist ) ) ;
                    new_train_vec( ii - ( bandwidth + 1 ) + match_max_val_ix( min_val_ix ) ) += train_vec( ii ) ;
                    assigned_cluster_centre_ix( ii ) = ii - ( bandwidth + 1 ) + match_max_val_ix( min_val_ix ) ;
                   else
                    error( 'There is an unknown error in instance of ==2 matching max_vals with unequal distances. Write the code to deal with this edge case.' ) ;
                   endif

           endif ## 

    elseif ( numel( match_max_val_ix ) > 2  ) ## There is an instance of >2 matching max_vals.
    ## There must be one max val closer than the others or two equally close
     centre_window_dist = abs( ( bandwidth + 1 ) .- match_max_val_ix ) ;
     centre_window_dist_min = min( centre_window_dist ) ;
     centre_window_dist_min_ix = find( centre_window_dist == centre_window_dist_min ) ;
     
       if ( numel( centre_window_dist_min_ix ) == 1 ) ## there is one closet ix
        assigned_ix = find( assigned_cluster_centre_ix == ii ) ;
        
            if ( !isempty( assigned_ix ) ) ## should always be true because at least the one original ii == ii
               new_train_vec( ii - ( bandwidth + 1 ) + centre_window_dist_min_ix ) += train_vec( ii ) ;  
               assigned_cluster_centre_ix( ii ) = ii - ( bandwidth + 1 ) + centre_window_dist_min_ix ;
               assigned_cluster_centre_ix( assigned_ix ) = ii - ( bandwidth + 1 ) + centre_window_dist_min_ix ;
            elseif ( isempty( assigned_ix ) ) ## but cheap insurance
               new_train_vec( ii - ( bandwidth + 1 ) + centre_window_dist_min_ix ) += train_vec( ii ) ;  
               assigned_cluster_centre_ix( ii ) = ii - ( bandwidth + 1 ) + centre_window_dist_min_ix ;
            endif
         
       elseif ( numel( centre_window_dist_min_ix ) == 2 ) ## there are 2 equally close ix
        assigned_ix = find( assigned_cluster_centre_ix == ii ) ;
        
            if ( !isempty( assigned_ix ) ) ## should always be true because at least the one original ii == ii
               ix_before = find( assigned_ix < ii ) ;
               ix_after = find( assigned_ix > ii ) ;
               new_train_vec( ii - ( bandwidth + 1 ) + centre_window_dist_min_ix( 1 ) ) += train_vec( ii ) / 2 ;
               new_train_vec( ii - ( bandwidth + 1 ) + centre_window_dist_min_ix( 2 ) ) += train_vec( ii ) / 2 ;
               assigned_cluster_centre_ix( assigned_ix( ix_before ) ) = ii - ( bandwidth + 1 ) + centre_window_dist_min_ix( 1 ) ;
               assigned_cluster_centre_ix( assigned_ix( ix_after ) ) = ii - ( bandwidth + 1 ) + centre_window_dist_min_ix( 2 ) ;
               assigned_cluster_centre_ix( ii ) = ii ; ## bit of a kluge             
            elseif ( isempty( assigned_ix ) ) ## but cheap insurance 
               ## no other assigned_cluster_centre_ix values to account for, so just split equally
               new_train_vec( ii - ( bandwidth + 1 ) + centre_window_dist_min_ix( 1 ) ) += train_vec( ii ) / 2 ;
               new_train_vec( ii - ( bandwidth + 1 ) + centre_window_dist_min_ix( 2 ) ) += train_vec( ii ) / 2 ;
               assigned_cluster_centre_ix( ii ) = ii ; ## bit of a kluge
            endif
        
       else
        error( 'Unknown error in numel( match_max_val_ix ) > 2.' ) ;
       endif
     ##error( 'There is an instance of >2 matching max_vals. Write the code to deal with this edge case.' ) ;
    else
     error( 'There is an unknown error in instance of >2 matching max_vals. Write the code to deal with this edge case.' ) ;
    endif

  endif ## main if end

endfor

if ( sum( ( train_vec == new_train_vec ) ) == length_train_vec )
 has_converged = 1 ;
else
 train_vec = new_train_vec ;
endif

endwhile

cluster_centre_ix = unique( assigned_cluster_centre_ix ) ;
cluster_centre_ix( cluster_centre_ix == 0 ) = [] ;

endfunction

The reason for this re-write was to accommodate a cross validation routine, which is described in the paper Cluster Validation by Prediction Strength, and a simple outline of which is given in this stackexchange.com answer.

My Octave code implementation of this is shown in the code box below. This is not exactly as described in the above paper because the number of clusters, K, is not exactly specified due to the above function automatically determining K based on the data. The routine below is perhaps more accurately described as being inspired by the original paper.

## create train and test data sets
########## UNCOMMENT AS NECESSARY #####
time_ix = [ 1 : 198 ] ; ## Monday
##time_ix = [ 115 : 342 ] ; ## Tuesday
##time_ix = [ 259 : 486 ] ; ## Wednesday
##time_ix = [ 403 : 630 ] ; ## Thursday
##time_ix = [ 547 : 720 ] ; ## Friday
##time_ix = [ 1 : 720 ] ; ## all data
#######################################

all_cv_solutions = zeros( size( data_matrix , 3 ) , size( data_matrix , 3 ) ) ;

n_iters = 1 ;
for iter = 1 : n_iters ## related to # of rand_ix sets generated 
 
rand_ix = randperm( size( data_matrix , 1 ) ) ; 
train_ix = rand_ix( 1 : round( numel( rand_ix ) * 0.5 ) ) ;
test_ix = rand_ix( round( numel( rand_ix ) * 0.5 ) + 1 : end ) ;
train_data_matrix = sum( data_matrix( train_ix , time_ix , : ) ) ;
test_data_matrix = sum( data_matrix( test_ix , time_ix , : ) ) ;

all_proportions_indicated = zeros( 1 , size( data_matrix , 3 ) ) ;

for cv_ix = 1 : size( data_matrix , 3 ) ; ## related to delta_turning_point_filter n_bar parameter

for bandwidth = 1 : size( data_matrix , 3 )

## train set clustering
if ( bandwidth == 1 )
[ train_out , cluster_centre_ix_train , assigned_cluster_centre_ix_train ] = blurred_maxshift_1d_linear_V2( train_data_matrix(:,:,cv_ix) , bandwidth ) ;
elseif( bandwidth > 1 )
[ train_out , cluster_centre_ix_train , assigned_cluster_centre_ix_train ] = blurred_maxshift_1d_linear_V2( train_out , bandwidth ) ;
endif
train_out_pairs_mat = zeros( numel( assigned_cluster_centre_ix_train ) , numel( assigned_cluster_centre_ix_train ) ) ;
 for ii = 1 : numel( cluster_centre_ix_train )
  cc_ix = find( assigned_cluster_centre_ix_train == cluster_centre_ix_train( ii ) ) ;
  train_out_pairs_mat( cc_ix , cc_ix ) = 1 ;
 endfor
train_out_pairs_mat = triu( train_out_pairs_mat , 1 ) ; ## get upper diagonal matrix

## test set clustering
if ( bandwidth == 1 )
[ test_out , cluster_centre_ix_test , assigned_cluster_centre_ix_test ] = blurred_maxshift_1d_linear_V2( test_data_matrix(:,:,cv_ix) , bandwidth ) ;
elseif( bandwidth > 1 )
[ test_out , cluster_centre_ix_test , assigned_cluster_centre_ix_test ] = blurred_maxshift_1d_linear_V2( test_out , bandwidth ) ;
endif

all_test_out_pairs_mat = zeros( numel( assigned_cluster_centre_ix_test ) , numel( assigned_cluster_centre_ix_test ) ) ;
test_out_pairs_clusters_proportions = ones( 1 , numel( cluster_centre_ix_test ) ) ;
 for ii = 1 : numel( cluster_centre_ix_test )
  cc_ix = find( assigned_cluster_centre_ix_test == cluster_centre_ix_test( ii ) ) ;
  all_test_out_pairs_mat( cc_ix , cc_ix ) = 1 ;
  test_out_pairs_mat = all_test_out_pairs_mat( cc_ix , cc_ix ) ;
  test_out_pairs_mat = triu( test_out_pairs_mat , 1 ) ; ## get upper diagonal matrix
  test_out_pairs_mat_sum = sum( sum( test_out_pairs_mat ) ) ;
   if ( test_out_pairs_mat_sum > 0 )
   test_out_pairs_clusters_proportions( ii ) = sum( sum( train_out_pairs_mat( cc_ix , cc_ix ) .* test_out_pairs_mat ) ) / ...
                                                test_out_pairs_mat_sum ;
   endif
 endfor

all_proportions_indicated( bandwidth ) = min( test_out_pairs_clusters_proportions ) ;
all_cv_solutions( bandwidth , cv_ix ) += all_proportions_indicated( bandwidth ) ; 

endfor ## bandwidth for loop

endfor ## end of cv_ix loop

endfor ## end of iter for

all_cv_solutions = all_cv_solutions ./ n_iters ;

surf( all_cv_solutions ) ; xlabel( 'BANDWIDTH' , 'fontsize' , 20 ) ; ylabel( 'CV IX' , 'fontsize' , 20 ) ;

I won't discuss the workings of the code any further as readers are free to read the original paper and my code interpretation of it. A typical surface plot of the output is shown below.

The "bandwidth" plotted along the front edge of the surface plot is one of the input parameters to the "blurred_maxshift_1d_linear" function, whilst the "lookback" is a parameter of the original data generating function which identifies local highs and lows in a price time series. There appears to be a distinct "elbow" at "lookback" = 6 which is more or less consistent for all values of "bandwidth." Since the underlying data for this is 10 minute OHLC bars, the ideal "lookback" would, therefore, appear to be on the hourly timeframe.

However, having spent some considerable time and effort to get the above working satisfactorily, I am now not so sure that I'll actually use the above code. The reason for this is shown in the following animated GIF.

This shows K segmentation of the exact same data used above, from K = 1 to 19 inclusive, using R and its Ckmeans.1d.dp package, with vignette tutorial here. I am particularly attracted to this because of its speed, compared to my code above, as well as its guarantees with regard to optimality and reproducability. If one stares at the GIF for long enough one can see possible, significant clusters at index values (x-axis) which correspond, approximately, to particularly significant times such as London and New York market opening and closing times: ix = 55, 85, 115 and 145.

More about this in my next post. 

A Temporal Clustering Function

Recently a reader contacted me with a view to collaborating on some work regarding the Delta phenomenon but after a brief exchange of e-mails this seems to have petered out. However, for my part, the work I have done has opened a few new avenues of investigation and this post is about one of them.

One of the problems I set out to solve was clustering in the time domain, or temporal clustering as I call it. Take a time series and record the time of occurance of an event by setting to 1, in an otherwise zero filled 1-dimensional vector the same length as the original time series, the value of the vector at time index tx and repeat for all occurances of the event. In my case the event(s) I am interested in are local highs and lows in the time series. This vector is then "chopped" into segments representing distinct periods of time, e.g. 1 day, 1 week etc. and stacked into a matrix where each row is one complete period and the columns represent the same time in each period, e.g. the first column is the first hour of the trading week, the second column the second hour etc. Sum down the columns to get a final 1-dimensional vector of counts of the timing of events happening within each period over the entire time series data record. A chart of such is shown below.

The coloured vertical lines show the opening and closing times of the London and New York sessions (7am to 5pm in their respective local times) for one complete week at a 10 minute bar time scale, in this case for the GBP_USD forex pair. This is what I want to cluster.

The solution I have come up with is the Octave function in the code box below

## Copyright (C) 2020 dekalog
## 
## This program is free software: you can redistribute it and/or modify it
## under the terms of the GNU General Public License as published by
## the Free Software Foundation, either version 3 of the License, or
## (at your option) any later version.
## 
## This program is distributed in the hope that it will be useful, but
## WITHOUT ANY WARRANTY; without even the implied warranty of
## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
## GNU General Public License for more details.
## 
## You should have received a copy of the GNU General Public License
## along with this program.  If not, see
## .

## -*- texinfo -*- 
## @deftypefn {} {@var{centre_ix} =} blurred_maxshift_1d_linear (@var{train_vec}, @var{bandwidth})
##
## Clusters the 1 dimensional vector TRAIN_VEC using a "centred" sliding window of length  2 * BANDWIDTH + 1.
##
## Based on the idea of the Blurred Meanshift Algorithm.
##
## The "centre ix" value of the odd length sliding window is assigned to the
## maximum value ix of the sliding window. The centre_ix, if it is not the 
## maximum value, is then set to zero. A pass through the whole length of
## TRAIN_VEC is completed before any assignments are made.
##
## @seealso{}
## @end deftypefn

## Author: dekalog 
## Created: 2020-10-17

function new_train_vec = blurred_maxshift_1d_linear ( train_vec , bandwidth )

if ( nargin < 2 )
 bandwidth = 1 ;
endif

if ( numel( train_vec ) < 2 * bandwidth + 1 )
 error( 'Bandwidth too wide for length of train_vec.' ) ;
endif

length_train_vec = numel( train_vec ) ;
assigned_cluster_centre_ix = zeros( size( train_vec ) ) ;

## initialise the while condition variable
has_converged = 0 ;

while ( has_converged < 1 )
 
new_train_vec = zeros( size( train_vec ) ) ;

## do the beginning and end of train_vec first
[ ~ , ix ] = max( train_vec( 1 : 2 * bandwidth + 1 ) ) ;
new_train_vec( ix ) = sum( train_vec( 1 : bandwidth ) ) ;

[ ~ , ix ] = max( train_vec( end - 2 * bandwidth : end ) ) ;
new_train_vec( end - 2 * bandwidth - 1 + ix ) = sum( train_vec( end - bandwidth + 1 : end ) ) ;

for ii = 2 * bandwidth + 1 : numel( train_vec ) - bandwidth
 [ ~ , ix ] = max( train_vec( ii - bandwidth : ii + bandwidth ) ) ;
 new_train_vec( ii - bandwidth  - 1 + ix ) += train_vec( ii ) ; 
endfor

if ( sum( ( train_vec == new_train_vec ) ) == length_train_vec )
 has_converged = 1 ;
else
 train_vec = new_train_vec ;
endif

endwhile

endfunction

I have named the function "blurred_maxshift_1d_linear" as it is inspired by the "blurred" version of the Mean shift algorithm, operates on a 1-dimensional vector and is "linear" in that there is no periodic wrapping of the data within the function code. The two function inputs are the above type of data, obviously, and an integer parameter "bandwidth" which controls the size of a moving window over the data in which the data is shifted according to a maximum value, hence maxshift rather than meanshift. I won't discuss the code further as it is pretty straightforward.

A chart of a typical clustering solution is (bandwidth setting == 2)

where the black line is the original count data and red the clustering solution. The bandwidth setting in this case is approximately equivalent to clustering with a 50 minute moving window. 

The following heatmap chart is a stacked version of the above where the bandwidth parameter is varied from 1 to 10 inclusive upwards, with the original data being at the lowest level per pane.

The intensity reflects the counts at each time tx index per bandwidth setting. The difference between the panes is that in the upper pane the raw data is the function input per bandwidth setting, whilst the lower pane shows hierarchical clustering whereby the output of the function is used as the input to the next function call with the next higher bandwidth parameter setting.

More in due course.

Candlestick Pattern Scanner Functions

Since my last currency strength candlestick chart post it seemed to make sense to be able to scan said charts for signals, so below is the code for two Octave functions which act as candlestick pattern scanners. The code is fully vectorised and self-contained, and on my machine they can scan more than 300,000 OHLC bars for 27/29 separate patterns in less than 0.5 seconds. Both functions have a self-explanatory help and within the body of the code there are ample comments which describe the patterns. Enjoy!

Bullish Reversal Indicator

## Copyright (C) 2020 dekalog
## 
## This program is free software: you can redistribute it and/or modify it
## under the terms of the GNU General Public License as published by
## the Free Software Foundation, either version 3 of the License, or
## (at your option) any later version.
## 
## This program is distributed in the hope that it will be useful, but
## WITHOUT ANY WARRANTY; without even the implied warranty of
## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
## GNU General Public License for more details.
## 
## You should have received a copy of the GNU General Public License
## along with this program.  If not, see
## .

## -*- texinfo -*- 
## @deftypefn {} {@var{bullish_signal_matrix} =} candle_bullish_reversal (@var{high}, @var{low}, @var{close}, @var{open}, @var{downtrend}, @var{avge_hi_lo_range})
##
## The HIGH, LOW, CLOSE and OPEN should be vectors of equal length, and are required.
##
## Internally the DOWNTREND is determined as a bar's CLOSE being lower than the CLOSE
## of the bar 5 bars previously. This can be over-ruled by an optional, user supplied
## vector DOWNTREND of zeros and ones, where one(s) represent a bar in a DOWNTREND.
## If a DOWNTREND vector consists solely of ones, this effectively turns off the
## discriminative ability of the DOWNTREND condition as all bars will be considered
## to be DOWNTREND bars. If the DOWNTREND vector consists solely of zeros, no bar
## will be classified as being in a DOWNTREND and those candle patterns that require a 
## DOWNTREND as a condition will not be indicated.
##
## The AVGE_HI_LO_RANGE is calculated as a rolling 5 bar simple moving average of
## the HI-LO range of bars. A bar is considered a "long line" if its HI-LO range
## is greater than the AVGE_HI_LO_RANGE. This can also be over-ruled, as above, by
## an optional, user supplied vector AVGE_HI_LO_RANGE, which consists of zeros and
## ones, with ones representing "long line" bars. An AVGE_HI_LO_RANGE vector of all
## ones or zeros will have a similar effect as described for the DOWNTREND vector.
##
## The candlestick pattern descriptions are predominantly taken from
##
## https://www.candlescanner.com
##
## and the bullish only patterns implemented are:
##
## Bullish Reversal High Reliability
##
##     01 - Bullish Abandoned Baby
##
##     02 - Concealing Baby Swallow
##
##     03 - Kicking
##
##     04 - Morning Doji Star
##
##     05 - Morning Star
##
##     06 - Piercing Line
##
##     07 - Three Inside Up
##
##     08 - Three Outside Up
##
##     09 - Three White Soldiers
##
## Bullish Reversal Moderate Reliability
##
##     10 - Breakaway
##
##     11 - Counter Attack
##
##     12 - Doji Star
##
##     13 - Dragonfly Doji
##
##     14 - Engulfing
##
##     15 - Gravestone Doji
##
##     16 - Harami Cross
##
##     17 - Homing Pigeon
##
##     18 - Ladder Bottom
##
##     19 - Long Legged Doji
##
##     20 - Matching Low
##
##     21 - Meeting Lines
##
##     22 - Stick Sandwich
##
##     23 - Three Stars in the South
##
##     24 - Tri Star
##
##     25 - Unique Three River Bottom
##
## Bullish Reversal Low Reliability
##
##     26 - Belt Hold
##
##     27 - Hammer
##
##     28 - Harami
##
##     29 - Inverted Hammer
##
## The output BULLISH_SIGNAL_MATRIX has a row length the same as the OHLC
## input vectors and 29 columns. The matrix is a zero filled matrix with the
## value 1 when a candle pattern is indicated. The row ix of the value is the
## row ix of the last candle in the pattern. The column ix corresponds to the
## pattern index given above.
##
## @seealso{candle_bearish_reversal}
## @end deftypefn

## Author: dekalog 
## Created: 2020-08-10

function bullish_signal_matrix = candle_bullish_reversal ( varargin )
 
if ( nargin < 4 || nargin > 6 )
  print_usage () ;
elseif ( nargin == 4 )
  high = varargin{1} ; low = varargin{2} ; close = varargin{3} ; open = varargin{4} ;
  downtrend = close < shift( close , 5 ) ;
  avge_hi_lo_range = filter( [0.2 ; 0.2 ; 0.2 ; 0.2 ; 0.2 ] , 1 , high .- low ) ; ## 5 bar simple moving average
elseif ( nargin == 5 )
  high = varargin{1} ; low = varargin{2} ; close = varargin{3} ; open = varargin{4} ;
  downtrend = varargin{5} ;
  avge_hi_lo_range = filter( [0.2 ; 0.2 ; 0.2 ; 0.2 ; 0.2 ] , 1 , high .- low ) ; ## 5 bar simple moving average
elseif ( nargin == 6 )
  high = varargin{1} ; low = varargin{2} ; close = varargin{3} ; open = varargin{4} ;
  downtrend = varargin{5} ; avge_hi_lo_range = varargin{6} ;
endif
 
bullish_signal_matrix = zeros( size( close , 1 ) , 29 ) ;

## pre-compute some basic vectors for re-use below, trading memory use for
## speed of execution
downtrend_1 = shift( downtrend , 1 ) ;
downtrend_2 = shift( downtrend , 2 ) ;
downtrend_3 = shift( downtrend , 3 ) ;
long_line = ( high .- low ) > shift( avge_hi_lo_range , 1 ) ;
long_line_1 = shift( long_line , 1 ) ;
long_line_2 = shift( long_line , 2 ) ;
long_line_3 = shift( long_line , 3 ) ;
long_line_4 = shift( long_line , 4 ) ;
open_1 = shift( open , 1 ) ;
open_2 = shift( open , 2 ) ;
open_3 = shift( open , 3 ) ;
open_4 = shift( open , 4 ) ;
open_5 = shift( open , 5 ) ;
high_1 = shift( high , 1 ) ;
high_3 = shift( high , 3 ) ;
low_1 = shift( low ,1 ) ;
low_2 = shift( low ,2 ) ;
low_3 = shift( low , 3 ) ;
close_1 = shift( close , 1 ) ;
close_2 = shift( close , 2 ) ;
close_3 = shift( close , 3 ) ;
close_4 = shift( close , 4 ) ;
black_body = close < open ;
black_body_1 = shift( black_body , 1 ) ;
black_body_2 = shift( black_body , 2 ) ;
black_body_3 = shift( black_body , 3 ) ;
black_body_4 = shift( black_body , 4 ) ;
white_body = close > open ;
white_body_1 = shift( white_body , 1 ) ;
white_body_2 = shift( white_body , 2 ) ;
doji = ( close == open ) ;
doji_1 = shift( doji , 1 ) ;
doji_2 = shift( doji , 2 ) ;
body_high = max( [ open , close ] , [] , 2 ) ;
body_high_1 = shift( body_high , 1 ) ;
body_high_2 = shift( body_high , 2 ) ;
body_low = min( [ open , close ] , [] , 2 ) ;
body_low_1 = shift( body_low , 1 ) ;
body_low_2 = shift( body_low , 2 ) ;
body_range = body_high .- body_low ;
body_range_1 = shift( body_range , 1 ) ;
body_range_2 = shift( body_range , 2 ) ;
body_midpoint = ( body_high .+ body_low ) ./ 2 ;
body_midpoint_1 = shift( body_midpoint , 1 ) ;
body_midpoint_2 = shift( body_midpoint , 2 ) ;
wick_gap_up = ( low > high_1 ) ;
wick_gap_down = ( high < low_1 ) ;
wick_gap_down_1 = shift( wick_gap_down , 1 ) ;
black_marubozu = ( high == open ) .* ( low == close ) ;
black_marubozu_1 = shift( black_marubozu , 1 ) ;
black_marubozu_2 = shift( black_marubozu , 2 ) ;
black_marubozu_3 = shift( black_marubozu , 3 ) ;
white_marubozu = ( high == close ) .* ( low == open ) ;

## 01 - Bullish Abandoned Baby
## First candle
##     a candle in a downtrend
##     black body
## Second candle
##     a doji candle
##     the high price below the prior low price
## Third candle
##     white body
##     the low price above the prior high price
bullish_signal_matrix(:,1) = downtrend_2 .* black_body_2 .* ...
                              doji_1 .* wick_gap_down_1 .* ...
                              white_body .* wick_gap_up ; 
## 02 - Concealing Baby Swallow
## First candle
##     a Black Marubozu candle in a downtrend
## Second candle
##     a Black Marubozu candle
##     candle opens within the prior candle's body
##     candle closes below the prior closing price
## Third candle
##     a High Wave basic candle with no lower shadow
##     candle opens below the prior closing price
##     upper shadow enters the prior candle's body
## Fourth candle
##     black body
##     candle’s body engulfs the prior candle’s body including the shadows
bullish_signal_matrix(:,2) = downtrend_3 .* black_marubozu_3 .* ...
                              black_marubozu_2 .* (open_2 < high_3) .* (open_2 > low_3) .* (close_2 < close_3) .* ...
                              (high_1 .- body_high_1 >= 3.*body_range_1) .* (low_1 == body_low_1) .* (open_1 < close_2) .* (high_1 > body_low_2) .* ...
                              black_body .* (open > high_1) .* (close < low_1) ;
## 03 - Kicking
## First candle
##     a Black Marubozu
##     appears on as a long line
## Second candle
##     a White Marubozu
##     price gaps upward
##     appears on as a long line
bullish_signal_matrix(:,3) = black_marubozu_1 .* long_line_1 .* ...
                              white_marubozu .* (low > high_1) .* long_line ;
## 04 - Morning Doji Star
## First candle
##     a candle in a downtrend
##     black body
## Second candle
##     a doji candle
##     a doji body below the previous candle body
##     the high price above the previous candle low price
## Third candle
##     white body
##     candle body above the previous candle body
##     the closing price above the midpoint of the first candle body
bullish_signal_matrix(:,4) = downtrend_2 .* black_body_2 .* ...
                              doji_1 .* (open_1 < body_low_2) .* (high_1 > low_2) .* ...
                              white_body .* (body_low > body_high_1) .* (close > body_midpoint_2) ;
## 05 - Morning Star
## First candle
##     a candle in a downtrend
##     black body
## Second candle
##     white or black body
##     the candle body is located below the prior body
## Third candle
##     white body
##     the candle body is located above the prior body
##     the candle closes at least halfway up the body of the first line
bullish_signal_matrix(:,5) = downtrend_2 .* black_body_2 .* ...
                              (body_high_1 < body_low_2) .* ...
                              white_body .* (body_low > body_high_1) .* (close >= body_midpoint_2) ;
## 06 - Piercing Line
## First candle
##     a candle in a downtrend
##     black body
## Second candle
##     white body
##     the opening below or equal of the prior low
##     the closing above the midpoint of the prior candle's body
##     the closing below the previous opening
bullish_signal_matrix(:,6) = downtrend_1 .* black_body_1 .* ...
                              white_body .* (open <= low_1) .* (close > body_midpoint_1) .* (close < open_1) ;
## 07 - Three Inside Up
## First candle
##     a candle in a downtrend
##     black body
## Second candle
##     white body
##     the candle body is engulfed by the prior candle body
## Third candle
##     the closing price is above the previous closing price
bullish_signal_matrix(:,7) = downtrend_2 .* black_body_2 .* ...
                              white_body_1 .* (body_high_1 < body_high_2) .* (body_low_1 > body_low_2) .* ...
                              (close > close_1) ;
## 08 - Three Outside Up
## First candle
##     a candle in a downtrend
##     black body
## Second candle
##     white body
##     candle’s body engulfs the prior (black) candle’s body
## Third candle
##     closing price above the previous closing price
##     white body
bullish_signal_matrix(:,8) = downtrend_2 .* black_body_2 .* ...
                              white_body_1 .* (body_high_1 > body_high_2) .* (body_low_1 < body_low_2) .* ...
                              (close > close_1) .* white_body ;
## 09 - Three White Soldiers
## First candle
##     a candle in a downtrend
##     white body
## Second candle
##     white body
##     the opening price within the previous body
##     the closing price above the previous closing price
## Third candle
##     white body
##     the opening price within the previous body
##     the closing price above the previous closing price
bullish_signal_matrix(:,9) = downtrend_2 .* white_body_2 .* ...
                              white_body_1 .* (open_1 < body_high_2) .* (open_1 > body_low_2) .* (close_1 > close_2) .* ...
                              white_body .* (open < body_high_1) .* (open > body_low_1) .* (close > close_1) ; 
## 10 - Breakaway
## First candle
##     a tall black candle
## Second candle
##     a black candle
##     candle opens below the previous closing price (downward price gap, shadows can overlap)
## Third candle
##     a white or black candle
##     candle opens below the previous opening price
## Fourth candle
##     a black candle
##     candle closes below the previous closing price
## Fifth candle
##     a tall white candle
##     candle opens above the previous closing price
##     candle closes above the second line's opening price and below the first line's opening price
bullish_signal_matrix(:,10) = black_body_4 .* long_line_4 .*...
                               black_body_3 .* (open_3 < close_4) .* ...
                               (open_2 < open_3) .* ...
                               black_body_1 .* (close_1 < close_2) .* ...
                               white_body .* long_line .* (open > close_1) .* (close > open_3) .* (close < open_4) ; 
## 11 - Counter Attack
## First candle
##     a candle in a downtrend
##     black body
## Second candle
##     white body
##     the opening price is lower than the previous closing price
##     the closing price is at or higher than the previous closing price
bullish_signal_matrix(:,11) = downtrend_1 .* black_body_1 .* ...
                               white_body .* (open < close_1) .* (close >= close_1) ;
## 12 - Doji Star
## First candle
##     a candle in a downtrend
##     black body
## Second candle
##     a doji candle
##     a body below the first candle's body
bullish_signal_matrix(:,12) = downtrend_1 .* black_body_1 .* ...
                               doji .* (open < body_low_1) ; 
## 13 - Dragonfly Doji
##     Opening, closing and maximum prices are the same or very similar
##     Long lower shadow
##     appears on as a long line
bullish_signal_matrix(:,13) = (open == close) .* (close == high) .* long_line ;

## 14 - Engulfing
## First candle
##     a candle in a downtrend
##     black body
## Second candle
##     white body
##     candle's body engulfs the prior (black) candle's body
bullish_signal_matrix(:,14) = downtrend_1 .* black_body_1 .* ...
                               white_body .* (body_high > body_high_1) .* (body_low < body_low_1) ; 
## 15 - Gravestone Doji
##     Opening, closing and minimum prices are the same or very similar
##     Long upper shadow
##     appears on as a long line
bullish_signal_matrix(:,15) = (open == close) .* (close == low) .* long_line ; 

## 16 - Harami Cross
## First candle
##     a candle in a downtrend
##     black body
##     appears on as a long line
## Second candle
##     a doji candle with two shadows
##     the candle (including shadows) is engulfed by the previous candle's body
bullish_signal_matrix(:,16) = downtrend_1 .* black_body_1 .* long_line_1 .* ...
                               doji .* (high > close) .* (low < close) .* (high < body_high_1) .* (low > body_low_1) ;
## 17 - Homing Pigeon
## First candle
##     a candle in a downtrend
##     black body
## Second candle
##     black body
##     candle’s body engulfed by the prior candle’s body
bullish_signal_matrix(:,17) = downtrend_1 .* black_body_1 .* ...
                               black_body .* (body_high < body_high_1) .* (body_low > body_low_1) ;
## 18 - Ladder Bottom
## First candle
##     a tall black candle
## Second candle
##     a tall black candle
##     candle opens below the previous opening price
## Third candle
##     a tall black candle
##     candle opens below the previous opening price
## Fourth candle
##     a black candle
##     candle closes below the previous closing price
##     candle has a long upper shadow
## Fifth candle
##     a tall white candle
##     candle opens above the previous opening price
##     candle closes at or above the third line's opening price and below the first line's opening price
bullish_signal_matrix(:,18) = black_body_4 .* long_line_4 .* ...
                               black_body_3 .* (open_3 < open_4) .* long_line_3 .* ...
                               black_body_2 .* (open_2 < open_3) .* long_line_2 .* ...
                               black_body_1 .* (close_1 < close_2) .* ((high_1 .- body_high_1) > body_range_1) .* ...
                               white_body .* (open > open_1) .* (close >= open_2) .* (close < open_5) ; 
## 19 - Long Legged Doji
##     a doji candle
##     opening and closing prices are the same or similar
##     upper and lower shadow are very long
##     body is located in the middle of the candle or nearly mid-range
##     appears on as a long line
bullish_signal_matrix(:,19) = doji .* ((high .- close) > 0) .* ((close .- low) > 0) .* long_line ;

## 20 - Matching Low
## First candle
##     a candle in a downtrend
##     black body
##     no lower shadow
##     appears as a long line
## Second candle
##     black body
##     the opening price is below the previous opening price
##     the closing price is at the level of the previous closing price
##     no lower shadow
bullish_signal_matrix(:,20) = downtrend_1 .* black_body_1 .* (close_1 == low_1) .* long_line_1 .* ...
                               black_body .* (open < open_1) .* (close == close_1) .* (close == low) ;
## 21 - Meeting Lines
## First candle
##     a candle in a downtrend
##     black body
##     appears as a long line
## Second candle
##     white body
##     the closing price is equal to the previous closing price
bullish_signal_matrix(:,21) = downtrend_1 .* black_body_1 .* long_line_1 .* ...
                               white_body .* (close == close_1) ;
## 22 - Stick Sandwich
## First candle
##     a candle in a downtrend
##     black body
## Second candle
##     white body
##     the opening price is higher than the previous closing price
##     the closing price is higher than the previous opening price
## Third candle
##     black_body
##     the opening price is higher than the previous closing price
##     the closing price is equal to or higher than first line close
bullish_signal_matrix(:,22) = downtrend_2 .* black_body_2 .* ...
                               white_body_1 .* (open_1 > close_2) .* (close_1 > open_2) .* ...
                               black_body .* (open > close_1) .* (close >= close_2) ;
## 23 - Three Stars in the South
## First candle
##     a candle in a downtrend
##     black body
##     long lower shadow
## Second candle
##     black body
##     the opening below the prior opening
##     the closing below or at the prior closing
##     the low above the prior low
## Third candle
##     a marubozu candle with black body
##     appears as a short line
##     a candle is located within the prior candle
bullish_signal_matrix(:,23) = downtrend_2 .* black_body_2 .* ((body_low_2 .- low_2) > body_range_2) .* ...
                               black_body_1 .* (open_1 < open_2) .* (close_1 <= close_2) .* (low_1 > low_2) .* ...
                               black_marubozu .* (high < high_1) .* (low > low_1) ;
## 24 - Tri Star
## First candle
##     a doji candle in a downtrend
## Second candle
##     a doji candle
##     a body below the prior body
## Third candle
##     a doji candle
##     a body above the prior body
bullish_signal_matrix(:,24) = downtrend_2 .* doji_2 .* ...
                               doji_1 .* (open_1 < close_2) .* ...
                               doji .* (open > close_1) ;
## 25 - Unique Three River Bottom
## First candle
##     black candle
##     a candle in a downtrend
## Second candle
##     black candle
##     a body within the prior body
##     the lower shadow is at least twice longer than the body
##     the low price below the prior low price
## Third candle
##     white candle
##     a body located below the prior body
##     the low price above the prior low price
bullish_signal_matrix(:,25) = downtrend_2 .* black_body_2 .* ...
                               black_body_1 .* (body_high_1 < body_high_2) .* (body_low_1 > body_low_2) .* ((body_low_1 .- low_1) >= body_range_1) .* (low_1 < low_2) .* ...
                               white_body .* (body_high < body_low_1) .* (low > low_1) ;
## 26 - Belt Hold
##     white body
##     no lower shadow
##     short upper shadow
##     appears as a long line
bullish_signal_matrix(:,26) = downtrend .* white_body .* (open == low) .* long_line ;

## 27 - Hammer
##     downtrend
##     white or black candle with a small body
##     no upper shadow or the shadow cannot be longer than the body
##     lower shadow between two to three times longer than the body
##     if the gap is created at the opening or at the closing, it makes the signal stronger
##     appears as a long line
bullish_signal_matrix(:,27) = downtrend .* (doji == 0) .* ((high .- body_high) < body_range) .* ((body_low .- low) > 2.*body_range) .* long_line ;

## 28 - Harami
## First candle
##     a candle in a downtrend
##     black body
## Second candle
##     white body
##     candle's body engulfed by the prior (black) candle's body
bullish_signal_matrix(:,28) = downtrend_1 .* black_body_1 .* ...
                               white_body .* (body_high < body_high_1) .* (body_low > body_low_1) ;
## 29 - Inverted Hammer 
## First candle
##     a candle in a downtrend
##     black body
## Second candle
##     a white or black candle with a small body
##     no lower shadow or the shadow cannot be longer than then body
##     upper shadow at least 2.5 times longer than the body
##     the open below or at the level of the previous candle's close
bullish_signal_matrix(:,29) = downtrend_1 .* black_body_1 .* ...
                               (doji == 0 ) .* ((body_low .- low) < body_range) .* ((high .- body_high) >= 2.5.*body_range) .* (open <= close_1) ;
                               
bullish_signal_matrix( 1 : 10 , : ) = 0 ;

endfunction

and Bearish Reversal Indicator

## Copyright (C) 2020 dekalog
## 
## This program is free software: you can redistribute it and/or modify it
## under the terms of the GNU General Public License as published by
## the Free Software Foundation, either version 3 of the License, or
## (at your option) any later version.
## 
## This program is distributed in the hope that it will be useful, but
## WITHOUT ANY WARRANTY; without even the implied warranty of
## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
## GNU General Public License for more details.
## 
## You should have received a copy of the GNU General Public License
## along with this program.  If not, see
## .

## -*- texinfo -*- 
## @deftypefn {} {@var{bearish_signal_matrix} =} candle_bearish_reversal (@var{high}, @var{low}, @var{close}, @var{open}, @var{uptrend}, @var{avge_hi_lo_range})
##
## The HIGH, LOW, CLOSE and OPEN should be vectors of equal length, and are required.
##
## Internally the UPTREND is determined as a bar's CLOSE being higher than the CLOSE
## of the bar 5 bars previously. This can be over-ruled by an optional, user supplied
## vector UPTREND of zeros and ones, where one(s) represent a bar in an UPTREND.
## If an UPTREND vector consists solely of ones, this effectively turns off the
## discriminative ability of the UPTREND condition as all bars will be considered
## to be UPTREND bars. If the UPTREND vector consists solely of zeros, no bar
## will be classified as being in an UPTREND and those candle patterns that require an 
## UPTREND as a condition will not be indicated.
##
## The AVGE_HI_LO_RANGE is calculated as a rolling 5 bar simple moving average of
## the HI-LO range of bars. A bar is considered a "long line" if its HI-LO range
## is greater than the AVGE_HI_LO_RANGE. This can also be over-ruled, as above, by
## an optional, user supplied vector AVGE_HI_LO_RANGE, which consists of zeros and
## ones, with ones representing "long line" bars. An AVGE_HI_LO_RANGE vector of all
## ones or zeros will have a similar effect as described for the UPTREND vector.
##
## The candlestick pattern descriptions are predominantly taken from
##
## https://www.candlescanner.com
##
## and the bearish only patterns implemented are:
##
## Bearish Reversal High Reliability
##
##     01 - Bearish Abandoned Baby
##
##     02 - Dark Cloud Cover
##
##     03 - Evening Doji Star
##
##     04 - Evening Star
##
##     05 - Kicking
##
##     06 - Three Black Crows
##
##     07 - Three Inside Down
##
##     08 - Three Outside Down
##
##     09 - Upside Gap Two Crows
##
## Bearish Reversal Moderate Reliability
##
##     10 - Advance Block
##
##     11 - Breakaway
##
##     12 - Counter Attack
##
##     13 - Deliberation
##
##     14 - Doji Star
##
##     15 - Dragonfly Doji
##
##     16 - Engulfing
##
##     17 - Gravestone Doji
##
##     18 - Harami Cross
##
##     19 - Identical Three Crows
##
##     20 - Long Legged Doji
##
##     21 - Meeting Lines
##
##     22 - Tri Star
##
##     23 - Two Crows
##
## Bearish Reversal Low Reliability
##
##     24 - Belt Hold 
##
##     25 - Hanging Man
##
##     26 - Harami
##
##     27 - Shooting Star
##
## The output BEARISH_SIGNAL_MATRIX has a row length the same as the OHLC
## input vectors and 27 columns. The matrix is a zero filled matrix with the
## value 1 when a candle pattern is indicated. The row ix of the value is the
## row ix of the last candle in the pattern. The column ix corresponds to the
## pattern index given above.
##
## @seealso{candle_bullish_reversal}
## @end deftypefn

## Author: dekalog 
## Created: 2020-08-10

function bearish_signal_matrix = candle_bearish_reversal ( varargin )

if ( nargin < 4 || nargin > 6 )
  print_usage () ;
elseif ( nargin == 4 )
  high = varargin{1} ; low = varargin{2} ; close = varargin{3} ; open = varargin{4} ;
  uptrend = close > shift( close , 5 ) ;
  avge_hi_lo_range = filter( [0.2 ; 0.2 ; 0.2 ; 0.2 ; 0.2 ] , 1 , high .- low ) ; ## 5 bar simple moving average
elseif ( nargin == 5 )
  high = varargin{1} ; low = varargin{2} ; close = varargin{3} ; open = varargin{4} ;
  uptrend = varargin{5} ;
  avge_hi_lo_range = filter( [0.2 ; 0.2 ; 0.2 ; 0.2 ; 0.2 ] , 1 , high .- low ) ; ## 5 bar simple moving average
elseif ( nargin == 6 )
  high = varargin{1} ; low = varargin{2} ; close = varargin{3} ; open = varargin{4} ;
  uptrend = varargin{5} ; avge_hi_lo_range = varargin{6} ;
endif

bearish_signal_matrix = zeros( size( close , 1 ) , 27 ) ;

## pre-compute some basic vectors for re-use below, trading memory use for
## speed of execution
uptrend_1 = shift( uptrend , 1 ) ;
uptrend_2 = shift( uptrend , 2 ) ;
long_line = ( high .- low ) > shift( avge_hi_lo_range , 1 ) ;
long_line_1 = shift( long_line , 1 ) ;
long_line_4 = shift( long_line , 4 ) ;
open_1 = shift( open , 1 ) ;
open_2 = shift( open , 2 ) ;
open_3 = shift( open , 3 ) ;
high_1 = shift( high , 1 ) ;
high_2 = shift( high , 2 ) ;
high_4 = shift( high , 4 ) ;
low_1 = shift( low , 1 ) ;
close_1 = shift( close , 1 ) ;
close_2 = shift( close , 2 ) ;
close_4 = shift( close , 4 ) ;
black_body = close < open ;
black_body_1 = shift( black_body , 1 ) ;
black_body_2 = shift( black_body , 2 ) ;
white_body = close > open ;
white_body_1 = shift( white_body , 1 ) ;
white_body_2 = shift( white_body , 2 ) ;
white_body_3 = shift( white_body , 3 ) ;
white_body_4 = shift( white_body , 4 ) ;
doji = ( close == open ) ;
doji_1 = shift( doji , 1 ) ;
doji_2 = shift( doji , 2 ) ;
body_high = max( [ open , close ] , [] , 2 ) ;
body_high_1 = shift( body_high , 1 ) ;
body_high_2 = shift( body_high , 2 ) ;
body_low = min( [ open , close ] , [] , 2 ) ;
body_low_1 = shift( body_low , 1 ) ;
body_low_2 = shift( body_low , 2 ) ;
body_midpoint = ( body_high .+ body_low ) ./ 2 ;
body_midpoint_1 = shift( body_midpoint , 1 ) ;
body_midpoint_2 = shift( body_midpoint , 2 ) ;
candle_midpoint = ( high .+ low ) ./ 2 ;
candle_midpoint_1 = shift( candle_midpoint , 1 ) ; 
body_range = body_high .- body_low ;
wick_gap_up = low > shift( high , 1 ) ;
wick_gap_up_1 = shift( wick_gap_up , 1 ) ;
wick_gap_down = high < shift( low , 1 ) ;
black_marubozu = ( high == open ) .* ( low == close ) ;
white_marubozu = ( high == close ) .* ( low == open ) ;
white_marubozu_1 = shift( white_marubozu , 1 ) ;

## 01 - Bearish Abandoned Baby
## First candle
##     a candle in an uptrend
##     white body
## Second candle
##     a doji candle
##     the low price above the prior high price
## Third candle
##     black body
##     the high price below the prior low price
bearish_signal_matrix(:,1) = uptrend_2 .* white_body_2 .* ...
                              doji_1 .* wick_gap_up_1 .* ...
                              black_body .* wick_gap_down ; 
## 02 - Dark Cloud Cover
## First candle
##     a candle in an uptrend
##     white body
## Second candle
##     black body
##     the opening above or equal of the prior high
##     the closing below the midpoint of the prior candle
##     the closing above the previous opening
bearish_signal_matrix(:,2) = uptrend_1 .* white_body_1 .* ...
                              black_body .* (open >= high_1) .* (close < candle_midpoint_1) .* (close > open_1) ;  
## 03 - Evening Doji Star
## First candle
##     a candle in an uptrend
##     white body
## Second candle
##     a doji candle
##     a doji body above the previous candle body
##     the low price below the previous candle high price
## Third candle
##     black body
##     candle body below the previous candle body
##     the closing price below the midpoint of the first candle body
bearish_signal_matrix(:,3) = uptrend_2 .* white_body_2 .* ...
                              doji_1 .* (open_1 > body_high_2) .* (low_1 < high_2) .* ...
                              black_body .* (body_high < close_1) .* (close < body_midpoint_2) ;
## 04 - Evening Star
## First candle
##     a candle in an uptrend
##     white body
## Second candle
##     white or black body
##     the candle body is located above the prior body
## Third candle
##     black body
##     the candle body is located below the prior body
##     the candle closes at least halfway down the body of the first line
bearish_signal_matrix(:,4) = uptrend_2 .* white_body_2 .* ...
                              (body_low_1 > body_high_2) .* ...
                              black_body .* (body_high < body_low_1) .* (close <= body_midpoint_2) ;
## 05 - Kicking
## First candle
##     a White Marubozu
##     appears on as a long line
## Second candle
##     a Black Marubozu
##     price gaps downward
##     appears on as a long line
bearish_signal_matrix(:,5) = white_marubozu_1 .* long_line_1 .* ...
                              black_marubozu .* (high < low_1) .* long_line ; 
## 06 - Three Black Crows
## First candle
##     a candle in an uptrend
##     black body
## Second candle
##     black body
##     the opening price within the previous body
##     the closing price below the previous closing price
## Third candle
##     black body
##     the opening price within the previous body
##     the closing price below the previous closing price
bearish_signal_matrix(:,6) = uptrend_2 .* black_body_2 .* ...
                              black_body_1 .* (open_1 < body_high_2) .* (open_1 > body_low_2) .* (close_1 < close_2) .* ...
                              black_body .* (open < body_high_1) .* (open > body_low_1) .* (close < close_1) ; 
## 07 - Three Inside Down
## First candle
##     a candle in an uptrend
##     white body
## Second candle
##     black body
##     the candle body is engulfed by the prior candle body
## Third candle
##     black_body
##     the closing price is below the previous closing price
bearish_signal_matrix(:,7) = uptrend_2 .* white_body_2 .* ...
                              black_body_1 .* (body_high_1 < body_high_2) .* (body_low_1 > body_low_2) .* ...
                              black_body .* (close < close_1) ;
## 08 - Three Outside Down
## First candle
##     a candle in an uptrend
##     white body
## Second candle
##     black body
##     candle’s body engulfs the prior (white) candle’s body
## Third candle
##     closing price below the previous closing price
##     black body
bearish_signal_matrix(:,8) = uptrend_2 .* white_body_2 .* ...
                              black_body_1 .* (body_high_1 > body_high_2) .* (body_low_1 < body_low_2) .* ...
                              black_body .* (close < close_1) ;
## 09 - Upside Gap Two Crows
## First candle
##     a candle in an uptrend
##     white body
## Second candle
##     black body
##     candle's body above the previous candle's body
## Third candle
##     black body
##     candle’s body engulfs the prior candle’s body
bearish_signal_matrix(:,9) = uptrend_2 .* white_body_2 .* ...
                              black_body_1 .* (body_low_1 > body_high_2) .* ...
                              black_body .* (body_high > body_high_1) .* (body_low < body_low_1) ;
## 10 - Advance Block
## First candle
##     a candle in an uptrend
##     white body
## Second candle
##     white body
##     the opening price is within the previous body
##     the closing price is above the previous closing price
## Third candle
##     white body
##     the opening price is within the previous body
##     the closing price is above the previous closing price
bearish_signal_matrix(:,10) = uptrend_2 .* white_body_2 .* ...
                               white_body_1 .* (open_1 < body_high_2) .* (open_1 > body_low_2) .* (close_1 > close_2) .* ...
                               white_body .* (open < body_high_1) .* (open > body_low_1) .* (close > close_1) ;
## 11 - Breakaway
## First candle
##     a tall white candle
## Second candle
##     a white candle
##     candle opens above the previous closing price (upward price gap, shadows can overlap)
## Third candle
##     a white or black candle
##     candle opens above the previous opening price
## Fourth candle
##     a white candle
##     candle closes above the previous closing price
## Fifth candle
##     a tall black candle
##     candle opens below the previous closing price
##     candle closes below the second line's opening price and above the first line's closing price
##     the price gap formed between the first and the second line is not closed
bearish_signal_matrix(:,11) = white_body_4 .* long_line_4 .* ...
                               white_body_3 .* (open_3 > close_4) .* ...
                               (open_2 > open_3 ) .* ...
                               white_body_1 .* (close_1 > close_2) .* ...
                               black_body .* long_line .* (open < close_1) .* (close < open_3) .* (close > close_4) .* (low > high_4) ;
## 12 - Counter Attack
## First candle
##     a candle in an uptrend
##     white body
## Second candle
##     black_body
##     the opening price is higher than the previous closing price
##     the closing price is at or lower than the previous closing price
bearish_signal_matrix(:,12) = uptrend_1 .* white_body_1 .* ...
                               black_body .* (open > close_1) .* (close <= close_1) ;
## 13 - Deliberation
## First candle
##     a candle in an uptrend
##     white body
## Second candle
##     white body
##     the opening price is above the previous opening price
##     the closing price is above the previous closing price
## Third candle
##     white body
##     the opening price is slightly lower or higher than the previous closing price
##     the closing price is above the previous closing price
bearish_signal_matrix(:,13) = uptrend_2 .* white_body_2 .* ...
                               white_body_1 .* (open_1 > open_2) .* (close_1 > close_2) .* ...
                               white_body .* (open > body_midpoint_1) .* (close > close_1) ;
## 14 - Doji Star
## First candle
##     a candle in an uptrend
##     white body
## Second candle
##     a doji candle
##     a body above the first candle's body
bearish_signal_matrix(:,14) = uptrend_1 .* white_body_1 .* ...
                               doji .* (open > close_1) ; 
## 15 - Dragonfly Doji
##     Opening, closing and maximum prices are the same or very similar
##     Long lower shadow
##     appears on as a long line
bearish_signal_matrix(:,15) = (open == close) .* (close == high) .* long_line ;

## 16 - Engulfing
## First candle
##     a candle in an uptrend
##     white body
## Second candle
##     black body
##     candle's body engulfs the prior (white) candle's body
bearish_signal_matrix(:,16) = uptrend_1 .* white_body_1 .* ...
                               black_body .* (body_high > body_high_1) .* (body_low < body_low_1) ;
## 17 - Gravestone Doji
##     Opening, closing and minimum prices are the same or very similar
##     Long upper shadow
##     appears on as a long line
bearish_signal_matrix(:,17) = (open == close) .* (close == low) .* long_line ;

## 18 - Harami Cross
## First candle
##     a candle in an uptrend
##     white body
## Second candle
##     a doji candle with two shadows
##     candle's body engulfed by the prior (white) candle's body
bearish_signal_matrix(:,18) = uptrend_1 .* white_body_1 .* ...
                               doji .* (high > close) .* (low < close) .* (high < body_high_1) .* (low > body_low_1) ;
## 19 - Identical Three Crows
## First candle
##     a candle in an uptrend
##     black body
## Second candle
##     black body
##     the opening price at or near the prior close
## Third candle
##     black body
##     the opening price at or near the prior close
bearish_signal_matrix(:,19) = uptrend_2 .* black_body_2 .* ...
                               black_body_1 .* (open_1 == close_2) .* ...
                               black_body .* (open == close_1) ;
## 20 - Long Legged Doji
##     a doji candle
##     opening and closing prices are the same or similar
##     upper and lower shadow are very long
##     body is located in the middle of the candle or nearly mid-range
##     appears on as a long line
bearish_signal_matrix(:,20) = doji .* ((high .- close) > 0) .* ((close .- low) > 0) .* long_line ;

## 21 - Meeting Lines
## First candle
##     a candle in an uptrend
##     white body
##     appears as a long line
## Second candle
##     black body
##     the closing price is equal to the previous closing price
bearish_signal_matrix(:,21) = uptrend_1 .* white_body_1 .* long_line_1 .* ...
                               black_body .* (close == close_1) ;
## 22 - Tri Star
## First candle
##     a doji candle in an uptrend
## Second candle
##     a doji candle
##     a body above the prior body
## Third candle
##     a doji candle
##     a body below the prior body
bearish_signal_matrix(:,22) = uptrend_2 .* doji_2 .* ...
                               doji_1 .* (open_1 > close_2) .* ...
                               doji .* (open < close_1) ;
## 23 - Two Crows
## First candle
##     a candle in an uptrend
##     white body
## Second candle
##     black body
##     the closing price above the prior closing price (gap between bodies)
## Third candle
##     black body
##     the opening price within the prior body
##     the closing price within the body of the first line (gap close)
bearish_signal_matrix(:,23) = uptrend_2 .* white_body_2 .* ...
                               black_body_1 .* (close_1 > close_2) .* ...
                               black_body .* (open < body_high_1) .* (open > body_low_1) .* (close < close_2) ;
## 24 - Belt Hold
##     black body
##     no upper shadow
##     short lower shadow
##     appears as a long line
bearish_signal_matrix(:,24) = black_body .* (open == high) .* (low < close) .* long_line ;

## 25 - Hanging Man
##     uptrend
##     white or black candle with a small body
##     no upper shadow or the shadow cannot be longer than the body
##     lower shadow at least two times longer than the body
##     if the gap is created at the opening or at the closing, it makes the signal stronger
##     appears as a long line
##     the body fully located above the trendline
bearish_signal_matrix(:,25) = uptrend .* (doji == 0) .* ((high .- body_high) < body_range) .* ((body_low .- low) >= 2.*body_range) ; 

## 26 - Harami
## First candle
##     a candle in an uptrend
##     white body
## Second candle
##     black body
##     candle's body engulfed by the prior (white) candle's body
bearish_signal_matrix(:,26) = uptrend_1 .* white_body_1 .* ...
                               black_body .* (body_high < body_high_1) .* (body_low > body_low_1) ;
## 27 - Shooting Star
##     white or black candle with a small body
##     no lower shadow or the shadow cannot be longer than the body
##     upper shadow at least two times longer than the body
##     if the gap is created at the opening or the closing, it makes the signal stronger
##     appears as a long line
bearish_signal_matrix(:,27) = (doji == 0) .* ((body_low .- low) < body_range) .* ((high .- body_high) >= 2.*body_range) ;

bearish_signal_matrix( 1 : 10 , : ) = 0 ;

endfunction