核密度估计是在概率论中用来估计未知的密度函数,属于非参数检验方法之一,由Rosenblatt (1955)和Emanuel Parzen(1962)提出,又名Parzen窗(Parzen window)。
假设我们有n个数X1-Xn,我们要计算某一个数X的概率密度有多大。核密度估计的方法是这样的:
其中K为核密度函数,
h为设定的窗宽。
库文档分享
核密度估计的原理其实是很简单的。在我们对某一事物的概率分布的情况下。
如果某一个数在观察中出现了,我们可以认为这个数的概率密度很大,和这个数比较近的数的概率密度也会比较大,而那些离这个数远的数的概率密度会比较小。
基于这种想法,针对观察中的第一个数,我们都可以f(x-xi)去拟合我们想象中的那个远小近大概率密度。
当然其实也可以用其他对称的函数。针对每一个观察中出现的数拟合出多个概率密度分布函数之后,取平均。
如果某些数是比较重要,某些数反之,则可以取加权平均。
库文档分享
但是核密度的估计并不是,也不能够找到真正的分布函数。我们可以举一个极端的例子:在R中输入:
plot(density(rep(0, 1000)))
可以看到它得到了正态分布的曲线,但实际上呢?从数据上判断,它更有可能是一个退化的单点分布。
库文档分享
但是这并不意味着核密度估计是不可取的,至少他可以解决许多模拟中存在的异方差问题。比如说我们要估计一下下面的一组数据:
(10)
dat<c(rgamma(300,shape=2,scale=2),rgamma(100,shape=10,scale=2))
库文档分享
库文档分享
可以看出它是由300个服从gamma(2,2)与100个gamma(10,2)的随机数构成的,他用参数统计的办法是没有办法得到一个好的估计的。那么我们尝试使用核密度估计:
plot(density(dat),ylim=c(0,))
库文档分享
将利用正态核密度与标准密度函数作对比
dfn<-function(x,a,alpha1,alpha2,theta){
a*dgamma(x,shape=alpha1,scale=theta)+(1-a)*dgamma(x,shape=alpha2,scale=theta)}
pfn<-function(x,a,alpha1,alpha2,theta){
a*pgamma(x,shape=alpha1,scale=theta)+(1-a)*pgamma(x,shape=alpha2,scale=theta)}
curve(dfn(x,,2,10,2),add=T,col="red")
库文档分享
得到下图:
(红色的曲线为真实密度曲线)
库文档分享
可以看出核密度与真实密度相比,得到大致的估计是不成问题的。至少趋势是得到了的。如果换用gamma分布的核效果无疑会更好,但是遗憾的是r中并没有提供那么多的核供我们挑选(其实我们知道核的选择远没有窗宽的选择来得重要),所以也无需介怀。
R中提供的核:kernel = c("gaussian", "epanechnikov", "rectangular", "triangular", "biweight","cosine", "optcosine")。
库文档分享
我们先来看看窗宽的选择对核密度估计的影响:
dfn1<-function(x){
*dnorm(x,3,1)+*dnorm(x,-3,1)}
par(mfrow=c(2,2))
curve(dfn1(x),from=-6,to=6)
data<-c(rnorm(200,3,1),rnorm(200,-3,1))
plot(density(data,bw=8))
plot(density(data,bw=))
plot(density(data,bw=))
库文档分享
R语言与核密度估计非参数统计 来自淘豆网m.daumloan.com转载请标明出处.