轉(zhuǎn)錄組PCR去重

原始數(shù)據(jù)查看:

查看行數(shù)

grep '^@' ctr.R1.fq | sort |uniq|wc -l

27787625 也就是乘以四后的行數(shù)遇八,即fq的name不會重

比對:

index=/media/pc/disk1/sun/refdata/ensembl_GRCm38/03_bowtie2_index/GRCm38

bowtie2 -p 20 -x $index -1 ctr.R1.fq -2 ctr.R2.fq -S ctr.sam &&

bowtie2 -p 20 -x $index -1 ko.R1.fq -2 ko.R2.fq -S ko.sam

比對結(jié)果:

查看sam行數(shù):

igv:

samtools sort --threads 10 -m 2G -o?ctr.bam?ctr.sam &&

samtools sort --threads 10 -m 2G -o??ko.bam? ko.sam

samtools index?ctr.bam &&

samtools index?ko.bam

提取mapped:

samtools view -h -@ 10 -F 4 ctr.sam > ctr.mapped.sam &&

samtools view -h -@ 10 -F 4 ko.sam > ko.mapped.sam

wc -l 查看mapped的sam行數(shù):

對mapped的進行igv:

samtools sort --threads 10 -m 3G -o?ctr.mapped.bam?ctr.mapped.sam &&

samtools sort --threads 10 -m 3G -o??ko.mapped.bam? ko.mapped.sam &&

samtools index?ctr.mappedbam &&

samtools index?ko.mapped.bam

#################################

trim:

trimmomatic PE -phred33 -threads 20 ctr.R1.fq ctr.R2.fq ctr.R1.paired.fq ctr.R1.unpaired.fq ctr.R2.paired.fq ctr.R2.unpaired.fq ILLUMINACLIP:/home/pc/miniconda3/pkgs/trimmomatic-0.38-1/share/trimmomatic-0.38-1/adapters/TruSeq3-PE.fa:2:30:10:8:true LEADING:3 TRAILING:3 SLIDINGWINDOW:4:15 AVGQUAL:20 MINLEN:36

trimmomatic PE -phred33 -threads 20?ko.R1.fq?ko.R2.fq?ko.R1.paired.fq?ko.R1.unpaired.fq?ko.R2.paired.fq ko.R2.unpaired.fq ILLUMINACLIP:/home/pc/miniconda3/pkgs/trimmomatic-0.38-1/share/trimmomatic-0.38-1/adapters/TruSeq3-PE.fa:2:30:10:8:true LEADING:3 TRAILING:3 SLIDINGWINDOW:4:15 AVGQUAL:20 MINLEN:36?

比對:

bowtie2 -p 20 -x $index -1 ctr.R1.paired.fq -2 ctr.R2.paired.fq -S ctr.paired.sam &&

bowtie2 -p 20 -x $index -1 ko.R1.paired.fq -2 ko.R2.paired.fq -S ko.paired.sam

經(jīng)過trim后比對情況并未改善:



1.fatsquniq:

網(wǎng)站https://sourceforge.net/projects/fastuniq/

用法:

list文件

去重命令:

fastuniq -i list -o fastuniq.ctr.R1.fq -p fastuniq.ctr.R2.fq

fastuniq -i list2 -o fastuniq.ko.R1.fq -p fastuniq.ko.R2.fq

去重后的行數(shù)

grep '^@' fastuniq.ctr.R1.fq | sort | uniq | wc -l

22295210 乘以四后也是行數(shù);

去重率:

ctrl:80.23%

ko:86.58%

比對:

bowtie2 -p 20 -x $index -1 fastuniq.ctr.R1.fq -2 fastuniq.ctr.R2.fq -S fastuniq.ctr.sam?&&

bowtie2 -p 20 -x $index -1?fastuniq.ko.R1.fq?-2?fastuniq.ko.R2.fq?-S?fastuniq.ko.sam

比對結(jié)果:

查看sam:

44590443 fastuniq.ctr.sam

44155981 fastuniq.ko.sam

igv:

