Использование метода наименьших квадратов

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

Предполагается, что процесс описывается одномерным уравнением 2-го порядка

W = a₀ + a₁x + a₂x², 0 £ x £ 6.

Считаем, что величина х измеряется точно, а W – с ошибкой e, имеющей нормальное распределение с нулевым математическим ожиданием и единичной дисперсией
М(e) = 0, s²(e) = 1.

 

 

Выборка десяти случайных пар ( ) представлена в табл. 3.1 в графах 2 и 3.

 

Таблица 3.1

№	x		Wm	e
	4,8608 4,2396 2,7792 0,5988 3,2136 4,5156 5,9340 1,5852 4,4880 4,0932	9,28 9,40 7,88 1,86 7,77 8,73 8,33 5,16 7,28 9,22	8,848 8,821 7,460 2,039 8,056 8,874 8,118 4,994 8,872 8,767	0,432 0,579 0,420 -0,179 -0,286 -0,144 0,212 0,166 -1,592 0,453

Метод наименьших квадратов заключается в том, что неизвестные (искомые) коэффициенты а₀ , а₁ , а₂ должны минимизировать функцию, представляющую собой сумму квадратов невязок e_j:

.

 

Минимум некоторой функции, как известно, находится в точке , где все частные производные этой функции по переменным а₀, а₁, а₂равны нулю.

Для определения частных производных, распишем функцию G через ее предполагаемый вид:

.

Возьмем от функции G производные по а₀, а₁, а₂ :

;

;

.

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

Решая полученную систему, получим а₀ = –0,161; а₁ = 3,929; а₂ = –0,427.

Таким образом, математическая модель будет иметь вид

Wm = –0,161 + 3,929 x –0,427x². (3.2)

Проверим адекватность модели методом Фишера. Для этого заполним четвертый и пятый столбцы таблицы 3.1, подставляя в математическую модель (3.2) и затем в формулу (3.1) значения x_j из первого столбца.

Определим число степеней свободы системы по формуле

f_s = n – m – 1,

где n = 10 – количество экспериментальных точек; m = 3 – количество неизвестных коэффициентов. То есть f_s = 6.

Выборочная дисперсия вычисляется по формуле

.

Критерий Фишера вычисляется по формуле

.

По статистическим таблицам при 5%-м уровне риска (a = 0,05) находим пороговое значение критерия Фишера

.

Так как полученное значение F меньше критического (порогового), гипотеза об адекватности модели реальному процессу принимается.

 

3.3. Статистические методы проверки адекватности
математических моделей

 

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

Математически задача проверки адекватности модели формулируется как задача проверки предположения о том, что значение отклика модели W_m отличается от реального отклика системы W не более чем на заданную величину e*:

. (3.3)

Однако, истинное значение отклика системы никогда неизвестно. Полученный в результате эксперимента отклик в силу неконтролируемого дрейфа системы, разброса характеристик ее элементов и, наконец, просто ошибок измерения представляет собой случайную величину, отличающуюся от W. Поэтому при сравнении результатов математического и физического экспериментов будет получена совокупность случайных величин {e_i}: , среди которых могут оказаться как величины, удовлетворяющие условию (3.3), так и не удовлетворяющие ему.

Можно ли считать, что полученные отклонения (e_i > e*) объясняются случайными причинами или их наличие должно быть признано существенным, что приводит к отказу от проверяемой модели. Для решения этого вопроса на основе выборки случайных величин {e_i} строят статистические критерии, по которым оценивают адекватность модели.

Гипотеза об адекватности модели действительности (гипотеза Н₀) может быть сформулирована как предположение о том, что полученная совокупность {e_i} не дает оснований отказаться от рассматриваемой модели. Иными словами, модель удовлетворяет заданной точности e*.

Альтернативная гипотеза Н₁ состоит в том, что модель не отвечает заданным требованиям (3.3) и, следовательно, должна быть отвергнута.

