# From http://www.stat.ucl.ac.be/ISpersonnel/lecoutre/stats/fichiers/_gallery.pdf # 3-D plots # # mu1<-0 # setting the expected value of x1 mu2<-0 # setting the expected value of x2 s11<-10 # setting the variance of x1 s12<-15 # setting the covariance between x1 and x2 s22<-10 # setting the variance of x2 rho<-0.5 # setting the correlation coefficient between x1 and x2 x1<-seq(-10,10,length=41) # generating the vector series x1 x2<-x1 # copying x1 to x2 # f<-function(x1,x2) { term1<-1/(2*pi*sqrt(s11*s22*(1-rho^2))) term2<--1/(2*(1-rho^2)) term3<-(x1-mu1)^2/s11 term4<-(x2-mu2)^2/s22 term5<--2*rho*((x1-mu1)*(x2-mu2))/(sqrt(s11)*sqrt(s22)) term1*exp(term2*(term3+term4-term5)) } # setting up the function of the multivariate normal density # z<-outer(x1,x2,f) # calculating the density values # persp(x1, x2, z, main="Two dimensional Normal Distribution", sub=expression(italic(f)~(bold(x))==frac(1,2~pi~sqrt(sigma[11]~ sigma[22]~(1-rho^2)))~phantom(0)^bold(.)~exp~bgroup("{", list(-frac(1,2(1-rho^2)), bgroup("[", frac((x[1]~-~mu[1])^2, sigma[11])~-~2~rho~frac(x[1]~-~mu[1], sqrt(sigma[11]))~ frac(x[2]~-~mu[2],sqrt(sigma[22]))~+~ frac((x[2]~-~mu[2])^2, sigma[22]),"]")),"}")), col="lightgreen", theta=30, phi=20, r=5.0, d=0.1, expand=0.5, ltheta=90, lphi=180, shade=0.75, ticktype="detailed", nticks=5) # produces the 3-D plot # mtext(expression(list(mu[1]==0,mu[2]==0,sigma[11]==10,sigma[22]==10,sigma[12 ]==15,rho==0.5)), side=3) # adding a text line to the graph