samtools sort --threads 10 -m 2G -o fastuniq.ctr.bam fastuniq.ctr.sam &&?

samtools sort --threads 10 -m 2G -o fastuniq.ko.bam fastuniq.ko.sam &&

samtools index?fastuniq.ctr.bam &&

samtools index??fastuniq.ko.bam

提取mapped:

samtools view -h -@ 10 -F 4?fastuniq.ctr.sam >?fastuniq.ctr.mapped.sam &&

samtools view -h -@ 10 -F 4?fastuniq.ko.sam >?fastuniq.ko.mapped.sam

查看mapped行數(shù):

提取mapped后不能用igv秩铆。原因是缺少header:

解決:提取mapped時加-h參數(shù):ref:[E::sam_parse1] missing SAM header



2.seqkit:

網(wǎng)站https://bioinf.shenwei.me/seqkit/

用法:

命令:

seqkit rmdup -s ctr.R1.fq -o seqkit.ctr.R1.fq.gz -d seqkit.ctr.R1.duplicated.fq.gz -D seqkit.ctr.R1.duplicated.txt -j 5

seqkit rmdup -s ctr.R2.fq -o seqkit.ctr.R2.fq.gz -d seqkit.ctr.R2.duplicated.fq.gz -D seqkit.ctr.R2.duplicated.txt -j 5

seqkit rmdup -s ko.R1.fq -o seqkit.ko.R1.fq.gz -d seqkit.ko.R1.duplicated.fq.gz -D seqkit.ko.R1.duplicated.txt -j 5

seqkit rmdup -s ko.R2.fq -o seqkit.ko.R2.fq.gz -d seqkit.ko.R2.duplicated.fq.gz -D seqkit.ko.R2.duplicated.txt -j 5

查看行數(shù):

zcat seqkit.ctr.R1.fq.gz | wc -l -> 57594872 51.82% [去重]

zcat seqkit.ctr.R2.fq.gz | wc -l -> 63321348 56.97%??[去重]

zcat seqkit.ko.R1.fq.gz | wc -l -> 60717652 54.63%??[去重]

zcat seqkit.ko.R2.fq.gz | wc -l -> 67585076 66.26%??[去重]

比對:

bowtie2 -p 20 -x $index -1?seqkit.ctr.R1.fq.gz?-2?seqkit.ctr.R2.fq.gz?-S?seqkit.ctr.sam?&&

bowtie2 -p 20 -x $index -1?seqkit.ko.R1.fq?-2?seqkit.ko.R2.fq?-S?seqkit.ko.sam

報錯了:雙端不對齊

結(jié)論:

seqkit去重率約50-70%榛做,但是使用后雙端不齊浅乔,無法比對赌厅。但是可以保留去重掉的fastq焕刮。



3.super_deduper:

官網(wǎng)https://github.com/dstreett/Super-Deduper

新地址https://github.com/ibest/HTStream

用法

命令:

super_deduper -1 ctr.R1.fq -2?ctr.R1.fq -g -p super_deduper.ctrl

super_deduper -1 ko.R1.fq -2?ko.R1.fq?-g -p?super_deduper.ko

查看行數(shù):

zcat super_deduper.ctrl_R1.fastq.gz |wc -l -> 3915060 3.52% [去重]

zcat super_deduper.ctrl_R2.fastq.gz |wc -l -> 3915060 3.52% [去重]

zcat super_deduper.ko_R1.fastq.gz |wc -l -> 3971360 3.89%?[去重]

zcat super_deduper.ko_R2.fastq.gz |wc -l -> 3971360 3.89%?[去重]

修改參數(shù):

super_deduper -1 ctr.R1.fq -2 ctr.R1.fq -g -p super_deduper.ctrl -q 30

super_deduper -1 ko.R1.fq -2 ko.R1.fq -g -p super_deduper.ko -q 20

查看again:

zcat super_deduper.ctrl_R1.fastq.gz | wc -l -> 3915060?3.52% [去重]

zcat super_deduper.ctrl_R2.fastq.gz | wc -l -> 3915060?

