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

Full Text

ВВЕДЕНИЕ

Для прогнозирования физико-химических свойств различных веществ по их химической структуре традиционно используют аддитивные методы [1–6]. В рамках аддитивных методов (методов групповых вкладов) выдвигается предположение о том, что каждый фрагмент химической структуры вещества вносит свой вклад в значение того или иного физико-химического свойства [2, 3]. Первые попытки применения аддитивных методов к полимерным объектам были сделаны в работах [7–10], а значительный прогресс был достигнут в работе Ван-Кревелена [11]: в ней рассматривались возможности расчета и прогнозирования таких свойств полимеров, как плотность, коэффициент термического расширения, теплоемкость, механические свойства и другие. Данное направление впоследствии было расширено в работах J. Bicerano [12], А.А. Аскадского [13, 14], M. Salame [15–17]. Еще в начале 1990-х годов точность прогнозирования транспортных свойств полимеров аддитивными методами была неудовлетворительной [18], однако в конце того же десятилетия появились методы групповых вкладов с достаточно высокой прогностической способностью [19–21]. Позднее вышел цикл работ, посвященных различным аддитивным подходам к прогнозированию параметров транспорта газов в полимерах [22, 23]. Точность прогнозирования газотранспортных свойств методами, представленными в публикациях [22, 23], развивали в дальнейшем как в работах отечественной группы Ю.П. Ямпольского [24–27], так и в работах зарубежной группы под руководством S.B. Tantekin-Ersolmaz [28].

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

В настоящее время все большую популярность обретают нейросетевые модели прогнозирования [29–33]. Они используются для решения как научных, так и прикладных задач. Ранее для прогнозирования ряда физико-химических свойств полимеров использовали методы машинного обучения, а именно, метод, использующий на входе уникальный набор структурных дескрипторов¹ для описания каждого полимера (вместо структурной формулы) с обучением посредством методов гребневой регрессии (KRR²) или гауссовой регрессии (GPR³) на основе радиальных базисных функций (RBF⁴) и результаты работы были реализованы в формате онлайн-платформы прогнозирования Polymer Genome [30, 31]. J.W. Barnett и соавт. [32] схожим образом решали задачу поиска структур полимеров, имеющих заданные характеристики газопроницаемости, а в работе [34] метод, использующий уникальный набор дескрипторов, был объединен с формированием синтетической базы полимерных структур (для увеличения количества данных) методом молекулярно-динамического моделирования для прогнозирования доли свободного объема в полимерах. В работе группы С.В. Люлина [29] применяли сверточную нейронную сеть для прогнозирования температуры стеклования полимеров на основе сгенерированной базы данных для полиимидов, сформированной для получения массива данных, достаточного для проведения машинного обучения. Детальный обзор современных нейросетевых архитектур, связанных с химическими приложениями, приведен в работе P. Reiser и соавт. [35].

Лимитирующим фактором для применения нейронных сетей для прогнозирования газотранспортных свойств полимеров является ограниченность массива доступных экспериментальных данных. Как было сказано выше, в ряде работ [29, 34] используют созданные тем или иным расчетным методом массивы данных, чтобы проведение машинного обучения с применением современных нейронных моделей было реализуемо и обосновано. Поскольку в линейных инкрементальных моделях связи структура–свойство обычно рассматривают систему один газ–различные полимеры, количество экспериментальных данных (уравнений) для разных газов может отличаться весьма существенно и, как правило, не превышает 1500, что крайне мало для построения любых нейросетевых моделей. Так как для системы один полимер–разные газы известны хорошо работающие линейные зависимости логарифма коэффициента диффузии от квадрата кинетического диаметра газа d²_eff и логарифма коэффициента растворимости от параметра потенциала Леннард-Джонса (ε / k)_eff [36–38], существует принципиальная возможность построения обобщенной модели для всех полимеров и газов, т.е. объединения инкрементов, получаемых аддитивными методами, со свойствами газов. Такая обобщенная модель позволит существенно расширить массив экспериментальных данных для прогнозирования.

