1.找到你所感兴趣的基因家族
番茄(Solanum lycopersicum),最喜爱的蔬菜水果之一。摘录维基百科最基本的介绍,详细了解番茄的起源,自行Google。小编还是喜欢Transporter gene family,就觉得特别有意思。植物对于各种营养元素的吸收,都需要其帮助,一旦缺少了,轻则营养不良,重则一命呜呼。本次流程,我选择了The natural resistance-associated macrophage protein (NRAMP)家族。
The tomato (see pronunciation) is the edible, often red, fruit of the plant Solanum lycopersicum, commonly known as a tomato plant. The plant belongs to the nightshade family, Solanaceae.
2.获取基因家族pfam number
pfam主页
KEYWORD SEARCH灰色状态
搜索结果
NRAMP家族信息界面
model界面
3. 利用hmmsearch进行基因家族初步筛选
hmmsearch Nramp.hmm protein.fa > out
,一般我只用到这么简单的语法。Nramp.hmm 是上一步下载到的文件 protein.fa是番茄全基因组蛋白序列文件 out是重定向的输出的文件
Nramp.hmm文件
out输出文件的内容
大致思路:首先从out输出文件的内容中,将其中的geneID截取下来,然后再根据ID号将蛋白序列从protein.fa文件中获取所有家族成员。
代码如下:
# 截取id号
vim out
# 获取id号所在的行号,然后再用sed命令截取行,再用grep命令将id号匹配并重定向。
在vim命令模式下,输入“:set nu”
# sed命令截取,并用管道符直接输入给grep,匹配重定向到id文件
sed -n '17,26p' out | grep -o "Sol.*\.1" > id
# 利用samtools工具来进行序列提取
# 首先建立索引文件
samtools faidx protein.fa
# 再将id好作为输入,之后在重定向
# 参考链接:https://www.biostars.org/p/49820/
xargs samtools faidx protein.fa < id > nramp_protein
less nramp_protein
# 得到的序列文件是含有回车符的,我利用一个perl单行命令将fasta格式的多行序列变成单行的fasta格式序列,链接:http://www.biotrainee.com/thread-291-1-1.html
perl -pe '/^>/ ? print "\n" : chomp' in.fasta | tail -n +2 > out.fasta
# 最后在samrt网站确认是否是该家族成员,进行最后的鉴定。链接:http://smart.embl.de/smart/set_mode.cgi?NORMAL=1
行号显示
会多这样一个fai后缀的索引文件
含有NRAMP结构域的基因