MetaTrader 5 / 例
English
preview
MQL5における二変量コピュラ（第2回）：MQL5でのアルキメデスコピュラの実装

MQL5における二変量コピュラ（第2回）：MQL5でのアルキメデスコピュラの実装

MetaTrader 5 — 例 |
15 0
Francis Dube
Francis Dube

はじめに

前回の記事では、二変量楕円型コピュラ、具体的には正規コピュラおよびtコピュラのMQL5実装を紹介しました。これらのコピュラを応用するにあたり、少なくとも2つの課題があります。第一に、代数的表現が複雑であること、第二に、対称性が高いためデータに含まれる非対称な依存構造を十分に表現できない場合があることです。これらの課題を背景として、より扱いやすい代数表現を持ち、かつ非対称な依存関係をモデル化できるコピュラの研究が進められてきました。その中でも特に広く利用されているのが、アルキメデスコピュラのクラスです。本記事では、代表的な二変量アルキメデスコピュラの特性を解説し、それらをMQL5で実装する方法について説明します。また、既存の二変量コピュラライブラリを拡張し、フランク、ジョー、グンベル、クレイトン、N13、N14コピュラの実装を追加していきます。


アルキメデスコピュラの紹介

二変量アルキメデスコピュラは、均一な周辺分布を持つ2つのランダム変数間の依存構造をモデル化するために統計で使用される特定のタイプのコピュラC(u,v)です。その本質的な特徴は生成表現にあり、これはジェネレータ関数と呼ばれる、単一の連続かつ厳密単調減少で凸な関数ϕによって表現されます。

ジェネレータ関数

生成元ϕは条件ϕ(1)=0を満たす必要があります。この構造により、アルキメデスコピュラは高い対称性(C(u,v)=C(v,u))を持ちます。また、ジェネレータ関数を変更するだけで多様な依存構造を簡潔に表現できるという利点があります。アルキメデスコピュラは、本質的には一般的なコピュラよりも柔軟性に優れていますが、実務上はその柔軟性は通常ある程度制限されます。理論的には、ジェネレータ関数の選択肢は無限に存在します。関数にわずかな変更を加えるだけで、新たに固有のコピュラが生成されます。しかし実際には、単一のスカラー入力変数によって定義される固有のジェネレータ関数により規定されるパラメトリックファミリーを選択するのが一般的です。これにより、パラメータ推定およびモデリングが大幅に容易になります。そのため、代表的なアルキメデスコピュラの多くは、ジェネレータ関数に組み込まれた単一のパラメータによって特徴付けられます。このパラメータが依存の強さを制御します。

ジェネレータ関数はアルキメデスコピュラを定義し、確率変数間の依存構造全体を捉えます。平易に言えば、その役割は周辺確率を、単純な加算によって容易に結合できる依存構造へと変換することにあります。この関数は周辺変数を変換します。ジェネレータは厳密単調減少であるため、低い確率は大きな値へ写像され、高い確率は0に写像されます。これは確率スケールの反転に相当します。ジェネレータ関数の具体的な形状およびパラメータが、結果として得られるコピュラファミリーを完全に決定し、ひいては2変数間の依存関係の正確なあり方を規定します。次のセクションでは、二変量フランクコピュラからアルキメデスコピュラの検討を開始します。


二変量フランクコピュラ

二変量フランクコピュラは、裾依存性を示さずに正の依存性と負の依存性の両方をモデル化できる対称アルキメデスコピュラです。これは以下のジェネレータ関数によって定義されます。

フランクジェネレータ関数

θは依存関係の強さと方向を決定します。フランクコピュラのジェネレータにおいて、θが0の場合は定義されません。θが正の値であれば正の依存関係を示し、負の値であれば負の依存関係を示します。θが0に近づくと、コピュラは独立コピュラに近づきます。フランクコピュラは放射対称であり、分布の両裾を同等に扱います。そのため、上側裾依存性や上側裾依存性は示しません。この特性により、分布の中心付近で対称的な中程度の依存構造をモデル化する場合に適しています。たとえば、極端な同時変動が起こりにくい連続変数の解析に向いています。

frank.mqh内で定義されるCFrankクラスは、CBivariateCopulaを継承し、フランクコピュラを表現します。

//+------------------------------------------------------------------+
//|                                                        frank.mqh |
//|                                  Copyright 2025, MetaQuotes Ltd. |
//|                                             https://www.mql5.com |
//+------------------------------------------------------------------+
#property copyright "Copyright 2025, MetaQuotes Ltd."
#property link      "https://www.mql5.com"
#include "base.mqh"
//+------------------------------------------------------------------+
//| Frank Copula.                                                    |
//+------------------------------------------------------------------+
class CFrank:public CBivariateCopula
  {
private:
   class CInt_Function_1_Func : public CIntegrator1_Func
     {
   private:
      double         m_th;
   public:
      //---
                     CInt_Function_1_Func(void) {}
                    ~CInt_Function_1_Func(void) {}
      void              set_theta(double theta)
        {
         m_th = theta;
        }
      virtual void      Int_Func(double x,double xminusa,double bminusx,double &y,CObject &obj)
        {
         y = x/m_th/(exp(x)-1.0);
        }
     };
   class CBrent1:public CBrentQ
     {
   private:
      double            debyel(double theta)
        {
         CInt_Function_1_Func fint;
         fint.set_theta(theta);
         CObject obj;
         CAutoGKStateShell s;
         double integral;
         CAutoGKReportShell rep;
         //---
         CAlglib::AutoGKSmooth(0.0,theta,s);
         //---
         CAlglib::AutoGKIntegrate(s,fint,obj);
         //---
         CAlglib::AutoGKResults(s,integral,rep);

         CAutoGKReport report = rep.GetInnerObj();

         if(report.m_terminationtype<0.0)
            Print(__FUNCTION__, " integration error ",report.m_terminationtype);

         return integral;
        }
   public:
                     CBrent1(void)
        {
        }
                    ~CBrent1(void)
        {
        }
      virtual double objective(double u, double v,double y)
        {
         return (1.0 - 4.0/u+4.0*debyel(u)/u) - y;
        }
     };

   virtual double    theta_hat(const double tau) override
     {
      CBrent1 fun;
      return fun.minimize(0,tau,-100.0,100.0,2.e-12,4.0*2.e-16,100);
     }
   vector            generate_pair(double v1, double v2)
     {
      vector out(2);
      out[0] = v1;
      out[1] = -1.0/m_theta*log(1.0+(v2*(1.0-exp(-m_theta)))/(v2*(exp(-m_theta*v1)-1.0)-exp(-m_theta*v1)));
      return out;
     }
protected:
   virtual double    pdf(double u,double v) override
     {
      double et,eut,evt,pd;
      et = exp(m_theta);
      eut = exp(u*m_theta);
      evt = exp(v*m_theta);
      pd = et*eut*evt*(et - 1.0)*m_theta/pow(et+eut*evt-et*eut-et*evt,2.0);
      return pd;
     }
   virtual double    cdf(double u, double v) override
     {
      return (-1.0 / m_theta * log(1.0 + (exp(-1.0 * m_theta * u) - 1.0) * (exp(-1.0 * m_theta * v) - 1.0)/ (exp(-1 * m_theta) - 1.0)));
     }
   virtual double    condi_cdf(double u, double v) override
     {
      double enut = exp(-u*m_theta);
      double envt = exp(-v*m_theta);
      double ent = exp(-m_theta);
      double denominator = ((ent - 1.0) + (enut - 1.0) * (envt - 1.0));
      if(denominator)
         return (envt * (enut - 1.0)/ (denominator));
      else
         return EMPTY_VALUE;
     }
   virtual vector    pdf(vector &u,vector &v) override
     {
      vector eut,evt,pd;
      double et = exp(m_theta);
      eut = exp(u*m_theta);
      evt = exp(v*m_theta);
      pd = et*eut*evt*(et - 1.0)*m_theta/pow(et+eut*evt-et*eut-et*evt,2.0);
      return pd;
     }
   virtual vector    cdf(vector &u, vector &v) override
     {
      return (-1.0 / m_theta * log(1.0 + (exp(-1.0 * m_theta * u) - 1.0) * (exp(-1.0 * m_theta * v) - 1.0)/ (exp(-1 * m_theta) - 1.0)));
     }
   virtual vector    condi_cdf(vector &u, vector &v) override
     {
      vector enut = exp(-1.0*u*m_theta);
      vector envt = exp(-1.0*v*m_theta);
      double ent = exp(-m_theta);

      return (envt * (enut - 1.0)/ ((ent - 1.0) + (enut - 1.0) * (envt - 1.0)));
     }
public:
                     CFrank(void)
     {
      m_threshold = 1.e-10;
      m_copula_type = FRANK_COPULA;
      m_invalid_theta = 0.0;
     }
                    ~CFrank(void)
     {
     }
   virtual matrix    Sample(ulong num_samples) override
     {
      double unf_v[],unf_c[];

      vector v,c,u;

      if(!MathRandomUniform(0.0,1.0,int(num_samples),unf_v) ||
         !MathRandomUniform(0.0,1.0,int(num_samples),unf_c) ||
         !v.Assign(unf_v) || !c.Assign(unf_c))
        {
         Print(__FUNCTION__, " failed to get uniform random numbers ", GetLastError());
         return matrix::Zeros(0,0);
        }

      matrix out(num_samples,2);

      for(ulong irow = 0; irow<num_samples; irow++)
         out.Row(generate_pair(v[irow],c[irow]),irow);

      return out;

     }

  };
//+------------------------------------------------------------------+

θパラメータは、データサンプルから推定したランク相関係数τ（ケンドールのタウ）に基づいて決定されます。

シータとケンドールのタウの関係


ここで、D₁​(θ)は1次のデバイ関数です。下図は、θパラメータが5のフランクコピュラからサンプリングした合成データの散布図です。

フランクコピュラのサンプル


二変量クレイトンコピュラ

次に、二変量クレイトンコピュラについて見ていきます。クレイトンコピュラは、強い下側裾依存性や、結合分布の尾部における非対称性をモデル化できる点が特徴です。これは以下のジェネレータ関数によって定義されます。

クレイトンジェネレータ関数

θの値が大きいほど、下側裾依存性は強くなります。クレイトンコピュラの下側裾依存係数は次の式で表されます。

クレイトン下側裾依存性

上側裾依存性は存在せず(λU​=0)、極端に低い値が同時に発生する状況をモデル化するのに適しています。θ>0の範囲でのみ正の依存関係を表現し、θが0に近づくと独立コピュラに収束します。θ = 0の場合、クレイトンコピュラは定義されません。クレイトンコピュラは非対称であるため、同時に値が下がるケースを表現するのに効果的です。一方、値が同時に上昇するケースは捉えません。clayton.mqhヘッダファイルには、クレイトンコピュラを表現するCClaytonクラスの定義が含まれています。