Однако этого все еще недостаточно для реализации машинного обучения на современных нейросетевых моделях. Тем не менее классические нейросетевые модели на таких массивах данных уже работоспособны. В случае, если на выходе нейронной сети присутствует только одно числовое значение, классические нейросетевые модели по своей сути ничем не отличаются от обычных нелинейных функций многих переменных, со всеми присущими им достоинствами (гибкость и возможность описать почти любые непротиворечивые данные) и недостатками (такими как сложность подбора параметров модели, а при отсутствии теоретического базиса, как правило, такие модели имеют очень низкую предсказательную способность или – в терминах нейронных сетей – являются переобученными). Хорошая предсказательная способность может быть достигнута только на огромных массивах экспериментальных данных (50000 экспериментальных точек и более [39]), что существенно ограничивает применимость нейронных сетей в современных научных изысканиях. В связи с этим в настоящей работе применяются классические нейронные сети как вспомогательный инструмент для уменьшения ошибки прогнозирования с применением обобщенной модели, о которой было сказано выше.

ЭКСПЕРИМЕНТАЛЬНАЯ ЧАСТЬ

Краткое рассмотрение применяемых нейросетевых моделей

Рассмотрим классические полносвязные нейронные сети (рис. 1). Любая нейронная сеть разделяется структурно на несколько слоев: входной слой, один или несколько внутренних слоев и выходной слой. Каждый слой содержит в себе один или несколько нейронов. Нейрон – это некая функция, которая суммирует все, что приходит на вход по межнейронным связям, полученное значение преобразует функцией активации f и передает на выход нейрона. Количество нейронов во входном слое обычно соответствует размерности подаваемого на вход вектора, хотя в некоторых случаях, добавляют еще один нейрон смещения (bias neuron). Нейроны смещения не имеют входа. Нейроны входного слоя обычно имеют рефлексивную функцию активации (т.е. никоим образом не преобразуют входное значение, а сразу передают его на выход нейрона). Каждый нейрон внутреннего слоя на вход принимает несколько значений от нейронов предыдущего слоя, которые вначале суммируются, а затем полученное значение подвергается преобразованию с помощью функции активации f (как правило нелинейной) и подается на выход нейрона. В полносвязных сетях каждый выход нейрона связан со входами всех нейронов следующего слоя. Каждая межнейронная связь (синапс) имеет вес w_i, на который умножается значение с выхода нейрона, перед поступлением на вход нейрона следующего слоя. Проще говоря, нейронная сеть попеременно чередует линейное и нелинейное преобразование. К настоящему времени описано использование более 20 различных нелинейных функций активации нейронов [39], однако наиболее часто применяют либо логистическую функцию активации (она удобна в теоретических исследованиях, но слишком медленно рассчитывается для использования в реальных задачах), либо различные модификации функции ReLU [40].

 

Рис. 1. Схематическое изображение полносвязной классической нейронной сети (MLP) и расчетной схемы нейрона (для простоты изображена нейронная сеть, получающая на вход всего три значения, в реальности их бывает значительно больше, иногда на несколько порядков). MLP – Multilayer Perceptron (в русскоязычной литературе: многослойный перцептрон Розенблатта). Цветные рисунки можно посмотреть в электронной версии.

 

Все описанное выше относится только к внутренней структуре нейронной сети. Классическая нейронная сеть – в нашем случае это нелинейная функция, так как на выходе у нее будет только одно числовое значение. Перед описанием остаточных нейронных сетей стоит рассмотреть, что собой представляет в математическом плане нейронная сеть, если все функции активации заменить на рефлексивные (см. описание к рис. 4). В этом нам поможет рис. 2. Вывод о линейности полученной модели не удивителен, потому что из математики давно известно, что суперпозиция линейных преобразований – это линейное преобразование.

 

