Аппроксимация законов распределения ортогональными полиномами Текст научной статьи по специальности «Математика»

С.А. Прохоров, И.А. Лезин, И.В. Солдатова

АППРОКСИМАЦИЯ ЗАКОНОВ РАСПРЕДЕЛЕНИЯ ОРТОГОНАЛЬНЫМИ ПОЛИНОМАМИ

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

При обработке статистического материала часто приходится решать вопрос о том, как подобрать для данного статистического ряда теоретическую кривую распределения, выражающую лишь существенные черты статистического материала. Такая задача, относящаяся к классу аппроксимативных задач, в данном случае называется задачей выравнивания (сглаживания) статистических рядов [1].

Авторами разработана автоматизированная система аппроксимативного анализа плотностей распределения вероятностей, являющаяся дальнейшим развитием разработанных ранее систем, используемая как при подготовке специалистов в области обработки информации и управления, так и при проведении научных исследований [2, 3].

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

первичной статистической обработки временных рядов;

идентификации законов распределения по фазовым портретам;

аппроксимации плотностей распределения ортогональными полиномами Лежандра, Чебышева, Лагерра, Эрмита [4];

сплайн-интерполяции; работы с базой данных.

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

Последовательность этапов решения задачи сглаживания представлена на рис. 1 [2].

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

нижеследующие этапы.

1. Определить диапазон изменения статистическо-

2. При неизвестном законе распределения определить ширину дифференциального коридора:

Р и с. 1. Блок-схема сглаживания статистических рядов

х — х

д X___ тах тгп

~ м

(1)

где М - количество дифференциальных коридоров [1].

го ряда хтгп тах

Из выражения (1) следует, что значение величины М не определено и, следовательно, исследователь должен задать её самостоятельно, ориентируясь на рекомендации специалистов. Как следует из литературных источников, значение М рекомендуется выбирать в диапазоне от 10 до 20 [1].

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

Ах =

' 86 Их!

I 4 '\г.

(2)

где 6 - погрешность аппроксимации функции распределения; "(х) тшх - максимум модуля вто-

рой производной функции распределения [2];

М = еп(

Ах

(3)

3. Определить частоту попадания анализируемой случайной величины в у-тый дифференциальный коридор:

- 1 N

р1 = N16

(4)

где

1, если еП 1

х, - х„

Ах

если xi = у Ах;

2

0, иначе

(5)

- индикатор состояния.

Следует отметить, что 6^ у+1 = 1/2, если х,=)Ах л х^хтшх, т.е. в этом случае в у и у+1 коридоры добавляется по 1/2.

4. Если частота попадания в какой-либо к-тый дифференциальный коридор мала (д<0,01*0,02), то для уменьшения влияния случайности его объединяют с к+1 коридором. Эта операция может быть применена неоднократно.

Исходя из условия нормировки

М А

Xк) • Ах = 1, (6)

)=1

вычисляем высоту столбцов:

н, = А

Ах

(7)

Далее на интервале [хшп,хтшх1

с ] строится линейная сплайн-модель. Для этого введем массив

пар {ху’Уу }, = 1’М , где х, - середина,-того интервала гистограммы, а у, = Н, - соответственно высота у-того столбца. Значения сплайн-модели рассчитываются следующим образом:

У1 +(У2 - У1 )/(х2 - х1 Мх - х1 )>х < х1>

/,.М (х) = 1 у, +(У,+1 - уу )/(х)+1 - ху Мх - ху) х) < х < х,+1’ (8)

.Ум-1 + (Ум - Ум-1)/(хм - хм-1) •(х - хм-1 )’х > хм-Здесь М (х) - линейный сплайн, построенный на М дифференциальных коридорах, который

представляет собой модель аппроксимируемой плотности вероятности.

Очевидно, что внешний вид сплайн-модели сильно зависит от количества дифференциальных коридоров, причем утверждать однозначно, при каком именно значении М сплайн-модель наилучшим образом соответствует истинной плотности вероятности, невозможно. Поэтому зачастую целесообразнее применять механизм “усреднения” модели плотности вероятности по нескольким сплайн-моделям.

