Взаимная информация как критерий для поэтапного отбора признаков

Введение

Взаимная информация является ценным инструментом для выявления эффективных предикторов, особенно когда речь идет о сложных, нелинейных взаимосвязях. Он может выявлять зависимости, которые могут быть упущены другими методами, что делает его особенно подходящим для моделей, которые могут использовать такие сложности. В настоящей статье рассматривается применение взаимной информации при выборе признаков, основное внимание уделяется алгоритму, предложенному Ханчуан Пенгом (Hanchuan Peng), Фухуэй Лонгом (Fuhui Long) и Крисом Дином (Chris Ding) в их исследовательской работе под названием "Выбор признаков на основе взаимной информации: Критерии максимальной зависимости, максимальной релевантности и минимальной избыточности» ("Feature Selection Based on Mutual Information: Criteria of Max-Dependency, Max-Relevance, and Min-Redundancy").

Начнем с обсуждения оценки взаимной информации для непрерывных переменных, а затем углубимся в сам процесс отбора признаков. Наконец, мы проиллюстрируем эффективность алгоритма на примерах, использующих как синтетические, так и реальные наборы данных.

Оценка взаимной информации непрерывных переменных

Взаимная информация, обозначаемая как I(X;Y), количественно определяет общую информацию между двумя случайными величинами, X и Y. Для дискретных переменных ее вычисление относительно просто и включает суммирование общих и предельных вероятностей.

Однако непрерывные переменные представляют собой серьезную проблему из-за их бесконечного диапазона возможных значений. Распространенный подход к решению этой проблемы заключается в дискретизации непрерывных переменных по ячейкам и последующем применении формулы дискретной взаимной информации. К сожалению, точность этого метода очень чувствительна к ширине диапазона или количеству выбранных диапазонов, что может привести к существенным различиям в расчетной взаимной информации.

Чтобы вычислить взаимную информацию для непрерывных переменных, мы должны перейти от дискретных сумм к непрерывным интегралам. Это требует знания функций общей и предельной плотности вероятности (PDF), которые часто недоступны из-за ограниченности данных. Метод Парзеновского окна обеспечивает решение путем аппроксимации PDF непосредственно на основе данных. Этот метод включает в себя размещение оконной функции над данными и вычисление взвешенной суммы точек данных внутри окна. Вес, присвоенный каждой точке данных, уменьшается с увеличением ее расстояния от центра окна. Перемещая это окно по пространству данных и вычисляя взвешенную сумму в каждой точке, мы можем построить оценку функции плотности вероятности.

Этот метод особенно полезен в сценариях, где базовое распределение неизвестно или является сложным. Адаптируя форму и размер окна, мы можем повысить плавность и точность оценки плотности. Однако важно отметить, что выбор параметров окна может существенно повлиять на качество оценки. Эффективность оценки плотности Парцена зависит от выбора подходящей весовой функции, которая в идеале должна быть равна единице и демонстрировать быстрое затухание по мере удаления аргумента от нуля. Популярным выбором для этой весовой функции является Гауссовская функция, характеризуемая масштабным параметром sigma.

Однако определение оптимального значения для sigma является сложной задачей, часто приводящей к значительным колебаниям расчетной плотности. Такая чувствительность к выбору сигмы является основным недостатком метода Парзеновского окна, что делает его далеко не идеальным для таких задач, как оценка предикторов-кандидатов, где необходимым условием является точная оценка плотности.

Альтернативой методу Парзеновского окна для оценки взаимной информации между непрерывными переменными является адаптивное разбиение. Адаптивное разбиение предлагает значительное улучшение по сравнению с наивным подходом фиксированного биннинга. Путем итеративного разделения областей с высоким содержанием информации адаптивное разбиение эффективно фокусирует вычислительные ресурсы на наиболее важных областях пространства данных. Такой целенаправленный подход приводит к более точным оценкам взаимной информации, в частности, при работе со сложными, неравномерными распределениями. Ключевое преимущество адаптивного разбиения заключается в его способности адаптироваться к базовой структуре данных. Рекурсивный процесс разбиения гарантирует, что области со значительными зависимостями будут дополнительно разделены, в то время как области с минимальной информацией останутся нетронутыми. Такой динамический подход позволяет избежать распространенных ошибок фиксированного разбиения, которые могут привести к чрезмерному сглаживанию или перегруппировке данных.

Кроме того, адаптивное разбиение включает в себя статистические тесты для управления процессом разбиения. Оценивая статистическую значимость потенциальных подразделений, адаптивное разбиение позволяет отличить подлинную информацию от случайного шума. Это помогает предотвратить переобучение, что в противном случае могло бы привести к завышенным оценкам взаимной информации. Основная проблема наивного метода заключается в том, что он уделяет чрезмерное внимание областям бивариатного домена с небольшим количеством случаев или вообще без них, в то же время потенциально пренебрегая плотными областями, где содержится большая часть информации.  

Адаптивное разбиение начинается с крупнозернистого разбиения пространства данных. Для каждого раздела вычисляется взаимная информация между переменными. Если это значение превышает заранее установленное пороговое значение, раздел далее подразделяется на более мелкие подразделы. Этот рекурсивный процесс разбиения продолжается до тех пор, пока не будет выполнен критерий остановки. Если процесс разбиения остановить слишком рано, результирующие разбиения могут быть слишком грубыми, что приведет к потере детализации и смещение в сторону понижения в предполагаемой взаимной информации. И наоборот, если процесс разбиения продолжается слишком долго, разделы могут стать слишком мелкими, что приведет к переобучению и завышению оценки взаимной информации из-за шума.

Чтобы определить, следует ли далее разбивать раздел на подразделы, используется тест хи-квадрат. Этот тест оценивает независимость между двумя переменными внутри раздела. Раздел разделен на четыре подраздела, создавая таблицу сопряженности 2x2. Подсчитывается количество точек данных, попадающих в каждый из четырех подразделов. В соответствии с нулевой гипотезой о независимости двух переменных ожидаемое количество точек данных в каждом подразделе рассчитывается на основе предельных итоговых значений. Затем вычисленная статистика хи-квадрат сравнивается с критическим значением из распределения хи-квадрат с 1 степенью свободы. Если вычисленное значение превышает критическое, нулевая гипотеза независимости отклоняется, что указывает на существенную взаимосвязь между двумя переменными в рамках раздела.

Если тест хи-квадрат является значительным, раздел далее подразделяется на более мелкие подразделы. Этот процесс продолжается рекурсивно до тех пор, пока разделы не станут достаточно однородными или пока не будет выполнен критерий остановки. Если тест хи-квадрат не является значительным, раздел считается однородным и дальнейшее подразделение не требуется.

Чтобы разобраться с возможностью возникновения более сложных взаимосвязей внутри раздела, которые могут быть пропущены при простом разделении 2x2, алгоритм включает дополнительную проверку. Если тест хи-квадрат 2х2 не выявляет существенной взаимосвязи, но раздел остается относительно большим, выполняется более точное разбиение 4х4. Затем к этому разделению 4х4 применяется тест хи-квадрат, чтобы оценить наличие неслучайного распределения. Если тест 4х4 также не выявляет существенной взаимосвязи, то раздел считается однородным и рассчитывается его вклад в общую взаимную информацию. В этом случае дальнейшее деление не требуется.

Если тест хи-квадрат 2х2 или 4х4 указывает на неравномерное распределение в пределах раздела, раздел делится на четыре меньших подраздела. Эти подразделы обрабатываются по-разному путем проверки их размера. Если подраздел слишком мал, он считается конечным. Вычисляется его вклад в общую взаимную информацию, и раздел больше не рассматривается для дальнейшего разделения. В противном случае, если подраздел достаточно велик, он помечается как кандидат на будущее дополнительное разделение.

Стандартная статистика по тесту хи-квадрат два на два приведена в указанном ниже уравнении.

Адаптивное разбиение на MQL5

Класс Capm в mutualinfo.mqh реализует алгоритм адаптивного разделения для оценки взаимной информации для непрерывных переменных. Для управления рекурсивным процессом разбиения он использует пользовательскую структуру IntStack. Эта структура содержит информацию о разделах и любых подразделах, содержащихся в ней. Структура состоит из шести членов, подробно описанных ниже:

• Xstart и Xstop: Они определяют диапазон индексов для одной переменной в пределах текущего охватывающего раздела.
• Ystart и Ystop: Они определяют диапазон индексов для второй переменной в пределах текущего охватывающего раздела.
• DataStart и DataStop: Они задают начальный и конечный индексы точек данных внутри охватывающего раздела, размечая подразделы.
  

