Результаты моделирования представлены на рис. 2. 7.
Полученный график оценки спектральной плотности шума с доверительным коридором показывает наличие равномерного распределения мощности сигнала по частотному спектру с последующим спадом мощности на частоте выше 1, 1 рад/с.
Далее необходимо выполнить параметрическое оценивание ТОУ.
2. 8. Параметрическое оценивание данных
Параметрическое оценивание экспериментальных данных проводится с целью определения параметров модели заданной структуры путем минимизации выбранного критерия качества модели (чаще всего – среднего квадрата рассогласования выходов объекта и его постулируемой модели).
Для проведения параметрического оценивания
массив экспериментальных данных необходимо разделить условно на две части (не обязательно равные) >> zdanv=zdan(1:500);
>> zdane=zdan(501:1000).
Первая часть массива данных будет использоваться для параметрического оценивания и построения модели ТОУ. Вторая часть необходима будет для верификации (проверки качества) модели, определения адекватности полученной модели и определения погрешностей идентификации. Необходимо отметить, что параметрическая идентификация в пакете System Identification Toolbox MATLAB выполняется в дискретном виде и полученные модели, являются дискретными.
В пакете System Identification Toolbox рассмотрены различные виды моделей ТОУ, которые с различной степенью достоверности описывают объект автоматизации. Для выбора наиболее приемлемой структуры и вида моделей при параметрическом оценивании экспериментальных данных в пакете System Identification Toolbox MATLAB имеются специальные функции
· параметрического оценивания,
· задания структуры модели,
· изменения и уточнения структуры модели и выбора структуры модели.
Функция оценивания ar оценивает параметры модели авторегресии:
A(z) y(t) = e(t),
где: A(z) = 1 + a1z – 1 + a2z – 2 +...+ a naz – na , т.е. коэффициенты полинома A(z), при моделировании скалярных временных последовательностей. Функция имеет синтаксис:
th = ar(y,n)
Или другое написание, позволяющее изменять параметры моделирования:
[th,refl]=ar(y,n,approach,win,maxsize,T)
где аргументы:
y – вектор-столбец данных, содержащий N элементов;
n – порядок модели (число оцениваемых коэффициентов);
approach – аргумент (строковая переменная) определяет метод оценивания:
• 'ls' – метод наименьших квадратов;
• 'yw' – метод Юла-Уокера;
• 'burg' – метод Бэрга (комбинация метода наименьших квадратов с минимизацией гармонического среднего);
• 'gl' – метод с использованием гармонического среднего;
Если любое из данных значений заканчивается нулем (например, burg0), то вычисление сопровождается оцениванием корреляционных функций.
Аргумент win (строковая переменная) используется в случае отсутствия части данных:
• win =‘now’ – используются только имеющиеся данные (используется по умолчанию, за исключением случая approach = ‘yw’);
• window = 'prw’ – отсутствующие начальные данные заменяются нулями, так что суммирование начинается с нулевого момента времени;
• window = 'pow’ – последующие отсутствующие данные заменяются нулями, так что суммирование расширяется до момента времени N + n;
• window = ‘ppw’ – и начальные, и последующие отсутствующие данные заменяются нулями (используется в алгоритме Юла-Уокера);
Арумент maxsize определяет максимальную размерность задачи; Т – интервал дискретизации.
Возвращаемые величины:
• th – полученная модель авторегрессии в тета-формате (внутреннем матричном формате представления параметрических моделей пакета System Identification);
• relf – информация о коэффициентах и функции потерь;
Для использования функция параметрического оценивания ar необходимо из массива экспериментальных данных, записанных в файле dan выделить выходную переменную у с помощью команды>> y=dan.y,
что равносильно команде
>> y=get(dan,'y')
>> th =ar(y,4)
Discrete-time IDPOLY model: A(q)y(t) = e(t)
A(q) = 1 - 1.348 q^-1 + 0.6695 q^-2 - 0.2531 q^-3 - 0.04431 q^-4
Estimated using AR ('fb'/'now') from data set y
Loss function 0.0140559 and FPE 0.0141588
Sampling interval: 1
Полная информация о модели авторегрессии th может быть получена с помощью команды:
>> present(th)
Discrete-time IDPOLY model: A(q)y(t) = e(t)
A(q) = 1 - 1.348 (+-0.03022) q^-1 + 0.6695 (+-0.05018) q^-2
- 0.2531 (+-0.05018) q^-3 - 0.04431 (+-0.03023) q^-4
Estimated using AR ('fb'/'now') from data set y
Loss function 0.0140559 and FPE 0.0141588
Sampling interval: 1
Created: 17-Dec-2009 10:51:00
Last modified: 17-Dec-2009 10:51:00
В информации приведены сведения о том, что модель является дискретной и для оценивания ее параметров используется прямой-обратный метод (разновидность метода наименьших квадратов), на что указывает строковая переменная 'fb' (используется по умолчанию); для построения модели используются только имеющиеся данные у, на что указывает строковая переменная 'now' (используется по умолчанию); определены: функция потерь Loss function, как остаточная сумма квадратов ошибки и так называемый теоретический информационный критерий Акейке (Akaike's Information Theoretic Criterion – AIC) FPE; интервал дискретизации Sampling interval.
Следующая функция arx оценивает параметры модели AR и ARX: параметры модели ARX, представленной зависимостью:
A(z)y(t) = B(z) u(t) + e(t) ,
или в развернутом виде:
y(t) + а1y(t-1) + …+ аna y(t-n) = b1 u(t) + b2 u(t - 1) + …+ bnb u(t - m) + e(t).
Здесь и ниже e(t) – дискретный белый шум.
B(z) = b1 + b2 z-1 + …+ bbn z-nb + 1
Функция имеет следующий синтаксис:
dar = arx(z,nn).
Или другое написание, позволяющее изменять параметры моделирования:
dar = arx(z,nn,maxsize,T),
где z – экспериментальные данные;
nn – задаваемые параметры модели (аргумент nn содержит три параметра: na – порядок ( число коэффициентов) полинома A(z); nb – порядок полинома B(z); nk – величина задержки;
maxsize - максимальная размерность задачи;
Т – интервал дискретизации.
Естественно задаться вопросом: Какую степень полинома выбрать? Известно, что с увеличением порядка полиномов улучшается степень адекватности модели реальному объекту. Однако при этом получаются громоздкие выражения, и увеличивается время моделирования. Поэтому для нахождения оптимального порядка полиномов можно воспользоваться функциями выбора структуры модели:
Функция arxstruc вычисляет функции потерь для ряда различных конкурирующих ARX моделей с одним выходом:
v = arxstruc(ze,zv,NN),
или v = arxstruc(ze,zv,NN, maxsize);
где: ze,zv – соответственно, матрицы экспериментальных данных для оценивания и верификации моделей; NN – матрица задания конкурирующих структур со строками вида nn = [na nb nk];
maxsize - максимальная размерность задачи.
Возвращаемая величина v – матрица, верхние элементы каждого столбца которой (кроме последнего) являются значениями функции потерь для ARX моделей, структура которых отображается последующими элементами столбцов (т. е. каждый столбец соответствует одной модели). Первый элемент последнего столбца – это число значений экспериментальных данных для верификации моделей.