хх

тшх тгп

С использованием (8) для расчета сплайн-модели плотности вероятности усредненная модель определяется по формуле

М +1

X ^, т (х)

т /-Л____ т=М-I

І,,М (х )= 21 +1 ■ (9)

В данном случае в точке х суммируются значения 21+1 сплайн-моделей, построенных на одних и тех же исходных выборках для различного числа коридоров М (от М-1 до М+1), затем

сумма делится на количество просуммированных значений 21+1, тем самым определяется

среднее значение для указанного ряда сплайнов.

Следующим этапом решения задачи сглаживания является идентификация законов распределения [2]. В том случае, когда закон распределения в виде функций определенного вида с неизвестными параметрами подобрать не удается, единственным способом является применение ортогональных разложений в том или ином базисе. Отметим, что плотность распределения вероятности может находиться в следующих диапазонах (- да, да), (- да, а], \Ь, да) (см. рис. 2):

-4

-3

2

x

Р и с. 2. Типовые законы распределения: а- нормальный закон; б- экспоненциальный закон

Таким образом, при построении модели плотности распределения вероятности необходимы ортогональные функции или полиномы, ортогональные как на конечном, так и на бесконечном (полубесконечном) интервалах: Лежандра, Чебышева, Лагерра, Эрмита [4].

Рассмотрим возможность аппроксимации плотности вероятности произвольного вида ортогональными полиномами. Дана функция вида f (х), определенная на интервале [L,R]. Она может быть разложена в абсолютно сходящийся ряд вида [4]:

А ¥

