//+------------------------------------------------------------------+ //| EconometricsM.mqh | //| Copyright 2000-2026, MetaQuotes Ltd. | //| www.mql5.com | //+------------------------------------------------------------------+ //--- #include "EconometricsA.mqh" //--- //+------------------------------------------------------------------+ //| struct for coefficients testing | //+------------------------------------------------------------------+ struct CoefficientStats { double estimate; // point estimate double std_error; // standard error double t_stat; // t-statistic double p_value; // one-sided p-value double conf_low; // low bound of confidence interval double conf_high; // high bound of confidence interval }; //+------------------------------------------------------------------+ //| struct for prognosed data | //+------------------------------------------------------------------+ struct SPrognose { vector xnew; // regressors for prognose (xnew[0]==1!) double point_progn; // point prognose (mean) double conf_low; // low bound of confidence interval double conf_high; // high bound of confidence interval double progn_low; // low bound of prognose interval double progn_high; // high bound of prognose interval void print() // print struct { PrintFormat("point prognose: %.3f",point_progn); PrintFormat("confidence interval, low bound: %.3f, high bound: %.3f",conf_low,conf_high); PrintFormat("prognose interval, low bound: %.3f, high bound: %.3f",progn_low,progn_high); }; }; //+------------------------------------------------------------------+ //| struct for prognosed data for TSLS | //+------------------------------------------------------------------+ struct SPrognose_TSLS { vector xnew_exo; // regressors for prognose (xnew_exo[0]==1!) vector xnew_endo; // regressors for prognose vector znew_iv; // instruments for prognose double point_progn; // point prognose (mean) double conf_low; // low bound of confidence interval double conf_high; // high bound of confidence interval double progn_low; // low bound of prognose interval double progn_high; // high bound of prognose interval void print() // print struct { PrintFormat("point prognose: %.3f",point_progn); PrintFormat("confidence interval, low bound: %.3f, high bound: %.3f",conf_low,conf_high); PrintFormat("prognose interval, low bound: %.3f, high bound: %.3f",progn_low,progn_high); }; }; //+------------------------------------------------------------------+ //| Computation of parameters and residuals for multiple regression | //| Inputs: y - target variable, X - regressors | //| Outputs: b - parameters, e - residuals, c - residuals' SD | //+------------------------------------------------------------------+ void regression(vector& y, matrix& X, vector& b, vector& e, double& c) { //--- Input data validation ulong k=X.Cols(); if(k<1) {Print("error: empty X"); return;} ulong n=X.Rows(); if(n 0) stats[i].t_stat = stats[i].estimate / stats[i].std_error; else stats[i].t_stat = 0.0; // One-sided p-value computation using Student's t-distribution stats[i].p_value = MathCumulativeDistributionT(MathAbs(stats[i].t_stat),n-k,false,false,ner); if(ner!=0) {Print("regression_Newey_West() error: MathCumulativeDistributionT() error: ",ner); return;} // 95% Confidence interval for the coefficient (critical t-value approx 1.96 for large n) // For precision, we fetch the exact critical value from Student's inverse CDF t_crit=MathQuantileT(0.025,n-k,false,false,ner); if(ner!=0) {Print("regression_Newey_West() error: MathQuantileT() error: ",ner); return;} stats[i].conf_low = stats[i].estimate - t_crit * stats[i].std_error; stats[i].conf_high = stats[i].estimate + t_crit * stats[i].std_error; } //--- Forecast and interval boundaries computation // 1. Point forecast (mean prediction) // Vector dot product: prog.xnew (k x 1) and b (k x 1) prog.point_progn = prog.xnew.Dot(b); // 2. Variance of the forecast mean (for Confidence Interval) // Quadratic form: x_new^T * Vb * x_new double var_ci = prog.xnew.Dot(Vb.MatMul(prog.xnew)); double std_ci = MathSqrt(var_ci); // 3. Variance of the individual prediction (for Prediction Interval) double var_pi = var_ci + se2; double std_pi = MathSqrt(var_pi); // Get critical t-value for 95% threshold (pre-calculated t_crit from coefficients statistics computing can be reused) t_crit=MathQuantileT(0.025,n-k,false,false,ner); if(ner!=0) {Print("regression_Newey_West() error: MathQuantileT() error: ",ner); return;} // 4. Compute final boundaries for both intervals prog.conf_low = prog.point_progn - t_crit * std_ci; prog.conf_high = prog.point_progn + t_crit * std_ci; prog.progn_low = prog.point_progn - t_crit * std_pi; prog.progn_high = prog.point_progn + t_crit * std_pi; } //+-------------------------------------------------------------------+ //| Conditional Least Squares (CLS) for ARMA(p,q) model | //| | //| ARGUMENTS: | //| [in] y - Dependent variable vector (size n x 1) | //| [in] p - Order of the AR component | //| [in] q - Order of the MA component | //| [in] y_ - Pre-sample history of dependent variable (p x 1) | //| [in] e_ - Pre-sample history of residuals/errors (q x 1) | //| [out] e - residuals/errors (n x 1) | //| [out] stats - Dynamic array to store coefficients statistics | //| [out] prog - Structure containing next step forecast results | //+-------------------------------------------------------------------+ void regression_CLS(const vector& y, const ulong p, const ulong q, const vector& y_, const vector& e_, vector& e, CoefficientStats& stats[], SPrognose& prog, const ulong max_iter=1000, const double min_change=1e-5) { //--- Input data validation ulong n=y.Size(),k=p+q+1; if(max_iter==0) {Print("regression_CLS() error: max_iter==0"); return;} if(min_change<=0.0) {Print("regression_CLS() error: min_change<=0.0"); return;} if(n<=k) {Print("regression_CLS() error: y too short"); return;} if(y_.Size()!=p) {Print("regression_CLS() error: wrong y_ size"); return;} if(e_.Size()!=q) {Print("regression_CLS() error: wrong e_ size"); return;} //-- AR-part coefficients initialization with OLS matrix X0=matrix::Ones(n,p+1); for(ulong i=0;i=j)?y[i-j]:y_[p+i-j]; if(X0.Rank()=i)?b[p+i]*J[t-i][0]:0.0; for(ulong j=1;j<=p;++j) { J[t][j]=(t>=j)?-y[t-j]:-y_[p+t-j]; for(ulong i=1;i<=q;++i) J[t][j]-=(t>=i)?b[p+i]*J[t-i][j]:0.0; } for(ulong j=1;j<=q;++j) { J[t][p+j]=(t>=j)?-e[t-j]:-e_[q+t-j]; for(ulong i=1;i<=q;++i) J[t][p+j]-=(t>=i)?b[p+i]*J[t-i][p+j]:0.0; } } //--- delta_b = -(J^T * J)^(-1) * J^T * e JT_J = J.Transpose().MatMul(J); // Add Lambda to the main diagonal (Levenberg-Marquardt modification) for(ulong i = 0; i < k; ++i) JT_J[i][i] += lambda; // Compute the damped step delta_b = -1.0 * JT_J.Inv().MatMul(J.Transpose()).MatMul(e); //--- 1. Apply the step temporarily b += delta_b; //--- 2. Recalculate residuals vector e for the current parameters for(ulong t = 0; t < n; ++t) { e[t] = y[t] - b[0]; for(ulong j = 1; j <= p; ++j) e[t] -= (t >= j) ? b[j] * y[t - j] : b[j] * y_[p + t - j]; for(ulong j = 1; j <= q; ++j) e[t] -= (t >= j) ? b[p + j] * e[t - j] : b[p + j] * e_[q + t - j]; } rss_new = e.Dot(e); //--- 3. Check if the step improved the model if(rss_new < rss_old) { // Successful iteration: accept updates lambda /= v; // Decrease damping rss_old = rss_new; b_old = b; // Update backups e_old = e; // Convergence check based on parameter changes change = delta_b.Norm(VECTOR_NORM_P, 2.0); if(change < min_change) { PrintFormat("L-M converged successfully after %d iterations.", iter + 1); break; } } else { // Unsuccessful iteration: rollback updates and increase damping lambda *= v; // Increase damping to behave more like gradient descent b = b_old; // Rollback coefficients e = e_old; // Rollback residuals to previous valid state // Guard against infinite loops with huge lambda if(lambda > 1e10) { Print("L-M stopped: lambda became too large (optimization stuck)."); return; } } } if(iter>=max_iter) {Print("L-M not converged"); return;} //--- Compute unbiased residual variance (sigma^2) double se2 = e.Dot(e) / double(n - k); //--- Coefficients covariance matrix for non-linear OLS //--- Calculated as sigma^2 * (J^T * J)^(-1) using the final iteration's Jacobian matrix matrix Vb = se2 * J.Transpose().MatMul(J).Inv(); //--- Extract standard errors and compute coefficients statistics ArrayResize(stats, (int)k); int ner = 0; // 95% Confidence interval for the coefficient double t_crit = MathQuantileT(0.025, double(n - k), false, false, ner); if(ner != 0) { Print("regression_CLS() error: MathQuantileT() failed with code ", ner); return; } for(ulong i = 0; i < k; ++i) { stats[i].estimate = b[i]; // Standard error is the square root of the diagonal element of Vb stats[i].std_error = MathSqrt(Vb[i][i]); // T-statistic computation (testing H0: b[i] == 0) if(stats[i].std_error > 0) stats[i].t_stat = stats[i].estimate / stats[i].std_error; else stats[i].t_stat = 0.0; // One-sided p-value computation using Student's t-distribution stats[i].p_value = MathCumulativeDistributionT(MathAbs(stats[i].t_stat), double(n - k), false, false, ner); if(ner != 0) { Print("regression_CLS() error: MathCumulativeDistributionT() failed with code ", ner); return; } stats[i].conf_low = stats[i].estimate - t_crit * stats[i].std_error; stats[i].conf_high = stats[i].estimate + t_crit * stats[i].std_error; } //--- Forecast and prediction intervals computation (One-Step-Ahead) // Create a full regressor vector for the forecast bar (size: p + q + 1) prog.xnew.Resize(k); // 1. Set Constant components prog.xnew[0] = 1.0; // 2. Fill AR lags from the very end of the sample vector y for(ulong j = 1; j <= p; ++j) prog.xnew[j] = (n >= j) ? y[n - j] : y_[p + (n - j)]; // 3. Fill MA lags from the very end of the calculated residuals vector e for(ulong j = 1; j <= q; ++j) prog.xnew[p + j] = (n >= j) ? e[n - j] : e_[q + (n - j)]; // 4. Compute Point Forecast (Mean prediction) prog.point_progn = prog.xnew.Dot(b); // 5. Variance of the forecast mean (Confidence Interval) // For non-linear OLS, the gradient vector at the forecast point is equal to prog.xnew double var_ci = prog.xnew.Dot(Vb.MatMul(prog.xnew)); double std_ci = MathSqrt(var_ci); // 6. Variance of the individual future value (Prediction Interval) double var_pi = var_ci + se2; double std_pi = MathSqrt(var_pi); // Pre-calculated negative t_crit from step 3 is reused here (with adjusted signs) prog.conf_low = prog.point_progn - t_crit * std_ci; prog.conf_high = prog.point_progn + t_crit * std_ci; prog.progn_low = prog.point_progn - t_crit * std_pi; prog.progn_high = prog.point_progn + t_crit * std_pi; } //+------------------------------------------------------------------+ //| Durbin-Wu-Hausman (DWH) Endogeneity Test (Regression-Based) | //| H0: Regressors are exogenous (OLS is consistent and efficient) | //| H1: Regressors are endogenous (OLS is biased, 2SLS is required) | //| Returns: p-value of the test (p-value < 0.05 means endogeneity) | //+------------------------------------------------------------------+ double DWH_test(const vector &y, const matrix &X_exo, const matrix &X_endo, const matrix &Z_iv) { ulong n = y.Size(); ulong p = X_exo.Cols(); ulong r = X_endo.Cols(); ulong m = Z_iv.Cols(); //--- 1. Input data validation if(n < 1 || X_exo.Rows() != n || X_endo.Rows() != n || Z_iv.Rows() != n) { Print("DWH Test error: matrix row dimensions mismatch or empty data."); return 1.0; } if(m < r) { Print("DWH Test error: too few instruments (m < r). Test cannot be performed."); return 1.0; } //--- 2. Construct the full instrument matrix Z = [X_exo | Z_iv] matrix Z = X_exo.Concat(Z_iv, 1); if(Z.Rank() < Z.Cols()) { Print("DWH Test error: multicollinearity in the instrument matrix Z."); return 1.0; } //--- 3. FIRST STAGE: Regress each endogenous variable on all instruments to get residuals matrix V_hat; V_hat.Resize(n, r); matrix Pl = Z.Transpose().MatMul(Z).Inv().MatMul(Z.Transpose()); vector bl, el, ytmp; for(ulong i = 0; i < r; ++i) { ytmp = X_endo.Col(i); bl = Pl.MatMul(ytmp); //--- Coefficients of the first stage el = ytmp - Z.MatMul(bl); //--- Residuals containing the "toxic" endogenous part V_hat.Col(el, i); //--- Save residuals as a column } //--- 4. SECOND STAGE: Fit the Unrestricted (Augmented) Model via OLS: Y ~ X_exo + X_endo + V_hat matrix X_orig = X_exo.Concat(X_endo, 1); matrix X_augmented = X_orig.Concat(V_hat, 1); ulong k_aug = X_augmented.Cols(); // Total variables in the augmented model (p + r + r) //--- Validate sample size after augmentation if(n <= k_aug) { Print("DWH Test error: too few observations for the augmented model."); return 1.0; } matrix X_aug_T = X_augmented.Transpose(); matrix XX_aug_inv = X_aug_T.MatMul(X_augmented).Inv(); vector b_aug = XX_aug_inv.MatMul(X_aug_T).MatMul(y); //--- Calculate Sum of Squares for the Unrestricted Model (RSS_unrestricted) vector e_aug = y - X_augmented.MatMul(b_aug); double rss_unrestricted = e_aug.MatMul(e_aug); //--- 5. Fit the Restricted Model (Standard OLS): Y ~ X_exo + X_endo matrix X_orig_T = X_orig.Transpose(); vector b_ols = X_orig_T.MatMul(X_orig).Inv().MatMul(X_orig_T).MatMul(y); //--- Calculate Sum of Squares for the Restricted Model (RSS_restricted) vector e_ols = y - X_orig.MatMul(b_ols); double rss_restricted = e_ols.MatMul(e_ols); //--- 6. Compute the F-statistic for the joint significance of V_hat coefficients //--- Degrees of freedom: numerator = r (number of restrictions), denominator = n - k_aug double h_numerator = (rss_restricted - rss_unrestricted) / (double)r; double h_denominator = rss_unrestricted / (double)(n - k_aug); if(h_denominator <= 0.0) { Print("DWH Test warning: perfect fit or zero residual variance in the augmented model."); return 1.0; } double hausman_F = h_numerator / h_denominator; int h_err = 0; //--- Compute the right-tailed p-value from the Fisher F-distribution double f_cdf = MathCumulativeDistributionF(hausman_F, (double)r, (double)(n - k_aug), h_err); double p_value = 1.0 - f_cdf; //--- 7. Print diagnostic results to the terminal log PrintFormat("DWH Test | F-stat: %.4f | p-value: %.6f", hausman_F, p_value); if(p_value < 0.05) Print("DWH Test Result: Reject H0. Endogeneity is significant. Use 2SLS."); else Print("DWH Test Result: Fail to reject H0. No significant endogeneity. OLS is preferred."); return p_value; } //+------------------------------------------------------------------+ //| Two stage least square algorithm (TSLS) | //+------------------------------------------------------------------+ void TSLS(vector& y, matrix& X_exo, matrix& X_endo, matrix& Z_iv, SPrognose_TSLS& prog, CoefficientStats& stats[], double confidence_level = 0.95) { //--- Extract data dimensions: n - observations, p - exogenous, r - endogenous, m - clean instruments ulong n=y.Size(), q=X_exo.Cols(), r=X_endo.Cols(), m=Z_iv.Cols(); //--- Input data validation if(n<1) {Print("TSLS error: empty y"); return;} if(X_exo.Rows()!=n) {Print("TSLS error: wrong X_exo rows number"); return;} if(X_endo.Rows()!=n) {Print("TSLS error: wrong X_endo rows number"); return;} if(Z_iv.Rows()!=n) {Print("TSLS error: wrong Z_iv rows number"); return;} if(n <= (q + r) || n <= (q + m)) {Print("TSLS error: too few observations (n) for the number of variables"); return;} if(m 0) { stats[i].t_stat = stats[i].estimate / stats[i].std_error; //--- MathCumulativeDistributionT returns the left-tail probability P(T <= t) double p_cumulative = MathCumulativeDistributionT(MathAbs(stats[i].t_stat), (double)df, err_code); stats[i].p_value = 1.0 - p_cumulative; } else { stats[i].t_stat = 0.0; stats[i].p_value = 1.0; } //--- Calculate the confidence interval bounds for the coefficient stats[i].conf_low = stats[i].estimate - (t_crit * stats[i].std_error); stats[i].conf_high = stats[i].estimate + (t_crit * stats[i].std_error); } //--- Out-of-sample forecasting stage using the instrumented approach //--- 1. Construct the complete out-of-sample instrument vector Z0 = [xnew_exo | znew_iv] vector z0 = prog.xnew_exo.Concat(prog.znew_iv); //--- 2. Clean the future: project the future endogenous variables using first-stage coefficients prog.xnew_endo = z0.MatMul(B_first); //--- 3. Form the final cleaned out-of-sample regressor vector for the second-stage equation vector x0_hat = prog.xnew_exo.Concat(prog.xnew_endo); //--- 4. Calculate the symmetric point forecast (mean prediction) using vector dot product prog.point_progn = x0_hat.MatMul(b_vec); //--- 5. Calculate the variance scale factor (g-factor) for the out-of-sample bar matrix X0_mat(1, x0_hat.Size()); X0_mat.Row(x0_hat, 0); matrix X0_T = X0_mat.Transpose(); matrix shift_mat = X0_mat.MatMul(XX_inv).MatMul(X0_T); //--- Extract the scalar value from the 1x1 matrix using explicit indexing double g_factor = shift_mat[0][0]; //--- 6. Calculate the standard errors for both the confidence interval and the prediction interval double se_conf = MathSqrt(sigma_sq * g_factor); double se_pred = MathSqrt(sigma_sq * (1.0 + g_factor)); //--- 7. Derive half-widths for the interval bounds using the critical t-value double conf_width = t_crit * se_conf; double pred_width = t_crit * se_pred; //--- 8. Write upper and lower boundaries into the SPrognose_TSLS structure fields prog.conf_low = prog.point_progn - conf_width; prog.conf_high = prog.point_progn + conf_width; prog.progn_low = prog.point_progn - pred_width; prog.progn_high = prog.point_progn + pred_width; } //+------------------------------------------------------------------+ //| Computation of the correlation matrix (assumed X1==const) | //| Inputs: y - target variable, X - regressors | //| Output: CM - correlation matrix | //+------------------------------------------------------------------+ void corr_matrix(vector& y, matrix& X, matrix& CM) { //--- Input data validation ulong k=X.Cols(); if(k<1) {Print("corr_matrix() error: empty X"); return;} ulong n=X.Rows(); if(y.Size()!=n) {Print("corr_matrix() error: wrong size of y"); return;} matrix yX=X; yX.Col(y,0); CM=yX.CorrCoef(false); } //+------------------------------------------------------------------+ //| Residuals plot vs. price bar index | //+------------------------------------------------------------------+ void t_residuals_plot(vector& residuals) { double e[]; vector2array(residuals,e); t_residuals_plot(e); } //+------------------------------------------------------------------+ //| EPDF of residuals vs. normal density with residual SD | //| nofx - number of points on plot | //+------------------------------------------------------------------+ void epdf_vs_normalpdf(vector& residuals, int nofx = 30) { double e[]; vector2array(residuals,e); epdf_vs_normalpdf(e,nofx); } //+------------------------------------------------------------------+ //| QQ-plot of residuals vs. normal distribution with residual SD | //+------------------------------------------------------------------+ void qq_plot(vector& residuals) { double e[]; vector2array(residuals,e); qq_plot(e); } //+------------------------------------------------------------------+ //| Correlogram of residuals | //+------------------------------------------------------------------+ void correlogram(vector& residuals) { double e[]; vector2array(residuals,e); correlogram(e); } //+------------------------------------------------------------------+ //| Scatter plot of (x, y) points with line y = a * x + b | //+------------------------------------------------------------------+ void scatter_plot(vector& x, vector& y, bool add_line = false, double a = 0.0, double b = 0.0) { double xa[],ya[]; vector2array(x,xa); vector2array(y,ya); scatter_plot(xa,ya,add_line,a,b); } //+------------------------------------------------------------------+ //| Scatter plot of (Xi, y) and fitted regression line | //| Inputs: X - regressors matrix, y - target variable, | //| i - regressor's index (X column index) | //+------------------------------------------------------------------+ void scatter_plot_Xi_y(matrix& X, ulong i, vector& y) { if(i>=X.Cols()) { Print("input regressor's index is out of range"); return; } double a,b,c; vector e; regression1(y,X.Col(i),a,b,c,e); scatter_plot(X.Col(i),y,true,a,b); } //+------------------------------------------------------------------+ //| Scatter plot of (Xi, Xj) and fitted regression line | //| Inputs: X - regressors matrix, | //| i, j - regressors' indexes (X columnes indexes) | //+------------------------------------------------------------------+ void scatter_plot_Xi_Xj(matrix& X, ulong i, ulong j) { if(i>=X.Cols()||j>=X.Cols()) { Print("input regressor's index is out of range"); return; } double a,b,c; vector e; regression1(X.Col(j),X.Col(i),a,b,c,e); scatter_plot(X.Col(i),X.Col(j),true,a,b); } //+------------------------------------------------------------------+ //| Partial regression plot for Xi | //| Inputs: X - regressors matrix, y - target variable, | //| i - regressor's index (X column index) | //+------------------------------------------------------------------+ void partial_regression_plot(matrix& X, ulong i, vector& y) { //--- Input data validation ulong k=X.Cols(); if(k<2) { Print("too few regressors"); return; } if(i>=k) { Print("input regressor's index is out of range"); return; } double a,b,c; vector x=X.Col(i),ey,ex,e1,bv; //--- matrix without Xi matrix X_1=X; if(i!=k-1) X_1.Col(X_1.Col(k-1),i); X_1.Resize(X_1.Rows(),k-1); regression(y,X_1,bv,ey,c); regression(x,X_1,bv,ex,c); regression1(ey,ex,a,b,c,e1); scatter_plot(ex,ey,true,a,b); } //+------------------------------------------------------------------+ //| Scatter plot of (Yfit, Yreal) and line y = x | //| Inputs: X - regressors matrix, b - parameters vector | //| y - real values of dependent variable | //+------------------------------------------------------------------+ void scatter_plot_Yfit_Yreal(matrix& X, vector& b, vector& y) { //--- Input data validation ulong k=X.Cols(); if(k<1) {Print("error: empty X"); return;} if(b.Size()!=k) {Print("error: wrong size of b"); return;} ulong n=X.Rows(); if(y.Size()!=n) {Print("error: wrong size of y"); return;} scatter_plot(X.MatMul(b),y,true,0.0,1.0); } //+------------------------------------------------------------------+ //| Scatter plot of (Yfit, residuals) | //| Inputs: X - regressors matrix, b - parameters vector | //| e - residuals | //+------------------------------------------------------------------+ void scatter_plot_Yfit_residuals(matrix& X, vector& b, vector& e) { //--- Input data validation ulong k=X.Cols(); if(k<1) {Print("error: empty X"); return;} if(b.Size()!=k) {Print("error: wrong size of b"); return;} ulong n=X.Rows(); if(e.Size()!=n) {Print("error: wrong size of e"); return;} scatter_plot(X.MatMul(b),e); } //+------------------------------------------------------------------+ //| R^2, coefficient of determination | //| Inputs: y - real values of dependent variable, e - residuals | //+------------------------------------------------------------------+ double R2(vector& y, vector& e) { //--- Input data validation ulong n=y.Size(); if(n<2) {Print("R2() error: y too short"); return 0.0;} if(e.Size()!=n) {Print("R2() error: different sizes of y and e"); return 0.0;} double Sy=y.Var(0); if(Sy<=DBL_MIN) {Print("R2() error: y = const"); return 0.0;} return 1.0-e.Var(0)/Sy; } //+------------------------------------------------------------------+ //| R^2_adj, adjusted coefficient of determination | //| Inputs: y - real values of dependent variable, e - residuals | //| k - number of regressors | //+------------------------------------------------------------------+ double R2_adj(vector& y, vector& e, ulong k) { //--- Input data validation ulong n=y.Size(); if(n<2) {Print("R2_adj() error: y too short"); return 0.0;} if(n<=k) {Print("R2_adj() error: n <= k"); return 0.0;} if(e.Size()!=n) {Print("R2_adj() error: different sizes of y and e"); return 0.0;} double Sy=y.Var(1); if(Sy<=DBL_MIN) {Print("R2_adj() error: y = const"); return 0.0;} return 1.0-e.Var((int)k)/Sy; } //+------------------------------------------------------------------+ //| VIF, Variance Inflation Factor for X2, ..., Xk (assumed X1=const)| //| Input: X - regressors | //| Output: vif - (k-1)-size VIF vector | //+------------------------------------------------------------------+ void VIF(matrix& X, vector& vif) { //--- Input data validation ulong k=X.Cols(), n=X.Rows(); if(k<2) {Print("VIF() error: not enough regressors"); return;} if(n=k) {Print("error: out of regressors range"); return;} for(int i=0;i=0;--i) { if(irs[i]to) return; if(to>=A.Cols()) return; if(A.Rows()<2) return; double m, s; for(ulong i=from;i<=to; ++i) {m=A.Col(i).Mean(); s=A.Col(i).Std(); A.Col((A.Col(i)-m)/s,i);} } //+------------------------------------------------------------------+ //| Helper functions to export vector, matrix and vector+matrix | //| to text file (for data analysis in other programs) | //+------------------------------------------------------------------+ string vector2string(vector& v) { ulong k=v.Size(); if(k<1) return ""; string res=DoubleToString(v[0]); for(ulong i=1;i