//+------------------------------------------------------------------+
//|                                                      clayton.mqh |
//|                                  Copyright 2025, MetaQuotes Ltd. |
//|                                             https://www.mql5.com |
//+------------------------------------------------------------------+
#property copyright "Copyright 2025, MetaQuotes Ltd."
#property link      "https://www.mql5.com"
#include "base.mqh"
//+------------------------------------------------------------------+
//| Clayton copula.                                                  |
//+------------------------------------------------------------------+
class CClayton:public CBivariateCopula
  {
private:
   virtual double    theta_hat(const double tau) override
     {
      return 2.0*tau/(1.0-tau);
     }
   vector            generate_pair(double v1, double v2)
     {
      vector out(2);
      double w = 0.0;

      out[0] = v1;
      out[1] = pow(pow(v1,-m_theta)*(pow(v2,(-m_theta / (1 + m_theta))) - 1.0)+1.0,-1.0/m_theta);
      return out;
     }
protected:
   virtual double    pdf(double u,double v) override
     {
      double u_part,v_part,pd;
      u_part = pow(u,(-1 - m_theta));
      v_part = pow(v,(-1 - m_theta));
      pd = ((1.0 + m_theta) * u_part * v_part*pow(-1.0 + u_part * u + v_part * v, (-2.0 - 1.0 / m_theta)));
      return pd;
     }
   virtual double    cdf(double u, double v) override
     {
      double expo = pow(MathMax(pow(u,-m_theta)+pow(v,-m_theta)-1.0,0.0),-1.0/m_theta);

      return expo;
     }
   virtual double    condi_cdf(double u, double v) override
     {
      double unt = pow(u,-m_theta);
      double vnt = pow(v,-m_theta);
      double tpow = 1.0/m_theta + 1.0;
      double denominator = pow(unt+vnt-1.0,tpow);
      if(denominator)
         return vnt/v/pow(unt+vnt-1.0,tpow);
      else
         return EMPTY_VALUE;
     }
   virtual double    inv_condi_cdf(double y, double V) override
     {
      if(m_theta < 0.0)
         return V;
      else
        {
         double a = pow(y,m_theta/(-1.0 - m_theta));
         double b = pow(V, m_theta);

         if(b==0)
            return 1.0;
         return pow((a + b - 1.0)/b, -1.0/m_theta);
        }
     }
   virtual vector    inv_condi_cdf(vector& y, vector& V) override
     {
      if(m_theta < 0.0)
         return V;
      else
        {
         vector a = pow(y,m_theta/(-1.0 - m_theta));
         vector b = pow(V, m_theta);

         if(b.CompareEqual(vector::Zeros(b.Size()))==0)
            return vector::Ones(V.Size());
         return pow((a + b - 1.0)/b, -1.0/m_theta);
        }
     }
   virtual vector    pdf(vector &u,vector &v) override
     {
      vector u_part,v_part,pd;
      u_part = pow(u,(-1 - m_theta));
      v_part = pow(v,(-1 - m_theta));
      pd = ((1.0 + m_theta) * u_part * v_part*pow(-1.0 + u_part * u + v_part * v, (-2.0 - 1.0 / m_theta)));
      return pd;
     }
   virtual vector    cdf(vector &u, vector &v) override
     {
      vector out(u.Size());
      for(ulong i = 0; i<u.Size(); ++i)
         out[i] = cdf(u[i],v[i]);
      return out;
     }
   virtual vector    condi_cdf(vector &u, vector &v) override
     {
      vector unt = pow(u,-m_theta);
      vector vnt = pow(v,-m_theta);
      double tpow = 1.0/m_theta + 1.0;
      return vnt/v/pow(unt+vnt-1.0,tpow);
     }
public:
                     CClayton(void)
     {
      m_threshold = 1.e-10;
      m_bounds[0] = -1.0;
      m_bounds[1] = DBL_MAX;
      m_invalid_theta = 0.0;
      m_copula_type = CLAYTON_COPULA;
     }
                    ~CClayton(void)
     {
     }
   virtual matrix    Sample(ulong num_samples) override
     {
      double unf_v[],unf_c[];

      vector v,c,u;

      if(!MathRandomUniform(0.0,1.0,int(num_samples),unf_v) ||
         !MathRandomUniform(0.0,1.0,int(num_samples),unf_c) ||
         !v.Assign(unf_v) || !c.Assign(unf_c))
        {
         Print(__FUNCTION__, " failed to get uniform random numbers ", GetLastError());
         return matrix::Zeros(0,0);
        }

      matrix out(num_samples,2);

      for(ulong irow = 0; irow<num_samples; irow++)
         out.Row(generate_pair(v[irow],c[irow]),irow);

      return out;

     }

  };
//+------------------------------------------------------------------+

ここでは、クレイトン二変量分布からサンプリングした合成データの散布図を示します。 

二変量クレイトンサンプルデータセット


二変量グンベルコピュラ

二変量グンベルコピュラはアルキメデスコピュラの一種で、上側裾依存性を効果的にモデル化し、上側裾と下側裾の間に非対称性を持ちます。これは以下のジェネレータ関数によって定義されます。

グンベルジェネレータ関数

θの値が大きいほど依存関係は強くなります。グンベルコピュラは上側裾依存係数を持ちます。

グンベルの上側裾依存性

また、下側裾依存性は存在しないため、ある変数の極端に高い値が他の変数の極端に高い値と同時に現れる傾向を捉えるのに有用です。グンベルコピュラは正の依存関係のみを表現でき、θ = 1の場合は独立コピュラに収束します。グンベルコピュラは、gumbel.mqhに定義されたCGumbelクラスとして実装されています。

//+------------------------------------------------------------------+
//|                                                       gumbel.mqh |
//|                                  Copyright 2025, MetaQuotes Ltd. |
//|                                             https://www.mql5.com |
//+------------------------------------------------------------------+
#property copyright "Copyright 2025, MetaQuotes Ltd."
#property link      "https://www.mql5.com"
#include "base.mqh"
//+------------------------------------------------------------------+
//| Gumbel Copula.                                                   |
//+------------------------------------------------------------------+
class CGumbel:public CBivariateCopula
  {
private:
   class CBrent1:public CBrentQ
     {
   public:
                     CBrent1(void)
        {
        }
                    ~CBrent1(void)
        {
        }
      virtual double objective(double u, double v,double y)
        {
         return (u*(1.0 - log(u)/v)) - y;
        }
     };
   virtual double    theta_hat(const double tau) override
     {
      return 1.0/(1.0-tau);
     }
   vector            generate_pair(double v1, double v2)
     {
      vector out(2);
      double w = 0.0;
      if(v2>m_threshold)
        {
         CBrent1 fun;
         w = fun.minimize(m_theta,v2,m_threshold,1.0,2.e-12,4.0*2.e-16,100);
        }
      else
         w = 1.e10;
      out[0] = exp(pow(v1,1.0/m_theta)*log(w));
      out[1] = exp(pow((1.0-v1),1.0/m_theta)*log(w));
      return out;
     }
protected:
   virtual double    pdf(double u,double v) override
     {
      double u_part,v_part,expo,pd;
      u_part = pow(-1.0*log(u),m_theta);
      v_part = pow(-1.0*log(v),m_theta);
      expo = pow(u_part+v_part,1.0/m_theta);
      pd = 1.0/(u*v) * (exp(-1.0*expo)* u_part/(-1.0*log(u)) * v_part/(-1.0*log(v))*(m_theta+expo-1.0)*pow(u_part+v_part,(1.0/m_theta-2.0)));
      return pd;
     }
   virtual double    cdf(double u, double v) override
     {
      double expo = pow(pow(-1.0*log(u),m_theta) + pow(-1.0*log(v),m_theta),(1. / m_theta));

      return exp(-1.0*expo);
     }
   virtual double    condi_cdf(double u, double v) override
     {
      double expo = pow(pow(-1.0*log(u),m_theta)+pow(-1.0*log(v),m_theta),((1 - m_theta) / m_theta));
      return (cdf(u,v) * expo * pow(-1.0*log(v),m_theta-1.0)/v);
     }
   virtual vector    pdf(vector &u,vector &v) override
     {
      vector u_part,v_part,expo,pd;
      u_part = pow(-1.0*log(u),m_theta);
      v_part = pow(-1.0*log(v),m_theta);
      expo = pow(u_part+v_part,1.0/m_theta);
      pd = 1.0/(u*v) * (exp(-1.0*expo)* u_part/(-1.0*log(u)) * v_part/(-1.0*log(v))*(m_theta+expo-1.0)*pow(u_part+v_part,(1.0/m_theta-2.0)));
      return pd;
     }
   virtual vector    cdf(vector &u, vector &v) override
     {

      vector expo = pow(pow(-1.0*log(u),m_theta) + pow(-1.0*log(v),m_theta),(1. / m_theta));
      return exp(-1.0*expo);
     }
   virtual vector    condi_cdf(vector &u, vector &v) override
     {

      vector expo = pow(pow(-1.0*log(u),m_theta)+pow(-1.0*log(v),m_theta),((1 - m_theta) / m_theta));
      return (cdf(u,v) * expo * pow(-1.0*log(v),m_theta-1.0)/v);
     }
public:
                     CGumbel(void)
     {
      m_threshold = 1.e-10;
      m_bounds[0] = 1.0;
      m_bounds[1] = DBL_MAX;
      m_copula_type = GUMBEL_COPULA;
     }
                    ~CGumbel(void)
     {
     }
   virtual matrix    Sample(ulong num_samples) override
     {
      double unf_v[],unf_c[];

      vector v,c,u;

      if(!MathRandomUniform(0.0,1.0,int(num_samples),unf_v) ||
         !MathRandomUniform(0.0,1.0,int(num_samples),unf_c) ||
         !v.Assign(unf_v) || !c.Assign(unf_c))
        {
         Print(__FUNCTION__, " failed to get uniform random numbers ", GetLastError());
         return matrix::Zeros(0,0);
        }

      matrix out(num_samples,2);

      for(ulong irow = 0; irow<num_samples; irow++)
         out.Row(generate_pair(v[irow],c[irow]),irow);

      return out;

     }

  };
//+------------------------------------------------------------------+

下図は、θパラメータが5のグンベルコピュラからサンプリングした合成データの散布図です。

二変量グンベルサンプルデータセット


二変量ジョーコピュラ

二変量ジョーコピュラはアルキメデスコピュラの一種で、上側裾依存性が強く、上側裾と下側裾の間に非対称性を持つことで知られています。ある変数の極端に高い値が、もう一方の変数の極端に高い値と同時に現れやすい状況を捉えますが、低い値には同様の依存はありません。このコピュラは以下のジェネレータ関数によって定義されます。

ジョージェネレータ関数

θの値が大きいほど依存関係は強くなります。ジョーコピュラは常に正の依存関係のみを表現し、非楕円型コピュラに属するため、非線形な関係のモデリングにも柔軟です。重要なのは、上側裾依存係数を持つということです。

ジョーの裾依存性

下側裾依存性は存在しません。ジョーコピュラはjoe.mqhに実装されています。

