Вычисление некоторых характеристик неразорения страховой компании с помощью модели Лундберга Текст научной статьи по специальности «Математика»

Сер. 10. 2011. Вып. 3

ВЕСТНИК САНКТ-ПЕТЕРБУРГСКОГО УНИВЕРСИТЕТА

УДК 519.95 В. Н. Иголкин

ВЫЧИСЛЕНИЕ НЕКОТОРЫХ ХАРАКТЕРИСТИК НЕРАЗОРЕНИЯ СТРАХОВОЙ КОМПАНИИ С ПОМОЩЬЮ МОДЕЛИ ЛУНДБЕРГА*)

1. Вычисление вероятности неразорения. Рассмотрим модель Лундберга [1] разорения страховой компании с рекуррентным потоком исков, в которой интервалы между исками ти имеют распределение д(Ь), а сами иски Хи - распределение /(х). Обозначим моменты прихода исков Ьи, ти = Ьи — Ьи-1, Ьо = 0. Тогда вероятности неразорения в соседние моменты времени Ьт и Ьт-1, Рт, Рт-1 связаны так [2, с. 33]:

СЮ

Рт (и) = ^ Рт-1(у + и)^фст—X (У): (1)

-и

_ , , I 0, если и < 0,

Р0(и) = 1 1 ^ П

I 1, если и > 0.

Здесь Фст-х (у) - функция распределения величины г/и = сти — Хи, с - интенсивность получения страховых премий. Тогда для Р(и) = Иш„^то Рп(и) получим уравнение

СЮ

Р(и) = I Р(у + и)йФст-х(у). (2)

-и

Для существования нетривиальных Р(и) должно быть выполнено неравенство

сМ[ти] — М [Хи] > 0, (3)

М[...] означает математическое ожидание. Если распределение интервалов является эрланговским порядка п <?(£) = е~л*, то в изображениях по Лапласу Р(и) —> ф(р),

/(х) ^ ф(р) получим решение уравнения (2) [2, с. 34]

£ (р -А)^£*(фш(р))\р=х

*<»> ---------------------------------------• (4)

где ф(А), ф'(А),..., фп(А) - неизвестные константы, подлежащие определению.

Иголкин Владимир Николаевич — доцент кафедры математической теории моделирования систем управления факультета прикладной математики-процессов управления Санкт-Петербургского государственного университета. Количество опубликованных работ: 45. Научное направление: задачи оптимизации при случайных воздействиях. E-mail: Vigolkin@pobox.spbu.ru.

+ ) Работа выполнена при финансовой поддержке Российского фонда фундаментальных исследований (грант № 09-01-00360).

© В. Н. Иголкин, 2011

В работе [2] доказано, что при выполнении неравенства (3) знаменатель ф(р) имеет в правой полуплоскости n корней.

Тогда, в силу аналитичности, ф(р) в Rep > 0 и числитель ф(р) обращается в нуль при значениях р, равных этим корням. Это позволяет определить ф(р) с точностью до постоянного множителя, который можно найти из известного соотношения lim^o рф(р) = P (ж) = 1.

Если и распределение исков является эрланговским f(x) = то в [2] до-

казано, что все корни знаменателя простые и

-I r+1

фь>)=1- + Е—■

Р k=1 Р - Рк

где pi,... ,pr+i - корни знаменателя ф(р) из левой полуплоскости.

Таким образом, если задано распределение f (x), нужно обратить преобразование Лапласа из формулы (4) и получить P(и).

2. Вычисление последовательности Pn(и), n = 1, 2,.... Знание самой последовательности Pn(u), а не только ее предельного значения P(и), дает значительно большую информацию о финансовом состоянии страховой компании. Для вычисления последовательности Pn(u) можно использовать производящую функцию P(и, z) =

£~о Pk (u)zk.

2.1. Вычисление производящей функции. Нетрудно видеть, что функция P(и, z) удовлетворяет следующему интегральному уравнению:

СЮ

P(u,z) = 1 + z J P(и + y,z)d§cr-X(y).

Если использовать преобразование Лапласа и положить P(и, z) ^ ф(р, z), f (x) ^ Ф(р),

то

j + * ^-а)~+" Ё(р-^)км^г(ф(р^)ф(р))1р=х

