这个是有机结合生物信息学的linux和数据格式的练习题:
下载bowtie2软件后拿到示例数据:
mkdir -p ~/biosoft
cd ~/biosoft
wget -c https://sourceforge.net/projects/bowtie-bio/files/bowtie2/2.3.4.3/bowtie2-2.3.4.3-linux-x86_64.zip
unzip bowtie2-2.3.4.3-linux-x86_64.zip
cd ~/biosoft/bowtie2-2.3.4.3-linux-x86_64/example/reads
1)统计reads_1.fq 文件中共有多少条序列信息
grep -c '^@r' reads_1.fq
共10000条序列信息。
2)输出所有的reads_1.fq文件中的标识符(即以@开头的那一行)
grep '^@r' reads_1.fq
3)输出reads_1.fq文件中的 所有序列信息(即每个序列的第二行)
cat reads_1.fq | paste - - - - | cut -f 2 | less -SN
# awk '{if (NR%4 == 2) print}' reads_1.fq | less -SN
4)输出以‘+’及其后面的描述信息(即每个序列的第三行)
cat reads_1.fq | paste - - - - | cut -f 3 | less -SN
5)输出质量值信息(即每个序列的第四行)
cat reads_1.fq | paste - - - - | cut -f 4 | less -SN
6)计算reads_1.fq 文件含有N碱基的reads个数
cat reads_1.fq | paste - - - - | cut -f 2 | grep -n 'N' | wc -l
有6429条reads中包含了N碱基。
7)统计文件中reads_1.fq文件里面的序列的碱基总数
cat reads_1.fq | paste - - - - | cut -f 2 | wc -m
碱基总数为1098399。
8)计算reads_1.fq 所有的reads中N碱基的总数
cat reads_1.fq | paste - - - - | cut -f 2 | grep -o 'N' | wc -l
N碱基总数为26001。
9)统计reads_1.fq 中测序碱基质量值恰好为Q20的个数
cat reads_1.fq | paste - - - - | cut -f 4 | grep -o '5' | wc -l
# 20 + 33 ---> 53 ---> 5
Q20的个数为21369。
10)统计reads_1.fq 中测序碱基质量值恰好为Q30的个数
cat reads_1.fq | paste - - - - | cut -f 4 | grep -o '?' | wc -l
# 30 + 33 ---> 63 ---> ?
Q30的个数为21574。
11)统计reads_1.fq 中所有序列的第一位碱基的ATCGNatcg分布情况
cat reads_1.fq | paste - - - - | cut -f 2 | cut -c 1 | sort | uniq -c
12)将reads_1.fq 转为reads_1.fa文件(即将fastq转化为fasta)
cat reads_1.fq | paste - - - - | cut -f 1,2 | tr '\t' '\n' | tr '@' '>' > reads_1.fa
# awk '{if (NR%4 == 1) {print ">" substr($0, 2)}} {if (NR%4 == 2) print}' reads_1.fq > reads_1.fa
13)统计上述reads_1.fa文件中共有多少条序列
cat reads_1.fa | grep -v '>' | wc -l
# cat reads_1.fa | paste - - | cut -f 2 | wc -l
共10000条序列信息。
14)计算reads_1.fa文件中总的碱基序列的GC数量
grep -v '>' reads_1.fa | grep -o '[CG]' | wc -l
# cat reads_1.fa | paste - - | awk '{print $2}' | grep -o -e G -e C | wc -l
GC数量为529983。
15)删除 reads_1.fa文件中的每条序列的N碱基
tr -d 'N' < reads_1.fa
# cat reads_1.fa | tr -d "N" | grep 'N' -c
16)删除 reads_1.fa文件中的含有N碱基的序列
cat reads_1.fa | paste - - | cut -f 2 | sed '/N/d' | wc -l
17)删除 reads_1.fa文件中的短于65bp的序列
cat reads_1.fa | paste - - | cut -f 2 | awk '{if (length($1) >= 65) print}' | wc
18) 删除 reads_1.fa文件每条序列的前后五个碱基
cat reads_1.fa | paste - - | cut -f 2 | sed 's/^.....//' | sed 's/.....$//' | head
19)删除 reads_1.fa文件中的长于125bp的序列
cat reads_1.fa | paste - - | cut -f 2 | awk '{if (length($1) <= 125) print}' | wc
20)查看reads_1.fq 中每条序列的第一位碱基的质量值的平均值
cat reads_1.fq | paste - - - - | cut -f 4 | cut -c 1 | perl -alne '{print ord($_)-33}' | paste -s -d + | bc
# 涉及到perl单行命令,还需加强学习
参考:
http://www.bio-info-trainee.com/3575.html
http://ascii.911cha.com/
提示
echo $PATH |grep -o ":"|wc
cat reads_1.fq |paste - - - - |less -SN
cat reads_1.fq |paste - - - - |cut -f 2|grep -o [ATCGNatcg]|wc
# http://ascii.911cha.com/
# 63 ?
# 53 5
# 97 a
# 107 k
网友评论