Метод наименьших квадратов Гаусса с весовыми коэффициентами отклонений / Хабр

В моей практике метод наименьших квадратов Гаусса используется в двух случаях.

Когда производится измерение, для корректировки полученной величины.

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

В качестве примера возьму измерение МДС (магнитодвижущая сила) срабатывания геркона.

Согласно ГОСТ 25810 «Контакты магнитоуправляемые герметизированные. Методы измерения электрических параметров» МДС срабатывания определяется по значению тока, протекающего через измерительную катушку в момент срабатывания геркона.

Если упрощенно, в катушке создается нарастающий со скоростью 1A/мc ток до срабатывания геркона. Величина этого тока и есть МДС срабатывания.

Графически корректировку результатов измерения можно представить следующим образом.

Ось X – то, что измерено, ось Y – то, что должно быть.

В метрологической службе берутся эталонные герконы. Эталонный геркон – это геркон с измеренной метрологами МДС.

,,,, – МДС 1, 2, 3, 4, 5-го эталонных герконов.

Измеряем 5 эталонных герконов.

,,,,– измеренная МДС 1, 2, 3, 4, 5-го геркона.

Для справки. Метод наименьших квадратов Гаусса позволяет заменить кривую результатов прямой (рис. 1).

Рисунок 1. Прямая, вычисленная по методу Гаусса

Уравнение прямой:

Коэффициенты a, b однозначно определяют прямую.

Формулы Гаусса для a, b:

$a=\frac{n\sum{Y_iX_i}-\sum{X_i}\sum{Y_i}}{n\sum{X_i^2}-(\sum{X_i})^2}$ $b=\frac{\sum{Y_i}-a\sum{X_i}}{n}$

В нашем примере n=5,– измеренная МДС,- эталонная.

Таким образом, измерив МДС геркона, вносим поправку, получаем скорректированное значение:

В другом случае, при задании, к примеру, тока, решаем обратную задачу. Если необходимо задать ток, тогда в ЦАП заносим код тока равного:

$X_k=\frac{Y_k-b}{a}$

Вернемся к измерению МДС.

Погрешность измерения за счет влияния внешних электрических и магнитных полей не должна превышать 0.5А, и не должна быть более 2%.

Практика показала, что чем меньше величина МДС, тем меньшая погрешность допускается.

Моя идея – это ввести весовые коэффициенты отклонений в формулу критерия Гаусса. При этом функция, весовых коэффициентов измерений – это тоже прямая.

Тогда получаем.

Критерий метода наименьших квадратов Гаусса: сумма квадратов отклонений точек кривой от точек прямой минимальна.

$F(a,b)=\sum^{\substack{n}}_{\substack{i=1} }(Y_i-\overline{Y_i})^2$

Где– эталонная величина, $\overline{Y_i}$ – точка на прямой Гаусса.

Если определить , как вес отклонения для, то уравнение Гаусса примет вид:

$F(a,b)=\sum^{\substack{n}}_{\substack{i=1} }K_i(Y_i-\overline{Y_i})^2$

Получаем уравнения для a, b прямой:

$a=\frac{\sum{K_i}\sum{K_iY_iX_i}-\sum{K_iX_i}\sum{K_iY_i}} {\sum{K_i}\sum{K_iX_i^2}-(\sum{K_iX_i})^2}$ $b=\frac{\sum{K_iY_i}-a\sum{K_iX_i}} {\sum{K_i}}$

Пример. Допустим- линейная функция:

Пусть для двух крайних точек прямой:

Весовые коэффициенты заданы как:

$K_1-100\%=1$ $K_5-20\%=0.2$

Тогда получаем c, d прямой весовых коэффициентов по двум точкам:

$c=\frac{K_5-K_1} {X_5-X_1}=\frac{0.2-1}{65-8}\approx-0.014$ $d=K_1-cX_1=1+0.014*8\approx1.112$

Калькулятор (рис. 2) и тексты макросов демонстрируют идею метода наименьших квадратов Гаусса с весовыми коэффициентами отклонений.

Private Sub CommandButton1_Click()

MNK

MNKk

End Sub

' Метод наименьших квадратов Гаусса

Sub MNK()

Dim i

Dim n

Dim ny, nx, nyx, nx2

Dim a, b

n = 5

ny = 0

nx = 0

nyx = 0

nx2 = 0

For i = 3 To n + 2

ny = ny + Cells(i, 2)

nx = nx + Cells(i, 3)

nyx = nyx + Cells(i, 2) * Cells(i, 3)

nx2 = nx2 + Cells(i, 3) * Cells(i, 3)

Next

a = (n nyx - nx ny) / (n nx2 - nx nx)

