Nothing
## ipeglim/inst/thesis/fig/o3ef.R
##
## o3ef = order 3 exponential family of distribution (i.e., normal prior)
library(ipeglim)
rm(list=ls())
lmu <- c(-3, 2, 1, 0, -1, -2, 3)
e0 <- seq(0,2,0.25)
e2 <- seq(0,5,0.3)
hcube <- expand.grid(e0=e0, e2=e2, lmu=lmu)
e1 <- apply(hcube, 1, function(x){
x <- as.vector(x)
rt.e1 <- tryCatch(uniroot(eta1fn, lower=-20, upper=30, tol=1e-4, eta0=x[1], eta2=x[2], logm=x[3])$root, error=function(e) return(NaN))
return(rt.e1)
})
hcube$e1 <- e1
z <- range(e1, na.rm=TRUE)
tobj <- wireframe(e1~e0*e2, data=hcube,
group=lmu,
screen=list(y=60, x=15, y=20),
# zoom=0.85,
pretty=TRUE,
xlim=c(0,2),
ylim=c(0,5),
zlim=c(z[1],z[2]),
main=expression("Level Sets of "~E~"("~theta~"|"~eta[0]~","~eta[1]~","~eta[2]~")"),
sub=expression("log"~"("~mu~")"~"="~-3 ~"to"~3~"by"~1~"from left to right"),
xlab=expression(eta[0]),
ylab=expression(eta[2]),
zlab=expression(eta[1]),
par.settings = list(axis.line = list(col = "transparent")),
par.box=list(col='black', lwd=0.25),
scales=list(arrows=FALSE, col="black"), # show axese with ticks (z.ticks=2),
type="l"
)
print(tobj)
if(FALSE){
# xyplot(e1~e2|e0, data=hcube[hcube$e0 %in% c(0,0.25),])
}
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.