Так как выборка {e_i} случайна, решение о выборе одной из гипотез Н₀ или Н₁ носит вероятностный характер. При этом может быть допущена ошибка первого рода, состоящая в отказе от правильной модели (принимается Н₁, когда верна Н₀), или ошибка второго рода, состоящая в принятии ошибочной модели (принимается Н₀, когда верна Н₁). Вероятность ошибки первого рода обозначают через a, второго рода – b. Принято называть a риском разработчика, b – риском потребителя. Разумеется, желательно минимизировать как a, так и b. Однако, при заданном объеме экспериментальной выборки уменьшение a влечет за собой увеличение b.

На практике a задается на определенном уровне (a = 0,05; 0,01; 0,005; 0,001), при этом в 100a% случаев правильная модель отвергается.

Величина 1– b характеризует вероятность отказа от ошибочной модели, называется мощностью критерия и является мерой его эффективности.

Выбор вероятностей ошибок a и b при проверке конкретной модели зависит от ответственности решений, принимаемых на основе моделирования.

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

Для оценки гипотезы об адекватности модели существует несколько критериев:

1) Критерий согласия c² Пирсона.

2) Критерий Смирнова-Колмогорова.

3) Критерий Фишера и др.

При использовании критерия c² проверке подлежит гипотеза о том, что рассматриваемая модель адекватна исследуемой системе с вероятностью р (например, р = 0,95). Это значит, что при n независимых испытаниях np значений e_i должно удовлетворять условию (3.3) и лишь в (1– р)п случаях это условие может быть нарушено.

В результате случайного эксперимента для этих событий будут получены частоты n₁ и n₂: n₁» рп; n₂ » (1– р)п; (n₁ + n₂ = п).

Частоты n₁ и n₂ отличаются от точных вероятностных оценок или из-за несоответствия модели действительности (заданная вероятность р не соблюдается), или из-за случайных отклонений.

Для оценки предположения о том, что отклонения n₁ и n₂ от соответствующих вероятностей случайны, строится функция

,

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

Полученное значение U ^* сравнивается с табличным значением при заданном уровне риска a . Если U ^* превышает пороговое значение , модель должна быть отвергнута, и принимается гипотеза Н₁. Если U ^*£ , экспериментальные данные не противоречат гипотезе об адекватности модели, и принимается гипотеза Н₀.

Необходимым условием использования критерия c² является многочисленность экспериментальных данных (не меньше 20).

Критерий Смирнова-Колмогорова основан на максимальном значении отклонений

.

Для заданной экспериментальной выборки строится вспомогательная функция

,

которая сравнивается с пороговым значением l_n,a , определенным по таблицам распределения функции Смирнова-Колмогорова.

При модель должна быть отвергнута, а при экспериментальные данные не противоречат гипотезе об адекватности модели.

Критерий Смирнова-Колмогорова целесообразно использовать при относительно малых выборках, когда критерий c² оказывается неэффективным.

Критерий Фишера осуществляется путем анализа дисперсий. Если дисперсия, характеризующая ошибку эксперимента s²(W), известна, вычисляется выборочная дисперсия S ²(e) и составляется F-отношение:

.

Полученную величину F-отношения сравнивают с пороговым значением критерия Фишера F_{f s,¥, a} при заданном уровне риска a.

При F_{f s,¥} £ F_{f s,¥, a} полученная величина S ²(e) может быть объяснена случайным разбросом экспериментальных данных и, следовательно, нет оснований для отказа от проверяемой модели.

Если F_{f s,¥} > F_{f s,¥, a} , полученное расхождение результатов моделирования и экспериментальных данных знáчимо и, следовательно, модель должна быть отвергнута как недостаточно точная.

3.4. Идентификация параметров математической модели
силы резания токарной операции

 

Построим математическую модель силы резания при обработке круглой детали на токарном станке (Рис. 3.6).

 

 

Сила резания Р описывается математической моделью в виде позинома

P = C S^a V^b t^g, (3.4)

где S – продольная подача; V – скорость резания; t – глубина резания; С, a, b, g – неизвестные параметры.

Формула (3.4) является справочной. Для определения неизвестных параметров воспользуемся методом наименьших квадратов.

 

Пусть проведено n экспериментов, результаты которых сведены в таблицу 3.2.

 

Таблица 3.2

№	S(мм/об)	V(мм/с)	t(мм)	P(Кг)
2   n	S1 S2   Sn	V1 V2   Vn	t1 t2   tn	P1 P2   Pn

 

