NetworkMomentum/Include/NetworkMomentum/NM_DDTW.mqh
2026-07-30 03:34:33 +00:00

229 lines
7.4 KiB
MQL5

//+------------------------------------------------------------------+
//| NM_DDTW.mqh |
//| MMQ — Muhammad Minhas Qamar |
//| www.mql5.com/en/articles/23763 |
//+------------------------------------------------------------------+
#property copyright "MMQ — Muhammad Minhas Qamar"
#property link "https://www.mql5.com/en/articles/23763"
#property version "1.00"
#property strict
#include <NetworkMomentum\NM_Matrix.mqh>
//+------------------------------------------------------------------+
//| Lead-lag detection by Derivative Dynamic Time Warping. |
//| |
//| Section 2.2 of the paper. For a T x M block of returns it |
//| builds the skew-symmetric M x M matrix V of integer lags: how |
//| many steps series j must shift to best align with series i. |
//| |
//| Three stages: |
//| |
//| 1. Standardise each column to zero mean, unit population std |
//| (a zero std is replaced by 1, not skipped). |
//| |
//| 2. Replace each column by its derivative (Equation 5): |
//| DX[i] = ((X[i]-X[i-1]) + (X[i+1]-X[i-1])/2) / 2 |
//| for interior i, with the two ends copied inward. |
//| |
//| 3. For every pair, run an UNCONSTRAINED DTW on the derivative |
//| series and read the lag as the mode of (j - i) along the |
//| warping path. |
//| |
//| The DTW tie-breaking is fixed so the integer lags are |
//| deterministic rather than merely optimal: |
//| |
//| - forward cost D[i][j] = (a-b)^2 + min(diag, up, left), the |
//| squared local cost with an inf border and D[0][0]=0; |
//| |
//| - the path is recovered by argmin over [diag, up, left] with |
//| the FIRST minimum winning, so ties resolve diag > up > left; |
//| |
//| - the lag mode breaks its own ties toward the smaller lag. |
//+------------------------------------------------------------------+
class CNMDDTW
{
private:
//--- Equation 5 derivative of one series into out[].
static void Derivative(const double &x[],const int n,double &out[]);
//--- warping-path lag between two equal-length series.
static int WarpingLag(const double &x[],const double &y[],const int n);
public:
//--- skew-symmetric M x M lead-lag matrix from a T x M return block.
static void LeadLagMatrix(const NMMatrix &returns,NMMatrix &V);
};
//+------------------------------------------------------------------+
//| Equation 5 derivative, ends copied inward. |
//+------------------------------------------------------------------+
void CNMDDTW::Derivative(const double &x[],const int n,double &out[])
{
ArrayResize(out,n);
if(n<3)
{
for(int i=0;i<n;i++)
out[i]=0.0;
return;
}
for(int i=1;i<n-1;i++)
out[i]=((x[i]-x[i-1])+(x[i+1]-x[i-1])/2.0)/2.0;
out[0]=out[1];
out[n-1]=out[n-2];
}
//+------------------------------------------------------------------+
//| DTW lag between two equal-length series. |
//| Fills the (n+1)x(n+1) accumulated-cost matrix, back-tracks the |
//| optimal path, then returns the mode of (j-i) over that path. |
//+------------------------------------------------------------------+
int CNMDDTW::WarpingLag(const double &x[],const double &y[],const int n)
{
const int W=n+1;
double D[];
ArrayResize(D,W*W);
for(int k=0;k<W*W;k++)
D[k]=DBL_MAX;
D[0]=0.0;
//--- forward pass: squared local cost plus the cheapest predecessor.
for(int i=1;i<W;i++)
{
for(int j=1;j<W;j++)
{
double diff=x[i-1]-y[j-1];
double c=diff*diff;
double diag=D[(i-1)*W+(j-1)];
double up =D[(i-1)*W+j];
double left=D[i*W+(j-1)];
double best=diag;
if(up<best)
best=up;
if(left<best)
best=left;
D[i*W+j]=c+best;
}
}
//--- back-track from the corner; first minimum wins (diag > up > left).
int lag_diffs[];
ArrayResize(lag_diffs,2*W);
int cnt=0;
int i=n,j=n;
lag_diffs[cnt++]=j-i;
while(i>0 && j>0)
{
double diag=D[(i-1)*W+(j-1)];
double up =D[(i-1)*W+j];
double left=D[i*W+(j-1)];
int c=0; // 0 diag, 1 up, 2 left
double m=diag;
if(up<m)
{
m=up;
c=1;
}
if(left<m)
{
m=left;
c=2;
}
if(c==0)
{
i--;
j--;
}
else
if(c==1)
i--;
else
j--;
if(i>0 && j>0)
lag_diffs[cnt++]=j-i;
}
//--- mode of the collected (j-i), ties resolved toward the smaller lag.
int best_lag=0;
int best_count=-1;
for(int a=0;a<cnt;a++)
{
int cand=lag_diffs[a];
int count=0;
for(int b=0;b<cnt;b++)
if(lag_diffs[b]==cand)
count++;
if(count>best_count || (count==best_count && cand<best_lag))
{
best_count=count;
best_lag=cand;
}
}
return best_lag;
}
//+------------------------------------------------------------------+
//| Skew-symmetric M x M lead-lag matrix (Section 2.2). |
//+------------------------------------------------------------------+
void CNMDDTW::LeadLagMatrix(const NMMatrix &returns,NMMatrix &V)
{
const int T=returns.rows;
const int M=returns.cols;
V.Init(M,M);
//--- standardise every column, zero std replaced by 1.
NMMatrix z;
z.Init(T,M);
for(int c=0;c<M;c++)
{
double s=0.0;
for(int t=0;t<T;t++)
s+=returns.Get(t,c);
double mean=s/T;
double v=0.0;
for(int t=0;t<T;t++)
{
double d=returns.Get(t,c)-mean;
v+=d*d;
}
double sd=MathSqrt(v/T);
if(sd==0.0)
sd=1.0;
for(int t=0;t<T;t++)
z.Set(t,c,(returns.Get(t,c)-mean)/sd);
}
//--- derivative of each standardised column (Equation 5).
NMMatrix dz;
dz.Init(T,M);
double col[],dcol[];
ArrayResize(col,T);
for(int c=0;c<M;c++)
{
for(int t=0;t<T;t++)
col[t]=z.Get(t,c);
Derivative(col,T,dcol);
for(int t=0;t<T;t++)
dz.Set(t,c,dcol[t]);
}
//--- pairwise DTW lag, written skew-symmetrically.
double xi[],yj[];
ArrayResize(xi,T);
ArrayResize(yj,T);
for(int a=0;a<M;a++)
{
for(int b=a+1;b<M;b++)
{
for(int t=0;t<T;t++)
{
xi[t]=dz.Get(t,a);
yj[t]=dz.Get(t,b);
}
int lag=WarpingLag(xi,yj,T);
V.Set(a,b,lag);
V.Set(b,a,-lag);
}
}
}
//+------------------------------------------------------------------+