In [1]:
# rm(list=ls()) 
options(OutDec = ",")
#==============================================================================
# Figura elipse
#==============================================================================
silence <- suppressPackageStartupMessages # Para omitir mensagens de alertas
silence(library(ellipse))
xbarra  <- c(0,0)
S       <- matrix(c(1,0.75,0.75,1),2,2,byrow=T)
IS      <- solve(S)
DES     <- eigen(S)
c2      <- 6
eval    <- DES$values
evec    <- DES$vectors
y0      <- sqrt(c2 * eval[1]) * evec[,1] + xbarra
y1      <- sqrt(c2 * eval[2]) * evec[,2] + xbarra
y2      <- xbarra - sqrt(c2 * eval[1]) * evec[,1]
y3      <- xbarra - sqrt(c2 * eval[2]) * evec[,2]
In [2]:
#==============================================================================
# Elipse constante
#==============================================================================
v       <- ellipse(S,t = sqrt(1/c2),npoints=2000)
v[,1]   <- xbarra[1] + v[,1]
v[,2]   <- xbarra[2] + v[,2]
In [3]:
#==============================================================================
# Grafico
#==============================================================================
par(mfrow=c(1,1),lwd=2,cex.lab=1.5,cex.axis=1.5,lab=c(10,10,5),
    mar=c(4.5,5,2,2))
plot(v,col="blue",bty="n",xaxt="n",yaxt="n",type="l",lwd=3,xlim=c(-0.55,0.55),
    ylim=c(-0.55,0.55),xlab=expression(X[1]),ylab=expression(X[2]))
abline(v=0,col=1,lwd=2)
abline(h=0,col=1,lwd=2)
arrows(xbarra[1],xbarra[2],y0[1],y0[2],col="red",lwd=3,code=0,lty=2)
arrows(xbarra[1],xbarra[2],y1[1],y1[2],col="red",lwd=3,code=0,lty=2)
arrows(xbarra[1],xbarra[2],y2[1],y2[2],col="red",lwd=3,code=0,lty=2)
arrows(xbarra[1],xbarra[2],y3[1],y3[2],col="red",lwd=3,code=0,lty=2)
text(0.12,0.05,expression(theta),cex=1.5)
xseq <- seq(0.07,0.1,0.001) 
r    <- 0.1
lines(xseq,sqrt(r^2-xseq^2),lwd=3)
text(0.14,-0.45,expression(mu==list("(0","0)")),cex=1.5)
text(0.14,-0.52,expression(list(rho==0,75)),cex=1.5)
text(-0.2,-0.45,expression(c^2==paste(x^T,Sigma^{-1},x,sep="")),cex=1.5)
text( 0.4, 0.55,expression(y[1]==paste(v[1]^T,x,sep="")),cex=1.5)
text(-0.4, 0.55,expression(y[2]==paste(v[2]^T,x,sep="")),cex=1.5)