Для упрощения решения поставленной задачи прологарифмируем выражение (3.4):

.

Введем обозначения

Тогда формула (3.4) преобразуется к линейному виду:

. (3.5)

Метод наименьших квадратов сведется к минимизации функции

,

где – логарифмы экспериментальных значений силы резания, взятых из табл. 3.2;
– логарифмы силы резания, предсказанные с помощью математической модели (3.5), – логарифмы экспериментальных значений подачи, скорости и глубины резания, взятых из той же табл. (3.2), – логарифм неизвестного параметра С.

Возьмем производные от функции G по и приравняем их к нулю:

Разделим обе части уравнений на –2; вынесем , a, b, g за знак суммы; перенесем члены, не зависящие от , a, b, g, в правую часть:

(3.6)

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

Для определения параметров исходной модели (3.4) необходимо для коэффициента С (только для него) проделать операцию, обратную логарифмированию – потенцирование: С = . Коэффициенты a, b, g получаются непосредственно из решения системы (3.6).

Если в распоряжении исследователя имеются экспериментальные данные, для проверки адекватности математической модели действительности можно использовать методы математической статистики. Рассматриваемый ниже метод пригоден при изучении любых математических моделей. Однако конкретный анализ проводится на примере построенной модели силы резания при точении с помощью критерия согласия c², предложенного Пирсоном.

Гипотеза Н₀ формулируется как предположение о том, что отклонение e экспериментальных данных от значений, предсказанных моделью (3.4), с вероятностью р (доверительная вероятность) укладываются в некоторый толерантный интервал ±e*. Если это предположение правильно, то в толерантный интервал (Р ± e*) должно укладываться np отклонений e_i = |Р_i – Р_mi|. Вне толерантного интервала должно оказаться (1– p)n отклонений. Для ограниченной случайной выборки из n наблюдений эти события будут наблюдаться с частотой n₁ и n₂ , лишь приближенно совпадающие с соответствующими вероятностями:

n₁ » pn; n₂ » (1– p)n; n₁ + n₂ = n.

Необходимо установить, можно ли объяснить эти отклонения случайными причинами (в этом случае можно принять гипотезу Н₀) или же они не случайны – статистически значимы (в этом случае нужно принять альтернативную гипотезу Н₁).

Для этого вычисляется некоторая величина U, называемая статистикой:

.

Эту величину нужно сравнить с пороговым значением c²-критерия (c²_1,a) при принятом уровне риска a. Если U £ c²_1,a , наблюдаемые отклонения частот от соответствующих вероятностей можно объяснить случайностью и нет оснований для отказа от нуль-гипотезы Н₀. Если U > c²_1,a , то или произошло маловероятное событие (1– р), или наблюдаемые отклонения не случайны. В этом случае принимается гипотеза Н₁.

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

1) Можно изменить исходные предположения с тем, чтобы увеличить толерантный интервал ±e* или уменьшить доверительную вероятность р. При этом умéньшатся отклонения n₁ и n₂ от соответствующих вероятностей, и проверка может привести к принятию гипотезы Н₀. В этом случае моделью можно пользоваться, но нужно признать, что ее точность оказалась ниже, чем первоначально предполагалось.

2) Можно уменьшить уровень риска a (то есть вероятность отказа от правильной модели в результате неудачного эксперимента). Это приводит к увеличению порогового значения c²_1,a. Это, в свою очередь, может изменить оценку значения U. Однако нужно помнить, что при этом увеличивается риск признать правильной ошибочную модель.

3) Можно потребовать увеличения объема выборки, что, разумеется, приведет к увеличению точности оценки модели и уменьшению риска ошибок.

При проверке адекватности моделей действительности всегда рассматривается случай, когда за пределами толерантного интервала оказалось больше точек, чем ожидалось (n₁ < pn; n₂ > (1– pn)n). В противном случае опасений за точность модели не возникает, однако можно предположить, что величина толерантного интервала задана необоснованно большой. Если в результате проверки по критерию c² в этом случае будет получена величина U > c²_1,a , то завышение толерантного интервала (или занижение доверительной вероятности р) статистически значимо, и необходимо уменьшить e* или увеличить р.
В обоих случаях нужно признать, что модель оказалась точнее, чем ожидалось.