Full Text

Введение В настоящее время во всем мире чрезвычайно актуальны проблемы анализа и эффективной обработки больших объемов статистических данных (Big Data) [1; 2], накопленных в разных сферах человеческой деятельности. К одной из техник Big Data относится регрессионный анализ [3; 4]. С помощью регрессионных моделей успешно решается множество прикладных задач. Например, в [5] с помощью регрессионного анализа определены наиболее значимые факторы риска развития микрососудистых осложнений сахарного диабета второго типа, в [6] - получено значение предельно допустимого состава смесевого топлива для применения в тракторных дизелях. Традиционно много регрессионных моделей строится на основе экономической статистики (см, например, [7]). Такие модели называются эконометрическими. При всем этом одной из главных проблем, связанных с построением регрессионной модели, считается выбор ее спецификации [8], т.е. общего вида модели, в том числе состава и формы входящих в нее связей. Выбор состава входящих в уравнение регрессии переменных формализуется в виде так называемой задачи отбора наиболее информативных регрессоров (ОИР) [9] на основе некоторого критерия качества. Для решения проблемы ОИР существует множество различных методов, описание большинства из которых можно найти в [10]. Единственным из них, который гарантирует оптимальное решение задачи, считается метод «всех регрессий», реализация которого предполагает перебор всех возможных комбинаций переменных в модели. Но этот метод и самый трудоемкий из всех существующих, поэтому неэффективен при решении задач Big Data при большом числе переменных. Кроме того, построенную в результате реализации метода «всех регрессий» регрессионную модель бывает затруднительно или вовсе невозможно интерпретировать. Например, знаки ее оценок могут противоречить содержательному смыслу факторов, в ней может присутствовать мультиколлинеарность [11] и т.д. На сегодняшний день проблеме построения интерпретируемых моделей машинного обучения [12; 13] в научной литературе уделяется значительное внимание. Интерпретация построенной модели способствует выявлению и устранению ее «уязвимых» мест, что повышает доверие к модели у экспертов из данной предметной области. К сожалению, точного определения интерпретируемости модели не существует. Успешная попытка дать такое определение для регрессионных моделей была предпринята в работе [14], в которой рассмотрено необходимое и достаточное условие однозначной интерпретируемости регрессии. Однако для построения таких моделей в [14] предложен все тот же малоэффективный метод «всех регрессий». За последние десятилетия была существенно развита технология решения задач частично-целочисленного программирования, поэтому появился новый более эффективный метод нахождения оптимального решения задачи ОИР. В зарубежной литературе задача ОИР при оценивании линейных регрессий с помощью метода наименьших квадратов (МНК) формализуется в виде задачи частично-булевого квадратичного программирования (ЧБКП) (см, например, [15-17]), число ограничений которой зависит от объема выборки. В работе [18] автору впервые удалось свести такую задачу к задаче частично-булевого линейного программирования (ЧБЛП), число ограничений которой не зависит от объема выборки. В дальнейшем эта задача эволюционировала, в ней появлялись новые ограничения. Например, в [19] была сформулирована задача ЧБЛП, решение которой приводит к построению линейной регрессии с оптимальным числом объясняющих переменных, абсолютные вклады которых в общую детерминацию не меньше, чем число . В той же работе [19] было предложено определение вполне интерпретируемой линейной регрессии (ВИЛинР). В исследовании [20] определение было уточнено за счет ограничений на коэффициенты интеркорреляций. И в той же работе [20] проведены многочисленные вычислительные эксперименты, подтверждающие высокую эффективность предложенного метода для обработки Big Data. Например, удалось эффективно обработать выборку из 515345 наблюдений. Целью данной работы является обобщение предложенного в работах [18-20] математического аппарата, эффективно справляющегося с обработкой Big Data, для построения вполне интерпретируемых квазилинейных регрессионных моделей. 1. Технология построения вполне интерпретируемых квазилинейных регрессий Введем в рассмотрение элементарную квазилинейную регрессию (КЛинР) [9]: , , (1) где - объем выборки; - число объясняющих переменных; - i-е значение объясняемой переменной; - i-е значение j-й объясняющей переменной; , , , - неизвестные параметры; - i-я ошибка аппроксимации; - k-e элементарное преобразование j-й переменной, выбранное из набора . Линейная регрессия является частным случаем КЛинР (1). Достоинство КЛинР (1) в том, что они линейны по параметрам, поэтому их можно оценивать с помощью МНК. А недостаток заключается в проблематичности их интерпретации. Например, если оценена модель , то затруднительно каким-либо образом объяснить оценку 17 при переменной . В некоторых источниках, например, в [12], предлагается интерпретировать такую модель следующим образом: если квадрат переменной увеличится на единицу, то значение переменной увеличится на 17 единиц. На наш взгляд, такой подход к интерпретации не совсем корректен, поскольку не объясняется прямое влияние именно переменной на . Введем определение вполне интерпретируемой квазилинейной регрессии (ВИКЛинР). Для удобства запишем КЛинР (1) в виде , , (2) где , , , . Определение. КЛинР (2), оцененная с помощью МНК, называется вполне интерпретируемой, если: 1) каждая преобразованная объясняющая переменная входит в модель не более одного раза; 2) знаки всех коэффициентов корреляции , , , удовлетворяют содержательному смыслу решаемой задачи; 3) знаки всех оценок согласуются со знаками соответствующих коэффициентов корреляции , т.е. , , ; 4) все абсолютные вклады переменных в общую детерминацию удовлетворяют неравенствам , , ; 5) все коэффициенты интеркорреляций , , , , ; 6) коэффициент детерминации модели . Кратко прокомментируем данное определение. В нем условие № 2 означает, что еще до оценивания КЛинР (2) требуется анализировать соответствие знаков коэффициентов корреляции преобразованных объясняющих переменных с y содержательному смыслу решаемой задачи. В случае выявления противоречий необходимо либо увеличить объем выборки, либо исключить соответствующую преобразованную переменную из рассмотрения. На этом этапе желательно привлекать экспертов в данной предметной области. Экспертам следует помнить, что элементарные преобразования переменных могу искажать направление влияния объясняющих переменных на y. Стандартизованная регрессия для (2) записывается в виде , , где , , , , ; , , - неизвестные стандартизованные коэффициенты; , - новые ошибки аппроксимации. Если выполняется условие № 3, то, как отмечено в [19], становятся справедливы формулы для абсолютных вкладов переменных в общую детерминацию : , , . По этим критериям можно делать выводы о степени влияния преобразованных переменных на y. Условие № 4 означает, что каждая такая переменная должна вносить вклад в детерминацию не менее 0,01. Условие № 5 означает ограничение на эффект мультиколлинеарности. Уровень 0,4, при котором модель считается приближенно интерпретируемой, взят из [14]. Также эту границу можно вычислить, например, по формуле , где - заданный уровень значимости. В этом случае все коэффициенты интеркорреляций будут незначимы по t-критерию Стьюдента для уровня . Если коэффициент детерминации регрессии меньше 0,8, то вряд ли можно отнести такую модель ко вполне интерпретируемой, поэтому введено условие № 6. Рассмотрим решение проблемы интерпретации квазилинейных регрессий. Предположим, что оцененная ВИКЛинР имеет вид: , (3) где - элемент вектора , указывающий номер выбранного преобразования для j-й объясняющей переменной. Если все функции из набора непрерывны и монотонны на отрезках , где , , , то становятся справедливы следующие формулы критериев нелинейности [21] для каждой объясняющей переменной: , . (4) Каждый из критериев нелинейности (4) принимает значения от 0 до 1. Если , то преобразование j-й объясняющей переменной является линейным. А если , то преобразование j-й объясняющей переменной в значительной степени нелинейно. Сформулируем алгоритм интерпретации влияния j-й объясняющей переменной на y в квазилинейной регрессии (3) с монотонными функциями. 1. Вычислить критерий нелинейности для j-й переменной по формуле (4). 2. Если , то вместо объясняется коэффициент 3. Если , то сначала по формуле Лагранжа находится точка , а потом объясняются коэффициенты и на отрезках и соответственно. По аналогии с задачей построения вполне интерпретируемой линейной регрессии, рассмотренной, например, в [20], сформулируем следующую задачу ЧБЛП построения ВИКЛинР (2): , (5) , , (6) , (7) , (8) (9) , , , (10) , , , , , (11) , (12) где символами обозначены коэффициенты корреляции между переменными и M - большое положительное число; - бинарная переменная, принимающая значение 1, если j-я объясняющая переменная с k-м преобразованием входит в регрессию, и 0, если не входит; - верхняя граница интеркорреляций; - нижняя граница абсолютных вкладов переменных. Точность коэффициентов корреляций в этой задаче должна быть не менее 12 знаков. В задаче (5)-(12) ограничения (12) обеспечивают выполнение условия № 1 в предложенном выше определении, ограничения (7), (8) - условия № 3, ограничения (10) при - условия № 4, ограничения (11) при - условия № 5. Условие № 2 проверяется до решения этой задачи, а условие № 6 - после решения. Ограничения (6) предназначены для включения/исключения уравнений в линейную систему, с помощью которой находятся МНК-оценки. Целевая функция (5) означает выбор регрессии с максимальным значением коэффициент детерминации . Таким образом, решение задачи ЧБЛП (5)-(12) приводит к построению ВИКЛинР с оптимальным по критерию количеством регрессоров, в которой , вклады , а интеркорреляции . Если исследователь желает контролировать в задаче (5)-(12) коэффициенты вздутия дисперсии VIF, то ее необходимо дополнить линейными ограничениями из работы [22], а если контролировать значимость оценок по t-критерию Стьюдента, то ограничениями из работы [23]. Выбирать допустимые границы параметра M в задаче (5)-(12) можно по аналогии с процедурой выбора таких границ в задаче построения вполне интерпретируемой линейной регрессии [19]. Для этого ограничения (6)-(8) следует заменить следующими выражениями: , , , (13) , , (14) , , (15) где , , , - значение коэффициента детерминации регрессии, построенной со всеми регрессорами. Для нахождения нужно решить серию из задач линейного программирования с целевыми функциями и линейными ограничениями , , , , , , Значения параметров берутся как значения целевых функций этих задач в точках оптимума. Для нахождения чисел нужно решить серию из задач линейного программирования с теми же линейными ограничениями, но с целевыми функциями . Оптимальные решения задач (5)-(12) и (5), (9)-(15) не отличаются. 2. Моделирование прочности бетона на сжатие Для демонстрации работоспособности предложенного математического аппарата решалась задача обработки данных из строительной области. Как известно, одно из лидирующих мест среди строительных материалов во всем мире занимает бетон. Ключевой проблемой при его изготовлении считается подбор состава компонент смеси и их пропорций так, чтобы обеспечить максимальную его долговечность. Для решения этой проблемы, как и многих других, могут применяться методы математического моделирования. Важный показатель, которым характеризуется бетон, - прочность на сжатие. Математическому моделированию прочности бетона посвящено множество научных работ. Из них хотелось бы выделить статью [24], в которой приведен весьма солидный обзор существующих математических моделей. Их анализ позволил сделать вывод, что во всех моделях фигурирует так называемое водоцементное отношение по объему, т.е. отношение массы воды к массе цемента. Так, например, более 100 лет назад Д. Абрамс сформулировал правило [25]: , где - прочность бетона, , - некоторые константы, - водоцементное отношение. Как видно, чем выше водоцементное отношение, тем ниже прочность на сжатие. А, например, в работе М. Боломея [26] предложена следующая формула прочности бетона: , где , - некоторые константы, - цементно-водное отношение. Как видно, оно влияет на в противоположном направлении. Однако на прочность бетона может влиять не только водоцементное отношение, но и качество заполнителей, добавление пластификатора, метод замеса, температура воздуха и др. В связи с этим вызывает интерес статистическое моделирование прочности бетона на сжатие. Построению статистических моделей для прогнозирования прочности бетона посвящена работа [24]. В ней использованы статистические данные из хранилища [27] по следующим переменным: - прочность бетона на сжатие (МПа); - цемент (кг/м3 смеси); - шлак доменной печи (кг/м3 смеси); - зола (кг/м3 смеси); - вода (кг/м3 смеси); - пластификатор (кг/м3 смеси); - щебень (кг/м3 смеси); - песок (кг/м3 смеси); - возраст бетона (дни). Как отмечено в [24], выборка из 920 наблюдений была разделена на обучающую, объемом 700, и контрольную, объемом 220. На сегодняшний день на сайте [27] объем этой выборки составляет уже 1030 наблюдений. В [24] авторами предложены следующие спецификации нелинейных регрессионных моделей: (16) (17) (18) Коэффициенты детерминации оцененные с помощью МНК регрессий (16)-(18), составили 0,721669, 0,772118 и 0, 778409 соответственно. С использованием эконометрического пакета Gretl модели (16)-(18), построенные по выборке объема 700, были переоценены по новой выборке объема 1030. Новые значения составили 0,713829, 0,754206 и 0,762094. Как видно, за счет новой информации качество аппроксимации моделей несколько ухудшилось. Лучшая из этих трех регрессий, для которой может быть использована для прогнозирования. Но ее никак нельзя считать вполне интерпретируемой, поскольку она не удовлетворяет ни одному из шести условий, рассмотренных выше. Так, во-первых, некоторые преобразованные переменные входят в наилучшую модель более одного раза, например, ; во-вторых, не проведен анализ согласованности знаков преобразованных переменных, например, , содержательному смыслу задачи; в-третьих, знаки некоторых оценок не согласуются со знаками соответствующих коэффициентов корреляции, например, оценка при переменной составляет 7,83, а коэффициент корреляции между и равен -0,18. В-четвертых, не работают формулы для абсолютных вкладов переменных в общую детерминацию; в-пятых, интеркорреляции между некоторыми преобразованными переменными близки к единице по абсолютной величине, например, между и такой коэффициент равен 0,941, что свидетельствует о наличии мультиколлинеарности; в-шестых, коэффициент детерминации модели оказался ниже 0,8. Была поставлена цель по данным, расположенным на сайте [27], построить ВИКЛинР прочности бетона на основе решения задачи ЧБЛП (5), (9)-(15). Для этого была использована программа «ВИнтер-2», подробное описание которой можно найти в [28]. «ВИнтер-2» позволяет в зависимости от выбранных пользователем начальных параметров автоматически формулировать для решателя LPSolve задачи ЧБЛП для построения, в частности, ВИКЛинР. Поскольку цементно-водное отношение фигурирует в большинстве известных математических моделей, было принято решение использовать одну переменную вместо двух переменных и . Предварительно был проведен анализ соответствия знаков коэффициентов корреляции объясняющих переменных с y содержательному смыслу задачи. Оказалось, что , , , , , , . Экспертами было принято решение о согласованности знаков абсолютно всех этих коэффициентов смыслу задачу. Поскольку переменные , , , содержат нулевые значения, что не позволяет использовать такие элементарные преобразования, как , и др., было решено увеличить их на 1. Для формирования задачи ЧБЛП (5), (9)-(15) для решателя LPSolve в «ВИнтер-2» были выбраны следующие начальные параметры: 1) объясняющие переменные - , , , , , , ; 2) элементарные преобразования - , , , , , , , , ; 3) точность вещественных чисел - 12 знаков; 4) наименьший абсолютный вклад ; 5) наибольшая интеркорреляция Запустив «ВИнтер-2», сначала автоматически были сформированы все возможные комбинации преобразованных переменных, общее число которых составило 63. Затем «ВИнтер-2» сама проверила согласованность знаков коэффициентов корреляции преобразованных переменных с y содержательному смыслу задачи. В итоге все 63 переменных проверку прошли. После чего была автоматически сформирована задача ЧБЛП (5), (9)-(15). Значения больших чисел M в ограничениях (13)-(15) программа также осуществила самостоятельно, для чего была решена необходимая серия задач линейного программирования. Далее сформированная задача ЧБЛП вручную была перенесена в решатель LPSolve, в котором была построена следующая КЛинР: . (19) В уравнении (19) в скобках под коэффициентами указаны наблюдаемые значения t-критерия Стьюдента, а над ними - абсолютные вклады переменных в общую детерминацию. Для КЛинР (19) , т.е. по этому критерию качество ее аппроксимации выше, чем у оцененной модели (18) на 0,0561. При этом спецификация регрессии (19), у которой 5 неизвестных параметров, проще спецификации (18) с восьмью параметрами. Модель (19) удовлетворяет всем шести вышеперечисленным требованиям, поэтому ее справедливо можно считать ВИКЛинР. Дополнительно к этому все ее оценки значимы по t-критерию Стьюдента для уровня значимости , а коэффициенты вздутия дисперсии VIF, равные 1,094, 1,085, 1,007 и 1,004 соответственно, говорят об отсутствии мультиколлинеарности. Для интерпретации ВИКЛинР (19) сначала были найдены значения критериев нелинейности преобразованных переменных по формулам (4): , , , Затем проводилась интерпретация влияния каждой объясняющей переменной на y на основе предложенного выше алгоритма. Для наглядности на рисунке сплошными линиями выделены входящие в уравнение (19) нелинейные функции, а пунктиром - заменяющие их при интерпретации прямые и ломаные. 1. Поскольку , то разница между кривой и прямой линией несущественна (рисунок). Тогда справедлива следующая интерпретация: с увеличением цементно-водного отношения на одну единицу (при неизменных значениях остальных переменных) прочность бетона на сжатие y увеличивается примерно на 16,368 МПа. а б в г Рис. Прямые и ломаные, заменяющие нелинейные функции (авторские результаты) 2. Так как , т.е. графики кривой и заменяющей ее прямой совпадают (рисунок, б), то можно сделать следующий вывод: с увеличением шлака доменной печи на 1 кг/м3 смеси (при неизменных значениях остальных переменных) прочность бетона на сжатие y увеличивается в среднем на 0,064 МПа. Полученный результат о влиянии шлака доменной печи на прочность бетона подтверждают, например, исследования [29; 30]. 3. Ввиду того, что , кривая заменена ломаной, у которой уравнение первого звена при , а уравнение второго звена при (рисунок, в). Тогда можно дать следующую интерпретацию: если , то с увеличением пластификатора на 1 кг/м3 смеси (при неизменных значениях остальных переменных) прочность бетона на сжатие y увеличивается в среднем на 2,904 МПа, а если , то на 0,0196 МПа. 4. На том основании, что кривая заменена ломаной, у которой уравнение первого звена при , а уравнение второго звена при (рисунок, г). Тогда справедлива следующая интерпретация: если то с увеличением возраста бетона на один день (при неизменных значениях остальных переменных) прочность бетона на сжатие y увеличивается в среднем на 0,57 МПа, а если , то на 0,0492 МПа. Заключение 1. Достоинство сформулированной задачи ЧБЛП (5)-(12) построения ВИКЛинР состоит в том, что число ограничений в ней не зависит от объема выборки. Так, если бы в данных, по которым моделировалась прочность бетона, было бы не 1030 наблюдений, а, например, более миллиона, то можно предположить, что время построения модели по-прежнему осталось бы в допустимых пределах. Тем самым предложенная в статье технология построения ВИКЛинР обладает потенциалом в области решения задач обработки Big Data. 2. Недостаток задачи ЧБЛП (5)-(12) состоит в том, что не до конца ясно, как выбирать большие числа M в ограничениях (6)-(8). Например, при моделировании прочности бетона на компьютере с процессором AMD Ryzen 3 4300U (2,70 GHz) и объемом оперативной памяти 16 ГБ на решение задачи при использовании ограничений (13)-(15) вместо (6)-(8) ушло 611 секунд. Было установлено, что если в этой задаче выбрать в ограничениях (6)-(8) число M = 500, то время ее решения составит уже 302 секунды, т.е. примерно в 2 раза меньше. Это классическая проблема, решение которой ищут как отечественные, так и зарубежные ученые. 3. Несмотря на то, что структура построенной ВИКЛинР (19) была выбрана автоматически, влияние объясняющих переменных на y в ней согласуется как с содержательным смыслом задачи, так и с другими существующими математическими моделями. К сожалению, пока в программе «ВИнтер-2» реализовано только 9 элементарных функций, поэтому за счет использования более широкого набора преобразований качество аппроксимации модели (19) может быть серьезно улучшено. Тем не менее построенную регрессию (19) можно использовать для прогнозирования прочности бетона в зависимости от входящих в нее факторов. 4. С помощью предложенной технологии можно успешно решать задачи выявления скрытых нелинейных зависимостей в больших наборах данных абсолютно в любой предметной области.

About the authors

M. P. Bazilevskiy

Irkutsk State Transport University

References

Statistics

Views

Abstract - 37

PDF (Russian) - 24

Refbacks

• There are currently no refbacks.