生成染色体长度文件chrom.sizes

染色体长度文件在多种分析中都需要用到，可以在R中生成，亦可在LINUX中成生，下面以hg19为例。

一、在LINUX中

参考：TwoBit Sequence Archives

# 第一步，将.fa转为.2bit，顾名思义，这是一个二进制文件
$ faToTwoBit hg19.fa hg19.2bit

-rwxrwxrwx 1 eric eric 779M Nov 25 14:24 hg19.2bit
-rwxrwxrwx 1 eric eric 3.0G Sep 27 19:32 hg19.fa

# 第二步，

$ twoBitInfo hg19.2bit stdout | sort -k1,1 -k2,2n > hg19.genome.chrom.sizes
$ wc -l hg19.genome.chrom.sizes
# 93 hg19.genome.chrom.sizes
$ head hg19.genome.chrom.sizes
chr1    249250621
chr10   135534747
chr11   135006516
chr11_gl000202_random   40103
chr12   133851895
chr13   115169878
chr14   107349540
chr15   102531392
chr16   90354753
chr17   81195210

二、在R中

> x <- Seqinfo(genome="hg19")
> x
Seqinfo object with 93 sequences (1 circular) from hg19 genome:
  seqnames       seqlengths isCircular genome
  chr1            249250621      FALSE   hg19
  chr2            243199373      FALSE   hg19
  chr3            198022430      FALSE   hg19
  chr4            191154276      FALSE   hg19
  chr5            180915260      FALSE   hg19
  ...                   ...        ...    ...
  chrUn_gl000245      36651      FALSE   hg19
  chrUn_gl000246      38154      FALSE   hg19
  chrUn_gl000247      36422      FALSE   hg19
  chrUn_gl000248      39786      FALSE   hg19
  chrUn_gl000249      38502      FALSE   hg19

> m <- as.data.frame(x)
> n <- data.frame(rownames(m),m$seqlengths)
> write.table(n,"./hg19.size.all",row.names=FALSE,col.names=FALSE,quote=FALSE) # 注意，去掉行名、列名，以及引号

# 回到LINUX下，对其进行sort，
$ sort -k1,1 -k2,2n hg19.size.all > hg19.all.size
$ wc -l hg19.all.size
# 93 hg19.all.size
$ head hg19.all.size
chr1 249250621
chr10 135534747
chr11 135006516
chr11_gl000202_random 40103
chr12 133851895
chr13 115169878
chr14 107349540
chr15 102531392
chr16 90354753
chr17 81195210

对比一下，这两个文件是一样的，当然，在LINUX中要相对简单一点，但是由于不涉及使用.fa文件，处理时间上会有优势，几分钟就能完成，在LINUX中耗时要长一些。