Рис. 2. Математическая тождественность линейной модели и нейронной сети с рефлексивными функциями активации (для простоты рассматривается нейронная сеть только с тремя входными значениями).

 

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

Понятие остаточных нейронных сетей обязано своим появлением современной сверточной архитектуре ResNet, появившейся в 2015 году благодаря работе K. He c сотрудниками [41]. Основная идея этой архитектуры заключалась в добавлении к каждому слою нейронов обходного пути, иными словами, при необходимости слой нейронов просто передает значения своих входов на выходы (по сути, используя рефлексивную функцию активации y (x) = x). Эта модификация позволила успешно обучать очень глубокие сети размером до 1000 слоев. До появления ResNet многие исследователи считали обучение нейронных сетей размером более 100 слоев практически невозможным из-за так называемой “проблемы исчезающих градиентов” [39]. В настоящее время почти все современные архитектуры нейронных сетей в том или ином виде реализуют идею обходного пути, однако обычно это довольно глубокие сверточные сети. Мы же попытались реализовать эту идею на небольшой классической полносвязной нейронной сети прямого распространения (многослойном перцептроне Розенблатта [42], модифицированном в соответствии современными представлениями).

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

$Y=\alpha F_{lin}(x_1, \dots , x_n)+(1 – \alpha )F_{nonlin}(x_1, \dots , x_n)$, (1)

где F_lin (x₁, …, x_n) – линейная функция, F_nonlin (x₁, …, x_n) – нелинейная функция, α ∈ [0, 1], причем, если линейная модель функционирует удовлетворительно, то следует ожидать того, что значение α будет близко к единице, иными словами, вклад нелинейной составляющей будет мал.

Описание обобщенной линейной модели

Построение обобщенной линейной модели проводили с использованием базы данных ИНХС “Газоразделительные параметры стеклообразных полимеров” [43]. Данные по экспериментальным логарифмам коэффициентов проницаемостей газов при 35 °С для каждой из 3052 систем полимер–газ были дополнены 49 структурными дескрипторами полимера, полученными по алгоритмам ранее использовавшимся в модели атомных вкладов [25] с помощью программы RIADA и двумя наиболее подходящими, по нашему мнению, свойствами газа. В качестве свойств газа использовали эффективный параметр потенциала Леннард-Джонса (ε/k)_eff и квадрат кинетического диаметра газа d²_eff, приведенные в работе [44]. Итак, для построения модели составляли переопределенную линейную систему уравнений, в каждом уравнении которой левая часть представляла собой сумму произведений числовых значений структурных дескрипторов A₁–A₄₉ на соответствующие переменные x₁–x₄₉ для полимера и сумму двух произведений A₅₀x₅₀ и A₅₁x₅₁, относящихся к системе полимер–газ (см. уравнения (2), (3)). Правая часть каждого уравнения системы представляла собой десятичный логарифм коэффициента проницаемости газа, деленный на квадрат кинетического диаметра газа, (log₁₀P35)/d_eff². Решив эту систему уравнений, получим одновременно и параметры модели x₁–x₄₉ для структурных дескрипторов, и оценку констант Теплякова–Мирса – (K₁ + K₃) и K₄ (которые также являются параметрами модели).

$\frac{{\log }_{10}P35}{d_{eff}^2}=\sum _{i=1}^{51}{x}_i{A}_i$, (2)

$\frac{{\log }_{10}P35}{d_{eff}^2}=\left(\sum _{i=1}^{49}{x}_i{A}_i\right)+K_4\frac{{\left(\epsilon /k\right)}_{eff}}{d_{eff}^2}+\frac{\left(K_1+K_3\right)}{d_{eff}^2}$, (3)

где [A₁, …, A₄₉] – структурные дескрипторы; $A_{50}=\frac{{\left(\epsilon /k\right)}_{eff}}{d_{eff}^2}$; $A_{51}=\frac{1}{d_{eff}^2}$; [x₁, … , x₅₁] – переменные (x₅₀ ≡ K₄, x₅₁ ≡ K₁ + K₃).