//+------------------------------------------------------------------+
//|                                                          joe.mqh |
//|                                  Copyright 2025, MetaQuotes Ltd. |
//|                                             https://www.mql5.com |
//+------------------------------------------------------------------+
#property copyright "Copyright 2025, MetaQuotes Ltd."
#property link      "https://www.mql5.com"
#include "base.mqh"
//+------------------------------------------------------------------+
//| Joe Copula.                                                      |
//+------------------------------------------------------------------+
class CJoe:public CBivariateCopula
  {
private:
   class CInt_Function_1_Func : public CIntegrator1_Func
     {
   private:
      double         m_th;
   public:
      //---
                     CInt_Function_1_Func(void) {}
                    ~CInt_Function_1_Func(void) {}
      void              set_theta(double theta)
        {
         m_th = theta;
        }
      virtual void      Int_Func(double x,double xminusa,double bminusx,double &y,CObject &obj)
        {
         y = (1.0 - pow(1.0 - x,m_th)) * pow(1.0 - x,1.0 - m_th) * log(1.0 - pow(1.0 - x,m_th)) / m_th;
        }
     };
   class CBrent1:public CBrentQ
     {
   public:
                     CBrent1(void)
        {
        }
                    ~CBrent1(void)
        {
        }
      virtual double objective(double u, double v,double y)
        {
         return (u-1.0/v*(log(1.0-pow(1.0-u,v))*(1.0-pow(1.0-u,v))/pow(1.0-u,(v-1.0)))) - y;
        }
     };
   class CBrent2:public CBrentQ
     {
   public:
                     CBrent2(void)
        {
        }
                    ~CBrent2(void)
        {
        }
      virtual double objective(double u, double v,double y)
        {
         CInt_Function_1_Func fint;
         fint.set_theta(u);
         CObject obj;
         CAutoGKStateShell s;
         double integral;
         CAutoGKReportShell rep;
         //---
         CAlglib::AutoGKSmooth(0.0,1.0,s);
         //---
         CAlglib::AutoGKIntegrate(s,fint,obj);
         //---
         CAlglib::AutoGKResults(s,integral,rep);

         CAutoGKReport report = rep.GetInnerObj();

         if(report.m_terminationtype<0.0)
            Print(__FUNCTION__, " integration error ",report.m_terminationtype);

         return (1.0+4.0*integral)-y;
        }
     };
   virtual double    theta_hat(const double tau) override
     {
      CBrent2 fun;
      return fun.minimize(0,tau,1.0,100.0,2.e-12,4.0*2.e-16,100);
     }
   vector            generate_pair(double v1, double v2)
     {
      vector out(2);
      double w = 0;
      if(v2>m_threshold)
        {
         CBrent1 fun;
         w = fun.minimize(m_theta,v2,m_threshold,1.0-m_threshold,2.e-12,4.0*2.e-16,100);
        }
      else
         w = m_threshold;
      out[0] = 1.0 - pow(1.0 - pow(1.0 - pow(1.0 - w, m_theta), v1),(1.0 / m_theta));
      out[1] = 1.0 - pow(1.0 - pow(1.0 - pow(1.0 - w, m_theta), (1.0-v1)),(1.0 / m_theta));
      return out;
     }
protected:
   virtual double    pdf(double u,double v) override
     {
      double up,vp,pd;
      up = pow(1.0-u,m_theta);
      vp = pow(1.0-v,m_theta);
      pd = (up/(1.0-u)*vp/(1.0-v)*pow(up+vp-up*vp,1.0/m_theta-2.0)*(m_theta-(up-1.0)*(vp-1.0)));
      return pd;
     }
   virtual double    cdf(double u, double v) override
     {
      double up = pow(1.0-u,m_theta);
      double vp = pow(1.0-v,m_theta);
      return 1.0 - pow(up+vp-up*vp,1.0/m_theta);
     }
   virtual double    condi_cdf(double u, double v) override
     {
      double up = pow(1.0-u,m_theta);
      double vp = pow(1.0-v,m_theta);
      return (-(-1.0+up)*pow(up+vp-up*vp,-1.0+1.0/m_theta)*vp/(1.0-v));
     }
   virtual vector    pdf(vector &u,vector &v) override
     {
      vector up,vp,pd;
      up = pow(1.0-u,m_theta);
      vp = pow(1.0-v,m_theta);
      pd = (up/(1.0-u)*vp/(1.0-v)*pow(up+vp-up*vp,1.0/m_theta-2.0)*(m_theta-(up-1.0)*(vp-1.0)));
      return pd;
     }
   virtual vector    cdf(vector &u, vector &v) override
     {
      vector up = pow(1.0-u,m_theta);
      vector vp = pow(1.0-v,m_theta);
      return 1.0 - pow(up+vp-up*vp,1.0/m_theta);
     }
   virtual vector    condi_cdf(vector &u, vector &v) override
     {

      vector up = pow(1.0-u,m_theta);
      vector vp = pow(1.0-v,m_theta);
      return (-1.0*(-1.0+up)*pow(up+vp-up*vp,-1.0+1.0/m_theta)*vp/(1.0-v));
     }
public:
                     CJoe(void)
     {
      m_threshold = 1.e-10;
      m_copula_type = JOE_COPULA;
      m_bounds[0] = 1.0;
      m_bounds[1] = DBL_MAX;
     }
                    ~CJoe(void)
     {
     }
   virtual matrix    Sample(ulong num_samples) override
     {
      double unf_v[],unf_c[];

      vector v,c,u;

      if(!MathRandomUniform(0.0,1.0,int(num_samples),unf_v) ||
         !MathRandomUniform(0.0,1.0,int(num_samples),unf_c) ||
         !v.Assign(unf_v) || !c.Assign(unf_c))
        {
         Print(__FUNCTION__, " failed to get uniform random numbers ", GetLastError());
         return matrix::Zeros(0,0);
        }

      matrix out(num_samples,2);

      for(ulong irow = 0; irow<num_samples; irow++)
         out.Row(generate_pair(v[irow],c[irow]),irow);

      return out;

     }

  };
//+------------------------------------------------------------------+

下図は、θパラメータが5のジョーコピュラからサンプリングした合成データセットで表される依存関係です。

二変量ジョーサンプルデータセット


二変量N13コピュラ

N13コピュラは、アルキメデス型の1パラメータ族で、正の依存関係のみをモデル化します。他のアルキメデスコピュラと異なり、漸近的な裾依存性は示しません。θの値が大きくなるほど、全体の一致度が強くなります。θ = 0の場合、独立コピュラに近づきます。裾係数が0であるため、極端な同時変動のクラスタリングではなく、分布の中間〜上中位での依存を重視したモデリングに適しています。これは、正規コピュラやフランクコピュラよりも、中間領域の依存を強く捉えたい場合に有効です。N13コピュラは以下のジェネレータ関数によって定義されます。

n13ジェネレータ関数

MQL5では、n13.mqh内でCN13クラスとして実装されています。 

//+------------------------------------------------------------------+
//|                                                          n13.mqh |
//|                                  Copyright 2025, MetaQuotes Ltd. |
//|                                             https://www.mql5.com |
//+------------------------------------------------------------------+
#property copyright "Copyright 2025, MetaQuotes Ltd."
#property link      "https://www.mql5.com"
#include "base.mqh"
//+------------------------------------------------------------------+
//| N13 Copula (Nelsen 13)                                           |
//+------------------------------------------------------------------+
class CN13:public CBivariateCopula
  {
private:
   class CInt_Function_1_Func : public CIntegrator1_Func
     {
   private:
      double         m_th;
   public:
      //---
                     CInt_Function_1_Func(void) {}
                    ~CInt_Function_1_Func(void) {}
      void              set_theta(double theta)
        {
         m_th = theta;
        }
      virtual void      Int_Func(double x,double xminusa,double bminusx,double &y,CObject &obj)
        {
         y = -((x - x * pow((1 - log(x)),1 - m_th) - x * log(x)) / m_th);
        }
     };
   class CBrent1:public CBrentQ
     {
   public:
                     CBrent1(void)
        {
        }
                    ~CBrent1(void)
        {
        }
      virtual double objective(double u, double v,double y)
        {
         return (u+1.0/v*(u-u*pow(1.0-1.0*log(u),1.0-v)-u*log(u))) - y;;
        }
     };
   class CBrent2:public CBrentQ
     {
   public:
                     CBrent2(void)
        {
        }
                    ~CBrent2(void)
        {
        }
      virtual double objective(double u, double v,double y)
        {
         CInt_Function_1_Func fint;
         fint.set_theta(u);
         CObject obj;
         CAutoGKStateShell s;
         double integral;
         CAutoGKReportShell rep;
         //---
         CAlglib::AutoGKSmooth(0.0,1.0,s);
         //---
         CAlglib::AutoGKIntegrate(s,fint,obj);
         //---
         CAlglib::AutoGKResults(s,integral,rep);

         CAutoGKReport report = rep.GetInnerObj();

         if(report.m_terminationtype<0.0)
            Print(__FUNCTION__, " integration error ",report.m_terminationtype);

         return (1.0 + 4.0*integral) - y;
        }
     };
   virtual double    theta_hat(const double tau) override
     {
      CBrent2 fun;
      return fun.minimize(0,tau,1.e-7,100.0,2.e-12,4.0*2.e-16,100);
     }
   vector            generate_pair(double v1, double v2)
     {
      vector out(2);
      double w = 0.0;
      if(v2>m_threshold)
        {
         CBrent1 fun;
         w = fun.minimize(m_theta,v2,m_threshold,1.0-m_threshold,2.e-12,4.0*2.e-16,100);
        }
      else
         w = m_threshold;
      out[0] = exp(1. - pow(v1 * (pow((1.0 - log(w)),m_theta) - 1.) + 1.,(1. / m_theta)));
      out[1] = exp(1. - pow((1.0-v1) * (pow((1.0 - log(w)),m_theta) - 1.) + 1.,(1. / m_theta)));
      return out;
     }
protected:
   virtual double    pdf(double u,double v) override
     {
      double Cuv,u_part,v_part,numerator,denominator;
      Cuv = cdf(u,v);
      u_part = pow(1.0 - log(u),m_theta);
      v_part = pow(1.0 - log(v),m_theta);
      numerator = Cuv * u_part *v_part * (-1.0+m_theta + pow(-1.0+u_part +v_part,1.0/m_theta))*pow(-1.0+u_part+v_part,1.0/m_theta);
      denominator = u*v*(1.0-1.0*log(u))*(1.0 - log(v))*pow(-1.0+u_part+v_part,2.0);
      return (numerator+1.e-8)/(denominator+1.e-8);
     }
   virtual double    cdf(double u, double v) override
     {
      double u_part = pow(1.0-1.0*log(u),m_theta);
      double v_part = pow(1.0-1.0*log(v),m_theta);
      double cdf = exp(1.0 - pow(-1.0+u_part+v_part,1.0/m_theta));
      return cdf;
     }
   virtual double    condi_cdf(double u, double v) override
     {
      if(!m_theta)
         return EMPTY_VALUE;

      double u_part = pow(1.0 - log(u),m_theta);
      double v_part = pow(1.0 - log(v),m_theta);
      double cuv = cdf(u,v);

      double numerator = cuv *pow(-1.0+u_part+v_part,1.0/m_theta)*v_part;
      double denominator = v*(-1.0+u_part+v_part)*(1.0-1.0*log(v));

      return (numerator+1.e-8)/(denominator+1.e-8);
     }
   virtual vector    pdf(vector &u,vector &v) override
     {
      vector Cuv,u_part,v_part,numerator,denominator;
      Cuv = cdf(u,v);
      u_part = pow(1.0 - log(u),m_theta);
      v_part = pow(1.0 - log(v),m_theta);
      numerator = Cuv * u_part *v_part * (-1.0+m_theta + pow(-1.0+u_part +v_part,1.0/m_theta))*pow(-1.0+u_part+v_part,1.0/m_theta);
      denominator = u*v*(1.0-1.0*log(u))*(1.0 - log(v))*pow(-1.0+u_part+v_part,2.0);
      return (numerator+1.e-8)/(denominator+1.e-8);
     }
   virtual vector    cdf(vector &u, vector &v) override
     {
      vector u_part = pow(1.0-1.0*log(u),m_theta);
      vector v_part = pow(1.0-1.0*log(v),m_theta);
      vector cdf = exp(1.0 - pow(-1.0+u_part+v_part,1.0/m_theta));
      return cdf;
     }
   virtual vector    condi_cdf(vector &u, vector &v) override
     {

      vector u_part = pow(1.0 - log(u),m_theta);
      vector v_part = pow(1.0 - log(v),m_theta);
      vector cuv = cdf(u,v);

      vector numerator = cuv *pow(-1.0+u_part+v_part,1.0/m_theta)*v_part;
      vector denominator = v*(-1.0+u_part+v_part)*(1.0-1.0*log(v));

      return (numerator+1.e-8)/(denominator+1.e-8);
     }
public:
                     CN13(void)
     {
      m_threshold = 1.e-10;
      m_bounds[0] = 0.0;
      m_bounds[1] = DBL_MAX;
      m_copula_type = N13_COPULA;
     }
                    ~CN13(void)
     {
     }
   virtual matrix    Sample(ulong num_samples) override
     {
      double unf_v[],unf_c[];

      vector v,c,u;

      if(!MathRandomUniform(0.0,1.0,int(num_samples),unf_v) ||
         !MathRandomUniform(0.0,1.0,int(num_samples),unf_c) ||
         !v.Assign(unf_v) || !c.Assign(unf_c))
        {
         Print(__FUNCTION__, " failed to get uniform random numbers ", GetLastError());
         return matrix::Zeros(0,0);
        }

      matrix out(num_samples,2);

      for(ulong irow = 0; irow<num_samples; irow++)
         out.Row(generate_pair(v[irow],c[irow]),irow);

      return out;

     }

  };