//+------------------------------------------------------------------+
//| int type stack structure                                         |
//+------------------------------------------------------------------+
struct IntStack
  {
   int Xstart ;     // X value (rank) at which this rectangle starts
   int Xstop ;      // And stops
   int Ystart ;     // Ditto for Y
   int Ystop ;
   int DataStart ;  // Starting index into indices for the cases in this
   int DataStop ;   // rectangle, and the (inclusive) ending index
  };

Конструктор Capm принимает три аргумента:

• no_split: Логический флаг, определяющий, как обрабатывать значения, которые приблизительно равны. Если установлено значение true, алгоритм гарантирует, что такие значения будут сгруппированы во время разбиения, предотвращая их разделение по разным разделам.
• dep_vals: Вектор, содержащий значения одной из анализируемых переменных. Эта переменная, как правило, является той, взаимная информация о которой будет рассчитываться на основе множества других переменных.
• chi_critical_value: Пороговое значение для тестов хи-квадрат, используемых в процессе разбиения. Это значение определяет уровень значимости для тестов и влияет на глубину разбиения. Как правило, подходят значения от 4,0 до 8,0, при этом 6,0 является обычным выбором с универсальной полезностью.
  

//+------------------------------------------------------------------+
//|  constructor                                                     |
//+------------------------------------------------------------------+
Capm::Capm(bool no_split,vector &dep_vals, double chi_critical_value = 6.0)
  {
   m_no_split = no_split;
   m_size = int(dep_vals.Size()) ;
   m_chi_crit = chi_critical_value;

   if(ArrayResize(m_indices,m_size)<0 || !m_tempvals.Resize(m_size) || ArrayResize(m_ranks,m_size)<0
      || (m_no_split && ArrayResize(m_tied_ranks,m_size)<0))
     {
      Print(__FUNCTION__, " Arrayresize error ", GetLastError());
      return;
     }

   for(int i=0 ; i<m_size ; i++)
     {
      m_tempvals[i] = dep_vals[i] ;
      m_indices[i] = i ;
     }

   qsortdsi(0, m_size-1, m_tempvals, m_indices) ;

   for(int i=0 ; i<m_size ; i++)
     {
      m_ranks[m_indices[i]] = i ;
      if(! m_no_split)
         continue ;
      if(i < m_size-1  &&
         m_tempvals[i+1] - m_tempvals[i] < 1.e-12 * (1.0+fabs(m_tempvals[i])+fabs(m_tempvals[i+1])))
         m_tied_ranks[i] = 1 ;
      else
         m_tied_ranks[i] = 0 ;
     }
  }

Конструктор класса Capm начинается с установки флага no_split, который определяет, следует ли учитывать связи в данных, и параметра chi_critical_value, который устанавливает пороговое значение для теста хи-квадрат, используемого для выявления значимых взаимосвязей между переменными. Затем конструктор выделяет память для массивов для хранения данных, рангов и связанных значений. Он копирует значения переменных во временный вектор и сортирует их. Исходные индексы точек данных сохраняются в массиве m_indices. Наконец, конструктор присваивает ранги отсортированным точкам данных и идентифицирует любые связанные значения, отмечая их в массиве m_tied_ranks, если установлен флаг no_split.

//+------------------------------------------------------------------+
//| Mutual information using the Adaptive partitioning method        |
//+------------------------------------------------------------------+
class Capm
  {

public:
                     Capm(bool no_split,vector &dep_vals, double chi_critical_value = 6.0) ;
                    ~Capm(void) ;
   double            fit(vector &raw) ;

private:
   bool               m_no_split;  // respect ties
   int m_size ;                        // Number of cases
   int m_ranks[] ;                     // 'Dependent' variable ranks
   int m_tied_ranks[] ;                // tied[i] != 0 if case with rank i == case with rank i+1
   double m_chi_crit ;                 // Chi-square test criterion
   int               m_indices[];      // indices
   vector            m_tempvals;       // temp values
  } ;

Метод fit() в классе Capm предназначен для оценки взаимной информации между одной переменной (предоставленной конструктору) и другой переменной, предоставленной в качестве аргумента. В качестве входных данных он принимает один вектор, представляя значения второй переменной. Этот вектор должен иметь тот же размер, что и вектор, предоставленный конструктору. Вызывая метод fit() множество раз с разными векторами, мы можем эффективно вычислять взаимную информацию между одной переменной и различными переменными.

