Обобщенный линейный образцовый рабочий процесс

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

Шаг 1. Загрузите данные.

Загрузите ирисовые данные Фишера. Извлеките строки, которые имеют классификацию versicolor или virginica. Это строки 51 - 150. Создайте логические переменные отклика, которые являются true для versicolor цветов.

load fisheriris
X = meas(51:end,:); % versicolor and virginica
y = strcmp('versicolor',species(51:end));

Шаг 2. Соответствуйте обобщенной линейной модели.

Соответствуйте бином обобщил линейную модель к данным.

mdl = fitglm(X,y,'linear',...
    'distr','binomial')

mdl = 
Generalized linear regression model:
    logit(y) ~ 1 + x1 + x2 + x3 + x4
    Distribution = Binomial

Estimated Coefficients:
                   Estimate      SE       tStat      pValue 
                   ________    ______    _______    ________

    (Intercept)     42.638     25.708     1.6586    0.097204
    x1              2.4652     2.3943     1.0296     0.30319
    x2              6.6809     4.4796     1.4914     0.13585
    x3             -9.4294     4.7372    -1.9905    0.046537
    x4             -18.286     9.7426    -1.8769    0.060529


100 observations, 95 error degrees of freedom
Dispersion: 1
Chi^2-statistic vs. constant model: 127, p-value = 1.95e-26

Шаг 3. Исследуйте результат, рассмотрите альтернативные модели.

Некоторые $p$ - значения в столбце pValue не являются очень маленькими. Возможно, модель может быть упрощена.

Смотрите, включают ли приблизительно 95% доверительных интервалов для коэффициентов 0. Если так, возможно, эти образцовые условия могли быть удалены.

confint = coefCI(mdl)

confint = 5×2

   -8.3984   93.6740
   -2.2881    7.2185
   -2.2122   15.5739
  -18.8339   -0.0248
  -37.6277    1.0554

Только два из предикторов имеют коэффициенты, доверительные интервалы которых не включают 0.

Коэффициенты 'x1' и 'x2' имеют самое большое $p$ Значения. Протестируйте, могли ли оба коэффициента быть нулем.

M = [0 1 0 0 0     % picks out coefficient for column 1
     0 0 1 0 0];   % picks out coefficient for column 2
p = coefTest(mdl,M)

p = 0.1442

$p$ - значение приблизительно 0,14 не является очень маленьким. Исключите те условия из модели.

mdl1 = removeTerms(mdl,'x1 + x2')

mdl1 = 
Generalized linear regression model:
    logit(y) ~ 1 + x3 + x4
    Distribution = Binomial

Estimated Coefficients:
                   Estimate      SE       tStat       pValue  
                   ________    ______    _______    __________

    (Intercept)     45.272     13.612      3.326    0.00088103
    x3             -5.7545     2.3059    -2.4956      0.012576
    x4             -10.447     3.7557    -2.7816     0.0054092


100 observations, 97 error degrees of freedom
Dispersion: 1
Chi^2-statistic vs. constant model: 118, p-value = 2.3e-26

Возможно, было бы лучше иметь stepwiseglm, идентифицируют модель первоначально.

mdl2 = stepwiseglm(X,y,...
    'constant','Distribution','binomial','upper','linear')

1. Adding x4, Deviance = 33.4208, Chi2Stat = 105.2086, PValue = 1.099298e-24
2. Adding x3, Deviance = 20.5635, Chi2Stat = 12.8573, PValue = 0.000336166
3. Adding x2, Deviance = 13.2658, Chi2Stat = 7.29767, PValue = 0.00690441

mdl2 = 
Generalized linear regression model:
    logit(y) ~ 1 + x2 + x3 + x4
    Distribution = Binomial

Estimated Coefficients:
                   Estimate      SE       tStat      pValue 
                   ________    ______    _______    ________

    (Intercept)     50.527     23.995     2.1057    0.035227
    x2              8.3761     4.7612     1.7592    0.078536
    x3             -7.8745     3.8407    -2.0503    0.040334
    x4              -21.43     10.707    -2.0014     0.04535


100 observations, 96 error degrees of freedom
Dispersion: 1
Chi^2-statistic vs. constant model: 125, p-value = 5.4e-27

stepwiseglm включал 'x2' в модель, потому что это не добавляет и не удаляет условия с $p$ - значения между 0,05 и 0.10.

Шаг 4. Ищите выбросы и исключите их.

Исследуйте график рычагов искать влиятельные выбросы.

plotDiagnostics(mdl2,'leverage')

Существует одно наблюдение с рычагами близко к одному. Используя Data Cursor, кликните по точке и найдите, что это имеет индекс 69.

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

oldCoeffs = mdl2.Coefficients.Estimate;
mdl3 = fitglm(X,y,'linear',...
    'distr','binomial','pred',2:4,'exclude',69);
newCoeffs = mdl3.Coefficients.Estimate;
disp([oldCoeffs newCoeffs])

   50.5268   50.5268
    8.3761    8.3761
   -7.8745   -7.8745
  -21.4296  -21.4296

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

Шаг 5. Предскажите вероятность, что новый цветок является versicolor.

Используйте mdl2, чтобы предсказать вероятность, что цветок со средними измерениями является versicolor. Сгенерируйте доверительные интервалы для своего прогноза.

[newf,newc] = predict(mdl2,mean(X))

newf = 0.5086

newc = 1×2

    0.1863    0.8239

Модель дает почти 50%-ю вероятность, что средний цветок является versicolor с широким доверительным интервалом об этой оценке.

Документация