f (х)=Z P Ук (хX (10)

к =0

где Рк - коэффициенты Фурье, ук (х) - семейство базисных функций, ортогональных на определенном интервале с весом т(х). Это семейство характеризуется интегралом:

R

У к (х \Vn (х МХ УХ

L

где [L,R] - интервал, на котором выполняется условие ортогональности указанных функций.

Следует подчеркнуть, что на практике приходится ограничиваться конечным числом членов ряда (10). Это приводит к появлению методической погрешности, значение которой зависит в немалой степени от способа оценки параметров модели.

Поэтому для модели аппроксимирующей функции

а m

f (х )=! рк •Ук (х X (12)

к = 0

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

Ik II2, k=n,

(її)

0, k Ф n,

б

а

Ч

L

определяются по формуле [3]

f (Х)-Ё ßk Ук (x)

m(x)dx = min, (13)

л

Рк =| /(х)У {х)/л{х)(1х . (14)

ь

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

л т

А = | /2 (х)м(х)х - £ Рк . (15)

Ь к=0

Так как выбранная система ортогональных полиномов определена на конкретном интервале, а интервал существования аппроксимируемой функции \хтп, хтах ] произволен, то они могут не совпадать. Поэтому приходится вводить коэффициенты линейного переноса, и формула (14) примет вид

, хтах

Рк = - |/(хук(х - Ь)/ а)( - Ь)/ аК. (16)

а

X,

Примеры аппроксимации различными системами полиномов и значения коэффициентов линейного переноса а и Ь представлены в табл. 1.

При численном нахождении значения интеграла (16) в силу того, что на границах интервала аппроксимации трудно добиться хорошего приближения, имеет смысл предварительно преобразовать аппроксимируемую функцию так, чтобы ее значение на одном из концов интервала (или сразу на обоих) стало равно нулю. В таком случае удается точнее посчитать коэффициенты и избежать на концах аппроксимируемой функции явления Гиббса [5].

Преобразование осуществляется по формуле

/(0)(х)= /(х)-й • х - с. (17)

Значения коэффициентов й и с рассчитываются в зависимости от того, сколько концов функции

требуется свести к нулю. Если к нулю сводится левый конец, то коэффициенты рассчитывают-

ся по следующим формулам:

\й = 0,

;( У (18)

С = / (хтт У

й = 0,

/( У (19)

С = / (Хтах У

Если сводится только правый конец, то

Для обоих концов сразу коэффициенты принимают значения:

Jd = (f (xmax )- f (xmin ))/(xmax - Xmin ) (20)

|c = f (x ■ )-d • X ■ . ( )

J vvmm ^ min

Далее вычисляются коэффициенты разложения для функции f ^о)(х):

1 Xmax

ßk0 = - J f (0)(x)^k ( -b)/a)m(x - b)/a)dx. (21)

а

xm

Полученные коэффициенты необходимо пересчитать в Рк , которые используются в выражении (12). Массивы коэффициентов отличаются друг от друга только первыми двумя элементами. Примеры расчета коэффициентов для различных ортогональных систем представлены в табл. 2.

Следующим приемом аппроксимации плотности распределения вероятностей, позволяющим уменьшить погрешность, является двухсторонняя аппроксимация. Суть этого метода заключается в том, что какая-то определенная точка х (например, точка экстремума) разбивает интервал аппроксимации функции \хШп,хтах] на два независимых интервала \хтп,х] и [х,хтах]. На каждом из них функцию аппроксимируют отдельно, что в ряде случаев позволяет значи-

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

/ ч Г0, х < 0,

1(хН (22)

[1, х > 0.

Т а б л и ц а 1

Примеры расчета коэффициентов

Полиномы Основные соотношения Коэффициенты

Лежандра [-1;1] mW=1; kdiv2 Pk (x)= E as ■ Xk~2s , s=0 ( ) V 2 • k+1 (2k—2s)/ s 42 2k • s/{k—s)/(k—2s)/ ßk _ - f f(x)Pk ((x - b)/a)dx a J xmin x — x ■ a max min ~ 2 ’ x + x ■ b max min 2

Чебышева [-1;1] m(x)= , ; V1 — x To (x ) = Vp kdiv2 ; Tk (x)= £ as • xk—2s s=0 a, =(— 1)' J2• 2k-’s—' •k(k(—s — V p s/-(k — 2s)/ ,xmax ßk _“ f Ax)Tk{{x-b)/a)/l—((x-b)/a)2 dx a J xmin x — x ■ a max min ~ 2 ’ x + x ■ b max min 2

Лагерра [0; ¥) ^(x) = e—x L0 (x)= 1 k Lk (x)= E as • xs ’ s=0 as =( 1)s • 2 k/ (s/)2 •(k — s)/ x max , ч ßk _ jf-)(-—b)/a)^--:—b)/ad, x min a _ 1, b _ xmin

Эрмита (— ¥ ;+да) k div 2 Hk (x)= E as • xk—2s , s=0 ( ,)s 1 k/ 2k—2s ^ =( 14 2k -P s/( — 2s)/ x ßk _ j fix)-x—b)/a)e--(---’da2dx, x min a _ 1, b _ i^max + xmin ~)l2

Выражение (12) с учетом (22) представим в виде

а тл

/ (х) = Е Р ,Л • У к Л ((х - ЬЛ )/ а Л ) •1(х - х)+

к=0 (23)

тП

+ Е Р П • У к П ( - ЬП )/ а П ) •1(х - Х )•

\=0

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

1 X

РкЛ = I-------------------------------------------------------1 [ /(х)ккЛ ( - ЬЛ )/ аЛ )Ил ( - ЬЛ )/ аЛ К ; (24)

\аЛ\ ^

I ^ I Хт;п

\аП\ х

для левой и правой ветви соответственно.

1 хтах

РкП = I---------1 I /(х )кП ( - ЬП )/ ап )мп ((х - ЬП )/ ап )йх (25)

Т а б л и ц а 2 Пересчет коэффициентов при сведении к нулю

Полиномы Пересчет коэффициентов

