非参数统计翻译23、24章.doc

STA333 第23讲 核密度估计23.1 存在的问题在统计学入门之初，你所运用的代表性的去描绘一组值的分布的工具就是直方图。例如：在3.2部分你看到了如何在uadata资料组获得大学新生的平均分数的直方图: site uadata hist(uadata$gpa.endyr1, main= , xlab=Freshman GPA)虽然这在很多情况下都足够用，但直方图实际上充满了潜在的问题。本节课的目的就是从什么是直方图这一基本思想入手，去寻找直方图想要做什么，然后采用一些“修饰”去改善它。23.2 估计一个密度直方图是一种估计概率密度或者叫做概率分布，并通过运用一个选择的样本来估计变量值群体的方法。还记得在统计学基础中，一个概率分布下的总面积等于1。因为面积代表着问题中变量所产生的概率。为了构造一个直方图，我们把所有数据覆盖的总区间分成若干相等的子区间，也就是我们所说的小方块。每次，一个数据都会落入特定的子区间，然后，一个宽度为1的方块会被放置在它上面。当我们构造一个直方图时，我们需要考虑以下两点：小方块的大小以及小方块的末端点。例: 翼展. 以下数据是源于建立在1956至1984年间的飞机的翼展日志，（完整的数据可在Bowman & Azzalini (1997) Applied Smoothing Techniques for Data Analysis中找到）我们用了这些数据的一部分，即就是下面这些观察值：2 22 42 62 82 102 122 142 162 182 202 222我只用了一部分数据用作说明，否则图形会变得很拥挤。在下面数据点用与x轴的交叉来描绘。下面是相同的样本用相同密度的两种表现。如果我们选择0作为开始断点并且令小方块宽度为0.5，我们的直方图就是左边的那个。看起来这个分布是单峰并且向右歪斜的。但是如果我们用相同的宽度但是将开始断点调至0.25，我们的直方图就是右边的那个。现在我们就得到了一个完全不同的密度估计：看起来它是双峰的(有两个峰)！ 这个例子阐明了用传统直方图可能会遇到的三个问题。即：直方图不光滑。在第一个例子中，平均分数是用一个连续的量度来刻画，但是直方图并没有绝对显示出连续。一个直方图更多的是连续分布的一种不连续“类属”表现。所以，它舍弃了某些东西。直方图依赖于小方块的边界点。在上面的翼展例子中能看到这个问题。直方图依赖于小方块的宽度。直方图中显示的细节水平取决于小方块边界所留的间隙有多么宽。23.3 核密度估计(KDE)我们可以通过运用核密度估计量来减轻上述的前两个问题。为了消除小方块中末点的的依赖性，我们选择集中数据点上的方块，而不是修改方块的末端点。以下是利用翼展案例对这两种的差别进行描述。在左边的直方图中，我们放置了一个宽度为1，高度为1/12的方块（星罗密布的盒子），因为共有12个数据点，我们将它们相加。该密度估计（实线的那个）比我们之前看到的任何一个块头都小，因为我们正开始提取数据中的一些好的结构。它显示：密度是双峰的。这就是我们所说的箱核密度估计。它仍然是不连续的，因为我们用了一个不连续的核作为我们的基础。问题：这种核密度估计下的总面积是什么？光滑的核如果我们用一些“光滑的”核作为我们的基础，我们将会得到一个光滑的核密度估计。从而也能消除直方图的第一个问题。 一个光滑的核函数是满足下面数学性质的一些函数：如果你不懂微积分学，其实它的目的就是让一个核区域下的面积为1。（我保证，我在竭力使复杂的数学难度最小化）通常，K(x) 是一个对称的，单峰的概率密度分布函数。最常见的选择实际上是用正态分布作为一个核。如果你这样做了，那么核密度估计程序可以被刻画如下：Thanks to Ricardo Gutierrez-Osuna of Wright State University for this picture.小的正态核被加起来去创造密度估计，就像是箱核密度估计。但是光滑的核是如何被利用去建立一个样本分布的平滑表现呢？这大概反映了实际的总体分布是连续的这一事实。例:大学进入许可。以下是如何利用R从uadata数据库中建立一个GPA数据的核密度估计图。将改图与本节课第一页的直方图作对比。 plot(density(uadata$gpa.endyr1)2223.4 带宽的选择 (AKA: 总像那么回事儿)所以，我们已经击破了传统直方图面临的三个问题中的两个。但是，第三个“小方块宽度”的问题呢？ 核密度估计中的这个带宽问题近似于传统直方图中选择小方块的宽度。实际上，在KDE问题中带宽的选择有许多潜在的主观因素。核函数K(x)的选择对最后的结果实际上影响很小。（这也是常量经常被用作核的原因）但是应该选择宽度为多少的核作为核呢？这是一个关键的问题，因为：用过宽的带宽将会过度平滑KDE，导致分布的重要特征被模糊化。用过窄的带宽将会不足平滑KDE，导致过多的捕获细节从而过分强调分布中的一些微不足道或者随机的人工误差。以下这四个图完美的总结了这些问题。这儿有相同数据的四个不同核密度估计，都运用了高斯核函数，只是带宽不同：右上角的KDE看起来是这儿最佳的图。但是，如何选择一个最佳的带宽呢？大多数时候，它取决于你如何定义“最佳”。很明显，所选择的最佳带宽是那个能绘制出尽可能接近真实分布模型的KDE的带宽。然而，它需要你首先知道真正的分布到底是什么嗯，这确实是个问题。但是，有一些经验法则很有用。其中被广泛使用的一个是去选择一个带宽，使其值为:其中 s 为数据样本的标准差。有时也建议，与其使用 s , 也可以用小s 和IQR/1.34.来做运算。如果你的数据有离群值或者一些很奇怪的东西，这种方法会更好。例:大学进入许可。R中系统默认的密度点（参考之前的例子）使用的带宽为0.1538. 你可以通过“bw=”选项来修改这一默认值。 plot(density(uadata$gpa.endyr1, bw=0.08)这个带宽（我自己果断选择的，没什么科学依据啦）显然对数据的平滑有一点不足。22第23讲 练习题一个关于300个顾客周内正午时分在市中心Columbus Chipotle餐厅的等候时间（以分钟来记）的随机样本被记录在此。数据在R数据库waittime.Rdata中可寻。1、利用数据和概率的频率化解释去估计P(X 1), 也就是估计随机选出的一个周内正午时分来就餐的顾客需要等待超过一分钟的概率。2、为周内正午时分在该市中心Columbus Chipotle餐厅就餐的等候时间的概率分布寻找一个核密度估计，用系统默认的带宽。3. 重复问题1，但是用问题2的KDE去估计P(X 1)。STA333 第24讲LOESS非参数回归24.1 推广简单的回归思想简单线性回归是用来模拟预测变量 X 和反应变量 Y 之间关系的一种广泛使用的工具，它可以被用来：1、评价 X 和 Y 之间线性关系的显著性和强度。2、当给定一个假设的X 值时预测Y 的未来平均值。在一开始（当然了，也是在复习），我们给出一个R中简单线性回归的一个简单例子。.例：肌肉质量。一个人的肌肉质量会随着年龄的增长而下降。为了探索这一关系在女性中的体现，一个营养学家在年龄从40岁至79岁，每10岁一组，共5组的的人中随机抽取出15个女性。数据可在R的musclemass中寻得。变量数据分别是肌肉质量和年龄。我们来做以下3点：1、绘制一个散点图形象化的研究这种关系。2、拟合一个回归模型将肌肉质量与年龄联系起来。3、预测65岁女性的平均肌肉质量，然后找到该预测置信度为95%的置信区间。.解：首先，绘制出散点图： plot(massage, data=musclemass, main=Muscle Mass vs Age)2、散点图就在左侧。（我们将在最后一节见到这个图）我们可以看到，随着年龄增长，肌肉质量趋于下降。应该同样注意到的是这种下降是一种粗略的线性样式，因此，一个简单的线性回归在这儿看起来似乎是可行的。我们现在在R中拟合这个模型，lm()函数需要被用到，随后summary()函数可以看到回归拟合的结果。 fit summary(fit)Coefficients:Estimate Std. Error t value Pr(|t|)(Intercept) 156.3466 5.5123 28.36 2e-16 *age -1.1900 0.0902 -13.19 predict(fit, newdata=data.frame(age=65), int=conf)fit lwr upr1, 78.99686 76.69872 81.295我们可以以95%的置信度确信，65岁女性的真正平均肌肉质量是在76.7和81.3之间。同时，线性回归线也在拟合完模型后通过下面的命令覆盖在散点图上。 abline(fit)简单线性回归（如果你有许多预测变量的话，叫它多重回归也OK）是数据分析中非常有用的一个工具：它提供了一个百宝盒，帮助我们在分析X 与Y的关系时提取有用的信息。所以，有什么问题吗？一切看起来似乎都很好。噢，等等还记得之前的回归假定吗？为了防止你忘记（或者其实就没记得过），我把它们列举如下：误差假定我们假设误差有三个性质：(1)它们是正态分布的(2) 它们的方差为常数(3)它们是独立的。用符号表示就是：i iid N(0, 2).线性假定我们假定对线性回归模型的假定是正确的。异常观察值偶尔，一些观察值并不符合这一模型。这些独立的观察值有潜力去戏剧性的改变我们的拟合结果，甚至改变我们对模型的选择。看起来有好多要求啊，但是任何好的统计学家都会在相信标准回归分析的结果之前检验这些假定以确保它们是可信的。虽然在违反假定时我们可以运用一些高科技（比如：数据转化），偶尔你会发现，运用简单线性回归在有些时候即使是想一下都“差太远”了。换而言之，简单线性回归在有些数据面前“太简单”啦。本节课我们处理的主要问题是：违反线性假定的情况。当你拟合一个简单线性回归时，你给你的X/Y之间的联系强加了一个线性结构。万一这种联系不是线性的呢？违反正态和常数方差假定的情况。如果误差（残差）不是正态分布或者我们看到一个图显示方差并不是常数，那么我们在之前几页所找到的置信区间就不正确了，直白的说，就是没用的。 给简单线性回归输入一个非参程序。有时，这也被叫做散点图平滑化。也叫做LOWESS，它是LOcally WEighted Scatterplot Smoothing（局部加权回归法）的简称。有时也被叫做LOESS （没有 “W”）。别问我为什么。R有两套程序(分别是lowess 和loess) 都适用于这个模型。Loess比较新，是旧版本lowess的一个改进，所以我们将使用loess。24.2 LOESS基础LOESS在20世纪70年代晚期才发展起来，它是众多建立在传统线性回归基础上的现代建模法之一。现代回归法被用来在传统参数程序拟合的不好或者必须经过很大的努力才能有效利用时发挥作用。LOESS将线性回归的简易性与非线性回归的灵活性结合起来。它是通过拟合简单的模型来集中数据的子集进而建立一个逐点描述数据间的关系的函数将这两种性质结合起来。实际上，这种方法最大的吸引力在于，你不需要列举数据的全部形式的函数来拟合模型，只需要拟合部分数据。例如，这儿有一个对musclemass数据的应急的LOESS拟合曲线，用R函数loess，下面是绘出的图。 plot(massage, data=musclemass, main=Muscle Mass vs Age) out curve(predict(out, newdata=data.frame(age = x), add=TRUE)因为它拟合的太密集了，LOESS一直不太可能在实际中使用它。直到最近，这些个方法被有意识的设计，通过将现在的计算能力的潜力发挥到最大，来实现我们的模型目标。传统的方法是难以实现的。LOESS具体的说明了一种方法，即就是我们所说的局部加权回归法。在数据的每一个点上，一个低水平的多项式（既非线性，也非二次）用来模拟数据，随着预测变量X的值逐渐接近被估计的点，模型就通过加权最小平方这一技术完成了。加权最小平方这种方法给接近估计值的点比较多的权重，远离估计值的点较少的权重。然后该点回归函数的值可以通过利用改点的预测变量值所估计出来的二项式求得。LOESS拟合就在利用回归函数计算完n个数据点的值后完成。这种方法的许多细节是很灵活的。特别地，你可以控制下列三个属性： 多项式模型的阶数。通常它是线性的（直线）或者是二次的（二阶多项式）。 平滑参数。这给一个特定X值的“临近值”了一个定义。也叫做平滑参数，因为它控制了LOESS回归函数的灵活性。大的值会对应最光滑的函数，随着数据的波动它摆动的最小。平滑参数越小，回归函数与数据贴合的越紧密。太小的回归参数值是不可取的，因为回归函数最终会获得数据中的随机误差。 在R函数loess中，光滑参数是由论据跨度给得。 一个权函数决定了在拟合LOESS曲线时在一个邻近的区间上的每一个数据在一个特殊的点扮演了多大的角色。权函数给了估计值的临近点最大的权重，给了估计值最远的点最小的权重。权重的使用基于预测变量空间上临近的点较离得较远的点更有可能以简单的方式联系起来这一思想。按照这种逻辑，遵循局部模型最紧密的点会对局部模型参数产生最大的影响。反之，不太贴合局部模型的点则对局部模型参数估计影响不大。 在数学上，LOESS使用的传统权函数是下面的三次方权函数：更多关于LOESS计算的“残酷”细节，参见我们课堂上将发给大家的讲义。 LOESS的优点：LOESS对描述关系的函数的规格没有要求。 LOESS非常灵活，这使得它对于模拟那些没有理论模型存在的复杂程序非常理想。 它可以被用来查证一个简单的模型（比如一个被强加了线性关系的模型）是否是合理的。 LOESS的缺点： LOESS要求相当大且密集取样的数据以便于拟合出好的模型。这并不奇怪，因为LOESS需要局部结构上正确的经验信息来做局部拟合。 LOESS无法将拟合出的回归函数以数学公式形式呈现，这将会使得分析出的结果之间的交流变得困难。 LOESS有点儿倾向于体现数据集中离群值的效果，就像其他方法一样。24.3 例子例: 小孩的呼吸率。一个高的呼吸率是小孩的呼吸道感染与否的潜在诊断指标。为判断是否真正“高”，一个医师必须对正常的呼吸率的分布有一个清楚的认识。为了这一目的，意大利研究人员测量了年龄从15天至3岁不等（以月份给出）的n = 618个小孩的呼吸率。数据存在于R的respiratory中，用数据做以下几点：1. 拟合一个数据的非参数LOESS回归曲线，并且2. 利用它去计算12个月大的小孩的真正平均呼吸率的95%的置信区间。 解决方案：以下是做这项工作的程序：par(mfrow=c(1,2) # set up two panels for our plotsattach(respiratory)bootpred - numeric(1000) # create vector to hold bootstrapped predictionshyp.age - 12 # set the hypothesized x value for predictionrate.lo - loess(rateage,respiratory)pred - predict(rate.lo, data.frame(age = hyp.age)plot(rateage, data=respiratory) # make scatterplot of the data pointslines(age,rate.lo$fit, col=2) # draw in the fitted loess curvefor (i in 1:1000) n - nrow(respiratory)index - sample(1:n, size=n, replace=TRUE)bootdata - data.frame(rate=respiratory$rateindex,age=respiratory$ageindex)bootfit - loess(rateage, bootdata)bootpredi slrfit predict(slrfit, newdata=data.frame(age=12), int=conf)fit lwr upr1, 38.7036 38.07893 39.32828简单线性回归输出的95%置信度的12个月大的小孩的真正平均呼吸率的区间是（38.08, 39.32 )，而LOESS模型输出的结果是(35.06, 37.46)。很明显，简单线性回归“丢失”了随着年龄增长呼吸率浅浅的减少现象。你可以想象一下，一条直线会通过上面的LOESS拟合线的上方。这种转化为简单线性回归的方法过高的估计了该点的值。你可以（为了好玩和学习）生成一个引导指令的LOESS曲线图。下面是以respiratory为例阐明这一程序。为了好玩哈，我选择将跨度值降低到0.3：跟默认跨度0.75比较起来这个跨度会使拟合变