sed 命令从文本或者标准输入中每次读入一行数据。
我们先从简单的实例出发,看下该命令怎么将一列中的 chrm12 , chrom2 等转换成 chr12 , chr2 的格式。
虽然示例文件处理仅仅只有三行,但我们可以将这种处理方式运用到上G甚至更大的数据文件中,而不用打开整个文件进行处理。并且,可以借助重导向实现对数据处理结果的输出。
sed 替换命令采用的格式是
sed 会自动搜索符合 pattern 的字符串,然后修改为 replacement (我们想要修改后的样子)。一般默认 sed 只替换第一个匹配的 pattern ,我们可以通过添加全局标识 g 将其应用到数据的所有行中。
如果我们想要忽略匹配的大小写,使用 i 标识
默认 sed 命令支持基本的POSIX正则表达式(BRE),可以通过 -E 选项进行拓展(ERE)。很多的Linux命令都这种方式,像常用的 grep 命令。
再看一个实例,如果我们想把 chr1:28647389-28659480 这样格式的文字转换为三列,可以使用:
我们聚焦在第二个命令 sed 上。初看杂乱无章,但是从最大的结构看依旧是
先看 pattern 部分,这是由几个简单正则表达式组成的复合体,几个 () 括起来的字符串可以单独看。第一个匹配 chr 加上一个非冒号的字符,第二个和第三个都是匹配多个数字。最开始的 ^ 表示以 chr 起始(前面没有字符),各个括号中间的是对应的字符。整体的 pattern 的目的就是为了找到文本中符合这种模式的字符串,如果只是想把这个模式找出来的话,几个括号可以不用加。显然这几个括号的作用就是将它们划分成多个域,帮助 sed 进行处理。可以看到 replacement 部分存在 \1 , \2 , \3 ,它恰好对应 () 的顺序。这样我们在中间插入 \t 制表符,就可以完成我们想要的功能:将原字符串转换为三列。
我本身对字符串并不是非常熟悉,懂一些元字符,可能讲解的不是很到位。不熟悉正则表达式的朋友,可以学习和参考下 学习正则表达式 ,是我从Github上Copy到的非常好的学习资料,有兴趣也可以Fork学习。
上山的路总是有很多条,我们下面看下其他实现该功能的办法:
这三种方式看起来都非常简单有效。它处理字符串的思路不是从匹配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 这样格式的文字转换为三列,可以使用:
我们聚焦在第二个命令 sed 上。初看杂乱无章,但是从最大的结构看依旧是
先看 pattern 部分,这是由几个简单正则表达式组成的复合体,几个 () 括起来的字符串可以单独看。第一个匹配 chr 加上一个非冒号的字符,第二个和第三个都是匹配多个数字。最开始的 ^ 表示以 chr 起始(前面没有字符),各个括号中间的是对应的字符。整体的 pattern 的目的就是为了找到文本中符合这种模式的字符串,如果只是想把这个模式找出来的话,几个括号可以不用加。显然这几个括号的作用就是将它们划分成多个域,帮助 sed 进行处理。可以看到 replacement 部分存在 \1 , \2 , \3 ,它恰好对应 () 的顺序。这样我们在中间插入 \t 制表符,就可以完成我们想要的功能:将原字符串转换为三列。
我本身对字符串并不是非常熟悉,懂一些元字符,可能讲解的不是很到位。不熟悉正则表达式的朋友,可以学习和参考下 学习正则表达式 ,是我从Github上Copy到的非常好的学习资料,有兴趣也可以Fork学习。
上山的路总是有很多条,我们下面看下其他实现该功能的办法:
这三种方式看起来都非常简单有效。它处理字符串的思路不是从匹配pattern然后替换入手,不对,应该说是不是从匹配所有pattern然后替换入手。处理的关键是只处理字符串中看似无用的连字符 : 与 - ,将其替换成制表符从而轻松完成分割。
sed 's/:/\t/' | sed 's/-/\t/' 可以通过 -e 选项写为 sed -e 's/:/\t/' -e 's/-/\t/' ,效果等价。
默认, sed 会输出每一行的结果,用 replacement 替换 pattern ,但实际中我们可能会因此得到不想要的结果。比如下面的这个例子。
如果我们想要抓出 gtf 文件第九列的转录名,可能会使用以下命令
我们可以发现一些没有转录名行的结果是输出整行,这可不是我们想要的。一种解决办法是在使用 sed 之前先抓出有 transcript_id 的行。其实 sed 命令本身也可以通过选项和参数设定解决这个问题,这里我们可以用 -n 选项关闭 sed 输出所有行,在最末的 / 后加 p 只输出匹配项。
注意方括号内 ^ 是非(取反)的意思。
解释如下:
+ 号的使用是一种非贪婪的方法。很多新手会用 * ,这是贪婪 *** 作,往往会得不偿失,需要注意喔。
使用 * 时它会尽量多地去匹配符合要求的模式。
我们也可以用 sed 命令来获取特定范围的行,比如说我要取出头10行,可以使用
20到50行
当然 sed 的功能特性远远不止这些,有待于大家更多地挖掘。不过需要注意的是,尽量让工具干它最擅长的事情。如果是复杂地大规模计算,还是最好写个Python脚本。
首先需要记住 连续 命令和 管道 命令的区别:前者是简单地一个一个按顺序运行程序(一般用 &&或者 );后者前一个程序的输出结果会直接传到下一个命令程序的输入中(这不就是流程化 *** 作么,用 | 分隔)。
子shell可以让我们在一个独立的shell进程中执行连续命令。
首先看个例子
发现仅仅加了个括号,结果就不同了。第二个命令就用了子shell,它把两个 echo 命令放进单独的空间执行后将结果传给下游。
子shell在对 gtf 文件进行 *** 作时有个非常有意思有用的用处。我们如果想对 gtf 文件排序,但是又想要保留文件头部注释信息,我们就能够用两次 grep *** 作分别抓出注释和非注释信息,然后又把它结合在一起。下面看看效果,用 less 进行检查:
可以看到,子shell确实能够给我们提供非常有用的 *** 作去组合命令实现想要的功能。
很多生信命令行工具需要提供多个输入和输出参数,这用在管道命令里可能会导致非常低效的情形(管道只接受一个标准输入和输出)。幸好,我们可以使用命令管道来解决此类问题。
命名管道 ,也成为FIFO(先入先出,额,这不是队列么:smile:)。它是一个特殊的排序文件,命名管道有点像文件,它可以永久保留在你的文件系统上(估计本质就是文件吧~)。
我们用 mkfifo 来生成它
可以它看它权限的第一个字符是p,指代是pipe。说明是个特殊文件。
我们像文件一样对它进行一些 *** 作
比如当使用一个生信命令行工具
in1.fq in2.fq 就可以上游输出数据到 processing_tool 的命名管道;同理 out1.fq out2.fq 可以是命名管道用来写进输出数据。
但这样我们每次都得不停地创建和删除这些文件,解决办法是使用匿名管道,也叫进程替换。
不能光说,看看例子就知道和理解了。
echo 命令运行后使用了进程替换,产生匿名文件,然后匿名文件被重导向 cat 命令。
把它用到工具上,就变成了(假定上游zcat下游执行grep命令)
关于Linux数据处理工具内容全部整理发布在我的博客上。 详情点击
一、格式介绍(一)gtf 文件。GTF 为General Transfer Format缩写,跟 GFF2格式类似。相信大家做转录组分析时候经常会看到Cufflinks或者Stringtie软件对转录组进行定量与组装会时产生一个gtf文件,里面包含的信息如下: 每列信息的含义如下:seqname - 序列的ID,可以是染色体的ID也可以是Scaffold或者Contig的ID。source - 产生此文件的软件,如Stringtie产生的则为Stringtie,CUfflinks产生的则为Cufflinks,不知道的使用点 “.” 表示。feature - 后面start和end之间区域代表的特征,如果此区域是基因,则此处为gene,如果是外显子,则为exon,如果是转录本,则为transcript,如果是非编码RNA则为lncRNA,如果是重复序列,则为TE,等等,主要表明这一块区域的特征。start -上述feature的在序列上的起始位置。end - 上述feature的在序列上的终止位置。score - 一个浮点数值,也可以为点 “.” 。有值的时候代表上述feature的可靠 性。因为无论是gene还是mRNA,都是基于预测差生的,因而必然会有一个值来衡量预测准确性。strand - + (forward)或者 - (reverse),代表上述feature是位于正链还是负链上。frame - 内含子相位,只能为'0', '1' or '2',或者为点 “.”。 '0' 代表feature起始碱基为三联体密码子的第一个碱基, '1' 代表三联体密码子的第2个碱基, 2代表第3个碱基。attribute -备注列。主要备注该feature的一些信息,常见的是gene或者transcript等的ID信息以及FPKM值等,多个备注信息之间通常用分号分隔。 (二)gff 格式。为general feature format缩写,目前采用的是version 3,即我们常说的gff3文件。该文件常用来对基因组进行注释,表示基因,外显子,CDS,UTR等在基因组上的位置。众多基因预测软件如Glean,EVM,AUGUSTUS等会产生此格式文件。 与gtf文件不同之处只是在第9列。此列格式为“标签=值”(tag=value),标签与值之间用“=”,不同的标签之间用“;”隔开,一个标签可以有多个值,不同值用“,”分割。二、gtf与gff转换以及对GFF文件进行过滤。常采用的软件是gffread,为Cufflinks自带的一个程序,他不仅可以实现GTF与GFF的互相转换,而且还可以对GFF文件进行过滤处理。以下是gffread的帮助信息: Usage: gffread <input_gff>[-g <genomic_seqs_fasta>| <dir>][-s <seq_info.fsize>] [-o <outfile.gff>] [-t <tname>] [-r [[<strand>]<chr>:]<start>..<end>[-R]] [-CTVNJMKQAFGUBHZWTOLE] [-w <exons.fa>] [-x <cds.fa>] [-y <tr_cds.fa>] [-i <maxintron>] <input_gffmatch>为一个GFF/GTF文件,必填的一个文件 常用参数介绍: -g 序列文件,即GFF/GTF文件第一列ID对应的序列文件。 -i 丢弃掉内含子大于的转录本(mRNA/transcript) -r 起始和终止位置,填写示例100.10000即为输出与100到10000有重叠的所有转录组,也可以限制序列ID及链,填写示例:+Chr1:100..10000。 -R 丢弃掉此范围的转录本,与-r相反。 -U 丢弃掉 single-exon的转录本 -C 丢低调无CDS的转录本。 -V 丢弃掉含有移码突变的转录本。 -H 如果使用了-V,则重新检查并调整内含子相位,避免由于翻译起始位点选择的位置不对导致移码突变的产生。 -B 如果使用了-V, 对于单外显子基因,则重新检查相反的链,是否存在移码突变。 -N 丢弃掉多外显子基因剪接位点不是常见的 GT-AG, GC-AG or AT-AC序列。 -J 丢弃掉没有起始密码子或者终止密码子的转录本,仅保留有完整编码框的转录本。 --no-pseudo:过滤掉含有 'pseudo' 的注释信息 -M/--merge : 合并完全相同的或者存在包含关系的转录本。 -d:使用 -M ,将合并信息输出到文件中 --cluster-only: 类似于 --merge 但是不合并转录本 -K 对于-M 选项:also collapse shorter, fully contained transcripts with fewer introns than the container -Q 对于-M 选项:移除包含关系的转录本的限制条件:多外显子转录本将会合并,如果他们内含子位置完全一样,单外显子转录本只需要有80%一样即可合并。 --force-exons: 使GFF features的最小层级为exon -E 对于重复的 ID或者 GFF/GTF 其他潜在的格式问题给出警告信息。 -Z 将内含子小于4 bp的邻近的两个外显子合并为一个。 -w 输出每个转录本的外显子序列 -x 输出CDS序列 -W 对于 -w 和 -x 选项,输出外显子位置坐标到FASTA序列的ID中 -y 输出蛋白质序列 -L 将Ensembl GTF 转换为 GFF3 conversion (implies -Fshould be used with -m) -o 输出"filtered" 后的GFF文件 。 -T -o 参数将输出 GTF格式。 示例命令: 1.GFF转换GTF gffread input.gff3 -T -o out.gtf‘ 2.GTF转换GFF3 gffread input.gtf -o out.gff3 3.根据GFF或者GTF提取蛋白质,CDS和外显子序列 gffread gene.gff3 -g genome.fa -x cds.fa -y pep.fa -w cdna.fa三、GFF文件比较主要采用gffcompare(https://github.com/gpertea/gffcompare),其主要具有三个功能:1)评估Cufflinks/Stringtie等转录本组装软件的准确性;2)合并多个GFF/GTF中重叠的部分(多个样本组装结果的合并)3)可以对一个或多个GTF/GFF文件的注释相对于参考的GTF/GFF文件进行分类(with "class codes" assigned to transcripts as per their relationship with the matching/overlapping reference transcript),如Pacbio预测的GTF与参考GFF比较,修正和评估参考的注释结果。 Usage: gffcompare [-r <reference_mrna.gtf>[-R]] [-G] [-T] [-V] [-s <seq_path>] [-o <outprefix>] [-p <cprefix>] {-i <input_gtf_list>| <input1.gtf>[<input2.gtf>.. <inputN.gtf>]} 常用参数介绍: -i 多个GTF 文件时,使用此选项较方便,将多个GTF文件写在一个文件中,通过此选项传入即可。 -r 参考的 GTF/GFF文件 -R 针对的是-r参数,仅考虑参考与任何输入的注释文件有重叠的 。 -Q 针对的是-r参数,仅考虑输入的注释文件与任何参考有重叠的 。 (警告,这将丢弃所有的新的注释位点) -M 丢弃(忽略)掉输入的注释文件和参考注释文件中单外显子转录本 -N 丢弃(忽略)掉参考注释文件中单外显子转录本 -s 基因组序列文件 -e 当评估外显子准确性时,离参考末端外显子最远的距离(默认100) -d 转录本聚类时起始位点相差的最大距离 (默认100) -C 在.combined.gtf文件中包含 "contained" 类型的转录本 -F 如果仅是3’端不同,则不丢弃输入的GTF文件中被参考包含的冗余的转录本注释信息。 -G 不丢弃输入的GTF文件中被参考包含的冗余的转录本注释信息,主要是鉴于可变剪接。 -T 对于每一个输入文件不产生 .tmap 和 .refmap文件 -V 给出 GFF 解析时的警告信息 参考命令: gffcompare -r refChr.gff3 -R -G -o combine input.gtf 输出结果中有以下几个文件: combine.combined.gtf combine.loci combine.stats combine.tracking 其中在combine.combined.gtf中有一个class_code 代表输入的注释文件与参考注释文件相似性信息,具体如下: #Transfrag class codes PriorityCodeDescription 1=Complete match of intron chain 2cContained 3jPotentially novel isoform (fragment): at least one splice junction is shared with a reference transcript 4eSingle exon transfrag overlapping a reference exon and at least 10 bp of a reference intron, indicating a possible pre-mRNA fragment. 5iA transfrag falling entirely within a reference intron 6oGeneric exonic overlap with a reference transcript 7pPossible polymerase run-on fragment (within 2Kbases of a reference transcript) 8rRepeat. Currently determined by looking at the soft-masked reference sequence and applied to transcripts where at least 50% of the bases are lower case 9uUnknown, intergenic transcript 10xExonic overlap with reference on the opposite strand 11sAn intron of the transfrag overlaps a reference intron on the opposite strand (likely due to read mapping errors) 12.(.tracking file only, indicates multiple classifications) 由于输出文件几乎跟cuffcompare格式几乎是一样的, 更详细输出介绍参见http://cole-trapnell-lab.github.io/cufflinks/cuffcompare/。转自:http://www.huoyunjn.com/wuliuxinwen/2/33709819.htmgtf.gz文件是一种压缩文件,以.gz或者.tar.gz(.tgz)为扩展名,它的处理是在Linux和OSX系统里面都可以直接解压使用这种压缩文件。在Windows下常用压缩软件WinRAR打开gz文件,它相当于常见的RAR和ZIP格式。拓展资料: 1、简单的说,就是经过压缩软件压缩的文件叫压缩文件,压缩的原理是把文件的二进制代码压缩,把相邻的0,1代码减少,比如有000000,可以把它变成6个0 的写法60,来减少该文件的空间。2、压缩文件的基本原理是查找文件内的重复字节,并建立一个相同字节的"词典"文件,并用一个代码表示,比如在文件里有几处有一个相同的词"中华人民共和国"用一个代码表示并写入"词典"文件,这样就可以达到缩小文件的目的。欢迎分享,转载请注明来源:内存溢出
评论列表(0条)