На этапе построения модели x_i – это переменные. После построения модели полученные фиксированные числовые значения будем называть параметрами модели

Расчеты параметров линейной модели (x₁–x₅₁) проводили с помощью программы RIADA методом многомерной линейной регрессии на основе сингулярного разложения (SVD⁵) [45].

Описание нейросетевой модели

Построение нейросетевой модели проводили с использованием базы данных ИНХС "Газоразделительные параметры стеклообразных полимеров" [43]. Обучающий, валидационный и тестовый наборы данных получены разбиением полного набора данных (3052 системы полимер–газ) с использованием дивергенции Дженсена–Шенона (Jensen–Shannon divergence) [46] так, чтобы функции плотности распределения значений логарифмов коэффициентов проницаемости газов для трех полученных подмножеств были как можно более близкими, а информационная энтропия при этом была как можно выше. Обучающий набор содержал данные для 1534 систем полимер–газ, валидационный – для 901 системы, тестовый – для 617 систем. Нейросетевая модель структурно представляла собой остаточную трехслойную полносвязную нейронную сеть прямого распространения (52-5-4-1), содержащую 62 нейрона и 227 весов (рис. 3).

 

Рис. 3. Структура нейросети (а) и графическое представление весов межнейронных связей (б) после обучения нейросетевой модели, в которой K₂ = Σx_iA_i.

 

Параметры этой нейросетевой модели приведены на рис. 3б. Видно, что параметры находятся в приемлемом диапазоне (от –1.0 до +1.0), т.е. модель способна к активному дальнейшему обучению, и при необходимости может быть дообучена. С другой стороны рис. 3б показывает, что модель близка к сбалансированной, так как количество положительных параметров примерно равно количеству отрицательных, и отсутствуют отдельные резко выделяющиеся (на несколько единиц) параметры, т.е. не существует явно доминирующих нейронов. Следовательно, количество незначимых для модели параметров является малым.

Внутренние слои сети, помимо нейронов с функцией активации Leaky Rectified Linear Unit (LReLU) [40] (изображены на рис. 3а в виде эллипсов красного цвета), содержали по одному нейрону с рефлексивной функцией активации (эллипсы с бирюзовым фоном). Графики функции активации LReLu и рефлексивной функции активации приведены на рис. 4а и 4б соответственно. Входной и первый внутренний слои содержали по одному нейрону смещения (bias neurons), изображенных в виде эллипсов с белым фоном. Обучение (около 10000 эпох) проводили методом оптимизации Adam [47]. Для ускорения обучения начальное приближение весов, связанных с линейной частью модели (нейроны 65, 71), взяли из предварительно построенной линейной модели, остальные веса были получены по методу Xavier Glorot [48, 49] в диапазоне от –0.003 до +0.003 (т.е. взяты близкими к нулю, но получены современным способом инициализации начального приближения).

 

Рис. 4. Графическое представление использованных функций активации нейронов LReLU, при a = 3/8 (в русскоязычной литературе часто называется, как ReLU с утечкой) (а) и рефлексивная функция активации f (x) = x (б).

 

Подготовка данных

Входной вектор для каждой системы полимер–газ содержал в себе структурные дескрипторы полимера, полученные с использованием программы RIADA (49 шт.) для полимера и значения, (ε / k)_eff / d_eff² и 1 / d_eff² для газа (уравнения (2), (3)). Каждую координату входного вектора стандартизовали по всему набору данных. Экспериментальные свойства (lg P35, десятичный логарифм коэффициента проницаемости газа при 35 °С) центрировали и шкалировали на диапазон [–0.85…+0.85]. При прогнозировании свойства переход к исходному диапазону осуществляли обратным преобразованием.

Обучение нейронной сети