//+------------------------------------------------------------------+

下図は、N13コピュラからサンプリングしたデータセットの可視化です。

二変量N13サンプルデータセット


二変量N14コピュラ

N14コピュラは、もう1つのアルキメデス型1パラメータ族で正の依存関係のみをモデル化しますが、N13とは異なり、上側裾と下側裾の両方への依存性を持つことが可能です。上側裾依存性は明確に現れ、下側側裾依存性はパラメータθが独立から離れるにつれて通常は弱く現れます。この特性により、極端に高い値の同時発生や、ある程度は極端に低い値の同時発生を捉えることができます。依存の強さはθの値が大きくなるほど増加し、下限では独立コピュラに近づきます。金融ペアの実データ分析では、裾の非対称により、N14コピュラが現実データに適合しやすいことが報告されています。N14コピュラには以下のジェネレータ関数があります。

n14ジェネレータ関数

MQL5では、n14.mqh内でCN14クラスとして実装されています。 

//+------------------------------------------------------------------+
//|                                                          n14.mqh |
//|                                  Copyright 2025, MetaQuotes Ltd. |
//|                                             https://www.mql5.com |
//+------------------------------------------------------------------+
#property copyright "Copyright 2025, MetaQuotes Ltd."
#property link      "https://www.mql5.com"
#include "base.mqh"
//+------------------------------------------------------------------+
//| N14 Copula (Nelsen 14)                                           |
//+------------------------------------------------------------------+
class CN14:public CBivariateCopula
  {
private:
   class CBrent:public CBrentQ
     {
   public:
                     CBrent(void)
        {
        }
                    ~CBrent(void)
        {
        }
      virtual double objective(double u, double v,double y) override
        {
         return ((-1.0*u) * (-2.0 + pow(u,1.0/v))) - y;
        }
     };
   virtual double    theta_hat(const double tau) override
     {
      return (1.0 + tau)/(2.0 - 2.0*tau);
     }
   vector            generate_pair(double v1, double v2)
     {
      vector out(2);
      double w = 0.0;
      if(v2>m_threshold)
        {
         CBrent brent;
         w = brent.minimize(m_theta,v2,m_threshold,1.0-m_threshold,2.e-12,4.0*2.e-16,100);
        }
      else
         w = m_threshold;
      out[0] = pow(1.0 + pow(v1 * pow(pow(w,-1.0/m_theta)-1.0,m_theta),1.0/m_theta),-1.0*m_theta);
      out[1] = pow(1.0 + pow((1.0-v1) * pow(pow(w,-1.0/m_theta)-1.0,m_theta),1.0/m_theta),-1.0*m_theta);
      return out;
     }
protected:
   virtual double    pdf(double u,double v) override
     {
      double u_ker,v_ker,u_part,v_part,cdf_ker,numerator,denominator;
      u_ker = -1.0+pow(u,1.0/m_theta);
      v_ker = -1.0+pow(v,1.0/m_theta);
      u_part = pow(-1.0 + pow(u,-1.0/m_theta),m_theta);
      v_part = pow(-1.0 + pow(v,-1.0/m_theta),m_theta);
      cdf_ker = 1.0 + pow(u_part+v_part,1.0/m_theta);
      numerator = (u_part * v_part *(cdf_ker - 1.0)*(-1.0 + m_theta + 2.0*m_theta*(cdf_ker-1.0)));
      denominator = pow(u_part + v_part,2.0) * pow(cdf_ker,(2.0 + m_theta))* u * v * u_ker * v_ker * m_theta;

      return (numerator+1.e-8)/(denominator+1.e-8);
     }
   virtual double    cdf(double u, double v) override
     {
      double u_part = pow(-1.0+pow(u,-1.0/m_theta),m_theta);
      double v_part = pow(-1.0+pow(v,-1.0/m_theta),m_theta);
      double cdf = pow(1.0 + pow(u_part+v_part,1.0/m_theta),-1*m_theta);
      return cdf;
     }
   virtual double    condi_cdf(double u, double v) override
     {
      double v_ker = -1.0 + pow(v, -1.0 / m_theta);
      double u_part = pow(-1.0 + pow(u, -1.0 / m_theta),m_theta);
      double v_part = pow(-1.0 + pow(v, -1.0 / m_theta),m_theta);
      double cdf_ker = 1.0 + pow(u_part + v_part, (1.0/ m_theta));

      double numerator = v_part * (cdf_ker - 1.0);
      double denominator = pow(v, (1.0 + 1.0 / m_theta)) * v_ker * (u_part + v_part) * pow(cdf_ker,1.0 + m_theta);

      return (numerator+1.e-8)/(denominator+1.e-8);
     }
   virtual vector    pdf(vector&u,vector&v) override
     {
      vector u_ker,v_ker,u_part,v_part,cdf_ker,numerator,denominator;
      u_ker = -1.0+pow(u,1.0/m_theta);
      v_ker = -1.0+pow(v,1.0/m_theta);
      u_part = pow(-1.0 + pow(u,-1.0/m_theta),m_theta);
      v_part = pow(-1.0 + pow(v,-1.0/m_theta),m_theta);
      cdf_ker = 1.0 + pow(u_part+v_part,1.0/m_theta);
      numerator = (u_part * v_part *(cdf_ker - 1.0)*(-1.0 + m_theta + 2.0*m_theta*(cdf_ker-1.0)));
      denominator = pow(u_part + v_part,2.0) * pow(cdf_ker,(2.0 + m_theta))* u * v * u_ker * v_ker * m_theta;

      return (numerator+1.e-8)/(denominator+1.e-8);
     }
   virtual vector    cdf(vector&u, vector&v) override
     {
      vector u_part = pow(-1.0+pow(u,-1.0/m_theta),m_theta);
      vector v_part = pow(-1.0+pow(v,-1.0/m_theta),m_theta);
      vector cdf = pow(1.0 + pow(u_part+v_part,1.0/m_theta),-1*m_theta);
      return cdf;
     }
   virtual vector    condi_cdf(vector&u, vector&v) override
     {
      vector v_ker = -1.0 + pow(v, -1.0 / m_theta);
      vector u_part = pow(-1.0 + pow(u, -1.0 / m_theta),m_theta);
      vector v_part = pow(-1.0 + pow(v, -1.0 / m_theta),m_theta);
      vector cdf_ker = 1.0 + pow(u_part + v_part, (1.0/ m_theta));

      vector numerator = v_part * (cdf_ker - 1.0);
      vector denominator = pow(v, (1.0 + 1.0 / m_theta)) * v_ker * (u_part + v_part) * pow(cdf_ker,1.0 + m_theta);

      return (numerator+1.e-8)/(denominator+1.e-8);
     }
public:
                     CN14(void)
     {
      m_threshold = 1.e-10;
      m_bounds[0] = 1.0;
      m_bounds[1] = DBL_MAX;
      m_copula_type = N14_COPULA;
     }
                    ~CN14(void)
     {
     }
   virtual matrix    Sample(ulong num_samples) override
     {
      double unf_v[],unf_c[];

      vector v,c,u;

      if(!MathRandomUniform(0.0,1.0,int(num_samples),unf_v) ||
         !MathRandomUniform(0.0,1.0,int(num_samples),unf_c) ||
         !v.Assign(unf_v) || !c.Assign(unf_c))
        {
         Print(__FUNCTION__, " failed to get uniform random numbers ", GetLastError());
         return matrix::Zeros(0,0);
        }

      matrix out(num_samples,2);

      for(ulong irow = 0; irow<num_samples; irow++)
         out.Row(generate_pair(v[irow],c[irow]),irow);

      return out;

     }
  };
//+------------------------------------------------------------------+

以下の散布図は、N14コピュラからサンプリングされたデータセットを示しています。

二変量N14サンプルデータセット

ここまでで実装したコピュラを数えると、二変量データセットにおけるさまざまな依存関係を捉えられる十分な数が揃いました。次のセクションでは、コピュラによって捉えられた依存関係モデルから何が推測できるかを解説します。これにより、コピュラをどのように応用すべきかをより深く理解できます。


コピュラ演算の解釈

コピュラは変数間の依存関係を捉えることが知られています。コピュラのパラメータは依存の強さを表し、これは直感的に理解しやすいです。さらに、コピュラは確率論に基づいているため、この依存構造を通して、各変数の特性も推測できます。これはコピュラの累積分布関数(CDF)、確率密度関数(PDF)、および条件付きCDFを評価することでおこなわれます。このような推論は、生の変数の分位数またはパーセンタイルを表す変数に対しておこなわれます。推論の対象となる変数は、生データの分位数や百分位数で表されます。分位数や百分位数は、ソートされたデータセットを同じ大きさの区間に分割する統計的指標であり、特定のデータ点の分布内での位置や相対的な立ち位置を示すのに役立ちます。

分位数uとvが与えられた場合、コピュラCDF C(u,v)は結合確率P(U≤u,V≤v)を表します。これは、「1つの変数がu分位数以下、もう1つがv分位数以下になる確率」を、両変数間の依存関係を考慮して定量化します。たとえば、u=0.8、v=0.9の場合、1つの変数が分布の下位80%に入り、もう1つの変数が下位90%に入る確率を示します。MQL5コード上では、コピュラオブジェクトのCopula_CDF()メソッドで計算できます。

