Файл: Математическое моделирование в экологии.doc

где τ — переменная интегрирования;

α — коэффициент;

у — зависимая переменная, которая называется состоянием системы;

t — независимая переменная, которая может быть не только временем;

у₀ — не зависящее от времени начальное условие;

x(τ) — детерминированная входная функция (возмущающая сила).

Для того чтобы получить наблюдаемую зависимую переменную Y(t), к функции у(τ) следует добавить ненаблюдаемую ошибку ε(t). Для дискретных наблюдений

Y(t_i) = y(t_i)+ ε(t_i),

а для непрерывных переменных

Y(_i) = y(t)+ ε(t)

Если параметр, α заменить его оценкой , то остаточная ошибка определяется

Е(t) = Y(t) - (t).

Целью оценивания параметров является получение значения параметра, а в процессе наблюдений Y(t_i) и Y(t). Чтобы сделать это, необходимо знать функцию х(t) и иметь некоторую информацию о характере ε(t).

Более общей является модель, содержащая систему линейных обыкновенных дифференциальных уравнений первого порядка с постоянными коэффициентами (рис. 12.6)

.

Рис. 12.6. Многомерный процесс с несколькими входами

В матричной форме модель имеет вид:

, .

; Y= ; Х= .

Тогда решение можно записать в форме

.

Модель, содержащую одно или несколько линейных дифференциальных уравнений более высокого порядка с постоянными коэффициентами, например

,

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

Тогда уравнение второго порядка будет представлено в виде двух уравнений первого порядка

, ,

однако функция, а в действительности является производной.

Общая нелинейная модель первого порядка имеет вид

,

где F(α,Y,t) представляет собой весьма общую нелинейную функцию. Уравнение, за редким исключением, не имеет аналитического решения и его следует решать численными методами.

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

эксперименты должны быть поставлены так, чтобы измерялся вектор производных dy/dt а не сам вектор Y. Если наблюдаемой переменной является производная, то процедура оценивания вообще не затрагивает дифференциального уравнения; параметры и начальные условия можно оценить методами регрессионного анализа.

Другой способ (менее удовлетворительный) позволяет избежать операций с производными при оценивании. Он состоит в использовании численных значений производных, полученных по наблюдениям величины Y.

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

Исследуем численную оценку. Численное дифференцирование детерминирова-нных переменных предусматривает вычисление dy/dt или высших производных при некотором произвольном значении независимой переменной t (напри-мер, t₀) по заданному ряду значений у в некотором интервале вблизи t₀. Детерминированная ошибка в производной становится меньше, если данные концентрируются около значения t₀, расположенного в середине интервала изменения t, чем когда значение t₀ попадает на тот или другой конец интервала. Вместо непрерывной производной можно использовать любой из интерполяционных многочленов (Лагранжа, Грама).

Большинство схем аппроксимации для производных можно записать в общей форме

)

где а_i— постоянные;

D — дифференциальный оператор;

k — порядок производной;

h — фиксированное приращение независимой переменной;

(m+1) — число используемых опорных точек.

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

.

Если принять дисперсии всех величин у_i, равными друг другу, то

.

Отсюда видно, что чем меньше интервал h и чем больше членов в формуле, тем больше ошибка в производной. Эта ошибка растет с увеличением порядка производной.

Пример. В табл. 12.4 представлены численные ошибки для не- скольких разностных формул, использованных для вычисления температурного градиента в стенке (dТ/dz)z₀ по измеренным значениям температуры, которая считалась детерминированной переменной.

Таблица 12.4

Численные ошибки использованных для вычисления данных

Здесь аппроксимация градиента d/Т/dz и «порядок ошибки» маскируют значительно более важный источник ошибки при вычислении dT/dz, а именно случайную ошибку, порождаемую измерением.

Предположим, что все величины Т_i в табл. 12.4 обладают одинаковым стандартным отклонением 0,01 (1% от Т₀) или дисперсией в 10^-4. Тогда для многочлена четвертого порядка

при измерениях: T_o =-1,000; Т₁=-0,588; Т₂=-0,295; Т₃=-0,259; T₄=-0,305; h = 0,1. Тогда , что составляет около 16% от .

Оценивание методом наименьших квадратов. Если наблюдения Y для откликов модели представляют собой непрерывные функции времени в интервале от t = 0 до t = t_i то МНК требует минимизировать величину

,

где Г — ковариационная матрица;

Ф — проинтегрированное по времени значение квадрата ошибки. Если наблюдения производились в дискретные моменты времени t_i, i= 1, 2, ..., n, то, согласно критерию Маркова, следует минимизировать величину

.

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

Величина ψ в общей форме

ψ (α,у₀,t_i) = Y(t_i) - ε(t_i),

ε(t_i) = Y(t_i) - ψ (t_i).

т.е.

Для минимизации дифференцируем функцию Ф по у₀ и α.

Приравниваем ее нулю

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

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

Пример. Пусть имеем модель

.

Тогда

.

Дифференцируя

по у₀ а затем по α и заменяя в получившихся выражениях у₀ и α на их оценки получим: