Создайте 1000 нормальных случайных чисел из нормального распределения со средним значением 3 и стандартным отклонением 5.

rng('default') % For reproducibility
x = normrnd(3,5,[1000,1]);

Найдите оценки параметров и 99% доверительные интервалы.

[muHat,sigmaHat,muCI,sigmaCI] = normfit(x,0.01)

muHat = 2.8368

sigmaHat = 4.9948

muCI = 2×1

    2.4292
    3.2445

sigmaCI = 2×1

    4.7218
    5.2989

muHat - среднее значение выборки, и sigmaHat - квадратный корень несмещенного оценщика дисперсии. muCI и sigmaCI содержат 99% доверительные интервалы параметров среднего и стандартного отклонения соответственно. Первая строка является нижней границей, а вторая - верхней границей.

Поиск MLE набора данных с цензурой с помощью normfit. Использовать statset чтобы указать параметры итеративного алгоритма, которые normfit используется для вычисления MLE для данных, подвергнутых цензуре, и последующего поиска MLE.

Загрузите образцы данных.

load lightbulb

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

Найдите индексы флуоресцентных лампочек.

idx = find(lightbulb(:,2) == 0);

Предположим, что время жизни соответствует нормальному распределению, и найдите MLE параметров нормального распределения. Второй входной аргумент normfit определяет уровень достоверности. Пройти внутрь [] для использования значения по умолчанию 0,05. Третий входной аргумент указывает цензурную информацию.

censoring = lightbulb(idx,3) == 1;
[muHat1,sigmaHat1] = normfit(lightbulb(idx,1),[],censoring)

muHat1 = 9.4966e+03

sigmaHat1 = 3.0640e+03

Отображение параметров алгоритма по умолчанию, normfit используется для оценки параметров нормального распределения.

statset('normfit')

ans = struct with fields:
          Display: 'off'
      MaxFunEvals: 200
          MaxIter: 100
           TolBnd: 1.0000e-06
           TolFun: 1.0000e-08
       TolTypeFun: []
             TolX: 1.0000e-08
         TolTypeX: []
          GradObj: []
         Jacobian: []
        DerivStep: []
      FunValCheck: []
           Robust: []
     RobustWgtFun: []
           WgtFun: []
             Tune: []
      UseParallel: []
    UseSubstreams: []
          Streams: {}
        OutputFcn: []

Сохраните параметры с другим именем. Изменение способа отображения результатов (Display) и допуск окончания для целевой функции (TolFun).

options = statset('normfit');
options.Display = 'final';
options.TolFun = 1e-10;

Кроме того, можно задать параметры алгоритма, используя аргументы пары имя-значение функции statset.

options = statset('Display','final','TolFun',1e-10);

Найдите MLE с новыми параметрами алгоритма.

[muHat2,sigmaHat2] = normfit(lightbulb(idx,1),[],censoring,[],options)

Successful convergence: Norm of gradient less than OPTIONS.TolFun

muHat2 = 9.4966e+03

sigmaHat2 = 3.0640e+03

normfit отображает отчет о конечной итерации.

Функция normfit находит среднее значение выборки и квадратный корень несмещенного оценщика дисперсии без цензуры. Среднее значение выборки равно MLE среднего параметра, но квадратный корень несмещенного оценщика дисперсии не равен MLE параметра стандартного отклонения.

Поиск параметров нормального распределения с помощью normfit, преобразовать их в MLE, а затем сравнить отрицательный log вероятностей оценок с помощью normlike.

Создайте 100 нормальных случайных чисел из стандартного нормального распределения.

rng('default') % For reproducibility
n = 100;
x = normrnd(0,1,[n,1]);

Найдите среднее значение выборки и квадратный корень несмещенного оценщика дисперсии.

[muHat,sigmaHat] = normfit(x)

muHat = 0.1231

sigmaHat = 1.1624

Преобразуйте квадратный корень несмещенного оценщика дисперсии в MLE параметра стандартного отклонения.

sigmaHat_MLE = sqrt((n-1)/n)*sigmaHat

sigmaHat_MLE = 1.1566

Разница между sigmaHat и sigmaHat_MLE ничтожно мало для больших n.

Кроме того, можно найти MLE с помощью функции. mle.

phat = mle(x)

phat = 1×2

    0.1231    1.1566

phat(1) и phat(2) - MLE среднего значения и параметра стандартного отклонения соответственно.

Подтвердите, что вероятность регистрации MLE (muHat и sigmaHat_MLE) больше логарифмической вероятности непредвзятых оценщиков (muHat и sigmaHat) с помощью normlike функция.

logL = -normlike([muHat,sigmaHat],x)

logL = -156.4424

logL_MLE = -normlike([muHat,sigmaHat_MLE],x)

logL_MLE = -156.4399