黄健
细节决定成败
2023-3-30 13:41
阅读:1896

在即将出版的《Perl生物信息学编程》一书中,有这样一个实例练习:用Perl单行命令来确定仅见于skin的多肽、仅见于muscle的多肽以及两者都有的多肽。两个数据文件skin、muscle见GitHub。书中给出的参考答案如下。

1、仅见于皮肤的多肽:

perl -ne 'if(!$#ARGV){$h{$_}=1} else{print if !$h{$_}}' muscle skin >skin.prl

2、仅见于肌肉的多肽:

perl -ne 'if(!$#ARGV){$h{$_}=1} else{print if !$h{$_}}' skin muscle >muscle.prl

3、两者都有的多肽:

perl -ne 'if(!$#ARGV){$h{$_}=1} else{print if $h{$_}}' skin muscle >both.prl

简要解释:通过-e进入单行模式,执行单引号之间的代码;-n选项相当于while (<>) { ... },会依次读取单引号中单行命令后列出的文件。在单行命令模式中,后面提供的文件名作为命令行参数,存入了@ARGV,而#ARGV则返回了@ARGV数组最后一个元素的下标。由于Perl中数组下标默认从0开始,上述命令行有两个文件名作为参数,因此刚开始时,$#ARGV为1(可在BEGIN中检验)。但是,一旦读入第一个文件,相当于$ARGV=shift@ARGV,这时@ARGV里面就只剩下一个文件了,对应的$#ARGV为0了。上述单行命令巧妙地利用了这一点,如第一个单行命令,当读入第一个文件muscle时,$#ARGV为0,那么!$#ARGV就为真,因此muscle中的每一条多肽被记录到哈希%h的键中,对应的值都是1;而读入第二个文件skin时,$#ARGV为-1,因此!$#ARGV为假,执行else语句,skin的每条多肽作为键代入哈希%h检验,如果不存在,!$h{$_}为真,打印的相应多肽就是skin里有但muscle里没有的多肽。余此类推。

现在来验证一下。用最简单的Linux命令来统计只见于skin的多肽数目,comm -23 skin muscle | wc -l,结果为644条;只见于muscle的多肽数目,comm -13 skin muscle | wc -l,结果为679条;二者都有的多肽,comm -12 skin muscle | wc -l,结果为465条。wc -l skin.prl muscle.prl both.prl,结果分别为644、679、464条。啥?二者都有的多肽,为啥Perl单行命令的结果比comm命令的结果少了一条。comm -12 skin muscle > both.comm;comm both.prl both.comm发现,Perl的结果中少了一条多肽YWRSGGM。tail skin muscle发现,在skin中YWRSGGM在倒数第二行,有换行符;在muscle中,YWRSGGM在倒数第一行,没有换行符,因此造成了上述问题。一个快捷的解决方案是给muscle最后一行加个换行符。或者如下,更可确保无误:

perl -nE 'chomp; if(!$#ARGV){$h{$_}=1} else{say if $h{$_}}' skin muscle | wc -l

转载本文请联系原作者获取授权,同时请注明本文来自黄健科学网博客。

链接地址:https://wap.sciencenet.cn/blog-204973-1382368.html?mobile=1

收藏

分享到:

当前推荐数:0
推荐到博客首页
网友评论0 条评论
确定删除指定的回复吗?
确定删除本博文吗?