6.4. Спектральный анализ

Спектральный анализ, предназначенный для определения передаточных функций линейных систем, возник при использовании статистических методов в спектральном оценивании. Этот метод хорошо изложен в книгах Дженкинса и Уоттса [193, гл. 10] и Бриллингера [63, гл. 6]; кроме того, метод широко обсуждается во многих других учебниках по анализу временных рядов. В данном разделе применяется несколько необычный подход, состоящий в получении стандартных процедур посредством сглаживания эмпирической оценки передаточной функции.

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

истинная передаточная функция является гладкой функцией

Далее, если частотный интервал мал в смысле незначительного изменения на нем то

являются некоррелированными и несмещенными оценками примерно одной и той же константы дисперсии которых равны по лемме 6.1

Здесь мыпренебрегли членами, стремящимися к нулю при

Если предположить, что постоянна на интервале

то, как хорошо известно, наилучшей (в смысле минимума дисперсии) оценкой этой константы является средне-взвешенное «измерений» (6.39) для частот (6.40), причем каждое измерение взвешивается в соответствии с обратной величиной ее дисперсии (ср. с задачей и леммой II.2, (11.58) Приложения II):

Для больших можно с хорошей степенью аппроксимации использовать интегралы, которые соответствуют суммам (по Риману) в (6.41а):

Если передаточная функция не постоянна на интервале (6.40), разумно использовать дополнительное взвешивание, которое выделяет близкие к частоты:

Здесь сосредоточенная вокруг функция, а - «параметр формы», кратко обсуждаемый ниже.

Очевидно, (6.42) соответствует функции

Далее, если спектр шума известен, оценка (6.43) может быть реализована. Если же спектр не известен, можно принять следующее допущение: предположим, что спектр шума не претерпевает значительных изменений на частотных интервалах, соответствующих "ширине” функции взвешивания

Тогда можно заменить на что приводит к сокращению константы в (6.43). При выполнении условия (6.45) оценка

является, таким образом, хорошим приближением оценки (6.43),

Следует заметить, что при невыполнении (6.45) может быть лучше добавить процедуру оценивания и использовать эту оценку в (6.43).

Связь с процедурой Блэкмена-Тьюки. Рассмотрим знаменатель в (6.46). Он представляет собой средне-взвешенное периодограммы Используя результат (2.72), получаем при

где спектр определенный соотношениями Более того, если

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

Аналогично, поскольку

числитель в (6.46)

является оценкой взаимного сиектра между входной и выходной переменными.

Таким образом, оценка передаточной функции (6.46) представляет собой отношение двух спектральных оценок

что имеет определенный смысл в свете (2.78). Спектральные оценки (6.48) и (6.50) являются стандартными, известными в литературе оценками спектра и взаимного спектра как сглаженных периодограмм; см. работы Блэкмена и Тьюки [53], Джен кинса и Уоттса [193] или Бриллингера [63].

Существует общий способ выражения этих оценок. Коэффициенты ряда Фурье для периодограммы равны

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

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

Поскольку интеграл (6.48) представляет собой свертку, его коэффициенты ряда Фурье равны произведеиию (6.52) и поэтому разложение в ряд Фурье функции (6.48) имеет вид

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

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

Весовая функция: частотное окно. Обсудим теперь весовую функцию В спектральном анализе ее часто называют частотным окном (аналогично, называют временным окном.) Если это окно "широкое”, то будет взвешиваться большое количество различных частот в интервале (6.40). Это должно привести к малой дисперсии оценки . В то же время широкое окно будет включать оценки для частот, существенно далеких от со средними значениями, которые могут значительно отличаться от Это вызовет большое смещение. Таким образом, ширина окна позволит управлять компромиссом между смещением и дисперсией. Чтобы сделать это управление компромиссом несколько более формальным, используем для описания ширины скаляр у, причем большая величина у соответствует узкому окну.

Таблица 6.1 . (см. скан) Некоторые частотные окна для спектрального анализа

Будем характеризовать окно следующими числами:

При возрастании у (и сужении частотного окна) значение уменьшается, растет.

Некоторые типичные окна представлены в таблице 6.1. (Другие примеры можно найти в [63, табл. 3.3.1].) Заметим, что в (6.55) масштабирующий параметр у был выбран так, что Частотные окна показаны на рис. 6.1. Для этих окон имеем

Бартлет:

Парзен:

Хэмминг:

Рис. 6.1. Некоторые простые частотные окна. Сплошная линии - окно Парнеза, штриховая линия - окно Хэмминга, точечная линия - окно Бартлетта;

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

Асимптотические свойства сглаженной оценки. Оценка (6.46) изучалась в различных работах но спектральному анализу. См., например, книги Дженкинса и Уоттса [193, гл. 10] и Бриллингера [63, гл. 60]. Результаты, являющиеся асимптотическими как по так и по 7, могут быть получены следующим образом (см. Приложение Рассмотрим оценку (6.64) и предположим, что истинная система удовлетворяет условиям леммы 6.1. Тогда имеем смещение

штрих и двойной штрих обозначают дифференцирование один раз и дважды соответственно, дисперсия

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

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

Здесь

Можно показать также справедливость следующих дополнительных результатов (см. [63, гл. 6] и задачи

1. Оценки асимптотически некоррелированы и каждая имеет дисперсию, равную половине величины (6.59). (6.62)

2. Оценки на разных частотах асимптотически некоррелированы. (6.63)

3. Оценки для произвольного набора частот имеют асимптотически нормальное совместное распределение со средними и дисперсиями, определяемыми

4. По поводу свойств см. задачу

Из (6.60) видно, что желаемое свойство окна состоит в малости как так и Можно также вычислить значение параметра ширины окна , минимизирующее Предположим, что как , так и неограниченно возрастают, а отношение стремится к нулю и, таким образом, можно использовать асимптотические выражения. Допустим также, что имеет место (6.57), т. е. Тогда из (6.60) следует

Рис. 6.2. Смоделированные данные системы (6.67)

Рис. 6.3. График амплитуды эмпирической оценки передаточной функции, основанной на данных, представленных на рис. 6.2. Ровная линия - истинный график амплитуды для системы (6.67)

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

Оптимальный выбор у приводит к среднеквадратичной ошибке, убывающей как

При практическом использовании результаты (6.65) и (6.66) формально не могут быть достигнуты. Вместо этого типичная процедура должна состоять в выборе начального значения (см. табл. 6.1) и последующем вычислении и нанесении на график соответствующих оценок для различных значений у. При возрастании у появляется все больше и больше характерных деталей оценки. Это происходит благодаря уменьшению смещения (более ясное проявление истинных резонансных пиков и т.п.), а также в связи с возрастанием дисперсии (ложные, случайные пики). Процедура должна быть остановлена, когда пользователь почувствует появление преобладающих ложных деталей.

Пример 6.1. Моделировалась система.

где белый шум с дисперсией 1, а в качестве входного воздействия использовалось более 1000 реализаций псевдослучайного двоичного сигнала (см. раздел 14.3). Часть записанных в результате данных показана на рис. 6.2. Соответствующие эмпирические оценки передаточной функции изображены на рис. 6.3. Оценка формировалась в соответствии с (6.46) для окна Парзеиа с различными значениями 7. На рис. 6.4 представлены результаты для . Здесь разумным выбором ширины окна оказывается

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

(кликните для просмотра скана)

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

Затем оценка может быть образована как среднее арифметическое

или как средневзвешенное в соответствии с обратными величинами дисперсий:

с

представляющей собой периодограмму подгруппы. Обратная величина дисперсии оценки равна но сокращается при образовании (6.70).

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