Scientific journal
Modern high technologies
ISSN 1812-7320
"Перечень" ВАК
ИФ РИНЦ = 0,899

DESCRIPTIVE SIGNAL SMOOTHING IN A SINGLE ALGORITHM NONPARAMETRIC IDENTIFICATION OF TECHNICAL SYSTEMS

Voskoboynikov Yu.E. 1, 2 Boeva V.A. 1
1 Novosibirsk State University of Architecture and Civil Engineering (Sibstrin)
2 Novosibirsk State Technical University
Mathematical models of many technical systems have the form of Volterra integral equation of the first kind with a difference kernel. For such systems, the non-parametric identification problem reduces to the estimation of pulse transition characteristic of the system from the registered noise-contaminated values of input and output signals. In actual practice, some identification schemes use the square-wave signal with constant amplitude at some instant as the input signal of identification system. For such an input signal, the pulse transition characteristic is an output signal first derivative. Except that, derivative calculation procedure is an ill-posed problem characterized by an essential singularity of estimated derivative non-stability concerning noise terms in output signal. To find a unique stable solution for the problem, for observation data different smoothing algorithms exist, and the most effective is cubic smoothing splines algorithm. There are boundary conditions for these splines coefficients unambiguous determination. Unfortunately, the equal to zero second spline derivatives that are traditional natural boundary conditions not allow taking properties of identification problem into account. What is more, smoothing cubic splines algorithms existed keep from using qualitative a priori or a posteriori information about identification problem required solution when coefficients calculation is. For these reasons, in this paper authors present smoothing cubic splines algorithm with identifying pulse transition characteristic qualitative data taking into account to improve identification procedure accuracy. All the studies investigated prove the confident efficiency of smoothing algorithm proposed and identification procedure as a whole.
identification problem
Volterra integral equation of the first kind
smoothing cubic splines
boundary conditions imposition
taking account of qualitative a priori information

Для моделирования различных динамических систем достаточно часто используются интегральные модели [1, с. 4–5, 163–164]. Если параметры системы не меняются во времени, то такая система называется стационарной, и ее поведение описывается интегральным уравнением Вольтерра I рода с разностным ядром [2, с. 109–128]:

voskob01.wmf (1)

где k(t) – импульсная переходная функция (ИПФ) системы; φ(τ), f(t) – ее входной и выходной сигналы. Для физически реализуемой системы выполняется условие k(t) = 0 при t < 0. Задача непараметрической идентификации системы (1) заключается в вычислении оценки для ИПФ k(t) по измеренным значениям входного и выходного сигналов идентифицируемой системы.

На практике для идентификации стационарных систем часто (в силу своей простоты) применяется схема, в которой на вход идентифицируемой системы подается (в момент t = 0) ступенчатый сигнал, амплитуда которого постоянна. Заметим, что если амплитуда равна 1, то такой сигнал называется функцией Хевисайда (выполнение этого условия предполагается в дальнейшем). Для такого входного сигнала можно показать [2, c. 163–164], что:

voskob02.wmf (2)

где fН(t) – выходной сигнал (реакция системы), если на вход подана функция Хевисайда. Несмотря на хорошо разработанные алгоритмы численного дифференцирования, использование (2) на практике связано с той же проблемой некорректности, что и при решении уравнения (1), так как операция дифференцирования является некорректно поставленной задачей [3]. Чаще всего это проявляется в неустойчивости операции дифференцирования, когда даже небольшой уровень шума регистрации выходного сигнала вызывает очень большие ошибки в полученной (на основе (2)) оценке ИПФ. Для получения устойчивой оценки ИПФ зарегистрированный (с ошибками измерения) сигнал voskob03.wmf первоначально сглаживают (фильтруют), а затем применяют операцию дифференцирования. Часто на практике для этих целей используют (из-за их простоты построения и наличия соответствующего программного обеспечения в ряде математических пакетов) сглаживающие кубические сплайны (СКС) [4; 5]. К сожалению, в получаемых оценках ИПФ при высоком уровне шума (10–15 %) присутствуют (даже при оптимальных параметрах сглаживания) значительные ошибки идентификации (особенно при малых и больших значениях времени t). Поэтому целью данной статьи является разработка способов уменьшения ошибки идентификации путем:

