Кирьянов. Самоучитель MathCad 11. Кирьянов д в
 Скачать 10.75 Mb.
|
|
Рис. 15.12. Линейная регрессия (листинг 15.7 или Для расчета линейной регрессии в Mathcad имеются два дублирующих друг друга способа. Правила их применения представлены в листингах и 15.8. Результат обоих листингов получается одинаковым (рис. 15.12). — вектор из двух элементов коэффициентов линейной регрессии 390 Часть III. Численные методы intercept — коэффициент ь линейной регрессии slope — коэффициента линейной регрессии; • х — вектор действительных данных аргумента; • у — вектор действительных данных значений того же размера. Листинг Линейная регрессиях ( 0 1 2 3 4 у ( 4.1 2.4 3 4.3 3.6 5.2 line (x , у 0.414 f := line у) + line (x, у) Листинг Другая форма записи линейной регрессии х := ( 0 1 2 3 4 у := { 4 . 1 2 . 4 3 4 . 3 3 . 6 5 . 2 i n t e r c e p t ( x , уху c e p t ( х , у ) + s l o p e ( x , у ) В Mathcad имеется альтернативный алгоритм, реализующий не минимизацию суммы квадратов ошибок, а энную линейную регрессию для расчета коэффициентов аи (листинг 15.9). — вектор из двух элементов (Ь,а) коэффициентов линейной регрессии • х,у — векторы действительных данных одинакового размера. Листинг 15.9. Построение линейной регрессии двумя разными методами (продолжение листинга 15.7) 2.517 (x, у) = 0.55 g it , у) + medf it (x , у) Различие результатов среднеквадратичной и медиан-медианной регрессии иллюстрируется 15.13. Глава 15. Обработка данных 1 5 . 1 3 . Линейная регрессия по методу наименьших квадратов и методу медиан (листинги 15.7 и 15.9) 15.2.2. Полиномиальная регрессия В реализована регрессия одним полиномом, отрезками нескольких полиномов, а также двумерная регрессия массива данных. Полиномиальная регрессия Полиномиальная регрессия означает приближение данных полиномом степени 15.14). При к=1 ном является прямой линией, при к — параболой, при к=з — кубической параболой и т. д. Как правило, на практике применяются Примечание Для построения регрессии полиномом к-й степени необходимо наличие по крайней мере точек данных. В Mathcad полиномиальная осуществляется комбинацией встроенной функции regress и полиномиальной интерполяции (см. разд. 15.1.2). • — вектор коэффициентов для построения полиномиальной регрессии данных interp(s,x,y, t) — результат полиномиальной регрессии; • • х — вектор действительных данных аргумента, элементы которого расположены в порядке возрастания; • у — вектор действительных данных значений того же размера Часть III. Численные методы к — степень полинома регрессии (целое положительное число — значение аргумента полинома регрессии. Внимание! Для построения полиномиальной регрессии после функции r e g r e s s Вы обязаны использовать функцию i n t e r p Рис. 15.14. Регрессия полиномами разной степени (коллаж результатов листинга 15.10 для разных к) Пример полиномиальной регрессии квадратичной параболой приведен в листинге 15.10. Листинг 15.10. Полиномиальная регрессиях тку, к) А ( t) i n t e r p ( s , х, у , Регрессия отрезками полиномов Помимо приближения массива данных одним полиномом имеется возможность осуществить регрессию сшивкой отрезков (точнее говоря, участков, т. кони имеют криволинейную форму) нескольких полиномов. Для этого имеется встроенная функция loess, применение которой аналогично функции И 15.15). Глава Обработка данных — вектор коэффициентов для построения регрессии данных отрезками полиномов t) — результат полиномиальной регрессии; • • х — вектор действительных данных аргумента, элементы которого расположены в порядке возрастания; • у — вектор действительных данных значений того же размера — параметр, определяющий размер отрезков полиномов (положительное число, хорошие результаты дает значение порядка Параметр span задает степень сглаженности данных. При больших значениях регрессия практически не отличается от регрессии одним полиномом (например, span=2 дает почти тот же результат, что и приближение точек параболой Листинг 15.11. Регрессия отрезками полиномов х : = ( О у : = { 4 . 1 2 . 4 3 4 . 3 s := l o e s s ( x , у 0 . 9 Ат. Совет одним полиномом эффективна, когда множество точек выглядит как полинома регрессия отрезками полиномов оказывается полезной в противоположном случае 5 4 О 2 О 3 о 4 / 5 G 15.15. Регрессия отрезками полиномов (листинг 15.11) 394 Часть III. Численные методы Двумерная полиномиальная регрессия По аналогии с одномерной полиномиальной регрессией и двумерной интерполяцией (см. разд. 15.1.5) Mathcad позволяет приблизить множество точек поверхностью, которая определяется многомерной полиномиальной зависимостью. В качестве аргументов встроенных функций для построения полиномиальной регрессии должны стоять в этом случае не векторы, а соответствующие матрицы regress — вектор коэффициентов для построения полиномиальной регрессии данных span) — вектор коэффициентов для построения регрессии данных отрезками полиномов — скалярная функция, аппроксимирующая данные выборки двумерного поля по координатам хи у кубическими сплайнами вектор вторых производных, созданный одной из сопутствующих функций loess ИЛИ х — матрица размерности определяющая пары значений аргумента (столбцы соответствуют метками у — вектор действительных данных размерности N; • span параметр, определяющий размер отрезков полиномов; • к — степень полинома регрессии (целое положительное число — вектор из двух элементов, содержащий значения аргументов хи у, для которых вычисляется интерполяция. Внимание! Для построения регрессии не предполагается никакого предварительного упо- рядочивания данных (как, например, для двумерной интерполяции, которая требует их представления в виде матрицы NXN). В связи с этим данные представляются как вектор. Двумерная полиномиальная регрессия иллюстрируется листингом ирис. 15.16. Сравните стиль представления данных для двумерной регрессии с представлением тех же данных для двумерной сплайн-интерполяции (см. листинги ее результаты с исходными данными (см. рис. 15.10) и их сплайн-интерполяцией (см. рис. Листинг Двумерная полиномиальная регрессия 0 0 0 0 1 1 1 1 1 2 2 2 2 2 3 3 3 3 3 4 4 4 4 0 1 0 2 0 3 0 4 0 0 1 0 2 0 3 0 4 0 0 1 0 2 0 3 0 4 0 0 1 0 2 0 3 0 4 0 0 1 0 2 0 3 0 4 0 Глава 15. Обработка данных : = ( l 1 0 1 . 1 1 . 2 1 2 3 2 . 1 1 . 5 1 . 3 3 . 3 5 1 . 7 2 1 . 3 3 3 . 7 2 . 1 2 . 9 1 . 5 2 S r e g r e s s ( X , Z , 3 ) i .. к О .. к := interp i — - к • 40 к Примечание Обратите внимание на знаки транспонирования в листинге. Они применены для корректного представления аргументов (например z в качестве вектора, а не строки). Рис. 15.16. Двумерная полиномиальная регрессия (листинг 15.12) 15.2.3. Регрессия специального вида Кроме рассмотренных, в Mathcad встроено еще несколько видов трехпара- метрической регрессии. Их реализация несколько отличается от приведенных выше вариантов регрессии тем, что для них, помимо массива данных, требуется задать некоторые начальные значения коэффициентов Используйте соответствующий вид регрессии, если хорошо представляете себе, какой зависимостью описывается Ваш массив данных. Когда тип регрессии плохо отражает последовательность данных, то ее результат часто бывает неудовлетворительными даже сильно различающимся в зависимости от вы 396 Часть Численные методы бора начальных значений. Каждая из функций выдает вектор уточненных параметров (x,y,g) — регрессия экспонентой х функцией f (x,y,g) регрессия синусоидой +c; pwfit(x,y,g) — регрессия степенной функцией (x,y,g) — регрессия логарифмической функцией f (x) +с — регрессия логарифмической функцией х вектор действительных данных аргумента; • у — вектор действительных значений того же размера — вектор из трех элементов, задающий начальные значения Примечание Правильность выбора начальных значений можно оценить по результату регрессии если функция, выданная Mathcad, хорошо приближает зависимость ух, значит они были подобраны удачно. Пример расчета одного из видов трехпараметрической регрессии (экспоненциальной) приведен в листинге 15.13 и на рис. 15.17. В предпоследней строке листинга выведены в виде вектора вычисленные коэффициенты а,ь,с, а в последней строке через эти коэффициенты определена искомая функция f (х). Листинг 15.13. Экспоненциальная регрессиях 1 2 3 4 у ( 2.4 3 4.3 3.6 5.2 С := expf it ( худ С 0.544 3.099 f (t) • • + Глава 15. Обработка данных 397 Примечание Многие задачи регрессии данных различными ми зависимостями ух) можно свести к более надежной, с вычислительной точки зрения, линейной регрессии. Делается это с помощью соответствующей замены пере- менных. б 5 4 3 - - - 0 о 1 О 2 О 3 с,х о - 5 Рис. 15.17. Экспоненциальная регрессия (листинг 15.13) Регрессия общего вида В Mathcad можно осуществить регрессию в виде линейной комбинации где — любые функции пользователя, a — подлежащие определению коэффициенты. Кроме того, имеется путь проведения регрессии более общего вида, когда комбинацию функций и искомых коэффициентов задает сам пользователь. Приведем встроенные функции для регрессии общего вида и примеры их использования (листинги 15.14 и надеясь, что читатель при необходимости найдет более подробную информацию об этих специальных возможностях в справочной системе и Mathcad О (x,y,F) — вектор параметров линейной комбинации функций пользователя, осуществляющей регрессию данных — вектор параметров, реализующих регрессию данных с помощью функций пользователя общего видах — вектор действительных данных аргумента, элементы которого расположены в порядке возрастания; • у — вектор действительных значений того же размера) — пользовательская векторная функция скалярного аргумента — вектор начальных значений параметров регрессии размерности N; Часть III. Численные методы — векторная функция размерности составленная из функции пользователя и ее N частных производных по каждому из параметров с. Листинг Регрессия линейной комбинацией пользователях бух +С l i n f i t ( х , ух х ) х +х +Листинг 15.15. Регрессия общего х:=(0 1 2 3 4 5 6 тух, С := х С := it (х, у , g , С = 0.057 f (х Сглаживание и фильтрация При анализе данных часто возникает задача их фильтрации, заключающаяся в устранении одной из составляющих зависимости Наиболее часто Глава 15. Обработка данных целью фильтрации является подавление быстрых вариаций которые чаще всего обусловлены шумом. В результате из быстроосциллирующей висимости получается другая, сглаженная зависимость, в которой доминирует более низкочастотная составляющая. Наиболее простыми и эффективными рецептами сглаживания (можно считать регрессию различного вида (см. разд Однако регрессия часто уничтожает информативную составляющую данных, оставляя лишь наперед заданную пользователем зависимость. Часто рассматривают противоположную задачу фильтрации — устранение медленно меняющихся вариаций в целях исследования высокочастотной составляющей. В этом случае говорят о задаче устранения тренда. Иногда интерес представляют смешанные задачи выделения вариаций путем подавления как более быстрых, таки более медленных вариаций. Одна из возможностей решения связана с применением полосовой фильтрации. Несколько примеров программной реализации различных вариантов фильтрации приведены в данном разделе. Встроенные функции для сглаживания В Mathcad имеется несколько встроенных функций, реализующих различные алгоритмы сглаживания данных — сглаживание алгоритмом "бегущих медиан — сглаживание на основе функции Гаусса — локальное сглаживание адаптивным алгоритмом, основанное на анализе ближайших соседей каждой пары данных вектор действительных данных аргумента (для его элементы должны быть расположены в порядке возрастания); • у — вектор действительных значений того же размера, что их ширина окна сглаживания. Все функции имеют в качестве аргумента векторы, составленные из массива данных, и выдают в качестве результата вектор сглаженных данных того же размера. Функция medsmooth предполагает, что данные расположены равномерно Примечание Подробную информацию об алгоритмах, заложенных в функции сглаживания, Вы найдете в справочной системе Mathcad в статье Smoothing (Сглаживание), находящейся в разделе Statistics (Статистика).. Часто бывает полезным совместить сглаживание с последующей интерполяцией или регрессией. Соответствующий пример приведен в листинге 15.16 Часть III. Численные методы для функции Результат работы листинга показан на рис. кружки обозначают исходные данные, крестики — сглаженные, пунктирная кривая — результат сплайн-интерполяции). Сглаживание тех же данных при помощи "бегущих медиан" и функции Гаусса с разным значением ширины окна пропускания показаны на 15.19 и 15.20, соответственно Листинг 15.16. Сглаживание с последующей сплайн-интерполяцией х : = < 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 у : = ( 4 . 1 2 . 4 3 4 . 3 3 . 6 5 . 2 5 . 9 5 4 . 7 4 3 . 5 3 . 9 3 2 . 7 3 . 7 4 . 8 z ; = s u p s m o o t h [ х , у ) s : = c s p l i n e ( x , z Ах Адаптивное сглаживание (листинг 15.16) 15.19. Сглаживание "бегущими медианами Глава 15. Обработка данных , у , 2 ООО г 14 16 15.20. Сглаживание при помощи функции k s m o o t h 15.3.2. Скользящее усреднение Помимо встроенных в Mathcad, существует несколько популярных алгоритмов сглаживания, на одном из которых хочется остановиться особо. Самый простой и очень эффективный метод — это скользящее усреднение. Его суть состоит в расчете для каждого значения аргумента среднего значения по соседним данным. Число w называют окном скользящего усреднения; чем оно больше, тем больше данных участвуют в расчете среднего, тем более сглаженная кривая получается. На 15.21 показан результат скользящего усреднения одних и тех же данных (кружки) с разным окном 4 3 2 - \ 0 0 0 0 0 s J •"o 0 о 10 •, 0 о 0_ o" Рис. Скользящее усреднение с разными w=3 листинг 15.17, коллаж трех графиков Часть III. Численные методы (пунктир), (штрихованная кривая) и (сплошная кривая. Видно, что при малых сглаженные кривые практически повторяют ход изменения данных, а при больших w — отражают лишь закономерность их медленных вариаций. Чтобы реализовать в Mathcad скользящее усреднение, достаточно очень простой программы, приведенной в листинге 15.17. Она использует только значения у, оформленные в виде вектора, неявно предполагая, что они соответствуют значениям аргументах, расположенным через одинаковые промежутки. Вектор х применялся лишь для построения графика результата (рис. Листинг 15.17. Сглаживание скользящим усреднением { Об у 2.4 3 4.3 3.6 5.2 5.9 5 4.7 4 3.5 3.9 3 2.7 3.7 4.8 5.4 15 N rows (у 17 i 0 .. N - 1 i i '.- i f w , X X j - 0 j - i + 1 Примечание Приведенная программная реализация скользящего усреднения самая простая, но не самая лучшая. Возможно, Вы обратили внимание, что все кривые скользящего среднего на рис. 15.21 слегка "обгоняют" исходные данные. Почему так происходит, понятно согласно алгоритму, заложенному в последнюю строку листинга 15.17, скользящее среднее для каждой точки вычисляется путем усреднения значений предыдущих w точек. Чтобы результат скользящего усреднения был более адекватным, лучше применить центрированный алгоритм расчета предыдущими последующим значениям. Он будет немного сложнее, поскольку придется учитывать недостаток точек не только вначале (как это сделано в программе с помощью функции условия i f ) , но ив конце массива исходных данных Устранение тренда Еще одна типичная задача возникает, когда интерес исследований заключается не в анализе медленных (или низкочастотных вариаций сигнала у(х) (для чего применяется сглаживание данных, а в анализе быстрых его изменений. Часто бывает, что быстрые (или высокочастотные вариации накладываются определенным образом медленные, которые обычно называют Глава 15. Обработка данных 403 трендом. Часто тренд имеет заранее предсказуемый вид, например линейный. Чтобы устранить тренд, можно предложить последовательность действий, реализованную в листинге 15.18. 1. Вычислить регрессию f{x}, например линейную, исходя из априорной информации о тренде (предпоследняя строка листинга. Вычесть изданных (х ffx) (последняя строка листинга Листинг Устранение тренда j 10 11 12 13 14 15 Т 4 . 4 5 . 7 5 . 3 5 . 6 5 . 2 5 . 9 7 . 7 6 . 7 7 6 . 5 8 . 1 8 . 7 9 . 7 9 . 8 1 0 . 4 ) ( у ) = 1 7 i 0 .. N - 1 f f t ) у ) + l i n e у ) • t Т На рис. 15.22 показаны исходные данные (кружками, выделенный с помощью регрессии линейный тренд (сплошной прямой линией) и результат устранения тренда (пунктир, соединяющий крестики 15.22. Устранение тренда (листинг 15.18) Полосовая фильтрация В предыдущих разделах была рассмотрена фильтрация быстрых вариаций сигнала (сглаживание) и его медленных вариаций (снятие тренда). Иногда требуется выделить среднемасштабную составляющую сигнала, уменьшив Часть III. Численные методы как более быстрые, таки более медленные его компоненты. Одна извоз- можностей решения этой задачи связана с применением полосовой фильтрации на основе последовательного скользящего усреднения. Рис. 15.23. Результат полосовой фильтрации (листинг Алгоритм полосовой фильтрации приведен в листинге 15.19, а результат его применения показан на рис. 15.23 сплошной кривой. Алгоритм реализует такую последовательность операций. Приведение массива данных у к нулевому среднему значению путем его вычитания из каждого элемента у (третья и четвертая строки листинга. Устранение из сигнала у высокочастотной составляющей, имеющее целью получить сглаженный сигнал middle, например с помощью скользящего усреднения с малым окном (в листинге. Выделение из сигнала middle низкочастотной составляющей slow, например, путем скользящего усреднения с большим окном w (в листинге с помощью снятия тренда (см. разд. 15.3.3). 4. Вычитание из сигнала middle тренда slow (последняя строка листинга), тем самым выделяя среднемасштабную составляющую исходного сигнала у. |