系列目录:
基因家族分析 | 番茄Nramp基因家族分析(一)
基因家族分析 | 番茄Nramp基因家族分析(二)
以下分析完全参照lxmic的简书内容(基于全基因组的基因家族分析),侵删~
1.准备数据
①下载番茄全基因组CDS序列、protein序列、gff注释文件、基因组fa文件
ITAG3.2_CDS.fa
ITAG3.2_proteins.fa
ITAG3.2_gene_models.gff
S_lycopersicum_chromosomes.3.00.fa
②选定目标基因家族
目标基因家族: The natural resistance-associated macrophage protein (NRAMP)
③下载NRAMP的hmm模型 (https://pfam.xfam.org/)
Nramp.hmm #PF01566
2.利用hmmsearch初步筛选基因家族成员
hmmsearch Nramp.hmm Solanum_lycopersicum.SL2.50.pep.all.fa> out.txt
less -S out.txt
#下图结果可以看到在番茄蛋白序列中找到了9个属于Nramp基因家族的成员
3.批量获取家族成员信息
# 截取id号
vim out.txt
# 获取id号所在的行号,然后再用sed命令截取行,再用grep命令将id号匹配并重定向。
在vim命令模式下,输入“:set nu”
# sed命令截取,并用管道符直接输入给grep,匹配重定向到id文件
sed -n '17,26p' out.txt | grep -o "Sol.*\.1" > id.txt
# 利用samtools工具来进行序列提取# 首先建立索引文件
/data/zhaoqingyuan/biosoft/samtools-1.3.1/samtools faidx Solanum_lycopersicum.SL2.50.pep.all.fa (最近镜像出现问题不能直接调用samtools)
# 再将id好作为输入,之后在重定向
xargs /data/zhaoqingyuan/biosoft/samtools-1.3.1/samtools faidx Solanum_lycopersicum.SL2.50.pep.all.fa < id.txt > nramp_protein_1.fasta &
less nramp_protein_1.fasta
# 得到的序列文件是含有回车符的,我利用一个perl单行命令将fasta格式的多行序列变成单行的fasta格式序列,链接:http://www.biotrainee.com/thread-291-1-1.html
perl -pe '/^>/ ? print "\n" : chomp' nramp_protein_1.fasta | tail -n +2 > nramp_protein.fasta
# 最后在samrt网站确认是否是该家族成员,进行最后的鉴定。
4.CDS序列获取
# Nramp_num文件含有id号
/data/zhaoqingyuan/biosoft/samtools-1.3.1/samtools faidx Solanum_lycopersicum.SL2.50.cds.all.fa
xargs /data/zhaoqingyuan/biosoft/samtools-1.3.1/samtools faidx Solanum_lycopersicum.SL2.50.cds.all.fa < id.txt > Nramp_cds_1.fasta
less Nramp_cds_1.fasta
# 发现CDS序列有回车,分为很多行,用Perl单行命令将其变成一行
perl -pe '/^>/ ? print "\n" : chomp' Nramp_cds_1.fasta | tail -n +2 > Nramp_cds.fasta
less Nramp_cds.fasta
5.Genomic DNA序列获取
# 根据geneid批量获取基因的染色体位置信息文件
grep -f id.txt ITAG3.2_gene_models.gff | cut -f 1,4,5,9 | grep "ID=mRNA" | sed 's/;Name.*//g' | sed 's/ID.*://g' > gene.bed
# 利用bedtools工具,将基因序列提取出来并重定向到新文件
/data/zhaoqingyuan/biosoft/samtools-1.3.1/samtools faidx Solanum_lycopersicum.SL2.50.dna.toplevel.fa
bedtools getfasta -fi S_lycopersicum_chromosomes.3.00.fa -bed gene.bed -name > gene_seq.fasta (bedtools直接安装在环境中)
通过以上分析我们得到了以下数据:
id.txt
Nramp_cds.fasta
nramp_protein.fasta
gene.bed
gene_seq.fasta
接下来利用获得的这些数据进行可视化。
参考:
01 基于全基因组的基因家族分析(1):数据准备
02 基于全基因组的基因家族分析(2):SlNRAMP家族基因成员鉴定
03 基于全基因组的基因家族分析(3):SlNRAMP家族基因CDS和Genomic DNA序列获取