【shell笔记>文本处理|专项】Linux数据文本处理工具(3)

用Sed进行流编辑

sed命令从文本或者标准输入中每次读入一行数据。

我们先从简单的实例出发,看下该命令怎么将一列中的chrm12,chrom2等转换成chr12chr2的格式。

wangsx@SC-201708020022:~/tmp$ cat chrms.txt
chrom1    3214482    3216968
chrom1    3214234    3216968
chrom1    3213425    3210653

wangsx@SC-201708020022:~/tmp$ sed 's/chrom/chr/' chrms.txt
chr1    3214482    3216968
chr1    3214234    3216968
chr1    3213425    3210653

虽然示例文件处理仅仅只有三行,但我们可以将这种处理方式运用到上G甚至更大的数据文件中,而不用打开整个文件进行处理。并且,可以借助重导向实现对数据处理结果的输出。

sed替换命令采用的格式是

s/pattern/replacement/

sed会自动搜索符合pattern的字符串,然后修改为replacement(我们想要修改后的样子)。一般默认sed只替换第一个匹配的pattern,我们可以通过添加全局标识g将其应用到数据的所有行中。

s/pattern/replacement/g

如果我们想要忽略匹配的大小写,使用i标识

s/pattern/replacement/i

默认sed命令支持基本的POSIX正则表达式(BRE),可以通过-E选项进行拓展(ERE)。很多的Linux命令都这种方式,像常用的grep命令。

再看一个实例,如果我们想把chr1:28647389-28659480这样格式的文字转换为三列,可以使用:

wangsx@SC-201708020022:~/tmp$ echo "chr1:28647389-28659480" | \
>  sed -E 's/^(chr[^:]+):([0-9]+)-([0-9]+)/\1\t\2\t\3/'
chr1    28647389        28659480

我们聚焦在第二个命令sed上。初看杂乱无章,但是从最大的结构看依旧是

s/pattern/replacement/

先看pattern部分,这是由几个简单正则表达式组成的复合体,几个()括起来的字符串可以单独看。第一个匹配chr加上一个非冒号的字符,第二个和第三个都是匹配多个数字。最开始的^表示以chr起始(前面没有字符),各个括号中间的是对应的字符。整体的pattern的目的就是为了找到文本中符合这种模式的字符串,如果只是想把这个模式找出来的话,几个括号可以不用加。显然这几个括号的作用就是将它们划分成多个域,帮助sed进行处理。可以看到replacement部分存在\1,\2,\3,它恰好对应()的顺序。这样我们在中间插入\t制表符,就可以完成我们想要的功能:将原字符串转换为三列。

我本身对字符串并不是非常熟悉,懂一些元字符,可能讲解的不是很到位。不熟悉正则表达式的朋友,可以学习和参考下学习正则表达式,是我从Github上Copy到的非常好的学习资料,有兴趣也可以Fork学习。

上山的路总是有很多条,我们下面看下其他实现该功能的办法:

wangsx@SC-201708020022:~/tmp$ echo "chr1:28647389-28659480" |  sed 's/[:-]/\t/g'
chr1    28647389        28659480
wangsx@SC-201708020022:~/tmp$ echo "chr1:28647389-28659480" |  sed 's/:/\t/' | sed 's/-/\t/'
chr1    28647389        28659480
wangsx@SC-201708020022:~/tmp$ echo "chr1:28647389-28659480" |  tr ':-' '\t'
chr1    28647389        28659480

这三种方式看起来都非常简单有效。它处理字符串的思路不是从匹配pattern然后替换入手,不对,应该说是不是从匹配所有pattern然后替换入手。处理的关键是只处理字符串中看似无用的连字符:-,将其替换成制表符从而轻松完成分割。

sed 's/:/\t/' | sed 's/-/\t/'可以通过-e选项写为sed -e 's/:/\t/' -e 's/-/\t/',效果等价。

默认sed命令支持基本的POSIX正则表达式(BRE),可以通过-E选项进行拓展(ERE)。很多的Linux命令都这种方式,像常用的grep命令。

再看一个实例,如果我们想把chr1:28647389-28659480这样格式的文字转换为三列,可以使用:

wangsx@SC-201708020022:~/tmp$ echo "chr1:28647389-28659480" | \
>  sed -E 's/^(chr[^:]+):([0-9]+)-([0-9]+)/\1\t\2\t\3/'
chr1    28647389        28659480