Ф(р’z) = ----------------"7 Art+1(-l)rt ,, ,---------------•

1 + Z (p-aV+i Ф(Р)

Рассмотрим простейший случай, когда оба распределения являются экспоненциальными с параметрами А и л и А > ц. Тогда

^ = = р^р ~ х^р + ^ + zXfJ'^

L(p, z) = (р + fj,)\p - А + zpA], А =

А + л

Корни знаменателя Q(p, z) есть

Ро = 0, pi.2 = 0.5

A“M±V1_ (ГТ^(л+м)

Числитель Ь(р, г) должен равняться нулю при р = р\. Из этого условия найдем константу А, А = Л2~^ . Тогда

Ф(р,г)= >Л(!> + ':) = -(^ + *»

р(р - Р2 )Р1 Р1 V Р Р - Р2 J

_ М _ 1 , м

сі —---------, с2 — 1 И----------,

Р2 Р2

Следовательно,

ф{р, г)

Р(и, г)

1 Л (Л + м)(1 - \/1 - (л+^2)

1 — г \ р

1 -

2 ц(р-р2)

1г

А + /х

2/л

11

4А^г

(аТТ^

(( \ — и, А + и 4Аиг

X ехр ----------------------л / 1 — —--;---777 I М

(А + и)2

Обозначим а = Ь = (х+іі)2 • Отсюда

Р(и, г)

1г

1 — ("1 — \/1 — 62^ 1 ехр(—(2\/1 — Ъги)

и

В работе [2, с. 56] получено следующее представление Р(и, г):

Р(и, г) =

1г

1-----е

и

53 гтЬт( — 1)т-1У ,х

т—3 (пи)К N

N=0 ' э = 0

Его можно переписать так:

р (и,г) = 1 + '^2

т=1

1 _ _е-м« и

І=1

3=1

І—3 (пи)я N

ЕЦ^гЕ(-1)8ВД:5

N=0

Таким образом,

т І І—3 ( ^ N

Рт(и) = 1 - -е-^ЬЧ-ІЕ 1^гЕ(-1)8С'^С'о7-

и І=1 3 = 1 N=0 ■ я = 0

(5)

Вычисление Рт(и) с помощью рекуррентного соотношения. Соотношение (1) при с = 1 можно переписать в виде

СЮ

Рт (и) ^ Рт—1(х')^фсг—X (х и')

О

1

1

1

т

г

X

сю

Ф'(х — и) = J д(х + 4 — и)/(£)Л.

х(х —и,0)

В случае, если распределения д(Ь) и /(х) являются эрланговскими или их смесями, то интегралы в (6) берутся и Рт(и) вычисляются в аналитической форме последовательно, начиная с Р\(и).

В дипломной работе студентки факультета прикладной математики-процессов управления (ПМ-ПУ) СПбГУ Е. В. Блиновой «Исследование последовательности вероятностей неразорения в модели Лундберга-Крамера» (2010 г.) представлены программы в среде Мар1е, вычисляющие последовательности Рт(и) по формулам (5) и (6) для экспоненциальных распределений.

Пример.

д(Ь) = Ав—М, / (х) = /лв—мх,

с = 1, А = 0.5, ц =1, и = 4.

В этом случае вероятность неразорения на бесконечном интервале Р(и) равна 1 —

ле(м-А)« = о 9з2 м

Значения Рт(и), сосчитанные по формулам (5) и (6), получаются достаточно близкими: Р1 = 0.994, Р5 = 0.960, Р10 = 0.942, Р25 = 0.933, Р30 = 0.93 26, Р55 = 0.93 23. В этом примере М[ст — X] = 1 и последовательность Рт(и) достаточно быстро сходится к Р(и) .

Если М[ст—X] > 0 и небольшое, то последовательность Рт (и) сходится значительно медленнее.

При А = 0.8, л =1, с =1, и = 4 величина Р(4) = 0.641. А члены последовательности Рт(4) равны: Р1 = 0.992, Р10 = 0.829, Р35 = 0.698, Р150 = 0.644, Рз25 = 0.641.

В случае, когда функции д(£), /(х) общего вида, Рп(и) можно вычислять в узлах ит = шк, используя формулу (6) и какую-нибудь формулу численного интегрирования.

Пусть имеется таблица Рп—1(ии), к = 1,...,М. Величина Рп—1 (им) ~ 1, если М достаточно велико. Тогда

им сю

Рп(ит) ~ ^ Рп—1(ФФСТ—X (х ит ) + ^ ^Фст—X (х ит)

0 им

им

^ Рп — 1 (х)^фст—X (х ит) + 1 Фст—X (иМ ит) .

0

Интеграл может быть вычислен, например, по формуле трапеций.

3. Вероятность неразорения до момента £. Будем использовать ту же модель эволюции капитала страховой компании, д(4), /(х) функции плотности интервалов между исками и величинами самих исков. Ранее мы рассматривали величины Рп(и) вероятности неразорения на шаге п при начальном капитале и и их предельную величину Р(и) .

Если же нас интересует Р((и) - вероятность неразорения до некоторого момента времени 4, введем в рассмотрение величины Рп^(и) - вероятности того, что до момента

4 поступило п исков и не наступило разорение. Тогда Рг(и) = £ю=0 Рп,г(и).

В работе [2, с. 32] для Рп,Ь (и) получено следующее рекуррентное соотношение:

Ь п+сЬ

Рп,ь(и) = ! / Рп-1,ь-т(и + ст — у)д(т)/(у)(1у&г, (7)

0 0

Ро,ь(и) = 1 — О(Ь).

Тогда РЬ(и) удовлетворяет интегральному уравнению

Ь п+с(Ь-т)

Рь(и)= Ро,ь(и)+ ! д(т)dт J Р— (и + с( — т) — у)/(у)<1у. (8)

00

Если распределения д(Ь) и /(х) эрланговские, то интегралы в (7) являются берущимися и Рп,Ь(и) вычисляются в замкнутой форме.

В дипломной работе студентки факультета ПМ-ПУ СПбГУ А. С. Ивановой «Исследование вероятности неразорения страховой компании до момента Ь в модели Лундберга-Крамера» (2010 г.) представлены программы в среде Маріє для расчета РЬ.

3.1. Вычисление РЬ(и) в случае экспоненциальных распределений. В этом случае сначала находятся Рп,Ь(и) с помощью формулы (7). Интегралы вычисляются аналитически. При этом быстро растет количество слагаемых в выражении для Рп,Ь(и). С ростом п величины Рп,Ь(и) быстро убывают, и, когда Рп,Ь(и) станет достаточно малым, итерационный процесс заканчивается. Имеем РЬ(и) = £п=0 Рк,Ь(и)+£°=п Рк,Ь(и), «хвост» £°=п Рк,г(и) можно оценить следующим образом: Рк,Ь(и) = РЬ(и\к)пк, где

Рі(и\к) - условная вероятность неразорения, если пришло к исков; тгк = Ц^-е-Л4.

Если ХЬ < 1,то Пк монотонно убывают; если ХЬ > 1,то Пк растут до некоторого пп-1, затем монотонно убывают. Пусть такое п найдено. Тогда для отброшенного «хвоста» получим

ОО ОО

0 <^3 Рк,ь(и) Пк ,

к=п к=п

что позволяет оценить погрешность вычисления РЬ(и) = £п=о Рк,Ь(и).

Примеры.

1. Пусть Х = 0.5, л = 1, Ь = 1, и = 0.375. Тогда Р1(0.375) « £к=0 Рк,1(0.375) =

0.798,

5

1 — 13 Пк = 0.14 ■ 10-4, Р10(0.375) = 0.5892.

к=0

Для этого примера нетрудно найти решение уравнения (2): Р(0.375) = 0.5855. Ясно, что при достаточно больших Ь РЬ(и) близко Р(и).

2. Пусть Х = 0.5, л =1, Ь = 10, и = 1. Тогда Р10(1) = 0.6938, а Р(1) = 0.6968.

3. Пусть Х = 0.5, л =1, Ь = 20, и = 1.25. Тогда Р20(1.25) = 0.7368, а Р(1.25) = 0.7324.

3.2. Вычисление РЬ(и) в случае произвольных распределений. В этом случае также можно использовать формулу (7), вычисляя интегралы с помощью формулы трапеций. Величины Рп,Ьк (и^) строятся последовательно в узлах Ьк = кН, и^ = іі, к = 0,1,... ,М, Ь = МН, I = сН, і = 0,1,..., М, М1 = и + сЬ. Шаг Н подбирается в процессе вычислений.

С помощью данного алгоритма была вычислена величина Р1(0.375) для примера 1. Получен примерно тот же результат.

Для вычисления Рг(п) можно использовать и формулу (8). Перепишем ее следующим образом:

г и+с(г—т)

Рг(п) = Ро,г(п) + I д( - т)г1т ^ Рт(у)1 (п + с( - т) - у)^у. (9)

о о

Вероятность Рг(п) можно вычислять последовательно в узлах = кН, п3 = jl. При расчете интегралов используется формула трапеций. В процессе вычислений подбирается шаг Н. Тогда из формулы (9) получим расчетную формулу

кН ]1+с{кН—т)

Ркнф) = 1 - О(кН) + I д(кН - т) ^ Рт (у)/(I + с(кН - т) - у)<1у<1т =

оо

з1+с(кН—т)

о(

Н Г Н

= 1 - С(Щ + -д(кк) I Ро(у)/(и + к)1-у)г1у+-д(0)х

о

3 к-1 (3 + к—г)1

РкН(у)/и1 - у^у + 2^2 д((к - 1)Н) J Рн(у)!(и + к - г)1 - у^,у

Н

1 — С(кК) + —

д(Щ^1(Р0(0)/(и + к)1) + Р0(и + к)1)т +

3 + к — 1

+ ^53 Р0(в1)/((и + к - з)1) +

в=1

3 — 1

+ д(0)0.51(Ркн(0)/(и1)+Ркн(и1)1 (0) + 2 ]ТРкнЫ)/ ((и - з)1)) +

в = 1

к—1

+ 2 53 д((к - ъ)Н)0-51(РгН(0)/((и + к - ъ)1) + Рн ((и + к - Ъ)1) +

1=1

3 + к—г+1

^гН

+ 2 Ргн(з^)/((и + к - ъ - 8)1))

в=1

(10)

Сначала заполняется строка Р0(3) = 1 - 0(0). Строка РН(и) заполняется последовательно по I = 0,1,..., М. Затем заполняется строка Р2н(и) и т. д. Заметим, что величины Ркн(Ц) находятся и в правой части формулы (10).

Для д(Ь) = Хв—хг, /(х) = /1,в—^х X = 0.5, л =1, с = 1 был сосчитан ряд примеров, используя расчетную формулу (10). Приведем некоторые результаты: Р1.25(0.125) = 0.723, Р125(0.375) = 0.772, Р125(0.75) = 0.830, Р125(1.5) = 0.906, Р1.25(5) = 0.995, Ро.125(0.75) = 0.972, Ро.з75(0.75) = 0.927, Рз.125(0.75) = 0.745.

В соответствии с содержанием задачи Рг(п) монотонно растет по п и убывает по Ь.

4. Случай ступенчатых распределений. Распределения д(Ь) и /(х), как правило, заранее не известны. Обычно известны конечные выборки из соответствующих генеральных совокупностей. Для аппроксимации неизвестной функции плотности часто пользуются выборочной плотностью или гистограммой. Если нас интересует вероятность неразорения страховой компании Р(и) = Ишп^Ю) Рп(и), то нужно решать интегральное уравнение (2)

СЮ

Р(и) = J Р(у + и)йФст-х(у),

—и

СЮ

фст-х(у) = I с <1Р{х).

— Ю

Если в качестве распределений имеем гистограммы, то Фст—х (у) будет ступенчатой функцией с разрывами первого рода в узлах гистограммы ст — X. Предположим, что таких точек у к, конечное число к = 1, 2,...,Ы, и Чк - величина скачка в точке у и, Чк = Фст—Х(ук+1) — Фст—Х(ук). Тогда интегральное уравнение становится разностным:

Ю N

Р(и) = I Р(у + и^Фст—х(у) = 53 ЧзР(и + уз). (И)

— и 3=1

Уравнение (11) является разностным однородным уравнением порядка N с постоянными коэффициентами. Для нахождения Р(и) составим характеристическое уравнение. Будем искать Р(и) в виде Р(и) = вХи. Из (11) имеем

N

вХи =53 ЧзеХ(и+Уз}

3=1

откуда получаем характеристическое уравнение для Л

N

1 = Т, ЧзеХУз. (и)

з=1

При построении гистограмм величины уз следует выбирать соизмеримыми. Тогда при решении уравнения (12) нужно будет искать корни полинома. Представляют интерес лишь значения Л из левой полуплоскости, так как Р(и) ограничена на бесконечности. Заметим, что среди корней характеристического уравнения (12) есть корень Л1 = 0, так как £3=1 3 = 1. Если Л1, Л2,..., Лт - корни уравнения (12) из левой полуплоскости, то решение Р(и) будем искать в виде Р(и) = £^=1 сиеХки, где си - неизвестные коэффициенты, со = 1, так как Р(и) ^ 1 при и ^ то.

Для нахождения си, к = 1, 2,...,т, нужно определить т значений Р(и). Это можно сделать с помощью рекуррентной формулы (1), которая в рассматриваемом случае будет иметь вид Рт(и) = Х3=1 ЧзРт—1(и + уз). Если гистограмма функции Фст—Х(у) построена с шагом Н = уз — уз —1, то будем искать Рт(кН) последовательно по т, начиная с Р1(кН) = 1 — Ф(—кН). Значения Рт(кН) монотонно не возрастают по т, и, когда

Pm—l(kh) и Pm(kh) достаточно близки, итерации заканчиваются и Pm(kh) « P(kh). Неизвестные коэффициенты Ck находятся из системы P(kh) = £^=1 cseЛskh, к =

1, 2,...,m.

В дипломной работе студентки факультета ПМ-ПУ СПбГУ О. Т. Батыровой «Вероятность неразорения страховой компании в случае ступенчатых распределений» (2010 г.) разработана программа в среде Maple 13, которая реализует данный алгоритм, и приведены численные примеры.

Примеры.

1. g(t) = Ae—Лі^, f (x) = fie—^x, A = 0.5, ц = c = 1.

Тогда точным решением уравнения (2) является

P(u) = 1----------e( = -**)“ = 1 - 0.5e-°'5

c^

Затем были получены с помощью датчика псевдослучайных чисел выборки из распределений д(Ь) и /(х), объем каждой выборки равен 100. Были построены гистограммы с шагом 1 и получено характеристическое уравнение

0Ш11е—3Л + 0.0678е—2Л + 0.1284е—Л — 0.6668 + 0.1986еЛ +

+ 0.1036е2Л + 0.0758е3Л + 0.0334е4Л + 0.0191е5Л + 0.0188е6Л = 0.

Его единственным корнем из левой полуплоскости является Л2 = —0.4886, Л1 =0 также является корнем. Тогда Р(и) = 1 — се—°-4886и.

Было найдено Р(0) « 0.5555. Тогда с = 0.4445 и Р(и) = 1 — 0.4445е—°-4886и. Полу-

ченное приближенное решение близко к точному.

2. д(Ь) = Л(ХЬ)е—Л1, /(х) = ц(рх)е—^х, Л = 1, ц = 2, с = 1, и = 0.1. Используя

формулу (4), можно найти точное решение уравнения (2) и Р(0.1) = 0.5905.

Далее для распределений д(Ь) и /(х) были получены выборки объемом 100, построены гистограммы с шагом 1, выведено характеристическое уравнение и определены корни из левой полуплоскости: Л2 = —0.8181, Л3 4 = —1.3590 ± 2.2090г. Тогда Р (и) = 1 + с1еЛ2и + с2 еЛзи + с3еЛ4и.

Для определения с1, с2, сз найдены приближенные значения Р(0) = 0.6418, Р(1) = 0.7668, Р(2) = 0.8546. Тогда

P(u) = 1 - 0.4707e—Q + 2Re

(0.0053 + 0.0004i)e(—l'359Q+2'2Q9Qi)u

Р(0.1) = 0.5752. Найденное приближенное значение Р(0.1) близко к точному.

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

Литература

1. Grundell G. Aspect of risk theory. New York: Springer-Verlag, 1992. 17G p.

2. Иголкин В. Н., Ковригин А. Б. Финансовые потоки и их флуктуации. СПб.: Изд-во С.-Петерб. ун-та, 2GG6. 134 с.

Статья рекомендована к печати проф. Л. А. Петросяном.

Статья принята к печати 1G марта 2G11 г.