The book, The Leverage Space Trading Model by RALPH VINCE

Page 131, binned data and joint probabilities

MktSysA=c(-150.00,-45.33,-45.33,13.00,13.00,13.00,13.00,13.00,
          79.67,79.67,79.67,136)
MktSysB=c(253.00,-1000.00,-64.43,-64.43,-64.43,253.00 ,253.00,
          448.00,-64.43,-64.43,-64.43,253)
MktSysC=c(533.00,220.14,220.14,-500.00,533.00,220.14,799.00,
          220.14,-325.00,220.14,533.00,220.14)
PL=cbind(MktSysA,MktSysB,MktSysC)

p=c(0.076923077,0.076923077,0.153846154,0.076923077,0.076923077,
    0.076923077,0.076923077,0.076923077,0.076923077,0.076923077,
    0.076923077,0.076923077)
cbind(PL,p)
##       MktSysA  MktSysB MktSysC          p
##  [1,] -150.00   253.00  533.00 0.07692308
##  [2,]  -45.33 -1000.00  220.14 0.07692308
##  [3,]  -45.33   -64.43  220.14 0.15384615
##  [4,]   13.00   -64.43 -500.00 0.07692308
##  [5,]   13.00   -64.43  533.00 0.07692308
##  [6,]   13.00   253.00  220.14 0.07692308
##  [7,]   13.00   253.00  799.00 0.07692308
##  [8,]   13.00   448.00  220.14 0.07692308
##  [9,]   79.67   -64.43 -325.00 0.07692308
## [10,]   79.67   -64.43  220.14 0.07692308
## [11,]   79.67   -64.43  533.00 0.07692308
## [12,]  136.00   253.00  220.14 0.07692308
f=c(.307,0,.693)
BL=apply(PL,2,min)

ret=f*t(PL)/abs(BL)
Net_HPR=1+colSums(ret)
Net_HPR_p=Net_HPR^p
GHPR=prod(Net_HPR_p)
GHPR
## [1] 1.248538

Page 137, generate sim dataset, compute GHPR with f. 

q=12
b=.8# 1 minus b is drawndown percentage
sample_size=60000 #6,000,000 samples need 10 minutes running
outcomes=1:12
raw=sample(outcomes,size=q*sample_size,replace = TRUE,prob =p)
data=matrix(raw,nrow=sample_size)

f=c(.191,.007,.165)
ret=f*t(PL)/abs(BL)
Net_HPR=1+colSums(ret)
Net_HPR_p=Net_HPR^p
GHPR=prod(Net_HPR_p)
GHPR#1.087
## [1] 1.0867

Page 137, compute beta for each branch in the generated dataset.

beta==1 means not a drawdown branch or actual drawdown is less than the percentage.

#minimum by 1, then prod next one.
cumprodmin<-function(x){
 reg=1
 res=x
 for(i in 1: length(x)){
  reg=min(1,reg*x[i])
  res[i]=reg
 }
 return(res)
}

beta<-function(brch){
 HPR=Net_HPR[brch]
 HPR_Products=cumprodmin(HPR) 
 HPR_Products_b=HPR_Products-b
 ratio=sum(HPR_Products_b)/sum(abs(HPR_Products_b))
 return(as.integer(ratio))
}

#This indicates a 9.8 percent chance of a 20 percent drawdown
#over the next 12 months.
betas=apply(data,1,beta)
RD=1-mean(betas)
RD #risk of drawdown
## [1] 0.09885