**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