環(huán)境變量的補(bǔ)充
PATH只是眾多環(huán)境變量中的一個變量,用于存儲可執(zhí)行文件所在的目錄,以便在用戶輸入命令時可以查詢的到。尤其是自己寫的腳本或安裝的程序,系統(tǒng)不會知道它們在哪個路徑下,需要我們?nèi)ヌ峁┙o系統(tǒng)這些新的路徑,學(xué)名叫設(shè)置環(huán)境變量。
此外常用到的環(huán)境變量還有LD_LIBARY_PATH: 指定動態(tài)鏈接庫 (so文件)的位置,一般在安裝軟件出錯時會用到;PYTHONPATH: 指定Python的安裝包的路徑;PERL5LIB: 指定perl的安裝包的路徑。
設(shè)置環(huán)境變量要注意2點(diǎn):1. 設(shè)置新的環(huán)境變量時一般要包含原始的環(huán)境變量,不能覆蓋;2. 注意自己的目錄和系統(tǒng)環(huán)境變量的目錄的順序,想讓哪個先被找到,就先放哪個。
文件排序
seq: 產(chǎn)生一系列的數(shù)字; man seq查看其具體使用。我們這使用seq產(chǎn)生下游分析所用到的輸入文件。
# 產(chǎn)生從1到10的數(shù),步長為1ct@ehbio:~$ seq 1 1012345678910# 產(chǎn)生從1到10的數(shù),步長為1,用空格分割ct@ehbio:~$ seq -s ‘ ‘ 1 101 2 3 4 5 6 7 8 9 10# 產(chǎn)生從1到10的數(shù),步長為2# 如果有3個數(shù),中間的數(shù)為步長,最后一個始終為最大值ct@ehbio:~$ seq -s ‘ ‘ 1 2 101 3 5 7 9# 還記得前面提到的標(biāo)準(zhǔn)輸入和標(biāo)準(zhǔn)輸出吧# 后臺回復(fù) 標(biāo)準(zhǔn)輸入 查看ct@ehbio:~$ cat <(seq 0 3 17) testct@ehbio:~$ cat test 036912153915
sort: 排序,默認(rèn)按字符編碼排序。如果想按數(shù)字大小排序,需添加-n參數(shù)。
# 可能不符合預(yù)期的排序,系統(tǒng)首先排0,然后排1, 3, 6, 9ct@ehbio:~$ sort test012151533699# 按數(shù)字大小排序ct@ehbio:~$ sort -n test033699121515
sort -u: 去除重復(fù)的行,等同于sort | uniq。
ct@ehbio:~$ sort -nu test03691215
sort file | uniq -d: 獲得重復(fù)的行。(d=duplication)
ct@ehbio:~$ sort -n test | uniq -d3915
sort file | uniq -c: 獲得每行重復(fù)的次數(shù)。
# 第一列為每行出現(xiàn)的次數(shù),第二列為原始的行ct@ehbio:~$ sort -n test | uniq -c ?1 0 ?2 3 ?1 6 ?2 9 ?1 12 ?2 15# 換一個文件看的更清楚ct@ehbio:~$ cat a> b> c> b> a> e> d> a> END# 第一列為每行出現(xiàn)的次數(shù),第二列為原始的行ct@ehbio:~$ sort test2 | uniq -c ? ? ?3 a ? ? ?2 b ? ? ?1 c ? ? ?1 d ? ? ?1 e# 在執(zhí)行uniq操作前,文件要先排序,不然結(jié)果很詭異ct@ehbio:~$ cat test2 | uniq -c ? ? ?1 a ? ? ?1 b ? ? ?1 c ? ? ?1 b ? ? ?1 a ? ? ?1 e ? ? ?1 d ? ? ?1 a
整理下uniq -c的結(jié)果,使得原始行在前,每行的計(jì)數(shù)在后。
awk是一個強(qiáng)大的文本處理工具,其處理數(shù)據(jù)模式為按行處理。每次讀入一行,進(jìn)行操作。OFS: 輸出文件的列分隔符 (output file column separtor);FS為輸入文件的列分隔符 (默認(rèn)為空白字符)。awk中的列從第1到n列,分別記錄為$1, $2 … $n。BEGIN表示在文件讀取前先設(shè)置基本參數(shù);與之相對應(yīng)的是END,只文件讀取完成之后進(jìn)行操作。不以BEGIN, END開頭的{}就是文件讀取、處理的部分。
# 管道符還記得吧,后臺回復(fù) 管道 可查看# awk的操作就是鍍金上一步的結(jié)果,去除多余的空白,然后調(diào)換2列ct@ehbio:~$ sort test2 | uniq -c | awk ‘BEGIN{OFS=””;}{print $2, $1}’a ? ?3b ? ?2c ? ?1d ? ?1e ? ?1
對兩列文件,安照第二列進(jìn)行排序, sort -k2,2n。
# 第二列按數(shù)值大小排序ct@ehbio:~$ sort test2 | uniq -c | awk ‘BEGIN{OFS=””;}{print $2, $1}’ | sort -k2, 2nc ? ?1d ? ?1e ? ?1b ? ?2a ? ?3# 第二列按數(shù)值大小排序# 第二列相同的再按第一列的字母順序的逆序排序 (-r)# 注意看前3行的順序與上一步結(jié)果的差異ct@ehbio:~$ sort test2 | uniq -c | awk ‘BEGIN{OFS=””;}{print $2,$1}’ | sort -k2,2n -k1,1re ? ?1d ? ?1c ? ?1b ? ?2a ? ?3
FASTA序列提取
生成單行序列FASTA文件,提取特定基因的序列,最簡單的是使用grep命令。
grep在前面也提到過,以后還會經(jīng)常提到,主要用途是匹配文件中的字符串,以此為基礎(chǔ),進(jìn)行一系列的操作。如果會使用正則表達(dá)式,將會非常強(qiáng)大。正則表達(dá)式版本很多,幾乎每種語言都有自己的規(guī)則,本文檔不會展開,用到哪個提哪個。
# 生成單行序列FASTA文件ct@ehbio:~$ cat >SOX2> ACGAGGGACGCATCGGACGACTGCAGGACTGTC> >POU5F1> ACGAGGGACGCATCGGACGACTGCAGGACTGTC> >NANOG> CGGAAGGTAGTCGTCAGTGCAGCGAGTCCGT> ENDct@ehbio:~$ cat test.fasta >SOX2ACGAGGGACGCATCGGACGACTGCAGGACTGTC>POU5F1ACGAGGGACGCATCGGACGACTGCAGGACTGTC>NANOGCGGAAGGTAGTCGTCAGTGCAGCGAGTCCGT# grep匹配含有SOX2的行# -A 1 表示輸出的行中,包含匹配行的下一行 (A: after)ct@ehbio:~$ grep -A 1 ‘SOX2’ test.fasta >SOX2ACGAGGGACGCATCGGACGACTGCAGGACTGTC# 也可以使用AWK# 先判斷當(dāng)前行是不是 > 開頭,如果是,表示是序列名字行,替換掉大于號,取出名字。# sub 替換, sub(被替換的部分,要替換成的,待替換字符串)# 如果不以大于號開頭,則為序列行,存儲起來。# seq[name]: 相當(dāng)于建一個字典,name為key,序列為值。然后就可以使用name調(diào)取序列。ct@ehbio:~$ awk ‘BEGIN{OFS=FS=””}{if($0~/>/) {name=$0; sub(“>”, “”, name);} else seq[name]=$0;}END{print “>SOX2”; print seq[“SOX2”]}’ test.fasta>SOX2ACGAGGGACGCATCGGACGACTGCAGGACTGTC
多行FASTA序列提取要麻煩些,一個辦法就是轉(zhuǎn)成單行序列,用上面的方式處理。
sed和tr都為最常用的字符替換工具。
ct@ehbio:~$ cat >SOX2> ACGAGGGACGCATCGGACGACTGCAGGACTGTC> ACGAGGGACGCATCGGACGACTGCAGGACTGTC> ACGAGGGACGCATCGGACGACTGCAGGAC> >POU5F1> CGGAAGGTAGTCGTCAGTGCAGCGAGTCCGT> CGGAAGGTAGTCGTCAGTGCAGCGAGTCC> >NANOG> ACGAGGGACGCATCGGACGACTGCAGGACTGTC> ACGAGGGACGCATCGGACGACTGCAGG> ACGAGGGACGCATCGGACGACTGCAGGACTGTC> ACGAGGGACGCATCGGACGACTGCAGGACTGT> END# 給>號開頭的行的行尾加個TAB鍵,以便隔開名字和序列# TAB鍵不可見,直接看看不大# ()表示記錄匹配的內(nèi)容,1則表示()中記錄的匹配的內(nèi)容# 后面我們專門講sedct@ehbio:~$ sed ‘s/^(>.*)/1/’ test.fasta >SOX2 ? ?ACGAGGGACGCATCGGACGACTGCAGGACTGTCACGAGGGACGCATCGGACGACTGCAGGACTGTCACGAGGGACGCATCGGACGACTGCAGGAC>POU5F1 ? ?CGGAAGGTAGTCGTCAGTGCAGCGAGTCCGTCGGAAGGTAGTCGTCAGTGCAGCGAGTCC>NANOG ? ?ACGAGGGACGCATCGGACGACTGCAGGACTGTCACGAGGGACGCATCGGACGACTGCAGGACGAGGGACGCATCGGACGACTGCAGGACTGTCACGAGGGACGCATCGGACGACTGCAGGACTGT#使用cat -A 可以顯示文件中所有的符號# ^I 表示tab鍵# $表示行尾ct@ehbio:~$ sed ‘s/^(>.*)/1/’ test.fasta | cat -A>SOX2^I$ACGAGGGACGCATCGGACGACTGCAGGACTGTC$ACGAGGGACGCATCGGACGACTGCAGGACTGTC$ACGAGGGACGCATCGGACGACTGCAGGAC$>POU5F1^I$CGGAAGGTAGTCGTCAGTGCAGCGAGTCCGT$CGGAAGGTAGTCGTCAGTGCAGCGAGTCC$>NANOG^I$ACGAGGGACGCATCGGACGACTGCAGGACTGTC$ACGAGGGACGCATCGGACGACTGCAGG$ACGAGGGACGCATCGGACGACTGCAGGACTGTC$ACGAGGGACGCATCGGACGACTGCAGGACTGT$# 把所有的換行符替換為空格# tr這個命令,前面提到過,若想不起來 `man tr`查看# 主意第二個參數(shù),引號內(nèi)為空格。ct@ehbio:~$ sed ‘s/^(>.*)/1/’ test.fasta | tr ” ‘ ‘>SOX2 ? ? ACGAGGGACGCATCGGACGACTGCAGGACTGTC ACGAGGGACGCATCGGACGACTGCAGGACTGTC ACGAGGGACGCATCGGACGACTGCAGGAC >POU5F1 ? ? CGGAAGGTAGTCGTCAGTGCAGCGAGTCCGT CGGAAGGTAGTCGTCAGTGCAGCGAGTCC >NANOG ? ? ACGAGGGACGCATCGGACGACTGCAGGACTGTC ACGAGGGACGCATCGGACGACTGCAGG ACGAGGGACGCATCGGACGACTGCAGGACTGTC ACGAGGGACGCATCGGACGACTGCAGGACTGT # 把最后一個空格替換為換行符ct@ehbio:~$ sed ‘s/^(>.*)/1/’ test.fasta | tr ” ‘ ‘ | sed -e ‘s/ $//’>SOX2 ? ? ACGAGGGACGCATCGGACGACTGCAGGACTGTC ACGAGGGACGCATCGGACGACTGCAGGACTGTC ACGAGGGACGCATCGGACGACTGCAGGAC >POU5F1 ? ? CGGAAGGTAGTCGTCAGTGCAGCGAGTCCGT CGGAAGGTAGTCGTCAGTGCAGCGAGTCC >NANOG ? ? ACGAGGGACGCATCGGACGACTGCAGGACTGTC ACGAGGGACGCATCGGACGACTGCAGG ACGAGGGACGCATCGGACGACTGCAGGACTGTC ACGAGGGACGCATCGGACGACTGCAGGACTGT# 把 ?‘ >’替換為換行符 注意被替換的是 空格+大于號# 當(dāng)連用多個替換命令時,使用-e 隔開ct@ehbio:~$ sed ‘s/^(>.*)/1/’ test.fasta | tr ” ‘ ‘ | sed -e ‘s/ $//’ -e ‘s/ >/>/g’>SOX2 ? ? ACGAGGGACGCATCGGACGACTGCAGGACTGTC ACGAGGGACGCATCGGACGACTGCAGGACTGTC ACGAGGGACGCATCGGACGACTGCAGGAC>POU5F1 ? ? CGGAAGGTAGTCGTCAGTGCAGCGAGTCCGT CGGAAGGTAGTCGTCAGTGCAGCGAGTCC>NANOG ? ? ACGAGGGACGCATCGGACGACTGCAGGACTGTC ACGAGGGACGCATCGGACGACTGCAGG ACGAGGGACGCATCGGACGACTGCAGGACTGTC ACGAGGGACGCATCGGACGACTGCAGGACTGT# 把所有的空格替換掉ct@ehbio:~$ sed ‘s/^(>.*)/1/’ test.fasta | tr ” ‘ ‘ | sed -e ‘s/ $//’ -e ‘s/ >/>/g’ -e ‘s/ //g’>SOX2 ? ?ACGAGGGACGCATCGGACGACTGCAGGACTGTCACGAGGGACGCATCGGACGACTGCAGGACTGTCACGAGGGACGCATCGGACGACTGCAGGAC>POU5F1 ? ?CGGAAGGTAGTCGTCAGTGCAGCGAGTCCGTCGGAAGGTAGTCGTCAGTGCAGCGAGTCC>NANOG ? ?ACGAGGGACGCATCGGACGACTGCAGGACTGTCACGAGGGACGCATCGGACGACTGCAGGACGAGGGACGCATCGGACGACTGCAGGACTGTCACGAGGGACGCATCGGACGACTGCAGGACTGT# 把TAB鍵轉(zhuǎn)換為換行符ct@ehbio:~$ sed ‘s/^(>.*)/1/’ test.fasta | tr ” ‘ ‘ | sed -e ‘s/ $//’ -e ‘s/ >/>/g’ -e ‘s/ //g’ -e ‘s///g’ >SOX2ACGAGGGACGCATCGGACGACTGCAGGACTGTCACGAGGGACGCATCGGACGACTGCAGGACTGTCACGAGGGACGCATCGGACGACTGCAGGAC>POU5F1CGGAAGGTAGTCGTCAGTGCAGCGAGTCCGTCGGAAGGTAGTCGTCAGTGCAGCGAGTCC>NANOGACGAGGGACGCATCGGACGACTGCAGGACTGTCACGAGGGACGCATCGGACGACTGCAGGACGAGGGACGCATCGGACGACTGCAGGACTGTCACGAGGGACGCATCGGACGACTGCAGGACTGT
或者簡單點(diǎn),直接用前面的awk略微做下修改。
# 差別只在一點(diǎn)# 對于單行fasta文件,只需要記錄一行,seq[name]=$0# 對于多好fasta文件,需要把每一行序列都加到前面的序列上,seq[name]=seq[name]$0ct@ehbio:~$ awk ‘BEGIN{OFS=FS=””}{if($0~/>/) {name=$0; sub(“>”, “”, name);} else seq[name]=seq[name]$0;}END{print “>SOX2”; print seq[“SOX2”]}’ test.fasta>SOX2ACGAGGGACGCATCGGACGACTGCAGGACTGTCACGAGGGACGCATCGGACGACTGCAGGACTGTCACGAGGGACGCATCGGACGACTGCAGGAC