– построения СКС с заданием левого и правого краевых условий в виде значений первых производных;

– использования при построении СКС качественной априорной информации о поведении выходного сигнала fН(t) или ИПФ идентифицируемой системы.

Материалы и методы исследования

Предположим, что зарегистрированный в узлах voskob04.wmf зашумленный выходной сигнал voskob05.wmf допускает представление:

voskob06.wmf, (3)

где ηi – шум с нулевым средним и дисперсией voskob07.wmf. Для сглаживания зашумленных значений voskob08.wmf обратимся к сглаживающим кубическим сплайнам (СКС). Напомним, что сглаживающий кубический сплайн Sf(x) на каждом отрезке voskob09.wmf представляет собой полином третьей степени вида [4; 5]:

voskob10.wmf (4)

и является дважды непрерывно дифференцируемым на всем интервале [0, T]. Для однозначного вычисления коэффициентов СКС ai, bi, ci, di задают левые и правые краевые условия. Например, естественные краевые условия определяют нулевые значения вторых производных в узлах t1, tN, т.е.

voskob11.wmf. (5)

Показано [4; 6], что СКС с естественными краевыми условиями доставляет минимальное значение функционалу:

voskob12.wmf (6)

где voskob13.wmf – весовые множители. Меняя параметр сглаживания α в интервале [0, ∞), можно изменять ошибку сглаживания. Величину параметра сглаживания, минимизирующего среднеквадратическую ошибку (СКО) сглаживания, назовем оптимальным параметром и обозначим αopt. Алгоритмы оценивания αopt рассмотрены в работах [6, с. 62–67; 7].

К сожалению, задание краевых условий voskob14.wmf может привести к значительным ошибкам идентификации (ошибкам дифференцирования СКС) при малых и больших значениях аргумента СКС. Для иллюстрации этого факта на рисунке, а показаны значения «точной» ИПФ (сплошная кривая) и оценки voskob15.wmf (штриховая кривая), вычисленной дифференцированием СКС, построенной при α = αopt и краевых условиях (5), при этом шум измерения имел относительный уровень 0.15.

voskob1a.tif voskob1b.tif

а) б)

Оценки импульсной переходной функции

Видна большая ошибка идентификации при малых значениях t: вместо нулевого значения при t = 0 оценка voskob16.wmf, но по мере увеличения аргумента влияние левого краевого условия уменьшается. Аналогичную ошибку можно наблюдать для значений времени, близких к правому концу интервала: voskob17.wmf. Для уменьшения таких ошибок необходимо задать краевые условия, исходя из специфики решаемой задачи. Например, для ряда динамических систем (колебательные звенья второго порядка и др.) известно, что k(0) = 0, k(T) = 0. Учитывая эту и подобную априорную информацию, целесообразно краевые условия задать значениями первой производной, что, несомненно, приведет к увеличению точности идентификации. Поэтому перейдем к построению алгоритма вычисления коэффициентов сплайна при таких краевых условиях.

Составим систему уравнений (используя результаты работы [6]) для вычисления значений voskob18.wmf, в случае, когда слева в качестве краевых условий будут задаваться величины первых производных:

voskob19.wmf. (7)

Систему уравнений можно представить в следующем виде:

voskob20.wmf (8)

Коэффициенты системы определяются соотношениями, в которых voskob21a.wmf voskob21b.wmf – шаги между соответствующими узлами сетки

voskob22.wmf

После вычисления решения этой системы коэффициенты сплайна определяются выражениями:

voskob23.wmf

Очевидно, что, вычислив коэффициенты СКС, можно найти первую производную для любого значения voskob26.wmf:

voskob27.wmf

На рисунке, а точечной кривой показаны значения ИПФ, вычисленной как производная СКС с краевыми условиями (7). Видно, что точность этой оценки (в граничных точках интервала) существенно выше, чем при использовании краевых условий (5).

