Mathematica Code for Chernoff's Distribution

(*   The following Mathematica code computes the density function,   *)
(*   distribution function, quantiles, and moments of the random     *)
(*   variable  Z =  argmax_t ( B(t) - t^2)  where   B(t)  is         *)
(*   two-sided Brownian motion starting from zero.  It is based on   *)
(*   the results of Groeneboom, 1985, 1989.                          *)

<= 1;
p[t_,m_]:=p2[t,m]/; t< 1;

(*  Next we calculate the function  h  *)	

h[x_]:= 2*(2/Pi)^(1/2) * 
    Integrate[((2*x+v^2)*v^2 + (1/2)*(2*x + v^2)^2)*Exp[-v^2*(2*x+v^2)^2/2],
       {v,0,Infinity}];

(*   The next step is to calculate  g1,  the part of g2 depending on  p  *)  
(*   this is the hard part of g for x >= -1                              *)

g1[x_,m_]:=
   (2*Pi)^(-1/2)*Integrate[p[t,m]*Exp[-t*(2*x+t)^2 /2], {t,0,Infinity}];


(*  Now we calculate the function g for x >= -1; and call it g2		 *)

g2[x_,m_]:= 2*x -g1[x,m] + h[x];

(*   Now we calculate the function g for x < -1; and call it g3:	 *)

g3term[x_,k_]:=Exp[2^(1/3)*(-Rai[[k]])*Abs[x]]/AiryAiPrime[-Rai[[k]]];
g3[x_,m_]:=4^(1/3)*Exp[-(2/3)*Abs[x]^3]Sum[g3term[x,k],{k,1,m}];



(*  The density f_Z  and distribution function  F_Z   are defined   *)
(*  in terms of  g2 and g3  as follows:                             *)
(*  f_Z[x] == fd2[x] ,  F_Z[x] == Fd2[x]                            *)

fd2[x_]:=g2[x,18]*g2[-x,18]/2 /; x<=1;
fd2[x_]:=g2[x,18]*g3[-x,100]/2/; x> 1;
Fd2[x_] := .5 + NIntegrate[fd2[y], {y,0,x}];


(*  Here is the normal approximation suggested      *) 
(*  by Dykstra and Carolan (1997)                   *)

c:= .52
ndist=NormalDistribution[0.0, .52]
fapprox[x_] := PDF[ndist,x];
Fapprox[x_] := CDF[ndist, x];
Qapprox[u_]:=  Quantile[ndist,u]

(*  Now we compute quantiles:                       *)

QQ[p_] := FindRoot[Fd2[x] ==p, {x,Qapprox[p]}]

(*  Now we calculate moments:                       *)

MM[k_] := Integrate[x^k*fd2[x],{x,0,Infinity}]