我们聚焦在第二个命令sed上。初看杂乱无章,但是从最大的结构看依旧是

s/pattern/replacement/

先看pattern部分,这是由几个简单正则表达式组成的复合体,几个()括起来的字符串可以单独看。第一个匹配chr加上一个非冒号的字符,第二个和第三个都是匹配多个数字。最开始的^表示以chr起始(前面没有字符),各个括号中间的是对应的字符。整体的pattern的目的就是为了找到文本中符合这种模式的字符串,如果只是想把这个模式找出来的话,几个括号可以不用加。显然这几个括号的作用就是将它们划分成多个域,帮助sed进行处理。可以看到replacement部分存在\1,\2,\3,它恰好对应()的顺序。这样我们在中间插入\t制表符,就可以完成我们想要的功能:将原字符串转换为三列。

我本身对字符串并不是非常熟悉,懂一些元字符,可能讲解的不是很到位。不熟悉正则表达式的朋友,可以学习和参考下学习正则表达式,是我从Github上Copy到的非常好的学习资料,有兴趣也可以Fork学习。

上山的路总是有很多条,我们下面看下其他实现该功能的办法:

wangsx@SC-201708020022:~/tmp$ echo "chr1:28647389-28659480" |  sed 's/[:-]/\t/g'
chr1    28647389        28659480
wangsx@SC-201708020022:~/tmp$ echo "chr1:28647389-28659480" |  sed 's/:/\t/' | sed 's/-/\t/'
chr1    28647389        28659480
wangsx@SC-201708020022:~/tmp$ echo "chr1:28647389-28659480" |  tr ':-' '\t'
chr1    28647389        28659480

这三种方式看起来都非常简单有效。它处理字符串的思路不是从匹配pattern然后替换入手,不对,应该说是不是从匹配所有pattern然后替换入手。处理的关键是只处理字符串中看似无用的连字符:-,将其替换成制表符从而轻松完成分割。

sed 's/:/\t/' | sed 's/-/\t/'可以通过-e选项写为sed -e 's/:/\t/' -e 's/-/\t/',效果等价。

默认,sed会输出每一行的结果,用replacement替换pattern,但实际中我们可能会因此得到不想要的结果。比如下面的这个例子。

如果我们想要抓出gtf文件第九列的转录名,可能会使用以下命令

wangsx@SC-201708020022:~/database$ zgrep -v "^#" gencode.v27lift37.annotation.gtf.gz | head -n 3 | \
> sed -E 's/.*transcript_id "([^"]+)".*/\1/'
chr1    HAVANA  gene    11869   14409   .       +       .       gene_id "ENSG00000223972.5_2"; gene_type "transcribed_unprocessed_pseudogene"; gene_name "DDX11L1"; level 2; havana_gene "OTTHUMG00000000961.2_2"; remap_status "full_contig"; remap_num_mappings 1; remap_target_status "overlap";
ENST00000456328.2_1
ENST00000456328.2_1

我们可以发现一些没有转录名行的结果是输出整行,这可不是我们想要的。一种解决办法是在使用sed之前先抓出有transcript_id的行。其实sed命令本身也可以通过选项和参数设定解决这个问题,这里我们可以用-n选项关闭sed输出所有行,在最末的/后加p只输出匹配项。

wangsx@SC-201708020022:~/database$ zgrep -v "^#" gencode.v27lift37.annotation.gtf.gz | head -n 3 | sed -E -n 's/.*transc
ript_id "([^"]+)".*/\1/p'
ENST00000456328.2_1
ENST00000456328.2_1

注意方括号内^是非(取反)的意思。

