Regiones intrinsecas creibles

Caso Exponencial

PDF[ExponentialDistribution[θ],x]

model[x_,θ_]:=(E^(-x θ)) θ;

k[θ1_,θ2_]:=NIntegrate[model[x,θ2]Log[((model[x,θ2])/(model[x,θ1]))],{x,0,∞}];
δ[θ1_,θ2_]:=Min[k[θ2,θ1],k[θ1,θ2]];
Timing[δ[1,2]]

θθ=2;Plot[δ[θ0,θθ],{θ0,0.01,10}]

Simplify[PDF[GammaDistribution[α,1/β],θ],{α>0,β>0,θ>0}]

ga[θ_,α_,β_]:=(((E^(-β θ)) (β^α) (θ^(-1+α)))/(Gamma[α]));
post[θ_,n_,t_]:=ga[θ,n,t];
nn=12; tt=4.88;
fig1=Plot[post[θ,nn,tt],{θ,0,10},Axes->False, Frame->True]

Solucion numérica

d[θ0_,n_,t_]:=nn NIntegrate[δ[θ0,θ] post[θ,n,t],{θ,0,∞}];
Timing[d[3,nn,tt]]

pts=Table[{θ0,d[θ0,nn,tt]},{θ0,0.1,6,0.25}]

dap=Interpolation[pts];
fig2=Plot[dap[θ0],{θ0,0,10},Axes->False, Frame->True]

fig3=Show[GraphicsArray[{{fig1},{fig2}}]]

FindMinimum[dap[θ0],{θ0,2}]

estθexact=θ0/.FindMinimum[dap[θ0],{θ0,2}][[2]]

hh=2.5;FindRoot[dap[θ0]==hh,{θ0,1,2}]

corte1[h_]:=θ0/.FindRoot[dap[θ0]==h,{θ0,1,2}];
corte1[2]

corte2[h_]:=θ0/.FindRoot[dap[θ0]==h,{θ0,2,5}];
corte2[2]

prob[h_]:=NIntegrate[post[θ,nn,tt],{θ,corte1[h],corte2[h]}];
prob[2]

pp=0.95;FindRoot[prob[h]==pp,{h,2,3}]

cut[p_]:=h/.FindRoot[prob[h]==p,{h,2,3}];
hh=cut[pp];
{hh,corte1[hh],corte2[hh]}

Solución aproximada

Integrate[(1/θ),θ]

phi[θ_]:=Log[θ];
theta[ϕ_]:=Exp[ϕ];

estϕ=NIntegrate[phi[θ]post[θ,nn,tt],{θ,0,∞}]

estθapp=theta[estϕ]

{estθexact, estθapp}

nq[p_]:=Quantile[NormalDistribution[0,1],(1+p)/2];
nq[0.95]

regϕ={estϕ-((nq[pp])/(Sqrt[nn])),estϕ+((nq[pp])/(Sqrt[nn]))}

regθ=theta[regϕ]

{{corte1[hh],corte2[hh]},regθ}

Caso Uniforme

PDF[UniformDistribution[0,θ],x]

model[x_,θ_]:=If[0≤x≤θ,(1/θ),0];

k[θ1_,θ2_]:=If[θ2<θ1,Log[((θ1)/(θ2))],∞];
δ[θ1_,θ2_]:=Min[k[θ2,θ1],k[θ1,θ2]];
Timing[δ[1,2]]

θθ=2;Plot[δ[θ0,θθ],{θ0,0.01,10}]

Simplify[PDF[ParetoDistribution[τ,α],θ],{α>0,τ>0,θ>τ}]

pa[θ_,α_,τ_]:=α (θ^(-1-α)) (τ^α);
post[θ_,n_,t_]:=pa[θ,n,t];
θθ=2;nn=12; tt=1.8068;
fig1=Plot[post[θ,nn,tt],{θ,tt,4},Axes->False, Frame->True,PlotRange->All]

NIntegrate[post[θ,nn,tt],{θ,tt,∞}]

Solucion numérica

d[θ0_,n_,t_]:=nn NIntegrate[δ[θ0,θ] post[θ,n,t],{θ,t,∞}];
Timing[d[3,nn,tt]]

pts=Table[{θ0,d[θ0,nn,tt]},{θ0,0.1,6,0.1}]

dap=Interpolation[pts];
fig2=Plot[dap[θ0],{θ0,0,6},Axes->False, Frame->True]

fig3=Show[GraphicsArray[{{fig1},{fig2}}]]

FindMinimum[dap[θ0],{θ0,2}]

estθexact=θ0/.FindMinimum[dap[θ0],{θ0,2}][[2]]

hh=2.5;FindRoot[dap[θ0]==hh,{θ0,1,2}]

corte1[h_]:=Max[tt,θ0/.FindRoot[dap[θ0]==h,{θ0,1,2}]];
corte1[2]

corte2[h_]:=θ0/.FindRoot[dap[θ0]==h,{θ0,2,5}];
corte2[2]

{{nn,tt},{corte1[2],corte2[2]}}

prob[h_]:=NIntegrate[post[θ,nn,tt],{θ,corte1[h],corte2[h]}];
prob[2]

pp=0.95;FindRoot[prob[h]==pp,{h,2,3}]

cut[p_]:=h/.FindRoot[prob[h]==p,{h,2,3}];
hh=cut[pp];
{hh,corte1[hh],corte2[hh]}

Solución aproximada

Integrate[(1/θ),θ]

phi[θ_]:=Log[θ];
theta[ϕ_]:=Exp[ϕ];

estϕ=NIntegrate[phi[θ]post[θ,nn,tt],{θ,tt,∞}]

estθapp=theta[estϕ]

{estθexact, estθapp}

nq[p_]:=Quantile[NormalDistribution[0,1],(1+p)/2];
nq[0.95]

regϕ={estϕ-((nq[pp])/(Sqrt[nn])),estϕ+((nq[pp])/(Sqrt[nn]))}

regθ=theta[regϕ]

{{corte1[hh],corte2[hh]},regθ}