感興趣的基因信息包含在bedGraph文件中，下面命令是對其文件格式進(jìn)行轉(zhuǎn)換，一般進(jìn)行到bam文件可視化的效果比較好。
1. bedGraph轉(zhuǎn)bed文件
BedGraph ，的數(shù)據(jù)和bed文件很類似，ChIPseq數(shù)據(jù)做完peak calling后的bed文件最短只有三列，染色體序號，染色體起始位置和結(jié)束位置。如下所示，前面的聲明和Wig類似，后面的四列分別表示染色體序號，起始位置，結(jié)束位置和value值。相當(dāng)于為bed文件的延伸格式。
track type=bedGraph name="BedGraph Format" description="BedGraph format" visibility=full color=200,100,0 altColor=0,100,200 priority=20
chr19 49302000 49302300 -1.0
chr19 49302300 49302600 -0.75
chr19 49302600 49302900 -0.50
chr19 49302900 49303200 -0.25
chr19 49303200 49303500 0.0
chr19 49303500 49303800 0.25
chr19 49303800 49304100 0.50
chr19 49304100 49304400 0.75
chr19 49304400 49304700 1.00
所以我們想要得到bed文件只需要提取bedGraph的前三列即可，同時注意不要第一行，利用grep -v命令
# Convert bedGraph to bed file
grep -v track GSM1252087_edm2-4_RNAseq.bedGraph | cut -f 1-3 > GSM1252087_edm2-4_RNAseq.bed
2. 建立genome size文件
genome size文件是為了最后一步轉(zhuǎn)化為bam文件所必須的，samtools可以很簡單的建立index文件
# Build genome index file
samtools faidx Arabidopsis_thaliana.TAIR10.dna.toplevel.fa
生成的索引文件genome.fasta.fai,是一個文本文件，分成了5列。第一列是子序列的名稱；第二列是子序列的長度；個人認(rèn)為“第三列是序列所在的位置”，因為該數(shù)字從上往下逐漸變大，最后的數(shù)字是genome.fasta文件的大??；第4和5列不知是啥意思。于是通過此文件，可以定位子序列在fasta文件在磁盤上的存放位置，直接快速調(diào)出子序列。
genome size文件包含index文件的前兩行，也就是chromosome信息和子序列的長度，所以我們可以提取作為genome size。
# Build the genome size file
awk {'print "Chr"$1,"\t",$2'} Arabidopsis_thaliana.TAIR10.dna.toplevel.fa.fai > Arabidopsis_genomeFile.txt
3. bed轉(zhuǎn)bam
bedtools工具提供的bedtobam命令
# Change the bed file to bam file
cat GSM1252087_edm2-4_RNAseq.bed | awk '{x++; printf "%s\tread%d\n",$0,x}' | bedtools bedtobam -g Arabidopsis_genomeFile.txt -i - > GSM1252087_edm2-4_RNAseq.bam
