劉小澤寫于2020.5.21
記錄一下踩的坑
上報錯
手上有一些ATAC-seq的雙端測序數(shù)據(jù)（25bp）進(jìn)行了fastp質(zhì)控過濾后準(zhǔn)備進(jìn)行比對，結(jié)果比對過程發(fā)生報錯，信息如下
# Bowtie2
Error, fewer reads in file specified with -1 than in file specified with -2
terminate called after throwing an instance of 'int'
(ERR): bowtie2-align died with signal 6 (ABRT) (core dumped)
貌似顯示雙端測序數(shù)據(jù)中兩個文件reads數(shù)不等
如何查看測序數(shù)據(jù)的信息
簡單使用腳本
zcat raw/test_1.fastq.gz | wc -l| awk '{print $1/4}'
# 4237351
zcat raw/2d_2.fastq.gz | wc -l| awk '{print $1/4}'
# 4237351
使用小工具
seqkit stat raw/test_1.fastq.gz raw/test_2.fastq.gz
# file format type num_seqs sum_len min_len avg_len max_len
# raw/test_1.fastq.gz FASTQ DNA 4,237,351 105,933,775 25 25 25
# raw/test_2.fastq.gz FASTQ DNA 4,237,351 105,933,775 25 25 25
看到原始數(shù)據(jù)的兩個文件中reads數(shù)相等
既然原始數(shù)據(jù)一致，那么質(zhì)控過濾后的呢？
果然，利用seqkit stat發(fā)現(xiàn)clean的一對PE數(shù)據(jù)中reads數(shù)不一致
看了一下如何進(jìn)行過濾的
sample=test
fq1=~/atac/raw/${sample}_1.fastq.gz
fastp -i $fq1 -o clean/${sample}.fq1.gz \
-I $fq2 -O clean/${sample}.fq2.gz \
--thread=4 --length_required=10 --n_base_limit=5 \
--compression=6 -R ${sample}_1 1>log/log.fastp.${sample}.1.txt 2>&1
fq2=~/atac/raw/${sample}_2.fastq.gz
fastp -i $fq2 -o clean/${sample}.fq2.gz \
--thread=4 --length_required=10 --n_base_limit=5 \
--compression=6 -R ${sample}_2 1>log/log.fastp.${sample}.2.txt 2>&1
忽然想起來，這個腳本是從ChIP-seq腳本中復(fù)制粘貼過來的，只改了對應(yīng)的fq1和fq2，但一般ChIP-seq都是單端的，所以這個fastp就是對單端數(shù)據(jù)的過濾
問題的根源在于：這里把雙端數(shù)據(jù)當(dāng)成了兩個單端分別進(jìn)行處理，那么結(jié)果就不一致啦
將代碼修改為雙端對應(yīng)的就好了
# 如果是雙端，就需要指定-i -o (for input 1)
# -I -O (for input 2)
sample=test
fq1=~/atac/raw/${sample}_1.fastq.gz
fq2=~/atac/raw/${sample}_2.fastq.gz
fastp -i $fq1 -o clean/${sample}.fq1.gz \
-I $fq2 -O clean/${sample}.fq2.gz \
--thread=4 --length_required=10 --n_base_limit=5 \
--compression=6 -R ${sample} 1>log/log.fastp.${sample}.txt 2>&1
收獲
一般來說，軟件的報錯信息都會顯示在前面，需要觀察加發(fā)散思維
例如
# bowtie2
Error: reads file does not look like a FASTQ file
terminate called after throwing an instance of 'int'
(ERR): bowtie2-align died with signal 6 (ABRT) (core dumped)
看到這里，就要去看原始的fq文件（主要信息就是第2行和第4行），可能存在長度不一致的情況
zcat fastq.gz | paste - - - - | awk -F"\t" '{ if (length($2) != length($4)) print $0 }'
# 如果要挑出來這些reads
zcat fastq.gz | paste - - - - | awk -F"\t" '{ if (length($2) != length($4)) print $0 }' | tr '\t' '\n' > error_reads.fastq
再例如，一個R語言的報錯
安裝rMATS時看著滿屏的報錯，眼花繚亂，但后面其實都在重復(fù)說找不到libgsl.so.0這個庫
于是關(guān)鍵詞搜索一下libgsl.so.0 rmats：https://www.biostars.org/p/336864/
或者使用：sudo apt-get install libgsl0-dev （https://cloud.tencent.com/developer/article/1366294）
歡迎關(guān)注我們的公眾號～_～
我們是兩個農(nóng)轉(zhuǎn)生信的小碩，打造生信星球，想讓它成為一個不拽術(shù)語、通俗易懂的生信知識平臺。需要幫助或提出意見請后臺留言或發(fā)送郵件到jieandze1314@gmail.com
