545.sagews

#Plot1: marginal likelihood - here the maximum is attained 
# over the parabola. Outside of the maximum set the 
# function is decreasing with n exponentially fast
var('u,v')
plot3d(exp(-15*(u-v^2)^2),(u,-1,1),(v,-1,1))

#Plot2A: regular case in 3d - if the maximum is unique 
# we can use the standard Laplace's approximation.
# In a small neighborhood of $\hat{theta}$ use
# quadratic approximation.
var('u,v')
p1=plot3d(exp(-2*(.5*(u-1/2)^2+2*(v-1/2)^2)),(u,0,1),(v,0,1),color=(0.7,0,0),opacity=0.3)
p2=plot3d(exp(-10*(.5*(u-1/2)^2+2*(v-1/2)^2)),(u,0,1),(v,0,1),color=(.7,.7,0),opacity=0.5)
p3=plot3d(exp(-50*(.5*(u-1/2)^2+2*(v-1/2)^2)),(u,0,1),(v,0,1),color=(0.1,.8,0.1),opacity=0.7)
p4=plot3d(exp(-200*(.5*(u-1/2)^2+2*(v-1/2)^2)),(u,0,1),(v,0,1),color=(.1,0.1,0.8),opacity=1)
show(p1+p2+p3+p4)

#Plot2B: regular case in 3d - the same as (2A) but we
# show a section of the plot in u=0.5 which contains
# the maximum
var('u,v')
p1=plot3d(exp(-2*(.5*(u-1/2)^2+2*(v-1/2)^2)),(u,0,.5),(v,0,1),color=(0.7,0,0))
p2=plot3d(exp(-10*(.5*(u-1/2)^2+2*(v-1/2)^2)),(u,0,.5),(v,0,1),color=(.7,.7,0))
p3=plot3d(exp(-50*(.5*(u-1/2)^2+2*(v-1/2)^2)),(u,0,.5),(v,0,1),color=(0.1,.8,0.1))
p4=plot3d(exp(-200*(.5*(u-1/2)^2+2*(v-1/2)^2)),(u,0,.5),(v,0,1),color=(.1,0.1,0.8))
show(p1+p2+p3+p4)

#Plot3: a concrete singular example - the normalized likelihood
# in the case when the maximum set is singular
var('u,v')
p1=plot3d(exp(-10*(u^2-v^3)^2),(u,-1,1),(v,-1,1),color=(0.9,0.9,0.9),opacity=0.3)
p2=plot3d(exp(-50*(u^2-v^3)^2),(u,-1,1),(v,-1,1),color=(0.9,0.9,0.5),opacity=0.6)
p3=plot3d(exp(-150*(u^2-v^3)^2),(u,-1,1),(v,-1,1),color=(0.1,0.5,0.1))
show(p1+p2+p3)

#Plot4: a fiber of $\hat{p} if covariance is nonzero -  
# the model is unidentifiable but the maximum 
# parameter set is non-singular 
u, v = var('u,v')
f1 = (u, 1/(10*u), v)
f2 = (u, 1/(10*u), v)
p1 = parametric_plot3d(f1, (u,-1,-0.1), (v,-1,1), texture="orange")
p2 = parametric_plot3d(f2, (u,0.1,1), (v,-1,1), texture="violet")
show(p1+p2)

#Plot5: a fiber of $\hat{p} if covariance is zero -  
# the maximum parameter set is singular 
u, v = var('u,v')
f1 = (u, 0, v)
f2 = (0, u, v)
p1 = parametric_plot3d(f1, (u,-1,1), (v,-1,1), texture="orange")
p2 = parametric_plot3d(f2, (u,-1,1), (v,-1,1), texture="violet")
show(p1+p2)

#Plot6: the marginal likelihood if the covariance is nonzero -  
# for large n the behavior is nice 
var('u,v')
p1=plot3d(exp(-20*(u*v-1/10)^2),(u,-1,1),(v,-1,1),color=(.8,.8,.8),opacity=0.3)
p2=plot3d(exp(-200*(u*v-1/10)^2),(u,-1,1),(v,-1,1),color=(.9,.9,.9),opacity=0.5)
p3=plot3d(exp(-500*(u*v-1/10)^2),(u,-1,1),(v,-1,1),color=(0.1,0.3,0.5))
show(p1+p2+p3)

#Plot7: the marginal likelihood if the covariance is zero -  
# the largest contribution to the integral comes from 
# a neighborhood of the singularity 
var('u,v')
plot3d(exp(-50*u^2*v^2),(u,-1,1),(v,-1,1),color=(0.1,0.3,0.5))