我们拿到bed文件,通常需要按染色体排序后,对于相同染色体的还需进一步按染色体坐标位置排序,要尽量使染色体的排列顺序与bam文件中的染色体顺序一致。这样在使用某些同时处理bam和bed文件的生信工具的处理两文件的时候就不会报错。
Linux自带sort命令,其默认排序顺序是基于字典序排序的,用以下命令,可按染色体的字典序排序,染色体一致的然后再按基因组坐标位置大小排序,命令行(为演示染色体排序效果,本人把染色体信息一列单独提取出来,并取unique,统计了频数)和结果如下:
[xuxiong@~ xiongxu]$ cut -c4- /data2/data/b37/exons.bed|sort -k1,1 -k2,2n -k3,3n |cut -f1|uniq -c
24921 1
10330 10
13883 11
13653 12
4351 13
7886 14
9001 15
11101 16
14487 17
303 17_ctg5_hap1
1 17_gl000204_random
1 17_gl000205_random
3616 18
14408 19
127 19_gl000209_random
12 1_gl000191_random
22 1_gl000192_random
18313 2
6037 20
2506 21
5394 22
13939 3
9339 4
65 4_ctg9_hap1
4 4_gl000193_random
2 4_gl000194_random
10671 5
12478 6
852 6_apd_hap1
2102 6_cox_hap2
2076 6_dbb_hap3
1539 6_mann_hap4
2003 6_mcf_hap5
2012 6_qbl_hap6
1945 6_ssto_hap7
11848 7
4 7_gl000195_random
8333 8
10184 9
10 9_gl000201_random
13 M
1 Un_gl000212
23 Un_gl000213
1 Un_gl000215
4 Un_gl000218
7 Un_gl000219
1 Un_gl000221
11 Un_gl000222
16 Un_gl000223
52 Un_gl000228
2 Un_gl000237
4 Un_gl000242
8837 X
1059 Y
然而我们通常需要将染色体按数字顺序排序,同时把与某条染色体相关的未能组装的contig排在其后面,这样才能尽可能与bam文件中的染色体顺序一致。本人用perl代码实现和结果如下:
[xuxiong@~ xiongxu]$ cut -c4- /data2/data/b37/exons.bed |perl -ne 'chomp;@u=split(/\t/,$_);push @array,[@u];if(eof){print join("\n",map {join("\t",@{$_})} sort {my ($c)=$a->[0]=~/^(\d+)/;my ($d)=$b->[0]=~/^(\d+)/; (defined($c) && defined($d) ? ($c <=> $d or $a->[0] cmp $b->[0]) : ((defined($d) || defined($c) || $a->[0]=~/[XY]/ && $b->[0]=~/[XY]/ || $a->[0]!~/[XY]/ && $b->[0]!~/[XY]/) ? $a->[0] cmp $b->[0] : $b->[0] cmp $a->[0] )) or $a->[1] <=> $b->[1] or $a->[2] <=> $b->[2] or $a->[3] <=> $b->[3]} @array),"\n";}' |cut -f1|uniq -c
24921 1
12 1_gl000191_random
22 1_gl000192_random
18313 2
13939 3
9339 4
65 4_ctg9_hap1
4 4_gl000193_random
2 4_gl000194_random
10671 5
12478 6
852 6_apd_hap1
2102 6_cox_hap2
2076 6_dbb_hap3
1539 6_mann_hap4
2003 6_mcf_hap5
2012 6_qbl_hap6
1945 6_ssto_hap7
11848 7
4 7_gl000195_random
8333 8
10184 9
10 9_gl000201_random
10330 10
13883 11
13653 12
4351 13
7886 14
9001 15
11101 16
14487 17
303 17_ctg5_hap1
1 17_gl000204_random
1 17_gl000205_random
3616 18
14408 19
127 19_gl000209_random
6037 20
2506 21
5394 22
8837 X
1059 Y
13 M
1 Un_gl000212
23 Un_gl000213
1 Un_gl000215
4 Un_gl000218
7 Un_gl000219
1 Un_gl000221
11 Un_gl000222
16 Un_gl000223
52 Un_gl000228
2 Un_gl000237
4 Un_gl000242
如有更好的方法,欢迎讨论
想学生物信息?