非参数估计
本文由R第六组共同完成
核函数估计概论密度函数
选取了Gaussian Kernel公式。构建一个名为kernel.function的函数,参数为x和h。以下为kernel.function的具体形式
kernel.function<- function(x,h)
{x=sort(x)
n=length(x)
s=0;t=0;y=0
for(i in 2:n)
s[i]=0
for(i in 1:n)
{for(j in 1:n)
s[i]=s[i]+exp(-((x[i]-x[j])^2)/(2hh))
t[i]=s[i]
}
for(i in 1:n)
y[i]=t[i]/(nhsqrt(2*pi))
z=complex(re=x,im=y)
hist(x,freq=FALSE)
lines(z)
}
lines(z)
把核密度的函数画出来
opar<-par(no.readonly=TRUE)
par(mfrow=c(2,2))
x=rnorm(100) 将x设置为100个服从标准正态分布的随机数
kernel.function(x,0.5)
kernel.function(x,1)
kernel.function(x,1.5)
kernel.function(x,2)
可以看出随着h的增加,拟合程度先变好,后变差。而且尾部的拟合程度越来越差。
求最优h:构造两个函数分别为fun1和fun2
fun1为CV(h)表达式减号前面部分,及对f(x)平方求积分一式。
在fun1中,先筛选出符合[0,1]范围内的x,再把所有符合条件的x代入f(x)中,最后对f(x)的平方求和。
fun1<-function(x,h)
{x1<-x[x>0]
x2<-x[x<1]
k=intersect(x1,x2)
k=sort(k)
n=length(k)
s=0;t=0;y=0
for(i in 2:n)
s[i]=0
for(i in 1:n)
{for(j in 1:n)
s[i]=s[i]+exp(-((k[i]-k[j])^2)/(2hh))
t[i]=s[i]
}
for(i in 1:n)
y[i]=t[i]/(nhsqrt(2pi))
sum=sum(yy)
print(sum)
}
fun2为CV(h)表达式减号后面部分。用循环结构即可实现。
fun2<-function(x,h){
sum2<-0
for(i in 1:n)
{for(j in 1:n)
{
if(j!=i)
sum2 <- sum2+exp(-((x[i]-x[j])^2)/(2hh))
}}
sum2 <- sum22/(n(n-2)hsqrt(2*pi))
}
构造函数f,f为fun1-fun2的差值。f中的参数为x和h
f<-function(x,h){(fun1(x,h)-fun2(x,h))}
for (h in seq(0.1,2,0.1)) {
f(x,h) }
局部多项式估计函数:
my_fun <- function(x,h,p){
n <- length(mtcarsKaTeX parse error: Expected '}', go