CDFと補完的に、コピュラのPDF、c(u,v)があります。MQL5ではCopula_PDF()メソッドで評価します。CDFがある点までの累積確率を示すのに対し、PDFはその点における確率密度の集中度を示します。PDFが高い場合、その組み合わせ(u,v)が生じやすくなります。具体的には、コピュラPDFは次の式で表されます。

コピュラPDF

PDFの値が高い場合、それは独立の場合と比べて変数間の局所的依存が強いことを意味します。一方、PDFが低い、あるいはほぼゼロの場合は、その組み合わせが同時に起こる可能性が低いことを示します。

条件付き確率は、コピュラのCDFを一方の変数で偏微分することで求められます。通常、C(v∣u)やC(u∣v)と表され、一方の分位数が固定されたときの他方の条件付き挙動を捉えます。リスク管理の文脈では、たとえば「資産Aのリターンが10パーセンタイル(u=0.10)のとき、資産Bのリターンが5パーセンタイル以下(V≤0.05)になる確率はどれくらいか？」のような質問に使われます。 MQL5では、これをConditional_Probability()メソッドで実装します。第一引数は固定する変数、第二引数は条件付き分位数を表します。さらにInv_Conditional_Probability()メソッドを使うと、条件付きCDFの逆を求めることができます。条件付き確率yと条件付き分位数vが与えられたとき、C(v∣u)=yを満たすuを計算します。次のセクションでは、このコピュラ推論をペアトレード戦略の設計にどのように応用するかを見ていきます。


コピュラを用いたペアトレード戦略

距離法や共和分といった従来のペアトレード手法では、資産価格間の非線形または非対称な関係を考慮できず、取引スプレッドが正規分布に従い線形関係であると仮定しているため、時としてうまく機能しないことがあります。理論上、コピュラを用いた手法はこの問題に対する解決策となり得ます。コピュラは非正規性や裾依存性を捉える柔軟なモデリングを可能にするためです。このフレームワークでは、コピュラを用いて条件付き確率を計算し、過小評価や過大評価のイベントを定義することで取引シグナルを生成します。一方の資産が過小評価されているとみなされるのは、ペアのリターンが与えられた場合に低リターンの条件付き確率が非常に高いときです。逆に、過大評価されているとみなされるのは、その条件付き確率が非常に低いときです。つまり、取引のエントリーは線形の閾値ではなく、結合分布から導かれる分位数ベースの信頼水準に基づいて判断されます。

コピュラシグナル

資産価格の依存関係のモデリングにおけるコピュラの理論上の適合性が示唆されている一方で、その限界にも注意する必要があります。まず第一に、コピュラは定常性を持つ多変量データセットに適用した場合に最も効果的です。データが非定常である場合、周辺分布やコピュラ関数の推定パラメータは、推定対象となった特定期間にのみ有効な関係を捉えることになり、将来の期間において性質が変化している場合には、推論が信頼できず、予測精度も低下します。非定常時系列は、単に共通のトレンドによって高い相関や依存関係を示すことがよくあります。もし周辺分布が定常でない場合、コピュラモデルは変化する周辺特性の影響を捉えてしまい、変数間の本質的で時間不変な関係ではなく、単なる人工的な依存構造を捉えてしまう可能性があります。

非定常時系列を扱う際、確率積分変換によって一様分布に変換しても、必ずしも独立同分布(IID)にはならず、これは通常コピュラ推定の入力条件となります。定常性があることで、周辺変数の挙動と依存構造を明確に分離でき、これがコピュラアプローチを採用する主な利点です。周辺分布が非定常の場合、推定されたコピュラは周辺の非定常ダイナミクスを暗黙的に捉えてしまい、両者の区別が曖昧になります。それにもかかわらず、学術文献にはコピュラモデルに基づくペアトレード戦略の事例が数多く報告されています。

なお、定常性の要件は標準的または静的なコピュラモデルに対してのものです。データが本質的に非定常である場合には、専門的な手法が存在します。具体的には、データを事前処理して定常系列に変換した後にコピュラを適用する方法や、時間変動型個ぴゅらや非定常コピュラのように依存パラメータが時間や他の共変量に応じて変化するモデルがあります。これらのモデルは非定常データに対応するために設計されています。この種のコピュラについては別記事で取り上げる予定です。現時点では、標準的または静的なコピュラを用いて、生の資産価格ペアの依存関係をモデル化する有効性を検証します。このタイプの戦略を効果的に適用するためには、まず取引対象として適した銘柄を選定する方法を定義する必要があります。


取引ペアの選定

コピュラを用いたペアトレードにおいて銘柄を選定する際には、従来のペアトレード戦略で用いられる方法も応用可能です。代表的な方法として共和分法と距離法があります。共和分法は一般的に最も厳密な方法とされています。目的は、2つの非定常価格系列の線形結合が定常となるペアを見つけることです。価格系列に対してはEngle-Granger二段階法やジョハンセン検定が用いられます。定常スプレッドを持つペアは、価格が長期的に従う安定した均衡関係を示します。距離法（残差平方和、SSD）は、よりシンプルで非パラメトリックなアプローチです。

残差平方和

正規化された価格系列または累積リターン間の歴史的な残差平方和(SSD)を最小化するペアを選定します。SSDが最も小さいペアは、過去に最も密接に連動して動いたペアと見なされます。

テキストでは、銘柄ペアの選定に一致度指標を使用します。一致度指標は非パラメトリックな手法で、資産価格系列間の依存の強さや性質を定量化するものです。これはペアトレード戦略の核心となる原則です。代表的な非パラメトリックの一致度指標は以下のとりです。

  • ケンドールのタウ：一致ペアと不一致ペアの出現確率の差に基づく指標で、依存の強さを評価する頑健な指標です。
  • スピアマンの順位相関係数：本質的に、元データではなく、ランク化したデータに基づき計算したピアソン相関係数です。単調関係（対応する変数が同方向または逆方向に動く傾向）の強さを測定します。

一致度が高いほど、銘柄間の共動性が強いことを示します。ケンドールのタウやスピアマンの順位相関係数の値が最も高いペアを選定します。値が1に近い場合は、銘柄価格が一緒に動く傾向が強く、ペアトレードで長期的な関係を活用できる可能性が高いことを示します。MetaTrader 5では、スクリプトConcordance.ex5が銘柄リストを入力として受け取り、最も一致度の高いペアを計算して選定します。

//+------------------------------------------------------------------+
//|                                                  Concordance.mq5 |
//|                                  Copyright 2025, MetaQuotes Ltd. |
//|                                             https://www.mql5.com |
//+------------------------------------------------------------------+
#property copyright "Copyright 2025, MetaQuotes Ltd."
#property link      "https://www.mql5.com"
#property version   "1.00"
#property script_show_inputs
#include<dependence.mqh>
input string   input_symbols = "AUDUSD,NZDUSD,XAUUSD,XAGUSD,EURUSD,XAUEUR";
input datetime input_start_date = D'2025.01.01 00:00';
input ulong    input_count_bars = 260;
input ENUM_TIMEFRAMES input_time_frame = PERIOD_D1;
input ENUM_DEPENDENCE_MEASURE dependence_measure = KENDAL_TAU;
input bool Use_Returns = false;
//+------------------------------------------------------------------+
//| Script program start function                                    |
//+------------------------------------------------------------------+
void OnStart()
  {
//---
   string test_symbols[];
   int num_symbols = StringSplit(input_symbols,StringGetCharacter(",",0),test_symbols);
   if(StringFind(input_symbols,",",StringLen(input_symbols)-1) == StringLen(input_symbols)-1)
     {
      --num_symbols;
      ArrayRemove(test_symbols,uint(num_symbols),1);
     }

   if(num_symbols<2)
     {
      Print(" 2 or more symbols expected ");
      return;
     }

   matrix prices = matrix::Zeros(input_count_bars,ulong(num_symbols));
   vector dwnload;
   for(ulong i = 0; i<prices.Cols(); ++i)
     {
      int try = 10;
      ResetLastError();
      while((!dwnload.CopyRates(test_symbols[i],input_time_frame,COPY_RATES_CLOSE,input_start_date,input_count_bars) ||
         !prices.Col(dwnload,i)) && try)
       {
         Sleep(10000);
         --try;
       }
      /*if(StringFind(test_symbols[i],"XAUEUR")>=0)
        {
         prices.Col(dwnload*prices.Col(i-1),i);
        }*/
      if(GetLastError())
       {
        Print("Could not download data");
        return;
       }
     }
     
   if(Use_Returns)
     {
      prices = log(prices);
      prices = np::diff(prices,1,false);
     }

   matrix dep_mat = dependence(prices,dependence_measure);
   Print(dep_mat);
   Print(EnumToString(dependence_measure), "\n","Using ", Use_Returns?"Returns":"Raw prices");
   if(test_symbols.Size()>2)
     {
      matrix mm = dep_mat.TriL(-1);
      int index = int(mm.ArgMax());
      Print("Pair with highest dependence: ", test_symbols[int(index/int(dep_mat.Cols()))], " and ", test_symbols[int(MathMod(index,dep_mat.Cols()))]);
     }
  }
//+------------------------------------------------------------------+

下図は、銘柄群AUDUSD、NZDUSD、XAUUSD、XAGUSD、EURUSD、XAUEURを対象にConcordance.ex5スクリプトを実行した際の出力例です。

Concordance.ex5スクリプト出力

統計的な測定値に加えて、選定したペアの依存関係は経済的なファンダメンタルズに裏付けられていることが望ましいです。Concordance.ex5スクリプトのテスト結果を見ると、XAUEURとXAUUSDが最も高い共動性を示しました。このペアの依存関係は、両方とも同じ基礎資産に由来しているため、容易に説明できます。これにより、取引に適した銘柄ペアが確定しました。次のステップは、この銘柄ペアの依存特性に最も適合するコピュラモデルを選定することです。


コピュラモデルの選定

適切なコピュラモデルを選定する際の主な方法には、モデル選択基準と適合度(Goodness-of-Fit, GoF)検定があります。モデル選択基準（情報量基準）は、モデルの適合度と複雑さ（パラメータ数）のバランスを評価します。一般的に、値が低いほど、より適切なモデルと判断されます。代表的なものには、赤池情報量規準(AIC)、ベイズ情報量基準(BIC)、Hannan-Quinn情報量基準(HQIC)があります。  手順は一般的に3段階です。まず、複数の候補コピュラモデルのパラメータを推定し、各モデルの対数尤度を計算します。次に、AIC、BIC、HQICのいずれかを用いて、各候補モデルの情報量基準を算出します。最後に、すべての候補モデルの情報量基準を比較し、最も低い値を持つモデルを選定します。

別の方法として、特定のコピュラモデルが観測データに適合しているかを統計的に判定する適合度検定があります。これらの検定は、特定のコピュラモデルによって捉えられた依存構造が、観測データと整合しているかどうかを統計的に判定するものです。これらは正式な仮説検定であり、帰無仮説は「データは仮定されたコピュラ族に従って生成されている」となります。使用できる検定統計量の例には以下があります。

  • Cramér-von Mises統計量：経験的コピュラとパラメトリックコピュラの二乗差を測定する
  • コルモゴロフ–スミルノフ統計量：経験的コピュラとパラメトリックコピュラの絶対差の最大値を測定する

