→_→用到的命令&脚本
grep、awk、sed、cut,faSomeRecords(此脚本可以网上下载)
> 任务
:从拟南芥总蛋白质序列中抽取指定ID的蛋白序列
> 问题
目标ID与蛋白序列fa文件中的ID不完全相符。
目标ID的list文件如下:
拟南芥总蛋白序列的格式如下:
#观察两个文件可以发现,list文件ID不含.1、.2等后缀,且list文件中ID存在小写格式。直接调用faSomeRecords脚本抽取蛋白序列将无法正确匹配。
> 解决
第一步
用grep模糊匹配,从拟南芥总蛋白序列fa文件中抓取含相应ID的行。
$cat 1.txt |xargs -i grep -i "{}" xx.fa > 1.names
#grep -i 表示忽略大小写
#cat显示list文件内容,xargs -i 则表示将上一个命令的内容逐行输入,{ }用以指代输入的内容。
执行第一步得到如下结果:
第二步:
观察文件,发现pep在每行均存在,可将pep定义为列分隔符,输出第一列。
$cut -d 'p' -f 1 1.names > 1.hits 或
$awk -F 'pep' '{print $1}' 1.names > 1.hits
#awk -F ' ' 定义分列符
执行第二步得到如下结果:
第三步:
去掉>符号
$sed -i 's/>//g' 1.hits
#将1.hits文件中的>替换为空。
@至此,list文件修改完毕,下面直接正常使用 faSomeRecords脚本。
第四步:
抽取出目标序列
$./faSomeRecords xx.fa 1.hits hits.fa
#faSomeRecords要赋予可执行权限(chmod u+x faSomeRecords)
#hits.fa为抽取得到的目标序列。
###
后记
不甚智能(foolish)的方法,但解了我一时之急,权且记下。供琢磨改进~
<<<新手录