Загрузите выборочные данные.

load('lightbulb.mat');

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

Соответствуйте Cox пропорциональная модель опасностей в течение времени жизни лампочек, также объясняя цензурирование. Переменный предиктор является типом лампы.

b = coxphfit(lightbulb(:,2),lightbulb(:,1), ...
'Censoring',lightbulb(:,3))

b = 4.7262

Оценка отношения опасности $$e^b$$ = 112.8646. Это означает, что опасность для ламп накаливания является 112.86 раз опасностью для флуоресцентных ламп.

Загрузите выборочные данные.

load('lightbulb.mat');

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

Соответствуйте Cox пропорциональная модель опасностей, также объясняя цензурирование. Переменный предиктор является типом лампы.

b = coxphfit(lightbulb(:,2),lightbulb(:,1),...
'Censoring',lightbulb(:,3))

b = 4.7262

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

statset('coxphfit')

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

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

coxphopt = statset('coxphfit');
coxphopt.Display = 'final';
coxphopt.MaxIter = 50;

Запустите coxphfit с новыми параметрами алгоритма.

b = coxphfit(lightbulb(:,2),lightbulb(:,1),...
'Censoring',lightbulb(:,3),'Options',coxphopt)

Successful convergence: Norm of gradient less than OPTIONS.TolFun

b = 4.7262

coxphfit отображает отчет относительно итоговой итерации. Изменение максимального количества итераций не влияло на содействующую оценку.

Сгенерируйте данные Weibull в зависимости от предиктора X.

rng('default') % for reproducibility
X = 4*rand(100,1);
A = 50*exp(-0.5*X); 
B = 2;
y = wblrnd(A,B);

Значения отклика сгенерированы от распределения Weibull параметром формы в зависимости от переменного предиктора X и масштабный коэффициент 2.

Соответствуйте Cox пропорциональная модель опасностей.

[b,logL,H,stats] = coxphfit(X,y);
[b logL]

ans = 1×2

    0.9409 -331.1479

Содействующая оценка 0.9409, и логарифмическое значение вероятности –331.1479.

Запросите статистику модели.

stats

stats = struct with fields:
                    covb: 0.0158
                    beta: 0.9409
                      se: 0.1256
                       z: 7.4889
                       p: 6.9462e-14
                   csres: [100x1 double]
                  devres: [100x1 double]
                 martres: [100x1 double]
                  schres: [100x1 double]
                 sschres: [100x1 double]
                  scores: [100x1 double]
                 sscores: [100x1 double]
    LikelihoodRatioTestP: 6.6613e-16

Ковариационная матрица содействующих оценок, covb, содержит только одно значение, которое равно отклонению содействующей оценки в этом примере. Содействующая оценка, beta, совпадает с b и равно 0,9409. Стандартная погрешность содействующей оценки, se, 0.1256, который является квадратным корнем из отклонения 0.0158. $$z$$- статистическая величина, z, beta/se = 0.9409/0.1256 = 7.4880. P-значение, p, указывает что эффект X является значительным.

Постройте оценку Cox базовой функции оставшегося в живых вместе с известной функцией Weibull.

stairs(H(:,1),exp(-H(:,2)),'LineWidth',2)
xx = linspace(0,100);
line(xx,1-wblcdf(xx,50*exp(-0.5*mean(X)),B),'color','r','LineWidth',2)
xlim([0,50])
legend('Estimated Survivor Function','Weibull Survivor Function')

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