多くの場合、パラメトリックブートストラップを用いて帰無分布をシミュレーションし、p値を算出します。p値が高い場合、帰無仮説は棄却されず、モデルはデータに対して妥当と判断されます。現状、MQL5にはコピュラ用のコルモゴロフ–スミルノフやCramér-von Misesのネイティブ実装は存在しません。そのため、モデル選択基準による方法に依存する必要があります。

//+------------------------------------------------------------------+
//|                                              CopulaSelection.mq5 |
//|                                  Copyright 2025, MetaQuotes Ltd. |
//|                                             https://www.mql5.com |
//+------------------------------------------------------------------+
#property copyright "Copyright 2025, MetaQuotes Ltd."
#property link      "https://www.mql5.com"
#property version   "1.00"
#property script_show_inputs
#include <Copulas\Bivariate\bivariate_copula.mqh>
#include <ECDF\linear_cdf.mqh>
//--- input parameters
input string   FirstSymbol="XAUUSD";
input string   SecondSymbol="XAUEUR";
input datetime TrainingDataStart=D'2025.01.01';
input ulong    HistorySize=260;
input bool     SaveModel = false;
input string   EcdfModel = "model.ecdf";
input string   CopulaModel = "model.copula";
//---
//string   NormalizingSymbol="EURUSD";
//+------------------------------------------------------------------+
//| Script program start function                                    |
//+------------------------------------------------------------------+
void OnStart()
  {
//---
   vector p_1,p_2,p_3,p_n;
   matrix pdata,pobs;
   if(!p_1.CopyRates(FirstSymbol,PERIOD_CURRENT,COPY_RATES_CLOSE,TrainingDataStart,HistorySize) ||
      !p_2.CopyRates(SecondSymbol,PERIOD_CURRENT,COPY_RATES_CLOSE,TrainingDataStart,HistorySize) ||
      p_1.Size()!=p_2.Size() ||
      !pdata.Resize(p_1.Size(),2) ||
      !pdata.Col(p_1,0) ||
      !pdata.Col(p_2,1))
     {
      Print(" failed to collect and initialize rates matrix ", GetLastError());
      return;
     }
//---
   CLinearCDF qt();
   if(!qt.fit(pdata))
      return;
//---
   pobs = qt.to_quantile(pdata);
//---
   if(SaveModel)
     {
      CFileBin file;
      file.Open(EcdfModel,FILE_WRITE|FILE_BIN|FILE_COMMON);
      if(!qt.save(file.Handle()))
         Print(" Failed to save ", EcdfModel);
      file.Close();
     }
//---
   vector lowest = vector::Zeros(8);
   CBivariateCopula *bcop[8];
//---
   bcop[0] = new CClayton();
   bcop[1] = new CFrank();
   bcop[2] = new CGumbel();
   bcop[3] = new CJoe();
   bcop[4] = new CN13();
   bcop[5] = new CN14();
   bcop[6] = new CGaussian();
   bcop[7] = new CStudent();
//---
   for(uint i = 0; i<bcop.Size(); ++i)
     {
      bool fitted = bcop[i].Fit(pobs.Col(0),pobs.Col(1));
      //---
      if(fitted)
        {
         double ll = bcop[i].Log_likelihood(pobs.Col(0),pobs.Col(1));
         lowest[i] = aic(ll,int(pobs.Rows()));
         Print(EnumToString((ENUM_COPULA_TYPE)bcop[i].Type()));
         Print(" sic ", sic(ll,int(pobs.Rows())));
         Print(" aic ", lowest[i]);
         Print(" hqic ", hqic(ll,int(pobs.Rows())));
        }
     }
//---
   if(SaveModel)
     {
      ulong shift = lowest.ArgMin();
      CFileBin file;
      file.Open(CopulaModel,FILE_WRITE|FILE_BIN|FILE_COMMON);
      if(!bcop[shift].Save(file.Handle()))
         Print("Failed to save ", CopulaModel);
      file.Close();
     }
//---
   for(uint i = 0; i<bcop.Size(); delete bcop[i], ++i);
  }
//+------------------------------------------------------------------+

スクリプトCopulaSelection.ex5では、ユーザーが銘柄ペア、開始日、そしてサンプルとして使用するバーの数を指定でき、これらのサンプルを基に候補となるコピュラモデルを構築できます。このスクリプトは、各候補コピュラモデルの対数尤度を計算した後、情報量基準を算出します。このスクリプトでは、選択されたコピュラモデルを、対応する経験累積分布関数(ECDF)モデルとともに保存することができます。結果はMetaTrader 5の[操作ログ]タブに出力され、ユーザーが確認できるようになっています。以下は、XAUUSDとXAUEURの銘柄ペアにスクリプトを適用した際の出力例です。

CopulaSelectionスクリプト出力

スクリプトの出力から、フランクコピュラが最も情報量基準の値が低いことが明確に示されており、サンプルデータに対して最適なモデルであることが示唆されます。保存されたモデルは、コピュラを用いたペアトレード戦略を実装するインジケーターやエキスパートアドバイザー(EA)で使用されます。


コピュラを用いたペアトレード戦略の例

ここで実装する戦略は、Liew, R.Q.とWu, Y.によるJournal of Derivatives & Hedge Fundsの論文 の手法を基に、HudsonとThamesのアプローチに沿って修正版として構築したものです。戦略の核は、サンプル価格データから学習したコピュラモデルによりペア間の依存関係を捉え、現在の価格の変位値から計算された条件付き確率を用いて売買シグナルを生成する点にあります。エントリーロジックは以下のとおりです。

  • 「C(u_1 | u_2) <= 0.05」かつ「C(u_2 | u_1) >= 0.95」の場合、銘柄1が過小評価され、銘柄2が過大評価されていると判断し、銘柄1に買いシグナル、銘柄2に売りシグナルを同時に発生させる
  • 「C(u_1 | u_2) >= 0.95」かつ「C(u_2 | u_1) <= 0.05」の場合、銘柄2が過小評価され、銘柄1が過大評価されていると判断し、銘柄2に買いシグナル、銘柄1に売りシグナルを同時に発生させる

変数u_1とu_2は、それぞれ銘柄1と銘柄2の分位数を表します。ポジションの決済方法は以下の2つがあります。

  • いずれかの条件付き確率が上限のエグジット閾値を超える、または下限の閾値を下回った場合に両方のポジションを決済する
  • 片方の条件付き確率が上限の閾値を超え、同時にもう片方の条件付き確率が下限の閾値を下回った場合に両方のポジションを決済する

この戦略は純粋に検証目的で実装しています。理論的には、価格が学習データの最大値や最小値を超えると、モデルはすぐに破綻し意味をなさなくなることは明らかです。実装の目的は、コピュラが実際に取引機会を示唆できるかを評価することです。そのため、コピュラモデルを学習した期間と同じ期間で戦略をテストします（インサンプルテスト）。そのため、GoldCopulaSignals.ex5というインジケーターが作成されました。このインジケーターは、CopulaSelection.ex5スクリプトで保存されたECDFおよびコピュラモデルを読み込み、XAUEURとXAUUSDの分位数から計算された条件付き確率を表示します。

//+------------------------------------------------------------------+
//|                                            GoldCopulaSignals.mq5 |
//|                                  Copyright 2025, MetaQuotes Ltd. |
//|                                             https://www.mql5.com |
//+------------------------------------------------------------------+
#property copyright "Copyright 2025, MetaQuotes Ltd."
#property link      "https://www.mql5.com"
#property version   "1.00"
#property indicator_separate_window
#property indicator_buffers 4
#property indicator_plots   4
#property indicator_maximum 1.1
#property indicator_minimum -0.1
#property indicator_level1  0.05
#property indicator_level2  0.5
#property indicator_level3  0.95
#include<Copulas\Bivariate\bivariate_copula.mqh>
#include<ECDF\linear_cdf.mqh>
//--- plot S1
#property indicator_label1  "S1"
#property indicator_type1   DRAW_LINE
#property indicator_color1  clrRed
#property indicator_style1  STYLE_SOLID
#property indicator_width1  1
//--- plot S2
#property indicator_label2  "S2"
#property indicator_type2   DRAW_LINE
#property indicator_color2  clrBlue
#property indicator_style2  STYLE_SOLID
#property indicator_width2  1
//---
#property indicator_label3  "EURGold"
#property indicator_type3   DRAW_NONE

//---
#property indicator_label4  "USDGold"
#property indicator_type4  DRAW_NONE

//--- input parameters
input ENUM_COPULA_TYPE Copulatype = FRANK_COPULA;
input string   CopulaModelName = "model.copula";
input string   ECDFModelName = "model.ecdf";
//---
string   FirstSymbol="XAUUSD";
string   SecondSymbol="XAUEUR";
//--- indicator buffers
double         S1Buffer[];
double         S2Buffer[];
double         S3Buffer[];
double         S4Buffer[];
matrix         mat_;
CLinearCDF* qt;
CBivariateCopula *bcop;
CFileBin savedcop,savedcdf;
//+------------------------------------------------------------------+
//| Custom indicator initialization function                         |
//+------------------------------------------------------------------+
int OnInit()
  {

//---
   qt = new CLinearCDF();
//---
   switch(Copulatype)
     {
      case CLAYTON_COPULA:
         bcop = new CClayton();
         break;
      case FRANK_COPULA:
         bcop = new CFrank();
         break;
      case GUMBEL_COPULA:
         bcop = new CGumbel();
         break;
      case JOE_COPULA:
         bcop = new CJoe();
         break;
      case N13_COPULA:
         bcop = new CN13();
         break;
      case N14_COPULA:
         bcop = new CN14();
         break;
      case GAUSSIAN_COPULA:
         bcop = new CGaussian();
         break;
      case STUDENT_COPULA:
         bcop = new CStudent();
         break;
     }
//---
   savedcop.Open(CopulaModelName,FILE_READ|FILE_BIN|FILE_COMMON);
   if(!bcop.Load(savedcop.Handle()))
     {
      Print("Failed to load copula model ", CopulaModelName, " ", GetLastError());
      return INIT_FAILED;
     }
//---
   if((ENUM_COPULA_TYPE)bcop.Type() != Copulatype)
     {
      Print("Invalid copula model. Instantiated copula is ", EnumToString((ENUM_COPULA_TYPE)bcop.Type()),"\n but specified copula in input settings is ", EnumToString(Copulatype));
      return INIT_FAILED;
     }
//---
   savedcop.Close();
//---
   savedcdf.Open(ECDFModelName,FILE_READ|FILE_BIN|FILE_COMMON);
   if(!qt.load(savedcdf.Handle()))
     {
      Print("Failed to load the ECDF model ", ECDFModelName, " ", GetLastError());
      return INIT_FAILED;
     }
//---
   savedcdf.Close();
//---
   mat_.Resize(1,2);
//--- indicator buffers mapping
   SetIndexBuffer(0,S1Buffer,INDICATOR_DATA);
   SetIndexBuffer(1,S2Buffer,INDICATOR_DATA);
   SetIndexBuffer(2,S3Buffer,INDICATOR_DATA);
   SetIndexBuffer(3,S4Buffer,INDICATOR_DATA);
//---
   ArraySetAsSeries(S1Buffer,true);
   ArraySetAsSeries(S2Buffer,true);
   ArraySetAsSeries(S3Buffer,true);
   ArraySetAsSeries(S4Buffer,true);
//---
   PlotIndexSetString(0,PLOT_LABEL,FirstSymbol);
   PlotIndexSetString(1,PLOT_LABEL,SecondSymbol);
//---
   return(INIT_SUCCEEDED);
  }
