R语言与非参数统计(核密度估计)
核密度估计是在概率论中用来估计未知的密度函数,属于非参数检验方法之一,由Rosenblatt (1955)和Emanuel Parzen(1962)提出,又名Parzen窗(Parzen window)。
假设我们有n个数X1-Xn,我们要计算某一个数X的概率密度有多大。核密度估计的方法是这样的:
其中K为核密度函数,h为设定的窗宽。
核密度估计的原理其实是很简单的。在我们对某一事物的概率分布的情况下。如果某一个数在观察中出现了,我们可以认为这个数的概率密度很大,和这个数比较近的数的概率密度也会比较大,而那些离这个数远的数的概率密度会比较小。基于这种想法,针对观察中的第一个数,我们都可以f(x-xi)去拟合我们想象中的那个远小近大概率密度。当然其实也可以用其他对称的函数。针对每一个观察中出现的数拟合出多个概率密度分布函数之后,取平均。如果某些数是比较重要,某些数反之,则可以取加权平均。
但是核密度的估计并不是,也不能够找到真正的分布函数。我们可以举一个极端的例子:在R中输入:
plot(density(rep(0, 1000)))
可以看到它得到了正态分布的曲线,但实际上呢?从数据上判断,它更有可能是一个退化的单点分布。
但是这并不意味着核密度估计是不可取的,至少他可以解决许多模拟中存在的异方差问题。比如说我们要估计一下下面的一组数据:
可以看出它是由300个服从gamma(2,2)与100个gamma(10,2)的随机数构成的,他用参数统计的办法是没有办法得到一个好的估计的。那么我们尝试使用核密度估计:
plot(density(dat),ylim=c(0,0.2))
将利用正态核密度与标准密度函数作对比
得到下图:
(红色的曲线为真实密度曲线)
可以看出核密度与真实密度相比,得到大致的估计是不成问题的。至少趋势是得到了的。如果换用gamma分布的核效果无疑会更好,但是遗憾的是r中并没有提供那么多的核供我们挑选(其实我们知道核的选择远没有窗宽的选择来得重要),所以也无需介怀。
R中提供的核:kernel = c("gaussian", "epanechnikov", "rectangular", "triangular", "biweight","cosine", "optcosine")。
我们先来看看窗宽的选择对核密度估计的影响:
得到下图,我们可以清楚的看到带宽为0.8恰好合适,其余的不是拟合不足便是过拟合。
窗宽究竟该如何选择呢?
我们这里不加证明的给出最佳窗宽选择公式:
(这个基于积分均方误差最小的角度得到的)
这里介绍两个可操作的窗宽估计办法:(这两种方法都比较容易导致过分光滑)
1、 Silverman大拇指法则
这里使用R(phi’’)/sigma^5估计R(f’’),phi代表标准正态密度函数,得到h的表达式:
h=(4/(3n))^(*1/5)*sigma
2、 极大光滑原则
h=3*(R(K)/(35n))^(1/5)*sigma
当然也有比较麻烦的窗宽估计办法,比如缺一交叉验证,插入法等,可以参阅《computational statistics》一书
我们用上面的两种办法得到的窗宽是多少,他的核密度估计效果好吗?
我们还是以上面的混合正态数据为例来看看效果。
使用大拇指法则,将数据n=400,sigma=3.030658,带入公式,h=0.9685291
使用极大光滑原则,假设K为正态核,R(K)=1/(sqrt(2*pi)),h=1.121023
可以看出他们都比我们认为的h=0.8要大一些,作图如下:
plot(density(data,bw=0.9685))
plot(density(data,bw=1.1210))
由我们给出的
以Gauss核为例做核密度估计
用Gauss核做核密度估计的R程序如下(还是使用我们的混合正态密度的例子):
作图如下:
最后说一个R的内置函数density()。其实我觉得如果不是为了简要介绍核密度估计的一些常识我们完全可以只学会这个函数
先看看函数的基本用法:
density(x, ...)
## Default S3 method:
density(x, bw = "nrd0", adjust = 1,
kernel = c("gaussian", "epanechnikov", "rectangular",
"triangular", "biweight",
"cosine", "optcosine"),
weights = NULL, window = kernel, width,
give.Rkern = FALSE,
n = 512, from, to, cut = 3, na.rm = FALSE, ...)
对重要参数做出较为详细的说明:
X:我们要进行核密度估计的数据
Bw:窗宽,这里可以由我们自己制定,也可以使用默认的办法nrd0: Bandwidth selectors for Gaussian kernels。我们还可以使用bw.SJ(x,nb = 1000, lower = 0.1 * hmax, upper = hmax, method = c("ste","dpi"), tol = 0.1 * lower),这里的method =”dpi”就是前面提到过的插入法,”ste”代表solve-the-equationplug-in,也是插入法的改进
Kernel:核的选择
Weights:对比较重要的数据采取加权处理
对于上述混合正态数据data,有
> density(data)
Call:
density.default(x = data)
Data: data (400 obs.); Bandwidth 'bw' = 0.8229
x y
Min. :-7.5040 Min. :0.0000191
1stQu.:-3.5076 1st Qu.:0.0064919
Median : 0.4889 Median :0.0438924
Mean :0.4889 Mean :0.0624940
3rdQu.: 4.4853 3rd Qu.:0.1172919
Max. :8.4817 Max. :0.1615015
知道带宽:h=0.8229(采取正态密度核)那么带入密度估计式就可以写出密度估计函数。
最后以faithful数据集为例说明density的用法:
R数据集faithful是old faithful火山爆发的数据,其中“eruption”是火山爆发的持续时间,waiting是时间间隔
对数据“eruption”做核密度估计
R程序:
知道h= 0.3348
作图:
数据分析咨询请扫描二维码
CDA数据分析师在中国航信高科技产业园进行了面向测试度量的数据分析培训课程,培训人数近2 ...
2024-05-01CDA数据分析师走进深圳迈瑞生物医疗电子股份有限公司,在迈瑞总部展开了为期两天的培训,本次课程参训人员线上及线下近百人, ...
2024-05-01CDA数据分析师在合肥市对合肥阳光新能源科技有限公司开展了为期8天的企业内训。 合肥阳光新能源科技 ...
2024-05-01CDA数据分析师走进海尔大学,进行了《数据治理与数据中台建设的道与术》专题培训,培训现场爆满,近百人参加了此次培训。 ...
2024-05-01在中国银行苏州分行培训中心开始数据分析师培训,此次培训课程共10天内容,包括Excel、MySQL、概率论与数理统计、SPSS等内容, ...
2024-05-01从实际的业务需求出发,结合行业的典型应用特点,围绕实际的商业问题,探讨数据挖掘、机器学习模型在金融领域的应用,包括获客、信用评分、细分画像、交叉销售、反欺诈、违规识别、时序预测、运筹优化、流程挖掘九个方面,形成 ...
2024-05-01本次培训课程为线上+线下的模式,由于学员编程能力不一、部分学员没有编程基础,故提供统计学、python基 ...
2024-05-01华夏银行信用卡中心-机器学习培训 1、课程亮点 取材于业界一流企业和顶级咨询公司的行业实践;已经被证明是人人 ...
2024-05-01主 题:数据中台建设及数据分析应用主题分享 1. 数据中台市场洞察 2. 主流数据中台产品比较 3. 某企业数据中 ...
2024-05-01围绕“数据驱动”战略,全力打造我行 300 人数字化人才梯队,着力培养数字化管理人才、大数据专业团队 ...
2024-05-01在当今数据驱动的商业环境中,数据分析成为了企业决策的重要依据。通过对大量数据的收集、处理和分析,企业能够更好地理解市场 ...
2024-04-29在人工智能(AI)的世界里,提示词(Prompt)是一种强大的工具,它能够引导AI按照用户的需求产生特定的输出。本文将深入探讨AI ...
2024-04-29CDA立足未来职场,拓展前沿视野——对外经贸大学保险学院举办“三全育人大讲堂”分享行业最新动态。 ...
2024-04-294月2日,CDA数据分析师创始发起人兼协会理事长赵坚毅博士受邀在浙江万里学院举办了一场以“数字化能力在职场中的作用” ...
2024-04-29随机森林(Random Forests)现在机器学习中比较火的一个算法,是一种基于Bagging的集成学习方法,能够很好地处理分类和回归的问 ...
2022-12-23方差分析是数据分析中常用的一种统计分析方法,接下来让我们简单了解一下方差分析的基本思想和原理吧。 方差分析(Analysis ...
2022-12-23来源:关于数据分析与可视化 关于streamlit-aggrid 数据排序 表格样式的调整 数据 ...
2022-08-03作者:麦叔 定义 「把上面晦涩的概念汇成一句话就是:」 ❝ 回调函数就是一个被作为参 ...
2022-08-03现今,高学历人群日益增多,物以稀为贵的高学历光环淡去。无论本科生还是研究生,甚至博士生,求职竞争力都大不如前,就业压力越来越大。
2022-06-01某家企业10个人面试,有9个本科生……如何脱颖而出,除得体的举止和良好的沟通力外,证书成重要筹码,这也是很多人考证的关键所在。
2022-04-14