Обучение нейросетевой модели (около 10000 эпох) проводили методом оптимизации Adam [47]. В связи с небольшим объемом исходных данных по сравнению с современными нейросетями обучение осуществляли в многопоточном режиме (48 потоков) на процессоре AMD Ryzen 9 5900X без использования разбиения обучающей выборки на минибатчи⁶, и оно заняло не более 40 мин. Для предотвращения сильного отклонения отдельных значений от требуемых, на начальном этапе обучения использовали метод штрафных функций (penalty function method) в виде использующем функцию с ограниченным сверху ε-нечувствительным штрафом: penalty (x) = p₁exp (–4 (p₂ / x)²) , где x – отклонение от желаемого значения, p₁ = 100, p₂ = 0.2.

Исключение предпосылок к переобучению

Для предотвращения переобучения использовали метод ранней остановки [39]. Дроп-аут [50] не использовали в связи с малым количеством весов (227).

РЕЗУЛЬТАТЫ И ОБСУЖДЕНИЕ

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

По модели диффузия–растворение проницаемость пенетранта является результатом двух процессов – диффузии через сплошную среду и его растворения в ней [51]. Первый процесс характеризуется коэффициентом диффузии D, второй – коэффициентом растворимости S. Коэффициент растворимости соответствует угловому коэффициенту изотермы сорбции пенетранта в этой среде [52]. Для итогового процесса при невысоких давлениях (активностях), где оба параметра не зависят от давления, можно записать

$P=DS$, (4)

где P – коэффициент проницаемости пенетранта, причем каждый из сомножителей (D и S) для фиксированного полимерного образца зависит от свойств газа.

Температурная зависимость коэффициента диффузии пенетранта зачастую линеаризуется в аррениусовских координатах и описывается как

$\ln D=\ln D_0–E_D/RT$. (5)

Здесь E_D – энергия активации диффузии, D₀ – предэкспоненциальный множитель, характеризующий энтропию активации [53]. P. Meares предположил, что E_D определяется энергией, затраченной на образование цилиндрической полости (работа против сил когезии):

$E_D=1/4\times \pi d^2\lambda CED$, (6)

где d – диаметр цилиндра, соответствующий диаметру частицы пенетранта, λ – длина цилиндра, соответствующего длине диффузионного скачка пенетранта, CED – плотность энергии когезии [53]. Комбинация уравнений (5) и (6) в итоге приводит к общей зависимости, связывающей коэффициент диффузии и диаметр частицы пентранта:

$lgD=K_1–K_2{d}^2$, (7)

 

(K₁ и K₂ – константы для фиксированного образца полимера). Для диаметра молекул пенетранта мы использовали шкалу Теплякова–Мирса, построенную на основе тщательного анализа коэффициентов диффузии молекул газов и паров в полимерах, находящихся в различных релаксационных состояниях [44].

Температурная зависимость коэффициента растворимости также линеаризуется в аррениусовских координатах

$\ln S=\ln S_0+\Delta H_S/RT$, (8)

где угловой коэффициент ΔH_S соответствует энтальпии сорбции пенетранта. Энтальпию сорбции можно представить как тепловой эффект суммы двух процессов – конденсации (ΔH_c) и смешения (ΔH_m) [52]:

$\Delta H_S=\Delta H_c+\Delta H_m$. (9)

Как правило, |ΔH_c| >> |ΔH_m|, и поэтому энтальпия сорбции в основном определяется энтальпией конденсации. Известно, что энтальпия конденсации (теплота испарения) хорошо коррелирует с такими свойствами неполярных пенетрантов, как температура кипения [11] и поверхностное натяжение [54, 55]. В литературных источниках представлено более 20 различных формул, характеризующих связь энтальпии испарения с критическими параметрами и температурой кипения [56, 57]. Однако следует подчеркнуть, что использование подобного рода сложных зависимостей направлено на обеспечение высокой точности расчета для широкого круга различных по химической природе соединений. В случае ограниченного набора газов и легких паров – при небольших требованиях к точности – можно показать, что теплота испарения коррелирует также с критической температурой, параметром потенциала Леннард-Джонса (ε/k) и инвариантом потенциала Леннард-Джонса (ε/k)σ³/V_m (рис. 5). С учетом этого из уравнений (8) и (9) можно получить приближенные выражения для оценки десятичного логарифма коэффициента растворимости:

