Кирьянов Д. - MathCad 11 (1077323), страница 63
Текст из файла (страница 63)
Численные методы398• Gfx,o — векторная функция размерности N+I, составленная изфункции пользователя и ее N частных производных по каждому изпараметров с.Листинг 15.14. Регрессия линейной комбинацией функций пользователях : = (01у := ( 4 . 123452 . 434 . 3б )j13 . 65 . 25 . 9х+F{x) :=С :- l i n f i t ( х , у , F)3.957Гч0 . 8 5 45 . 6 0 5 х 10*f ( х ) :=х + С2х+Листинг 15.15. Регрессия общего видах:=(0 1 2 3 4 5 6 )у := ( 4 . 12.434.33.65.2т5.9/с,х\Со • х0.1 "G (х, С} :=хCrxС := genf it (х, у , g , G)3.163С =0.057С,-хf (х) :=С 0 -х15.3. Сглаживание и фильтрацияПри анализе данных часто возникает задача их фильтрации, заключающаясяв устранении одной из составляющих зависимости y ( x i ) .
Наиболее частоГлава 15. Обработка данных399целью фильтрации является подавление быстрых вариаций y(xi), которыечаще всего обусловлены шумом. В результате из быстроосциллирующей зависимости у(хА) получается другая, сглаженная зависимость, в которой доминирует более низкочастотная составляющая.Наиболее простыми и эффективными рецептами сглаживания (smoothing)можно считать регрессию различного вида (см. разд.
15.2). Однако регрессиячасто уничтожает информативную составляющую данных, оставляя лишьнаперед заданную пользователем зависимость.Часто рассматривают противоположную задачу фильтрации — устранениемедленно меняющихся вариаций в целях исследования высокочастотнойсоставляющей. В этом случае говорят о задаче устранения тренда. Иногдаинтерес представляют смешанные задачи выделения сред немасштабных вариаций путем подавления как более быстрых, так и более медленных вариаций. Одна из возможностей решения связана с применением полосовойфильтрации.Несколько примеров программной реализации различных вариантов фильтрации приведены в данном разделе.15.3.1.
Встроенные функции для сглаживанияВ Mathcad имеется несколько встроенных функций, реализующих различные алгоритмы сглаживания данных.• medsmooth(y,b) — сглаживание алгоритмом "бегущих медиан";• ksmooth(x,y,b) — сглаживание на основе функции Гаусса;•I supsmooth(x,y) — локальное сглаживание адаптивным алгоритмом, основанное на анализе ближайших соседей каждой пары данных;• х — вектор действительных данных аргумента (для supsmooth его элементы должны быть расположены в порядке возрастания);• у — вектор действительных значений того же размера, что и х;• ь — ширина окна сглаживания.Все функции имеют в качестве аргумента векторы, составленные из массиваданных, и выдают в качестве результата вектор сглаженных данных того жеразмера. Функция medsmooth предполагает, что данные расположены равномерно.(Примечание)Подробную информацию об алгоритмах, заложенных в функции сглаживания,Вы найдете в справочной системе Mathcad в статье Smoothing (Сглаживание),находящейся в разделе Statistics (Статистика)..Часто бывает полезным совместить сглаживание с последующей интерполяцией или регрессией.
Соответствующий пример приведен в листинге 15.16Часть III. Численные методы400для функции supsmooth. Результат работы листинга показан на рис. 15.18(кружки обозначают исходные данные, крестики — сглаженные, пунктирная кривая — результат сплайн-интерполяции). Сглаживание тех же данныхпри помощи "бегущих медиан" и функции Гаусса с разным значением ширины окна пропускания показаны на рис. 15.19 и 15.20, соответственно.I Листинг 15.16. Сглаживание с последующей сплайн-интерполяциейх:=<0у:=(4.1122.43344.3563.675.289105.954.711412133.53.9143152.7z ;= s u p s m o o t h [х , у)s := c s p l i n e ( x , z)А (t);= i n t e r p ( s , х , z , t )Рис. 15.18.
Адаптивное сглаживание (листинг 15.16)Рис. 15.19. Сглаживание "бегущими медианами"|16 )3.7Т4.85.4)ТГлава 15. Обработка данных401и1гksnoottifx ,у ,2)О О Оfcsnoothlx ,у,5)10121416Рис. 1 5 . 2 0 . Сглаживание при помощи функции k s m o o t h15.3.2. Скользящее усреднениеПомимо встроенных в Mathcad, существует несколько популярных алгоритмов сглаживания, на одном из которых хочется остановиться особо. Самыйпростой и очень эффективный метод — это скользящее усреднение. Егосуть состоит в расчете для каждого значения аргумента среднего значенияпо соседним w данным. Число w называют окном скользящего усреднения;чем оно больше, тем больше данных участвуют в расчете среднего, тем более сглаженная кривая получается.
На рис. 15.21 показан результат скользящего усреднения одних и тех же данных (кружки) с разным окном w=36105 -(J1•"o40It4oooo"o-v-—0\о0310s1100_•,0021о1ISРис. 1 5 . 2 1 . Скользящее усреднение с разными w=3 , 5 , 1 5(листинг 15.17, коллаж трех графиков)Часть III. Численные методы402(пунктир), w=5 (штрихованная кривая) и w=i5 (сплошная кривая). Видно,что при малых w сглаженные кривые практически повторяют ход измененияданных, а при больших w — отражают лишь закономерность их медленныхвариаций.Чтобы реализовать в Mathcad скользящее усреднение, достаточно оченьпростой программы, приведенной в листинге 15.17.
Она использует толькозначения у, оформленные в виде вектора, неявно предполагая, что они соответствуют значениям аргумента х, расположенным через одинаковыепромежутки. Вектор х применялся лишь для построения графика результата(рис. 15.21).Листинг 15.17. Сглаживание скользящим усреднениемх:= { О12345б7891011121314у:={4.1 2.4 3 4.3 3.6 5.2 5.9 5 4.7 4 3.5 3.9 31516)2.7 3.7 4.8 5.4w:= 15N := rows (у)Ы = 17i := 0 .. N - 1iК\А'.- i fw ,iлX *Xj - 0 j - i-w+li + 1ПримечаниеПриведенная программная реализация скользящего усреднения самая простая,но не самая лучшая.
Возможно, Вы обратили внимание, что все кривые скользящего среднего на рис. 15.21 слегка "обгоняют" исходные данные. Почему такпроисходит, понятно: согласно алгоритму, заложенному в последнюю строкулистинга 15.17, скользящее среднее для каждой точки вычисляется путем усреднения значений предыдущих w точек. Чтобы результат скользящего усреднения был более адекватным, лучше применить центрированный алгоритмрасчета по w/2 предыдущим и w/2 последующим значениям. Он будет немногосложнее, поскольку придется учитывать недостаток точек не тольков начале (как это сделано в программе с помощью функции условия i f ) , но ив конце массива исходных данных.15.3.3. Устранение трендаЕще одна типичная задача возникает, когда интерес исследований заключается не в анализе медленных (или низкочастотных) вариаций сигнала у(х)(для чего применяется сглаживание данных), а в анализе быстрых его изменений.
Часто бывает, что быстрые (или высокочастотные) вариации накладываются определенным образом на медленные, которые обычно называютГлава 15. Обработка данных403трендом. Часто тренд имеет заранее предсказуемый вид, например линейный. Чтобы устранить тренд, можно предложить последовательность действий, реализованную в листинге 15.18.1. Вычислить регрессию f{x}, например линейную, исходя из априорнойинформации о тренде (предпоследняя строка листинга).2. Вычесть из данных у (х) тренд ffx) (последняя строка листинга).! Листинг 15.18.
Устранение трендаj10у:-(5.14.4N:=rows (у)i : = 0 ..fft)5.75.35.65.25.97.76.71171213146.55.98.1158.7169.7Т9.810.4)ТN = 17N - 1: - l i n e f x , у) о + l i n e ( х , у) 1 • tНа рис. 15.22 показаны исходные данные (кружками), выделенный с помощью регрессии линейный тренд (сплошной прямой линией) и результатустранения тренда (пунктир, соединяющий крестики).Рис. 15.22. Устранение тренда (листинг 15.18)15.3.4.
Полосовая фильтрацияВ предыдущих разделах была рассмотрена фильтрация быстрых вариацийсигнала (сглаживание) и его медленных вариаций (снятие тренда). Иногдатребуется выделить среднемасштабную составляющую сигнала, уменьшивЧасть III. Численные методы404как более быстрые, так и более медленные его компоненты. Одна из возможностей решения этой задачи связана с применением полосовой фильтрации на основе последовательного скользящего усреднения.Рис. 15.23.
Результат полосовой фильтрации (листинг 15.19)Алгоритм полосовой фильтрации приведен в листинге 15.19, а результат егоприменения показан на рис. 15.23 сплошной кривой. Алгоритм реализуеттакую последовательность операций:1. Приведение массива данных у к нулевому среднему значению путем еговычитания из каждого элемента у (третья и четвертая строки листинга).2. Устранение из сигнала у высокочастотной составляющей, имеющее цельюполучить сглаженный сигнал middle, например с помощью скользящегоусреднения с малым окном w (в листинге 15.19 w=3).3.
Выделение из сигнала middle низкочастотной составляющей slow, например, путем скользящего усреднения с большим окном w (в листинге 15.19 w=7) либо, с помощью снятия тренда (см. разд. 15.3.3).4. Вычитание из сигнала middle тренда slow (последняя строка листинга),тем самым выделяя среднемасштабную составляющую исходного сигнала у.Листинг 15.19. Полосовая фильтрациях:= (Оу : = ( 5 . 1meanY74.4:=mean55.3(у)5.65.285.97.71 01 11 21 31 41 51 6т8.76.77.57.98.18.79.78.81 0 . 4 )тГлава 15. Обработка данных405у := у - meanYN := rows ( у)i :=0 ..N=17N- 1middle^ := ifi < w,j = 0j = i-w+1i + 1w:=7middlei := ifi <w,=Оmiddle•j = i-w+1i + 1m := m i d d l e - slow15.4.
Интегральные преобразованияИнтегральные преобразования массива сигнала у(х) ставят в соответствиевсей совокупности данных у(х) некоторую функцию другой координатыF ( V ) . Рассмотрим встроенные функции для расчета интегральных преобразований, реализованных в Mathcad.15.4.1. Преобразование ФурьеМатематический смысл преобразования Фурье состоит в представлениисигнала у(х) в виде бесконечной суммы синусоид вида Ffv) sin{vx).Функция Ffv) называется преобразованием Фурье или интегралом Фурье, илиФурье-спектром сигнала.
Ее аргумент v имеет смысл частоты соответствующей составляющей сигнала. Обратное преобразование Фурье переводитспектр F(V) в исходный сигнал у(х}. Согласно определению,F{v)у(х) • e x p (-ivx) d x .Как видно, преобразование Фурье является существенно комплексной величиной, даже если сигнал действительный.Преобразование Фурье действительных данныхПреобразование Фурье имеет огромное значение для различных математических приложений, и для него разработан очень эффективный алгоритм,Часть III. Численные методы406называемый БПФ (быстрым преобразованием Фурье). Этот алгоритм реализован в нескольких встроенных функциях Mathcad, различающихся нормировками.Оfft(y) — вектор прямого преобразования Фурье;•FFT(y) — вектор прямого преобразования Фурье в другой нормировке;Оifft(v) — вектор обратного преобразования Фурье;•ке;IFFT(V)— вектор обратного преобразования Фурье в другой нормиров-у — вектор действительных данных, взятых через равные промежуткизначений аргумента;v — вектор действительных данных Фурье-спектра, взятых через равные промежутки значений частоты.Внимание!Аргумент прямого Фурье-преобразования, т.
е. вектор у, должен иметь ровно 2"элементов ( п — целое число). Результатом является вектор с 1+2"' 1 элементами. И наоборот, аргумент обратного Фурье-преобразования должен иметь1+2" 1 элементов, а его результатом будет вектор из 2" элементов. Если числоданных не совпадает со степенью 2, то необходимо дополнить недостающиеэлементы нулями.Рис. 15.24. Исходные данные и обратное преобразование Фурье(листинг 15.20}Пример расчета Фурье-спектра для суммы трех синусоидальных сигналовразной амплитуды (показанных в виде сплошной кривой на рис.