Выбрав вид функции регрессии, т.е. вид рассматриваемой модели зависимости Y от Х (или Х от У), например, линейную модель y x =a+bx, необходимо определить конкретные значения коэффициентов модели.

При различных значениях а и b можно построить бесконечное число зависимостей вида y x =a+bx т.е на координатной плоскости имеется бесконечное количество прямых, нам же необходима такая зависимость, которая соответствует наблюдаемым значениям наилучшим образом. Таким образом, задача сводится к подбору наилучших коэффициентов.

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

Обозначим: Y i - значение, вычисленное по уравнению Y i =a+bx i . y i - измеренное значение, ε i =y i -Y i - разность между измеренными и вычисленными по уравнению значениям, ε i =y i -a-bx i .

В методе наименьших квадратов требуется, чтобы ε i , разность между измеренными y i и вычисленными по уравнению значениям Y i , была минимальной. Следовательно, находим коэффициенты а и b так, чтобы сумма квадратов отклонений наблюдаемых значений от значений на прямой линии регрессии оказалась наименьшей:

Исследуя на экстремум эту функцию аргументов а и с помощью производных, можно доказать, что функция принимает минимальное значение, если коэффициенты а и b являются решениями системы:

(2)

Если разделить обе части нормальных уравнений на n, то получим:

Учитывая, что (3)

Получим , отсюда , подставляя значение a в первое уравнение, получим:

При этом b называют коэффициентом регрессии; a называют свободным членом уравнения регрессии и вычисляют по формуле:

Полученная прямая является оценкой для теоретической линии регрессии. Имеем:

Итак, является уравнением линейной регрессии.