zcat super_deduper.ko_R1.fastq.gz | wc -l?-> 3971360?3.89%?[去重]

zcat super_deduper.ko_R2.fastq.gz | wc -l?-> 3971360? 改了參數(shù)后竟然沒變

比對:

bowtie2 -p 20 -x $index -1?super_deduper.ctrl_R1.fastq.gz?-2?super_deduper.ctrl_R2.fastq.gz?-S?super_deduper.ctrl.sam?&&

bowtie2 -p 20 -x $index -1?super_deduper.ko_R1.fastq.gz?-2?super_deduper.ko_R2.fastq.gz?-S?super_deduper.ko.sam

比對結(jié)果:

查看sam行數(shù):

igv:

samtools sort --threads 10 -m 2G -o super_deduper.ctrl.bam super_deduper.ctrl.sam && samtools sort --threads 10 -m 2G -o super_deduper.ko.bam super_deduper.ko.sam

samtools index?super_deduper.ctrl.bam &&

samtools index??super_deduper.ko.bam

提取mapped:

samtools view -h -@ 10 -F 4?super_deduper.ctrl.sam >?super_deduper.ctrl.mapped.sam &&

samtools view -h -@ 10 -F 4?super_deduper.ko.sam >?super_deduper.ko.mapped.sam

查看mapped的行數(shù)

結(jié)論:

super_deduper率約3-4%注益,去重去的太多碴巾。原理是什么呢?



4.fastp:

官網(wǎng)https://github.com/OpenGene/fastp

用法:

命令:

fastp -i ctr.R1.fq -I ctr.R2.fq -o fastp.ctr.R1.fq -O fastp.ctr.R2.fq

fastp -i ko.R1.fq -I ko.R2.fq -o fastp.ko.R1.fq -O fastp.ko.R2.fq

查看:

ctrl
ko

查看:

fastp.ctr.R1.fq ->?27082408*4 97.46%?[去重]

fastp.ctr.R2.fq?->?27082408*4?97.46%?[去重]

fastp.ko.R1.fq?->?24873179*4 97.54% [去重]

fastp.ko.R2.fq?->?24873179*4?97.54%?[去重]

比對:

bowtie2 -p 20 -x $index -1?fastp.ctr.R1.fq?-2?fastp.ctr.R2.fq?-S?fastp.ctr.sam?&&

bowtie2 -p 20 -x $index -1?fastp.ko.R1.fq?-2?fastp.ko.R2.fq?-S?fastp.ko.sam

比對情況:

提取mapped:

samtools view -h -@ 10 -F 4?fastp.ctr.sam >?fastp.ctr.mapped.sam &&

samtools view -h -@ 10 -F 4?fastp.ko.sam >?fastp.ko.mapped.sam

查看行:

35,942,791 fastp.ctr.mapped.sam

32,536,623 fastp.ko.mapped.sam

igv:

samtools sort --threads 10 -m 2G -o?fastp.ctr.mapped.bam?fastp.ctr.mapped.sam&& samtools sort --threads 10 -m 2G -o?fastp.ko.mapped.bamfastp.ko.mapped.sam &&

samtools index?fastp.ctr.mapped.bam??&&

samtools index?fastp.ko.mapped.bam



5.samtools:

官網(wǎng):http://www.htslib.org/

samtools sort --threads 10 -n? ctr.bam?-o ctr.sort.bam &&

samtools fixmate?--threads 10 -m ctr.sort.bam ctr.fixmate.bam &&

samtools sort?--threads 10 ctr.fixmate.bam?-o ctr.positionsort.bam &&

samtools markdup --threads 10 -r ctr.positionsort.bam ctr.markdup.bam

samtools sort --threads 10 -n??ko.bam?-o?ko.sort.bam &&

samtools fixmate?--threads 10?-m?ko.sort.bam?ko.fixmate.bam &&

samtools sort?--threads 10?ko.fixmate.bam?-o?ko.positionsort.bam &&

samtools markdup?--threads 10?-r?ko.positionsort.bam ko.markdup.bam?

