№>1(9)2008
С.А.Айвазян
Баиесовскии подход в эконометрическом анализе
Предлагаемая в данном номере журнала консультация посвящена так называемому байесовскому подходу в эконометрическом анализе, основанному на субъективно-вероятностном способе операционализации принципа максимального использования (наряду с исходными статистическими данными) априорной информации об исследуемом процессе.
Байесовские методы широко распространены в теории и практике эконометриче-ского анализа и являются обязательной составной частью современных учебных программ магистерского уровня по эконометрике в ведущих университетах мира. Особенно заметные преимущества (по сравнению с классическими методами) с точки зрения точности получаемых статистических выводов они имеют в условиях относительно малых выборок, что весьма характерно для эконометрического моделирования.
1. «Философия» и общая логическая схема байесовского подхода
Пусть в описании рассматриваемой эконометрической модели (закона распределения анализируемой случайной величины, функции регрессии, временного ряда, системы одновременных уравнений и т. п.) участвует 5-мерный параметр 0 = (0,,92,...,05)т и нашей задачей является построение наилучшей, в определенном смысле, статистической оценки 0 этого параметра по имеющимся к-мерным наблюдениям X, - (х,(1),х,(2),...,х,(к) )т, / - 1,2,...,п. Верхний индекс Т здесь и в дальнейшем означает операцию транспонирования вектора или матрицы, прописными буквами будут обозначаться векторные величины (записываемые как векторы-столбцы), а строчные буквы будут использоваться для обозначения одномерных (возможных или наблюденных) значений анализируемых случайных величин.
Байесовский подход является одним из возможных способов формализации и операционализации тезиса, в справедливости которого нет видимых причин сомневаться: степень нашей разумной уверенности в некотором утверждении (касающемся, например, неизвестного численного значения интересующего нас параметра) возрастает и корректируется по мере пополнения имеющейся у нас информации относительно исследуемого явления. Могут быть различные формы интерпретации и подтверждения этого тезиса, в том числе не имеющие отношения к байесовскому подходу. Одна из них выражена, например, в свойстве состоятельности оценки 0п неизвестного параметра 0: чем больше объем выборки п, на основании которой мы строим свою оценку 0п, тем большей информацией об этом параметре мы располагаем и тем ближе (в смысле сходимости (0 п к 0 по вероятности) к истине наше заключение.
Специфика именно байесовского способа операционализации этого тезиса основана на двух положениях.
93
Не 1(9) 2008
1) Во-первых, «степень нашей разумной уверенности» в справедливости некоторого утверждения численно выражается в виде вероятности. Это означает, что вероятность в байесовском подходе выходит за рамки ее интерпретации в терминах условий статистического ансамбля (см. п. В.2.1 в [Айвазян, Мхитарян (2001а)]), но относится к одной из категорий субъективной школы теории вероятностей.
2) Во-вторых, статистик при принятии решения использует в качестве исходной информации одновременно информацию двух типов: априорную и содержащуюся в исходных статистических данных (см. п. В.3.2 в [Айвазян, Мхитарян (2001а)]). При этом априорная информация предоставлена ему в виде некоторого априорного распределения вероятностей анализируемого неизвестного параметра, которое описывает степень его уверенности в том, что этот параметр примет то или иное значение, еще до начала сбора исходных статистических данных. По мере же поступления исходных статистических данных статистик уточняет (пересчитывает) это распределение, переходя от априорного распределения к апостериорному, используя для этого известную формулу Байеса:
Р (Л |В} =
Р (Л} -Р (В| А,-}
N ,
(В| А,}-Р (А,}
1=1
(1)
которая определяет правило вычисления условной вероятности события А, (при условии, что событие В уже имело место) по безусловной вероятности события А, и условным вероятностям Р (В |А,}, 1 = 1, 2, ...,М При этом предполагается, что А,, А2,..., Ам образуют полную систему событий, а событие В имеет ненулевую вероятность (т. е. Р(В} >0).
Общая логическая схема байесовского метода оценивания значений параметров представлена на рис. 1.
<ъ §
§
¡В
I
и <ъ
8-<ъ 5
0
1
т оз
! с
! и
8 и <ъ >5 <в
Щ
Рис. 1. Общая логическая схема байесовского подхода в статистическом оценивании
Рассмотрим реализацию схемы байесовского оценивания неизвестного параметра.
Априорные сведения о параметре 0 основаны на предыстории функционирования анализируемого процесса (если таковая имеется) и на профессиональных теоретических соображениях о его сущности, специфике, особенностях и т. п. В конечном итоге эти априорные сведения должны быть представлены в виде функции р( 0), задающей априорное распределение параметра и интерпретируемой как вероятность того, что параметр примет значение, равное 0, если параметр дискретен, или как функция плотности распределения в точке 0, если параметр непрерывен по своей природе.
Исходные статистические данные Хи X 2,..., Хп порождаются в соответствии с законом распределения вероятностей ((X10), где под ((Х| 0) понимается значение функции плотно-
94
Не 1(9)2008
-V
ф Консультации
95
сти наблюдаемой случайной величины 2 = (2(1),2(2),...,2 1к))т вточкеX, если 2 — непрерыв- Ц на, или вероятность Р{2 = X10}, если 2 дискретна (при условии, что значение неизвестного § параметра равно 0). По умолчанию предполагается, что наблюдения
Х„Х2.....Хп (2) 3
при фиксированном 0 являются статистически взаимонезависимыми, т. е. образуют случайную выборку из анализируемой генеральной совокупности.Так что, получая исходные статистические данные (2), мы к имеющейся априорной информации о параметре (в виде функции р(0)) присоединяем соответствующую выборочную (эмпирическую) информацию.
Соответственно, функция правдоподобия ¿(X,,...,Хп|0) (условная, при данном 0) имеющихся наблюдений (2) определится (с учетом их условной взаимонезависимости) соотношением
¿(X „ X 2,..., Хп |0) = Д X110) • Г (X 210) •... • Г(Хп |0). (3)
Вычисление апостериорного распределения р(0|X1,...,Xn) осуществляется с помощью формулы Байеса (1) (или ее непрерывного аналога), в которой роль события А играет событие, заключающееся в том, что значение оцениваемого параметра равно 0, а роль условия В — событие, заключающееся в том, что значения п наблюдений, произведенных в анализируемой генеральной совокупности, зафиксированы на уровнях X,, X2,..., Xn .Соответственно, имеем:
р^......^) = Р(0)ЦX......^ |0) . (4)
Построение байесовских точечных и интервальных оценок основано на использовании знания апостериорного распределения р(0|X,,...,Xn), задаваемого соотношением (4). В частности, в качестве байесовских точечных оценок (0(Б) используют среднее или модальное значение этого распределения, т. е.:
0 ¡Бср, = Е(0| X......^) = ^0p(0|X ......^ )С0, (5)
0^од = агд пахр^^,,...,Xn).
Отметим, что для определения общего вида апостериорной плотности ~ (0^,..., Xn )нам достаточно знать только числитель правой части (4), так как знаменатель этого выражения играет роль нормирующего множителя и от 0 не зависит (это существенно упрощает процесс практического построения оценок 0сБр) и 0^д).
Отметим также одно важное оптимальное свойство оценки 0 сБр). Пусть 0(X,,...,Xn) — любая оценка параметра 0. Оказывается, если качество любой оценки 0 (X,,..., Xn) измерять так называемым апостериорным байесовским риском
Я(Б) (X,,..., Xn) = Е{(0 (X,,..., ^) -0)2 X,,..., ^} = |(0( X,,..., Xn) -0)2 рВД,..., ^ )С0
или его средним (усреднение — по всем возможным выборкам (2)) значением Яср), то байесовская оценка (5) является наилучшей и в том и в другом смысле.
Нв 1(9) 2008
Для построения байесовского доверительного интервала для параметра 0 необходимо вычислить по формуле (4) функцию р(0|Хь...Хп) апостериорного закона распределения
1 + Р0
параметра 0, а затем по заданной доверительной вероятности Р0 определить 100-и
1— Р 2
100—-—%-ные точки этого закона, которые и дают соответственно левый и правый концы
искомой интервальной оценки.
Заметим, что байесовский способ оценивания может давать весьма ощутимый выигрыш в точности при ограниченных объемах выборок по сравнению с традиционным «частотным» подходом. В процессе же неограниченного роста объема выборки п оба подхода будут давать, в силу их состоятельности, все более похожие результаты.
«Узкие места» или три главных вопроса, возникающие при практической реализации байесовского подхода:
О как выбрать общий вид (т. е. параметрическое семейство р(0; Э)) априорного распределения оцениваемого параметра?
П) как подобрать численные значения Э0 параметров Э, определяющие конкретный вид априорного распределения при уже сделанном выборе общего вида р(0; Э)?
¡¡¡) как преодолеваются трудности реализации формулы (4) при вычислении апостериорного распределения ~(0|Х1,...,Хп)?
2. Априорные распределения, сопряженные с наблюдаемой генеральной совокупностью (определение и условие существования)
В решении сформулированных выше трех главных вопросов практической реализации байесовского подхода существенную роль играют распределения, сопряженные с наблюдаемой генеральной совокупностью (или, что то же, — распределения, сопряженные с функцией правдоподобия ЦХ-|,... Хп 10)).
¡| Определение 1. Семейство априорных распределений 6 = (р(0;Э)} называется сопря-
| женным по отношению к наблюдаемой генеральной совокупности ((Х10) (или по отноше-
| нию к функции правдоподобия ЦХЬ...,Хп|0), если и апостериорное распределение
<5 р(01Хт,...,Хп), вычисленное по формуле (4), снова принадлежит этому же семейству
Другими словами, семейство распределений 6 сопряжено с ^(Хп,..., Хп 10), если оно замк-
Ф I л \
§ нуто относительно операции (4) пересчета априорного распределения в апостериорное.
I
т оз
Таким образом, использование в качестве априорных законов распределения вероятностей (з.р.в.) сопряженных по отношению к I плотностей «расшивает» узкое место (¡¡¡): по-§: скольку общий вид апостериорного з.р.в. в этом случае известен, остается лишь уметь переев считывать значения его параметров Э при переходе от априорного распределения к апо-с
! и
<§ орных оказывается в широком классе случаев вполне естественным и оправданным, что позволяет получить ответ и на вопрос (¡).
Но всегда ли существует сопряженное по отношению к заданной функции 1(Х 1, Х2,... Хп 10) распределение, и если оно существует, то как его найти?
стериорному.
Как мы увидим позже (см. ниже п. 3), использование сопряженных з.р.в. в качестве апри-
96
Не 1(9)2008
Условие существования сопряженного семейства априорных распределений: |
если функция правдоподобия ¿(X,,...,Xn |0) представима в форме
¿(X,,...Xn |0) = у(Г, (X,,...Xn), ...,Гт (X,,...Xn); 0) • X,,...Xn), (6)
гдеТ](X,,...Xn)(} = ,, 2,...,т) и,...Xn) — некоторые функции от наблюденийX,,...,Xп, не зависящие от параметров 0, то существует семейство6 = р(0; С) априорных распределений, сопряженное с¿(X,,...,Xn|0).
Проверка условия существования сопряженного априорного распределения на ряде примеров.
Пример 1.
Анализируемая (наблюдаемая) генеральная совокупность нормальна с неизвестным значением среднего Е2 = 9 и известной дисперсией 02 = а2 (будем обозначать в дальнейшем подобный факт в форме 2е N(9; а0), где 2 — наблюдаемая случайная величина, а нижний индекс ц определяет ее размерность;так что, если 2 = (2ш,...,2(к))т — вектор, то 2еМк(0;) означает, что многомерная случайная величина размерности к распределена нормально с вектором средних значений 0 = (9,,92,...,9к)т и ковариационной матрицей ). В данном примере
Лу (X;-9)2 —(X-9)2 ( , Л п -^-У (X ,-X)2
Пример 2.
2 е N |9,;— I, где и среднее значение 9, = Е2, и 92— параметр точности
V
,
12=02.
Л
-V
ф Консультации
97
8 <в
О
¿( х,,..., х п |9) = П (X; |9) = 1 ■— I е 2а2» = е 2а° •[ I е 2а° = . (7)
Мы видим, что роль функции V(Т(х,...хп);9) из правой части (6) играет первый сомножи-
, П
тель в правой части (7), причем т = ,, Т(х,,...,хп) = X= —Ух, (достаточная статистика),
_ п ^
а следующие за V(X;9) сомножители правой части (7) от9 не зависят. Следовательно, семейство априорных, сопряженных с I распределений существует.
неиз-
вестны (т. е. 0 = (9,,92)). Воспользовавшись тем же представлением (7) для функции правдо-
, п
подобия ¿, убеждаемся, чтоT1(x1,...,xn) = X, Т2иь...,xn) = 52 = — У (X,- -X)2 (достаточные ста-
п ^
тистики), так что в данном случае т = 2 и семейство априорных, сопряженных (по отношению к ¿) распределений существует.
1 Функции T^(X1,...,Xn), участвующие в представлении (6) (если таковое существует), называются достаточными статистиками в задаче статистического оценивания параметров0 = (9, ,92,... Д)т. Размерность т векторной достаточной статистики Т(X1,...,Xn) = (T1(X1,...,Xn),...,Tm(X,,...,X,)) конечна при п и зависит от специфики функции и размерности $ оцениваемого параметра ©.Достаточные статистики играют важную роль в теории и приложениях математической статистики. В частности, они используются в задаче построения наилучших несмещенных оценок в следующей схеме: пусть0 — некоторая несмещенная оценка параметров 0 и 1:г тогда 0, = Е(01Т) будет снова несмещенной оценкой параметров 0, причем 1:г ^ < 1:г Е0 (под Е0 понимается ковариационная матрица вектора 0).
п
п
2
Нв 1(9) 2008
Пример 3.
Анализируется биномиально распределенная случайная величина £,е (М) — число «успехов» в серии из М испытаний Бернулли, где е — неизвестная вероятность «успеха» в одном таком испытании, а М — общее число (известное) испытаний Бернулли в рассматриваемой серии, так что
((х|е) = р(м) = х|е} = смех(1—е)м—х, х = 0,1,2,...,м.
Наблюдаются п таких серий. Тогда
п Ух, пМ—у х, п
цх 1,х2,...,хп|е) = Псмех(1—е)м—х = е» (1—е) = Псм,
=1 =1
где х — число «успехов» в -й серии.
п
Поэтому в рамках общего представления (6) в данном случае имеем: т = 1, Т (х1,..., хп )=У х —
1=1
достаточная статистика, что подтверждает существование априорного, сопряженного с Ь распределения параметра е.
Пример 4.
Анализируется отрицательнобиномиально распределенная случайная величина £,(е; К) — число испытаний в схеме Бернулли до К-го появления интересующего нас события, гдее — неизвестная вероятность появления этого события при одном испытании, а К — некоторое заданное целое положительное число. Тогда
((х|е) = р ш к) = х|е} = с—леК (1—е)х—К, х = к, к+1,...,
так что
п у х,—Кп п
их ......хп |е) =ПсКт—1еК (1—е)х—К = еКп (1—е) « .^—1.
, =1 , =1 ^ п
| Поэтому в рамках общего представления (6) в данном случае имеем: т = 1, Т (х1,..., хп )=У х, —
;=1
| достаточная статистика, что подтверждает существование априорного, сопряженного с I § распределения параметра е.
и <ъ
Ц Пример 5.
| В данном примере речь идет об оценивании параметра е пуассоновского з.р.в., т. е. о
5
о е х
<5 ((х|е) = Р (^ = х|е} = — е-9, х = 0,1,2,...,
в х!
:
8 так что
Ч п
О п ех; У х, п ( 1 л
Цх 1,х......хп |е) = П~ е-9= е-пе -еУ П| —|.
§ !=! х,! !=! ^ х,!
и
03
8 Сравнивая с общим представлением (6), в данном случае имеем:т = 1, Т(х1,...,хп)=ух, —
'=1
щ достаточная статистика, что подтверждает существование априорного, сопряженного с I распределения параметра е.
98
На 1(B) 2008
Пример 6. Ц
Анализируется экспоненциально распределенная (без сдвига) случайная величина с неиз- ¡3 вестным значением параметра масштаба е, т. е.
ч
|ее~9х при х > 0 о
f( х|9) = ,п
10 при х < 0
Соответственно:
L(х......хп |9) =
j х, I 9
9ne Ul 1 при х, > 0
0 при х, < 0
п
В рамках общего представления (6) в данном случае имеем: т = 1; Т(х1,...,хп) = у х, —
,=1
достаточная статистика, что подтверждает существование априорного, сопряженного с Ь распределения параметра е.
Пример 7.
Анализируется случайная величина, распределенная равномерно на отрезке [0;е] при неизвестном значении параметра е, т.е.
((х|е) = {епри0 -х-е.
[0 при х ( [0; е]
Соответственно:
Цххп |9) = | -Ч при 9 >хmax(n) = maxх,.
1^9 ) i< ¡< п
Следовательно, в рамках общего представления (6) имеем: m = 1; T(х1,...,хп) = хтах(п)— достаточная статистика, что подтверждает существование априорного, сопряженного с L распределения параметра 9.
Пример 8.
Анализируется модель распределения Парето с неизвестным значением параметра формы 9, т. е.
[9 х 9 >
f( х|9) = ^ при х > х о, [0 при х < х0
где пороговое значение х0 считается заданным. Соответственно:
-п9
\п | gп I п
Цх1.....хп|е) = епхп9-Щх,) =еп-дп >
гдедп = |Цх,^ —среднее геометрическое значение наблюденийх1,х2,...,хп анализируемой случайной величины. Следовательно, обращаясь к (6), имеем т = 1, Т(х1...хп) = дп —
99
Не 1(9) 2008
достаточная статистика, что подтверждает существование априорного, сопряженного с I распределения параметра 9.
Пример 9.
Рассмотрим классическую линейную модель множественной регрессии (КЛММР, см., например, [Айвазян (200,), §2.2]) с нормальными, в среднем нулевыми, взаимонезависимыми и гомоскедастичными остатками е,,е2,...,еп:
У = Х0 + Б,
где
У = (У,, У 2,., У,
Х=
, X ™ , X ™
, X п
Л к) Л
Л к)
Л к)
(8)
наблюденные значения соответственно зависимой (у) и объясняющих (X = (, X(,),...,X(к) )т) переменных, е = (е,,е2,...,еп)т — случайные регрессионные остатки, а 0 = (90,9,,...,9к)т и Ь = (йе ,•)- — неизвестные значения параметров модели. Напомним, что значения X, в соответствии с требованиями КЛММР, являются неслучайными и что упомянутые выше свойства регрессионных остатков формулируются в форме условий:
ее N,10;^I,
(9)
,
<ъ §
§
¡В
I
и <ъ
! <ъ 5
0
1
т оз
! с
1 и
8 и <ъ >5 <в
Щ
гдеIп — единичная матрица размерности п, ковариационная матрица остатков2е = — Iп,апа-раметр Ь = (йе,)-1 обычно называют параметром точности. Ь
С учетом (8)-(9) функция правдоподобия наблюдений (X, У) может быть представлена в форме:
п
Ь 2
¿(Х,У|0; Ь) =
-(У-Х0)' (У-Х0)
п
(2 я )2
(Ю)
Но (У - Х0)т(У - Х0) = (У - Х0 + хе? - Х0)т(У - Х0Э + Х0Э - Х0) = [(У - Х0Э) + Х((Э -0)]т X х [(У-X(©) + Х((© -0)], где 0 = (ХтХ)ХтУ — оценка метода наименьших квадратов параметров регрессии 0. Поэтому
(У - Х0)т( У - Х0) = (У - Х0)т(У - Х0) + (0Э-0)т Хт X((©-0), (И)
так как (У- Х0)тX((© - 0) = [X((© - 0)]т(У- Х0Э) = (УтХ - 0т ХтХ)(0 - 0) = [ УтХ - ((ХтХГ1 Х тУ )т х х (X т X)]((0 -0) = [У т X - Ут X(Xт X)Xт X)]((0 -0) = 0.
Возвращаясь к (Ю) и выражая в (И) сумму квадратов МНК-оцененных остатков
(У - Х0)т(У - Х0 )через оценку остаточной дисперсии а2 =
,
п - к -!
(У- Х(()т(У-Х0), имеем:
¿(Х;У|0; Ь) =-
-• Ь2 е
^^ а 2| ь --(((-0 )т Xт Х( 0-0) 2 ) 2
(,2)
(2я):
и
е
п
100
-^ № 1(9) 2008
Отметим, чтост2 и 0, в конечном счете, определяются по YT Y, XT Y и XT X, так что ивдан- ц ном случае функция правдоподобия L представима в форме (6), в которой набор достаточ- § ных статистикT(X; Y) конечен (при п ^ да) и определяется статистиками YT Y и XT Y. Следовательно, существует априорное распределение параметров 0 и h, сопряженное с L. ^
3. Генезис априорных сопряженных распределений
Оказывается, для широкого класса наблюдаемых генеральных совокупностей, функции правдоподобия которых допускают представление (6) (т.е. эти генеральные совокупности располагают сопряженным с L априорным распределением своих параметров), справедливо следующее утверждение о генезисе сопряженных априорных распределений:
если в байесовском подходе стартовать с априорного распределения, не несущего никакой дополнительной по отношению к имеющимся статистическим данным полезной информации об оцениваемых параметрах, то первый же переход от нее по формуле (4) к апостериорному распределению приведет нас к семейству распределений, сопряженному с наблюдаемой генеральной совокупностью2.
Именно этот прием поиска априорного распределения, сопряженного с анализируемой функцией правдоподобия, представимой в формуле (6), и предлагается использовать в байесовском подходе.
3.1. Априорные распределения, отражающие «скудость априорных знаний» (САЗ-априорные распределения)
Для математической формализации ситуаций, в которых исследователь не располагает никакой полезной априорной информацией о значениях оцениваемого параметра, Джеф-фрис (см. [Jeffreys (1957)]) предложил следующие два правила выбора соответствующего априорного распределения:
(а) если оцениваемый скалярный параметр 9 может (теоретически) принимать значения на конечном интервале [9min, 9max ] или на бесконечном интервале от -да до +да, то априорную функцию плотности p(9) следует считать постоянной на соответствующем интервале;
(б) если же из смысла оцениваемого параметра вытекает, что он может принимать любые положительные значения, то следует считать постоянной на всей числовой прямой (-да; + да) функцию плотности распределения логарифма от значения параметра, т.е. p(ln 9) = const при 9 £ (0; + да).
Будем называть такие априорные распределения «распределениями, отражающими скудость априорных знаний» или коротко — «САЗ-априорными распределениями». Соответственно, их одномерные функции плотности будем обозначать pСАЗ (9), а многомерные —
p САЗ (0).
Тот факт, что для определенных таким образом на бесконечной прямой (полупрямой) априорных распределений нарушается известное правило нормировки функции плотности
2 Строгое доказательство этого утверждения для однопараметрического экспоненциального семейства наблюдаемых генеральных совокупностей см. в [Ghosh et al. (2006), п. 5.1.5]. Однако справедливость этого утверждения подтверждается (непосредственной проверкой) и для весьма широкого класса наблюдаемых генеральных совокупностей, не принадлежащих экспоненциальному семейству.
101
Не 1(9) 2008
вероятности (поскольку при этом JpСА3 (9) d9 ф 1, но JpСА3 (9) d9 = да, где интегрирование проводится по всем возможным значениям 9), не доставляет «технических неудобств»: во-первых, пересчет такой «несобственной» априорной функции плотности pСА3 (9) в апостериорную по формуле (4) дает уже обычную (собственную) функцию плотности pСА3 (9|X1 ,...,Xn), а во-вторых, при любых сколь угодно больших значениях C плотность
/а\ I— при 9е [-C; +C] Р саз (9) = Ьс
[0 при 9g [-C; +C]
с
минимизируетэнтропийную меру информации— H = Jp(9)lnp(9) d9, содержащейся в плот- с
ности p(9) относительно параметра 9 (см., например, [Зельнер (1980), с. 59]). Последнее обстоятельство подтверждает обоснованность использования равномерных распределений pСА3 (9) = const или pСА3 (ln9) = const в качестве априорных распределений, отражающих скудость априорных знаний (или САЗ-априорныхраспределений).
Замечание 1. Общий вид апостериорного распределения p(9|X1,...,Xn), вычисляемого по формуле (4), определяется, с точностью до нормирующей константы, лишь числителем правой части этой формулы. Поэтому в дальнейшем при анализе равенств, справедливых с точностью до нормирующей константы, мы будем использовать знак Следуя этому правилу, сама формула (4) может быть представлена в виде:
Р(®|Х......Xn )p p(©)-L(X 1,...,Xn |0). (4')
Замечание 2. При анализе многомерных параметров 0 = (9—,...,9S)T априорные, в том числе САЗ-априорные, распределения обычно предполагают статистическую независимость компонент9i,...9s, т.е.
p(0) = p( 9i) • p(9 2) •... • p( 9 s). (13)
ф
g И, наконец, в заключение этого пункта определим вид априорной плотности p(Q)дляслу-§ чая p(ln9) = const, т. е. в ситуации, когда параметр9может принимать любые, нотолько поло-| жительные значения.
§ Пусть F9(y) = P{9< y} — функция распределения параметра 9.Тогда
<и
i[ F9 (y) = P{9 < y} = P{ln9 < lny} = Fln9 (lny).
Ег
<ъ
Соответственно, функция плотности распределения 9 будет:
о
is ft ( y) = ?!М = 5^n 9 (ln y) •50M = ^n 9 (ln y) • Ipl,
dy 5(ln y) dy y y
так как по условию fln 9 (ln y) = p(ln 9) = const. Так что в сокращенной записи имеем для положи-
<а
Л
с л
>s тельнозначных параметров 9:
I 1
8 p саз(9)~-, (14а)
¡3 9
щ
'(I а для параметров 9 с возможными значениями, заполняющими всю числовую прямую, щ
p саз (9) = const. (14б)
102
На 1(B) 2008
3.2. Общий подход к выводу семейства априорных распределений, Ц
сопряженных с наблюдаемой генеральной совокупностью §
'S
Общий подход к выводу семейства априорных распределений, сопряженных с наблюдаемой генеральной совокупностью, основан на утверждении об их генезисе, сформулиро- ^ ванном в начале пункта 3. Из этого утверждения вытекает, в частности, следующая общая схема определения такого семейства.
Шаг 1: проверка условия (6) существования семейства априорных распределений, сопряженных с функцией правдоподобия L для наблюдаемой генеральной совокупности.
Шаг 2: если функция правдоподобия L допускает представление (6) (т. е. если существует семейство сопряженных априорных распределений p(0;Л)), то осуществляется вывод САЗ-апостериорного распределения ~САЗ (0|Х1,...,Хп) по формуле (4'), т. е.
~САЗ (01 Xi.....Хп pСАЗ (0) -L(Xi.....Хп |0). (15)
Правая часть соотношения (15) и будет определять общий вид семейства априорных распределений, сопряженных с наблюдаемой генеральной совокупностью, характеризуемой функцией правдоподобия L(Х1,X2,...,Хп 10).
Продемонстрируем реализацию этой общей схемы на рассмотренных выше примерах 1-9. Очевидно, нам остается реализовать лишь шаг 2 из этой схемы, так как шаг 1 уже был реализован выше (см. п. 2).
Пример 1 (продолжение).
2 £ N1(9;ст2), где 9 = E2, — оцениваемый (неизвестный) параметр, а ст0 = D2, — известное (заданное) значениедисперсии наблюдаемой случайной величины. Ранее было установлено (см. выше, пример 1, формулу (7)), что в этом случае существует семейство сопряженных априорных распределений параметра 9. __(х-9)2
Определим p САЗ (9) = const и с учетом того, что L( х1,..., хп 19) ~ e 2 ст 2 (см. выше, форму-
лу (7)), имеем:
--^ (х-9)2
pcA3(9|х 1,...,хп)~ pсаз(9)-L(х 1,...,хп |9)~ e 2ст0
Но правая часть этого соотношения представляет собой (с точностью до нормирующего множителя, не зависящего от9) плотность нормального распределения со средним значением х и дисперсией ст0/п. Следовательно, семейство сопряженных априорных распределений неизвестного среднего значения 9 нормально распределенной генеральной совокупности (при известной дисперсии ст 0 = D2,) само принадлежит классу нормальных законов распределения.
Пример 2 (продолжение).
2 £ N 1 где и среднее значение 9, и параметр точности h = 1/D2, являются неизвестными (т. е. 0 = (9,h)). Ранее было установлено (см. выше, пример 2), что в этом случае существует семейство двумерных сопряженных априорных распределений параметра 0 = (9, h).
1
Определим pСАЗ (9) = const и pСАЗ (h) — и с учетом (13) и того, что
h
103
Не 1(9) 2008
n -h У (x,-в )2 n -h\n( в- x )2 +jr (x¡ -i)2 I 1 _nh_, e_x)2 n_1 -lly (x, -x)2 I h
L(x,..,xn | в, h) ~ h 2 e 2 « = h 2 e 2[ " J ~(nh)2 e 2 . h 22 e J , (16)
имеем:
i -nh(в-x)2 n-1-1 -fe(x,-x)2J h ~ слз(в, h|x 1,...,xn)~ p саз (в) ■ pCA3( h) • L( x 1,..., xn |в, h)~( nh)2 e 2 . h 2 e J . (17)
Но правая часть (17) представляет собой (с точностью до нормирующего множителя, не зависящего от в и h, см. Приложение 2) плотность двумерного гамма-нормального распределения
1 -ХоЛ(в-в )2
p(в, h)~( X 0 h)2 e 2 0 . h a-1e "ph (18)
с параметрами X0 = n,в0 = x,a = -—1 и p = — У(x, -x)2.
2 2 /=1
Следовательно, семейство сопряженных априорных распределений двумерного параметра 0 = (в,h), где в и h соответственно среднее значение и параметр точности наблюдаемой нормальной генеральной совокупности, принадлежит классу двумерных гамма-нормальных распределений (18).
Пример 3 (продолжение).
Наблюдаемая случайная величина £,в(M) подчиняется биномиальному з.р.в. с неизвестным значением вероятности «успеха» в и заданным числом испытаний Бернулли М.
Ранее было установлено (см. выше, пример 3), что существует семейство сопряженных априорных распределений параметра в. yx, nM yx.
Определим pCA3 (в) = 1для ве (0; 1) и с учетом того, что L(x,,..., xn |в)~ в = • (1-в) н , имеем:
У x, nM-у x ,
~ САз(в^ ......xn )~ Р САЗ (в) • L(x 1, ..., xn |в)~ в=1 • (1-в) = . (19)
Но правая часть соотношения (19) представляет собой (с точностью до нормирующего множителя, не зависящего от в) плотность бета-распределения
щ §
¡S
¡S
<5 Г(a + b) вa-Ul в b-1
u <ъ
s
¡2 с параметрами a = У x ¡ +1 и b = nM -У x ¡ +1 (участвующая в правой части (20) в выражени
о ,=1 '=1
g нормирующего множителя функция r(z) — это известная гамма-функция Эйлера, т.е.
I Р(в)6 (20)
Г(a) .Г(Ь)
и
2 Г(z) =f xz-1e-xdx).
Cf j
0
¡3. Следовательно, семейство сопряженных априорных распределений параметра 9 (веро-с ятности «успеха») наблюдаемой биномиально распределенной генеральной совокупности § принадлежит классу бета-распределений (20).
и 8
8 Пример 4 (продолжение).
иа Ранее было установлено (см. выше, пример 4), что отрицательно биномиально распределенная случайная величина 2,(9;К) имеет сопряженное априорное распределение пара-
104
Не 1(9)2008
метра 9 — вероятности «успеха» в одном испытании Бернулли. Как и в предыдущем приме-
У х, - Кп
ре, определяем рСАЗ (9) = 1 (для 9 е (0;1)).Тогда с учетом того, чтоЦх,,...,хп 19) ~ 9Кп (1-9) н
(см. выше, пример 4), имеем: ^
о
У х,-Кп
~саз(9|х ......Хпр саз(9)-Ц х 1,..., Хп |9) ~ 9Кп (1-9) и . (21)
Правая часть (21) представляет собой (с точностью до нормирующего множителя, незави-
п
сящегоот9) плотность бета-распределения (20) с параметрами а = Кп +1 и Ь = Ух , - Кп +1.
I=1
Так что семейство сопряженных априорных распределений параметра 9 (вероятности «успеха») наблюдаемой отрицательно биномиально распределенной случайной величины 2(9; К) принадлежит классу бета-распределений (20).
Пример 5 (продолжение).
Ранее было установлено (см. выше, пример 5), что параметр 9 пуассоновскогоз.р.в. имеет сопряженное априорное распределение. Из смысла параметра9следует, что он может при-
1
нимать только положительные значения, поэтому определяем р САЗ (9)—.Тогда с учетом
V 9
У х,
того, чтоЦхь...,хп|9)~ 9 = -е-п9(см. выше, пример 5), имеем:
Х> -1
~ саз (9|х ......хпр саз (9)-Ц х 1,..., хп |9) ~ 9 =1 - е-п9. (22)
Правая часть (22) представляет собой (с точностью до нормирующего множителя, независящего от 9) плотность гамма-распределения
р(9) —9а-1 е-рв, 9 > 0, (23)
Г(а)
п
с параметрами а = Ух, и р = п.
¡=1
Следовательно, семейство сопряженных априорных распределений параметра9 наблюдаемой генеральной совокупности принадлежит классу гамма-распределений (23).
Пример 6 (продолжение).
Ранее было установлено (см. выше, пример 6), что параметр масштаба 9 экспоненциального распределения имеет сопряженное априорное распределение. Поскольку 9 >0, опре-
1, „„ 4ь
деляем рСАЗ (9)—. Тогда с учетом того, чтои(х1,...,хп|9)~ 9п-е > (см. выше, пример 6),
имеем: 9
-I У х, |9
~саз(9|х 1,..., хп)~ р саз (9) -Ц х 1,..., х п |9)~ 9п-1- е ^ ^ . (24)
Правая часть (24) определяет (сточностью до нормирующего множителя, не зависящего
п
от9) плотность гамма-распределения (23) с параметрами а = п и р = У х, .Так что семейство
¡=1
105
Hb 1(9) 2008
сопряженных априорных распределений параметра масштаба 9 экспоненциально распределенной генеральной совокупности принадлежит классу гамма-распределений (23).
Пример 7 (продолжение).
Как мы видели (см. выше, пример 7), и при равномерно распределенной на отрезке [0;9] случайной величине неизвестный параметр 9 имеет сопряженное априорное распределение. Поскольку параметр 9 может принимать любые положительные значения,
1 f 1 ^n определяем pСАЗ (9)—. Тогда с учетом того, чтоL(xb...,xn 19) = 1 — 1 (и xmax(n) = maxXi <9),
9 v 9 ) i<i<n
имеем: v y
~ САЗ(0|Х 1,..., Xn
1 \ n+1
p саз (9)' L(x 1,..., Xn |0)~ при -> x max( n) (25)
0 при e< X max ( n)
Но правая часть соотношения (25) представляет собой (с точностью до нормирующего множителя, не зависящего от 9) плотность распределения Парето вида
(аВа
р(9) =|"9^т- ПРи 9 > 9т1п (26)
[0 при 9<9т1п
с параметром формы а = п и некоторым параметром сдвига 9т,п > хтах (п). Следовательно, семейство сопряженных априорных распределений параметра 9 равномерно (на [0;9]) распределенной случайной величины принадлежит классу распределений Парето вида (26).
Пример 8 (продолжение).
В данном примере речь идет о наблюдаемой генеральной совокупности, подчиняющей-щ ся распределению Парето с неизвестным значением параметра формы 9 и некоторым задан-| ным значением параметра сдвига х0 (см. выше, пример 8), так что
5 ^ ^п9 -Г п 1п Г Ми
1 L(x......х„ |0) =0| • g-nn ~ 0n • е
<S i Iх о
Щ S " -ч —
Л где gn = П I . Но тогда САЗ-апостериорная функция плотности распределения парамет-
I V '=1 ) ( 1 I
g ра 0 будет иметь вид | с учетом того, что pСАЗ (0) — I:
§ V 0)
" ~ ,л| ^ ,, л„ , 4nln(SlH.0
~ п.„ 1
" ~ саз(0|х 1,..., Xn)~ p саз (0) • L(x i,..., х„ |0) ~ 0n n1 • е L VXo )j . (27) g
g Мы видим, что правая часть соотношения (27) определяет (с точностью до нормирующе-
'! го множителя, не зависящего от параметра 0) плотность гамма-распределения (23) с пара-
§ метром а = n и параметром ß = n ln| — I, так что сопряженные априорные распределения
£ V х о )
иа параметра формы 0 наблюдаемой Парето-распределенной генеральной совокупности принадлежат семейству гамма-распределений.
106
Не 1(9)2008
Пример 9 (продолжение).
Выше при рассмотрении нормальной классической линейной модели множественной
I
<в
т >Х
регрессии с неизвестными значениями коэффициентов регрессии 0 - (90,9,,...,9к)' и па- ч
1 2 ч
раметра точности Ь - —- (где а2 - Ое,) мы убедились в том, что у параметров 90,9,,...,9к, Ь о а2
существуют сопряженные априорные распределения. Определим теперь общий вид сопряженного априорного распределения р(0; Ь) параметров 0 и Ь. С учетом «Замечания 2» (см. выше) и положительных значений параметра Ь имеем:
1
рСАЗ (90,9,, ...,9к; Ь) - рСАЗ (90 ) • рСАЗ (9, ) • ... • РСАЗ (9к ) • РСАЗ (Ь) ~ -.
Ь
Используя полученное ранее выражение (12) для функции правдоподобия ЦX; У10;Ь), имеем:
, 1 п -Г^а21л -—(00-0 )т (хт х)(0-0)
~саз(0;Ь|X;У) ~ рсаз(0;Ь) • иX;У|0;Ь) ~ -• Ь 2 • е 1 2 > • е 2 -
Ь
П=*=2_л -Г^а21 ь А 0-0 )т (х т х)(0-0)
- Ь 2 • е 1 2 ^ • Ь 2 • е 2 . (28)
Но правая часть соотношения (28) определяет (с точностью до нормирующего множителя, не зависящего от 0 и Ь) так называемое многомерное гамма-нормальное распределение
а - п - к -1 „ п - к - К,
с параметром сдвига (0, матрицей точности (X 1 X)и параметрами а-—-— и р-—-—а2
(подробнее о многомерном гамма-нормальном распределении и его свойствах см. в Приложении 2). Напомним, что 0 - (XтX)-1 XтУ — это МНК-оценка параметров регрессии 0, а 1
а2 --;-(У- X0)Т(У - X0) — оценка остаточной дисперсии а2.
п - к -1
Таким образом, сопряженные априорные распределения параметров (0;Ь) нормальной классической линейной множественной регрессии имеют общий вид:
р(0;Ь) ~ Ь^|А0|2 е0-00>тЛ0<0-00> • Ьа-1е^, (29)
в котором конкретное задание векторного параметра сдвига 00, (к +1)х(к + 1)-матрицы точности Л0 и скалярных параметров а и р однозначно определяет априорный закон распределения параметров 0 и Ь (напомним, что к +1 — это общее число объясняющих переменных, включая свободный член, в анализируемой модели регрессии).
Очевидно, семейство многомерных гамма-нормальных распределений (29) является многомерным обобщением двумерного гамма-нормального распределения (18).
3.3. Рекомендации по подбору конкретных значений параметров в сопряженных априорных распределениях
Использование в качестве априорных законов распределения вероятностей (з.р.в.), сопряженных с наблюдаемой генеральной совокупностью (в ситуациях, когда они существуют), позволяет нам определить их общий вид, т. е. задает целое семейство априорных распределений {р(0;Э)}. Однако при реализации байесовского подхода мы должны оперировать конкретным априорным распределением, что требует знания числовых значений Э0 параметров Э, от которых наш априорный з.р.в. зависит. Как же подбирать эти значения Э0 в каж-
107
Не 1(9) 2008
дом конкретном случае? Ниже описывается один из возможных подходов к решению данной задачи.
В широком классе ситуаций можно исходить из того, что нам известны априорные средние значения оцениваемого параметра 00 = Е0 = (Е91,Б92, ...,Е95)т и их среднеквадратиче-ские ошибки Д, = 7097, д 2 = ->/о97,..., Д 5 = л/оё" . тогда параметры априорного распределения, как правило, могут быть определены методом моментов (в случае многомерного параметра 0 — с учетом «Замечания 2» о статистической независимости компонент вектора 0 в априорном распределении, см. (13)).
Продемонстрируем реализацию этого подхода на рассмотренных выше примерах.
1) Определение параметров априорного гамма-распределения (см. формулу (23) и примеры 5, 6, 8). Как известно (см., например, [Айвазян, Мхитарян (2001)]), среднее значение (Е9) и дисперсия (09) гамма-распределения выражаются через параметры а и р этого распределения по формулам:
Е9 = а, 09 = —.
Р р2
Подставляя в эти соотношения вместо Е9 и 09 соответственно заданные значения 90 и Д2, получаем в качестве решений системы из двух уравнений (относительно а и р):
а = е|, р = е^.
Д2 Д2
(30)
2) Определение параметров априорного бета-распределения (см. формулу (20) и примеры 3 и 4). Используя выражения для среднего и дисперсии бета-рапсределения (см., например, [Айвазян, Мхитарян (2001)]) и решая систему из двух уравнений
щ § §
¡В
I
и <ъ
8-<ъ 5
0
1
т оз
! с
! и
8 и <ъ >5 <в
Щ
Е9 = -
а + Ь
= 9 0
09 =
аЬ
(31)
(а + Ь )2( а + Ь +1)
= Д2
относительно а и Ь, получаем:
9 0(1-9 0
а=
-90, Ь =
9 0(1-9 0
-90
1-9 0 9 0
3) Определение параметров априорного распределения Парето (см. формулу (26) и пример 7). В данном случае параметр формы а и параметр сдвига 9т,п определяются поза-данным значениям 90 = Е9 и Д2 = 09 из системы уравнений
Е9 =
а9т
а-1
= 9 0
(32)
09 =
= Д2
(а- 1)2(а- 2)
Решение этой системы относительно а и 9тп дает:
а = 1 + ,Ц + Д-, 9тш =^90 • (а-1).
1
—I
а
(32')
2
а9т
108
Не 1(9)2008
4) Определение параметров двумерного гамма-нормального распределения (см. формулу (18) в примере 2). Из свойств двумерного гамма-нормального распределения следует (см. Приложение 2), что частное априорное распределение параметра Ь есть гамма-распределение с параметрами а и р. Поэтому, воспользовавшись заданными значениями Ь0 = ЕЬ и Д2Ь = ОЬ, составляем систему из двух уравнений относительно а и Р:
ЕЬ=1=Ьо .
ОЬ = —= Д2Ь'
Получаем решение:
Ь2о п Ьо а = —- и р = —.
Д2Ь Д2Ь
(33)
(
дения о t
а
р
( аЛ
ленияследует, чтоEt 2а|90;X0 —
V р
= 9 0 и Оt
/
(а 2а|9 0; X 0 р
р
2а
X 0 - а 2а - 2
(см. Приложение 1),
так что при заданных значениях 90 = Е9 и Д29 = О9 имеем:
I
<в £
ч ч
о
Для определения параметра X0 и параметра сдвига 90 воспользуемся тем, что частное априорное распределение параметра сдвига 9 есть обобщенное распределение Стьюдента
с 2а степенями свободы, параметром сдвига 90 и параметром точности, равным X0 а (све-
2а| 9
0; X0 ~ -распределении см. в Приложении 1). Из свойств этого распреде-
• Значение параметра сдвига в распределении (18) равно90;
Д2 = р
Д29 =
а л 1 Р
;, откудаX0 =7— --
X 0 - а а -1 " Д29 а -1
(напомним, что а и р уже определены соотношениями (33)).
(34)
5) Определение параметров многомерного гамма-нормального распределения
(см. формулу (29) в примере 9). Воспользуемся свойствами многомерного гамма-нормального распределения (см. Приложение 2). В соответствии с ними:
О) частное распределение числового сомножителя Ь является гамма-распределением с параметрами а и р;
(и) частное распределение параметра 0 есть обобщенное (к + 1)-мерное распределение
Стьюдента с 2а числом степеней свободы, параметром сдвига 00 и матрицей точности а
В =рЛ0 (мы обозначаем его как t (2а| 0 0; В ^распределение).
Свойство 0) позволяет (при заданных значенияхЬ0 = ЕЬ и Д2Ь = ОЬ) определить значения параметров а и р по той же формуле (33).
Свойство (П), дополненное «Замечанием 2» и правилами вычисления вектора средних
(а
значений и ковариационной матрицы (к + 1)-мерной случайной величины t 2а|00;— Л(
р
109
Не 1(9) 2008
(см. Приложение 1), позволяет определить остальные параметры распределения (29) — параметр сдвига 00 и элементы матрицы Л0. Действительно:
Е0 = ЕГ
2а|00; —Л0 |=00 (задано!)
2Э = 2
2а
г(2а|0 о;-—Л о> 2а- 2
а
-Лг
( Д2 Д 0
(35)
где Ду — заданные значения априорных дисперсий компонент вектора 0 = (90,91, ...,9к), } = 0,1,..., к.
Таким образом, векторный параметр сдвига в распределении (29) определяется заданным вектором априорных средних значений 00, а диагональные элементы Я(0) (} = 0,1,...,к) матрицы Л0 определяются из уравнений (35) по формуле:
Я.'}) =— •
1 —
Д2 а-1
(36)
щ § §
¡В
I
и <ъ
!
<ъ §
0
1
п
03
! с
1 и
8 и <ъ >5 <в
Щ
где значения а и — определены соотношениями (33).
4. Пересчет значений параметров при переходе от априорного сопряженного распределения к апостериорному
Поскольку, по определению, семейство сопряженных априорных распределений {р(0; Э)}замкнуто относительно операции (4) пересчета априорного распределения в апостериорное, то общий вид апостериорного распределения ~(0|ХЬ...,Хп) при использовании сопряженных априорных распределений нам известен, и нам лишь надо уметь пересчитывать параметры 0(Х,,...,Хп) этого апостериорного распределения по заданным параметрам Э0 априорного распределения и имеющимся наблюдениямX,,Х2,...,Хп.
Общая схема такого пересчета следующая. Пусть {р(0; Э)} — семейство априорных распределений, сопряженных с функцией правдоподобия Ь(х],...,хп 10) имеющихся у нас наблюдений (Э = (<,..., Сц) — вектор параметров, от которых зависит сопряженное априорное распределение р(0; Э)), и пусть Э0 — заданные (известные) значения параметров Э в анализируемом случае.Тогда с помощью ряда тождественных преобразований правая часть соотношения
~(0|Х 1,...,Хпр(0; Э0)Р(Х 1,...,Хп |0)
(37)
приводится, с точностью до множителей, не зависящих от 0, к виду р(0; Э(ХЬ..., Хп)), где последняя функция принадлежит семейству р(0; Э), а каждая из компонент с1}(Х^...,Хп) (} = 1, 2,...,ц) вектора параметров Э(Х1,...,Хп)является функцией от Э0 и {Х1,Х2,...,Хп}.
Продемонстрируем реализацию этой общей схемы на наших примерах (с разной степенью подробности).
0
Д21
0
110
Не 1(9)2008
Пример 1 (продолжение)
х |9)~ е 2(см (7)) п(9- П) =
42жД
(9-4)2__^ V _й) 2 1
~(9 |х 1,...,хп)~ е~ 2^2 - е 2о2/'
где
N С/-01 \ ° 0/П Д2)
^(х 1,..., Хп) = Е(9|х 1,..., Хп) =-
ё 2 (х 1,..., Хп) = О(9|х 1,..., Хп) =
1 1
-н--
о2/ п Д20
( 1 1 1-1
(38)
2/ п Д20 )
X 0 = X 0 + п, ~0 =
пх + X 0 9 0
п + X 0
~ п а = а +—,
2
Р=Р + ^(х( -х)2
2 ' = 2|- + — V п X 0
(39)
Пример 3 (продолжение).
Непосредственная реализация соотношения (37) в данном случае дает:
111
г
2 о0(хх-9 )',... 1 .- п 5
ё2 = Д2°), такчто ^2лД 0 |
Вданном примере Цх,, х2,..., хп |9)~ е 2о° (см. (7)), р(9; П) = 1__— е 2Д2° (т.е. ё, = 0°, ¡2
■(х-9)2 9-~ )2 О
2ё2
1 _ 1 — - х + — -9 о
Необходимые промежуточные выкладки нацелены на выделение полного квадрата раз-
~ 1 1 _ ности (9-ёт)2 из выражения-(9-ё1)2 +—^—(9-х)2 и не представляют принципиаль-
аё2 2 о 0/п
ных трудностей.
Мы видим, что среднее () и дисперсия (ё2) апостериорного нормального распределения являются определенным образом средневзвешенными значениями априорных и выборочных соответственно средних и дисперсий.
Пример 2 (продолжение).
При реализации общей схемы пересчета априорных параметров в апостериорные вданном случае следует учесть представление функции правдоподобия I в форме (16) (см. выше, пример 2), вид (18) априорной плотности двумерного гамма-нормального распределения (в котором вектор параметров П° = (X°; 9°; а; р)), а также справедливость тождества
п(9-х)2 +Xо(9-9°)2 = (Xо + п)^°°9° + пх 12 (9° -х)2.
V X ° + п ) X ° + п
Тогда вычисление р(9;Ь) по схеме (37) приводит нас снова к двумерному гамма-нормальному распределению вида (18), но с параметрами
и
Нв 1(9) 2008
У х, пМ-у*, а+ух, -1 Ь+пМ-ух,-1
О(9|х 1,...,Хп)~ еа-1(1 -9)Ь-1 -9« (1-9) = =9 = (1-9) « . (40)
Но правая часть (40) определяет (с точностью до нормирующего множителя) снова бета-распределение с параметрами:
п
а = а +у х,,
(41)
п
~ = Ь + пМ-у х,.
,=1
Пример 4 (продолжение).
Подставляя в правую часть соотношения (37)
р(9; О о) = р(9; а,Ь)~ 9а-1(1-9)Ь-1,
Ух, -Кп
иX1,...,Хп |9)~ 9Кп (1-9) " ,
имеем:
~(9|х 1,...,Хп)~ 9а+Кп-1(1 -9)
Ь+у х, - Кп-1
Мы видим, что апостериорное распределение параметра (вероятности «успеха») отрицательно-биномиального закона, так же как и априорное, является бета-распределением и что его параметры О = (а,Ь) определяются соотношениями:
п
а = а + Кп, ~ = Ь +у х, - Кп. (42)
=1
Пример 5 (продолжение).
Как мы видели ранее, функция правдоподобия наблюдений пуассоновской генеральной совокупности имеет вид:
8 Ух,
| Цх 1,...,хп |9)~ 9=1 е-п9. £
^ Так что, используя в качестве априорного распределения р(9;О) = р(9;а,р) параметра 9
§ гамма-распределение (23), имеем:
щ
^ п п
5 У х, а+у х, -1
¡2 ~(9|х 1,...,хп)~ 9а-1е-1м-9» е-п9 =9 « е"
§ о
| Тем самым подтверждается сопряженность априорного гамма-распределения, причем
^ апостериорное гамма-распределение определяется параметрами О = (а,р), где
! п ~
§ а = а +у х ,; р=р + п. (43)
I; ,=1
§ Пример 6 (продолжение).
§ Функция правдоподобия экспоненциально распределенных наблюдений (с параметром щ масштаба 9) имеет вид:
Г п 1
иа -[ух,^
Цх 1,...,хп|9) = 9пе ) .
112
Так что при априорном гамма-распределении параметра 9 имеем:
Не 1(9)2008
I
<в
-IVх, |9 -|Р+Е"I |9
~(9|хх п)~ 9а-1е-р9 -9 пе ^ = 9а+п -1е 1 " ^ . ч:
Ч
О
Мы видим, что апостериорное распределение параметра 9 снова подчиняется закону гамма-распределения (23), но с параметрами:
а = а + п,
Р=Р+Х х,. (44)
Пример 7 (продолжение).
Подставляя в правую часть соотношения (37) функцию правдоподобия равномерно распределенных (на отрезке [0;9]) наблюдений и функцию плотности распределения Парето (26) в качестве априорного распределения р(9; О) = р(9; а; 9т,п), имеем:
а9а 1
~(9|х 1,...,Хп) ~ —т^- — (при 9 > тэх(9т1п; х 1, х2,...,Хп}).
9 9
Отсюда следует, что апостериорное распределение параметра9 описывается, также как и априорное, законом Парето (26), но с параметрами:
а = а + п, ~т1п = тах{9тп;х 1,х2,...,хп}. (45)
Пример 8 (продолжение).
Как мы видели (см. выше, пример 8), функция правдоподобия Парето-распределенных наблюдений имеет вид:
п I п | —1]-9
Цх 1,...,хп|9)~ 9п -е
Подставляя ее в правую часть соотношения (37), а также, в качестве априорного распределения р(9;О0), плотность гамма-распределения (23), имеем:
п1п Г ^ Не
~(9|х 1,...,хп)~ 9а-1е-9пе что определяет гамма-распределение с параметрами:
а = а + п, Р=Р + п 1п| ^ |, (46)
где дп = [Цх1 — среднее геометрическое наблюдений х,,...,хп, а х0— параметр сдвига в анализируемом распределении Парето (его значение считается заданным).
Пример 9 (продолжение).
Байесовское оценивание коэффициентов регрессии 0 = (90,9Ь...,9к) и параметра Ь в нормальной классической модели множественной регрессии (8)-(9) предполагает использование апостериорного распределения ~(0;Ь|X,У) этих параметров, определяемого по схеме (37). Подставляя в правую часть соотношения (37) в качестве априорного многомер-
х
0
113
Не 1(9) 2008
ное гамма-нормальное распределение (29), а также функцию правдоподобия ЦX,У10,Ь) (12), преобразованную к виду
ЦХ,У|0;Ь)~ Ь
■ Ь 2 е 2
(3-0)1 (X 1X)(0-0)
получаем после ряда тождественных преобразований (см. [Де Грот (1974)]) апостериорную плотность р(0;Ь | X,У) в форме многомерного гамма-нормального распределения (29), параметры которого определяются по параметрам 00,Л0,а и р априорного распределения и наблюдениям (X,У) следующими соотношениями:
0О = (Л0 + XтX) 1(Л000 + XтУ) — параметр сдвига; Л0 =Л0 + XтX — матрица точности;
~ п а =а + —;
2
1,
Р =Р + ^[(У - X0o)т У + (0 0 -0о)т Л 0 0 о]
параметры частного апостериорного гамма-распределения параметра точности Ь.
(47)
к- -[п-^а21 ь
5. Примеры задач на точечное и интервальное байесовское оценивание
параметров модели
Задача 1. Анализ закона распределения домашних хозяйств определенной социально-экономической страты в заданном регионе страны по величине среднедушевого дохода г|.
Мы располагаем следующей информацией об анализируемой генеральной совокупности:
(а) логарифм (натуральный) величины среднедушевого дохода (т. е. 2 = 1п | домашних хозяйств рассматриваемой страты данного региона распределен нормально с неизвестным | средним значением 9 и известной дисперсией ст2 = 0,28;
§ (б) имеются результаты обследования п = 10 случайно отобранных от анализируемой | страты домашних хозяйств по среднедушевому доходу у, (в нижеследующей таблице даны § значения х, = 1п у -,):
I 1 2 3 4 5 6 7 8 9 10
XI 0,54 1,20 0,36 0,80 0,42 2,10 0,70 0,25 0.90 0,48
(в) из предыстории и опыта обследования домашних хозяйств той же страты в других ре-
§: гионах страны получены априорные значения среднего Е9 = 90 = 0,60 и дисперсии
¡1 Е9 = Д20 = 0,03.
о
с
| Требуется:
<§ Используя сопряженное априорное распределение параметра 9, получить байесовские точечную и интервальную (с уровнем доверия Р0 = 0,95) оценки средней величины логарифма среднедушевого дохода и сравнить их с соответствующими оценками метода максимального правдоподобия.
114
Не 1(9)2008
Решение. Мы уже знаем (см. выше, пример 1), что сопряженное априорное распределение в данном случае существует и является нормальным, причем параметры этого распределения непосредственно заданы (Б0 = 00 = 0,60 и 00 = Д 20 = 0,03). В соответствии с выведенными выше формулами пересчета (см. п. 4, формулы (38)) имеем:
00 = Е(0|х 1,...,хп) =
ст2/ п Д20
ст 0/ п
-0 0
= 0,691,
8 <в
Ч Ч
О
Д20 = О(0|х 1,...,Хп ) =
СТ 0/ п
Д20
= 0,015.
Соответственно:
0ю = Е(0|х 1,...,хп) = 0,691
и с вероятностью Р0 = 0,95 можем утверждать, что 0(Ь) - и 0,025 • Д 0 <0< 0(Ь) + и 0,025 - Д.
С учетом того, что 2,5%-ная точка стандартного нормального распределения и0 025 = 1,96 и Д 0 = д/ДГ = 0,120, имеем:
0е [0,451; 0,931] с вероятностью Р0 = 0,95.
Решение этих же задач, основанное на методе максимального правдоподобия, дает:
0мп = Х = 0,775 и 0е [0,447; 1,103] с вероятностью Р0 = 0,95
(концы последнего доверительного интервала вычислены по формулам 0мп ±и0 025 •СТ^).
' л/п
Мы видим, что использование априорной информации о неизвестном значении параметра 0 = Е(1п и применение, соответственно, байесовского подхода в данной задаче позволили уточнить оценку и, в частности, сузить интервальную оценку по сравнению с классическим подходом почти в полтора раза.
Задача 2. Оценка интенсивности вызовов, поступающих на пункт «Скорой помощи» в час.
Число вызовов 2,, поступающих на пункт «Скорой помощи» в час, описывается распределением Пуассона с неизвестным значением параметра 0 = Е2, (см. выше, пример 5). Результаты регистрации числа вызововх, (в час), зафиксированные втечение одной смены (длящейся 8 часов), приведены в следующей таблице:
I 1 2 3 4 5 6 7 8
XI 3 1 4 2 6 3 3 2
Из опыта работы аналогичных пунктов определено априорное среднее значение 00 = Е0 = 3,6, причем случайный разброс значений этого параметра характеризуется дисперсией Д20 = 00 = 0,09.
115
Д2
0
Нв 1(9) 2008
Требуется:
Используя сопряженное априорное распределение параметра 9, получить байесовские точечную и интервальную (с уровнем доверия Р0 = 0,95) оценки средней интенсивности 9 = Е2, вызовов, поступающих на пункт «Скорой помощи», и сравнить их с соответствующими оценками метода максимального правдоподобия.
Решение. Как было установлено выше (см. пример 5), сопряженное априорное распределение параметра 9 в этом случае существует и описывается гамма-законом, параметры а и р которого определяются (в соответствии с рекомендациями п. 3.3) из системы уравнений:
Е9 = а = 3,60
Р
Э9 = — = 0,09
. Р2
Отсюда а = 144 и р = 40.
В соответствии с выведенными выше (см. п.4) формулами пересчета (43) параметров апостериорного гамма-распределения имеем:
8
а=а+£ х, = 144 + 24 = 168,
! =1
р = Р + п = 40 + 8 = 48.
таким образом:
9(5) = Е(9|х 1,...,Хп) = а = 3,5,
Р
и можно утверждать, что с вероятностью Р0 = 0,95 справедливы неравенства
§ У 0,975 (а, Р ) <9<у 0 ,025 (а, р),
§
¡В
<5 шись известными формулами (см. приложение 2):
и
щ 1
| у „ (~;р) = х 2 (2а)
Ег <ъ 5
0
1
т аз
I 2
§ гдехч (т) и ич — 100д%-ные точки «хи-квадрат»-распределения и стандартного нормального >¡5 распределения, соответственно, имеем:
I
со 9е [2,97; 4,03] с вероятностью Р0 = 0,95.
и
>5 Решение этих же задач, основанное на методе максимального правдоподобия, дает: иа
9 мп = х = 3,0;
где уч(а, р) — 100д°%-ная точка гамма-распределения с параметрами а и р. Воспользовав-
2р
и, при т > 100:
X2(т) « т + ия ^2т,
116
Не 1(9)2008
9 е [1,80; 4,20;] с вероятностью Р0 = 0,953
Мы видим, что в данном случае использование априорной информации о параметре 9 в рамках байесовского подхода позволило сузить размах интервальной оценки более чем в два раза!
Задача 3. Оценка «необходимой доли брака»9 в продукции, производимой автоматической линией.
Предприятие приобрело новую автоматическую линию (АЛ). Для оценки так называемой «необходимой доли брака»9 — вероятности того, что произведенное этой АЛ в режиме стационарного функционирования изделие окажется некондиционным, — было проконтролировано п = 5 партий по М = 80 изделий в каждой партии. Число дефектных изделий £,, обнаруженных в партии изделий объема М, адекватно описывается биномиальным з.р.в. с параметрами 9 и М. Результаты контроля представлены в таблице (х, — это число дефектных изделий, обнаруженных в ¡-й проконтролированной партии):
/ 1 2 3 4 5
XI 2 0 3 1 2
8
<в £
ч ч
о
Кроме того, проведенный анализ работы аналогичных АЛ, установленных на других предприятиях, показал, что «необходимая доля брака» в среднем равна 0,01 и имеет разброс, характеризуемый среднеквадратическим отклонением 0,003.
Требуется:
Используя сопряженное априорное распределение параметра 9, получить байесовские точечную и интервальную (с уровнем доверия Р0 = 0,90) оценки «необходимой доли брака»9 и сравнить их с соответствующими оценками метода максимального правдоподобия.
Решение. Выше (см. п. 3) было установлено, что сопряженное априорное распределение параметра 9 в данном случае существует и описывается бета-распределением, параметры а и Ь которого определяются из системы (см. п.3.3):
а + Ь
= 0,01
аЬ
[а + Ь )2( а + Ь +1)
= (0,003)2
Решение этой системы дает а = 10 и Ь = 990. Воспользовавшись формулами пересчета (41), получаем значения параметров а и Ь апостериорного распределения 9:
а = а х, = 10 + 8 = 18,
¡ =1
5
~ = Ь + 5 • 80 "У х, = 1382.
3 Данная интервальная оценка основана на асимптотической правдоподобия 9 мп.
9;л(— -нормальности оценки максимального
117
а
Не 1(9) 2008
Таким образом:
9(5) = Е(9|х 1,...,хп) = —= 0,01286, — + —
и можно утверждать, что с вероятностью Р0 = 0,90 справедливы неравенства
р 0,95 (а, — ) < 9 < р 0,05 (—, —),
гдерч (а,—)— 100д°%-ная точка бета-распределения с параметрами а и — . Воспользовавшись известными равенствами
— — (2—;2—) 1 р„(а;— = —-^ — , (У1,V2) -
— + — (2—;2—) Ря (V2, V1)
и таблицами 100д-процентныхточек (V,, V 2) распределения Г с числами степеней свободы числителя V и знаменателя V2, имеем:
9е [0,0083; 0,0182] с вероятностью Р0 = 0,90.
Решение этих же задач, основанное на методе максимального правдоподобия (см., например, [Айвазян, Мхитарян (2001б), задача 1.18]), дает:
9 мп =Л & = 0,02,
пМ Т"1
9е [0,0085; 0,0315] с вероятностью Р0 = 0,90.
Размах этой интервальной оценки в 2,3 раза п ревосходит ширину байесовской интервальной оценки!
Задача 4. Оценка интервала движения автобуса.
Приходящий в случайные моменты времени на остановку п ассажир втечение п яти своих | п оездок фиксировал время ожидания автобуса (в минутах): х^ = 1,2; х2 = 2,5; х3 = 0,5; х4 = 3,2; § х5 = 2,9. Известно, что автобус ходит строго п о рас п исанию с интервалом в 9 минут, так что | время ожидания автобуса п ассажиром можно считать случайной величиной £,, п одчиненной § [0;9]-равномерному з.р.в. (см. выше, п ример 7). Пытаясь оценить интервал движения автобу-£ са, п ассажир сумел п олучить до п олнительную информацию о п араметре 9: из анализа о п ыта Ц работы различных автобусных маршрутов города, функционирующих в едином регламент-I ном режиме, следовало, что среднее значение этого п араметра равно 5,38 мин., а случай-§ ный разброс в его значениях характеризуется средним квадратическим отклонением, рав-* ным 1,39 мин.
03
| Требуется:
| Ис п ользуя со п ряженное а п риорное рас п ределение п араметра 9, п олучить байесовские точечную и интервальную (с уровнем доверия Р0 = 0,95) оценки для неизвестного интервала | движения автобуса и сравнить их с соответствующими оценками, основанными на методе 8 максимального п равдо п одобия.
>5
иа Решение. В п . 3 (см. п ример 7) было установлено, что со п ряженное а п риорное рас п ределение п араметра 9 в данном случае существует и о п исывается рас п ределением Парето с п а-
118
Не 1(9)2008
раметром формы а и параметром сдвига 9™п, которые определяются из системы уравнений (32), т.е. по формулам (32'). В нашем случае имеем:
90 , I., 5,382
а = 1 + ,1 + -2 = 1 + .1-. 2 А2 V 1,392
= 5,00,
<в
ч
о
9т|п = -90 • (а-1) = -• 5,38 • 4 = 4,30 (мин.). а5
Параметры апостериорного распределения Парето определяются формулами пересчета (45):
а = а + п = 5 + 5 = 10,
9тп = тах {9тп; х 1,..., х 5} = 4,30 (мин.
Соответственно:
9(Б) = Е(9|х 1,..., х 5) =
а • 9т1п
а-1
■ = 4,78 (мин.)
и 9е [90,975 (а; 9гтп);90,025 (а; 9т[п)] с вероятностью Р0 = 0,95, где 9ч(а; 9т|п) — это 100д%-ная точка распределения Парето с параметрами (а; 9т,п). Поскольку функция распределения Парето определяется соотношением
Г(9) = Р {^(а;9т,п) <9} = 1-
!"9 ^
^тип
то значения 90,975 (а; 9т,п) и 90,025 (а; 9т,п) определяется из уравнений соответственно:
( ~ V
9т
(
90,975 (а; 9г
9
тт /у
Л
90,025 (а; 9т
= 0,975,
= 0,025.
V
Решение этих уравнений относительно90,975 (а; 9т,п) и 90,025 (а; 9т,п) при а = 10 и 9т,п = 4,3
дает: ~ ~
90,975 (Й; 9т1п ) = 4,31 и 90,025 (а; 9т1п ) = 6,22,
так что
9е [4,31; 6,22] с вероятностью Р0 = 0,95.
Решение тех же задач, основанное на методе максимального правдоподобия, дает (см. [Айвазян, Мхитарян (2001б), задача 1.22]):
9 мп = 3,84 (мин.),
9е [3,22; 6,69] с вероятностью Р0 = 0,95
(здесь дается оценка максимального правдоподобия 9мп, подправленная на несмещенность). Ивданномслучае байесовский подход позволил сузить ширину доверительного интервала почти в 2 раза (точнее, в 1,82 раза).
119
Не 1(9) 2008
Задача 5. Оценка параметров модели зависимости душевых доходов от объема автономных инвестиций.
В нижеследующей таблице приведены макроэкономические данные по США, характеризующие среднедушевой доход yt и автономные инвестиции xt (в долларах, в дефлирован-ных, с помощью индекса стоимости жизни, ценах) за 1922-1941 годы. Инвестиции определены приближенно как разность между среднедушевым доходом и душевым расходом на личное потребление (данные заимствованы из работы Haavelmo T. Methods of Measuring the Marginal Propensity to Consume. — JASA, vol. 42 (1947), pp 105-122).
t 1 (1922) 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 (1941)
Xt 39 60 42 52 47 51 45 60 39 41 22 17 27 33 48 51 33 46 54 100
yt 433 483 479 486 494 498 511 534 478 440 372 381 419 449 511 520 477 517 548 629
Анализируется нормальная классическая линейная модель парной регрессии (см. пример 9 в пункте 2 при к = 1):
ух = 90 +9х + еГ, Г = 1,2, ...,20. (48)
Анализ предыстории и экспертных оценок модели позволил получить следующую априорную информацию о значениях параметров 90,9! и Ь = (йеГ)_1:
90 = Е90 = 330; 90 = Е91 = 2,85; Ь0 = ЕЬ = 0,002;
Д20 = D9 0 = 225;
Д2 = D91 = 0,01;
Д2 = Dh = 25 • 10 -
<ъ §
¡S
¡S
I
и <ъ
!
<ъ §
0
1
п ta
! С
! U
8 U
<ъ >s <в
Щ
Требуется:
Используя сопряженное априорное распределение параметров(90; 9!; Ь), получить байесовские точечные и интервальные (с уровнем доверия Р0 = 0,90) оценки этих параметров и сравнить их с соответствующими оценками метода максимального правдоподобия.
Решение. Проведенный в пунктах 2и3 анализ примера 9 показал, что в данном случае существует сопряженное с наблюдаемой генеральной совокупностью распределение параметров (90; 9Ь' Ь) и что оно описывается трехмерным гамма-нормальным распределением (29) с параметрами 00 = (90;9° )т,Л0,а и р, определяемыми в соответствии с рекоменда-
h0
h0
циями (33) и (36), т. е.: 00 = (330; 2,85)T; а = -2- = 16; P = fe = 8000;
Л 0 =
Д h
2,37 0 0 53333,3
Д2
Параметры апостериорного гамма-нормального распределения вычисляются в соответствии с формулами пересчета (47):
0 0 =
349,0 2,9
а = 26, P = 14578, Л0 =
22,37 907
907 100 176
Точечные байесовские оценки параметров (90; 9!; Ь) определяются средними значениями соответствующих частных апостериорных распределений. С учетом свойств (/) и (¡¡) многомерного гамма-нормального распределения (см. выше, п. 5 раздела 3.3) имеем:
120
Не 1(9)2008
©(Б) = Е(©|ХУ) =©0 = (349,0; 2,90)т, £
«
~ §
й (Б> = Е( Ь|ХУ) = а = 0,00178.
Р Ч
При выводе интервальных байесовских оценок также используются свойства (/) и (//) ° многомерного гамма-нормального распределения, а также факт г (2а)-распределенно-сти случайных величин (0 (Б) -0 1 )Л/В~ где параметр точности В вычисляется поблочным компонентам матрицы точности ~ частного апостериорного обобщенного многомерного Г(2а| © (Б>;В>-распределения по формуле
В = ~, -~~(]> • (49)
(см. Приложение 1б>. Участвующие в этом соотношении число ~, 1х(к - 1)-матрица В., (к -1)х1-матрица В и (к -1)х(к - 1>-матрица ~(]> определяются следующим блочным представлением матрицы
К1)
а В (0,040 1,618 В нашем случае к = 2, 1 = 0 или 1, матрица В = ~Л0 =
Р
(50)
так что
1,618 178,665
В0 = 0,040-(1,618)7178,665 = 0,0254 и В, = 178,665-(1,618)70,040 = 113,217. Следовательно, с вероятностью Р0 = 0,90 мы можем утверждать, что 100Б)-000,0254 < Г0,05 (52) и 10(Б) -0^7113,217 < Г0,05 (52), так что (с учетом того, что Г0,05(52) = 1,676) имеем:
00 е [338,5; 359,5] с вероятностью Р0 = 0,90,
01 е [2,743; 3,057] с вероятностью Р0 = 0,90.
Поскольку параметр й подчиняется апостериорному гамма-распределению с параметрами а и р, то
й е[у 0 95(а;р); у 0 05(а;р)] с вероятностью Р0 = 0,90. ~ 1
Используя соотношениеуч(а;р) =—%2(2а), имеем (сучетом%2,95(52)«36,4 и %2,05(52)« «69,8): 2р
й е [0,00125; 0,00239] с вероятностью Р0 = 0,90.
Оценивание модели (48) с помощью метода максимального правдоподобия (дающего в данном случае те же результаты, что и метод наименьших квадратов) приводит к следующим точечным и интервальным оценкам:
т ^ -1 уТ\/ _ Ю ЛЛ 1. О П1^т
й =
' ' мп
= (Xт X)-1 Xт Г = (344,7; 3,05)
1 -1 —(У - Х©мп)т(^ - Х©мп) = 0,0015,
1
316,1 < 00 < 373,3; 2,45 <01 < 3,64 с вероятностью Р0 = 0,90, й е [0,00093; 0,00285] с вероятностью Р0 = 0,90.
121
Не 1(9) 2008
Мы видим, что байесовский подход позволяет сузить доверительный интервал для 90 в 2,6 раза, для 9, — в 3,7 раза и для Ь — в 1,4 раза по сравнению с подходом, основанным на методе максимального правдоподобия.
6. Байесовский прогноз зависимой переменной, основанный на нормальной классической линейной модели множественной регрессии
Мы продолжаем рассматривать нормальную КЛММР
yt =6 о
i=1
,< i)
+ sf, t = 1,2,..., n,
или, в матричной записи, модель (8)-(9) (см. выше), в которой остатки е,• = е(Х,-) нормальны, гомоскедастичны и взаимнонекоррелированы при любом (а не только наблюденном) наборе значений объясняющих переменных.
Введем в рассмотрение, наряду с наблюденными значениями X и У анализируемых переменных X = (1;х(1),х(2),...,х(к')т и у, их прогнозные (на — тактов времени вперед) значения:
X =
1 х n 1 X n
1 х n
-га
Л 2)
Л2)
Лk) \
Л k)
Л k) ^ n+q у
и Y =
' yn+Л
У n+2
V n+q у
<ъ §
Ss ¡S
I
u <ъ
!
<ъ §
0
1
п ta
! С
1 и
8 и <ъ >s <в
Щ
а также соответствующие остатки ~ = (s„+,,sn+2,...,sn+q)T.
Тогда в соответствии с (8)-(9):
Г~ = X0 + Y
[Ye Nq (0; h-1 •Iq).
Для того чтобы строить точечные и интервальные оценки для Y по заданным значениям X,X и Y, очевидно, надо располагать плотностью условного распределения p(Y\X;X;Y), которую обычно называют «прогнозной функцие~ плотности вероятности». Но поскольку из (8а)-(9а) следует, что распределение вектора Y зависит также от параметров 0 и h, а они в байесовском подходе интерпретируются как случайные величины, имеющие соответствующее апостериорноераспределе~ие, то реализуется следующая схема определения прогнозной функции плотности p(Y\X;X;Y):
p(~| X;X; Y) = J Jp(~;0; h| X;X; Y)d0dh = J Jp(~|0; h; X;X; Y) • p(0; h| X;X; Y) d0dh.
(51)
Правая часть (51) получена с использованием формулы произведения условных вероятностей
~ ~ ~ ~ — — — (~ — X©)Т (~ — X©)
Р{Л8 |С) = Р{Л\8, С) • Р(В |С). С учетом того, что р(~|©; Ь; Х;Х; У) = р(~|©; Ь X)~ Ь2 • е 2 ~ , а р(©;Ь|X;X;У) = р(©;Ь|X;У)— гамма-нормальное распределение с параметрами ©0, А0, а
и р, определяемыми по параметрам ©0,Л0, а и р априорного гамма-нормального распределения р(©,Ь) по формулам (47), интегрирование в правой части (51) дает:
0 h
0 h
122
Не 1(9)2008
р(©| X; X; У )•
1 + -(© - X© 0)т В '(У - X©
V+
(52)
Уп-
11 Уп+т - Г1-Ро (п - к -1) ■—=; Уп+т + Г1-Ро (п - к -1)
, т = 1,2,..., д (55)
• байесовская прогнозная доверительная область ДУ для вектора У = (уп+1,...,уп+д)т состоит, с заданной вероятностью Р0, из всехтех© = (уп+,,..., уп+д )т, которые удовлетворя ют неравенству
1
-(У - X©(Б))тЕ©1(© -X©(Б>) < Г1_Р0(д; п - к -1), (56)
д
г
2 —
<в
О
где V = п - к -1 и В ' = - [I д - X(Л0 + Xт X + X X Г^т] (53)
Р
(подробное доказательство этого факта читатель найдет, например, в [Зельнер (19)80)]). Таким образом, мы пришли к тому, что условное распределение д-мерного вектора У при заданных значениях X,У и X описывается обобщенным многомерным Г-распределением с п - к - 1степенями свободы, параметром сдвига X© 0 и матрицей точности В', определенной соотношением (53) (т. е. (©| X; X; У) = Г(п - к -1^© 0;В'), см. Приложение 1б).
Используя известные свойства обобщенного Г-распределения Стьюдента (см. Приложение 1б), получаем следующие байесовские прогнозы для У:
• точечныйбайесовскийпрогноздлякомпонент вектора© определяется соотношением:
у п+т (прогнозное) = (©(Б>)т ■ Хп+т, т = 1,2,...,д; (54)
• интервальный байесовский прогноз для компонент вектора © с вероятностью Р0 определяется соотношением:
где ГЕ (V) — 100е%-ная точка стандартного Г(V©распределения Стьюдента, а величины с'т вычисляются по схеме (49)—(50) с заменой (к х к ©матрицы В на д х д-матрицу В', определенную соотношением (53);
где ГЕ (V,, V2) — 100е%-ная точка Г(V,, V2 )-распределения, (©(Б) — байесовская точечная оцен-
п - к -1 ©
ка параметров регрессии ©,а2© =-(В ) 1 — ковариационная матрица вектора У. Мож-
п - к - 3
но показать, что для модели (8а)-(9а) эта область имеет форму д-мерного эллипсоида.
Рассмотрим реализацию описанной выше схемы построения точечных и интервальных байесовских прогнозов значений зависимой переменной в нормальной КЛММР на нашем примере, проанализированном в задаче 5.
Задача 5 (продолжение).
В условиях примера, рассмотренного выше в задаче 5, требуется: по заданным (планируемым) значениям автономных инвестиций х21 = 120 и х22 = 140 по-строитьточечные и интервальные байесовские прогнозы для среднедушевыхдоходов населения у 21 и у 22, а также прогнозную доверительную область Д У для этих значений с уровнем
V
2
123
Нв 1(9) 2008
доверия Р0 - 0,90. Сравнить полученные решения с решениями, основанными на методе максимального правдоподобия.
Решение. Итак, в нашем случае:
X =Р 12°1; ~ -Н- ?
[1 140^ [у 22)
В соответствии с (52) плотность условного распределения вектора ~ при заданных X, X и Уописывается обобщенным многомерным Г-распределением с числом степеней свободы
(1 120 У349,0Л
V - 20 — 1 — 1 - 18, параметром сдвига соотношением (53).
Ч1 140)
2,9
и матрицей точности В', определенной
Произведя необходимые вычисления по формулам (53)-(56) и используя известные свойства обобщенного многомерного Г-распределения (см. Приложение 1б), имеем:
у 21 (прогн.) - 349 + 2,9 • 120 - 697,4,
у22(прогн.) - 349 + 2,9• 140 - 755,5,
, (0,00159 —0,00022^1 Г721,8 106,8^
В "[—0,00022 0,00152 ); ^ ~ -[ю6,8 757,4^,
у21 е [653,6; 741,2] с вероятностью Р0 - 0,90,
у 22 е [710,7; 800,3] с вероятностью Р0 - 0,90,
\(у 21 ^ 1 (у 21 — 697,4^1т ( 0,00159 —0,00022Уу 21 — 697,4^ 1 ~ -\[у 22J '2 [у 22 — 755,5) [—0,00022 0,00152 )[у 22 — 755,5, |.
Сравним эти результаты с соответствующими прогнозами, основанными на оценках ме-$ тода максимального правдоподобия:
• точечный прогноз § у2мп (прогнозное) - 344,7 + 3,05 • 120 - 710,7,
£ у 2м2п (прогнозное) - 344,7 + 3,05 • 140 - 771,7; |
щ • интервальный прогноз строится на основе Г(п — 2)-распределенности случайных вели-
§ чин
§ ум+пт (прогн.) — уп+т
п ---, т -1,2;
1 п 1 п (х ппт х)
03
! С
I 1 1
^ г, , , , ./-м, , лт п 1-.Л о _ ~>1Л V _ 1 1П „ —1/1П £.2
X (X/ — х)2
5 в нашем случае п - 20, хпп1 - 120, хпп2 - 140, стмп -7—- „ ,_ -666,67; стмп - 25,82;
и
(хпп1 — х )2 5572,^ (хпп 2 — х )2 8958,6
§ , , Ьмп 0,0015
- 0,976;-— --- - 1,569 и Г0,05 (18) - 1,734,так что:
щ ^ —\2 5711 ' -о 5711
X (х> — х )2 X (х> — х)
124
1=1
1 -1
Не 1(9)2008
у21 е [647,1; 774,3] с вероятностью P0 = 0,90, -Е
м
у22 е [699,2; 844,2] с вероятностью P0 = 0,90. >§
Ч
о
Мы видим, что и в прогнозе байесовский подход позволяет сузить ширину прогнозной интервальной оценки для у22 в 1,45 раза, а для у22 — в 1,62 раза!
Приложения
Некоторые сведения об одномерных и многомерных законах распределения вероятностей, используемые в байесовском подходе
Приложение 1а
Обобщенное одномерное распределение Стьюдента с V степенями свободы, параметром сдвига 90 и параметром точности с (^|90;с)-распределение)
Как известно (см., например, [Айвазян, Мхитарян (2001), гл.3]), стандартный з.р.в. Стьюдента с V степенями свободы (ст. св.) описывает распределение случайной величины
г(V) = , 20 , (П.1)
1 1
1Е22
где 20,2ь---,^ — статистически взаимонезависимые (0;а2)-нормально распределенные случайные величины. Значение соответствующей функции плотности вероятности (((|)(х) в точке хзадается соотношением:
г(—1 2
V)(х) = ^ 2 Д —1 2, (П.2)
12 у
V
причем Et(V) = 0 и Dt(V) =-(V > 2).
V- 2
Введем в рассмотрение случайную величину t(VI90; с), являющуюся линейной функцией от а именно:
^^ с) = -1= t(V) +90. (П.3)
Ыс
Легко показать, что функция плотности вероятности ^^с)(х) случайной величины t(VI90; с) в точке х имеет вид:
41 -гМ, 2Л -
1 2 Я, , с(х-90)21 2
0;с) (х) =--— | , (П.4)
2
причем Et(V |90; с) = 90 и Dt(V |90; с) = -•—(V > 2).
с V-2
125
Нв 1(9) 2008
Распределение, задаваемое плотностью (П.4), называют обобщенным распределением Стьюдента (или Г(у|90; с)-распределением) с параметром сдвига 90 и параметром точности с. Отметим, что в данном случае параметр точности с не есть величина, обратная к дисперсии случайной величины Г(у|90; с): дисперсия может и не существовать (при V < 2).
Приложение 1б
Обобщенное к-мерное (к >2) распределение Стьюдента с V степенями свободы, параметром сдвига ©0 = (9° ,92,...,9т и (кх к)-матрицей точности В (или так называемое IСу|©0; В)-распределение)
Стандартный к-мерный з.р.в. Стьюдента с V ст. св. описывает распределение к-мерной случайной величины
г (V) = (г (1)( v),t и( V), ...,г (к)(v))т, (П.5)
где каждая из компонентtj(V) — стандартная стьюдентовская случайная величина (П.1), и все компоненты tj(у)(j = 1,2,...,к) взаимнонекоррелированы. Функция плотности вероятности
Г (V)
^(v)(X) в точке X = (х(1),х(2),...,х(к) )т задается соотношением
^(Х) = \2 Л ■|1+-Хт ■ Х| 2, (П.6)
Г (v ^ к -VI
(™)2 ■г^^ ]
- V
причем ЕГ(V) = 0к и ковариационная матрица 2Г(v) =--1к, где 0к обозначает к-мерный
V-2
вектор-столбец из нулей, а Iк — единичная матрица размерности к.
Обобщенное ^-мерное распределение Стьюдента (или Г(V| ©0;В)-распределение) с V ст. св., параметром сдвига ©0 и матрицей точности В описывает распределение случайной величины
Г( V© 0; В) = С ■ Г (V) +© 0, (П.7)
<ъ §
§ -
§ где С — некоторая невырожденная к х к-матрица, В = (ССт)-1, © 0 = (90,9 2,..., 90), а Г (V) —
I
и <ъ
Ц
стандартная стьюдентовская к-мерная случайная величина (П.5), подчиняющаяся з.р.в. с плотностью (П.6). Значение функции плотности вероятности ^00.В)(Х) случайной величины Г(VI ©0;В) в точке Х = (х(1),х(2),...,х(к) )т задается соотношением
У , ч
I гГ^+к] 1
| 7|© 0; В) (Х) = 1 к 2Д В"2 [/1+""( Х-© 0)т В (Х-© 0)1 2, (П.8)
| ("V)2 г( V) 1 V J
& - V
<5 причем ЕГ (v| © 0; В) = © 0 и ковариационная матрица 2 В) =--В.
; v-2
§ Именно этим з.р.в. описывается сопряженное априорное (а следовательно, и апостери-§ орное) частное распределение вектора © = (90,9,,.,9к)т коэффициентов регрессии в нормальной КЛММР, а также условное (при фиксированных X,У и X) апостериорное рас-иа пределение вектора © = (уп+!,уп+2,...,уп+д )т прогнозных значений зависимой переменной в этой модели (см. в тексте консультации пример 9 и задачу 5).
126
Не 1(9)2008
Г(V!© 0; В) -
(г (1)( V© 0;
I
<в
При построении байесовских интервальных оценок и доверительных областей для параметров ©, также как и при построении байесовских интервальных оценок и доверительных областей для прогнозных значений используются следующие свойства Г©0;В)-рас- ч пределения. ^
Свойство (А). Пусть анализируемая к-мерная случайная величина ГСу| ©0; В) разбита на два подвектора Г(1)Су| © 0; В) и Г(2)^| ©0; В) соответственно размерностей к1 и к2 (к 1 п к2 - к), т. е.
Соответственно этому разобьются на блоки вектор средних значений ©0 и матрица точности В:
(©1(1) ^ © 0 -I ч| и В -0 [©1(2) |
Сформулируем свойство (А).
Частное (маржинальное) распределение вектора Г(1) © 0; В) является камерным обобщенным распределением Стьюдента Г(V© 0 (1); В (1)) с параметром сдвига © 0 (1) - (9°,..., 90, )т и матрицей точности В (1) - В11 — В12 В — В 21.
При построении интервальных оценок нас интересует частное распределение отдельной (/-й) компоненты Г©0;В) анализируемой к-мерной обобщенной стьюдентовской случайной величины ГСу| © 0; В). Соответственно, при этом используется частный случай данного свойства, когда в роли Г(1)^|©0;В) выступает компонента Г(/)^|©0;В).Тогда:
к 1 - 1; ©0(1)-90; В(1) -Ьц — В, • В(/)• В,, (П.9)
где Ьц — /-й диагональный элемент матрицы точности В, В/. - (Ьл,...,Ь./—1,Ь./п1,...,Ь/к) — (к — 1)-мерная строка, В./ - (Ь1 /,...,Ь// ,Ь/п1./,...,Ьк/ )т — (к — 1)-мерный столбец матрицы В(/), а В(/) — это (к — 1) х(к — 1)-матрица, получающаяся из матрицы В вычеркиванием из нее/-й строки и /-го столбца.
Заметим, что в данном случае В(1) — это числовой параметр точности в частном обобщенном одномерном стьюдентовском распределении компоненты Г(/)Су|©0;В).
Свойство (В) используется при построении доверительных областей для неизвестных значений параметров КЛММР или одновременно для нескольких прогнозных значений зависимой переменной и заключается в том, что статистика
у - -(г(V©0; В) —©0)тВ(Г(V©0; В) —©0)
к
асимптотически (по да) подчиняетсяГ(к;v)-распределению.
Поэтому, определяя из таблиц по заданной доверительной вероятности Р0 значение 100(1-Р0)%-ной точки Г1—Р0 (к;V) соответствующего Г-распределения, мы можем с помощью неравенства
127
Hb 1(9) 2008
1
-(г(V©0; В) —©0)тВ(г(V©0;В) —©0) < Г—Р0 (к;V) (П.10)
к
определить к-мерную область, в которую попадает 100Р0 % наблюдений случайной величины Г М © 0; В).
Приложение 2а
Двумерное гамма-нормальное распределение и его свойства
Совместное двумерное распределение случайной величины (9;Ь) называется гамма-нормальным, если его функция плотности вероятности р(9;Ь) задается (с точностью до нормирующего множителя) соотношением
2 X0Ь ^9_9 )2
р(9; Ь)~(X0 Ь)2 е 2 0 • Ьа—1е-рЬ, (П.11)
где Х0,90, а и р — некоторые числовые значения параметров этого семейства распределений.
Свойства двумерного гамма-нормального распределения
(I) Частное распределение параметра 9 есть одномерное обобщенное распределение Стьюдента (П.4) с 2а ст. св., параметром сдвига 90 и параметром точности с - X0 .а/р, т. е.
9- г 2а|9 0; X 0
[ р) 1-
а
Отсюда, в частности, следует, что случайная величина X0 — (9 — 90) подчиняется стандартному распределению Стьюдента с 2а ст. св. V р
(II) Частное распределение параметра Ь есть гамма-распределение (23) с параметрами
аа
(а, р) и, следовательно, ЕЬ - р ОЬ - — и Ь е [у (а;р);уе (а;р)] с вероятностью Р0 - 1 — 2е, где у—(а,р) — это 100—%-ная точка гамма-распределения с параметрами а и р. Отметим, что при а, кратном 0,5, справедлива формула:
(U
Ss , 1
¡g у„ (a,ß) = —%2 (2а), (П.12)
2ß
где %q (m) — это 100д%-ная точка «хи-квадрат»-распределения с m ст. св.
I
и <ъ
Л (III) Условное распределение параметра 9 (при условии заданности значения параметра h,
| т. е. при h - h0, где h0 —заданное число) является (90;1/X0h0 )-нормальным распределени-
| ем (вытекает из (П.11) при подстановке в правую часть этого соотношения заданного значе-
* ния h - h0).
g: Приложение 26 Многомерное (k + 1-мерное, k> 1) гамма-нормальное распределение и его свойства
с
'! Совместное k + 1-мерное распределение параметров 0 - (9,,92,...,9k)T и h называется
g многомерным гамма-нормальным, если его функция плотности вероятности p(0;h) за-
^ дается (с точностью до нормирующего множителя) соотношением:
'S
ГС k 1 h
Щ p(0;h)~ h2\A0\2 e~2<0-00>tЛ0<0-00> .ha~1e-ßh, (П.13)
128
Не 1(9)2008
где заданные численные значения векторного параметра сдвига 00 = (9° ,92 ,•••,90)', элементов (к хк)-матрицы точности Л0, а также параметров а и р однозначно определяют з.р.в. § параметров 0 и Ь.
ч
Свойства многомерного гамма-нормального распределения °
(I) Частное распределение векторного параметра 0 = (9,, 9 2,...,9 к) есть многомерное обобщенное распределение Стьюдента (П.8) с 2аст. св., параметром сдвига 0 = (9? ,92,... ,90 )т и матрицей точности В - —Л0.
(II) Частное распределение скалярного параметра Ь есть гамма-распределение (23) с параметрами а и р.
(III) Условное распределение векторного параметра 0 (при условии заданности значения параметра Ь, т.е. при Ь - Ь0, гдеЬ0 — заданное значение) является к-мерным (00;(Ь0Л0)- )-нормальным распределением.
Приложение 3
Некоторые сведения об априорных з.р.в., сопряженных по отношению к наблюдаемым генеральным совокупностям, зависящим от единственного неизвестного параметра
№ пп З.р.в. наблюдаемой генеральной совокупности Сопряженный априорный з.р.в.p(6), выражения для E6 иР6 Апостериорный з.р.в. p (6| х,, x2,..., xn), выражения для его параметров
1 (9; ст 2)-нормальный, , (х-9)2 /Г(X|9) - ' ^ 2ст2 л/2яст (значение дисперсии ст2 известно) (90; ст 2)-нормальный; Е9-9с; 09-ст 2 (90 и ст0 — заданы) (00;ст'02)-нормальный, где 0„- X + у00 и ст02 -а2/n(1 + у), 1 + у а у-а2/ па0
2 Экспоненциальный |9е~9х при х > 0 ^ (х |9) [0 при х < 0 р(9) - ^ 9а-1е-р9 (9>0) Г( а) гамма-распределение; Е9-а / р; 09-а / р2 (а и р — заданы) Гамма-распределение с параметрами а' - а + n, Р'-Р+£х, i-1
3 [0; 9]-равномерный: Г(х|9) -(9 для0 *х [0 для х й [0;9] |а9а . „ р(9)-Ь^ при х >0;(9>0) — (0 при х < 0 распределение Парето; Е9- а90 ; 09- а92 а- 1 (а- 1)2( а- 2) (а и 90 — заданы) Распределение Парето с параметрами а' - а + n, 00 -max{0o;Xi,х2,...,Xnl
4 Распределение Пуассона: 9х Р х} - — е-9 х! х - 0,1,2,... р(9) - ^ 90-1е-р9 (9>0) Г( а) гамма-распределение; Е9-а; 09 - — р р2 (а и р — заданы) Гамма-распределение с параметрами n а'-а + ^ х,-, i -1 Р'-Р + n
129
Не 1(9) 2008
Окончание
№ пп З.р.в. наблюдаемой генеральной совокупности Сопряженный априорный з.р.в.p(6), выражения для E6 иР6 Апостериорный з.р.в. p (6| x,, x2,..., xn), выражения для его параметров
5 Биномиальное распределение: Р{х} - СН 9х(1 — 9)Н—х (значение параметра N известно) р(9) -г(a+b) ea-1(i e)b-1 г( a )Г(Ь) (0 < e < 1) — бета-распределение; Ee- a ; a + b De- ab (a + b )2( a + b + 1) (a и b — заданы) Бета-распределение с параметрами п аа пХ х/, /-1 п Ь'-Ь п пН — Xх-, /-1
6 Отрицательное биномиальное распределение Р {х} - Сх—19к (1—9)х — к (значение параметра к известно) х - к, к п 1,. -г(a+b) ea-i(i e)b-1 Г( a )r(b) (0 <e< 1) — бета-распределение; Ee- a ; a + b De- ab (a + b )2( a + b + 1) (a и b — заданы) Бета-распределение с параметрами аа п кп, п Ь' -Ь пХх, — кп /-1
7 Распределение Парето \9х 9 . ах|9) -Ы при х * х0 \0 при х < х0 (значение параметра х0 — известно) p(e) - p" ea-1e-pe (e>0) Г( a) гамма-распределение; Ee-a; De- — p p2 (a и p — заданы) Гамма-распределение с параметрами а'-ап п, р'-рп п 1п[д-где дп х, )
щ
Ss iB
I
u <ъ
J
!
<ъ §
0
1
n <a
! С
1 u
8 u <u >s
<S U5
Список литературы
Айвазян С. А. (2001). Прикладная статистика и основы эконометрики. Том 2: Основы эконометрики. Издание 2-е. Юнити. §§ 2.1-2.3.
Айвазян С. А., Мхитарян В. С. (2001а). Прикладная статистика и основы эконометрики. Том 1: Теория вероятностей и прикладная статистика. Издание 2-е. Юнити. § 7.6.
Айвазян С. А., Мхитарян В. С. (20016). Прикладная статистика в задачах и упражнениях. М.: Юнити.
Де Гроот М. (1974). Оптимальные статистические решения. Пер. с англ. М.: Мир. Гл. 4, 5, 9 и §§11.10-11.12.
ЗельнерА. (1980). Байесовские методы в эконометрике. Пер. с англ. М.: Статистика. Главы 1-3. Ghosh J. K., Delampady M., Samanta T. (2006). An Introduction to Bayesian Analysis. Theory and Methods. Springer.
Jeffreys H. (1957). Scientific Inference. 2nd ed. Cambridge University Press. Lancaster A. (2004). An Introduction to Modern Bayesian Econometrics. Blackwell Publ.
130