Лежандра Р0 = р00 ) + (с + й • Ъ\^2, Рг = #+ й • а •^ Рк = Р(0 \ к = 2, ..., т

Чебышева Ро = р0 )+(с+й • Ъ\4Р, Р, = )+ йа •^ Рк = Рк\к = 2,...,т

Лагерра р0 = ь0о)+ с + й •(Ъ + а), Р, = р^-й • с, Рк = Р(\к = 2,...,т

Эрмита Р0 = Р0о)+(с + й • Ъ)- 0ж, Р, = Р,)+ й • а •^ Рк = Р(\к = 2,...,т

Коэффициенты линейного переноса также изменяются соответственно изменению интервалов интегрирования.

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

Для этого требуется выполнение условий

А тл

/ (Х) = Е рк л • У к л ( - Ьл Х а Л Х

к=0 (26)

А тП

/ (х ) = Е Рк ,П ■ Ук ,П ((Х - ЪП )/а П Х

к=0

Так как при конечных значениях тл и тП условие (26) не выполняется, для его обеспечения необходимо искать аналитическое выражение для / (х) в виде

А тл тП

/(х) = Е Ьк,Л • УкЛ ( - ЪЛ /аЛ ) •- х) + Е Ьк,П • УкП ( - ЪП /аП ) •1( - х). (27)

к=0 к=0

Для того, чтобы квадратическая погрешность аппроксимации была минимальной при до-

полнительных условиях (26), необходимо минимизировать Ал и АП по ЪпЛ и ЪпП соответственно:

А л =|

А п = 1

"‘Л

/ (х) - Е Ьк ,Л • У к Л ((х - ЬЛ ) / а Л )

к=0 _

тП

/ (х) - Е Ьк П • У к П ((х - Ьп )/ а Л )

к =0

ЭЛ Л ЭЛ

Л

^Х + 1Л Е Ьк,Л ■ Ук,Л ((Х - ЬЛ )/ аЛ ) ;

к=0

)

тЛ

^Х + 1П Е Ьк,П ■ Ук,П ((Х - ЬП )/ аП ) •

Найдем частные производные ЭЛ

Л и1лП

ЭЬп,Л ’ ЭЬп,П

и приравняем их к нулю:

дЬ.

ЭЬ.

= -2

= -2

/ (х) - Е Ьк ,Л • У к, Л ((х - ЬЛ )/ а Л )

У к, л ((х - Ьл )/а л К +

ПлУк,Л ((х ЬЛ )/аЛ )= 0;

П

/ (Х) - Е Ьк, П • У к, П (х - ЬП )/ а Л )

(29)

У к, п (х - Ьп )/а л К +

+ ПП Ук,П (х ЬП )/ аП /

Выполнив промежуточные преобразования, получим

Ьп,Л = Рп,Л 2 Уп,л (х ЬЛ )/аЛ /

Ьп,П = Рп,П П Уп,П (Х - ЬП )аП /

(30)

С учетом (26) определим

(31)

"V/

Е Рк л • у к л( - Ъ л )/а л) - /(х)

к=0

тл

Е Ук,л (х - Ъл )/ал )

к=0

тП '

Е Рк ,П ' у к ,П (х - Ъ П )/а П ) - /(х)

к=0

тП

Е ук,П (х - ЪП )/аП )

к =0

Подставив это выражение в выражение для оценки коэффициентов разложения (30), окончательно получим

Л