//+------------------------------------------------------------------+
//| On deinitialization                                              |
//+------------------------------------------------------------------+
void  OnDeinit(const int  reason)
  {
//---
   if(CheckPointer(qt) == POINTER_DYNAMIC)
      delete qt;
//---
   if(CheckPointer(bcop) == POINTER_DYNAMIC)
      delete bcop;
  }
//+------------------------------------------------------------------+
//| Custom indicator iteration function                              |
//+------------------------------------------------------------------+
int OnCalculate(const int32_t rates_total,
                const int32_t prev_calculated,
                const datetime &time[],
                const double &open[],
                const double &high[],
                const double &low[],
                const double &close[],
                const long &tick_volume[],
                const long &volume[],
                const int32_t &spread[])
  {
//---
   int32_t limit;
   if(prev_calculated<=0)
      limit= rates_total-2;
   else
      limit=rates_total-prev_calculated;
//---
   for(int32_t i =limit; i >=0 ; --i)
     {
      mat_[0,0] = iClose(FirstSymbol,PERIOD_CURRENT,i);
      mat_[0,1] = iClose(SecondSymbol,PERIOD_CURRENT,i);
      S3Buffer[i] = mat_[0,1];
      S4Buffer[i] = mat_[0,0];
      mat_ = qt.to_quantile(mat_);
      S1Buffer[i] = bcop.Conditional_Probability(mat_[0][0],mat_[0][1]);
      S2Buffer[i] = bcop.Conditional_Probability(mat_[0][1],mat_[0][0]);
     }
//--- return value of prev_calculated for next call
   return(rates_total);
  }
//+------------------------------------------------------------------+

以下のスクリーンショットは、インジケーターの表示（サブウィンドウ）を示しています。

GoldCopulaSignalsインジケーター

青い線はXAUEURの条件付き確率を示し、赤い線はXAUUSDの条件付き確率を示しています。戦略自体はSimpleCopulaPairsStrategy.ex5というEAとして実装されています。

//+------------------------------------------------------------------+
//|                                    SimpleCopulaPairsStrategy.mq5 |
//|                                  Copyright 2025, MetaQuotes Ltd. |
//|                                             https://www.mql5.com |
//+------------------------------------------------------------------+
#property copyright "Copyright 2025, MetaQuotes Ltd."
#property link      "https://www.mql5.com"
#property version   "1.00"
#include <Trade\Trade.mqh>
#include <Trade\SymbolInfo.mqh>
#include <Trade\AccountInfo.mqh>
#include <TimeOptions.mqh>
#include<Copulas\Bivariate\base.mqh>
#resource "\\Indicators\\GoldCopulaSignals.ex5"
//---
#define USDGOLD "XAUUSD"
#define EURGOLD "XAUEUR"
//--- input parameters
enum ENUM_EXIT_RULE
  {
   AND_EXIT=0,//"And exit"
   OR_EXIT//Or exit
  };
//--- input parameters
input ENUM_COPULA_TYPE Copulatype = FRANK_COPULA;
input string   CopulaModelName = "model.copula";
input string   ECDFModelName = "model.ecdf";
input double   OpenThresholdProbability_First = 0.05;
input double   OpenThresholdProbability_Second = 0.95;
input double   UpperExitThresholdProbability = 0.5;
input double   LowerExitThresholdProbability = 0.5;
input ENUM_EXIT_RULE ExitRule = OR_EXIT;
input ulong    SlippagePoints = 10;
input double   TradingLots = 0.01;
input ulong    MagicNumber = 18973984;
double   StopLossPoints = 0.0;
double   TakeProfitPoints = 0.0;
int  OpenShift = 1;
int  CloseShift = 1;
//---
ENUM_DAY_TIME_MINUTES TradeSessionOpen = 6;
ENUM_DAY_TIME_MINUTES TradeSessionClose = 1435;
//--- indicator buffers
CSymbolInfo usd,eur;
//--
double usd_sigs[2],eur_sigs[2];
int indicator_handle;
//+------------------------------------------------------------------+
//| Expert initialization function                                   |
//+------------------------------------------------------------------+
int OnInit()
  {
//---
   if(!usd.Name(USDGOLD) || !eur.Name(EURGOLD))
     {
      Print("Symbol info config error, ");
      return INIT_FAILED;
     }
//---
   indicator_handle = INVALID_HANDLE;
   int try
         = 10;
   while(indicator_handle == INVALID_HANDLE && try
            >0)
        {
         indicator_handle = iCustom(NULL,PERIOD_CURRENT,"::Indicators\\GoldCopulaSignals.ex5",Copulatype,CopulaModelName,ECDFModelName);
         --try;
        }
//---
   if(indicator_handle == INVALID_HANDLE)
     {
      Print(" failed to initialize the indicator ", GetLastError());
      return INIT_FAILED;
     }
//---
   return(INIT_SUCCEEDED);
  }
//+------------------------------------------------------------------+
//| Expert deinitialization function                                 |
//+------------------------------------------------------------------+
void OnDeinit(const int reason)
  {
//---
  }
//+------------------------------------------------------------------+
//| Expert tick function                                             |
//+------------------------------------------------------------------+
void OnTick()
  {
//---
   datetime sessionopen = iTime(NULL,PERIOD_D1,0) + int(TradeSessionOpen*60);
   datetime sessionclose = iTime(NULL,PERIOD_D1,0) + int(TradeSessionClose*60);
//---
   if(TimeCurrent()<sessionopen || TimeCurrent()>=sessionclose)
      return;
//---
   if(SumMarketOrders(MagicNumber))
     {
      long exit = GetExitSignal(CloseShift);
      if(exit)
        {
         if(CloseAll(MagicNumber,NULL,WRONG_VALUE,SlippagePoints)!=2)
            Print("Failed to close all opened orders ");
        }
     }
   else
     {
      long entry = GetEntrySignal(OpenShift);
      if(entry>0)
        {
         if(!OpenOrder(usd,ORDER_TYPE_BUY,StopLossPoints,TakeProfitPoints,TradingLots,SlippagePoints,MagicNumber))
            Print(" Failed long entry for ", usd.Name());
         if(!OpenOrder(eur,ORDER_TYPE_SELL,StopLossPoints,TakeProfitPoints,TradingLots,SlippagePoints,MagicNumber))
            Print(" Failed short entry for ", eur.Name());
        }
      else
         if(entry<0)
           {
            if(!OpenOrder(eur,ORDER_TYPE_BUY,StopLossPoints,TakeProfitPoints,TradingLots,SlippagePoints,MagicNumber))
               Print(" Failed long entry for ", eur.Name());
            if(!OpenOrder(usd,ORDER_TYPE_SELL,StopLossPoints,TakeProfitPoints,TradingLots,SlippagePoints,MagicNumber))
               Print(" Failed short entry for ", usd.Name());
           }
     }
  }
//+------------------------------------------------------------------+
//| get the entry signal                                             |
//+------------------------------------------------------------------+
long GetEntrySignal(int shift = 1)
  {
   static datetime last;
   if(shift>0)
     {
      if(iTime(NULL,PERIOD_CURRENT,0)<=last)
         return 0;
      else
         last = iTime(NULL,PERIOD_CURRENT,0);
     }
   else
     {
      datetime last = (datetime)MathMax(eur.Time(),usd.Time());
      if(TimeCurrent()<=last)
         return 0;
      //---
      eur.RefreshRates();
      usd.RefreshRates();
     }
//---
   if(CopyBuffer(indicator_handle,0,shift,2,usd_sigs)<2 ||
      CopyBuffer(indicator_handle,1,shift,2,eur_sigs)<2)
     {
      Print(__FUNCTION__, " error copying from indicator buffers ");
      return 0;
     }
   double current_prob_first,current_prob_second;
//---
   current_prob_first = usd_sigs[1];
   current_prob_second = eur_sigs[1];
//---
   if(current_prob_first<=OpenThresholdProbability_First && current_prob_second>=OpenThresholdProbability_Second)
     {
      //Print(__FUNCTION__, " go long ");
      return 1;
     }
   else
      if(current_prob_first>=OpenThresholdProbability_Second && current_prob_second<=OpenThresholdProbability_First)
        {
         //Print(__FUNCTION__, " go short ");
         return -1;
        }
//---
   return 0;
  }
//+------------------------------------------------------------------+
//| get the exit signal                                              |
//+------------------------------------------------------------------+
long GetExitSignal(int shift = 1)
  {
   datetime last = (datetime)MathMax(eur.Time(),usd.Time());

   if(shift>0)
     {
      if(iTime(NULL,PERIOD_CURRENT,0)<=last)
         return 0;
     }
   else
     {
      if(TimeCurrent()<=last)
         return 0;
     }
//---
   eur.RefreshRates();
   usd.RefreshRates();

//---
   if(CopyBuffer(indicator_handle,0,shift,2,usd_sigs)<2 ||
      CopyBuffer(indicator_handle,1,shift,2,eur_sigs)<2)
     {
      Print(__FUNCTION__, " error copying from indicator buffers ");
      return 0;
     }
//---
   double prev_prob_first,prev_prob_second;
   double current_prob_first,current_prob_second;
//---
   prev_prob_first = usd_sigs[0];
   prev_prob_second = eur_sigs[0];
//---
   current_prob_first = usd_sigs[1];
   current_prob_second = eur_sigs[1];
//---
   bool u1_up,u2_dwn,u2_up,u1_dwn;
   u1_up = (prev_prob_first<= LowerExitThresholdProbability && current_prob_first>=UpperExitThresholdProbability);
   u1_dwn = (prev_prob_first>=UpperExitThresholdProbability && current_prob_first<=LowerExitThresholdProbability);
   u2_up = (prev_prob_second<=LowerExitThresholdProbability && current_prob_second>=UpperExitThresholdProbability);
   u2_dwn = (prev_prob_second>=UpperExitThresholdProbability && current_prob_second<=LowerExitThresholdProbability);

   if(ExitRule == AND_EXIT)
     {
      if(u1_up && u2_dwn)
         return -1;
      else
         if(u2_up && u1_dwn)
            return 1;
     }
   else
     {
      if(u1_up || u2_dwn)
         return -1;
      else
         if(u2_up || u1_dwn)
            return 1;
     }
   return 0;
  }
//+------------------------------------------------------------------+
//|Enumerate all market orders                                       |
//+------------------------------------------------------------------+
uint SumMarketOrders(ulong magic_number=ULONG_MAX,const string sym=NULL,const ENUM_POSITION_TYPE ordertype=WRONG_VALUE)
  {
   CPositionInfo posinfo;
//---
   uint count=0;
//---
   for(int i = PositionsTotal()-1; i>-1;--i)
     {
      if(!posinfo.SelectByIndex(i))
         continue;
      if(magic_number<ULONG_MAX && posinfo.Magic()!=long(magic_number))
         continue;
      if(sym!=NULL && posinfo.Symbol()!=sym)
         continue;
      if(ordertype!=WRONG_VALUE && posinfo.PositionType()!=ordertype)
         continue;
      count++;
     }
//---
   return count;
  }
