主要内容

广义线性模型工作流

这个例子说明了如何拟合一个广义线性模型并分析结果。典型的工作流包括导入数据、拟合广义线性模型、测试模型质量、修改模型以提高模型质量以及基于模型进行预测等步骤。在这个例子中,您使用Fisher虹膜数据来计算花在两个类之一中的概率。

加载数据

载入费雪虹膜数据。

负载fisheriris

提取第51至150行,其中分类为versicolor或virginica。

X =量(51:最终,);

创建逻辑响应变量真正的杂色的,virginica。

y = strcmp (“多色的”、物种(51:结束));

拟合广义线性模型

对数据拟合一个二项广义线性模型。

mdl = fitglm (X, y,“线性”“分布”“二”
mdl =广义线性回归模型:logit(y) ~ 1 + x1 + x2 + x3 + x4Estimate 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个观测值,95个误差自由度127, p-value = 1.95e-26

根据模型显示,有些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

只有第四个预测者x3具有置信区间不包含0的系数。

的系数x1x2有大量p-值及其95%置信区间包括0。测试两个系数是否都为零。指定一个假设矩阵来选择的系数x1x2

M = [0 1 0 0 0 0];p = coefTest (mdl)
p = 0.1442

p-value大约是0.14,这个值不小。删除x1x2从模型。

mdl1 = removeTerms (mdl,“x1 + x2”
mdl1 =广义线性回归模型:logit(y) ~ 1 + x3 + x4分布=二项估计系数:估计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个观测值,97个误差自由度1 Chi^2-statistic vs. constant model: 118, p-value = 2.3e-26

或者,您可以使用stepwiseglm

mdl2 = stepwiseglm (X, y,“不变”“分布”“二”“上”“线性”
1.加x4, Deviance = 33.4208, Chi2Stat = 105.2086, PValue = 1.099298e-24加上x3, Deviance = 20.5635, Chi2Stat = 12.8573, PValue = 0.000336166加上x2, Deviance = 13.2658, Chi2Stat = 7.29767, PValue = 0.00690441
mdl2 =广义线性回归模型:logit(y) ~ 1 + x2 + x3 + x4估计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个观测值,96个误差自由度

p值(pValue)x2在系数表中大于0.05,但是stepwiseglm包括x2在模型中,因为p值(PValue)添加x2小于0.05。的stepwiseglm函数计算PValue使用适合的和没有x2,而函数计算pValue基于仅从最终模型计算出的近似标准误差。因此,PValuepValue

识别异常值

检查杠杆图以寻找有影响力的异常值。

plotDiagnostics (mdl2“杠杆”

图中包含一个轴对象。具有标题的axis对象杠杆的Case order plot包含两个类型为line的对象。这些对象表示杠杆、参考线。

如果一个观察值的杠杆大大超过,那么它可以被认为是一个离群值p / n,在那里p系数的个数是和吗n为观察次数。虚线参考线是建议阈值,由2 * p / n,对应于本图中的0.08。一些观察的杠杆值大于10 * p / n(即0.40)。确定这些观察点。

idxOutliers =找到(mdl2.Diagnostics。利用> 10 * mdl2.NumCoefficients / mdl2.NumObservations)
idxOutliers =4×119 21 57 85

看看当你拟合一个排除这些点的模型时,模型系数是否会改变。

oldCoeffs = mdl2.Coefficients.Estimate;mdl3 = fitglm (X, y,“线性”“分布”“二”...“PredictorVars”2:4,“排除”, idxOutliers);newCoeffs = mdl3.Coefficients.Estimate;disp ([oldCoeffs newCoeffs])
50.5268 44.0085 8.3761 5.6361 -7.8745 -6.1145 -21.4296 -18.1236

的模型系数mdl3不同于mdl2。这一结果表明,高杠杆点的响应与简化模型的预测值不一致。

预测花斑的概率

使用mdl3用平均尺寸来预测一朵花是花斑的概率。生成预测的置信区间。

[newf, newc] =预测(mdl3,意味着(X))
newf = 0.4558
newc =1×20.1234 - 0.8329

该模型给出了几乎46%的概率,平均花是花斑的,置信区间很宽。

另请参阅

||||||

相关的话题