HMRs

Calculate the amount overlapped between H1 and IMR90 HMRs.

cd ~/

bsub -q 16G -o stdout -e stderr "cat /work3/NRPB1219/hg18_h1_hmr.bed /work3/NRPB1219/hg18_imr90_hmr.bed | sortBed | mergeBed -i - | intersectBed -a - -b /work3/NRPB1219/hg18_h1_hmr.bed -wao | intersectBed -a - -b /work3/NRPB1219/hg18_imr90_hmr.bed -wao | awk -F $'\t' 'BEGIN { OFS=FS } { print \$1,\$2,\$3,\$6-\$5,sprintf(\"%.2f\",\$10/(\$3-\$2)),\$13-\$12,sprintf(\"%.2f\",\$17/(\$3-\$2)) }' > Output/h1_imr90_hmr_coverage.bed"

The first part of the full command line above produces the merged HMR BED output. We use cat to join the HMRs from H1 and IMR90, sort with sortBed and merge feature regions using mergeBed.

cat /work3/NRPB1219/hg18_h1_hmr.bed /work3/NRPB1219/hg18_imr90_hmr.bed | sortBed | mergeBed -i - | head

# Console output

chr1    552721  558206
chr1    558242  562180
chr1    703351  704946
chr1    714082  714363
chr1    739503  744016
chr1    751866  753677
chr1    783310  784282
chr1    794854  798791
chr1    800933  803634
chr1    829069  831267

Then, we pipe the output from the first part of the command to two consecutive intersectBed to identify the overlapped HMRs between H1 and merged output, and also IMR90 and merged output.

# Console output

# intersectBed
# col1~3: merged HMR; col4~10: H1; col11~17: IMR90
chr1    552721  558206  chr1    554333  558188  HYPO0   21      +       3855    chr1    552721  558206  HYPO0   16      +       5485
chr1    558242  562180  chr1    558271  560164  HYPO1   24      +       1893    chr1    558242  562180  HYPO1   13      +       3938
chr1    703351  704946  chr1    703534  704946  HYPO2   23      +       1412    chr1    703351  704626  HYPO2   39      +       1275
chr1    714082  714363  chr1    714082  714363  HYPO3   39      +       281     .       -1      -1      .       -1      .       0
chr1    739503  744016  .       -1      -1      .       -1      .       0       chr1    739503  744016  HYPO3   7       +       4513
chr1    751866  753677  chr1    751975  753029  HYPO4   59      +       1054    chr1    751866  753677  HYPO4   36      +       1811
chr1    783310  784282  .       -1      -1      .       -1      .       0       chr1    783310  784282  HYPO5   23      +       972
chr1    794854  798791  chr1    794854  795385  HYPO5   69      +       531     chr1    794854  798791  HYPO6   17      +       3937
chr1    800933  803634  chr1    801137  802542  HYPO6   23      +       1405    chr1    800933  803634  HYPO7   15      +       2701
chr1    829069  831267  chr1    829069  831267  HYPO7   76      +       2198    chr1    829437  830716  HYPO8   106     +       1279

Use head to print the first 10 lines of the file to standard output..

head ~/Output/h1_imr90_hmr_coverage.bed

# Console output

# col4 = H1 HMR; col5 = H1 HMR coverage; col6 = IMR90 HMR; col7 = IMR90 HMR coverage
chr1    552721  558206  3855    0.70    5485    1.00
chr1    558242  562180  1893    0.48    3938    1.00
chr1    703351  704946  1412    0.89    1275    0.80
chr1    714082  714363  281     1.00    0       0.00
chr1    739503  744016  0       0.00    4513    1.00
chr1    751866  753677  1054    0.58    1811    1.00
chr1    783310  784282  0       0.00    972     1.00
chr1    794854  798791  531     0.13    3937    1.00
chr1    800933  803634  1405    0.52    2701    1.00
chr1    829069  831267  2198    1.00    1279    0.58