查看行:

samtools view??ctr.markdup.bam | wc -l?39,334,011

samtools?view ko.markdup.bam | wc -l?39,170,695

提取mapped:

samtools view -h -@ 10 -F 4??ctr.markdup.bam > ctr.markdup.mapped.sam &&

samtools view -h -@ 10 -F 4??ko.markdup.bam >?ko.markdup.mapped.sam

查看mapped的行數(shù):

wc -l?ctr.markdup.mapped.sam? ->?20,069,023

wc -l?ko.markdup.mapped.sam ->?21,054,225

igv:

samtools sort --threads 10 -m 2G -o ctr.markdup.mapped.bam ctr.markdup.mapped.sam && samtools sort --threads 10 -m 2G -o ko.markdup.mapped.bam ko.markdup.mapped.sam &&

samtools index?ctr.markdup.mapped.bam &&

samtools index??ko.markdup.mapped.bam


6.picard:

官網(wǎng):https://broadinstitute.github.io/picard/

用法

直接刪除冗余:

java -jar /home/pc/biosoft/picard.jar MarkDuplicates REMOVE_DUPLICATES=true I=ctr.sort.bam O=picard.ctr.sorted.markdup.bam M=picard.ctr.markdup.txt

****************************************************

java -jar /home/pc/biosoft/picard.jar MarkDuplicates REMOVE_DUPLICATES=true I=ko.sort.bam O=picard.ko.sorted.markdup.bam M=picard.ko.markdup.txt

查看行:

samtools?view picard.ctr.sorted.markdup.bam | wc -l ->?39,302,611

samtools?view picard.ko.sorted.markdup.bam | wc -l ->?39,151,069

igv:

samtools index picard.ctr.sorted.markdup.bam

samtools index picard.ko.sorted.markdup.bam

提取mapped:

samtools view -h -@ 10 -F 4??picard.ctr.sorted.markdup.bam >?picard.ctr.sorted.markdup.mapped.sam &&?

samtools view -h -@ 10 -F 4??picard.ko.sorted.markdup.bam >?picard.ko.sorted.markdup.mapped.sam

查看行數(shù):

20,037,624 picard.ctr.sorted.markdup.mapped.sam

21,034,600 picard.ko.sorted.markdup.mapped.sam

igv:

samtools sort --threads 10 -m 2G -o?picard.ctr.sorted.markdup.mapped.bam?picard.ctr.sorted.markdup.mapped.sam?&&

samtools sort --threads 10 -m 2G -o?picard.ko.sorted.markdup.mapped.bam?picard.ko.sorted.markdup.mapped.sam &&

samtools index?picard.ctr.sorted.markdup.mapped.bam?&&

samtools index?picard.ko.sorted.markdup.mapped.bam

# 比對效率低丑搔?

用hisat2比對看竟然比對效率提高這么多