b = (ny - a * nx) / n

Cells(3, 5) = a

Cells(3, 6) = b

End Sub

' Метод наименьших квадратов Гаусса с весовыми коэффициентами отклонений

Sub MNKk()

Dim i

Dim n

Dim nk, nky, nkx, nkyx, nkx2

Dim a, b

n = 5

nk = 0

nky = 0

nkx = 0

nkyx = 0

nkx2 = 0

For i = 3 To n + 2

nk = nk + Cells(i, 15)

nky = nky + Cells(i, 15) * Cells(i, 13)

nkx = nkx + Cells(i, 15) * Cells(i, 14)

nkyx = nkyx + Cells(i, 15) Cells(i, 13) Cells(i, 14)

nkx2 = nkx2 + Cells(i, 15) Cells(i, 14) Cells(i, 14)

Next

a = (nk nkyx - nkx nky) / (nk nkx2 - nkx nkx)

b = (nky - a * nkx) / nk

Cells(3, 17) = a

Cells(3, 18) = b

End Sub

Приложение №1. Метод наименьших квадратов Гаусса

$F(a,b)=\sum^{\substack{n}}_{\substack{i=1} }(Y_i-\overline{Y_i})^2$ $=\sum^{\substack{n}}_{\substack{i=1} }(Y_i-(aX_i+b))^2$ $=\sum^{\substack{n}}_{\substack{i=1} } (Y_i^2-2Y_i(aX_i+b)+(aX_i+b)^2)$ $=\sum^{\substack{n}}_{\substack{i=1} } (Y_i^2-2Y_iaX_i-2Y_i b+a^2X_i^2+2aX_ib+b^2)$ $\frac {\partial F(a,b)} {\partial a} =\sum^{\substack{n}}_{\substack{i=1} } (-2Y_iX_i + 2aX_i^2 + 2X_ib) = 0$ $- \sum^{\substack{n}}_{\substack{i=1} } Y_iX_i + \sum^{\substack{n}}_{\substack{i=1} } aX_i^2 + \sum^{\substack{n}}_{\substack{i=1} } X_ib = 0$ $\frac {\partial F(a,b)} {\partial b} = \sum^{\substack{n}}_{\substack{i=1} } (-2Y_i + 2aX_i + 2b) = 0$ $\sum^{\substack{n}}_{\substack{i=1} } (-Y_i + aX_i + b) = 0$ $\sum^{\substack{n}}_{\substack{i=1} } -Y_i + \sum^{\substack{n}}_{\substack{i=1} } aX_i + \sum^{\substack{n}}_{\substack{i=1} } b = 0$ $\sum^{\substack{n}}_{\substack{i=1} } -Y_i + a\sum^{\substack{n}}_{\substack{i=1} } X_i + nb = 0$ $b = \frac {\sum Y_i -a\sum X_i } {n}$ $na\sum X_i^2 – a(\sum X_i)^2 = n\sum Y_iX_i - \sum X_i \sum Y_i$ $a = \frac {n\sum Y_i X_i - \sum X_i \sum Y_i} {n\sum X_i^2 – (\sum X_i)^2}$

Приложение №2. Метод наименьших квадратов Гаусса с весовыми коэффициентами отклонений

$F(a,b)=\sum^{\substack{n}}_{\substack{i=1} } K_i(Y_i-\overline{Y_i})^2$

– вес отклонения для

$=\sum^{\substack{n}}_{\substack{i=1} } K_i(Y_i – (aX_i+b))^2$ $=\sum^{\substack{n}}_{\substack{i=1} } (K_iY_i^2 –2K_iY_iaX_i-2K_iY_ib+K(i)a^2X_i^2+2K_iaX_ib+K_ib^2 )$ $\frac {\partial F(a,b)} {\partial A} = \sum^{\substack{n}}_{\substack{i=1} } (-2K_iY_iX_i + 2K_iaX_i^2 + 2K_iX_ib) = 0$ $-\sum K_iY_iX_i + a\sum K_iX_i^2 + b\sum K_iX_i = 0$ $\frac {\partial F(a,b)} {\partial b} = \sum^{\substack{n}}_{\substack{i=1} } (-2K_iY_i + 2K_iaX_i + 2K_ib) = 0$ $\sum –K_iY_i + a\sum K_iX_i + b\sum K_i = 0$ $b = \frac {\sum K_iY_i - a\sum K_iX_i } {\sum K_i}$ $a = \frac {\sum K_i \sum K_iY_iX_i - \sum K_iX_i \sum K_iY_i} {\sum K_i \sum K_iX_i^2 – (\sum K_iX_I)^2}$