//+------------------------------------------------------------------+
//|  mutual information for continuous variables                     |
//+------------------------------------------------------------------+
double Capm::fit(vector &raw)
  {
   int i, k, ix, iy, nstack, splittable ;
   int current_indices[], x[], fullXstart, fullXstop, fullYstart, fullYstop, ipos ;
   int trialXstart[4], trialXstop[4], trialYstart[4], trialYstop[4] ;
   int ipx, ipy, xcut[4], ycut[4], iSubRec, x_tied[], ioff ;
   int X_AllTied, Y_AllTied ;
   int centerX, centerY, currentDataStart, currentDataStop ;
   int actual[4], actual44[16] ;
   vector expected(16), xfrac(4), yfrac(4) ;
   double diff, testval;
   double px, py, pxy, MI ;

   IntStack stack[];

   if(ArrayResize(stack,1,256)<0)
     {
      Print(__FUNCTION__," arrayresize error ", GetLastError());
      return EMPTY_VALUE;
     }

   if(ArrayResize(current_indices,m_size)<0 || ArrayResize(x,m_size)<0
      || (m_no_split && ArrayResize(x_tied,m_size)<0))
     {
      Print(__FUNCTION__, " Arrayresize error ", GetLastError());
      return EMPTY_VALUE;
     }

   for(i=0 ; i<m_size ; i++)
     {
      m_tempvals[i] = raw[i] ;
      m_indices[i] = i ;
     }

   if(!qsortdsi(0, m_size-1, m_tempvals, m_indices))
     {
      Print(__FUNCTION__, " Arrayresize error ", GetLastError());
      return EMPTY_VALUE;
     }

   for(i=0 ; i<m_size ; i++)
     {
      x[m_indices[i]] = i ;
      if(! m_no_split)
         continue ;
      if(i < m_size-1  &&
         m_tempvals[i+1] - m_tempvals[i] < 1.e-12 * (1.0+fabs(m_tempvals[i])+fabs(m_tempvals[i+1])))
         x_tied[i] = 1 ;
      else
         x_tied[i] = 0 ;
     }

   for(i=0 ; i<m_size ; i++)
     {
      m_indices[i] = i ;
     }

   stack[0].Xstart = 0 ;
   stack[0].Xstop = m_size-1 ;
   stack[0].Ystart = 0 ;
   stack[0].Ystop = m_size-1 ;
   stack[0].DataStart = 0 ;
   stack[0].DataStop = m_size-1 ;
   nstack = 1 ;

   MI = 0.0 ;
   while(nstack > 0)
     {

      --nstack ;
      fullXstart = stack[nstack].Xstart ;
      fullXstop = stack[nstack].Xstop ;
      fullYstart = stack[nstack].Ystart ;
      fullYstop = stack[nstack].Ystop ;
      currentDataStart = stack[nstack].DataStart ;
      currentDataStop = stack[nstack].DataStop ;


      centerX = (fullXstart + fullXstop) / 2 ;
      X_AllTied = (x_tied.Size())  && (x_tied[centerX] != 0) ;
      if(X_AllTied)
        {
         for(ioff=1 ; centerX-ioff >= fullXstart ; ioff++)
           {
            if(! x_tied[centerX-ioff])
              {
               X_AllTied = 0 ;
               centerX -= ioff ;
               break ;
              }
            if(centerX + ioff == fullXstop)
               break ;
            if(! x_tied[centerX+ioff])
              {
               X_AllTied = 0 ;
               centerX += ioff ;
               break ;
              }
           }
        }

      centerY = (fullYstart + fullYstop) / 2 ;
      Y_AllTied = (m_tied_ranks.Size())  && (m_tied_ranks[centerY] != 0) ;
      if(Y_AllTied)
        {
         for(ioff=1 ; centerY-ioff >= fullYstart ; ioff++)
           {
            if(! m_tied_ranks[centerY-ioff])
              {
               Y_AllTied = 0 ;
               centerY -= ioff ;
               break ;
              }
            if(centerY + ioff == fullYstop)
               break ;
            if(! m_tied_ranks[centerY+ioff])
              {
               Y_AllTied = 0 ;
               centerY += ioff ;
               break ;
              }
           }
        }

      if(X_AllTied  ||  Y_AllTied)
         splittable = 0 ;
      else
        {
         trialXstart[0] = trialXstart[1] = fullXstart ;
         trialXstop[0] = trialXstop[1] = centerX ;
         trialXstart[2] = trialXstart[3] = centerX+1 ;
         trialXstop[2] = trialXstop[3] = fullXstop ;
         trialYstart[0] = trialYstart[2] = fullYstart ;
         trialYstop[0] = trialYstop[2] = centerY ;
         trialYstart[1] = trialYstart[3] = centerY+1 ;
         trialYstop[1] = trialYstop[3] = fullYstop ;

         for(i=0 ; i<4 ; i++)
            expected[i] = (currentDataStop - currentDataStart + 1) *
                          (trialXstop[i]-trialXstart[i]+1.0) / (fullXstop-fullXstart+1.0) *
                          (trialYstop[i]-trialYstart[i]+1.0) / (fullYstop-fullYstart+1.0) ;

         actual[0] = actual[1] = actual[2] = actual[3] = 0 ;
         for(i=currentDataStart ; i<=currentDataStop ; i++)
           {
            k = m_indices[i] ;
            if(x[k] <= centerX)
              {
               if(m_ranks[k] <= centerY)
                  ++actual[0] ;
               else
                  ++actual[1] ;
              }
            else
              {
               if(m_ranks[k] <= centerY)
                  ++actual[2] ;
               else
                  ++actual[3] ;
              }
           }

         testval = 0.0 ;
         for(i=0 ; i<4 ; i++)
           {
            diff = fabs(actual[i] - expected[i]) - 0.5 ;
            testval += diff * diff / expected[i] ;
           }

         splittable = (testval > m_chi_crit)  ?  1 : 0 ;

Метод начинается с инициализации стека для управления процессом разбиения. Весь набор данных изначально помещается в стек в качестве корневого раздела. Затем алгоритм итеративно извлекает раздел из стека и оценивает его однородность. Для оценки статистической значимости любых потенциальных взаимосвязей между переменными в рамках раздела, применяется тест хи-квадрат. Если тест указывает на существенную взаимосвязь, секция разбивается на четыре подраздела, и каждый подраздел помещается обратно в стек для дальнейшего рассмотрения. Этот процесс продолжается до тех пор, пока все разделы не будут соответствовать критерию остановки или не станут однородными.

if(! splittable && fullXstop-fullXstart > 30 && fullYstop-fullYstart > 30)
           {
            ipx = fullXstart - 1 ;
            ipy = fullYstart - 1 ;
            for(i=0 ; i<4 ; i++)
              {
               xcut[i] = (fullXstop - fullXstart + 1) * (i+1) / 4 + fullXstart - 1 ;
               xfrac[i] = (xcut[i] - ipx) / (fullXstop - fullXstart + 1.0) ;
               ipx = xcut[i] ;
               ycut[i] = (fullYstop - fullYstart + 1) * (i+1) / 4 + fullYstart - 1 ;
               yfrac[i] = (ycut[i] - ipy) / (fullYstop - fullYstart + 1.0) ;
               ipy = ycut[i] ;
              }
            for(ix=0 ; ix<4 ; ix++)
              {
               for(iy=0 ; iy<4 ; iy++)
                 {
                  expected[ix*4+iy] = xfrac[ix] * yfrac[iy] *
                                      (currentDataStop - currentDataStart + 1) ;
                  actual44[ix*4+iy] = 0 ;
                 }
              }
            for(i=currentDataStart ; i<=currentDataStop ; i++)
              {
               k = m_indices[i] ;
               for(ix=0 ; ix<3 ; ix++)
                 {
                  if(x[k] <= xcut[ix])
                     break ;
                 }
               for(iy=0 ; iy<3 ; iy++)
                 {
                  if(m_ranks[k] <= ycut[iy])
                     break ;
                 }
               ++actual44[ix*4+iy] ;
              }
            testval = 0.0 ;
            for(ix=0 ; ix<4 ; ix++)
              {
               for(iy=0 ; iy<4 ; iy++)
                 {
                  diff = fabs(actual44[ix*4+iy] - expected[ix*4+iy]) - 0.5 ;
                  testval += diff * diff / expected[ix*4+iy] ;
                 }
              }
            splittable = (testval > 3 * m_chi_crit)  ?  1 : 0 ;
           }
        }

      if(splittable)
        {

         for(i=currentDataStart ; i<=currentDataStop ; i++)
            current_indices[i] = m_indices[i] ;

         ipos = currentDataStart ;
         for(iSubRec=0 ; iSubRec<4 ; iSubRec++)
           {

            if(actual[iSubRec] >= 3)
              {
               stack[nstack].Xstart = trialXstart[iSubRec] ;
               stack[nstack].Xstop = trialXstop[iSubRec] ;
               stack[nstack].Ystart = trialYstart[iSubRec] ;
               stack[nstack].Ystop = trialYstop[iSubRec] ;
               stack[nstack].DataStart = ipos ;
               stack[nstack].DataStop = ipos + actual[iSubRec] - 1 ;
               ++nstack ;
               if(ArrayResize(stack,nstack+1,256)<0)
                 {
                   Print(__FUNCTION__," arrayresize error ", GetLastError());
                   return EMPTY_VALUE;
                 }

               if(iSubRec == 0)
                 {
                  for(i=currentDataStart ; i<=currentDataStop ; i++)
                    {
                     k = current_indices[i] ;
                     if(x[k] <= centerX  &&  m_ranks[k] <= centerY)
                        m_indices[ipos++] = current_indices[i] ;
                    }
                 }
               else
                  if(iSubRec == 1)
                    {
                     for(i=currentDataStart ; i<=currentDataStop ; i++)
                       {
                        k = current_indices[i] ;
                        if(x[k] <= centerX  &&  m_ranks[k] > centerY)
                           m_indices[ipos++] = current_indices[i] ;
                       }
                    }
                  else
                     if(iSubRec == 2)
                       {
                        for(i=currentDataStart ; i<=currentDataStop ; i++)
                          {
                           k = current_indices[i] ;
                           if(x[k] > centerX  &&  m_ranks[k] <= centerY)
                              m_indices[ipos++] = current_indices[i] ;
                          }
                       }
                     else
                       {
                        for(i=currentDataStart ; i<=currentDataStop ; i++)
                          {
                           k = current_indices[i] ;
                           if(x[k] > centerX  &&  m_ranks[k] > centerY)
                              m_indices[ipos++] = current_indices[i] ;
                          }
                       }
              }

            else
              {
               if(actual[iSubRec] > 0)
                 {
                  px = (trialXstop[iSubRec] - trialXstart[iSubRec] + 1.0) / m_size ;
                  py = (trialYstop[iSubRec] - trialYstart[iSubRec] + 1.0) / m_size ;
                  pxy = (double) actual[iSubRec] / m_size ;
                  MI += pxy * log(pxy / (px * py)) ;
                 }
              }
           }
        }
      else
        {
         px = (fullXstop - fullXstart + 1.0) / m_size ;
         py = (fullYstop - fullYstart + 1.0) / m_size ;
         pxy = (currentDataStop - currentDataStart + 1.0) / m_size ;
         MI += pxy * log(pxy / (px * py)) ;
        }
     }

   return MI;

Для разделов, которые считаются однородными либо из-за незначительного теста хи-квадрат, либо из-за минимального количества точек данных, вычисляется взаимная информация. Этот расчет включает в себя оценку общего и предельного распределений вероятностей внутри раздела и применение приведенной ниже формулы.

Процесс продолжается рекурсивно, при этом разделы разбиваются и оцениваются до тех пор, пока дальнейшее значительное разделение станет невозможным. Окончательная оценка взаимной информации получается путем суммирования вкладов всех разделов терминала, что позволяет получить комплексную оценку, основанную на структуре данных. Если метод обнаруживает ошибку, он возвращает значение, эквивалентное встроенной в MQL5 константе EMPTY_VALUE. В следующем разделе мы будем использовать класс Capm для реализации алгоритма поэтапного отбора признаков.

Отбор предиктора с использованием взаимной информации

При применении взаимной информации к задаче отбора признаков наша цель состоит в том, чтобы идентифицировать подмножество предикторов из большего набора потенциальных предикторов, что максимизирует совместную зависимость с целевой переменной. Эта совместная зависимость представляет собой обобщение взаимной информации, где одной из величин является набор переменных, а не одна переменная. Вычисление совместной зависимости для больших подмножеств предикторов становится трудноразрешимым с точки зрения вычислений из-за проклятия размерности. Это происходит потому, что требуемое многомерное разбиение данных приводит к чрезвычайно малому количеству ячеек, особенно по мере увеличения числа измерений (предикторов).

Кроме того, огромное количество возможных комбинаций признаков может быть ошеломляющим даже для наборов признаков среднего размера. Этот комбинаторный взрыв делает непрактичным исчерпывающий поиск по всем возможным подмножествам, делая необходимым использование эффективных стратегий поиска.

Распространенным подходом к отбору признаков является прямой поэтапный отбор. Этот метод начинается с пустой модели и итеративно добавляет признак, обеспечивающий наибольшее улучшение эффективности модели, измеряемое выбранным критерием. Хотя этот подход эффективен, он страдает от ограничения «жадности». В нем могут отсутствовать оптимальные комбинации признаков, например, в тех случаях, когда два признака в сочетании обеспечивают значительно лучшую прогностическую способность, чем любой из них по отдельности. Это связано с тем, что алгоритм фокусируется на постепенных улучшениях на каждом этапе, потенциально упуская из виду синергетические взаимосвязи между признаками.

В то время как существуют другие методы, такие как методы высшего порядка и обратного отбора, прямой поэтапный отбор остается наиболее практичным выбором из-за вычислительной эффективности, особенно при работе с большими наборами данных и ограниченными вычислительными ресурсами. Кроме того, при применении к конкретному контексту совместной зависимости прямой поэтапный отбор демонстрирует удивительное свойство: он может эффективно аппроксимировать оптимальное подмножество признаков, даже если он явно не оптимизирует меру совместной зависимости.

Пэнг (Peng), Лонг (Long) и Дин (Ding) предложили алгоритм отбора признаков, основанный на взаимной информации, известный как критерий максимальной релевантности и минимальной избыточности (MRMR). Этот алгоритм эффективно аппроксимирует оптимальное подмножество признаков без явного вычисления совместной зависимости, что является трудноразрешимым с точки зрения вычислений. Критерий MRMR основан на концепции релевантности и избыточности. Соответствие набора признаков S целевой переменной Y определяется как среднее значение взаимной информации между Y и каждым признаком в S.

Здесь ∣S∣ - количество признаков в S, а I(Xi,Y) - взаимная информация между признаками Xi и Y.

Хотя максимизация релевантности интуитивно понятна, она может привести к неоптимальному набору признаков из-за избыточности. Если мы просто отберем признаки на основе их индивидуальной релевантности целевой переменной, то в итоге можем получить набор сильно коррелированных признаков, которые предоставляют мало дополнительной информации. Чтобы решить эту проблему, нам необходимо учитывать не только релевантность признака, но и его избыточность по отношению к уже отобранным признакам.

Избыточность потенциального предиктора Xi по отношению к уже отобранному набору признаков S определяется как средняя взаимная информация между Xi и каждым признаком в S:

Где I(Xi,Xj) - это взаимная информация между потенциальным признаком Xi и признаком Xj, уже имеющимся в S.

Важно отметить, что математическое выражение для избыточности идентично выражению для релевантности. Единственное различие заключается в интерпретации. При вычислении взаимной информации между потенциальным признаком и целевой переменной мы называем это релевантностью. Однако при вычислении взаимной информации между потенциальным признаком и признаком, уже входящим в отобранный набор, мы называем это избыточностью.

По сути, обе концепции предполагают измерение объема информации, совместно используемой двумя переменными, но контекст определяет, будет ли она считаться релевантной или избыточной. В случае релевантности мы обращаем внимание на взаимосвязь между признаком и целевой переменной, в то время как в случае избыточности мы фокусируемся на том, насколько сильно существует совпадение между потенциальным признаком и уже отобранными признаками.

Подводя итог, можно сказать, что алгоритм MRMR работает следующим образом:

1. В качестве первого признака отбирается признак с наибольшей взаимной информацией с целевой переменной Y.
2. На каждом последующем этапе алгоритм отбирает признак, максимизирующий следующий критерий.
   

Этот критерий уравновешивает релевантность признака целевой переменной Y (т.е. I(X_i, Y)) с его избыточностью по отношению к уже отобранному признаку. Замечательным аспектом алгоритма MRMR является то, что он эффективно аппроксимирует оптимальное подмножество признаков, даже несмотря на то, что прямая оптимизация совместной зависимости является трудноразрешимой с точки зрения вычислений. Этот результат, хотя и не очевиден сразу, математически доказан в оригинальной статье.

Для оценки статистической значимости отобранных признаков можно использовать два теста на перестановку Монте-Карло (MCP):

• Значение индивидуального признака: Этот тест оценивает значимость каждого отдельного признака путем сравнения его соответствия целевой переменной с нулевым распределением, полученным путем перестановки. Переставляя значения признаков, мы создаем нулевое распределение, предполагающее, что признак не имеет отношения к цели. Если наблюдаемая релевантность значительно выше, чем измененные значения, можно заключить, что признак, скорее всего, действительно информативен.
• Общая значимость модели: Этот тест оценивает значимость всего набора отобранных признаков путем сравнения совокупной суммы релевантности отдельных признаков с нулевым распределением, полученным путем перестановки. Переставляя значения целевых переменных, мы создаем нулевое распределение, предполагающее, что признаки не имеют связи с целью. Если наблюдаемая совокупная релевантность значительно выше, чем измененные значения, можно заключить, что отобранный набор признаков, скорее всего, будет предсказывать целевую переменную.
  

Реализация поэтапного отбора на MQL5 на основе взаимной информации

Класс Cmrmr в mutualinfo.mqh реализует алгоритм MRMR и принимает в своем конструкторе следующие параметры:

• num_reps: Количество повторений для теста на перестановку Монте-Карло. Установка этого параметра на значение, меньшее или равное 1, отключает тесты.
• max_preds: Максимальное количество отбираемых признаков.
• chisquare_thresh: Пороговое значение теста хи-квадрат для метода адаптивного разбиения, используемого для оценки взаимной информации.
• m_verbose: Логический флаг, определяющий, следует ли отображать подробные выходные данные во время выполнения алгоритма.
  

//+-----------------------------------------------------------------------+
//|Relevance minus redundancy for building an optimal subset of predictors|
//+-----------------------------------------------------------------------+
class Cmrmr
  {
private:
   int               m_max_preds;
   int               m_reps;
   bool              m_verbose;
   int               m_samples;
   int               m_vars;
   matrix            m_preds;
   vector            m_target;
   vector            m_crits;
   double            m_chisquare_thresh;
   ulong             m_selected_vars[];
   matrix            m_critical_values;
   vector            mutualinfo(int which_size,int &which[],vector &targs);
public:
                     Cmrmr(int num_reps,int max_preds, double chisquare_thresh,bool verbose);
                    ~Cmrmr(void);
   bool              StepWise(matrix &predictors, vector &targets);
   matrix            GetCriticalValues(void);
   bool              GetSelectedVars(ulong &output[]);
  };

Чтобы инициировать процесс отбора признаков с помощью класса Cmrmr, необходимо инстанцировать объект этого класса, а затем вызвать его метод StepWise(). Метод StepWise() использует два основных входных сигнала:

• Матрица потенциальных предикторных переменных: Эта матрица должна содержать данные по всем потенциальным признакам, которые должен рассмотреть практик.
• Вектор целевой переменной: Этот вектор должен содержать значения целевой переменной.

Метод возвращает логическое значение, указывающее на успех или неудачу процесса отбора признаков. Метод StepWise() в классе Cmrmr инициализирует необходимые структуры данных, а затем переходит к циклу тестирования MCP. Если MCP-тесты не запрашиваются (`num_reps` меньше или равно 1), цикл выполняется только один раз. В рамках цикла метод вычисляет взаимную информацию между каждым потенциальным предиктором и целевой переменной, используя экземпляр класса Capm.

//+------------------------------------------------------------------+
//|  Stepwise feature selection based on mutual information          |
//+------------------------------------------------------------------+
bool Cmrmr::StepWise(matrix &predictors,vector &targets)
  {
   if(m_selected_vars.Size())
      ArrayFree(m_selected_vars);

   m_preds = predictors;
   m_samples = int(m_preds.Rows());
   m_vars = int(m_preds.Cols());
   m_max_preds = m_max_preds>=m_vars?m_vars:m_max_preds;
   m_target = targets;
   
   
   int i, j, k, ivar, irep;
   int index[], stepwise_mcpt_count[], solo_mcpt_count[], stepwise_ivar[], which_preds[],original_stepwise_ivar[] ;

   int nkept,best_ivar;
   vector casework(m_samples), sorted(m_samples), mutual(m_vars);
   vector crit(m_vars), relevance(m_vars), original_relevance(m_vars), current_crits(m_vars), sorted_crits(m_vars);
   double best_crit, dtemp, group_pvalue,solo_pvalue;
   vector stepwise_crit(m_vars), original_stepwise_crit(m_vars);
   double sum_relevance;
   vector original_sum_relevance(m_vars), sum_redundancy(m_vars);
   vector target = m_target;

   best_crit = -DBL_MAX;
   best_ivar = -1;
   nkept = m_max_preds;

   if(ArrayResize(index,m_vars)<0 || ArrayResize(stepwise_mcpt_count,m_vars)<0 ||
      ArrayResize(solo_mcpt_count,m_vars)<0 || ArrayResize(which_preds,m_vars)<0 || ArrayResize(stepwise_ivar,m_vars)<0 ||
      ArrayResize(original_stepwise_ivar,m_vars)<0)
     {
      Print(__FUNCTION__," array resize error ", GetLastError());
      return false;
     }

   int unif_error;
   for(irep=0 ; irep<m_reps ; irep++)
     {
      if(irep)
        {
         i = m_samples ;
         while(i > 1)
           {
            j = (int)(MathRandomUniform(0.0,1.0,unif_error) * i);
            if(unif_error)
              {
               Print(__FUNCTION__," Mathrandomuniform error ", GetLastError());
               return false;
              }
            if(j >= i)
               j = i - 1 ;
            dtemp = target[--i] ;
            target[i] = target[j] ;
            target[j] = dtemp ;
           }
        }

      for(i=0 ; i<m_vars ; i++)
         which_preds[i] = i ;

      crit = mutualinfo(m_vars,which_preds,target);

      for(ivar=0 ; ivar<m_vars ; ivar++)
        {
         relevance[ivar] = crit[ivar] ;
         if(ivar == 0  ||  crit[ivar] > best_crit)
           {
            best_crit = crit[ivar] ;
            best_ivar = ivar ;
           }
        }

      stepwise_crit[0] = best_crit ;
      stepwise_ivar[0] = best_ivar ;
      sum_relevance = best_crit ;
      if(irep == 0)
        {

         original_stepwise_crit[0] = best_crit ;
         original_stepwise_ivar[0] = best_ivar ;
         original_sum_relevance[0] = sum_relevance ;
         stepwise_mcpt_count[0] = 1 ;

         for(ivar=0 ; ivar<m_vars ; ivar++)
           {
            index[ivar] = ivar ;
            original_relevance[ivar] = sorted_crits[ivar] = current_crits[ivar] = crit[ivar] ;
            solo_mcpt_count[ivar] = 1 ;
           }

         if(!qsortdsi(0, m_vars-1, sorted_crits, index))
           {
            Print(__FUNCTION__, " failed qsort ");
            return false;
           }

         if(m_verbose)
            Print("       Variable         MI") ;

         for(i=m_vars-1 ; i>=0 ; i--)
           {
            k = index[i] ;

            if(m_verbose)
               PrintFormat("%15s %12.4lf",string(k), current_crits[k]) ;
           }
        }

      else
        {
         if(sum_relevance >= original_sum_relevance[0])
            ++stepwise_mcpt_count[0] ;
         for(ivar=0 ; ivar<m_vars ; ivar++)
           {
            if(relevance[ivar] >= original_relevance[ivar])
               ++solo_mcpt_count[ivar] ;
           }
        }

      for(i=0 ; i<m_vars ; i++)
         sum_redundancy[i] = 0.0 ;

      for(nkept=1 ; nkept<m_max_preds ; nkept++)
        {

         if(irep == 0 && m_verbose)
           {
            Print("Predictors so far   Relevance   Redundancy   Criterion") ;
            for(i=0 ; i<nkept ; i++)
              {
               k = stepwise_ivar[i] ;
               PrintFormat("%15s %12.4lf %12.4lf %12.4lf",string(k), relevance[k], relevance[k] - stepwise_crit[i], stepwise_crit[i]) ;
              }
           }

         k = 0 ;
         for(i=0 ; i<m_vars ; i++)
           {
            for(j=0 ; j<nkept ; j++)
              {
               if(stepwise_ivar[j] == i)
                  break ;
              }
            if(j == nkept)
               which_preds[k++] = i ;
           }
         if(k != (m_vars - nkept))
           {
            Print(__FUNCTION__, " failed assertion ", __LINE__);
            return false;
           }

         k = stepwise_ivar[nkept-1] ;

         casework = m_preds.Col(k);

         crit = mutualinfo(m_vars-nkept,which_preds,casework);

         for(i=0 ; i<(m_vars-nkept) ; i++)
           {
            k = which_preds[i] ;
            sum_redundancy[k] += crit[i] ;
            index[i] = k ;
            sorted_crits[i] = current_crits[i] = ((relevance[k] - sum_redundancy[k]) / double(nkept)) ;
            if(i == 0  ||  current_crits[i] > best_crit)
              {
               best_crit = current_crits[i] ;
               best_ivar = k ;
              }
           }

         stepwise_crit[nkept] = best_crit ;
         stepwise_ivar[nkept] = best_ivar ;
         sum_relevance += relevance[best_ivar] ;

         if(irep == 0)
           {
            original_stepwise_crit[nkept] = best_crit ;
            original_stepwise_ivar[nkept] = best_ivar ;
            original_sum_relevance[nkept] = sum_relevance ;
            stepwise_mcpt_count[nkept] = 1 ;

            if(!qsortdsi(0, m_vars-nkept-1, sorted_crits, index))
              {
               Print(__FUNCTION__, " failed qsort ");
               return false;
              }

            if(m_verbose)
              {
               Print("Additional candidates, in order of decreasing relevance minus redundancy") ;
               Print("       Variable     Relevance   Redundancy   Criterion") ;

               for(i=m_vars-nkept-1 ; i>=0 ; i--)
                 {
                  k = index[i] ;
                  PrintFormat("%15s %12.4lf %12.4lf %12.4lf",
                              string(k), relevance[k], sum_redundancy[k] / nkept,
                              relevance[k] - sum_redundancy[k] / nkept) ;
                 }
              }
           }

         else
           {
            if(sum_relevance >= original_sum_relevance[nkept])
               ++stepwise_mcpt_count[nkept] ;
           }

        }

     }
//---
   m_critical_values = matrix::Zeros(nkept,m_reps>1?5:3);
//---
   if(ArrayResize(m_selected_vars,nkept)<0)
     {
      Print(__FUNCTION__, " failed array resize ", GetLastError());
      return false;
     }
//---
   if(m_verbose)
     {
      if(m_reps > 1)
         Print("Final predictors ||  Relevance ||  Redundancy || Criterion  ||  Solo pval || Group pval") ;
      else
         Print("Final predictors  ||  Relevance ||  Redundancy ||  Criterion") ;
     }
//---
   for(i=0 ; i<nkept ; i++)
     {
      k = original_stepwise_ivar[i] ;
      m_selected_vars[i] = ulong(k);

      m_critical_values[i][0] = original_relevance[k];
      m_critical_values[i][1] = original_relevance[k] - original_stepwise_crit[i];
      m_critical_values[i][2] = original_stepwise_crit[i];

      if(m_critical_values.Cols()>3)
        {
         group_pvalue = (double) stepwise_mcpt_count[i] / (double) m_reps;
         solo_pvalue = (double) solo_mcpt_count[k] / (double) m_reps;
         m_critical_values[i][3] = solo_pvalue;
         m_critical_values[i][4] = group_pvalue;
        }

      if(m_verbose)
        {
         if(m_reps > 1)
            PrintFormat("%15s %12.4lf %12.4lf %12.4lf    %8.3lf    %8.3lf",string(k), m_critical_values[i][0], m_critical_values[i][1], m_critical_values[i][2],solo_pvalue,group_pvalue) ;
         else
            PrintFormat("%15s %12.4lf %12.4lf %12.4lf",string(k), m_critical_values[i][0], m_critical_values[i][1], m_critical_values[i][2]);
        }
     }

   return true;

  }

Метод mutualinfo(), который является частным методом в классе Cmrmr, отвечает за оценку взаимной информации между данным предиктором и целевой переменной с использованием метода адаптивного разбиения. Он возвращает вектор взаимных информационных приближений для выбранных предикторов, указанных в качестве индексов столбцов, в mutualinfo().

//+------------------------------------------------------------------+
//| continuous mutual infomation                                     |
//+------------------------------------------------------------------+
vector Cmrmr::mutualinfo(int which_size,int &which[],vector &targs)
  {
   vector res = vector::Zeros(m_vars);

   res.Fill(-DBL_MAX);

   vector vars;

   Capm mia(true,targs,m_chisquare_thresh);

   for(int i = 0; i<which_size; i++)
     {
      vars = m_preds.Col(which[i]);
      res[i] = mia.fit(vars);
     }

   return res;

  }

Метод StepWise() сохраняет вычисленные значения взаимной информации в векторе "релевантности", поскольку эти значения понадобятся в последующих итерациях для вычисления коэффициента избыточности. Признак с наибольшей взаимной информацией отбирается в качестве первого признака, а информация о нем сохраняется в переменных stepwise_crit и stepwise_ivar. Кроме того, инициализируется переменная sum_relevance для отслеживания совокупной релевантности выбранных признаков, которая позже будет использоваться для тестирования MCP.

Для начального повторения без перестановок алгоритм сохраняет исходные значения релевантности и критерия для первого отобранного признака. Он также инициализирует счетчики для общего теста модели на значимость и тестов на значимость отдельных признаков. Кроме того, выводится таблица, отображающая взаимную информацию каждого потенциального признака с целевой переменной, чтобы получить представление о первоначальном ранжировании признаков. Эта таблица помогает визуализировать первоначальный процесс отбора и служит референсом для сравнения с измененными значениями на этапе тестирования на значимость.

Если мы не находимся в начальном, неизмененном повторении, алгоритм переходит к выполнению двух тестов MCP, предназначенных для оценки статистической значимости набора признаков, отобранного алгоритмом MRMR. Для эффективного управления растущей избыточностью по мере добавления новых признаков к выбранному набору алгоритм использует массив `sum_redundancy`. Массив инициализируется нулями. Этот массив используется для хранения совокупной избыточности для каждого потенциального признака относительно уже отобранных признаков.

В начале, поскольку никакие признаки не выбраны, все значения в этом массиве обнуляются. Когда к отобранному набору добавляется новый признак, необходимо обновить ее избыточность с каждым из оставшихся потенциальных признаков. Избыточность признака относительно растущего отобранного набора определяется как средняя взаимная информация между потенциальным признаком и уже отобранными признаками. Вычисляется взаимная информация между новым признаком и каждым из оставшихся потенциальных признаков. Это количественно определяет, каким объемом информации новый признак делится с каждым из отобранных признаков. Значение избыточности для каждого оставшегося признака обновляется путем добавления значения взаимной информации между недавно добавленным признаком и существующим признаком.

После обновления значений избыточности алгоритм вычисляет критерий MRMR (Максимальная релевантность, минимальная избыточность) для каждого из оставшихся потенциальных признаков. Критерий MRMR рассчитывается путем вычитания средней избыточности каждого признака из его релевантности. Как только критерий MRMR будет рассчитан для всех оставшихся потенциальных признаков, будет выбран признак с наивысшим показателем MRMR. Алгоритм повторяет этот процесс итеративно, обновляя массив `sum_redundancy` и пересчитывая критерий MRMR каждый раз при добавлении нового признака. По мере того как отбирается все больше признаков, избыточность оставшихся признаков отобранным набором постоянно обновляется, гарантируя, что в процессе отбора признаков будут учитываться развивающиеся взаимосвязи между признаками.

Если это первоначальная, неизмененная репликация, исходные значения этих величин сохраняются для последующего сравнения во время тестов на перестановку. В противном случае счетчики для "групповых" и "одиночных" тестов на перестановку увеличиваются. После выполнения требуемого количества итераций алгоритм завершает работу, выводя на экран таблицу с краткой информацией об отобранных признаках, их релевантности и p-значениях, полученных в результате тестов на перестановку. После вызова и успешного выполнения `StepWise()` можно вызвать `GetCriticalValues()`, чтобы получить матрицу таблицы значений, отображаемую в подробном режиме. В то время как `GetSelectedVars()` извлекает индексы отобранных переменных (столбцов), соответствующие свойству `m_max_preds` класса.

//+------------------------------------------------------------------+
//| get the matrix of decision variables                             |
//| matrix rows arranged in terms of descending relevance of each    |
//| candidate predictor. Most relevant predictors critical values are|
//| in the first row, etc.                                           |
//| Column 0 is the calculated relevance                             |
//| Column 1 is the calculated redundancy                            |
//| Column 2 is the calculated Stepwise Mutual information           |
//| Column 3 and 4 will exist only if MCP test is specified          |
//| Column 3 is the Solo probability value                           |
//| Column 4 are the Group probability values                        |
//+------------------------------------------------------------------+
matrix Cmrmr::GetCriticalValues(void)
{
  return m_critical_values;
}
//+--------------------------------------------------------------------------+
//|get the indices of the most relevant predictors as determined by algorithm|
//+--------------------------------------------------------------------------+
bool Cmrmr::GetSelectedVars(ulong &output[])
{
 return (ArrayCopy(output,m_selected_vars)>0);
}

In the next section, we will see how the Cmrmr class is applied using synthetic and real-world datasets.

Примеры релевантности минус избыточность

Чтобы эффективно проиллюстрировать применение алгоритма MRMR, мы начнем с его применения к синтетическому набору данных. Набор данных состоит из 100 образцов и 10 потенциальных предикторных переменных (признаков). Целевая переменная Y вычисляется путем суммирования первых четырех столбцов матрицы, поэтому связь между целевой величиной и предикторами является детерминистской для этих столбцов. Однако столбцы с 4 по 7 не содержат никакой реальной информации о целевом объекте и служат несущественными признаками. Переменные в столбцах 8 и 9 представляют собой комбинации переменных с индексами столбцов {1,3} и {0,2} соответственно.

//---
   MathSrand(125);

   matrix rdata(100,10);

   rdata.Random(0.0,1.0);

   vector dep = rdata.Col(0) + rdata.Col(1) + rdata.Col(2) + rdata.Col(3);

   vector sum02 = rdata.Col(0) + rdata.Col(2);

   vector sum13 = rdata.Col(1) + rdata.Col(3);

   if(!rdata.Col(sum13,8) || !rdata.Col(sum02,9))
     {
      Print(" Failed column insertion ", GetLastError());
      return;
     }

Цель состоит в том, чтобы использовать алгоритм MRMR для определения столбцов как наиболее релевантных предикторов. Этот эксперимент реализован в скрипте RelevanceMinusRedundancy.mq5, который позволяет настраивать различные гиперпараметры, включая количество перестановок по методу Монте-Карло и максимальное количество отбираемых признаков. Скрипт настроен на запуск в подробном режиме, предоставляя подробные выходные данные в процессе отбора признаков.

//inputs
input int NumReplications = 100;
input int MaxPredictors = 10;
input double ChiSquareThreshold = 6.0;

Исходными данными скрипта RelevanceMinusRedundancy.mq5 являются индивидуальные значения взаимной информации между каждым потенциальным предиктором и целевой переменной. Взаимная информация расположена в порядке убывания, что позволяет легко идентифицировать наиболее релевантные признаки.

Variable         MI
      9       0.2675
      8       0.1696
      0       0.1662
      3       0.1336
      2       0.0823
      1       0.0645
      6       0.0000
      7       0.0000
      4       0.0000
      5       0.0000

Как и ожидалось, первым признаком, отобранным алгоритмом MRMR, является тот, который обладает наибольшей индивидуальной взаимной информацией с целевой переменной. В данном случае это тот, который указан в столбце с индексом 9.

Predictors so far   Relevance   Redundancy   Criterion
            9       0.2675       0.0000       0.2675

Вторая выбираемая функция - это переменная в столбце с индексом 8. Судя по результатам анализа, эта переменная имеет высокую релевантность, а также минимальную избыточность по сравнению с первой отобранной переменной.

Predictors so far   Relevance   Redundancy   Criterion
            9       0.2675       0.0000       0.2675
            8       0.1696       0.0000       0.1696

Третий отбираемый признак демонстрирует, как избыточность влияет на процесс отбора. Несмотря на то, что признаки в столбцах 0, 1, 2 и 3 непосредственно связаны с целевой переменной, их высокая избыточность по отношению к уже отобранным признакам (столбцы 8 и 9) снижает общий критерий MRMR. Следовательно, алгоритм отбирает менее непосредственно связанные функции, которые имеют меньшую избыточность по сравнению с существующим набором. Только после отбора нескольких несвязанных признаков алгоритм начнет определять приоритетность признаков в столбцах 0, 1, 2 и 3, поскольку их избыточность с отобранным набором уменьшается.

Additional candidates, in order of decreasing relevance minus redundancy
Variable     Relevance   Redundancy   Criterion
    5       0.0000       0.0000       0.0000
    4       0.0000       0.0000       0.0000
    7       0.0000       0.0000       0.0000
    6       0.0000       0.0000       0.0000
    1       0.0645       0.0600       0.0044
    0       0.1662       0.1351       0.0311
    2       0.0823       0.1426      -0.0603
    3       0.1336       0.1781      -0.0445

В то время как нерелевантные переменные в этом синтетическом примере имели нулевую релевантность, важно отметить, что в реальных сценариях нерелевантные признаки, тем не менее, могут демонстрировать некоторую степень взаимосвязи с целевой переменной, что приводит к ненулевой избыточности. В таких случаях алгоритм MRMR эффективно уравновешивает релевантность и избыточность, отдавая приоритет признакам, предоставляющим уникальную информацию.

Final predictors ||  Relevance ||  Redundancy || Criterion  || Solo pval|| Group pval
9                     0.2675        0.0000       0.2675       0.010       0.010
8                     0.1696        0.0000       0.1696       0.010       0.010
4                     0.0000        0.0000       0.0000       1.000       0.010
5                     0.0000        0.0000       0.0000       1.000       0.010
6                     0.0000        0.0000       0.0000       1.000       0.010
7                     0.0000        0.0000       0.0000       1.000       0.010
1                     0.0645        0.0738      -0.0093       0.010       0.010
0                     0.1662        0.1811      -0.0149       0.010       0.010
2                     0.0823        0.1077      -0.0254       0.010       0.010
3                     0.1336        0.1584      -0.0247       0.010       0.010

P-значения, полученные в результате тестов на перестановку Монте-Карло, дают ценную информацию о статистической значимости отобранных признаков. Более низкие p-значения указывают на более веские доказательства того, что эта функция действительно информативна, а не случайна. Для определения статистической значимости обычно используется пороговое p-значение, равное 0,05.

Далее приводится демонстрация алгоритма на более реалистичном наборе данных. Скрипт StepWiseFeatureSelectionByMutualInformation.mq5 обеспечивает практическое применение алгоритма MRMR по финансовым данным. Он основан на предпосылке прогнозирования будущей доходности определенного инструмента с использованием набора технических индикаторов. Скрипт рассчитывает четыре индикатора: Индекс денежного потока (MFI), Скользящая средняя (MA), индекс силы медведей и индекс силы быков за различные периоды ретроспективы. Эти индикаторы служат в качестве потенциальных признаков для процесса отбора признаков. Применяя алгоритм MRMR, скрипт стремится определить наиболее информативную комбинацию индикаторов, которые могут эффективно прогнозировать будущую доходность. Сюда включен отбор признаков, которые соответствуют целевой переменной (будущим доходам), при минимизации избыточности с другими отобранными признаками.

//+------------------------------------------------------------------+
//|                  StepWiseFeatureSelectionByMutualInformation.mq5 |
//|                                  Copyright 2024, MetaQuotes Ltd. |
//|                                             https://www.mql5.com |
//+------------------------------------------------------------------+
#property copyright "Copyright 2024, MetaQuotes Ltd."
#property link      "https://www.mql5.com"
#property version   "1.00"
#property script_show_inputs
#resource "\\Indicators\\LogReturns.ex5"
#include<mutualinfo.mqh>
#include<ErrorDescription.mqh>
//+------------------------------------------------------------------+
//|indicator type                                                    |
//+------------------------------------------------------------------+
enum SELECT_INDICATOR
  {
   MFI=0,//MFI
   MA,//MA
   BEARS,//BEARS
   BULLS//BULLS
  };
//--- input parameters
input int NumReplications = 100;
input int MaxPredictors = 10;
input double ChiSquareThreshold = 6.0;
input bool VerboseOutPut = false;
input uint     period_inc=2;//lookback increment
input uint     max_lookback=6;
input ENUM_MA_METHOD         AppliedMA = MODE_SMA;
input ENUM_APPLIED_PRICE     AppliedPrice = PRICE_CLOSE;
input datetime SampleStartDate=D'2019.12.31';
input datetime SampleStopDate=D'2022.12.31';
input string   SetSymbol="BTCUSD";
input ENUM_TIMEFRAMES SetTF = PERIOD_D1;
//----
string predictor_names[];                      // variable names
int size_sample,                      //training set size
    size_observations,                //size of of both training and testing sets combined
    maxperiod,                        //maximum lookback
    indicator_handle=INVALID_HANDLE;  //long moving average indicator handle
//---
vector indicator[];                   //indicator indicator values;
//---
matrix feature_matrix;          //full matrix of features;
//+------------------------------------------------------------------+
//| Script program start function                                    |
//+------------------------------------------------------------------+
void OnStart()
  {
//---get relative shift of sample set
   int samplestart,samplestop,num_features;
   samplestart=iBarShift(SetSymbol!=""?SetSymbol:NULL,SetTF,SampleStartDate);
   samplestop=iBarShift(SetSymbol!=""?SetSymbol:NULL,SetTF,SampleStopDate);
   num_features = int((max_lookback/period_inc)*4);
//---check for errors from ibarshift calls
   if(samplestart<0 || samplestop<0)
     {
      Print(ErrorDescription(GetLastError()));
      return;
     }
//---set the size of the sample sets
   size_observations=(samplestart - samplestop) + 1 ;
   maxperiod=int(max_lookback);
//---check for input errors
   if(size_observations<=0 || maxperiod<=0)
     {
      Print("Invalid inputs ");
      return;
     }
//---allocate memory
   if(ArrayResize(indicator,num_features+1)<0)
     {
      Print(ErrorDescription(GetLastError()));
      return;
     }
//----get the full collection of indicator values
   int period_len;
   int k=0;
//---
   for(SELECT_INDICATOR select_indicator = 0; select_indicator<4; select_indicator++)
     {
      for(int iperiod=0; iperiod<int((indicator.Size()-1)/4); iperiod++)
        {
         period_len=int((iperiod+1) * period_inc);
         int try
               =10;
         predictor_names.Push(EnumToString(select_indicator)+"_"+string(period_len));
         while(try)
           {
            switch(select_indicator)
              {
               case MFI:
                  indicator_handle=iMFI(SetSymbol!=""?SetSymbol:NULL,SetTF,period_len,VOLUME_TICK);
                  break;
               case MA:
                  indicator_handle=iMA(SetSymbol!=""?SetSymbol:NULL,SetTF,period_len,0,AppliedMA,AppliedPrice);
                  break;
               case BEARS:
                  indicator_handle=iBearsPower(SetSymbol!=""?SetSymbol:NULL,SetTF,period_len);
                  break;
               case BULLS:
                  indicator_handle=iBullsPower(SetSymbol!=""?SetSymbol:NULL,SetTF,period_len);
                  break;
              }

            if(indicator_handle==INVALID_HANDLE)
               try
                  --;
            else
               break;
           }

         if(indicator_handle==INVALID_HANDLE)
           {
            Print("Invalid indicator handle ",EnumToString(select_indicator)," ", GetLastError());
            return;
           }

         Comment("copying data to buffer for indicator ",period_len);
         try
               = 0;
         while(!indicator[k].CopyIndicatorBuffer(indicator_handle,0,samplestop,size_observations) && try
                  <10)
              {
               try
                  ++;
               Sleep(5000);
              }

         if(try
               <10)
               ++k;
         else
           {
            Print("error copying to indicator buffers ",GetLastError());
            Comment("");
            return;
           }

         if(indicator_handle!=INVALID_HANDLE && IndicatorRelease(indicator_handle))
            indicator_handle=INVALID_HANDLE;
        }
     }

//---resize matrix
   if(!feature_matrix.Resize(size_observations,indicator.Size()-1))
     {
      Print(__LINE__);
      Print(ErrorDescription(GetLastError()));
      Comment("");
      return;
     }
//---copy collected data to matrix
   for(ulong i = 0; i<feature_matrix.Cols(); i++)
      if(!feature_matrix.Col(indicator[i],i))
        {
         Print(__LINE__);
         Print(ErrorDescription(GetLastError()));
         Comment("");
         return;
        }
//---
   int try
         = 10;
   while(try
            >-1 && indicator_handle == INVALID_HANDLE)
        {
         indicator_handle=iCustom(SetSymbol!=""?SetSymbol:NULL,SetTF,"\\Indicators\\LogReturns",0,1,1);
         try
            --;
        }
//---
   if(try
         <0)
        {
         Print("Could not initialize returns indicator ");
         Print(ErrorDescription(GetLastError()));
         Comment("");
         return;
        }
   else
     {
      try
            = 10;
     }
//---
   while(try
            >-1 && !indicator[indicator.Size()-1].CopyIndicatorBuffer(indicator_handle,0,samplestop-1,size_observations))
        {
         Sleep(2000);
         try
            --;
        }
//---
   if(try
         <0)
        {
         Print("Could not collect returns indicator info ");
         Print(ErrorDescription(GetLastError()));
         Comment("");
         return;
        }
   else
     {
      IndicatorRelease(indicator_handle);
      Comment("");
     }
//--- display the dataset
   string console;
   for(uint i = 0; i<predictor_names.Size(); i++)
      console+=StringFormat(" %12s ",predictor_names[i]);
   console+=" NextBar Returns (Target) ";
   Print(console);
   for(ulong i = 0; i<feature_matrix.Rows(); i++)
     {
      console = "";
      for(ulong j = 0; j<feature_matrix.Cols(); j++)
         console += StringFormat(" %12.6lf ",feature_matrix[i][j]);
      console+=StringFormat(" %12.6lf ",indicator[indicator.Size()-1][i]);
      Print(console);
     }
//---
   Cmrmr mstep(NumReplications,MaxPredictors,ChiSquareThreshold,VerboseOutPut);
//---
   if(!mstep.StepWise(feature_matrix,indicator[indicator.Size()-1]))
      return;
//---
   Print(" Final predictor labels ");
   ulong variables[];
   if(mstep.GetSelectedVars(variables))
     {
      for(uint i = 0; i<variables.Size(); i++)
         Print(predictor_names[variables[i]]);
     }
   return;
  }
//+------------------------------------------------------------------+

Ниже приведен фрагмент набора данных индикаторов, собранных из терминала MetaTrader 5, а также целевая переменная в последнем столбце. Всего существует 12 потенциальных предикторов (индикаторов).

MFI_2         MFI_4         MFI_6          MA_2          MA_4          MA_6       BEARS_2       BEARS_4       BEARS_6       BULLS_2       BULLS_4       BULLS_6  NextBar Returns(Target)
0.000000     53.151442     46.921608   7213.654000   7279.552750   7260.007333    -76.371267   -101.867797   -107.113146     87.265733     61.769203     56.523854     -0.003564
0.000000     26.316595     34.451328   7180.962000   7244.558000   7255.420333    -41.795089    -70.889078    -83.462247     42.344911     13.250922      0.677753     -0.032447
0.000000      0.000000     36.720977   7053.740000   7133.697000   7204.281833   -127.731696   -210.813447   -251.244462    153.948304     70.866553     30.435538      0.054562

Ниже приведены результаты анализа предиктора с использованием параметров скрипта по умолчанию. Сначала посмотрим на таблицу значений взаимной информации. Здесь мы видим, что некоторые индикаторы имеют нулевую взаимную информацию с целью. Индикатор скользящей средней с длительностью периода равной 6, обладает наибольшей взаимной информацией по отношению к цели.

PS          Variable         MI
ND                 5       0.0308
LG                 4       0.0293
MJ                 3       0.0279
GM                 6       0.0227
JP                 9       0.0182
IS                 1       0.0165
MF                 8       0.0135
JI                 7       0.0126
HL                 0       0.0000
JO                 2       0.0000
GS                10       0.0000
HD                11       0.0000

Естественно, MA_6 - это первый отбор, сделанный алгоритмом. Поскольку первый отбор сделан, будет интересно посмотреть, есть ли какие-либо другие потенциальные предикторы, которые делятся с ней (MA_6) информацией. Это можно сделать, изучив приведенные ниже значения избыточности.

ME   Predictors so far   Relevance   Redundancy   Criterion
GH                 5       0.0308       0.0000       0.0308
II   Additional candidates, in order of decreasing relevance minus redundancy
FE          Variable     Relevance   Redundancy   Criterion
LE                 0       0.0000       0.0095      -0.0095
ED                 1       0.0165       0.0869      -0.0704
FD                 2       0.0000       0.1149      -0.1149
PE                11       0.0000       0.2768      -0.2768
ID                 8       0.0135       0.3251      -0.3116
NS                 7       0.0126       0.3492      -0.3366
CR                10       0.0000       0.3484      -0.3484
ES                 6       0.0227       0.4030      -0.3802
IS                 9       0.0182       0.4285      -0.4103
CR                 3       0.0279       2.5363      -2.5084
DR                 4       0.0293       3.0096      -2.9803

Следующим отбором для алгоритма будет MFI_2, поскольку он обладает наименьшей избыточностью среди всех оставшихся потенциальных предикторов, несмотря на нулевую релевантность к целевой переменной. Также обратите внимание, что предикторы в индексах столбцов 3 и 4 (MA_2 и MA_4 соответственно) имеют самые высокие уровни избыточности с MA_6. Это понятно, поскольку это один и тот же индикатор, рассчитанный с меньшими длинами окон. Ниже приведен предиктор, установленный после еще двух отборов, оба из которых также имеют нулевую релевантность к целевой переменной. Следует отметить главный предиктор в наборе кандидатов, которые не были отобраны. Алгоритм начинает рассматривать релевантные предикторы из-за разбавляющего эффекта, который последние два отбора оказали на отобранные предикторы в целом.

PQ   Predictors so far   Relevance   Redundancy   Criterion
FL                 5       0.0308       0.0000       0.0308
JK                 0       0.0000       0.0095      -0.0095
DK                 2       0.0000       0.1796      -0.1796
EE  Additional candidates, in order of decreasing relevance minus redundancy
JH         Variable     Relevance   Redundancy   Criterion
CH                6       0.0227       0.1469      -0.1241
PH                9       0.0182       0.1465      -0.1283
GI               10       0.0000       0.2024      -0.2024
PG                7       0.0126       0.2133      -0.2007
PF               11       0.0000       0.2308      -0.2308
GG                8       0.0135       0.2664      -0.2529
QG                1       0.0165       0.4679      -0.4514
KF                3       0.0279       0.8840      -0.8561
LF                4       0.0293       1.0372      -1.0079

Мы завершаем наш комментарий к результатам, забегая вперед и показывая окончательный набор предикторов, выбранных алгоритмом.

OQ  Final predictors ||  Relevance ||  Redundancy || Criterion  ||  Solo pval || Group pval
HN                5       0.0308       0.0000       0.0308       0.010       0.010
DJ                0       0.0000       0.0095      -0.0095       1.000       0.010
JG                2       0.0000       0.1796      -0.1796       1.000       0.010
HD                6       0.0227       0.1620      -0.1393       0.010       0.010
FQ                9       0.0182       0.2023      -0.1842       0.010       0.010
DL               11       0.0000       0.2798      -0.2798       1.000       0.010
KJ                1       0.0165       0.2932      -0.2767       0.010       0.010
OG                7       0.0126       0.3298      -0.3172       0.010       0.010
OE               10       0.0000       0.4151      -0.4151       1.000       0.010
JP                8       0.0135       0.4585      -0.4450       0.010       0.010
RN   Final predictor labels
NO  MA_6
DE  MFI_2
NN  MFI_6
KP  BEARS_2
DJ  BULLS_2
FL  BULLS_6
PQ  MFI_4
MK  BEARS_4
FM  BULLS_4
OF  BEARS_6

Заключение

В настоящей статье рассмотрено применение взаимной информации при отборе признаков. Основное внимание уделяется алгоритму MRMR (максимальная зависимость, максимальная релевантность, минимальная избыточность), предложенному авторами Пенгом (Peng), Лонгом (Long) и Дином (Ding). Мы начали с обсуждения оценки взаимной информации для непрерывных переменных, подчеркнув ее способность фиксировать сложные нелинейные зависимости между признаками и целевой переменной. В ходе нашего исследования алгоритма MRMR мы продемонстрировали, как он эффективно уравновешивает конкурирующие задачи по максимизации релевантности целевой переменной при минимизации избыточности уже отобранных признаков.

Путем итеративного добавления признаков на основе их показателей MRMR и оценки их статистической значимости с помощью тестов на перестановки Монте-Карло алгоритм гарантирует, что отобранный набор признаков обеспечивает наиболее информативные и уникальные предикторы для модели. Синтетические и реальные примеры продемонстрировали практическую полезность алгоритма MRMR в реальных сценариях обработки данных, где нерелевантные признаки и мультиколлинеарность часто усложняют отбор признаков. Способность алгоритма MRMR определять приоритетность релевантных признаков и избегать избыточных гарантирует, что он остается ценным инструментом в процессе отбора признаков, особенно в тех случаях, когда взаимосвязи между переменными являются сложными и нелинейными.

Весь код, на который мы ссылаемся в настоящей статье, прилагается ниже. В следующей таблице описаны все файлы исходного кода, которые прилагаются к статье.

Имя файла	Описание
MQL5/include/np.mqh	Заголовочный файл с матричными и векторными служебными функциями.
MQL5/include/mutualinfo.mqh	Заголовочный файл, содержащий определение класса Capm, который реализует адаптивный метод разбиения для оценки взаимной информации непрерывных переменных. И определение класса Cmrmr, реализующего поэтапный отбор признаков на основе взаимной информации.
MQL5/scripts/RelevanceMinusRedundancy.mq5	Скрипт, демонстрирующий использование класса Cmrmr с помощью синтетических данных
MQL5/scripts/StepWiseFeatureSelectionByMutualInformation.mq5	Другой сценарий также демонстрирует применение класса Cmrmr на более реалистичном наборе данных.