**A**s a continuation of the Bayesian resolution of last week riddle, I looked at a numeric resolution of the four secretaries problem, while in the train back from Rouen (and trying to block the chatter of my neighbours, a nuisance I find myself more and more sensitive to!). The target function is defined as

gainz=function(b,c,T=1e4,type="raw"){ x=matrix(runif(4*T),ncol=4) maz=t(apply(x,1,cummax)) zam=t(apply(x[,4:1],1,cummax)) if (type=="raw"){return(mean( ((x[,2]>b*x[,1])*x[,2]+ (x[,2]<b*x[,1])*((x[,3]>c*maz[,2])*x[,3]+ (x[,3]<c*maz[,2])*x[,4]))/maz[,4]))} if (type=="global"){return(mean( ((x[,2]>b*x[,1])*(x[,2]==maz[,4])+ (x[,2]<b*x[,1])*((x[,3]>c*maz[,2])*(x[,3]==maz[,4])+ (x[,3]<c*maz[,2])*(x[,4]==maz[,4])))))} if (type=="remain"){return(mean( ((x[,2]>b*x[,1])*(x[,2]==zam[,3])+ (x[,2]<b*x[,1])*((x[,3]>c*maz[,2])*(x[,3]==zam[,2])+ (x[,3]<c*maz[,2])*(x[,4]==zam[,2])))))}}

where the data is generated from a U(0,1) distribution as the loss functions are made scaled free by deciding to always sacrifice the first draw, x¹. This function is to be optimised in (b,c) and hence I used a plain vanilla simulated annealing R code:

avemale=function(T=3e4,type){ b=c=.5 maxtar=targe=gainz(b,c,T=1e4,type) temp=0.1 for (t in 1:T){ bp=b+runif(1,-temp,temp) cp=c+runif(1,-temp,temp) parge=(bp>0)*(cp>0)*gainz(bp,cp,T=1e4,type) if (parge>maxtar){ b=bs=bp;c=cs=cp;maxtar=targe=parge}else{ if (runif(1)<exp((parge-targe)/temp)){ b=bp;c=cp;targe=parge}} temp=.9999*temp} return(list(bs=bs,cs=cs,max=maxtar))}

with outcomes

- b=1, c=.5, and optimum 0.8 for the raw type
- b=c=1 and optimum 0.45 for the global type
- b undefined, c=2/3 and optimum 0.75 for the remain type