下载FATSQ，读取10条序列并计算每条序列的长度和GC含量

FASTQ

这是目前存储测序数据最普遍、最公认的一个数据格式，另一个是uBam格式，但这篇文章中不打算对其进行介绍。上面所讲的FASTA文件，它所存的都是已经排列好的序列（如参考序列），FASTQ存的则是产生自测序仪的原始测序数据，它由测序的图像数据转换过来，也是文本文件，文件大小依照不同的测序量（或测序深度）而有很大差异，小的可能只有几M，大的则常常有几十G上百G，文件后缀通常都是.fastq，.fq或者.fq.gz（gz压缩）

独特的格式：每四行成为一个独立的单元，称之为read。具体的格式描述如下：

第一行：以‘@’开头，是这一条read的名字，这个字符串是根据测序时的状态信息转换过来的，中间不会有空格，它是每一条read的唯一标识符，同一份FASTQ文件中不会重复出现，甚至不同的FASTQ文件里也不会有重复；
第二行：测序read的序列，由A，C，G，T和N这五种字母构成，这也是我们真正关心的DNA序列，N代表的是测序时那些无法被识别出来的碱基；
第三行：以‘+’开头，在旧版的FASTQ文件中会直接重复第一行的信息，但现在一般什么也不加（节省存储空间）；
第四行：测序read的质量值，这个和第二行的碱基信息一样重要，它描述的是每个测序碱基的可靠程度，用ASCII码表示。

1、NCBI上的SRA数据库下载sra格式的数据

wget -c https://sra-pub-run-odp.s3.amazonaws.com/sra/SRR2031888/SRR2031888   
#下载SRR2031888.sra
astq-dump SRR2031888  #sra文件转fastq文件

通过filezilla将linux系统中的fastq文件移动保存到电脑桌面上（也可以直接在linux系统中下载python来完成）

2、文件读取10条序列

打开python软件

file=open('/Users/andao/Desktop/SRR2031888.fastq',mode='r+')
#在python中以用于读写的模式打开fastq文件
line_num = 0
f = open("/Users/andao/Desktop/10read.txt", mode="w")
for line in file:
  if line_num < 40:
    print (line.strip(),file=f)
    line_num += 1
  else:
    break
file.close() # 关闭打开的文件
f.close() #此时前10条序列已经写入10read.txt文件中

print (line.strip())打印在屏幕上会发现是前10条序列

3、计算每条序列的长度和GC含量

s = open("/Users/andao/Desktop/length.GC.txt", mode="w")
import re
with open('/Users/andao/Desktop/10read.txt') as fq:
    fq2fa_dict = {}
    i = 0
    percent_GC = 0
    for line in fq:
        i += 1
        if i % 4 == 1:
            line = line.replace('\n', '')
            seq_name = '>' + line[1:]
            fq2fa_dict[seq_name] = ''
        elif i % 4 == 2:
            seq = line.replace('\n', '')
            fq2fa_dict[seq_name] = seq
            size = len(seq)
            nG = seq.count("g") + seq.count("G")
            nC = seq.count("c") + seq.count("C")
            percent_GC = (nG + nC) / size
            fq2fa_dict[seq_name] = size, percent_GC

print(fq2fa_dict,file=s) #打印内容输出到文件length.GC.txt中
s.close()