hNAD = read.table('data/hNAD',stringsAsFactors = F)
hNAD$distance = 0
acen = read.table('data/acen.bed',stringsAsFactors = F)
for(i in 1:dim(hNAD)[1]){
  chr = hNAD[i,1]
  t = acen[acen$V1==chr,]
  hNAD[i,4] = min(abs(hNAD[i,2]-t[1,2]),abs(hNAD[i,2]-t[1,3]),
                  abs(hNAD[i,3]-t[1,2]),abs(hNAD[i,3]-t[1,3]))
}
B = read.table('data/B_compartment.bed',stringsAsFactors = F)
B = B[,1:3]
B$distance = 0
for(i in 1:dim(B)[1]){
  chr = B[i,1]
  t = acen[acen$V1==chr,]
  B[i,4] = min(abs(B[i,2]-t[1,2]),abs(B[i,2]-t[1,3]),
                  abs(B[i,3]-t[1,2]),abs(B[i,3]-t[1,3]))
}
# Density estimations
denx <- density(hNAD$distance)
deny <- density(B$distance)
# Plot
plot(denx,
     ylim = c(0, max(c(denx$y, deny$y))),
     xlim = c(min(c(denx$x, deny$x)),
              max(c(denx$x, deny$x))),main='')
lines(deny)
polygon(denx, col = rgb(43/255, 168/255, 9/255, alpha = 0.6))
polygon(deny, col = rgb(16/255, 137/255, 114/255, alpha = 0.6))