第18題和第20題不會做咽块,求高手指教啊~~
該練習使用的數(shù)據(jù)為bowtie2中的示例數(shù)據(jù)龙助,獲取方法如下:
# 下載bowtie2軟件先
wget -c https://sourceforge.net/projects/bowtie-bio/files/bowtie2/2.3.5.1/bowtie2-2.3.5.1-linux-x86_64.zip
unzip bowtie2-2.3.5.1-linux-x86_64.zip
# 找到示例數(shù)據(jù)
cd bowtie2-2.3.5.1-linux-x86_64/example/reads
- 統(tǒng)計
reads_1.fq
文件中共有多少條序列信息
# fq文件格式有4行徽缚,均以@打頭实辑,統(tǒng)計@字符的數(shù)量即可 (此處有坑霜幼,字符@需要再行首)
less -S reads_1.fq | grep "^@" | wc
- 輸出所有的
reads_1.fq
文件中的標識符(即以@開頭的那一行)
# 篩選出來后重定向
less -S reads_1.fq | grep "^@" > head.txt
- 輸出
reads_1.fq
文件中的 所有序列信息(即每個序列的第二行)
# fq格式中的第二行為堿基序列瘾婿,觀察數(shù)據(jù)格式
less reads_1.fq | paste -d " " - - - - | less -SN | awk '{print$2}' > sequence.txt
- 輸出以‘+’及其后面的描述信息(即每個序列的第三行)
# 每個序列的第三行不是除了+都是空的咩
less reads_1.fq | paste -d " " - - - - | less -SN | awk '{print$3}' > line3.txt
- 輸出質(zhì)量值信息(即每個序列的第四行)
# 命令跟上面差不多
less reads_1.fq | paste -d " " - - - - | less -SN | awk '{print$4}' > line4.txt
- 計算
reads_1.fq
文件含有N
堿基的reads
個數(shù)
# 要用到上面生成的sequence.txt文件
less -SN sequence.txt | grep N | wc
- 統(tǒng)計文件中
reads_1.fq
文件里面的序列的堿基總數(shù)
# 還是要用到sequence.txt文件
cat sequence.txt | grep -o -i "[atcgn]" | wc
- 計算
reads_1.fq
所有的reads
中N
堿基的總數(shù)
# 很有必要再復習一下grep的各個參數(shù)的功能
cat sequence.txt | grep -o N | wc
- 統(tǒng)計
reads_1.fq
中測序堿基質(zhì)量值恰好為Q20
的個數(shù)
# Q20對應的ASCII值為5恬惯,質(zhì)量值在第4行
cat line4.txt | grep -o 5 | wc
- 統(tǒng)計
reads_1.fq
中測序堿基質(zhì)量值恰好為Q30
的個數(shù)
# Q30對應的ASCII值為?向拆,質(zhì)量值在第4行
cat line4.txt | grep -o ? | wc
- 統(tǒng)計
reads_1.fq
中所有序列的第一位堿基的ATCGNatcg
分布情況
# 篩選出來然后計數(shù)
less -S sequence.txt | cut -c-1 | sort | uniq -c
- 將
reads_1.fq
轉為reads_1.fa
文件(即將fastq
轉化為fasta
)
# 需明確.fq文件和.fa文件的異同±叶可知只需提取.fq文件前兩行浓恳,然后把"@"替換為">"即可,如下:
less reads_1.fq | paste - - - - | awk '{print $1,$2}' | tr " " "\n" | tr "@" ">" > reads_1.fa
- 統(tǒng)計上述
reads_1.fa
文件中共有多少條序列
# 相對.fq文件來說就更簡單了
cat reads_1.fa | grep "^>" | wc
- 計算
reads_1.fa
文件中總的堿基序列的GC
數(shù)量
cat reads_1.fa | grep -o GC | wc
- 刪除
reads_1.fa
文件中的每條序列的N
堿基
# 終于sed出場
cat reads_1.fa | sed 's/[nN]/""/g' > reads_1_non.fa
- 刪除
reads_1.fa
文件中的含有N
堿基的序列
# 先找到含有N的序列碗暗,然后刪除
cat reads_1.fa | grep -i -v n > reads_1_qual.fa
- 刪除
reads_1.fa
文件中的短于65bp
的序列
# 得復習awk語法了啊, 刪除短于65bp的颈将,那就選擇長于65bp的
cat reads_1.fa | paste - - | awk '{if(length($2)>65)print}' | tr "\t " "\n" > reads_1_65.fa
- 刪除
reads_1.fa
文件每條序列的前后五個堿基
# 看網(wǎng)上有人用下面的命令,實測是錯誤的
cat reads_1.fa | paste - - | cut -f 2 | cut -c 5- | cut -c -5 > reads_1_trim.fa
# 正確的是什么樣的言疗?
晴圾??噪奄???拨扶?
- 刪除
reads_1.fa
文件中的長于125bp
的序列
# 跟17題原理一樣纹因,篩選長度小于125bp的序列
cat reads_1.fa | paste - - | awk '{if(length($2)<125) print}' | tr "\t " "\n" > reads_1_125.fa
- 查看
reads_1.fq
中每條序列的第一位堿基的質(zhì)量值的平均值
# 思路應該是提取出條序列第一個堿基對應的質(zhì)量值逻悠,然后換算成可計算的數(shù)值,然后加和求平均數(shù)(吧)?可是怎么轉換成數(shù)字呢```
cat reads_1.fa | paste - - - - | awk '{print $4}' | cut -c -1 > reads_1.l4.fa
每天都有新的收獲的感覺很棒~