Регрессия может быть прямой (b>0) и обратной (b Пример 1. Результаты измерения величин X и Y даны в таблице:

x i	-2	0	1	2	4
y i	0.5	1	1.5	2	3

Предполагая, что между X и Y существует линейная зависимость y=a+bx, способом наименьших квадратов определить коэффициенты a и b.

Решение. Здесь n=5
x i =-2+0+1+2+4=5;
x i 2 =4+0+1+4+16=25
x i y i =-2 0.5+0 1+1 1.5+2 2+4 3=16.5
y i =0.5+1+1.5+2+3=8

и нормальная система (2) имеет вид

Решая эту систему, получим: b=0.425, a=1.175. Поэтому y=1.175+0.425x.

Пример 2. Имеется выборка из 10 наблюдений экономических показателей (X) и (Y).

x i	180	172	173	169	175	170	179	170	167	174
y i	186	180	176	171	182	166	182	172	169	177

Требуется найти выборочное уравнение регрессии Y на X. Построить выборочную линию регрессии Y на X.

Решение. 1. Проведем упорядочивание данных по значениям x i и y i . Получаем новую таблицу:

x i	167	169	170	170	172	173	174	175	179	180
y i	169	171	166	172	180	176	177	182	182	186

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

x i	y i	x i 2	x i y i
167	169	27889	28223
169	171	28561	28899
170	166	28900	28220
170	172	28900	29240
172	180	29584	30960
173	176	29929	30448
174	177	30276	30798
175	182	30625	31850
179	182	32041	32578
180	186	32400	33480
∑x i =1729	∑y i =1761	∑x i 2 299105	∑x i y i =304696
x=172.9	y=176.1	x i 2 =29910.5	xy=30469.6

Согласно формуле (4), вычисляем коэффициента регрессии

а по формуле (5)

Таким образом, выборочное уравнение регрессии имеет вид y=-59.34+1.3804x.
Нанесем на координатной плоскости точки (x i ; y i) и отметим прямую регрессии.

Рис 4

На рис.4 видно, как располагаются наблюдаемые значения относительно линии регрессии. Для численной оценки отклонений y i от Y i , где y i наблюдаемые, а Y i определяемые регрессией значения, составим таблицу:

x i	y i	Y i	Y i -y i
167	169	168.055	-0.945
169	171	170.778	-0.222
170	166	172.140	6.140
170	172	172.140	0.140
172	180	174.863	-5.137
173	176	176.225	0.225
174	177	177.587	0.587
175	182	178.949	-3.051
179	182	184.395	2.395
180	186	185.757	-0.243

Значения Y i вычислены согласно уравнению регрессии.

Заметное отклонение некоторых наблюдаемых значений от линии регрессии объясняется малым числом наблюдений. При исследовании степени линейной зависимости Y от X число наблюдений учитывается. Сила зависимости определяется величиной коэффициента корреляции.

После выравнивания получим функцию следующего вида: g (x) = x + 1 3 + 1 .

Мы можем аппроксимировать эти данные с помощью линейной зависимости y = a x + b , вычислив соответствующие параметры. Для этого нам нужно будет применить так называемый метод наименьших квадратов. Также потребуется сделать чертеж, чтобы проверить, какая линия будет лучше выравнивать экспериментальные данные.

Yandex.RTB R-A-339285-1

В чем именно заключается МНК (метод наименьших квадратов)

Главное, что нам нужно сделать, – это найти такие коэффициенты линейной зависимости, при которых значение функции двух переменных F (a , b) = ∑ i = 1 n (y i - (a x i + b)) 2 будет наименьшим. Иначе говоря, при определенных значениях a и b сумма квадратов отклонений представленных данных от получившейся прямой будет иметь минимальное значение. В этом и состоит смысл метода наименьших квадратов. Все, что нам надо сделать для решения примера – это найти экстремум функции двух переменных.

Как вывести формулы для вычисления коэффициентов

Для того чтобы вывести формулы для вычисления коэффициентов, нужно составить и решить систему уравнений с двумя переменными. Для этого мы вычисляем частные производные выражения F (a , b) = ∑ i = 1 n (y i - (a x i + b)) 2 по a и b и приравниваем их к 0 .

δ F (a , b) δ a = 0 δ F (a , b) δ b = 0 ⇔ - 2 ∑ i = 1 n (y i - (a x i + b)) x i = 0 - 2 ∑ i = 1 n (y i - (a x i + b)) = 0 ⇔ a ∑ i = 1 n x i 2 + b ∑ i = 1 n x i = ∑ i = 1 n x i y i a ∑ i = 1 n x i + ∑ i = 1 n b = ∑ i = 1 n y i ⇔ a ∑ i = 1 n x i 2 + b ∑ i = 1 n x i = ∑ i = 1 n x i y i a ∑ i = 1 n x i + n b = ∑ i = 1 n y i

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

n ∑ i = 1 n x i y i - ∑ i = 1 n x i ∑ i = 1 n y i n ∑ i = 1 n - ∑ i = 1 n x i 2 b = ∑ i = 1 n y i - a ∑ i = 1 n x i n

Мы вычислили значения переменных, при который функция
F (a , b) = ∑ i = 1 n (y i - (a x i + b)) 2 примет минимальное значение. В третьем пункте мы докажем, почему оно является именно таким.

Это и есть применение метода наименьших квадратов на практике. Его формула, которая применяется для поиска параметра a , включает в себя ∑ i = 1 n x i , ∑ i = 1 n y i , ∑ i = 1 n x i y i , ∑ i = 1 n x i 2 , а также параметр
n – им обозначено количество экспериментальных данных. Советуем вам вычислять каждую сумму отдельно. Значение коэффициента b вычисляется сразу после a .

Обратимся вновь к исходному примеру.

Пример 1

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

	i = 1	i = 2	i = 3	i = 4	i = 5	∑ i = 1 5
x i	0	1	2	4	5	12
y i	2 , 1	2 , 4	2 , 6	2 , 8	3	12 , 9
x i y i	0	2 , 4	5 , 2	11 , 2	15	33 , 8
x i 2	0	1	4	16	25	46

Решение

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

Воспользуемся методом наименьших квадратов, чтобы вычислить нужные нам коэффициенты a и b . Для этого подставим нужные значения из последнего столбца и подсчитаем суммы:

n ∑ i = 1 n x i y i - ∑ i = 1 n x i ∑ i = 1 n y i n ∑ i = 1 n - ∑ i = 1 n x i 2 b = ∑ i = 1 n y i - a ∑ i = 1 n x i n ⇒ a = 5 · 33 , 8 - 12 · 12 , 9 5 · 46 - 12 2 b = 12 , 9 - a · 12 5 ⇒ a ≈ 0 , 165 b ≈ 2 , 184

У нас получилось, что нужная аппроксимирующая прямая будет выглядеть как y = 0 , 165 x + 2 , 184 . Теперь нам надо определить, какая линия будет лучше аппроксимировать данные – g (x) = x + 1 3 + 1 или 0 , 165 x + 2 , 184 . Произведем оценку с помощью метода наименьших квадратов.

Чтобы вычислить погрешность, нам надо найти суммы квадратов отклонений данных от прямых σ 1 = ∑ i = 1 n (y i - (a x i + b i)) 2 и σ 2 = ∑ i = 1 n (y i - g (x i)) 2 , минимальное значение будет соответствовать более подходящей линии.

σ 1 = ∑ i = 1 n (y i - (a x i + b i)) 2 = = ∑ i = 1 5 (y i - (0 , 165 x i + 2 , 184)) 2 ≈ 0 , 019 σ 2 = ∑ i = 1 n (y i - g (x i)) 2 = = ∑ i = 1 5 (y i - (x i + 1 3 + 1)) 2 ≈ 0 , 096

Ответ: поскольку σ 1 < σ 2 , то прямой, наилучшим образом аппроксимирующей исходные данные, будет
y = 0 , 165 x + 2 , 184 .

Метод наименьших квадратов наглядно показан на графической иллюстрации. С помощью красной линии отмечена прямая g (x) = x + 1 3 + 1 , синей – y = 0 , 165 x + 2 , 184 . Исходные данные обозначены розовыми точками.

Поясним, для чего именно нужны приближения подобного вида.

Они могут быть использованы в задачах, требующих сглаживания данных, а также в тех, где данные надо интерполировать или экстраполировать. Например, в задаче, разобранной выше, можно было бы найти значение наблюдаемой величины y при x = 3 или при x = 6 . Таким примерам мы посвятили отдельную статью.

Доказательство метода МНК

Чтобы функция приняла минимальное значение при вычисленных a и b , нужно, чтобы в данной точке матрица квадратичной формы дифференциала функции вида F (a , b) = ∑ i = 1 n (y i - (a x i + b)) 2 была положительно определенной. Покажем, как это должно выглядеть.

Пример 2

У нас есть дифференциал второго порядка следующего вида:

d 2 F (a ; b) = δ 2 F (a ; b) δ a 2 d 2 a + 2 δ 2 F (a ; b) δ a δ b d a d b + δ 2 F (a ; b) δ b 2 d 2 b

Решение

δ 2 F (a ; b) δ a 2 = δ δ F (a ; b) δ a δ a = = δ - 2 ∑ i = 1 n (y i - (a x i + b)) x i δ a = 2 ∑ i = 1 n (x i) 2 δ 2 F (a ; b) δ a δ b = δ δ F (a ; b) δ a δ b = = δ - 2 ∑ i = 1 n (y i - (a x i + b)) x i δ b = 2 ∑ i = 1 n x i δ 2 F (a ; b) δ b 2 = δ δ F (a ; b) δ b δ b = δ - 2 ∑ i = 1 n (y i - (a x i + b)) δ b = 2 ∑ i = 1 n (1) = 2 n

Иначе говоря, можно записать так: d 2 F (a ; b) = 2 ∑ i = 1 n (x i) 2 d 2 a + 2 · 2 ∑ x i i = 1 n d a d b + (2 n) d 2 b .

Мы получили матрицу квадратичной формы вида M = 2 ∑ i = 1 n (x i) 2 2 ∑ i = 1 n x i 2 ∑ i = 1 n x i 2 n .

В этом случае значения отдельных элементов не будут меняться в зависимости от a и b . Является ли эта матрица положительно определенной? Чтобы ответить на этот вопрос, проверим, являются ли ее угловые миноры положительными.

Вычисляем угловой минор первого порядка: 2 ∑ i = 1 n (x i) 2 > 0 . Поскольку точки x i не совпадают, то неравенство является строгим. Будем иметь это в виду при дальнейших расчетах.

Вычисляем угловой минор второго порядка:

d e t (M) = 2 ∑ i = 1 n (x i) 2 2 ∑ i = 1 n x i 2 ∑ i = 1 n x i 2 n = 4 n ∑ i = 1 n (x i) 2 - ∑ i = 1 n x i 2

После этого переходим к доказательству неравенства n ∑ i = 1 n (x i) 2 - ∑ i = 1 n x i 2 > 0 с помощью математической индукции.

1. Проверим, будет ли данное неравенство справедливым при произвольном n . Возьмем 2 и подсчитаем:

2 ∑ i = 1 2 (x i) 2 - ∑ i = 1 2 x i 2 = 2 x 1 2 + x 2 2 - x 1 + x 2 2 = = x 1 2 - 2 x 1 x 2 + x 2 2 = x 1 + x 2 2 > 0

У нас получилось верное равенство (если значения x 1 и x 2 не будут совпадать).

1. Сделаем предположение, что данное неравенство будет верным для n , т.е. n ∑ i = 1 n (x i) 2 - ∑ i = 1 n x i 2 > 0 – справедливо.
2. Теперь докажем справедливость при n + 1 , т.е. что (n + 1) ∑ i = 1 n + 1 (x i) 2 - ∑ i = 1 n + 1 x i 2 > 0 , если верно n ∑ i = 1 n (x i) 2 - ∑ i = 1 n x i 2 > 0 .

Вычисляем:

(n + 1) ∑ i = 1 n + 1 (x i) 2 - ∑ i = 1 n + 1 x i 2 = = (n + 1) ∑ i = 1 n (x i) 2 + x n + 1 2 - ∑ i = 1 n x i + x n + 1 2 = = n ∑ i = 1 n (x i) 2 + n · x n + 1 2 + ∑ i = 1 n (x i) 2 + x n + 1 2 - - ∑ i = 1 n x i 2 + 2 x n + 1 ∑ i = 1 n x i + x n + 1 2 = = ∑ i = 1 n (x i) 2 - ∑ i = 1 n x i 2 + n · x n + 1 2 - x n + 1 ∑ i = 1 n x i + ∑ i = 1 n (x i) 2 = = ∑ i = 1 n (x i) 2 - ∑ i = 1 n x i 2 + x n + 1 2 - 2 x n + 1 x 1 + x 1 2 + + x n + 1 2 - 2 x n + 1 x 2 + x 2 2 + . . . + x n + 1 2 - 2 x n + 1 x 1 + x n 2 = = n ∑ i = 1 n (x i) 2 - ∑ i = 1 n x i 2 + + (x n + 1 - x 1) 2 + (x n + 1 - x 2) 2 + . . . + (x n - 1 - x n) 2 > 0

Выражение, заключенное в фигурные скобки, будет больше 0 (исходя из того, что мы предполагали в пункте 2), и остальные слагаемые будут больше 0 , поскольку все они являются квадратами чисел. Мы доказали неравенство.

Ответ: найденные a и b будут соответствовать наименьшему значению функции F (a , b) = ∑ i = 1 n (y i - (a x i + b)) 2 , значит, они являются искомыми параметрами метода наименьших квадратов (МНК).

Если вы заметили ошибку в тексте, пожалуйста, выделите её и нажмите Ctrl+Enter

• Tutorial

Введение

Я математик-программист. Самый большой скачок в своей карьере я совершил, когда научился говорить:«Я ничего не понимаю!» Сейчас мне не стыдно сказать светилу науки, что мне читает лекцию, что я не понимаю, о чём оно, светило, мне говорит. И это очень сложно. Да, признаться в своём неведении сложно и стыдно. Кому понравится признаваться в том, что он не знает азов чего-то-там. В силу своей профессии я должен присутствовать на большом количестве презентаций и лекций, где, признаюсь, в подавляющем большинстве случаев мне хочется спать, потому что я ничего не понимаю. А не понимаю я потому, что огромная проблема текущей ситуации в науке кроется в математике. Она предполагает, что все слушатели знакомы с абсолютно всеми областями математики (что абсурдно). Признаться в том, что вы не знаете, что такое производная (о том, что это - чуть позже) - стыдно.

Но я научился говорить, что я не знаю, что такое умножение. Да, я не знаю, что такое подалгебра над алгеброй Ли. Да, я не знаю, зачем нужны в жизни квадратные уравнения. К слову, если вы уверены, что вы знаете, то нам есть над чем поговорить! Математика - это серия фокусов. Математики стараются запутать и запугать публику; там, где нет замешательства, нет репутации, нет авторитета. Да, это престижно говорить как можно более абстрактным языком, что есть по себе полная чушь.

Знаете ли вы, что такое производная? Вероятнее всего вы мне скажете про предел разностного отношения. На первом курсе матмеха СПбГУ Виктор Петрович Хавин мне определил производную как коэффициент первого члена ряда Тейлора функции в точке (это была отдельная гимнастика, чтобы определить ряд Тейлора без производных). Я долго смеялся над таким определением, покуда в итоге не понял, о чём оно. Производная не что иное, как просто мера того, насколько функция, которую мы дифференцируем, похожа на функцию y=x, y=x^2, y=x^3.

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

Задача на ближайшее время: я поручил своим студентам понять, что такое линейно-квадратичный регулятор . Не постесняйтесь, потратьте три минуты своей жизни, сходите по ссылке. Если вы ничего не поняли, то нам с вами по пути. Я (профессиональный математик-программист) тоже ничего не понял. И я уверяю, в этом можно разобраться «на пальцах». На данный момент я не знаю, что это такое, но я уверяю, что мы сумеем разобраться.

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

Итак, даны две точки (x0, y0), (x1, y1), например, (1,1) и (3,2), задача найти уравнение прямой, проходящей через эти две точки:

иллюстрация

Эта прямая должна иметь уравнение типа следующего:

Здесь альфа и бета нам неизвестны, но известны две точки этой прямой:

Можно записать это уравнение в матричном виде:

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

Давайте заменим конкретные матрицы на их символьное представление:

Тогда (alpha, beta) может быть легко найдено:

Более конкретно для наших предыдущих данных:

Что ведёт к следующему уравнению прямой, проходящей через точки (1,1) и (3,2):

Окей, тут всё понятно. А давайте найдём уравнение прямой, проходящей через три точки: (x0,y0), (x1,y1) и (x2,y2):

Ой-ой-ой, а ведь у нас три уравнения на две неизвестных! Стандартный математик скажет, что решения не существует. А что скажет программист? А он для начала перепишет предыдующую систему уравнений в следующем виде:

В нашем случае векторы i,j,b трёхмерны, следовательно, (в общем случае) решения этой системы не существует. Любой вектор (alpha\*i + beta\*j) лежит в плоскости, натянутой на векторы (i, j). Если b не принадлежит этой плоскости, то решения не существует (равенства в уравнении не достичь). Что делать? Давайте искать компромисс. Давайте обозначим через e(alpha, beta) насколько именно мы не достигли равенства:

И будем стараться минимизировать эту ошибку:

Почему квадрат?

Мы ищем не просто минимум нормы, а минимум квадрата нормы. Почему? Сама точка минимума совпадает, а квадрат даёт гладкую функцию (квадратичную функцию от агрументов (alpha,beta)), в то время как просто длина даёт функцию в виде конуса, недифференцируемую в точке минимума. Брр. Квадрат удобнее.

Очевидно, что ошибка минимизируется, когда вектор e ортогонален плоскости, натянутой на векторы i и j .

Иллюстрация

Иными словами: мы ищем такую прямую, что сумма квадратов длин расстояний от всех точек до этой прямой минимальна:

UPDATE: тут у меня косяк, расстояние до прямой должно измеряться по вертикали, а не ортогональной проекцией. Вот этот комментатор прав.

Иллюстрация

Совсеми иными словами (осторожно, плохо формализовано, но на пальцах должно быть ясно): мы берём все возможные прямые между всеми парами точек и ищем среднюю прямую между всеми:

Иллюстрация

Иное объяснение на пальцах: мы прикрепляем пружинку между всеми точками данных (тут у нас три) и прямой, что мы ищем, и прямая равновесного состояния есть именно то, что мы ищем.

Минимум квадратичной формы

Итак, имея данный вектор b и плоскость, натянутую на столбцы-векторы матрицы A (в данном случае (x0,x1,x2) и (1,1,1)), мы ищем вектор e с минимум квадрата длины. Очевидно, что минимум достижим только для вектора e , ортогонального плоскости, натянутой на столбцы-векторы матрицы A :

Иначе говоря, мы ищем такой вектор x=(alpha, beta), что:

Напоминаю, что этот вектор x=(alpha, beta) является минимумом квадратичной функции ||e(alpha, beta)||^2:

Тут нелишним будет вспомнить, что матрицу можно интерпретирвать в том числе как и квадратичную форму, например, единичная матрица ((1,0),(0,1)) может быть интерпретирована как функция x^2 + y^2:

квадратичная форма

Вся эта гимнастика известна под именем линейной регрессии .

Уравнение Лапласа с граничным условием Дирихле

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

Изначальный коммит доступен . Для минимизации внешних зависимостей я взял код своего софтверного рендерера, уже на хабре. Для решения линейной системы я пользуюсь OpenNL , это отличный солвер, который, правда, очень сложно установить: нужно скопировать два файла (.h+.c) в папку с вашим проектом. Всё сглаживание делается следующим кодом:

For (int d=0; d<3; d++) { nlNewContext(); nlSolverParameteri(NL_NB_VARIABLES, verts.size()); nlSolverParameteri(NL_LEAST_SQUARES, NL_TRUE); nlBegin(NL_SYSTEM); nlBegin(NL_MATRIX); for (int i=0; i<(int)verts.size(); i++) { nlBegin(NL_ROW); nlCoefficient(i, 1); nlRightHandSide(verts[i][d]); nlEnd(NL_ROW); } for (unsigned int i=0; i &face = faces[i]; for (int j=0; j<3; j++) { nlBegin(NL_ROW); nlCoefficient(face[ j ], 1); nlCoefficient(face[(j+1)%3], -1); nlEnd(NL_ROW); } } nlEnd(NL_MATRIX); nlEnd(NL_SYSTEM); nlSolve(); for (int i=0; i<(int)verts.size(); i++) { verts[i][d] = nlGetVariable(i); } }

X, Y и Z координаты отделимы, я их сглаживаю по отдельности. То есть, я решаю три системы линейных уравнений, каждое имеет количество переменных равным количеству вершин в моей модели. Первые n строк матрицы A имеют только одну единицу на строку, а первые n строк вектора b имеют оригинальные координаты модели. То есть, я привязываю по пружинке между новым положением вершины и старым положением вершины - новые не должны слишком далеко уходить от старых.

Все последующие строки матрицы A (faces.size()*3 = количеству рёбер всех треугольников в сетке) имеют одно вхождение 1 и одно вхождение -1, причём вектор b имеет нулевые компоненты напротив. Это значит, я вешаю пружинку на каждое ребро нашей треугольной сетки: все рёбра стараются получить одну и ту же вершину в качестве отправной и финальной точки.

Ещё раз: переменными являются все вершины, причём они не могут далеко отходить от изначального положения, но при этом стараются стать похожими друг на друга.

Вот результат:

Всё бы было хорошо, модель действительно сглажена, но она отошла от своего изначального края. Давайте чуть-чуть изменим код:

For (int i=0; i<(int)verts.size(); i++) { float scale = border[i] ? 1000: 1; nlBegin(NL_ROW); nlCoefficient(i, scale); nlRightHandSide(scale*verts[i][d]); nlEnd(NL_ROW); }

В нашей матрице A я для вершин, что находятся на краю, добавляю не строку из разряда v_i = verts[i][d], а 1000*v_i = 1000*verts[i][d]. Что это меняет? А меняет это нашу квадратичную форму ошибки. Теперь единичное отклонение от вершины на краю будет стоить не одну единицу, как раньше, а 1000*1000 единиц. То есть, мы повесили более сильную пружинку на крайние вершины, решение предпочтёт сильнее растянуть другие. Вот результат:

Давайте вдвое усилим пружинки между вершинами:
nlCoefficient(face[ j ], 2); nlCoefficient(face[(j+1)%3], -2);

Логично, что поверхность стала более гладкой:

А теперь ещё в сто раз сильнее:

Что это? Представьте, что мы обмакнули проволочное кольцо в мыльную воду. В итоге образовавшаяся мыльная плёнка будет стараться иметь наименьшую кривизну, насколько это возможно, касаясь-таки границы - нашего проволочного кольца. Именно это мы и получили, зафиксировав границу и попросив получить гладкую поверхность внутри. Поздравляю вас, мы только что решили уравнение Лапласа с граничными условиями Дирихле. Круто звучит? А на деле всего-навсего одну систему линейных уравнений решить.

Уравнение Пуассона

Давайте ещё крутое имя вспомним.

Предположим, что у меня есть такая картинка:

Всем хороша, только стул мне не нравится.

Разрежу картинку пополам:

И выделю руками стул:

Затем всё, что белое в маске, притяну к левой части картинки, а заодно по всей картинке скажу, что разница между двумя соседними пикселями должна равняться разнице между двумя соседними пикселями правой картинки:

For (int i=0; i

Вот результат:

Код и картинки доступны

Сущность метода наименьших квадратов заключается в отыскании параметров модели тренда, которая лучше всего описывает тенденцию развития какого-либо случайного явления во времени или в пространстве (тренд – это линия, которая и характеризует тенденцию этого развития). Задача метода наименьших квадратов (МНК) сводится к нахождению не просто какой-то модели тренда, а к нахождению лучшей или оптимальной модели. Эта модель будет оптимальной, если сумма квадратических отклонений между наблюдаемыми фактическими величинами и соответствующими им расчетными величинами тренда будет минимальной (наименьшей):

где - квадратичное отклонение между наблюдаемой фактической величиной

и соответствующей ей расчетной величиной тренда,

Фактическое (наблюдаемое) значение изучаемого явления,

Расчетное значение модели тренда,

Число наблюдений за изучаемым явлением.

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

Инструментарий МНК сводится к следующим процедурам:

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

Вторая процедура. Определяется, какая линия (траектория) способна лучше всего описать или охарактеризовать эту тенденцию.

Третья процедура.

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

Таблица 9.1

Номер наблюдения
Урожайность, ц/га

Поскольку уровень технологии при производстве подсолнечника в нашей стране за последние 10 лет практически не изменился, значит, по всей видимости, колебания урожайности в анализируемый период очень сильно зависели от колебания погодно-климатических условий. Действительно ли это так?

Первая процедура МНК. Проверяется гипотеза о существовании тенденции изменения урожайности подсолнечника в зависимости от изменения погодно-климатических условий за анализируемые 10 лет.

В данном примере за «y » целесообразно принять урожайность подсолнечника, а за «x » – номер наблюдаемого года в анализируемом периоде. Проверку гипотезы о существовании какой-либо взаимосвязи между «x » и «y » можно выполнить двумя способами: вручную и при помощи компьютерных программ. Конечно, при наличии компьютерной техники данная проблема решается сама собой. Но, чтобы лучше понять инструментарий МНК целесообразно выполнить проверку гипотезы о существовании связи между «x » и «y » вручную, когда под рукой находятся только ручка и обыкновенный калькулятор. В таких случаях гипотезу о существовании тенденции лучше всего проверить визуальным способом по расположению графического изображения анализируемого ряда динамики - корреляционного поля:

Корреляционное поле в нашем примере расположено вокруг медленно возрастающей линии. Это уже само по себе говорит о существовании определенной тенденции в изменении урожайности подсолнечника. Нельзя говорить о наличии какой-либо тенденции лишь тогда, когда корреляционное поле похоже на круг, окружность, строго вертикальное или строго горизонтальное облако, или же состоит из хаотично разбросанных точек. Во всех остальных случаях следует подтвердить гипотезу о существовании взаимосвязи между «x » и «y », и продолжить исследования.

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

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

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

		Гипербола:

Парабола второго порядка: :

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

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

Нахождение значений параметров уравнения регрессии, в нашем случае параметров и , является сердцевиной МНК. Данный процесс сводится к решению системы нормальных уравнений.

(9.2)

Эта система уравнений довольно легко решается методом Гаусса. Напомним, что в результате решения, в нашем примере, находятся значения параметров и . Таким образом, найденное уравнение регрессии будет иметь следующий вид:

Предыдущая статья: Углерод — характеристика элемента и химические свойства Следующая статья: Углерод — характеристика элемента и химические свойства