Russian Federation
Angarskiy gosudarstvennyy tehnicheskiy universitet (Vychislitel'nye mashiny i kompleksy, Professor)
Russian Federation
The implementation of a computer program in C for calculating the parameters of a multiple linear regression equation using the ridge regression method is considered. The calculations were performed in matrix form using L2 regularization to minimize the risk of overfitting the model and reduce the im-pact of multicollinearity on the stability of the estimates obtained. The integration of the developed program with the technological monitoring system and the creation of a virtual heavy reformat octane analyzer based on it is demonstrated
multiple linear regression, ridge regression method, program development, virtual analyzer
Нахождение параметров в уравнении множественной линейной регрессии (МЛР) может быть затруднено при аналитическом решении, особенно в условиях большого объема данных. Применение численных методов позволяет решить эту задачу. Однако в ряде случаев требуется периодический пересчёт параметров для обеспечения корректной работы модели в условиях изменчивости входных данных. Для автоматизации данного процесса предлагается разработка компьютерной программы, которая будет периодически находить искомые параметры и подставлять их в уравнение МЛР.
Для решения поставленной задачи, классический метод наименьших квадратов для поиска параметров уравнения МЛР не подходит. Поскольку этот метод основан на минимизации суммы квадратов отклонений между наблюдаемыми данными и значениями, предсказанными моделью, он позволяет найти наилучшие аппроксимации для набора экспериментальных данных, но по причине того, что метод наименьших квадратов (МНК) является несмещённой оценкой с минимальной дисперсией, это имеет место быть только при отсутствии мультиколлинеарности - сильной линейной зависимости между объясняющими переменными. В случае, при наличии мультиколлинеарности, МНК оказывает деструктивное воздействие на свойства оценок, полученных с его помощью, поскольку при наличии сильно коррелирующих признаков, матрица Грама становится близкой к вырожденной, что приводит к неустойчивости оценок коэффициентов, огромной дисперсии параметров и плохой обобщающей способности математической модели [3].
|
|
(1) |
где 
– зависимая переменная, 
– свободный член, 
– параметры регрессии для каждой независимой переменной, 
– независимые переменные, которые влияют на , ε – случайная ошибка модели, случайные вариации или нелинейные эффекты, которые не объясняются независимыми переменными.
|
|
(2) |
где Q – сумма квадратов остатков, 
– реальное значение зависимой переменной.
Для исключения этих недостатков МНК в решении поставленной задачи планируется ввести штраф за сложность модели (L2-регуляризация). Введение ограничения на норму вектора весов позволит стабилизировать математическую модель и сделать её более устойчивой к шуму в данных и мультиколлинеарности [3].
Данное ограничение реализуется через модификацию целевой функции МНК, путём добавления к исходному функционалу штрафного слагаемого, пропорционального квадрату евклидовой нормы вектора весов. Этот переход эквивалентен решению задачи условной минимизации (нахождение экстремума при наличии ограничений) [4].
|
|
(3) |
где n – количество строк наблюдений (строк в данных), m – количество признаков, yi – реальное значение целевой переменной для i-го объекта, xij – значение j-го признака для i-го объекта, βj – искомые параметры модели, λ – коэффициент регуляризации (λ>0).
Для нахождения оптимального вектора параметров 
в матричном виде, необходимо минимизировать уравнение 3. Дифференцируя выражение по β и приравнивая результат к нулю, мы получим:
|
|
(4) |
где T - знак транспонирования; 
– единичная матрица размерности (m+1)x(m+1).
После приведения задачи к виду безусловной оптимизации и используя метод гребневой регрессии мы сознательно отказываемся от несмещенности, из-за чего модель показывает небольшую систематическую ошибку и «недоучивается» на обучающей выборке, но зато она становится более устойчивой. То есть при изменении данных, параметры модели почти не изменяются, так как регуляризация удерживает их в узких рамках. Благодаря этому, модель также игнорирует случайные выбросы (шум) в данных и мультиколлинеарность, концентрируясь на главных признаках [5].
Ключевым аспектом реализации данного метода является выбор оптимального значения параметра регуляризации λ. Поскольку теоретически предсказать его величину невозможно, на практике применяется метод перекрёстной проверки (кросс-валидации). В ходе этого метода данные необходимо разбить на обучающую и валидационную выборки, после чего выбирается такое значение λ, которое минимизирует ошибку на данных, не участвующих в обучении. Такой подход позволяет гарантировать высокую обобщающую способность итоговой модели. Графически это представлено на рисунке 1.
Рисунок 1 – График анализа зависимости, составляющих среднеквадратичной ошибки от параметра регуляризации λ
При λ=0 (что соответствует МНК), модель обладает минимальным смещением, но максимальной дисперсией, что характерно для переобучения в условиях мультиколлинеарности. По мере увеличения λ, на рисунке 1 наблюдается монотонное снижение дисперсии оценок, что сопровождается умеренным ростом смещения. Оптимальное значение λ находится в точке минимума функции общей ошибки по формуле 5. В этой точке достигается наилучшая обобщающая способность модели: она становится достаточно устойчивой к шуму в данных, сохраняя при этом высокую точность аппроксимации.
|
|
(5) |
где 
– это среднеквадратичная ошибка, рассчитанная для модели, к которой применена регуляризация с определённым коэффициентом регуляризации, n- количество строк данных.
В отличие от МНК, который инвариантен к масштабу переменных, результат гребневой регрессии напрямую зависит от единиц измерения признаков. Поэтому важно перед применением метода гребневой регрессии центрировать данные. Это обусловлено тем, что штрафное слагаемое 
применяется ко всем параметрам модели одинаково. Если один признак измерен в миллионах, а второй в единицах, то параметр регуляризации λ будет «штрафовать» признаки с разными весами, что приведёт к искажению модели. Для недопущения этой проблемы, перед обучением модели выполняется преобразование всех признаков по формуле:
|
|
(6) |
где 
– среднее значение j-го признака, 
– среднее квадратичное отклонение этого признака.
После этого все признаки центрируются относительно нуля и приводятся к единичной дисперсии. Это гарантирует, что штраф за сложность модели распределится равномерно между всеми предикторами, позволяя параметру регуляризации эффективно выявлять реальную значимость каждого признака независимо от его исходного масштаба и гарантирует достижение наилучшего компромисса между смещением и дисперсией модели.
После нахождения оптимальных параметров модели, необходимо вычислить свободный член β0 по формуле:
|
|
(7) |
где 
– среднее значение y до центрирования данных
Нахождение свободного члена после всех параметров обусловлено тем, что если результат равен , то часть этого результата объясняется средними значениями факторов x с параметрами β [4]. Всё, что осталось необъяснённым будет описано свободным членом β0. Если же попытаться найти β0 сразу и добавить его в качестве вектора единиц в выборку исходных данных, это приведёт к тому, что штраф применится и к свободному члену. То есть алгоритм искусственно уменьшит β0 пытаясь минимизировать общую сумму квадратов, что влечёт за собой систематическую ошибку, которая будет занижать или завышать конечный результат даже если все зависимости найдены верно.
На основе описанного метода, создана программа, разработка которой велась на языке программирования C(Си), стандарта C17 (ISO/IEC 9899:2018).
Выбор данного языка обусловлен широкой поддержкой компилятора GCC, а также необходимостью прямого управления вычислительными ресурсами при работе с большими матрицами. Кроме этого, важным фактором является потребность в независимости от сторонних библиотек, что обеспечивает кроссплатформенность программы и стабильную работу на различных конфигурации компьютеров обеспечивая максимальное быстродействие.
Программная реализация последовательно воплощает этапы: от центрирования данных и управления динамической памятью до подпрограмм инверсии и регуляризации матрицы Грама. Для хранения данных использовались двумерные массивы типа double, а операции (умножение, инверсия) реализованы с помощью стандартной библиотеки math.h и метода Гаусса. Транспонирование выполнено в виде вложенного цикла с обменом индексов. Сложение параметра регуляризации λ проводилось только по диагональным элементам. Алгоритмическая структура этой последовательности действий выглядит следующим образом:
НАЧАЛО ПРОГРАММЫ
1. ЧТЕНИЕ ИСХОДНЫХ ДАННЫХ
| Открыть файлы "X.txt" и "Y.txt"
| Определить размерность: N (столбцы), Rows (строки)
| Считать данные в матрицы matX и matY
| Закрыть файлы
2. ПОДГОТОВКА ИСХОДНЫХ ДАННЫХ (Центрирование)
| Для каждого признака j:
| Вычислить среднее meanX[j]
| Вычислить среднее ответов meanY
| Вычесть средние значения из каждого элемента matX и matY
3. ПОИСК ОПТИМАЛЬНОЙ РЕГУЛЯРИЗАЦИИ (Grid Search)
| Установить min_mse = +∞
|
| Для каждого λ из набора {0.0001, 0.001, ..., 100.0}:
| |
| | А. Расчет коэффициентов:
| | 1. Вычислить матрицу Грама: G = (Xᵀ * X)
| | 2. Регуляризация: прибавить λ к диагонали G
| | 3. Решить уравнение: beta = inv(G) * (Xᵀ * Y)
| |
| | Б. Валидация (MSE):
| | 1. Вычислить предсказания: Y_pred = X * beta
| | 2. Ошибка: MSE = средний квадрат разности (matY - Y_pred)
| |
| | В. Сохранение результата:
| | Если MSE < min_mse:
| | min_mse = MSE, best_λ = λ, best_beta = beta
| |
| Конец цикла по λ
4. ФИНАЛЬНЫЙ ЭТАП
| Рассчитать свободный член: b0 = meanY - Σ(best_beta[j] * meanX[j])
| Сформировать итоговое уравнение: y = b0 + b1*x1 + ... + bn*xn
| Освободить память
КОНЕЦ ПРОГРАММЫ
Далее рассматриваются результаты экспериментальной апробации разработанной программы и интеграции её в систему технологического мониторинга.
Система технологического мониторинга (СТМ) представляет собой программно-аппаратный комплекс, предназначенный для непрерывного сбора, обработки и анализа данных о состоянии промышленного оборудования, инженерных систем и параметров технологических процессов в реальном времени. Для проверки работоспособности разработанной программы, на базе СТМ реализован виртуальный анализатор для определения октанового числа тяжёлого риформата по двум методам: моторному и исследовательскому [2]. В качестве входных данных используются доступные поточные параметры. Алгоритм выполняет их центрирование, определяет параметр регуляризации и вычисляет коэффициенты уравнения МЛР (множественной линейной регрессии). Полученные значения записываются в соответствующие теги, что позволяет использовать их в расчетах виртуального анализатора или располагать на мнемосхемах и визуализировать работу в СТМ.
На рисунке 2 представлена схема архитектуры интеграции разработанной программы в СТМ. Необходимые данные с полевых датчиков, установленных на линии тяжелого риформата, передаются в блок виртуального анализатора. Математическое ядро системы, реализованное на языке C, в режиме реального времени обрабатывает входные векторы признаков с использованием алгоритма гребневой регрессии. Результат вычислений (прогнозное октановое число) транслируется на мнемосхему, представленную на рисунке 3.
Рисунок 2 - Схема интеграции виртуального анализатора на базе алгоритма гребневой регрессии в СТМ
Условия проведения эксперимента следующие: тестирование работы программы и проведение эксперимента проводилось на аппаратном обеспечении компьютера, представленном в таблице 1 и программном обеспечении, представленном в таблице 2.
Таблица 1 – Аппаратное обеспечение компьютера для проведения эксперимента с тестированием разработанной программы
|
Компонент |
Характеристики |
|
Процессор (CPU) |
Intel(R) Core (TM) i5-4570 3.2ГГц |
|
Архитектура |
X86-64 |
|
Оперативная память (RAM) |
8 Гб |
|
Видеоадаптер (GPU) |
Intel (R) HD Graphics 4600 (113 Мб) |
|
Накопитель |
SSD 240 Гб |
Таблица 2 - Программное обеспечение компьютера для проведения эксперимента с тестированием разработанной программы
|
Параметр |
Значение |
|
Операционная система |
Ubuntu 22.04.3 LTS (Kernel 5.15) |
|
Стандарт языка |
C17 (ISO/IEC 9899:2018) |
|
Компилятор |
GCC version 11.4.0 |
|
Оптимизация |
-02 (стандартная) |
|
Библиотеки |
glibc 2.35, OpenSSL 3.0.2 |
Для замера времени выполнения программы использовалась встроенная библиотека time.h. Тест проводился при условии, что данные необходимо считать с твердотельного накопителя и загрузить их в кэш оперативной памяти компьютера. Процессор свободен. Результаты замеров представлены в таблице 3.
Таблица 3 - Результаты измерения времени выполнения алгоритма (в мс)
|
Объём входных данных (N) |
Время выполнения, мс |
|
10000 |
2 |
По результатам проведённых замеров, можно сделать вывод, что реализованную программу на языке C можно использовать в системах реального времени для задач с умеренным количеством признаков. Однако для работы с очень большим количеством входных данных (больше 100000) и высокой размерностью целесообразно использовать итерационные численные методы (например, градиентный спуск), которые не требуют ресурсо-затратной операции инверсии матрицы.
Для наглядной демонстрации работоспособности разработанной программы, была создана мнемосхема, представленная на рисунке 3. На ней отражены расчётные виртуальным анализатором показания, которые сравниваются с лабораторными показаниями.
Рисунок 3 – Экспериментальная мнемосхема, отражающая сравнение вычислений, выполненных виртуальным анализатором, с полученными в лаборатории
Сравнительный анализ работы виртуального анализатора на основе разработанной программы для различных методов определения октанового числа выявил характерные особенности. Для моторного метода зафиксирована меньшая среднеквадратичная ошибка чем у исследовательского, что объясняется более высокой стабильностью параметров детонационной стойкости в жестких условиях испытаний. В то же время, для исследовательского метода наблюдается более высокий коэффициент детерминации при большем значении абсолютной ошибки. Данный факт указывает на то, что разработанная математическая модель более чувствительна к вариациям состава компонентов, определяющих исследовательское октановое число, и успешно адаптируется к более широкому диапазону входных данных.
Кроме этого, анализ графиков с мнемосхемы на рисунке 2 указывает на соответствие вычислений, выполненных моделью реальным. Для моторного метода достигнуто значение среднеквадратичной ошибки (MSE) равное 0,294, при коэффициенте детерминации (R2) 0,61, а для исследовательского метода MSE = 0,404, при R2 = 0,63. Полученные значения среднеквадратичной ошибки (MSE = 0,294...0,404) соответствуют отклонению в пределах 0,5-0,6 единиц октанового числа, что сопоставимо с допустимой прецизионностью стандартных лабораторных методов определения детонационной стойкости (ГОСТ 8226 – 2022 (для моторного метода) и ГОСТ 511 – 2022 (для исследовательского метода)). Этого достаточно, чтобы технолог мог оперативно подстроить режим установки, не дожидаясь анализов из лаборатории (которые делаются 2 - 4 часа и являются дорогостоящими). Коэффициент детерминации свыше 0,6 также подтверждает высокую прогностическую способность модели в условиях промышленной неопределенности состава сырья.
Благодаря тому, что разработка программы велась на языке программирования C(Си), она обладает большим потенциалом для масштабируемости. Программа может быть адаптирована для использования во встраиваемых системах, распределённых системах управления (РСУ), для проектирования аналогичных виртуальных анализаторов, и использована в других задачах, требующих автоматического расчета параметров моделей МЛР.
1. Demidchenko E.A., Istomin A.L. Programma dlya rascheta parametrov v modelyah mnozhestvennoy lineynoy regressii. // Sbornik nauchnyh trudov mo-lodyh uchenyh AnGTU. Angarsk, 2026.
2. Demidchenko E.A., Istomin A.L. Virtual'nyy analizator oktanovogo chisla tyazhelogo riformata // Vestnik AnGTU. Angarsk, 2025. – s. 131 – 133.
3. Hoerl, A. E., & Kennard, R. W. Ridge Regression: Biased Estimation for Nonorthogonal Problems. //Technometrics, 1970.
4. Hastie, T., Tibshirani, R., & Friedman, J. The Elements of Statistical Learning: Data Mining, Inference, and Prediction. //Springer, 2009.
5. James, G., Witten, D., Hastie, T., & Tibshirani, R. An Introduction to Sta-tistical Learning (ISL). //Springer, 2013.
