MCScanX共线性分析
程序员文章站
2024-03-02 12:52:22
...
Phytozome数据处理
下载gene.gff3和transcript_primaryTranscriptOnly.fa两个文件,接下来数据整理需要保证gff文件中的基因名与fasta文件中的基因名相同。代码对于拟南芥试用,其它植物因为序列命名问题,稍有差异。
GFF数据
Phytozome下载的gene.gff3文件中,mRNA的会注释出多个转录本,gene注释是唯一的,所以我们取出gene的注释进行整理。后来发现有些物种在改名字的时候可能会稍有不同,需要修改一下。
grep '\bgene\b' Athaliana_167_TAIR10.gene.gff3 | awk '{print $1"\t"$9"\t"$4"\t"$5}' | sed 's/ID=.*;Name=//g' > Ath.gff
整理结果如下:
Chr1 AT1G01010 3631 5899
Chr1 AT1G01020 5928 8737
Chr1 AT1G01030 11649 13714
Chr1 AT1G01040 23146 31227
Chr1 AT1G01050 31170 33153
Chr1 AT1G01060 33379 37871
Chr1 AT1G01070 38752 40944
Chr1 AT1G01073 44677 44787
Chr1 AT1G01080 45296 47019
FASTA数据处理
下载transcript_primaryTranscriptOnly.fa文件为每个基因的主要转录本序列,命名时也会带着转录本的编号,比如ATCG00020.1,需要把“.”及后面的内容去掉。
cut -d '.' -f 1 Athaliana_167_TAIR10.transcript_primaryTranscriptOnly.fa > Ath.fa
下面验证Ath.gff3和Ath.fa中序列的名字是否一一对应,没有输出结果表示一一对应
cut -f 2 Ath.gff > lst1; grep '>' Ath.fa | sed 's/>//g' > lst2; cat lst1 lst2 | sort | uniq -c | sed 's/^ *//g' | grep -v '^2'