原题出自生信技能树-fasta和fastq格式文件的shell小练习

该习题利用bowtie2软件自带的reads_1.fq文件

##1.统计reads_1.fq文件中一共有多少条序列信息
答：每个fq文件有4行，wc -reads_1.fq 除以4，就得到共有10000条序列

2.输出所有的reads_1.fq文件中的标识符（以@开头的那一行）
答：cat reads_1.fq | paste - - - -|cut -f1
读取reads_1.fq，然后每4行paste一下，-表示占位符.最后cut -f1就是第一行，-f2及时第二行，以此类推。

另外可以用awk '{if(NR%4==0)print}'  reads_1.fq      NR%4==0表示除以4=0，也就是第4行。NR%4==1，就是第二行。NR是awk里边的一个变量，就是行数。awk的变量可以记住

3.输出所有的reads_1.fq文件中的标识符（第2行）
4.输出所有的reads_1.fq文件中的标识符（第3行）
5.输出所有的reads_1.fq文件中的标识符（第四行）

6.计算reads_1.fq文件中含有N碱基的reads个数
答：awk '{if(NR%4==2)print}' reads_1.fq|grep N|wc 

先抽出fq文件的第二行，然后grep N ，计数。这里也可以用 grep -c N     就不用wc了

7.统计文件中reads_1.fq文件里边的序列的N碱基总数。
awk '{if(NR%4==2)print}' reads_1.fq|grep  -o N|wc
这里grep -o输出所有的N碱基，后边加上|wc就行了

8.统计文件中reads_1.fq文件里边的序列的碱基总数。
awk '{if(NR%4==2)print}' reads_1.fq|grep  -o [AGCTN] |wc
直接在fq文件第二行中寻找AGCTN（这就是所有碱基）。切记：加上[ ]
 
另一种方法：把fq文件每一条序列的第二行有多长计算出来，相加就可以了
awk '{if(NR%4==2)print length}' reads_1.fq |paste -s -d + |bc
length是awk里的一个参数，就是对就出每一行的长度，paste -s -d +   把这些lenght相加（会把算式paste出来，所有要加上bc输出最后结果）

9.统计reads_1.fq中测序碱基质量值恰好为Q20的个数
10.统计reads_1.fq中测序碱基质量值恰好为Q30的个数
答：确定fq文件第4行中寻找（第4行是质控信息），根据第4行的质控信息找到测序平台（这里是sanger -33），然后根据不同测序平台的质控值（有个表）找到对应的符号（Q20+33=53，53对应符号为5，Q30对应符号为？），最后在fq文件第4行进行寻找统计。
$awk '{if(NR%4==0)print}' reads_1.fq|grep -o "?"|wc      Q 20的把？换成5就行了。

11.统计reads_1.fq文件中所有序列的第一位碱基ATCGNatcg分布情况
答：先找到第一位碱基cut -c 1，找分布
$awk '{if(NR%4==2)print}' reads_1.fq|cut -c 1|sort |uniq -c

12.将reads_1.fq文件转为reads_1.fa文件（即fastq转fasta）
答：第一步:fa文件只需要fq文件的前两行就可以了，用第一题中的paste命令把fq文件前两行打印出来（此时只有一行），然后换行（把空格换成换行）
cat reads_1.fq | paste - - - -|cut -f1,2 |tr '\t' '\n'
第二步：fa文件以>开头，把@替换掉，加上tr '@ '>'

最终命令：cat reads_1.fq | paste - - - -|cut -f1,2|tr '\t' '\n'|tr '@' '>' >reads_1.fa

13.统计上述reads_1.fa文件中共有多少条序列
wc reads_1.fa

14.计算reads_1.fa文件中总的碱基序列的GC数量
答：按照fq相同的操作：awk '{if(NR%2==0)print}' reads_1.fa|grep  -o [GC] |wc
或者：cat reads_1.fa |paste -- |cut -f 2 |grep -o [GC]|wc
再或者直接grep -o [GC]  reads_1.fa|wc 因为fa文件中只有序列含有ATCG的

15.删除reads_1.fa文件中每条序列的N碱基
答：用tr命令-d就是删除  tr -d
cat reads_1.fa |tr -d 'N'

16.删除reads_1.fa文件中含N碱基的序列
答：grep抓取不含N 的序列 grep -v N    然后再转为fa格式
cat reads_1.fa|paste - - |grep -v N|tr '\t' '\n' 

17.删除reads_1.fa文件中短于65bp的序列
答：把大于65bp的揪出来就可以了，然后再转为fa格式
cat reads_1.fa|paste - - |awk '{if(length($2)>65)print}'|tr '\t' '\n' 

18.删除reads_1.fa文件中每条序列的前5个碱基
cat reads_1.fa|paste -- |cut -f2|cut -c 5-

后5个怎么删除？？
cat reads_1.fa|paste -- |cut -f2|cut -c -5

19.删除reads_1.fa文件中长于125bp的序列
同17题

20.查看reads_1.fq中每条序列的第一位碱基的质量的平均值
cat reads_1.fq |paste - - - - |cut -f4|cut -c1|sort|uniq -c
也可以用awk '{if(NR%4==0)print}' reads_1.fq |cut -c1|sort q -c

感谢生信技能树，jmmy大神！