[R] about different bandwidths in one graph

Peter Ehlers ehlers at ucalgary.ca
Wed Jul 18 02:29:32 CEST 2012


On 2012-07-17 11:15, chester123 wrote:
> Thank you in advance.
>
> Now I want to make comparison of the different bandwidth h in a normal
> distribution graph.
>
> This is the table of bandwidth h: thumb rule (normal)--0.00205; thumb
> rule(Epanech.)--0.00452; Plug-in (normal)--0.0009;
> Plug-in(Epanech.)--0.002.
>
> this is the condition: N=1010 data sample is from normal distribution
> N(0,0.0077^2). The grid points are taken to be [-0.05,0.05] and increment is
> 10. Bandwidth is taken the above h value r respectively and the kernel can
> be Epanechnikov kernel or Gaussian kernel.
>
> The following is my code:
> #########################################################
> # Define the Epanechnikov kernel function
> kernel<-function(x){0.75*(1-x^2)*(abs(x)<=1)}
> ############################################################### # Define the
> kernel density estimator
> kernden=function(x,z,h,ker){
> # parameters: x=variable; h=bandwidth; z=grid point; ker=kernel
> nz<-length(z)
> nx<-length(x)
> x0=rep(1,nx*nz)
> dim(x0)=c(nx,nz)
> x1=t(x0)
> x0=x*x0
> x1=z*x1
> x0=x0-t(x1)
> if(ker==1){x1=kernel(x0/h)}        # Epanechnikov kernel
> if(ker==0){x1=dnorm(x0/h)}       # normal kernel
> f1=apply(x1,2,mean)/h
> return(f1)
> }
> #################################################################################################################################
> Simulation for different bandiwidths and different kernels
> n=1010	# n=1010
> ker=1         # ker=1=>Epan; ker=0=> Gaussian
> h0=c(0.00452,0.001984)    # set initial bandwidths
> z=seq(-0.05,0.05,by=10)     # grid points
> nz=length(z)                    # number of grid points
> x=rnorm(1010, mean=0, sd=0.0077)      # simulate x-N(0,0.0077^2)
> if(ker==1){h_o=2.34*n^{-0.2}}            # bandwidth for Epanechnikov kernel
> if(ker==0){h_o=1.06*n^{-0.2}}            # bandwidth for normal kernel
> f1=kernden(x,z,h0[1],ker)
> f2=kernden(x,z,h0[2],ker)
> f3=kernden(x,z,h0[3],ker)
> f4=kernden(x,z,h0[4],ker)
> text1=c("True","h=0.0025","h=0.00452","h=0.0009","h=0.002")
> data=cbind(dnorm(z),f1,f2,f3,f4)	# combine them as a matrix win.graph()
> matplot(z,data,type="l",lty=1:5,col=1:5,xlab="",ylab="")
> legend(-1,0.2,text1,lty=1:5,col=1:5)
> ################################################################
>
> But the error message is "Error in strwidth(legend, units = "user", cex =
> cex, font = text.font) :
>    plot.new has not been called yet".
>
> I know something is wrong in the code but don't know where.

I would guess that you somehow de-activated the graphics device
before issuing the legend() call. I can't comment on the rest of
your code - haven't examined it.

Peter Ehlers



More information about the R-help mailing list