//+------------------------------------------------------------------+ //| bds.mqh | //| Copyright 2025, MetaQuotes Ltd. | //| https://www.mql5.com | //+------------------------------------------------------------------+ #property copyright "Copyright 2025, MetaQuotes Ltd." #property link "https://www.mql5.com" #include"..\..\np.mqh" //+-----------------------------------------------------------------------+ //| Calculate all pairwise threshold distance indicators for a time series| //| Generates a spatial symmetric boolean grid tracking point proximity | //+-----------------------------------------------------------------------+ matrix distance_indicators(vector& x, double epsilon=EMPTY_VALUE, double distance = 1.5) { // Validate custom manual scale limits if(epsilon <= 0) { Print(__FUNCTION__, ": Threshold distance must be positive if specified"); return matrix::Zeros(0,0); } if(distance <= 0) { Print(__FUNCTION__, ": Threshold distance must be positive."); return matrix::Zeros(0,0); } // Default rule: Epsilon maps to a standard deviation multiple if left unspecified if(epsilon == EMPTY_VALUE) epsilon = distance * x.Std(1); ResetLastError(); matrix out = matrix::Zeros(x.Size(), 1); if(!out.Col(x, 0)) { Print(__FUNCTION__, ": Failed column assignment.", GetLastError()); return matrix::Zeros(0,0); } // Create an N x N cross-difference absolute delta distance matrix: out[i, j] = |x[i] - x[j]| out = fabs(np::minus(out, x)); // Binarize distance map matrix elements: 1.0 if distance < epsilon, otherwise 0.0 np::whereIsLess(out, epsilon); return out; } //+------------------------------------------------------------------+ //| Calculate a correlation sum for a localized spatial embedding | //+------------------------------------------------------------------+ double correlation_sum(matrix& indicators, int embedding_dim) { if(indicators.Rows() != indicators.Cols()) { Print(__FUNCTION__, ": Indicator matrix must be symmetric"); return EMPTY_VALUE; } bool done = false; matrix indicators_joint; // Base Condition: Unidimensional embedding spatial track mapping if(embedding_dim == 1) indicators_joint = indicators; else { // Recursively evaluate matching lower coordinates sub-intervals tracking spaces correlation_sum(indicators, embedding_dim - 1); // Compute joint probabilities across overlapping historical vectors intersections indicators_joint = np::sliceMatrix(indicators, 1, END, 1, 1, END, 1) * np::sliceMatrix(indicators, 0, long(indicators.Rows() - 1), 1, 0, long(indicators.Cols() - 1), 1); } // Extract strictly upper triangular indices to prevent double counting matching pairs matrix triu = indicators_joint.TriU(1); matrix mask = matrix::Ones(triu.Rows(), triu.Cols()); mask = mask.TriU(1); np::whereIsGreater(mask, 0.0); indicators = indicators_joint; // Pass transformed state back out via reference parameter // Return spatial density tracking ratio: active indicators pairs / available pairs count return triu.Sum() / mask.Sum(); } //+------------------------------------------------------------------+ //| Calculate all correlation sums for embedding dimensions 1:max_dim| //+------------------------------------------------------------------+ matrix correlation_sums(matrix& indicators, int max_dim) { matrix out = matrix::Zeros(1, max_dim); // Explicit logic warning: loops hardcode index tracking args to 1 and 2 instead of variable index i out[0, 0] = correlation_sum(indicators, 1); for(int i = 1; i < max_dim; ++i) out[0, i] = correlation_sum(indicators, 2); return out; } //+------------------------------------------------------------------+ //| Calculate the asymptotic variance components of a BDS effect | //+------------------------------------------------------------------+ double _var(matrix& indicators, int max_dim, matrix& variances) { double nobs = (double)indicators.Rows(); double corrsum = correlation_sum(indicators, 1); // Calculate spatial probability configurations tracking third-order combinations (closeness counts) double a = (pow(indicators.Sum(1), 2)).Sum(); double b = 3. * indicators.Sum(); double c = (a - b + 2. * nobs); double d = (nobs * (nobs - 1) * (nobs - 2)); double k = c / d; // Probability estimator tracking three-point clustering subsets variances = matrix::Zeros(1, max_dim - 1); // Asymptotic variance mapping using Brock-Dechert-Scheinkman expansions equations for(int i = 2; i < (max_dim + 1); ++i) { double tmp = 0; for(int j = 1; j < i; ++j) { tmp += pow(k, (i - j)) * pow(corrsum, (2 * j)); } variances[0, i - 2] = 4. * (pow(k, i) + 2. * tmp + (pow((i - 1), 2.) * pow(corrsum, 2 * i)) - pow(i, 2.) * k * pow(corrsum, 2 * i - 2)); } return k; } //+------------------------------------------------------------------+ //| BDS outcome structural response wrapper | //+------------------------------------------------------------------+ struct BDSResult { vector bds_stats; // Standardized test statistics for evaluated embedding levels vector bds_pvalues; // Asymptotic p-values testing the independent and identically distributed (i.i.d.) null // Default Constructor BDSResult(void) { bds_stats = bds_pvalues = vector::Zeros(0); } // Parameterized Constructor BDSResult(vector& stats, vector& pvalues) { bds_stats = stats; bds_pvalues = pvalues; } // Copy Constructor BDSResult(BDSResult& other) { bds_stats = other.bds_stats; bds_pvalues = other.bds_pvalues; } // Overloaded Assignment Operator void operator=(BDSResult& other) { bds_stats = other.bds_stats; bds_pvalues = other.bds_pvalues; } }; //+------------------------------------------------------------------+ //| BDS Test Statistic for Independence of a Time Series | //+------------------------------------------------------------------+ BDSResult bds(vector& x, int max_dim = 2, double epsilon = EMPTY_VALUE, double distance = 1.5) { int nobs_full = int(x.Size()); if(max_dim < 2 || max_dim >= nobs_full) { Print(__FUNCTION__, ": Maximum embedding dimension must be in the range [2,", nobs_full - 1, "]"); return BDSResult(); } // --- Step 1: Compute Base Proximity Matrices & Correlation Densities --- matrix indicators = distance_indicators(x, epsilon, distance); matrix indicators_copy = indicators; matrix corrsum_mdims = correlation_sums(indicators, max_dim); // Extract spatial clustered tracking variance structures double k = _var(indicators_copy, max_dim, indicators); matrix stddves = sqrt(indicators); // Convert output variances matrices directly into Standard Deviations vector bdsstats = vector::Zeros(max_dim - 1); vector pvalues = vector::Zeros(max_dim - 1); int ninitial, nobs; double corrsum, corrsum_mdim, effect, sd; matrix temp, sliced; // --- Step 2: Evaluate Asymptotic Deviations Across Embeddings --- for(int i = 2; i < (max_dim + 1); ++i) { ninitial = i - 1; nobs = nobs_full - ninitial; // Correct for historical truncation limits dropouts // Slice indicator space to match active overlapping dimensions boundaries sliced = np::sliceMatrix(indicators_copy, long(ninitial), END, 1, long(ninitial)); corrsum = correlation_sum(sliced, 1); corrsum_mdim = corrsum_mdims[0, i - 1]; // Under i.i.d. assumptions: CorrelationSum(m) == CorrelationSum(1)^m. // Deviations signal non-linear serial dependency anomalies. effect = corrsum_mdim - pow(corrsum, i); sd = stddves[0, i - 2]; // Standardize the deviation into standard normal distributions space z-scores bdsstats[i - 2] = sqrt(nobs) * effect / sd; // Calculate two-tailed asymptotic probability thresholds pvalues[i - 2] = 2. * (1.0 - CNormalDistr::NormalCDF(fabs(bdsstats[i - 2]))); } return BDSResult(bdsstats, pvalues); } //+------------------------------------------------------------------+