Задание 3. Построение оценок максимального правдоподобия и их доверительной области для параметров модели Гомпертца распределения длительности события

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

,

где a – «начальная смертность», b – «темп старения». Плотность вероятности имеет вид

.

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

, (3.7)

где , если наблюдение нецензурировано, , если наблюдение цензурировано. Необходимым условием максимума выражения (3.7) является равенство нулю первых частных производных LnL по a и по b

,

откуда получаем, что оценка максимального правдоподобия параметра a равна

, (3.5)

m – число нецензурированных наблюдений (случаев).

Для получения оценки максимального правдоподобия параметра b необходимо решить другое уравнение. Запишем

.

Подставив выражение для оценки , получим уравнение для определения оценки b, которое решается численными методами

(3.6).

Подставляя найденное значение b в (3.5) получим оценку .

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

Для построения доверительной области значений параметров модели Гомпертца используется метод построения поверхности логарифма правдоподобия. Этот метод основан на том, что при большом числе наблюдений случайная величина имеет распределение с двумя степенями свободы. Из свойств распределения следует, что с вероятностью 0.95 справедливо неравенство .

Вычисления проводятся в следующем порядке:

1. строится поверхность функции ;

2. максимальное значение поверхности достигается при ;

3. вычисляется критическое значение и определяется область, задаваемая неравенством ;

найденная область является 95% доверительной областью для параметров модели Гомпертца.

3.1 Получите у преподавателя значения параметров a, b и n. Запишите оператор для генерации независимой выборки значений случайной величины, имеющей распределение, описываемое моделью Гомпертца с параметрами a и b. Для этого сначала генерируется независимая выборка из n чисел, имеющих экспоненциальное распределение со средним значением 1. Элементы выборки вычисляются по формуле .

z=exprnd(1,1,n);

x=log(b/a*z+1)/b;

3.2 Считая, что цензурирования нет, построить по формуле (3.7) поверхность логарифма правдоподобия,

% значения переменных a1, a2, b1 и b2 могут зависеть от сгенерированных данных и параметров модели Гомпертца

a1=0.00001;a2=1.5; b1=1; b2=12;

at=a1:0.01:a2;

bt=b1:0.1:b2;

LnL=zeros(length(at),length(bt));

for i=1:length(at)

for j=1:length(bt)

LnL(i,j)=n*log(at(i))+bt(j)*sum(x)-at(i)/bt(j)*(sum(exp(bt(j)*x))-n);

end

end

3.3 Вычислить критическое значение t, нарисовать область, ограничивающую множество решений неравенства

maxLnL=max(LnL(:));

t=maxLnL-3.0;

contour(a,b,LnL',[lev lev]);

xlabel(‘a’); ylabel(‘b’);

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

ah=a(max(LnL,[],2)==maxLnL); % оценка a

bh=b(max(LnL,[],1)==maxLnL); % оценка b

disp([‘a=’,num2str(ah),’ b=’,num2str(bh)])

3.5 Нанести на плоскости, построенной в пункте 3.3, символ * зелёного цвета с координатами a и b, а также символ О красного цвета с координатами ah и bh

hold on;

plot(a,b,’*g’);

plot(ah,bh,’Or’);

3.6 С помощью команды run Lab3 выполните записанные операторы. Значения параметров, построенный график и выведенные на экран числа внесите в отчёт о выполнении лабораторной работы. Сделайте вывод о близости оценок параметров модели Гомпертца, построенных методом максимального правдоподобия, к их истинным значениям.

Контрольные вопросы.

1. Чем отличаются модели Вейбулла и Гомпертца?
2. По каким формулам вычисляется логарифм правдоподобия для экспоненциального распределения, модели Вейбулла и модели Гомпертца?
3. Что такое цензурирование справа и как оно учитывается при построении правдоподобия?
4. По каким формулам вычисляются оценки параметров модели в методе максимального правдоподобия?
5. В чём заключается сущность метода определения 95% доверительной области через поверхность логарифма правоподобия?
6. Как в MATLAB построить 3D поверхность?
7. Как в MATLAB нарисовать контурное представление поверхности?
8. Опишите результат выполнния оператора at=at(LnL>=t). Здесь at и LnL – векторы одинаковой длины, t – число.