$lgS=K_3+K_4(\epsilon /k)$, (10)

$lgS=h_3+h_4\left(\Delta H_{vap}\right)$, (11)

$lgS=t_3+t_4\left(T_c\right)$, (12)

$lgS=g_3+g_4\left({\gamma }_{BP}\right)$, (13)

$lgS=K_3^{\#}+K_4^{\#}(\epsilon /k){\sigma }^3/V_m$. (14)

 

Рис. 5. Корреляция энтальпии испарения ΔH_vap c критической температурой T_cr (а), параметром потенциала Леннард-Джонса (ε / k), поверхностным натяжением и инвариантом потенциала Леннард-Джонса (ε / k) σ³ / V_m для константных газов и легких паров. Для построения рисунков использованы (ε / k) рассчитанные по методике, приведенной в работе [58].

 

Здесь K₃, h₃, t₃, g₃, K₃^# и K₄, h₄, t₄, g₄, K₄^# – константы для фиксированного полимера, ΔH_vap – энтальпия испарения пенетранта при температуре кипения, γ_BP – поверхностное натяжение жидкой фазы пенетранта при температуре кипения, [(ε / k) σ³ / V_m] – инвариант потенциала Леннард-Джонса.

Комбинирование уравнений (4), (7) и (10) приводит к обобщенной линейной зависимости, которая должна выполняться при постановке эксперимента в варианте один полимерный образец–разные газы:

$lgР=(K_1+K_3)–K_2{d}^2+K_4(\epsilon /k)$. (15)

Анализ значений K₁, K₂, K₃, K₄ для 30 полимеров трех различных типов, проведенный на основе данных из работы [44], позволяет заключить, что K₁ изменяется примерно на 23%, K₃ на 46%, K₄ на 90% (в 1.9 раза), K₂ на 850% (в 9.5 раз) [59]. Несмотря на то, что оценка диапазонов изменения данных параметров весьма приближенная, наилучшим кандидатом для моделирования методом атомных вкладов [25], по-видимому, является параметр K₂, так как он в относительном плане обладает наиболее широким диапазоном изменения значений. Математическое описание этой модели представлено уравнениями (2) и (3), а результаты отражены на рис. 6а.

 

Рис. 6. Сравнение экспериментальных и расчетных значений lgP35 / d_eff² для линейной (K₂ = Σx_iA_i) (а) и нейросетевой (б) моделей. Единицы измерения: P35 – Баррер, d_eff – ангстремы.

 

Помимо данной модели, аналогичным образом была построена модель, в которой параметр K₄ был представлен в виде Σ (x_iA_i), причем модель показала крайне низкую сходимость для Не, и все 89 точек пришлось отбросить. При этом коэффициент детерминации R² составил 0.7445. Для определения причин такого поведения были отдельно построены обобщенные модели для десятичных логарифмов коэффициентов диффузии (lg D35) и десятичных логарифмов коэффициентов растворимости (lg S35). Модели для логарифмов коэффициентов диффузии показали удовлетворительные результаты для всех газов:

$(8+lgD35)/\left(d_{eff}^2\right)=\Sigma (x_i{A}_i/d_{eff}^2)–K_2, R^2= 0.7856, получено K_2=0.501$; (16)

$(8+lgD35)/\left(d_{eff}^2\right)=K_1\left(d_{eff}^2\right)+\left(\Sigma x_i{A}_i\right), R^2= 0.9198, получено K_1=4.454$. (17)

Одна из моделей для десятичных логарифмов коэффициентов растворимости оказалась приемлемой для всех газов:

$(2+lgS35)/(\epsilon /k)=\Sigma (x_i{A}_i/(\epsilon /k\left)\right)+K_4, R^2=0.8260, K_4=0.0108$. (18)

Вторая модель, несмотря на хороший коэффициент детерминации, плохо предсказывала коэффициент растворимости гелия:

$(2+lgS35)/(\epsilon /k) =K_3/(\epsilon /k)+\left(\Sigma x_i{A}_i\right), R^2=0.8853, K_3=–1.492 (c He)$. (19)

Причины таких отклонений для гелия, с одной стороны, могут быть вызваны большими ошибками определения коэффициента диффузии и, следовательно, коэффициента растворимости гелия, особенно, для средне- и высокопроницаемых полимеров. С другой стороны, причина отклонений может заключаться в некорректном разделении параметров потенциала Леннард-Джонса, связанном, по мнению авторов книги [60], с самой математической формой потенциала Леннард-Джонса [61], не позволяющей отдельно получить параметры (ε / k) и σ:

$U\left(r\right)/k=4\left(\epsilon /k\right)\left[\frac{{\sigma }^{12}}{r^{12}}-\frac{{\sigma }^6}{r^6}\right]$.

Там же рекомендовано использовать оценки (ε / k) и σ из критических параметров газов по методике, предложенной в работе [58]:

$\sigma {\left(\frac{P_c}{T_c}\right)}^{1/3}=2.3551-0.0874\omega$, (21)

$\frac{\epsilon }{k{T}_c}=0.7915+0.1693\omega$, (22)

где ω – фактор ацентричности; σ – параметр потенциала Леннард-Джонса, Å; P_c – критическое давление, атм; T_c – критическая температура, K; ε – параметр потенциала Леннард-Джонса, эрг; k = 1.3805 × 10^–16 эрг/K – постоянная Больцмана.

Использование оценочных значения (ε / k) для гелия 3.775 и 30.62 K не привело к существенному улучшению ситуации. Первая оценка была выполнена на основе критической температуры и фактора ацентричности по методике [58], вторая оценка сделана на основе парахора и ψ-фактора в соответствии с формулами приведенными в работе [62]. Таким образом, предположение о некорректном разделении параметров потенциала не нашло своего подтверждения.

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

 

Рис. 7. Сравнение экспериментальных и расчетных значений lgP35 / (ε / k) для линейной (K₄ = Σx_iA_i) (а) и нейросетевой (б) моделей без He. Единицы измерения: P35 – Баррер, (ε / k) – градусы Кельвина.

 

Далее была построена линейная модель, в которой параметр (K₁ + K₃) был представлен в виде Σ (x_iA_i). Эта линейная модель (рис. 8а) оказалась хуже по сравнению с линейной моделью, представленной на рис. 6а, по коэффициенту детерминации (R² = 0.7627), поэтому она тоже была исключена из дальнейшего рассмотрения.

 

Рис. 8. Сравнение экспериментальных и расчетных значений lgP35 [Баррер] для линейной (а) и нейросетевой (б) моделей.

 

Таким образом, из всех моделей, построенных на основе выражения (13), наиболее удачной оказалась линейная модель, в которой K₂ для каждого полимера был представлен суммой атомных вкладов в виде Σ (x_iA_i). Построив модель, на выходе были получены инкременты (x_i) для расчета K₂ и фиксированные значения (K₁ + K₃) = 3.293 и K₄ = 0.01166, которые, как и инкременты, являются параметрами модели.

Далее была предпринята попытка улучшить линейную модель, путем введения слабой нелинейной составляющей при помощи остаточной нейронной сети. На вход этой нейронной сети подавали те же параметры, что и на вход линейной модели (т.е. структурные дескрипторы A₁–A₄₉ для полимера и два дескриптора A₅₀ = (ε / k)_eff / (d_eff²) и A₅₁ = 1 / (d_eff²) для газа). По сравнению с линейной моделью (рис. 6а), нейросетевая модель (рис. 6б), полученная в результате обучения, имеет меньшую погрешность прогноза логарифма коэффициента проницаемости газов и паров при 35 °С на тестовом наборе (в обучении не используется) и более близкое к нормальному распределение остатков.

Аналогичным образом были построены и обучены нейросетевые модели для линейной модели, в которой параметры K₂ (рис. 7б) и (K₁ + K₃) (рис. 8б) были представлены в виде Σ (x_iA_i). Обе модели также позволили улучшить результаты линейной модели (рис. 7а и 8а), однако следует подчеркнуть, что улучшение оказалось не таким значительным, как на рис. 6б. Обучение потребовало более 500000 эпох, несмотря на это, и линейная, и нейросетевая модели оказались хуже модели, представленной на рис. 6б.

Для линейной модели, в которой параметры K₂ были представлены в виде Σ (x_iA_i) и усовершенствованы нейросетевой моделью, было показано, что коэффициенты детерминации в случае конкретных газов (не менее 50 точек для каждой выборки (исходная, обучающая, тестовая, валидационная) имеют удовлетворительные значения (R² ≥ 0.80).

ЗАКЛЮЧЕНИЕ

Важным результатом данной работы является создание обобщенной линейной модели для прогнозирования десятичного логарифма коэффициента проницаемости газов и паров. Если ранее для расчета коэффициента проницаемости каждого газа использовали свои параметры линейной модели, то по итогам текущей работы в обобщенной линейной модели параметры x₁–x₄₉ могут применяться к любой системе полимер–газ, в том числе для прогнозирования коэффициента проницаемости неизученного газа, для которого известны параметры потенциала Леннард-Джонса и кинетический диаметр газа (или их оценки). Таким образом, обобщенная линейная модель на вход получает простейшие структурные дескрипторы для полимера и легкодоступные свойства газа, а на выходе выдает десятичный логарифм коэффициента проницаемости газа при 35 °С. Также показана возможность улучшения обобщенных линейных моделей путем введения слабой нелинейной составляющей, с помощью малых остаточных нейронных сетей. При использовании нейронных сетей для улучшения наиболее точной обобщенной модели (с K₂ = Σ x_iA_i) коэффициент детерминации возрастает с R² = 0.84 до R² = 0.91.

Работа выполнена в рамках Госзадания ИНХС РАН.

¹ В оригинальной работе понятие “уникальный набор структурных дескрипторов” описывается одним словом: fingerprints.

² KRR – Kernel Ridge Regression (в русскоязычной литературе: метод гребневой регрессии на ядрах Мерсера).

³ GPR – Gauss Process Regression (в русскоязычной литературе: метод гауссовой регрессии).

⁴ RBF – Radial Basis Funсtion (в русскоязычной литературе: радиальная базисная функция).

⁵ SVD – Singular Value Decomposition (в русскоязычной литературе: разложение по сингулярным числам).

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

About the authors

Д. А. Царев

Инстинут нефтехимического синтеза им. А.В. Топчиева Российской академии наук

Email: alentiev@ips.ac.ru
Russian Federation, 119991, Москва, Ленинский пр., 29

В. Е. Рыжих

Инстинут нефтехимического синтеза им. А.В. Топчиева Российской академии наук

Email: alentiev@ips.ac.ru
Russian Federation, 119991, Москва, Ленинский пр., 29

Н. А. Белов

Инстинут нефтехимического синтеза им. А.В. Топчиева Российской академии наук

Email: alentiev@ips.ac.ru
Russian Federation, 119991, Москва, Ленинский пр., 29

А. Ю. Алентьев

Инстинут нефтехимического синтеза им. А.В. Топчиева Российской академии наук

Author for correspondence.
Email: alentiev@ips.ac.ru
Russian Federation, 119991, Москва, Ленинский пр., 29

References

Supplementary files