Заметим, что и для краевых условий (7) параметр сглаживания можно эффективно выбирать, используя статистический алгоритм, построенный на основе критерия оптимальности и изложенный в работе [6, с. 62–67].

Результаты исследования и их обсуждение

Дескриптивный алгоритм сглаживания. На практике у экспериментатора имеется определенная информация (априорная или апостериорная) о поведении на отдельных интервалах времени либо самой ИПФ, либо выходного сигнала fН(t) идентифицируемой системы. Очевидно, что учет такой нетривиальной достоверной информации может также привести к уменьшению ошибки идентификации ИПФ. К сожалению, изложенный выше алгоритм вычисления коэффициентов сплайна не позволяет учитывать такую информацию. Изложим подход, позволяющий при построении сплайна учесть некоторые формы задания априорной информации.

Запишем функционал (6) и имеющуюся априорную информацию в дискретном виде через значения сглаживающего сплайна voskob29.wmf, которые представим в виде вектора sα с проекциями voskob30.wmf. Можно показать, что слагаемое voskob31.wmf в выражении (6) можно представить выражением:

voskob32.wmf (9)

где матрица Q имеет размеры N×N и определяется как voskob33.wmf. Элементы матрицы H размером voskob34.wmf определяются выражениями:

voskob35.wmf

а элементы трехдиагональной матрицы A размером voskob36.wmf определяются выражениями:

voskob37.wmf;

voskob38.wmf.

Введя диагональную матрицу весов P с элементами voskob39.wmf, функционал (6) (с учетом (9)) можно записать в виде:

voskob40.wmf (10)

Перейдем к построению системы ограничений, позволяющей задать априорную информацию об идентифицируемой ИПФ. Наиболее универсальной формой задания (хорошо согласующейся с функционалом (10)) является система линейных неравенств вида:

voskob41.wmf (11)

где размеры матрицы D и вектора d зависят от вида априорной информации. Например, информацию о знаках ИПФ или ее производных можно задать с помощью конечных разностей соответствующего порядка [8]. Пример такого задания приводится ниже.

Тогда значения дескриптивного (от глагола «описывать») сглаживающего сплайна определяются как решение задачи условной минимизации:

voskob42.wmf (12)

при ограничениях:

voskob43.wmf (13)

После решения этой вариационной задачи по полученным значениям voskob44.wmf дескриптивного сплайна строится интерполяционный сплайн, производная от которого является оценкой для идентифицируемой ИПФ динамической системы.

Заметим, что в состав многих математических пакетов входят модули (например, в пакете MathCAD это функция Minimize и блок решений Given), позволяющие решать задачи условной оптимизации большой размерности [9].

Результаты численных исследований. Для исследования изложенного подхода к построению дескриптивного сплайна вернемся к ранее рассмотренному примеру. Из рисунка, а видно, что на интервале [c, d] построенная оценка ИПФ с условиями (7) (точечная кривая кое-где совпадает со штриховой кривой – краевые условия (5)) должна быть неубывающей функцией. Колебания этой оценки обусловлены влиянием остаточного (после сглаживания) шума измерения. Поэтому зададим эту информацию (по сути дела уже апостериорную) в виде системы ограничений (13) для проекций вектора ИПФ с номерами voskob45.wmf, соответствующих интервалу [c, d]. Тогда ненулевые элементы матрицы D с размерами voskob46.wmf задаются через конечные разности второго порядка следующими соотношениями:

voskob47.wmf

Остальные элементы матрицы D равны 0. Вектор d размерности N – 2 является нулевым вектором. Очевидно, что описанные ограничения соответствуют разностям второго порядка voskob48a.wmf voskob48b.wmf, т.е. значения сплайна для индексов voskob49.wmf не убывают.

На рисунке, б в увеличенном масштабе точечной кривой показаны значения оценки ИПФ, полученной дифференцированием интерполяционного сплайна, построенного по значениям voskob50.wmf, которые являются решением вариационной задачи (12), (13). Видно, что новая оценка ИПФ на интервале [c, d] является неубывающей функцией, что соответствует физическим представлениям об идентифицируемой ИПФ, и эта оценка имеет меньшую ошибку идентификации.

Заключение

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