411 lines
14 KiB
MQL5
411 lines
14 KiB
MQL5
//+------------------------------------------------------------------+
|
|
//| recursions.mqh |
|
|
//| Copyright 2025, MetaQuotes Ltd. |
|
|
//| https://www.mql5.com |
|
|
//+------------------------------------------------------------------+
|
|
#property copyright "Copyright 2025, MetaQuotes Ltd."
|
|
#property link "https://www.mql5.com"
|
|
#include"base.mqh"
|
|
//---
|
|
#define SQRT2_OV_PI 0.79788456080286541
|
|
const double LNSIGMA_MAX = log(DBL_MAX) - 0.1;
|
|
//+------------------------------------------------------------------+
|
|
//| helper function |
|
|
//+------------------------------------------------------------------+
|
|
double bounds_check(double _sigma2, vector &varbounds)
|
|
{
|
|
if(_sigma2<varbounds[0])
|
|
_sigma2 = varbounds[0];
|
|
else
|
|
if(_sigma2>varbounds[1])
|
|
{
|
|
if(MathClassify(_sigma2) == FP_NORMAL)
|
|
_sigma2 = varbounds[1] +log(_sigma2/varbounds[1]);
|
|
else
|
|
_sigma2 = varbounds[1] + 1000;
|
|
}
|
|
return _sigma2;
|
|
}
|
|
//+------------------------------------------------------------------+
|
|
//|arch recursion utility |
|
|
//+------------------------------------------------------------------+
|
|
vector arch_recursion(vector& parameters, vector& resids, vector& sigma2, ulong p, ulong nobs, double backcast, matrix &var_bounds)
|
|
{
|
|
for(ulong t =0; t<nobs; ++t)
|
|
{
|
|
sigma2[t] = parameters[0];
|
|
for(ulong i = 0; i<p; ++i)
|
|
{
|
|
if((t-i-1)<0)
|
|
sigma2[t] += parameters[i + 1] * backcast;
|
|
else
|
|
sigma2[t] += parameters[i + 1] * pow(resids[t - i - 1], 2);
|
|
}
|
|
sigma2[t] = bounds_check(sigma2[t], var_bounds.Row(t));
|
|
}
|
|
return sigma2;
|
|
}
|
|
//+------------------------------------------------------------------+
|
|
//|Compute variance recursion for EGARCH models |
|
|
//+------------------------------------------------------------------+
|
|
vector egarch_recursion(vector& parameters, vector& resids, vector& sigma2, ulong p, ulong o, ulong q,ulong nobs, double backcast,matrix& var_bounds, vector& lnsigma2, vector& std_resids, vector& abs_std_resids)
|
|
{
|
|
for(ulong t = 0; t<nobs; ++t)
|
|
{
|
|
ulong loc = 0;
|
|
lnsigma2[t] = parameters[loc];
|
|
loc += 1;
|
|
for(ulong j = 0; j<p; ++j)
|
|
{
|
|
if((t - 1 - j) >= 0)
|
|
lnsigma2[t] += parameters[loc] * (abs_std_resids[t - 1 - j] - SQRT2_OV_PI);
|
|
loc += 1;
|
|
}
|
|
for(ulong j = 0; j<o; ++j)
|
|
{
|
|
if((t - 1 - j) >= 0)
|
|
lnsigma2[t] += parameters[loc] * std_resids[t - 1 - j];
|
|
loc += 1;
|
|
}
|
|
for(ulong j = 0; j<q; ++j)
|
|
{
|
|
if((t - 1 - j) < 0)
|
|
lnsigma2[t] += parameters[loc] * backcast;
|
|
else
|
|
lnsigma2[t] += parameters[loc] * lnsigma2[t - 1 - j];
|
|
loc += 1;
|
|
}
|
|
if(lnsigma2[t] > LNSIGMA_MAX)
|
|
lnsigma2[t] = LNSIGMA_MAX;
|
|
sigma2[t] = exp(lnsigma2[t]);
|
|
if(sigma2[t] < var_bounds[t, 0])
|
|
{
|
|
sigma2[t] = var_bounds[t, 0];
|
|
lnsigma2[t] = log(sigma2[t]);
|
|
}
|
|
else
|
|
if(sigma2[t] > var_bounds[t, 1])
|
|
{
|
|
sigma2[t] = var_bounds[t, 1] + log(sigma2[t]) - log(var_bounds[t, 1]);
|
|
lnsigma2[t] = log(sigma2[t]);
|
|
}
|
|
std_resids[t] = resids[t] / sqrt(sigma2[t]);
|
|
abs_std_resids[t] = MathAbs(std_resids[t]);
|
|
}
|
|
return sigma2;
|
|
}
|
|
|
|
//+------------------------------------------------------------------+
|
|
//|Compute variance recursion for GARCH and related models |
|
|
//+------------------------------------------------------------------+
|
|
vector garch_recursion(vector ¶meters, vector& fresids,vector& sresids, vector& sigma2, ulong p, ulong o, ulong q, ulong nobs, double backcast,matrix &var_bounds)
|
|
{
|
|
for(long t = 0; t<long(nobs); ++t)
|
|
{
|
|
long loc = 0;
|
|
sigma2[t] = parameters[loc];
|
|
loc+=1;
|
|
for(long j = 0; j<long(p); ++j)
|
|
{
|
|
if((t-1-j)<0)
|
|
sigma2[t] += parameters[loc] * backcast;
|
|
else
|
|
sigma2[t] += parameters[loc] * fresids[t-1-j];
|
|
loc+=1;
|
|
}
|
|
for(long j = 0; j<long(o); ++j)
|
|
{
|
|
if((t-1-j)<0)
|
|
sigma2[t] += parameters[loc] *0.5*backcast;
|
|
else
|
|
sigma2[t] += (parameters[loc] * fresids[t - 1 - j] * double((sresids[t - 1 - j] < 0)));
|
|
loc+=1;
|
|
}
|
|
for(long j = 0; j<long(q); ++j)
|
|
{
|
|
if((t - 1 - j) < 0)
|
|
sigma2[t] += parameters[loc] * backcast;
|
|
else
|
|
sigma2[t] += parameters[loc] * sigma2[t - 1 - j];
|
|
loc += 1;
|
|
}
|
|
sigma2[t] = bounds_check(sigma2[t], var_bounds.Row(t));
|
|
}
|
|
return sigma2;
|
|
}
|
|
|
|
//+------------------------------------------------------------------+
|
|
//|Compute variance recursion for EWMA/RiskMetrics Variance |
|
|
//+------------------------------------------------------------------+
|
|
vector ewma_recursion(double lam, vector& resids, vector& sigma2, ulong nobs, double _backcast)
|
|
{
|
|
matrix varbs = matrix::Ones(nobs,2);
|
|
vector temp = {-1.,1.7e308};
|
|
varbs = np::multiply(varbs,temp);
|
|
temp.Resize(3,1);
|
|
temp[0] = 0.0;
|
|
temp[1] = 1.0 - lam;
|
|
temp[2] = lam;
|
|
vector resids2 = pow(resids,2.);
|
|
sigma2 = garch_recursion(temp,resids2,resids,sigma2,1,0,1,nobs,_backcast,varbs);
|
|
return sigma2;
|
|
}
|
|
//+------------------------------------------------------------------+
|
|
//| Base class that all volatility updaters must inherit from. |
|
|
//+------------------------------------------------------------------+
|
|
class VolatilityUpdater
|
|
{
|
|
protected:
|
|
bool m_initialized;
|
|
virtual void update_tester(ulong t, vector& parameters, vector& resids, vector& _sigma2,matrix& varbounds)
|
|
{
|
|
update(t,parameters,resids,_sigma2,varbounds);
|
|
}
|
|
public:
|
|
VolatilityUpdater(void):m_initialized(false)
|
|
{}
|
|
VolatilityUpdater(VolatilityUpdater& other)
|
|
{
|
|
m_initialized = other.is_initialized();
|
|
}
|
|
void operator=(VolatilityUpdater& other)
|
|
{
|
|
m_initialized = other.is_initialized();
|
|
}
|
|
~VolatilityUpdater(void)
|
|
{}
|
|
virtual bool is_initialized(void) { return m_initialized; }
|
|
virtual void initialize_update(vector ¶meters, vector& backcast, ulong _nobs);
|
|
virtual void update(ulong t, vector ¶meters, vector& resids, vector& _sigma2, matrix& varbounds);
|
|
};
|
|
//+------------------------------------------------------------------+
|
|
//| garch updater |
|
|
//+------------------------------------------------------------------+
|
|
class GarchUpdater:public VolatilityUpdater
|
|
{
|
|
private:
|
|
ulong m_p, m_o, m_q;
|
|
double m_power;
|
|
double m_backcast;
|
|
public:
|
|
GarchUpdater(void)
|
|
{
|
|
}
|
|
GarchUpdater(ulong p, ulong o, ulong q, double power)
|
|
{
|
|
initialize(p,o,q,power);
|
|
}
|
|
GarchUpdater(GarchUpdater& other)
|
|
{
|
|
m_p = get_p();
|
|
m_o = get_o();
|
|
m_q = get_q();
|
|
m_power = get_power();
|
|
m_backcast = get_backcast();
|
|
m_initialized = other.is_initialized();
|
|
}
|
|
|
|
void operator=(GarchUpdater& other)
|
|
{
|
|
m_p = get_p();
|
|
m_o = get_o();
|
|
m_q = get_q();
|
|
m_power = get_power();
|
|
m_backcast = get_backcast();
|
|
m_initialized = other.is_initialized();
|
|
}
|
|
~GarchUpdater(void)
|
|
{
|
|
}
|
|
ulong get_p(void) { return m_p; }
|
|
ulong get_o(void) { return m_o; }
|
|
ulong get_q(void) { return m_q; }
|
|
double get_power(void) { return m_power; }
|
|
double get_backcast(void) { return m_backcast; }
|
|
bool initialize(ulong p, ulong o, ulong q, double power)
|
|
{
|
|
m_p = p;
|
|
m_o = o;
|
|
m_backcast = - 1.;
|
|
m_q = q;
|
|
m_power = power;
|
|
m_initialized=(m_p||m_o);
|
|
return m_initialized;
|
|
}
|
|
virtual void initialize_update(vector& parameters, vector& backcast, ulong nobs) override
|
|
{
|
|
m_backcast = backcast[0];
|
|
}
|
|
virtual void update(ulong t, vector ¶meters, vector& resids, vector& sigma2, matrix& var_bounds) override
|
|
{
|
|
ulong loc = 0;
|
|
sigma2[t] = parameters[loc];
|
|
loc += 1;
|
|
|
|
for(ulong j = 0; j<m_p; ++j)
|
|
{
|
|
if((t - 1 - j) < 0)
|
|
sigma2[t] += parameters[loc] * m_backcast;
|
|
else
|
|
sigma2[t] += parameters[loc] * pow(MathAbs(resids[t - 1 - j]),m_power);
|
|
loc += 1;
|
|
}
|
|
for(ulong j = 0; j<m_o; ++j)
|
|
{
|
|
if((t - 1 - j) < 0)
|
|
sigma2[t] += parameters[loc] * 0.5 * m_backcast;
|
|
else
|
|
sigma2[t] += (
|
|
parameters[loc]
|
|
* pow(MathAbs(resids[t - 1 - j]), m_power)
|
|
* (resids[t - 1 - j] < 0)
|
|
);
|
|
loc += 1;
|
|
}
|
|
for(ulong j = 0; j<m_q; ++j)
|
|
{
|
|
if((t - 1 - j) < 0)
|
|
sigma2[t] += parameters[loc] * m_backcast;
|
|
else
|
|
sigma2[t] += parameters[loc] * sigma2[t - 1 - j];
|
|
loc += 1;
|
|
}
|
|
sigma2[t] = bounds_check(sigma2[t], var_bounds.Row(t));
|
|
}
|
|
};
|
|
/*+------------------------------------------------------------------+
|
|
//| EWMAUpdater |
|
|
//+------------------------------------------------------------------+
|
|
class EWMAUpdater: public VolatilityUpdater
|
|
{
|
|
private:
|
|
bool m_estimate_lam;
|
|
double m_backcast;
|
|
vector m_params;
|
|
EWMAUpdater(void)
|
|
{
|
|
|
|
}
|
|
~EWMAUpdater(void)
|
|
{
|
|
|
|
}
|
|
virtual bool initialize(ulong p,ulong o, ulong q,double power) override
|
|
{
|
|
Print(__FUNCTION__, " wrong function call ");
|
|
return false;
|
|
}
|
|
virtual bool initialize(double lam) override
|
|
{
|
|
m_estimate_lam = (lam == EMPTY_VALUE);
|
|
m_params = vector::Zeros(3);
|
|
if(!m_estimate_lam)
|
|
{
|
|
m_params[1] = 1.-lam;
|
|
m_params[2] = lam;
|
|
}
|
|
|
|
m_initialized = true;
|
|
return m_initialized;
|
|
}
|
|
virtual void initialize_update(vector ¶meters, vector& backcast, ulong _nobs) override
|
|
{
|
|
if(m_estimate_lam)
|
|
{
|
|
m_params[1] = 1.0 - parameters[0];
|
|
m_params[2] = parameters[0];
|
|
}
|
|
m_backcast = backcast[0];
|
|
}
|
|
virtual void update(ulong t, vector ¶meters, vector& resids, vector& sigma2, matrix& var_bounds) override
|
|
{
|
|
sigma2[t] = m_params[0];
|
|
if(t == 0)
|
|
sigma2[t] += m_backcast;
|
|
else
|
|
sigma2[t] += (
|
|
m_params[1] * resids[t - 1] * resids[t - 1]
|
|
+ m_params[2] * sigma2[t - 1]
|
|
);
|
|
sigma2[t] = bounds_check(sigma2[t], var_bounds.Row(t));
|
|
}
|
|
};
|
|
//+------------------------------------------------------------------+
|
|
//|EGARCHUpdater |
|
|
//+------------------------------------------------------------------+
|
|
class EGARCHUpdater: public VolatilityUpdater
|
|
{
|
|
private:
|
|
ulong m_p, m_o, m_q;
|
|
double m_backcast;
|
|
vector m_insigma2, m_std_resids, m_abs_std_resids;
|
|
void _resize(ulong nobs)
|
|
{
|
|
if(m_insigma2.Size()<nobs)
|
|
m_insigma2 = m_std_resids = m_abs_std_resids = vector::Zeros(nobs);
|
|
}
|
|
public:
|
|
EGARCHUpdater(ulong p, ulong o, ulong q)
|
|
{
|
|
m_p = p;
|
|
m_o = o;
|
|
m_q = q;
|
|
m_backcast = 9999.99;
|
|
m_insigma2 = m_std_resids = m_abs_std_resids = vector::Zeros(0);
|
|
}
|
|
~EGARCHUpdater(void)
|
|
{
|
|
}
|
|
virtual void initialize_update(vector ¶meters, vector& backcast, ulong _nobs) override
|
|
{
|
|
m_backcast = backcast[0];
|
|
_resize(_nobs);
|
|
}
|
|
virtual void update(ulong t, vector ¶meters, vector& resids, vector& sigma2, matrix& var_bounds) override
|
|
{
|
|
if(t)
|
|
{
|
|
m_std_resids[t - 1] = resids[t - 1] / sqrt(sigma2[t - 1]);
|
|
m_abs_std_resids[t - 1] = MathAbs(m_std_resids[t - 1]);
|
|
}
|
|
|
|
m_insigma2[t] = parameters[0];
|
|
ulong loc = 1;
|
|
for(ulong j = 0; j<m_p; ++j)
|
|
{
|
|
if((t - 1 - j) >= 0)
|
|
m_insigma2[t] += parameters[loc] * (m_abs_std_resids[t - 1 - j] - SQRT2_OV_PI);
|
|
loc += 1;
|
|
}
|
|
for(ulong j = 0; j<m_o; ++j)
|
|
{
|
|
if((t - 1 - j) >= 0)
|
|
m_insigma2[t] += parameters[loc] * m_std_resids[t - 1 - j];
|
|
loc += 1;
|
|
}
|
|
for(ulong j = 0; j<m_q; ++j)
|
|
{
|
|
if((t - 1 - j) < 0)
|
|
m_insigma2[t] += parameters[loc] * m_backcast;
|
|
else
|
|
m_insigma2[t] += parameters[loc] * m_insigma2[t - 1 - j];
|
|
loc += 1;
|
|
}
|
|
if(m_insigma2[t] > LNSIGMA_MAX)
|
|
m_insigma2[t] = LNSIGMA_MAX;
|
|
sigma2[t] = exp(m_insigma2[t]);
|
|
if(sigma2[t] < var_bounds[t, 0])
|
|
{
|
|
sigma2[t] = var_bounds[t, 0];
|
|
m_insigma2[t] = log(sigma2[t]);
|
|
}
|
|
else
|
|
if(sigma2[t] > var_bounds[t, 1])
|
|
{
|
|
sigma2[t] = var_bounds[t, 1] + log(sigma2[t]) - log(var_bounds[t, 1]);
|
|
m_insigma2[t] = log(sigma2[t]);
|
|
}
|
|
}
|
|
};*/
|
|
|
|
//+------------------------------------------------------------------+
|