Article-22235-Multiple-Regr.../EconometricsA.mqh
alexeynikolaev2 2ba284c25f new files added
2026-07-09 07:53:01 +03:00

640 lines
22 KiB
MQL5

//+------------------------------------------------------------------+
//| EconometricsA.mqh |
//| Copyright 2000-2026, MetaQuotes Ltd. |
//| www.mql5.com |
//+------------------------------------------------------------------+
#include <Math\Stat\T.mqh>
#include <Math\Stat\ChiSquare.mqh>
#include <Math\Stat\F.mqh>
#include <Math\Stat\Math.mqh>
#include <Graphics\Graphic.mqh>
//---
//+------------------------------------------------------------------+
//| Computation of parameters and residuals |
//| for regression on a constant |
//+------------------------------------------------------------------+
void regression0(double& y[], double& a, double& c, double& e[])
{
//--- Initialization and input data validation
a = c = 0.0;
ArrayResize(e, 0);
int ny = ArraySize(y);
if(ny < 1)
{
Print("no data for regression0()");
return;
}
ArrayResize(e, ny);
//--- Parameters computation
a = MathMean(y);
c = MathStandardDeviation(y);
//--- Residuals computation
for(int i = 0; i < ny; ++i)
e[i] = y[i] - a;
}
//+------------------------------------------------------------------+
//| Calculation of confidence interval for parameter A |
//| in regression on a constant |
//+------------------------------------------------------------------+
void Aconf0(double a, double c, int n, double conf_level)
{
//--- Initialization and input data validation
if(c <= 0.0)
{
Print("c must be positive in Aconf0()");
return;
}
if(n < 2)
{
Print("n must be greater 1 in Aconf0()");
return;
}
if(conf_level <= 0.0 || conf_level >= 1.0)
{
Print("conf_level must be between 0 and 1 in Aconf0()");
return;
}
int err = 0;
//--- t-statistic calculation
//--- double t = MathQuantileT((1.0 - conf_level) / 2.0, n - 1, err); // faulty library function
double t = MathQuantileT_TMP((1.0 - conf_level) / 2.0, n - 1, err); // custom implementation
if(err != 0)
{
Print("in MathQuantileT() error ", err);
return;
}
t *= c / MathSqrt(n);
//--- Print result
PrintFormat("A between %.4f and %.4f with confidence level %.3f", a + t, a - t, conf_level);
}
//+------------------------------------------------------------------+
//| Test of hypothesis H₀: A = 0 vs H₁: A > 0 |
//| for regression on a constant |
//+------------------------------------------------------------------+
void Aeq0_test0(double a, double c, int n)
{
Print("test H0: A = 0 vs H1: A > 0 (for regession on constant) result:");
//--- Input data validation
if(c <= 0.0)
{
Print("c must be positive in Aeq0_test0()");
return;
}
if(n < 2)
{
Print("n must be greater 1 in Aeq0_test0()");
return;
}
//--- t-statistic calculation
double t = a * MathSqrt(n) / c;
int err = 0;
//--- p-value calculation
double p_value = MathCumulativeDistributionT(t, n - 1, false, false, err);
if(err != 0)
{
Print("in MathCumulativeDistributionT() error ", err);
return;
}
//--- Print result if success
PrintFormat("t = %.3f, p-value = %.3f", t, p_value);
}
//+------------------------------------------------------------------+
//| Calculation of prediction interval for Y |
//| in regression on a constant |
//+------------------------------------------------------------------+
void prognose_interval0(double a, double c, int n, double conf_level)
{
//--- Input data validation
if(c <= 0.0)
{
Print("c must be positive in prognose_interval0()");
return;
}
if(n < 2)
{
Print("n must be greater 1 in prognose_interval0()");
return;
}
if(conf_level <= 0.0 || conf_level >= 1.0)
{
Print("conf_level must be between 0 and 1 in prognose_interval0()");
return;
}
int err = 0;
//--- t-statistic calculation
//--- double t = MathQuantileT((1.0 - conf_level) / 2.0, n - 1, err); // faulty library function
double t = MathQuantileT_TMP((1.0 - conf_level) / 2.0, n - 1, err); // custom implementation
if(err != 0)
{
Print("in MathQuantileT() error ", err);
return;
}
t *= c * MathSqrt(1.0 + 1.0 / n);
//--- Print result if success
PrintFormat("Y between %.4f and %.4f with confidence level %.3f", a + t, a - t, conf_level);
}
//+------------------------------------------------------------------+
//| Computation of parameters and residuals for |
//| simple linear regression |
//+------------------------------------------------------------------+
void regression1(double& y[], double& x[], double& a, double& b, double& c, double& e[])
{
//--- Initialization and input data validation
a = b = c = 0.0;
ArrayResize(e, 0);
int ny = ArraySize(y);
if(ny < 3)
{
Print("no data for regression1()");
return;
}
if(ny != ArraySize(x))
{
Print("arrays y[] and x[] have different lengths in regression1()");
return;
}
//--- Check x[] is not constant
double Sx = MathStandardDeviation(x);
if(Sx <= DBL_MIN)
{
Print("x[]==const, use regression0()");
return;
}
//--- Parameters A and B estimation
if(!MathCorrelationPearson(y, x, b))
{
Print("MathCorrelationPearson() error in regression1()");
return;
}
ArrayResize(e, ny);
b *= MathStandardDeviation(y) / Sx;
a = MathMean(y) - b * MathMean(x);
//--- Residuals computation
for(int i = 0; i < ny; ++i)
e[i] = y[i] - a - b * x[i];
//--- Parameter C estimation
c = MathStandardDeviation(e) * MathSqrt(1.0 + 1.0 / (ny - 2.0));
}
//+------------------------------------------------------------------+
//| Test of hypothesis H₀: B = 0 vs H₁: B ≠ 0 for |
//| simple linear regression |
//+------------------------------------------------------------------+
void Beq0_test1(double b, double c, double var_x, int n)
{
Print("test H0: B = 0 vs H1: B != 0 (for pair regession) result:");
//--- Input data validation
if(c <= 0.0)
{
Print("c must be positive in Beq0_test1()");
return;
}
if(var_x <= 0.0)
{
Print("var_x must be positive in Beq0_test1()");
return;
}
if(n < 3)
{
Print("n must be 3 or greater in Beq0_test1()");
return;
}
//--- t-statistic calculation
double t = MathAbs(b * MathSqrt(var_x * (n - 1)) / c);
int err = 0;
//--- p-value computation
double p_value = 2.0 * MathCumulativeDistributionT(t, n - 2, false, false, err);
if(err != 0)
{
Print("in MathCumulativeDistributionT() error ", err);
return;
}
//--- Print result if success
PrintFormat("t = %.3f, p-value = %.3f", t, p_value);
}
//+------------------------------------------------------------------+
//| Calculation of prediction interval for Y in |
//| simple linear regression |
//+------------------------------------------------------------------+
void prognose_interval1(double x_new, double a, double b, double c, double var_x, double mean_x, int n, double conf_level)
{
//--- Input data validation
if(c <= 0.0)
{
Print("c must be positive in prognose_interval1()");
return;
}
if(n < 2)
{
Print("n must be greater 2 in prognose_interval1()");
return;
}
if(conf_level <= 0.0 || conf_level >= 1.0)
{
Print("conf_level must be between 0 and 1 in prognose_interval1()");
return;
}
int err = 0;
//--- t-statistic calculation
//--- double t = MathQuantileT((1.0 - conf_level) / 2.0, n - 1, err); // faulty library function
double t = MathQuantileT_TMP((1.0 - conf_level) / 2.0, n - 1, err); // custom implementation
if(err != 0)
{
Print("in MathQuantileT() error ", err);
return;
}
t *= c * MathSqrt(1.0 + 1.0 / n + MathPow((x_new - mean_x), 2) * (n - 1) / var_x);
//--- Print result if success
PrintFormat("Y between %.4f and %.4f with confidence level %.3f", a + b * x_new + t, a + b * x_new - t, conf_level);
}
//+------------------------------------------------------------------+
//| Residuals plot vs. price bar index |
//+------------------------------------------------------------------+
void t_residuals_plot(double& residuals[])
{
int n_residuals = ArraySize(residuals);
//--- Input data validation
if(n_residuals < 1)
{
Print("no data for t_residuals_plot()");
return;
}
//--- Array T[] as x on plot
double T[];
if(!MathSequenceByCount(1, n_residuals, n_residuals, T))
{
Print("MathSequenceByCount() error");
return;
}
//--- Plot construction
ChartSetInteger(0, CHART_SHOW, false);
CGraphic graphic;
graphic.Create(0, "G", 0, 0, 0, 750, 300);
graphic.CurveAdd(T, residuals, CURVE_LINES, "Residuals");
graphic.CurvePlotAll();
graphic.Update();
MessageBox("Press OK to close graph", "Close graph");
ChartSetInteger(0, CHART_SHOW, true);
graphic.Destroy();
}
//+------------------------------------------------------------------+
//| EPDF of residuals vs. normal density with residual SD |
//| nofx - number of points on plot |
//+------------------------------------------------------------------+
void epdf_vs_normalpdf(double& residuals[], int nofx = 30)
{
int n_residuals = ArraySize(residuals);
//--- Input data validation
if(n_residuals < 1)
{
Print("no data for epdf_vs_normalpdf()");
return;
}
//--- Array for x
double x[];
if(!MathSequenceByCount(MathMin(residuals), MathMax(residuals), nofx, x))
{
Print("sequence error");
return;
}
//--- Array for EPDF
double epdf[];
if(!MathProbabilityDensityEmpirical(residuals, nofx, x, epdf))
{
Print("MathProbabilityDensityEmpirical() error");
return;
}
//--- Array for theoretical PDF
double normal_pdf[];
if(!MathProbabilityDensityNormal(x, 0.0, MathStandardDeviation(residuals), normal_pdf))
{
Print("MathProbabilityDensityNormal() error");
return;
}
//--- Plot construction
ChartSetInteger(0, CHART_SHOW, false);
CGraphic graphic;
graphic.Create(0, "G", 0, 0, 0, 750, 400);
graphic.CurveAdd(x, epdf, CURVE_LINES, "Residuals");
graphic.CurveAdd(x, normal_pdf, CURVE_LINES, "Normal");
graphic.CurvePlotAll();
graphic.Update();
MessageBox("Press OK to close graph", "Close graph");
ChartSetInteger(0, CHART_SHOW, true);
graphic.Destroy();
}
//+------------------------------------------------------------------+
//| QQ-plot of residuals vs. normal distribution with residual SD |
//+------------------------------------------------------------------+
void qq_plot(double& residuals[])
{
int n_residuals = ArraySize(residuals);
//--- Input data validation
if(n_residuals < 1)
{
Print("no data for qq_plot()");
return;
}
//--- Copy input array and sort this copy
double residuals_copy[];
ArrayResize(residuals_copy, n_residuals);
ArrayCopy(residuals_copy, residuals);
ArraySort(residuals_copy);
//--- Probability sequence
double pseq[];
if(!MathSequenceByCount(0.5 / n_residuals, 1.0 - 0.5 / n_residuals, n_residuals, pseq))
{
Print("MathSequenceByCount() error");
return;
}
//--- Theoretical quantiles
double normal_qs[];
if(!MathQuantileNormal(pseq, 0.0, MathStandardDeviation(residuals_copy), normal_qs))
{
Print("MathQuantileNormal() error");
return;
}
//--- Plot construction
ChartSetInteger(0, CHART_SHOW, false);
CGraphic graphic;
graphic.Create(0, "G", 0, 0, 0, 750, 550);
graphic.CurveAdd(normal_qs, residuals_copy, CURVE_LINES, "QQ-plot");
graphic.CurveAdd(residuals_copy, residuals_copy, CURVE_LINES, "Y=X");
graphic.CurvePlotAll();
graphic.Update();
MessageBox("Press OK to close graph", "Close graph");
ChartSetInteger(0, CHART_SHOW, true);
graphic.Destroy();
}
//+------------------------------------------------------------------+
//| Correlogram of residuals |
//+------------------------------------------------------------------+
void correlogram(double& residuals[])
{
int n_residuals = ArraySize(residuals);
//--- Input data validation
if(n_residuals < 2)
{
Print("no data for correlogram()");
return;
}
//--- Number of lags
int ncg = (int)MathSqrt(n_residuals);
//--- Significance levels
double bound = 1.96 / MathSqrt(n_residuals);
double cg0 = 0.0, xb[2] = {0, ncg + 1}, upb[2] = {bound, bound}, lowb[2] = {-bound, -bound};
//--- Lags array creation
double lags[];
if(!MathSequenceByCount(1, ncg, ncg, lags))
{
Print("MathSequenceByCount() error");
return;
}
//--- ACF array calculation
double cg[];
ArrayResize(cg, ncg);
ArrayFill(cg, 0, ncg, 0.0);
for(int i = 0; i < n_residuals; ++i)
cg0 += residuals[i] * residuals[i];
for(int k = 1; k <= ncg; ++k)
{
for(int i = 0; i < n_residuals - k; ++i)
cg[k - 1] += residuals[i] * residuals[i + k];
cg[k - 1] /= cg0;
}
//--- Plot construction
ChartSetInteger(0, CHART_SHOW, false);
CGraphic graphic;
graphic.Create(0, "G", 0, 0, 0, 750, 300);
graphic.CurveAdd(lags, cg, CURVE_HISTOGRAM, "Correlogram");
graphic.CurveAdd(xb, upb, CURVE_LINES, "Upper bound");
graphic.CurveAdd(xb, lowb, CURVE_LINES, "Lower bound");
graphic.CurvePlotAll();
graphic.Update();
MessageBox("Press OK to close graph", "Close graph");
ChartSetInteger(0, CHART_SHOW, true);
graphic.Destroy();
}
//+------------------------------------------------------------------+
//| Scatter plot of (x, y) points with line y = a * x + b |
//+------------------------------------------------------------------+
void scatter_plot(double& x[], double& y[], bool add_line = false, double a = 0.0, double b = 0.0)
{
int ny = ArraySize(y);
//--- Input data validation
if(ny < 1)
{
Print("no data for regression1()");
return;
}
if(ny != ArraySize(x))
{
Print("arrays y[] and x[] have different lengths in scatter_plot()");
return;
}
//--- Regression line points calculation
double xreg[2] = {MathMin(x), MathMax(x)};
double yreg[2] = {a + b*xreg[0], a + b*xreg[1]};
//--- Plot construction
ChartSetInteger(0, CHART_SHOW, false);
CGraphic graphic;
graphic.Create(0, "G", 0, 0, 0, 750, 550);
graphic.CurveAdd(x, y, CURVE_POINTS, "scatter plot");
if(add_line)
graphic.CurveAdd(xreg, yreg, CURVE_LINES, "y = a * x + b");
graphic.CurvePlotAll();
graphic.Update();
MessageBox("Press OK to close graph", "Close graph");
ChartSetInteger(0, CHART_SHOW, true);
graphic.Destroy();
}
//+------------------------------------------------------------------+
//| Ljung-Box test for autocorrelation in residuals |
//| m - lags number, pq - sum of p and q in ARMA(p,q) model |
//+------------------------------------------------------------------+
void Ljung_Box_test(double& residuals[], int m = 10, int pq = 0)
{
Print("Ljung-Box test result:");
//--- Input data validation
if(pq < 0)
{
Print("pq must be nonnegative");
return;
}
if(m <= pq)
{
Print("m must be greater pq");
return;
}
int n_residuals = ArraySize(residuals);
if(n_residuals < m + 1)
{
Print("not enough data");
return;
}
//--- Initialization
double Q = 0.0, cov0 = 0.0, cov, ro;
//--- Test statistic calculation
for(int i = 0; i < n_residuals; ++i)
cov0 += residuals[i] * residuals[i];
for(int k = 1; k <= m; ++k)
{
cov = 0.0;
for(int i = 0; i < n_residuals - k; ++i)
cov += residuals[i] * residuals[i + k];
ro = cov / cov0;
Q += ro * ro / (n_residuals - k);
}
Q *= n_residuals * (n_residuals + 2);
//--- p-value computation
int err = 0;
double p_value = MathCumulativeDistributionChiSquare(Q, m - pq, false, false, err);
if(err != 0)
{
Print("MathCumulativeDistributionChiSquare() error: ", err);
return;
}
//--- Print result
PrintFormat("Q = %.3f, p-value = %.3f", Q, p_value);
}
//+------------------------------------------------------------------+
//| Jarque-Bera test for normality of distribution |
//+------------------------------------------------------------------+
void Jarque_Bera_test(double& residuals[])
{
Print("Jarque-Bera test result:");
//--- Input data validation
int n_residuals = ArraySize(residuals);
if(n_residuals < 1)
{
Print("not enough data");
return;
}
//--- Skewness calculation
double S = MathSkewness(residuals);
//--- Kurtosis calculation
double K = MathKurtosis(residuals);
//--- Test statistic calculation
double JB = n_residuals * (S * S + K * K / 4.0) / 6.0;
//--- p-value computation
double p_value = MathExp(-JB / 2.0);
//--- Print result
PrintFormat("JB = %.3f, p-value = %.3f", JB, p_value);
}
//+------------------------------------------------------------------+
//| Pettitt test for structural break (abrupt change in mean) |
//+------------------------------------------------------------------+
void Pettitt_test(double& residuals[])
{
Print("Pettitt test result:");
//--- Input data validation
int n_residuals = ArraySize(residuals);
if(n_residuals < 2)
{
Print("not enough data");
return;
}
//--- Initialization
double K = 0.0, U, W = 0.0;
int t_change = 0;
//--- Sample ranking computation
double rank[];
MathRank(residuals, rank);
//--- Test statistic calculation
for(int t = 1; t < n_residuals; ++t)
{
W += rank[t - 1];
U = MathAbs(2.0 * W - t * (n_residuals + 1.0));
if(U > K)
{
K = U;
t_change = t;
}
}
//--- p-value computation
double p_value = MathMin(2.0 * MathExp(-6.0 * K * K / (n_residuals * n_residuals * (n_residuals + 1.0))), 1.0);
//--- Print result
PrintFormat("K = %.3f, p-value = %.3f, change-point = %d", K, p_value, t_change);
}
//+------------------------------------------------------------------+
//| Student's t-distribution quantile |
//| (custom implementation replacing faulty library function) |
//| JUST A TEMPORARY REPLACEMENT! (Based on library variant) |
//+------------------------------------------------------------------+
double MathQuantileT_TMP(const double prob, const double nu, int &error_code)
{
if(!MathIsValidNumber(prob) || !MathIsValidNumber(nu))
{
error_code = ERR_ARGUMENTS_NAN;
return QNaN;
}
if(nu <= 0.0 || nu != MathRound(nu))
{
error_code = ERR_ARGUMENTS_INVALID;
return QNaN;
}
if(prob < 0.0 || prob > 1.0)
{
error_code = ERR_ARGUMENTS_INVALID;
return QNaN;
}
error_code = ERR_RESULT_INFINITE;
if(prob == 0.0)
return QNEGINF;
if(prob == 1.0)
return QPOSINF;
error_code = ERR_OK;
if(prob == 0.5)
return 0.0;
if(nu == 1.0)
return MathTan(M_PI * (prob - 0.5));
int err_code = 0;
double x_low = -1.0;
double x_high = 1.0;
while(MathCumulativeDistributionT(x_low, nu, err_code) > prob && x_low > -1e6)
x_low *= 2.0;
while(MathCumulativeDistributionT(x_high, nu, err_code) < prob && x_high < 1e6)
x_high *= 2.0;
double x = 0.0;
int max_iterations = 100;
double precision = 1e-12; // Высокая точность
for(int i = 0; i < max_iterations; i++)
{
x = (x_low + x_high) / 2.0;
double cdf = MathCumulativeDistributionT(x, nu, err_code);
if(MathAbs(x_high - x_low) < precision)
break;
if(cdf < prob)
x_low = x;
else
x_high = x;
if(i == max_iterations - 1)
error_code = ERR_NON_CONVERGENCE;
}
return x;
}
//+------------------------------------------------------------------+
//| Helper function to export array to text file |
//| (for data analysis in other programs) |
//+------------------------------------------------------------------+
void array2csv(double& a[], string fldr, string fnm)
{
//--- Initialization and input data validation
int ftxt = FileOpen(fldr + "\\" + fnm, FILE_WRITE | FILE_TXT | FILE_ANSI | FILE_COMMON);
if(ftxt == INVALID_HANDLE)
{
Print("FileOpen() error");
return;
}
int na = ArraySize(a);
if(na < 1)
{
Print("No data for array2csv()");
return;
}
//--- Writing array to file
FileWriteString(ftxt, DoubleToString(a[0]));
for(int i = 1; i < na; ++i)
FileWriteString(ftxt, "\n" + DoubleToString(a[i]));
FileClose(ftxt);
}
//+------------------------------------------------------------------+