Лекция 8. Нелинейная регрессия

Краткое содержание

1. Нелинейная регрессия: постановка задачи 2. Методы Ньютона, Гаусса-Ньютона, Левенберга-Марквардта 3. Функции в языке GNU Octave 4. Практическая реализация нелинейной регрессии в GNU Octave

Виды регрессии

1. Линейная регрессия

𝑦 = 𝑏 0 + 𝑏 1 𝑥 1 𝑋𝛽 = 𝑦; + ⋯ + 𝑏 𝑛 𝑥 𝑛 ⊤ 𝑋 −1 𝑋 ⊤ 𝑦

2. Нелинейная регрессия, поддающаяся линеаризации

Примеры: 𝑦 = 𝑏 𝑘 = 𝑘 0 0 + 𝑏 1 𝑓 1 (𝑥 1 , … , 𝑥 𝑛 ) + ⋯ + 𝑏 𝑛 𝑓 𝑛 (𝑥 1 , … , 𝑥 𝑛 ) exp − 𝐸 𝑎 𝑅𝑇 ⇒ ln 𝑘 = ln 𝑘 0 − 𝐸 𝑎 𝑅𝑇 ; 𝑐 𝑝 𝑇 = 𝑎 + 𝑏𝑇 + 𝑐𝑇 2 ;

3. Нелинейная регрессия, не поддающаяся линеаризации

𝑦 = 𝑓 𝛽 1 , … , 𝛽 𝑚 , 𝑥 1 , … , 𝑥 𝑛 1 Пример: 𝑦 = 1+exp −𝛽 1 −𝛽 2 𝑥 𝛽, ) Для решения таких задач используют численные методы: Ньютона, Гаусса-Ньютона, градиентного спуска, Левенберга-Марквардта

Численные методы и МНК

Нелинейная система уравнений

𝜑 𝛽 1 , … , 𝛽 𝑘 , 𝑥 11 , 𝑥 21 , … 𝑥 𝑚1 𝜑 𝛽 1 , … , 𝛽 𝑘 , 𝑥 12 , 𝑥 22 , … , 𝑥 𝑚2 … 𝜑 𝛽 1 , … , 𝛽 𝑘 , 𝑥 1𝑛 , 𝑥 2𝑛 , … , 𝑥 𝑚𝑛 = 𝑦 1 = 𝑦 2 = 𝑦 𝑛

Как правило, полученную систему уравнений нельзя решить аналитическими методами

Нужно использовать численные методы

Минимизируемая сумма квадратов

𝑖 − 𝑦 𝑖 2 = 𝛽, 𝑖 𝑖 𝑖

Система уравнений для поиска минимума (ср. с линейным МНК)

− 𝑦 𝑖 2 = 𝑓 𝑖 2 𝑖 𝜕𝛽 1 ⋯ = 0 = 0

Рассматриваемые численные методы

1. Метод Ньютона 2. Метод Гаусса-Ньютона 3. Метод Левенберга-Марквардта 𝜕𝛽 𝑚

Метод Ньютона

Одномерный случай (нелинейное уравнения)

𝑓 𝑥 ≈ 𝑓 𝑥 0 + 𝑓 ′ 𝑥 = 𝑥 0 − 𝑥 0 𝑓(𝑥 0 ) 𝑓 ′ 𝑥 𝑥 − 𝑥 0 0 = 0

Многомерный случай (система нелинейных уравнений)

𝛽 ≈ 𝐹 0 + 𝜕𝐹 0 𝜕𝛽 𝑖 𝑖 𝛽 𝑖 − 𝛽 𝑖 0 = 0 ⇒ 𝐹 0 + 𝛻𝐹 0 1.

2.

3.

Особенности метода Ньютона

Быстро сходится (квадратичная сходимость) В случае нескольких корней может быть найден любой из них (зависит от начального приближения) Скорость сходимости (и сама сходимость) может зависеть от начального приближения

Метод Ньютона: практическая реализация

𝑓 𝑓 𝑥 = cos 𝑥 − 𝑥 ′ 𝑥

Уравнение

3 = − sin 𝑥 − 3𝑥 2 = 0 = 0 Касательная Прибл.

решение

Программа format long; delta = 1e-15; x0 = 0.5; x = x0 + 2*delta; while abs(x-x0) >= delta x0 = x; f = cos(x0)-x0^3; df = -sin(x0)-3*x0^2; x = x0 – f/df; disp(x); end

Функция Начальное прибл.

2.000000000000000

1.348055393079852

1.001262924310922

0.880591138294078

0.865691456412747

0.865474078978736

0.865474033101617

0.865474033101614

МНК и метод Гаусса-Ньютона

Система уравнений и метод Ньютона

𝜕𝛽 1 ⋯ = 0 ⇒ 𝜕𝛽 𝑖 𝜕𝛽 𝑚 = 0 ≈ 𝜕𝛽 𝑖 0 0 + 𝑗 𝜕𝛽 𝑖 𝜕𝛽 𝑗 𝛽 𝑗 − 𝛽 𝑗 0 = 0 0

Матричная запись и метод Гаусса-Ньютона