解释如下:

  1. 首先,匹配字符串"transcript_id"之前0或多个任意字符(.表示除换行键的任意字符)。
  2. 然后,匹配和捕获一个或多个不是引号的字符,用的是[^"]

+号的使用是一种非贪婪的方法。很多新手会用*,这是贪婪操作,往往会得不偿失,需要注意喔。

wangsx@SC-201708020022:~/database$ sed -E 's/transcript_id "([^"]+)".*/\1/' greedy_example.txt
ENSMUST00000160944
wangsx@SC-201708020022:~/database$ sed -E 's/transcript_id "(.*)".*/\1/' greedy_example.txt
ENSMUST00000160944"; gene_name "Gm16088

使用*时它会尽量多地去匹配符合要求的模式。

我们也可以用sed命令来获取特定范围的行,比如说我要取出头10行,可以使用

sed -n '1,10p' filename

20到50行

sed -n '20,50p' filename

当然sed的功能特性远远不止这些,有待于大家更多地挖掘。不过需要注意的是,尽量让工具干它最擅长的事情。如果是复杂地大规模计算,还是最好写个Python脚本。

KISS原则:

Keep Incredible Sed Simple

高级Shell用法

子shell

首先需要记住连续命令和管道命令的区别:前者是简单地一个一个按顺序运行程序(一般用&&或者;);后者前一个程序的输出结果会直接传到下一个命令程序的输入中(这不就是流程化操作么,用|分隔)。

子shell可以让我们在一个独立的shell进程中执行连续命令。

首先看个例子

wangsx@SC-201708020022:~/database$ echo "this command"; echo "that command" | sed 's/command/step/'
this command
that step
wangsx@SC-201708020022:~/database$ (echo "this command"; echo "that command") | sed 's/command/step/'
this step
that step

发现仅仅加了个括号,结果就不同了。第二个命令就用了子shell,它把两个echo命令放进单独的空间执行后将结果传给下游。

子shell在对gtf文件进行操作时有个非常有意思有用的用处。我们如果想对gtf文件排序,但是又想要保留文件头部注释信息,我们就能够用两次grep操作分别抓出注释和非注释信息,然后又把它结合在一起。下面看看效果,用less进行检查:

wangsx@SC-201708020022:~/database$ (zgrep "^#" gencode.v27lift37.annotation.gtf.gz; \
>  zgrep -v "^#" gencode.v27lift37.annotation.gtf.gz | \
>  sort -k1,1 -k4,4n) | less

可以看到,子shell确实能够给我们提供非常有用的操作去组合命令实现想要的功能。

命令管道及进程替换

很多生信命令行工具需要提供多个输入和输出参数,这用在管道命令里可能会导致非常低效的情形(管道只接受一个标准输入和输出)。幸好,我们可以使用命令管道来解决此类问题。

命名管道,也成为FIFO(先入先出,额,这不是队列么:smile:)。它是一个特殊的排序文件,命名管道有点像文件,它可以永久保留在你的文件系统上(估计本质就是文件吧~)。

我们用mkfifo来生成它

wangsx@SC-201708020022:~/tmp$ ls -l fqin
prw-rw-rw- 1 wangsx wangsx 0 9月   3 20:45 fqin

可以它看它权限的第一个字符是p,指代是pipe。说明是个特殊文件。

我们像文件一样对它进行一些操作

wangsx@SC-201708020022:~/tmp$ cat fqin
hello, named pipes
[1]+  已完成               echo "hello, named pipes" > fqin
wangsx@SC-201708020022:~/tmp$ rm fqin

比如当使用一个生信命令行工具

processing_tool --in1 in1.fq --in2 in2.fq --out1 out1.fq --out2 out2.fq

in1.fq in2.fq就可以上游输出数据到processing_tool的命名管道;同理out1.fq out2.fq可以是命名管道用来写进输出数据。

但这样我们每次都得不停地创建和删除这些文件,解决办法是使用匿名管道,也叫进程替换。

不能光说,看看例子就知道和理解了。

wangsx@SC-201708020022:~/tmp$ cat <(echo "hello, process substitution")
hello, process substitution

echo命令运行后使用了进程替换,产生匿名文件,然后匿名文件被重导向cat命令。

把它用到工具上,就变成了(假定上游zcat下游执行grep命令)

processing_tool --in1 < (zcat file1) --in2 < (zcat file2) --out1 (gzip > outfile1) --out2 (gzip > outfile2)

关于Linux数据处理工具内容全部整理发布在我的博客上。详情点击

最后编辑于
©著作权归作者所有,转载或内容合作请联系作者
  • 序言:七十年代末,一起剥皮案震惊了整个滨河市,随后出现的几起案子,更是在滨河造成了极大的恐慌,老刑警刘岩,带你破解...
    沈念sama阅读 199,711评论 5 468
  • 序言:滨河连续发生了三起死亡事件,死亡现场离奇诡异,居然都是意外死亡,警方通过查阅死者的电脑和手机,发现死者居然都...
    沈念sama阅读 83,932评论 2 376
  • 文/潘晓璐 我一进店门,熙熙楼的掌柜王于贵愁眉苦脸地迎上来,“玉大人,你说我怎么就摊上这事。” “怎么了?”我有些...
    开封第一讲书人阅读 146,770评论 0 330
  • 文/不坏的土叔 我叫张陵,是天一观的道长。 经常有香客问我,道长,这世上最难降的妖魔是什么? 我笑而不...
    开封第一讲书人阅读 53,799评论 1 271
  • 正文 为了忘掉前任,我火速办了婚礼,结果婚礼上,老公的妹妹穿的比我还像新娘。我一直安慰自己,他们只是感情好,可当我...
    茶点故事阅读 62,697评论 5 359
  • 文/花漫 我一把揭开白布。 她就那样静静地躺着,像睡着了一般。 火红的嫁衣衬着肌肤如雪。 梳的纹丝不乱的头发上,一...
    开封第一讲书人阅读 48,069评论 1 276
  • 那天,我揣着相机与录音,去河边找鬼。 笑死,一个胖子当着我的面吹牛,可吹牛的内容都是我干的。 我是一名探鬼主播,决...
    沈念sama阅读 37,535评论 3 390
  • 文/苍兰香墨 我猛地睁开眼,长吁一口气:“原来是场噩梦啊……” “哼!你这毒妇竟也来了?” 一声冷哼从身侧响起,我...
    开封第一讲书人阅读 36,200评论 0 254
  • 序言:老挝万荣一对情侣失踪,失踪者是张志新(化名)和其女友刘颖,没想到半个月后,有当地人在树林里发现了一具尸体,经...
    沈念sama阅读 40,353评论 1 294
  • 正文 独居荒郊野岭守林人离奇死亡,尸身上长有42处带血的脓包…… 初始之章·张勋 以下内容为张勋视角 年9月15日...
    茶点故事阅读 35,290评论 2 317
  • 正文 我和宋清朗相恋三年,在试婚纱的时候发现自己被绿了。 大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
    茶点故事阅读 37,331评论 1 329
  • 序言:一个原本活蹦乱跳的男人离奇死亡,死状恐怖,灵堂内的尸体忽然破棺而出,到底是诈尸还是另有隐情,我是刑警宁泽,带...
    沈念sama阅读 33,020评论 3 315
  • 正文 年R本政府宣布,位于F岛的核电站,受9级特大地震影响,放射性物质发生泄漏。R本人自食恶果不足惜,却给世界环境...
    茶点故事阅读 38,610评论 3 303
  • 文/蒙蒙 一、第九天 我趴在偏房一处隐蔽的房顶上张望。 院中可真热闹,春花似锦、人声如沸。这庄子的主人今日做“春日...
    开封第一讲书人阅读 29,694评论 0 19
  • 文/苍兰香墨 我抬头看了看天上的太阳。三九已至,却和暖如春,着一层夹袄步出监牢的瞬间,已是汗流浃背。 一阵脚步声响...
    开封第一讲书人阅读 30,927评论 1 255
  • 我被黑心中介骗来泰国打工, 没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留,地道东北人。 一个月前我还...
    沈念sama阅读 42,330评论 2 346
  • 正文 我出身青楼,却偏偏与公主长得像,于是被迫代替她去往敌国和亲。 传闻我的和亲对象是个残疾皇子,可洞房花烛夜当晚...
    茶点故事阅读 41,904评论 2 341

推荐阅读更多精彩内容

  • linux资料总章2.1 1.0写的不好抱歉 但是2.0已经改了很多 但是错误还是无法避免 以后资料会慢慢更新 大...
    数据革命阅读 12,127评论 2 34
  • 本文承接之前写的三十分钟学会AWK一文,在学习完AWK之后,趁热打铁又学习了一下SED,不得不说这两个工具真的堪称...
    mylxsw阅读 4,379评论 3 74
  • 本文笔记源自这里——[实验楼]欢迎大家在下面交流其中有问题的地方喜欢请点收藏,每日更新(全部已亲自实践). 一. ...
    东皇Amrzs阅读 3,953评论 7 54
  • Spring Cloud为开发人员提供了快速构建分布式系统中一些常见模式的工具(例如配置管理,服务发现,断路器,智...
    卡卡罗2017阅读 134,564评论 18 139
  • 知识点 sort uniq cut wc sed命令 awk命令 crontab定时器 sort sort 命令对...