基因组测序模拟

杨利霞

基因组测序模拟
, @( m) Z6 R, B9 y基因组测序模拟

一、摘要

通过熟悉已有的基因组测序模拟和评估程序，加深全基因组鸟枪法测序原理的理解，并且能够编写程序模拟全基因组鸟枪法测序，理解覆盖度、测序深度、拷贝数等概念，设置测序相关参数，生成单端/双端测序结果文件
+ ~3 m, O7 B% E/ V( m6 }# ^
二、材料和方法
# a4 U+ r7 S: Y
1、硬件平台

处理器:Intel（R） Core（TM）i7-4710MQ CPU @ 2.50GHz
$ s* i7 V( I+ M% t7 N4 X安装内存(RAM)：16.0GB
! @) q+ i, q* q8 B9 F& p2 a  j3 w
2、系统平台
5 m2 n# c& D1 }5 Z* m& \$ x, }Windows 8.1，Ubuntu

3、软件平台

art_454
3 o% V" v2 U$ PGenomeABC http://crdd.osdd.net/raghava/genomeabc/
Python3.5
Biopython
4、数据库资源

  h& K0 p1 [* [3 U5 YNCBI数据库：https://www.ncbi.nlm.nih.gov/

' `* W# M7 w1 w/ N  d, Z& |* N5、研究对象
2 L) I" X- w- X% @
$ V& J9 d9 u( @8 n; Z- k0 [酵母基因组Saccharomyces cerevisiae S288c (assembly R64)
( q( D; z7 L3 L& v  vftp://ftp.ncbi.nlm.nih.gov/genomes/all/GCF/000/146/045/GCF_000146045.2_R64/GCF_000146045.2_R64_genomic.fna.gz
3 G& g& d$ k0 M8 w
$ }, T; q6 h9 D% ]6、方法

8 ]3 _, |" [) W% h, Hart_454的使用
/ l, _7 O1 L- v' y0 G首先至art系列软件的官网，下载软件，在ubuntu系统安装，然后阅读相关参数设置的帮助文档，运行程序。
) X! _' i1 S; m8 fGenomeABC
进入GenomeABC(http://crdd.osdd.net/raghava/genomeabc/)，输入参数，获得模拟测序结果。
4 q; z( B+ W( w7 t6 X3 H编程模拟测序
下载安装python，并且安装biopython扩展模块，编写程序，模拟单端/双端测序。
# ^" L, m1 i4 w' o  v; g三、结果
: T6 A7 d3 P5 C+ I) M9 ^+ v
1、art_454的运行结果
, {& z; m+ {. Y$ q, n% C
1 B9 w+ f0 `! l无参数art_454运行，阅读帮助文档

, a% ^$ @( e7 }4 b; [图表 1无参数art_454运行
% K* ~1 V+ p' S8 _对酵母基因组进行基因组单端测序模拟，FOLD_COVERAGE设为20，即覆盖度为20.
下图为模拟单端测序，程序运行过程及结果
+ j$ r# M& F. O; G5 f: Z
! q+ ]) b" C  ~5 |! d5 {& E图表 2 art454单端测序
: ?& I* W" ~1 N% K* T, u  ~6 z, Q
$ h& U5 |" u; ?# m5 |% @图表 3 art454单端模拟结果
双端测序模拟，FOLD_COVERAGE设为20，即覆盖度为20；MEAN_FRAG_LEN设为1500，即平均片段长度为1500；STD_DEV设为20，即长度的标准差为20
4 d5 ]# Y9 T8 t) G1 S下图为模拟双端测序，程序运行过程及结果
, T$ H3 @/ U0 t4 q/ ~4 a5 D
  G. d* W* N: p: \  P; @8 s图表 4 art454双端测序
. q% ~! P" c. \3 T: i" \
" v8 Q3 q$ [# C8 C8 o$ I, W图表 5 art454双端模拟结果
2、GenomeABC
! ?# a( Z2 c/ ^5 W下图为设置参数页面

5 A' N- y5 B; X; n9 B) b2 M下图为结果下载页面

8 Z3 g& m! d: k$ Q' w2 d! b) Q图表 6 结果下载页面
3、编程模拟测序结果
8 ?7 B! @' f( X4 O' L* M拷贝数是这里的N值；覆盖度是m，测序深度是宏观的量，在这里与覆盖度意思相同，就是测序仪10X，20X。
) ]- Y9 r4 d0 x9 h- G' {7 j单端测序

图表 7 程序模拟单端测序
. X  B7 l; J: G& D1 Y# j5 ]双端测序

图表 8 程序模拟双端测序
) t: K7 T" f! w8 T测序结果

5 O; m2 k( E0 ^. o, X- \/ z& N图表 9 结果文件

因为期望片段长度是600bp，在片段长度区间200-1000bp内，所以大部分的片段都没有删除。
测序结果统计表
* J6 y6 B2 u6 f2 {  `
1 K' N1 K7 |1 s8 P9 ~) B( h测序方式        基因组大小(bp)        片段长度区间 (bp)        N值        期望片段长度        克隆保留率        片段数量        Reads长度范围（bp）        Reads总数量        Reads总长度        覆盖度(m值)        理论丢失率(e-m)        覆盖率(1-e-m)
单端        12157kb        200-1000        10        600        0.95        107378        50-100        101968        7645.541kb        0.62889        0.53318        0.46682
单端        12157kb        200-1000        20        600        0.95        213722        50-100        202996        15227.882kb        1.25259        0.28576        0.71424
  |( m, M& P% D双端        12157kb        200-1000        10        600        0.95        106704        50-100        202770        15212.662kb        1.25134        0.28612        0.71388
双端        12157kb        200-1000        20        600        0.95        214212        50-100        407186        30534.265kb        2.51164        0.08114        0.91886
3 C7 q2 q6 w4 `/ X  v四、讨论和结论

程序运行方法

; r9 A9 e8 R4 `& \8 d- I在类的构造方法init()中，调整参数。
Averagefragmentlength为片段平均的长度；
' L# O! e8 O4 y5 h/ c1 N, Pminfragmentlength和maxfragmentlength是保留片段的范围；
cloneRetainprobability是克隆的保留率；
  y& ^3 {/ T8 ?; X4 U- X, V- B' w1 z, dminreadslength和maxreadslength是测序reads的长度范围

模拟测序的诸多方法都封装成了Sequencing类，只需要创建类，并调用singlereadsequencing()和pairreadsequencing()方法，传入文件名的参数即可。
0 O& t5 R% S6 t/ s( T7 I% E
附录

from Bio import SeqIO
from math import exp
9 x' |$ r" ^) N7 }7 yimport random

% H' \* P1 K3 A3 h3 ?  v- P2 r2 iclass Sequencing:
    # N代表拷贝份数
+ b9 T3 S. a, @7 s    def __init__(self)
        self.fragmentList = []
        self.readsID = 1
        self.readsList = []
        self.averagefragmentlength = 650
        self.minfragmentlength = 500
6 E* H* T* F* g- M( C        self.maxfragmentlength = 800
        self.cloneRetainprobability = 1
        self.minreadslength = 50
/ ^2 m6 m7 k  l        self.maxreadslength = 150
        self.N = 10
  Z8 n3 `0 {6 X/ |        self.genomeLength = 0
        self.allreadslength = 0
! T" k: p; N! u# G4 o' G1 P
* r9 ?+ n/ V7 R8 s' F& }    # 生成断裂点
% i$ e' v; `1 W8 p# J2 z, ?6 {$ q1 x    def generatebreakpoint(self, seqlen, averageLength):
* Q3 J1 X' `! G% J% K6 r9 J- n% I2 p* T        # 假设平均每500bp 产生一个断裂点(averageLength = 500),通过随机函数生成seqlen/500个随机断裂点(1到seqlen之间的随机整数)
        breakpoint = [random.randint(0, seqlen) for _ in range(int(seqlen / averageLength))]
        breakpoint.append(seqlen)
* r) I" y! C4 T" `6 {        breakpoint.append(0)
* t; |5 u5 B7 A0 t, `( P' _5 Y( T        # 把随机断裂点从小到大排序
  T5 G" q7 G) D3 p  y# K0 R+ l        breakpoint.sort()
, k7 V/ ^+ e/ \1 q0 D( n% @        return breakpoint
0 X7 ^8 E4 Z: }  ~4 \0 @
; N! O. q: ~1 [) p8 a- p3 d    # 沿断裂点打断基因组,并删除不符合长度要求的序列片段,定义片段范围:200-1000bp
* Z. w& j2 `! v9 n/ Y6 T) e    def breakgenome(self, seq, breakpoint, minfragmentlength, maxfragmentlength):
' v- u; j! [. @  M: P1 ~, m, n# Y        for i in range(len(breakpoint) - 1):
. ?1 E' U' i( j" Z& a5 B            fragment = seq[breakpoint:breakpoint[i + 1]]
2 n* A) j$ X! d0 L, J0 X/ r            if maxfragmentlength > len(fragment) > minfragmentlength:
                self.fragmentList.append(fragment)
        return self.fragmentList

    # 模拟克隆时的随机丢失情况,random.random()生成0-1的保留概率
+ h* U. G& x7 s9 O, ~5 s& _0 S    def clonefragment(self, fragmentList, cloneRetainprobability):
        clonedfragmentList = []
( U% r- K  a; S; M% [# m9 J        Lossprobability = [random.random() for _ in range(len(fragmentList))]
& `; a# m! ]7 ?" M) w        for i in range(len(fragmentList)):
( Y9 P5 I" @+ L6 Q# `* B# s            if Lossprobability <= cloneRetainprobability:
+ H# q' o: f6 G7 ^- D                clonedfragmentList.append(fragmentList)
        return clonedfragmentList

    # 模拟单端测序,并修改reads的ID号
+ F, b/ k3 u+ _    def singleread(self, clonedfragmentList):
2 K0 e7 p* `& J        for fragment in clonedfragmentList:
4 A' u0 B$ q; E) F2 R' `' V            fragment.id = ""
  _' L9 W# ^4 X9 M            fragment.name = ""
            fragment.description = fragment.description[12:].split(",")[0]
2 i$ E( a! j8 c, h$ }0 M            fragment.description = str(self.readsID) + "." + fragment.description
0 \4 d- p9 }/ ~4 e: G) ?            self.readsID += 1
2 O; N$ H$ }$ I% U, D6 C            readslength = random.randint(self.minreadslength, self.maxreadslength)
            self.allreadslength += readslength
            self.readsList.append(fragment[:readslength])

    def singlereadsequencing(self, genomedata, sequencingResult):
. W/ C1 T% @+ O' u# f7 d. v        for seq_record in SeqIO.parse(genomedata, "fasta"):
            seqlen = len(seq_record)
1 c: |4 O( ?4 z% A( Z8 l1 T            self.genomeLength += seqlen
% m. D1 I6 ]; |. K            for i in range(self.N):
2 e: `( W- x" c0 O+ L                # 生成断裂点
                breakpoint = self.generatebreakpoint(seqlen, self.averagefragmentlength)
' x6 l7 c5 }2 H4 A1 c. W) b                # 沿断裂点打断基因组
9 N6 w& l$ ^) m                self.breakgenome(seq_record, breakpoint, self.minfragmentlength, self.maxfragmentlength)
+ w: a) q  [  N+ ^        # 模拟克隆时的随机丢失情况
" Q, S5 I/ @/ Y- [0 `8 N! J, G- f        clonedfragmentList = self.clonefragment(self.fragmentList, self.cloneRetainprobability)
% T4 j, g0 j  h/ v" p        # 模拟单端测序
' Z% x( @/ m6 v' x  i        self.singleread(clonedfragmentList)
        SeqIO.write(self.readsList, sequencingResult, "fasta")

    def pairread(self, clonedfragmentList):
        for fragment in clonedfragmentList:
            fragment.id = ""
            fragment.name = ""
* @$ {4 D& _2 @2 S            description = fragment.description[12:].split(",")[0]
            fragment.description = str(self.readsID) + "." + description
            readslength = random.randint(self.minreadslength, self.maxreadslength)
, E# _; `7 m' D5 ]' o9 X2 e            self.allreadslength += readslength
            self.readsList.append(fragment[:readslength])
; g- A( T3 t. J/ J: D
/ C; v6 c# M$ Q% c2 t! O            readslength = random.randint(self.minreadslength, self.maxreadslength)
            self.allreadslength += readslength

- K+ V. |  x8 ?) n- Q% n; @) Z; k+ d            fragmentcomplement = fragment.reverse_complement()
* R8 K; e) f! N9 F2 m8 P            fragmentcomplement.id = ""
% A8 ~7 x* z! _" f) y            fragmentcomplement.name = ""
6 t- ?/ M( \5 M/ g5 s9 |            fragmentcomplement.description = str(self.readsID) + "." + description
            self.readsList.append(fragmentcomplement[:readslength])
. p6 o% G3 {+ x: s# K! P
) \; `9 o. l8 `% x! g. e' {            self.readsID += 1

    def pairreadsequencing(self,genomedata, sequencingResult_1, sequencingResult_2):
        for seq_record in SeqIO.parse(genomedata, "fasta"):
: e% {8 Y2 P1 L  e7 E+ g            seqlen = len(seq_record)
            self.genomeLength += seqlen
            for i in range(self.N):
                # 生成断裂点
                breakpoint = self.generatebreakpoint(seqlen, self.averagefragmentlength)
                # 沿断裂点打断基因组
                self.breakgenome(seq_record, breakpoint, self.minfragmentlength, self.maxfragmentlength)
5 Z' t2 a) e2 _. ^' o! |2 e        # 模拟克隆时的随机丢失情况
2 S, f9 h7 X' E$ y+ L        clonedfragmentList = self.clonefragment(self.fragmentList, self.cloneRetainprobability)
7 a3 [+ @$ \2 y        # 模拟双端测序
" H( `' B7 d. o; [! k        self.pairread(clonedfragmentList)
% q6 g/ @1 z! v' R$ d        readsList_1 = [self.readsList for i in range(len(self.readsList)) if i % 2 == 0]
        readsList_2 = [self.readsList for i in range(len(self.readsList)) if i % 2 == 1]
        SeqIO.write(readsList_1, sequencingResult_1, "fasta")
        SeqIO.write(readsList_2, sequencingResult_2, "fasta")
& s3 \* N2 a" |) R
* B6 K0 R  C$ f" g$ m% R    def resultsummary(self):
& j: n  M2 I1 b- l        print("基因组长度:" + str(self.genomeLength / 1000) + "kb")
        print("片段长度区间:" + str(self.minfragmentlength) + "-" + str(self.maxfragmentlength))
        print("N值:" + str(self.N))
        print("期望片段长度:" + str(self.averagefragmentlength))
        print("克隆保留率:" + str(self.cloneRetainprobability))
" o: P/ a) r9 f9 T9 o3 d        print("片段数量:" + str(len(self.fragmentList)))
        print("reads长度:" + str(self.minreadslength) + "-" + str(self.maxreadslength))
0 f4 P! O# A5 O        print("reads总数量:" + str(len(self.readsList)))
        print("reads总长度:" + str(self.allreadslength / 1000) + "kb")
        m = self.allreadslength / self.genomeLength
        print("覆盖度(m值):" + str(round(m, 5)))
        print("理论丢失率(e^-m):" + str(round(exp(-m), 5)))
        print("覆盖率(1-e^-m):" + str(round(1 - exp(-m), 5)))
- }- J( U, A& k' O- l# -------------------------------------------主程序-------------------------------------------
+ r) ~- a( J1 S/ }( r: e7 M2 J# 模拟单端测序
! {' f( V  g- xsequencingObj = Sequencing()
sequencingObj.singlereadsequencing("data/NC_025452.fasta", "result/virusSingleRead.fa")
sequencingObj.resultsummary()
! S* Q" T& l1 s# F! U* g' |
( q) N, V0 N& S# 模拟双端测序
sequencingObj = Sequencing()
8 c1 f- h* S1 @# W- |sequencingObj.pairreadsequencing("data/GCF_000146045.2_R64_genomic.fna", "result/yeastPairRead_1.fa", "result/yeastPairRead_2.fa")
' w9 i2 u- e% o) H; rsequencingObj.resultsummary()
from Bio import SeqIO
from math import exp
* r3 T9 a2 X: C1 p9 Jimport random

# a6 r* F% \' K3 K$ jclass Sequencing:
    # N代表拷贝份数
% a7 B# X# M# y5 C    def __init__(self):
" `2 }, t3 ?+ t) S4 r        self.fragmentList = []
  f! I7 O' r8 I) d2 h# H        self.readsID = 1
        self.readsList = []
9 H5 o$ Y5 _; L: ]  ^2 c5 j        self.averagefragmentlength = 650
        self.minfragmentlength = 500
        self.maxfragmentlength = 800
        self.cloneRetainprobability = 1
        self.minreadslength = 50
        self.maxreadslength = 150
        self.N = 10
        self.genomeLength = 0
        self.allreadslength = 0

. P+ y& u9 N8 {* {7 H! X; T    # 生成断裂点
    def generatebreakpoint(self, seqlen, averageLength):
9 `8 q- B' Z/ t% `* U        # 假设平均每500bp 产生一个断裂点(averageLength = 500),通过随机函数生成seqlen/500个随机断裂点(1到seqlen之间的随机整数)
        breakpoint = [random.randint(0, seqlen) for _ in range(int(seqlen / averageLength))]
8 ^1 [: r) h- w1 ^3 N8 v4 |        breakpoint.append(seqlen)
( v2 ]) B$ B& L3 I- w9 v        breakpoint.append(0)
        # 把随机断裂点从小到大排序
        breakpoint.sort()
        return breakpoint
! A. t- T# z5 d& b6 z7 b0 o
# N! [! e4 ~9 w1 y" y    # 沿断裂点打断基因组,并删除不符合长度要求的序列片段,定义片段范围:200-1000bp
    def breakgenome(self, seq, breakpoint, minfragmentlength, maxfragmentlength):
$ z. c0 H% e# I1 T- J8 T" _        for i in range(len(breakpoint) - 1):
1 ?8 y* D+ s: X' G3 l/ p. |7 a8 w            fragment = seq[breakpoint:breakpoint[i + 1]]
            if maxfragmentlength > len(fragment) > minfragmentlength:
                self.fragmentList.append(fragment)
        return self.fragmentList

    # 模拟克隆时的随机丢失情况,random.random()生成0-1的保留概率
    def clonefragment(self, fragmentList, cloneRetainprobability):
        clonedfragmentList = []
        Lossprobability = [random.random() for _ in range(len(fragmentList))]
        for i in range(len(fragmentList)):
) q( Z% g- \5 N: x: }            if Lossprobability <= cloneRetainprobability:
+ B$ _" v: @  H# p$ x                clonedfragmentList.append(fragmentList)
        return clonedfragmentList
( {/ z; V5 t5 D
+ p. I+ X" c. O( q6 A' Z    # 模拟单端测序,并修改reads的ID号
$ E# W2 i0 n8 S9 ~    def singleread(self, clonedfragmentList):
        for fragment in clonedfragmentList:
            fragment.id = ""
            fragment.name = ""
            fragment.description = fragment.description[12:].split(",")[0]
            fragment.description = str(self.readsID) + "." + fragment.description
5 `: \2 E$ K: R! W            self.readsID += 1
            readslength = random.randint(self.minreadslength, self.maxreadslength)
- P0 L4 q! h* L8 o2 @            self.allreadslength += readslength
            self.readsList.append(fragment[:readslength])

    def singlereadsequencing(self, genomedata, sequencingResult):
1 s* U* z7 C8 C# @6 C9 w- X        for seq_record in SeqIO.parse(genomedata, "fasta"):
            seqlen = len(seq_record)
            self.genomeLength += seqlen
            for i in range(self.N):
                # 生成断裂点
                breakpoint = self.generatebreakpoint(seqlen, self.averagefragmentlength)
' E. w& V0 t+ ?                # 沿断裂点打断基因组
                self.breakgenome(seq_record, breakpoint, self.minfragmentlength, self.maxfragmentlength)
        # 模拟克隆时的随机丢失情况
! Q4 K7 a( Z) [5 W        clonedfragmentList = self.clonefragment(self.fragmentList, self.cloneRetainprobability)
% ?3 `# X$ X0 u$ C% s0 {' p        # 模拟单端测序
        self.singleread(clonedfragmentList)
9 s) `" ?+ v0 K' i5 }  x4 G, q' S3 M        SeqIO.write(self.readsList, sequencingResult, "fasta")
2 ^; b! v4 \* ?0 _9 t, A& q  C# i
/ \- s5 e0 h0 U  g) s- Y    def pairread(self, clonedfragmentList):
        for fragment in clonedfragmentList:
2 q9 C, ^; v, N. ~( ?            fragment.id = ""
            fragment.name = ""
            description = fragment.description[12:].split(",")[0]
            fragment.description = str(self.readsID) + "." + description
            readslength = random.randint(self.minreadslength, self.maxreadslength)
            self.allreadslength += readslength
            self.readsList.append(fragment[:readslength])
; R9 c3 u: a% T6 N+ Y, P5 b
            readslength = random.randint(self.minreadslength, self.maxreadslength)
            self.allreadslength += readslength
0 W8 z' x& h, m$ R, Q5 i
            fragmentcomplement = fragment.reverse_complement()
            fragmentcomplement.id = ""
            fragmentcomplement.name = ""
- D3 G7 t% E" B            fragmentcomplement.description = str(self.readsID) + "." + description
            self.readsList.append(fragmentcomplement[:readslength])

7 b4 h. ]% e6 d$ h8 f: j% i            self.readsID += 1
, c' a4 m; E7 P4 q
    def pairreadsequencing(self,genomedata, sequencingResult_1, sequencingResult_2):
        for seq_record in SeqIO.parse(genomedata, "fasta"):
3 n) T% f( c& _            seqlen = len(seq_record)
& m! ~1 ~1 ?( m            self.genomeLength += seqlen
            for i in range(self.N):
8 {# h9 i/ w( ?                # 生成断裂点
                breakpoint = self.generatebreakpoint(seqlen, self.averagefragmentlength)
                # 沿断裂点打断基因组
                self.breakgenome(seq_record, breakpoint, self.minfragmentlength, self.maxfragmentlength)
- u: U  v5 p+ N/ i1 w% J        # 模拟克隆时的随机丢失情况
        clonedfragmentList = self.clonefragment(self.fragmentList, self.cloneRetainprobability)
        # 模拟双端测序
; Y+ ~) c0 A! M+ c1 q& Z& @        self.pairread(clonedfragmentList)
4 j8 Q/ G6 D/ E" ^        readsList_1 = [self.readsList for i in range(len(self.readsList)) if i % 2 == 0]
        readsList_2 = [self.readsList for i in range(len(self.readsList)) if i % 2 == 1]
8 _) S* h( ]- R' _        SeqIO.write(readsList_1, sequencingResult_1, "fasta")
  l) Y& n% {7 s9 I        SeqIO.write(readsList_2, sequencingResult_2, "fasta")
1 j4 [$ H/ ^" \0 S8 {  n
    def resultsummary(self):
        print("基因组长度:" + str(self.genomeLength / 1000) + "kb")
9 S6 R7 H$ O2 p        print("片段长度区间:" + str(self.minfragmentlength) + "-" + str(self.maxfragmentlength))
& N4 }2 x; l) |& C; Z) f        print("N值:" + str(self.N))
- G, W' ]9 q7 ~0 ~& m+ ?: M        print("期望片段长度:" + str(self.averagefragmentlength))
        print("克隆保留率:" + str(self.cloneRetainprobability))
  ~& k+ v" k  q' d        print("片段数量:" + str(len(self.fragmentList)))
7 c( _8 f1 |* R0 f9 ~3 K% b        print("reads长度:" + str(self.minreadslength) + "-" + str(self.maxreadslength))
" j0 R4 K9 {: i" I6 Z( Y        print("reads总数量:" + str(len(self.readsList)))
7 S1 E% Y6 v( y        print("reads总长度:" + str(self.allreadslength / 1000) + "kb")
        m = self.allreadslength / self.genomeLength
        print("覆盖度(m值):" + str(round(m, 5)))
        print("理论丢失率(e^-m):" + str(round(exp(-m), 5)))
3 b& j+ a- a& `, i8 z, ~" R$ V( C        print("覆盖率(1-e^-m):" + str(round(1 - exp(-m), 5)))
8 S; H9 n! }) v! O" F! \# -------------------------------------------主程序-------------------------------------------
# 模拟单端测序
sequencingObj = Sequencing()
sequencingObj.singlereadsequencing("data/NC_025452.fasta", "result/virusSingleRead.fa")
sequencingObj.resultsummary()

' M! S3 e( W; A9 N- ?$ q1 h: |0 F# 模拟双端测序
. ?" G" U1 _4 ^# P7 E, t$ MsequencingObj = Sequencing()
sequencingObj.pairreadsequencing("data/GCF_000146045.2_R64_genomic.fna", "result/yeastPairRead_1.fa", "result/yeastPairRead_2.fa")
sequencingObj.resultsummary()

- b/ p9 |/ a1 E" x0 a  R* _

* _, }. r6 k" O8 j, b/ [# h4 I  O
: R& H. z% E; E, ]0 d: a. ]) J

3477959497

不错。。。。。。。。。。。。。。。。。