0 ⊤ 𝛽 0 0 + 𝐽 ⊤ 𝛽 0 0 𝛻𝐹 𝑖 𝜕𝐹 = 𝜕𝛽 𝑖 𝜕 = 𝜕𝛽 𝑖

Градиент, якобиан и гессиан

𝑓 𝑗 2 𝑗 𝜕𝑓 𝑗 = 2 𝑓 𝑗 𝜕𝛽 𝑖 𝑗 = 2𝐽 ⊤ 𝑓; 𝐽 𝑖𝑗 = 𝜕𝑓 𝑖 𝜕𝛽 𝑗 = 𝜕𝜑( 𝑖 ) 𝜕𝛽 𝑗 𝐻 𝑖𝑗 𝜕𝐹 𝛽 = 𝜕𝛽 𝑖 𝜕𝛽 𝑗 𝜕 = 𝜕𝛽 𝑗 2 𝑓 𝑘 𝜕𝑓 𝑘 𝜕𝛽 𝑖 𝑘 = 2 𝑘 𝜕𝑓 𝑘 𝜕𝑓 𝑘 𝜕𝛽 𝑖 𝜕𝛽 𝑗 + 𝑓 𝑘 𝜕𝑓 𝑘 𝜕𝛽 𝑖 𝜕𝛽 𝑗 ≈ 2 𝜕𝑓 𝑘 𝜕𝛽 𝑖 𝑘 𝜕𝑓 𝑘 𝜕𝛽 𝑗 = 2𝐽 ⊤ 𝐽

Метод Левенберга-Марквардта

1.

2.

3.

4.

𝐽 ⊤ 𝐽 + 𝜆 ⊤ 𝑓

Каждая итерация включает в себя подбор шага

𝝀

:

Взять начальное (очень малое) значение 𝜆 Выполнить итерацию, найдя 𝑝 Если значение 𝐹 возросло, то увеличить шаг 𝜆 и вернуться к п.1

Если значение 𝐹 уменьшилось, то принять полученное значение перейти к следующей итерации и

Сходится медленнее метода Гаусса-Ньютона, но менее требователен к начальному приближению. Широко применяется на практике.

Метод Л.-М. – комбинация методов Гаусса-Ньютона и градиентного спуска

Обозначения: а) метод градиентного спуска: 𝛽 б) метод Гаусса-Ньютона: 𝐽 ⊤ 𝐽 𝛽 − (0) , 𝛽 (0) 0 – начальное приближение; итерации, 𝜆 – шаг, 𝐹 – минимизируемая функция 𝛽 – результат

О функциях в GNU Octave

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

Именованные функции Анонимные функции

Хранятся в файлах. Синтаксис:

% Справочная информация function [o1,…,om]=funcname(i1,…,in) % Тело функции end

Хранятся в переменных. Синтаксис:

func = @(arg1,…,argn) expression

Примеры:

sqr = @(x) x.^2 len = @(x,y) sqrt(x.^2 + y.^2)

Пример:

% SQREQ_ROOTS Finds the roots % of square equation function [x1,x2] = sqreq_roots(a,b,c) D = b.^2 - 4*a.*c; x1 = (-b + sqrt(D)) ./ (2*a); x2 = (-b - sqrt(D)) ./ (2*a); end

См. также: вложенные функции (nested function), MEX-файлы, объектно ориентированное программирование

Нелинейная регрессия: практическая реализация

Шаг 1. Задание аппроксимирующей функции

𝑦 = 𝛽 1 + 𝛽 2 exp(−𝛽 3 𝑥)

Шаг 3. Запись на языке Octave b0 = [10 10 10]; [b, res, J] = ...

lsqfit_lm(X,Y,@func, b0); [db, b_lb, b_ub, sb] = ...

lsqfit_ci(b, res, J); function [F, J] = func(b, x) F = b(1) + b(2)*exp(-b(3)*x); if nargout == 2 df_db1 = ones(size(x)); df_db2 = exp(-b(3)*x); df_db3 = -b(2)*exp( b(3)*x).*x; J = [df_db1 df_db2 df_db3]; end Шаг 2. Задание якобиана

𝜕𝑦 𝑖 = 1 𝜕𝛽 1 𝜕𝑦 𝑖 𝜕𝛽 2 𝜕𝛽 3 𝜕𝑦 𝑖 = exp −𝛽 = −𝑥𝛽 2 3 𝑥 exp(−𝛽 𝑖 3 𝑥 𝑖 )

Шаг 4. Подбор начального приближение и визуализация результатов

Нелинейная регрессия: доверительные интервалы

Линейная регрессия

𝑠 2 = 2 𝑋 ⊤ 𝑋 −1

Нелинейная регрессия

𝑠 2 = 2 𝐽 ⊤ 𝐽 −1 При этом якобиан и отклонения рассчитываются в точке

Исходная система уравнений

𝜑 𝛽, 𝑥 1 ⋯ 𝜑 𝛽, 𝑥 2 = 𝑦 1 = 𝑦 𝑛

Результат линеаризации в векторной форме Разложение в ряд Тейлора

𝛽 𝑖 − 𝑖 𝑖 𝜕 𝑖