最后編輯于
?著作權(quán)歸作者所有,轉(zhuǎn)載或內(nèi)容合作請聯(lián)系作者
  • 序言:七十年代末厦瓢,一起剝皮案震驚了整個濱河市提揍,隨后出現(xiàn)的幾起案子,更是在濱河造成了極大的恐慌旷痕,老刑警劉巖碳锈,帶你破解...
    沈念sama閱讀 222,183評論 6 516
  • 序言:濱河連續(xù)發(fā)生了三起死亡事件,死亡現(xiàn)場離奇詭異欺抗,居然都是意外死亡售碳,警方通過查閱死者的電腦和手機,發(fā)現(xiàn)死者居然都...
    沈念sama閱讀 94,850評論 3 399
  • 文/潘曉璐 我一進店門绞呈,熙熙樓的掌柜王于貴愁眉苦臉地迎上來贸人,“玉大人,你說我怎么就攤上這事佃声∫罩牵” “怎么了?”我有些...
    開封第一講書人閱讀 168,766評論 0 361
  • 文/不壞的土叔 我叫張陵圾亏,是天一觀的道長十拣。 經(jīng)常有香客問我,道長志鹃,這世上最難降的妖魔是什么夭问? 我笑而不...
    開封第一講書人閱讀 59,854評論 1 299
  • 正文 為了忘掉前任,我火速辦了婚禮曹铃,結(jié)果婚禮上缰趋,老公的妹妹穿的比我還像新娘。我一直安慰自己陕见,他們只是感情好秘血,可當(dāng)我...
    茶點故事閱讀 68,871評論 6 398
  • 文/花漫 我一把揭開白布。 她就那樣靜靜地躺著评甜,像睡著了一般灰粮。 火紅的嫁衣襯著肌膚如雪。 梳的紋絲不亂的頭發(fā)上忍坷,一...
    開封第一講書人閱讀 52,457評論 1 311
  • 那天谋竖,我揣著相機與錄音,去河邊找鬼承匣。 笑死蓖乘,一個胖子當(dāng)著我的面吹牛,可吹牛的內(nèi)容都是我干的韧骗。 我是一名探鬼主播嘉抒,決...
    沈念sama閱讀 40,999評論 3 422
  • 文/蒼蘭香墨 我猛地睜開眼,長吁一口氣:“原來是場噩夢啊……” “哼袍暴!你這毒婦竟也來了些侍?” 一聲冷哼從身側(cè)響起隶症,我...
    開封第一講書人閱讀 39,914評論 0 277
  • 序言:老撾萬榮一對情侶失蹤,失蹤者是張志新(化名)和其女友劉穎岗宣,沒想到半個月后蚂会,有當(dāng)?shù)厝嗽跇淞掷锇l(fā)現(xiàn)了一具尸體,經(jīng)...
    沈念sama閱讀 46,465評論 1 319
  • 正文 獨居荒郊野嶺守林人離奇死亡耗式,尸身上長有42處帶血的膿包…… 初始之章·張勛 以下內(nèi)容為張勛視角 年9月15日...
    茶點故事閱讀 38,543評論 3 342
  • 正文 我和宋清朗相戀三年胁住,在試婚紗的時候發(fā)現(xiàn)自己被綠了。 大學(xué)時的朋友給我發(fā)了我未婚夫和他白月光在一起吃飯的照片刊咳。...
    茶點故事閱讀 40,675評論 1 353
  • 序言:一個原本活蹦亂跳的男人離奇死亡彪见,死狀恐怖,靈堂內(nèi)的尸體忽然破棺而出娱挨,到底是詐尸還是另有隱情余指,我是刑警寧澤,帶...
    沈念sama閱讀 36,354評論 5 351
  • 正文 年R本政府宣布跷坝,位于F島的核電站酵镜,受9級特大地震影響,放射性物質(zhì)發(fā)生泄漏柴钻。R本人自食惡果不足惜淮韭,卻給世界環(huán)境...
    茶點故事閱讀 42,029評論 3 335
  • 文/蒙蒙 一、第九天 我趴在偏房一處隱蔽的房頂上張望顿颅。 院中可真熱鬧缸濒,春花似錦足丢、人聲如沸粱腻。這莊子的主人今日做“春日...
    開封第一講書人閱讀 32,514評論 0 25
  • 文/蒼蘭香墨 我抬頭看了看天上的太陽绍些。三九已至,卻和暖如春耀鸦,著一層夾襖步出監(jiān)牢的瞬間柬批,已是汗流浹背。 一陣腳步聲響...
    開封第一講書人閱讀 33,616評論 1 274
  • 我被黑心中介騙來泰國打工袖订, 沒想到剛下飛機就差點兒被人妖公主榨干…… 1. 我叫王不留氮帐,地道東北人。 一個月前我還...
    沈念sama閱讀 49,091評論 3 378
  • 正文 我出身青樓洛姑,卻偏偏與公主長得像上沐,于是被迫代替她去往敵國和親。 傳聞我的和親對象是個殘疾皇子楞艾,可洞房花燭夜當(dāng)晚...
    茶點故事閱讀 45,685評論 2 360

推薦閱讀更多精彩內(nèi)容