/ (х) - Е Р, Л • Ук, Л ((х - ЬЛ V а Л )

Ьп,Л = Рп,Л +-----------------тт----------------Уп,Л ((х - ЬЛ )/а Л )’

Л

ЕУк2, Л ((х - ЬЛ )/а Л )

к=0

(32)

/ (х )- Е Рк, П ■ У к, П ((х - Ь П )1а П )

Ьп,П = Рп,П + т Уп,П ((х - ЬП )/а П /

ЕУк2, П ((х - ЬП )/а П )

к=0

В табл, 3 приведены значения коэффициентов ЬпЛ и ЬпП для различных систем ортогональных полиномов.

Так как полиномы являются степенными функциями, то от представления (12) можно перейти к степенному представлению аппроксимирующей функции, пересчитав коэффициенты. В таком случае общий вид аппроксимирующей функции выглядит так:

/(Х) = Е Рк,/ • (Х - Ь )/ а)к > Х Є [Хтт , Хтах ] •

(33)

к=0

В табл. 4 приведены формулы расчета коэффициентов Рк / в формуле (32) для различных систем полиномов. Применяемые в табл. 4 обозначения а5 х+2т - это коэффициенты ах полинома степени s+2m, используемые в табл. 1.

Т а б л и ц а 3

Коэффициенты разложения при двусторонней аппроксимации

Полиномы Коэффициенты

Лежандра тл /(х )-Е Рк, л • рк (1) Ъ„л л = Рпл л + тг рп (1/ л Е р2 (1) к = 0 тП / (х )-Е Рк, п • рк (-1) Ъп,П = Рп,П + т Рп( 1 Е Р2 (-1) к = 0

Чебышева / (х)- Ъ Р 1 (1 о [2 тл Р 1 1— Р0,л + 1 Е Рк,л ^ ^(1)

Ъп,л Рп, л 1 / (х)- а а I 1 + 2 • тл п(/ п 1 11 тп ^ По,п + , -Е(-1)к • Рк,п

-п,п Гпм • 1 + 2• тп *пУ 4 п

Лагерра тл / (х ) - Е Рк, л Ъп, л = Рп, л + , , тл + 1 тП / (х ) - Е Рк, П Ъ = Р + к=0 ип,П Рп,П ^ Л тП +1

Таблица4 Коэффициенты степенного ряда

Полиномы Коэффициенты

Лежандра, Чебышева, Эрмита (т-к )&у2 Рк,/ = Е^ •а^-, к=0,..., т 5=0

Лагерра п Рк,/ = ЕР ■ ак,*, к = 0,...,т Б=к

Заметим, что по полученному аналитическому выражению плотности распределения вероятностей легко находится аналитическое выражение и функции распределения:

0, х < хшт,

I/((Х - Ъ)/ а^ Хшт < х < Хш

Хт;п

1, X > Хшах .

(34)

Коэффициенты восстановленной функции распределения ¡3^р определяются следующим образом:

РгЛ = Pi-\J ' а / І = 1 И + 1

я+1 (35)

Ьо,— =“Х А- ■((( - ь)/a)i •

i=l

Соответственно выражение для функции распределения с учетом (35) запишется в виде

А П + 1

- (х)=-(х - ьУа У • (36)

i=0

На рис. 3 приведен пример реализации случайного процесса - оборотов газотурбинного

двигателя при пуске. В качестве примера для выборки, приведенной на рис. 3, построена ап-

проксимационная аналитическая модель плотности вероятностей, график которой представлен на рис. 4.

50 100 150 200 250 300 350 400 450 500 550 600

Отсчеты

Р и с. 3. Обороты газотурбинного двигателя при пуске

Р и с. 4. Аппроксимированная плотность вероятностей

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

0

X

БИБЛИОГРАФИЧЕСКИЙ СПИСОК

1. ВентцельЕ.С. Теория вероятностей. М.: Наука, 1969. 576 с.

2. Прохоров С.А. Аппроксимативный анализ случайных процессов / 2-е изд., перераб. и доп. Самара: Самар. гос. аэрокосм. ун-т. 2001. 380 с.

3. Прохоров С. А. Моделирование и анализ случайных процессов. Лабораторный практикум / 2-е изд., перераб. и доп. Самара: СНЦ РАН, 2001. 380 с.

4. Д. Джексон. Ряды Фурье и ортогональные полиномы. М.: ГИИЛ, 1948. 260 с.

5. АндреАнго. Математика для электро- и радиоинженеров. М.: Наука, 1965. 780 с.

Поступила 19.05.2004 г.