//+------------------------------------------------------------------+
//|Enumerate all market orders                                       |
//+------------------------------------------------------------------+
ENUM_POSITION_TYPE LastOrderType(ulong magic_number=ULONG_MAX)
  {
   CPositionInfo posinfo;
//---
   ENUM_POSITION_TYPE direction_type = WRONG_VALUE;
   string symb = NULL;
//---
   for(int i = PositionsTotal()-1; i>-1;--i)
     {
      if(!posinfo.SelectByIndex(i))
         continue;
      if(magic_number<ULONG_MAX && posinfo.Magic()!=long(magic_number))
         continue;
      symb = posinfo.Symbol();
      direction_type = posinfo.PositionType();
      //---
      switch(direction_type)
        {
         case POSITION_TYPE_BUY:
           {
            if(symb == EURGOLD)
               return POSITION_TYPE_SELL;
            else
               if(symb == USDGOLD)
                  return POSITION_TYPE_BUY;
           }
         case POSITION_TYPE_SELL:
           {
            if(symb == EURGOLD)
               return POSITION_TYPE_BUY;
            else
               if(symb == USDGOLD)
                  return POSITION_TYPE_SELL;
           }
         default:
            return WRONG_VALUE;
        }
     }
//---
   return WRONG_VALUE;
  }
//+------------------------------------------------------------------+
//| Profit                                                           |
//+------------------------------------------------------------------+
double Profit(ulong magic_number=ULONG_MAX)
  {
   CPositionInfo posinfo;
//---
   double profit = 0.0;
//---
   for(int i = PositionsTotal()-1; i>-1;--i)
     {
      if(!posinfo.SelectByIndex(i))
         continue;
      if(magic_number<ULONG_MAX && posinfo.Magic()!=long(magic_number))
         continue;
      profit += (posinfo.Profit() + posinfo.Commission() + posinfo.Swap());
     }
//---
   return profit;
  }
//+------------------------------------------------------------------+
//|Open trade                                                        |
//+------------------------------------------------------------------+
bool OpenOrder(CSymbolInfo &syminfo, ENUM_ORDER_TYPE type,double stoplosspoints=0.0, double takeprofitpoints = 0.0,double lotsize=0.1,ulong slippage=10,ulong magic=12345678,string tradecomment=NULL)
  {
   CTrade trade;
//---
   CAccountInfo accinfo;
//---
   trade.SetExpertMagicNumber(magic); // magic
   trade.SetMarginMode();
   trade.SetDeviationInPoints(slippage);
//---
   if(!trade.SetTypeFillingBySymbol(syminfo.Name()))
     {
      Print(__FUNCTION__," Unknown Type filling mode for ", syminfo.Name());
      return false;
     }
//---
   bool ans=false;
   syminfo.RefreshRates();
   string sy = syminfo.Name();
//---
   double price = (type == ORDER_TYPE_BUY)?syminfo.Ask():syminfo.Bid();
//---
   price = NormalizeDouble(MathAbs(price),syminfo.Digits());
//---
   bool result = false;
//---
   switch(type)
     {
      case ORDER_TYPE_BUY:
         if(accinfo.FreeMarginCheck(sy,type,lotsize,price)<0.0)
           {
            Print("Insufficient funds to open long order ("+sy+"). Free Margin = ", accinfo.FreeMargin());
            return false;
           }
         result = trade.Buy(lotsize,sy,price,(stoplosspoints)?NormalizeDouble(price - MathAbs(stoplosspoints*syminfo.Point()),syminfo.Digits()):0.0,(takeprofitpoints)?NormalizeDouble(price + MathAbs(takeprofitpoints*syminfo.Point()),syminfo.Digits()):0.0,tradecomment);
         break;
      case ORDER_TYPE_SELL:
         if(accinfo.FreeMarginCheck(sy,type,lotsize,price)<0.0)
           {
            Print("Insufficient funds to open short order ("+sy+"). Free Margin = ", accinfo.FreeMargin());
            return false;
           }
         result = trade.Sell(lotsize,sy,price,(stoplosspoints)?NormalizeDouble(price + MathAbs(stoplosspoints*syminfo.Point()),syminfo.Digits()):0.0,(takeprofitpoints)?NormalizeDouble(price - MathAbs(takeprofitpoints*syminfo.Point()),syminfo.Digits()):0.0,tradecomment);
         break;
     }
//---
   if(!result)
      Print(__FUNCTION__, " ", trade.CheckResultRetcodeDescription());
//---
   return result;
  }
//+------------------------------------------------------------------+
//| Close market order                                               |
//+------------------------------------------------------------------+
bool CloseOrder(ulong ticket,ulong magic,ulong slp=10)
  {
//---
   CTrade trade;
//---
   trade.SetExpertMagicNumber(magic);
//---
   bool  result = trade.PositionClose(ticket,slp);
//---
   if(!result)
      Print(trade.CheckResultRetcodeDescription());
//---
   return result;
  }
//+------------------------------------------------------------------+
//|  Close all orders                                                |
//+------------------------------------------------------------------+
uint CloseAll(ulong magic_number=ULONG_MAX,const string sym=NULL,const ENUM_POSITION_TYPE ordertype=WRONG_VALUE,const ulong mslip=10)
  {
   CPositionInfo posinfo;
//---
   uint countclosed=0;
//---
   for(int i = PositionsTotal()-1; i>-1;--i)
     {
      if(!posinfo.SelectByIndex(i))
         continue;
      if(magic_number<ULONG_MAX && posinfo.Magic()!=long(magic_number))
         continue;
      if(sym!=NULL && posinfo.Symbol()!=sym)
         continue;
      if(ordertype!=WRONG_VALUE && posinfo.PositionType()!=ordertype)
         continue;
      if(CloseOrder(posinfo.Ticket(),mslip))
         countclosed++;
     }
//---
   return countclosed;
  }
//+------------------------------------------------------------------+

エントリーおよびエグジットの閾値はユーザーが調整可能であり、これらのパラメータの評価や最適化をおこなうことができます。また、EAではユーザー自身が任意のコピュラモデルやECDFモデルを指定してテストすることも可能です。なお、ここで示される結果は、銘柄レートのわずかな違いによって再現テスト時の結果と異なる場合があることに注意してください。

結果は次のとおりです。

EAレポート

エクイティグラフ

グラフを見ると、結果は概ね好意的であり、この種の戦略には一定の有効性があることが示唆されます。ただし、これらの結果はインサンプルテストに基づくものである点に留意する必要があります。価格が学習範囲を超えると、モデルは破綻し、インジケーターは無効なシグナルを出すようになります。この現象は、以下のスクリーンショットに示されているように、価格が新高値に向かって変動した際のインジケーターの挙動から確認できます。

GoldCopulsSignalsインジケーターのモデル破綻例

それでも、得られた結果は励みになるものであり、さらなる研究の動機付けとなります。


結論

本記事では、純粋なMQL5での代表的な二変量アルキメデスコピュラの実装を紹介しました。対象となったコピュラはフランク、ジョー、グンベル、クレイトン、N13、N14です。それぞれが捉えられる依存関係のタイプについても解説しました。記事の後半では、これらのコピュラをペアトレード戦略の構築に応用する方法を示し、実際に簡単な戦略を実装してその結果を観察しました。今後の記事では、コピュラについてさらに掘り下げていく予定です。以下に、記事内で参照したすべてのコードをまとめます。

ファイルまたはフォルダ   説明
MQL5/include/Copulas すべてのコピュラモデルを実装したヘッダファイルを格納したフォルダ
MQL5/include/Brent コピュラ実装で使用されるヘッダファイルを格納したフォルダ 
MQL5/include/ECDF ECDFの実装に関連するヘッダファイルを格納したフォルダ
MQL5/include/dependence.mqh 一致指標の実装を含むヘッダファイル
MQL5/include/info_criteria.mqh AIC、HQIC、BICを計算する関数を定義するファイル 
MQL5/include/np.mqh ベクトルと行列計算用のユーティリティ関数を定義するファイル 
MQL5/scripts/Concordance.mq5 候補銘柄の中から最も共動性の高いペアを見つけるスクリプト
MQL5/scripts/CopulaSelection.mq5 実装済みコピュラすべてのモデル選択基準を計算するスクリプト
MQL5/indicators/GoldCopulaSignals.mq5 EAのエントリーとエグジットのシグナルを生成するインジケーター 
MQL5/experts/SimpleCopulaPairsStrategy.mq5 コピュラを用いたペアトレード戦略を実装したEA 

MetaQuotes Ltdにより英語から翻訳されました。
元の記事: https://www.mql5.com/en/articles/19931

添付されたファイル |
ZIPをダウンロード
BivariateCopulae_Pt2.zip (53.39 KB)

警告: これらの資料についてのすべての権利はMetaQuotes Ltd.が保有しています。これらの資料の全部または一部の複製や再プリントは禁じられています。

この記事はサイトのユーザーによって執筆されたものであり、著者の個人的な見解を反映しています。MetaQuotes Ltdは、提示された情報の正確性や、記載されているソリューション、戦略、または推奨事項の使用によって生じたいかなる結果についても責任を負いません。

この著者による他の記事

ディスカッションに移動
MQL5における市場ポジショニング戦略の体系（第2回）： Nvidia向けマルチパターンのビット単位学習 MQL5における市場ポジショニング戦略の体系（第2回）： Nvidia向けマルチパターンのビット単位学習
管理可能なテスト期間において、特定の資産を特定の取引方向で検証する市場ポジショニングに関する新連載を継続します。前回の記事では、Nvidia Corp (NVDA)の株を対象に、RSIとDeMarkerオシレーターの組み合わせから5つのシグナルパターンを検証しました。本記事では残りの5パターンを取り上げ、さらに複数パターンの組み合わせにも踏み込みます。これには、10パターンすべての自由な組み合わせや、特定のペアのみを組み合わせる特殊パターンも含まれます。
プライスアクション分析ツールキットの開発（第49回）：トレンド系、モメンタム系、ボラティリティ系インジケーターを1つのMQL5システムに統合する プライスアクション分析ツールキットの開発（第49回）：トレンド系、モメンタム系、ボラティリティ系インジケーターを1つのMQL5システムに統合する
Multi Indicator Handler EAでMetaTrader 5のチャートをシンプルにしましょう。このインタラクティブなダッシュボードは、トレンド系、モメンタム系、ボラティリティ系インジケーターを1つのリアルタイムパネルに統合します。用途に応じてプロファイルを瞬時に切り替え、ワンクリックで表示と非表示を切り替えてチャートを整理し、プライスアクションに集中できます。本記事では、これをMQL5で自作してカスタマイズする手順をステップバイステップで解説します。
エラー 146 (「トレードコンテキスト ビジー」) と、その対処方法 エラー 146 (「トレードコンテキスト ビジー」) と、その対処方法
この記事では、MT4において複数のEAの衝突をさける方法を扱います。ターミナルの操作、MQL4の基本的な使い方がわかる人にとって、役に立つでしょう。
取引戦略の開発：バタフライオシレーター法 取引戦略の開発：バタフライオシレーター法
魅力的な数学概念であるバタフライ曲線を、実践的な取引ツールへと応用する方法を紹介します。バタフライオシレーターを構築し、それを基盤とした基本的な取引戦略を開発します。この戦略は、オシレーター特有の周期的シグナルと移動平均による従来型のトレンド確認を効果的に組み合わせることで、潜在的な市場エントリーポイントを特定するための体系的なアプローチを実現します。