寫(xiě)在前面：我是想著給課題組師弟師妹們培訓(xùn)出的一些例題。哈哈哈颅湘，沒(méi)想到也難到了自己戳稽，因此花上了一些時(shí)間思考著這些平時(shí)我們用python或者其他現(xiàn)有的軟件實(shí)現(xiàn)的功能灸芳，如何用三劍客解決察署？<br />只能說(shuō)三劍客真的比我想象中還有強(qiáng)大吧闷游，我對(duì)他們了解還不夠，日后加深了解贴汪。
Shell命令與生物信息相關(guān)的示例用法
- CDS文件內(nèi)提取對(duì)應(yīng)ID的命令- grep&sed
- 檢測(cè)CDS可編碼力 - grep & awk
- fastq轉(zhuǎn)fasta -awk
- 序列反向互補(bǔ)- awk&sed
CDS文件內(nèi)提取對(duì)應(yīng)ID的命令
ID：Acyan01G0000500.1 <br />file : final.gene.cds.fa
# 了解文件格式
less -S final.gene.cds.fa
#格式屬于一行ID一行序列
grep Acyan01G0000500.1 final.gene.cds.fa -A1
- 如果要搜索多個(gè)ID： Acyan01G0000500.1 Acyan01G0000800.1 Acyan03G0000900.1 的序列脐往？
vi idfile
grep -f idfile final.gene.cds.fa -A1
grep -f idfile final.gene.cds.fa -A1 |grep -v "^-"
了解 -f, -v的含義 "^" "$" 的含義
#grep參數(shù)學(xué)習(xí)：
-A n 指定輸出搜索目標(biāo)的后n行
-B n 指定輸出搜索目標(biāo)的前n行
-f <file> 以文件內(nèi)容為搜索目標(biāo)
-v 排除搜索到的行
"^" 行頭符號(hào)
"$" 行尾符號(hào)
思考題：對(duì)于不是一行ID+一行序列的如何提取扳埂？
sed -n '/LITCHI014205.m1/, /^>/p' Lchinesis_genome.Chr.cds| sed '$d'
事實(shí)上钙勃，這個(gè)方法我也想不到，全靠百度找答案聂喇。<br />思路：我們需要打印ID匹配到的一行到下一個(gè)基因的一行。參數(shù) -n 只輸出符合要求的行蔚携，其他行都不打酉Ｌ； 匹配具有ID和起始為">"字符之間的行內(nèi)容都實(shí)現(xiàn)打印酝蜒，最后我們的內(nèi)容是會(huì)包括了下一個(gè)基因ID行誊辉，所以刪除最后一行就可以實(shí)現(xiàn)了。
檢測(cè)CDS序列的編碼性
判斷ATG 起始位點(diǎn)
grep -v ">" final.gene.cds.fa | grep -v "^ATG" |head
grep -v ">" final.gene.cds.fa | grep -v "^ATG" > tmp
grep -f tmp final.gene.cds.fa -B1 |grep ">"
如何更加優(yōu)雅亡脑？
grep -f tmp final.gene.cds.fa -B1 |grep ">" > tmp && grep -f tmp final.gene.cds.fa -B1 |grep ">" |head
grep -f tmp final.gene.cds.fa -B1 |grep ">" |grep -f - final.gene.cds.fa -B1 |grep ">" > fail.id
如果遇到非一行名字一行序列格式的文件堕澄？
File : Lchinesis_genome.Chr.cds
grep ">" Lchinesis_genome.Chr.cds -A1|head
grep ">" Lchinesis_genome.Chr.cds -A1|grep -v ">"|grep -v "^-" | grep -v '^ATG'
grep ">" Lchinesis_genome.Chr.cds -A1|grep -v ">"|grep -v "^-" | grep -v '^ATG' |grep -f - Lchinesis_genome.Chr.cds -B1 |grep ">"
判斷是否滿(mǎn)足三的倍數(shù)邀跃？
awk -F "" '{if($1~/>/){print}else{print NF}}' final.gene.cds.fa
awk -F "" '{if($1~/>/){print}else{print NF}}' final.gene.cds.fa |cut -f1 -d" " |head
#整理可看格式
awk -F "" '{if($1~/>/){print}else{print NF%3}}' tmp|cut -f1 -d" " |sed ':a;N; s/\n/\t/'|sed 's/>//' | sort -k2,2nr
- 參數(shù)解析：
-F 指定分隔符，上述例子制定了分割符號(hào)為無(wú)
NF 是一個(gè)變量指代【域的數(shù)量】蛙紫，指代以分隔符為分割的域的總數(shù)拍屑；很復(fù)雜，如果分隔符默認(rèn)為"\t"坑傅，那么NF就是列數(shù)僵驰，但是由于我們分割符已經(jīng)設(shè)定了為無(wú)，所以是指代一行的總數(shù)唁毒；
$1 第一域蒜茴， ~ 一般是用于正則表達(dá)式上是否符合的符號(hào) // 用于正則式的包裝
所以第一個(gè)命令指代的是，如果行頭是以">"為符號(hào)的話浆西，我們就直接打臃鬯健；如果不是近零，我們就統(tǒng)計(jì)行的總字符數(shù)（核酸長(zhǎng)度）诺核；
NF%3 就是核酸長(zhǎng)度是否能整除3，如果有余數(shù)就是非3的倍數(shù)秒赤；
最后再用sed來(lái)去替換換行符為tab分隔符猪瞬；
關(guān)于:a;N; 這個(gè)陌生的符號(hào)，我下面有一個(gè)圖入篮，我感覺(jué)我也說(shuō)的不太清楚陈瘦。
最后的最后排序查看是否有余數(shù)不為0的ID
fastq-fasta轉(zhuǎn)換
#行首匹配@字符
awk '{if($1~/^@/){print}}' test1_R1.fq|head
#捕獲測(cè)序序列
awk '{if($1~/^[ATCG]/){print}}' test1_R1.fq|head
#捕獲測(cè)序序列為帶@和ATCG序列
awk '{if($1~/^@/){print}if($1~/^[ATCG]/){print}}' test1_R1.fq |head
#整理格式
awk '{if($1~/^@/){print}if($1~/^[ATCG]/){print}}' test1_R1.fq |sed 's/@/>/'|head
#重定向到新文件
awk '{if($1~/^@/){print}if($1~/^[ATCG]/){print}}' test1_R1.fq |sed 's/@/>/' > test_R1.fastq
相對(duì)上面的一些例子，fastq-fasta轉(zhuǎn)換就變得更加簡(jiǎn)單了潮售。
序列反向互補(bǔ)配對(duì)
實(shí)現(xiàn)反向
- file=test.seq
awk -F "" '{if($1~/>/){print}else{{for(i=NF;i>0;i--)printf $i}printf"\n"}}' test.seq|head -n2
awk -F "" '{if($1~/>/){print}else{{for(i=NF;i>0;i--)printf $i}printf"\n"}}' test_rev.seq
同樣設(shè)置分割符為無(wú)痊项，說(shuō)明逐個(gè)字符查看；
用到一個(gè)awk內(nèi)置判斷酥诽，if(){} else{}
如果以">"起始鞍泉，我們直接打印，否則我們就實(shí)現(xiàn)反向
實(shí)現(xiàn)反向用到了一個(gè)內(nèi)置的for()循環(huán)肮帐，i=NF 定義變量i為行數(shù)咖驮；i>0 必須大于0；i-- 實(shí)現(xiàn)每循環(huán)減i=i-1
然后逐一輸出域值（堿基）;
注意這里用到的是printf不是print训枢。其兩個(gè)區(qū)別在于托修，printf只輸出字符，不自帶換行符恒界，而print會(huì)在字符后加換行符睦刃。
當(dāng)完成所有else里的內(nèi)容再額外補(bǔ)充一次換行符，進(jìn)入下一個(gè)基因十酣；
實(shí)現(xiàn)互補(bǔ)
sed '/^[ATCG]/y/ATCG/TAGC/' test_rev.seq > test_rev_comp.seq
- 解析： 針對(duì)開(kāi)頭異ATCG的行進(jìn)行逐一替換涩拙；
sed 匹配條件/^[ATCG]/ 滿(mǎn)足以A/T/C/G開(kāi)頭的行际长，就可以排除了ID行；
然后再進(jìn)行y動(dòng)作的替換兴泥，這個(gè)替換在很多教程都沒(méi)怎么提及工育；也是一個(gè)替換動(dòng)作，但跟s有點(diǎn)不一樣郁轻；
它是一個(gè)逐一替換的動(dòng)作翅娶，A-T,T-A,C-G,G-C 同時(shí)實(shí)現(xiàn)替換，不會(huì)重復(fù)替換好唯；
進(jìn)一步實(shí)現(xiàn)了互補(bǔ)動(dòng)作
因?yàn)榧磳⒋蛩阒v解的內(nèi)容是三劍客竭沫，所以以上都以三劍客進(jìn)行整理命令；事實(shí)上同等替換的骑篙，容易理解的命令還是很多蜕提； 反向可以用rev， 替換可以用tr
ps： 真實(shí)數(shù)據(jù)是很多已經(jīng)造好的車(chē)子靶端，如：seqkit谎势。也不需要我們折騰著造輪子，但是知識(shí)往往都是在造輪子的過(guò)程中積累杨名。
