基因组测序模拟

杨利霞

基因组测序模拟
5 P+ [5 W* h! s/ o1 L# o' j# C基因组测序模拟

一、摘要
1 O9 P+ F" g; H! ]# T$ f/ C; h( {
( @$ W2 p3 g. X- F; l4 ?5 ~9 s通过熟悉已有的基因组测序模拟和评估程序，加深全基因组鸟枪法测序原理的理解，并且能够编写程序模拟全基因组鸟枪法测序，理解覆盖度、测序深度、拷贝数等概念，设置测序相关参数，生成单端/双端测序结果文件
- o1 O; [  t, _4 c/ f
二、材料和方法

' D1 S3 y$ C( c1、硬件平台

处理器:Intel（R） Core（TM）i7-4710MQ CPU @ 2.50GHz
; ~" I$ @/ x, d' r. O7 o) |安装内存(RAM)：16.0GB

2、系统平台
; F# N) o3 n( v4 U4 k( XWindows 8.1，Ubuntu
2 K% t0 z, v' N' G# x
/ P4 U0 }- L* U0 k- ~3、软件平台

6 z% e- e: k% }. [% A: uart_454
GenomeABC http://crdd.osdd.net/raghava/genomeabc/
Python3.5
& p. c) U- `0 V6 @3 I6 t4 nBiopython
4、数据库资源

: f! J: ]- S1 e. n' A, L5 V* qNCBI数据库：https://www.ncbi.nlm.nih.gov/
2 `. b- [* T# B: Z4 o
# T9 \6 u! _$ B4 _6 {5、研究对象
2 v# f; C6 B0 V% s
酵母基因组Saccharomyces cerevisiae S288c (assembly R64)
ftp://ftp.ncbi.nlm.nih.gov/genomes/all/GCF/000/146/045/GCF_000146045.2_R64/GCF_000146045.2_R64_genomic.fna.gz
: p/ P# Y4 u- H
6、方法
" G% }: }9 F0 L- x2 ?( b2 H" s% n
art_454的使用
首先至art系列软件的官网，下载软件，在ubuntu系统安装，然后阅读相关参数设置的帮助文档，运行程序。
GenomeABC
进入GenomeABC(http://crdd.osdd.net/raghava/genomeabc/)，输入参数，获得模拟测序结果。
: A3 ^/ p3 ~/ }( W/ c8 T/ g5 F! r$ T编程模拟测序
下载安装python，并且安装biopython扩展模块，编写程序，模拟单端/双端测序。
% T# O- m: ?$ b3 {7 Q' k三、结果

1、art_454的运行结果

无参数art_454运行，阅读帮助文档

图表 1无参数art_454运行
对酵母基因组进行基因组单端测序模拟，FOLD_COVERAGE设为20，即覆盖度为20.
下图为模拟单端测序，程序运行过程及结果

图表 2 art454单端测序
5 N3 h" D3 @8 s6 o
图表 3 art454单端模拟结果
双端测序模拟，FOLD_COVERAGE设为20，即覆盖度为20；MEAN_FRAG_LEN设为1500，即平均片段长度为1500；STD_DEV设为20，即长度的标准差为20
下图为模拟双端测序，程序运行过程及结果
2 @, E' B6 R- R# c$ L& D
# A/ D) y4 g. ^" Q图表 4 art454双端测序

图表 5 art454双端模拟结果
8 R8 J- j2 r7 j0 `; K4 u2、GenomeABC
; Q1 j7 ]6 S, D& o. d下图为设置参数页面

下图为结果下载页面
7 \5 O/ y9 Z; G/ W) I
图表 6 结果下载页面
3、编程模拟测序结果
拷贝数是这里的N值；覆盖度是m，测序深度是宏观的量，在这里与覆盖度意思相同，就是测序仪10X，20X。
单端测序
" J7 D$ @% x% d% T2 \
图表 7 程序模拟单端测序
% n' \. \9 _% V# U2 U双端测序

1 ]# B  I2 O! }  F$ ~图表 8 程序模拟双端测序
测序结果

图表 9 结果文件
" u: E) {/ N+ I  A# Z; {: \, ?0 ?  W/ O
因为期望片段长度是600bp，在片段长度区间200-1000bp内，所以大部分的片段都没有删除。
% U2 V# d; {) x$ V0 n, k1 i测序结果统计表
1 W/ R# p8 n2 e2 t& T# p+ m+ Z" E
$ S! i- a* g  r  c测序方式        基因组大小(bp)        片段长度区间 (bp)        N值        期望片段长度        克隆保留率        片段数量        Reads长度范围（bp）        Reads总数量        Reads总长度        覆盖度(m值)        理论丢失率(e-m)        覆盖率(1-e-m)
5 M5 ]# U( n* F: }; c6 }单端        12157kb        200-1000        10        600        0.95        107378        50-100        101968        7645.541kb        0.62889        0.53318        0.46682
+ q0 T$ t! }7 d; ]+ I; b& z单端        12157kb        200-1000        20        600        0.95        213722        50-100        202996        15227.882kb        1.25259        0.28576        0.71424
* j9 U9 o# I: A! E" Q9 y/ L双端        12157kb        200-1000        10        600        0.95        106704        50-100        202770        15212.662kb        1.25134        0.28612        0.71388
, z) k: T. m# l: V3 \$ F0 E% I双端        12157kb        200-1000        20        600        0.95        214212        50-100        407186        30534.265kb        2.51164        0.08114        0.91886
四、讨论和结论

程序运行方法

5 }/ _! T7 U0 O; d( I& O  |在类的构造方法init()中，调整参数。
& G& J) I% F: eAveragefragmentlength为片段平均的长度；
minfragmentlength和maxfragmentlength是保留片段的范围；
5 b1 b1 R/ k! b/ qcloneRetainprobability是克隆的保留率；
minreadslength和maxreadslength是测序reads的长度范围
* X8 P  A+ W' D) C  ~- @: Y
% c/ {, D$ o( Y模拟测序的诸多方法都封装成了Sequencing类，只需要创建类，并调用singlereadsequencing()和pairreadsequencing()方法，传入文件名的参数即可。

附录

, V2 y5 y- T2 N3 w1 a: }from Bio import SeqIO
9 X+ {5 ^$ c' L- w2 J8 |2 Jfrom math import exp
import random

class Sequencing:
    # N代表拷贝份数
9 W1 e- K, L* ~7 C    def __init__(self)
        self.fragmentList = []
        self.readsID = 1
        self.readsList = []
        self.averagefragmentlength = 650
        self.minfragmentlength = 500
8 Z# O* g6 ~, k        self.maxfragmentlength = 800
        self.cloneRetainprobability = 1
2 R) `' G0 y7 J2 I$ p  P        self.minreadslength = 50
7 i: e6 H1 `# X3 u  y        self.maxreadslength = 150
        self.N = 10
& S4 T% Z) b, a+ n% X        self.genomeLength = 0
4 F7 D4 J/ C7 {0 K' I; r" M5 C        self.allreadslength = 0
0 L/ p/ q5 V) a
& s) b! H1 P  _# P    # 生成断裂点
+ |. F1 \8 {! E( L    def generatebreakpoint(self, seqlen, averageLength):
* G- D2 Q' J" c' y4 ]2 b. H$ Y, H* C        # 假设平均每500bp 产生一个断裂点(averageLength = 500),通过随机函数生成seqlen/500个随机断裂点(1到seqlen之间的随机整数)
        breakpoint = [random.randint(0, seqlen) for _ in range(int(seqlen / averageLength))]
        breakpoint.append(seqlen)
        breakpoint.append(0)
        # 把随机断裂点从小到大排序
        breakpoint.sort()
* M) m7 Q9 s' A+ c/ D- u        return breakpoint
2 @7 a  S/ \- s& i, g. n7 m- A
# j7 }- A5 y2 o% R! O7 [# }+ i    # 沿断裂点打断基因组,并删除不符合长度要求的序列片段,定义片段范围:200-1000bp
    def breakgenome(self, seq, breakpoint, minfragmentlength, maxfragmentlength):
        for i in range(len(breakpoint) - 1):
            fragment = seq[breakpoint:breakpoint[i + 1]]
  g7 v% |' n6 L% L6 B            if maxfragmentlength > len(fragment) > minfragmentlength:
& @. a' x! J* x$ F1 k  w( {                self.fragmentList.append(fragment)
3 a( @" Z& Z8 O2 T' ]. @7 p9 L4 L! e        return self.fragmentList
3 A' F: Q( {$ b: x9 m: \% z
/ Y# T# [3 ]& `- @7 [: a    # 模拟克隆时的随机丢失情况,random.random()生成0-1的保留概率
    def clonefragment(self, fragmentList, cloneRetainprobability):
        clonedfragmentList = []
; Y! P  a$ E) K# U8 f; ~        Lossprobability = [random.random() for _ in range(len(fragmentList))]
7 K; |0 e$ T5 ]7 [# N  G# v$ h        for i in range(len(fragmentList)):
2 D9 o7 l" G8 S            if Lossprobability <= cloneRetainprobability:
                clonedfragmentList.append(fragmentList)
) O( o6 J* b, Y        return clonedfragmentList
6 m+ x& _! p0 S4 s  E: F' K& A
    # 模拟单端测序,并修改reads的ID号
    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
% N1 X: m* J- D+ Q            self.readsID += 1
            readslength = random.randint(self.minreadslength, self.maxreadslength)
) x! j$ g1 Y# b) y( b8 n            self.allreadslength += readslength
/ p3 E! w& W4 V# N# c( J            self.readsList.append(fragment[:readslength])
( r- X, n/ T" N/ S1 _  a' U
    def singlereadsequencing(self, genomedata, sequencingResult):
        for seq_record in SeqIO.parse(genomedata, "fasta"):
/ W: U9 s- {5 }( R& L& ]( {            seqlen = len(seq_record)
+ Z3 q/ U" K: H4 T            self.genomeLength += seqlen
            for i in range(self.N):
                # 生成断裂点
5 X7 h- L5 X: J8 c& m8 s% e! ^7 Y                breakpoint = self.generatebreakpoint(seqlen, self.averagefragmentlength)
                # 沿断裂点打断基因组
                self.breakgenome(seq_record, breakpoint, self.minfragmentlength, self.maxfragmentlength)
3 F# A1 _1 r6 y! }6 b  Y        # 模拟克隆时的随机丢失情况
' e. m1 D' U8 w        clonedfragmentList = self.clonefragment(self.fragmentList, self.cloneRetainprobability)
' C: ^. B) a5 W8 ?        # 模拟单端测序
% N0 C: j& l! u( t        self.singleread(clonedfragmentList)
        SeqIO.write(self.readsList, sequencingResult, "fasta")

    def pairread(self, clonedfragmentList):
        for fragment in clonedfragmentList:
  B0 |4 i, `# C' X7 v! F, k            fragment.id = ""
  v, V6 w5 b' i/ @            fragment.name = ""
5 ?$ @. r8 S1 F            description = fragment.description[12:].split(",")[0]
            fragment.description = str(self.readsID) + "." + description
( T" r- C3 C/ n8 j, c/ r# t& v8 ]- V            readslength = random.randint(self.minreadslength, self.maxreadslength)
            self.allreadslength += readslength
            self.readsList.append(fragment[:readslength])
+ V! S+ @' h4 O. r. A
! \( y1 p( S( q- N* r            readslength = random.randint(self.minreadslength, self.maxreadslength)
            self.allreadslength += readslength
7 ~  E( E( t' L) a6 i  {( {
/ p3 I6 I* h8 P+ A) \            fragmentcomplement = fragment.reverse_complement()
            fragmentcomplement.id = ""
# a8 g% k6 \  x8 m$ \' E: q            fragmentcomplement.name = ""
            fragmentcomplement.description = str(self.readsID) + "." + description
            self.readsList.append(fragmentcomplement[:readslength])

            self.readsID += 1

$ W: c# R3 v+ Y" f+ y    def pairreadsequencing(self,genomedata, sequencingResult_1, sequencingResult_2):
        for seq_record in SeqIO.parse(genomedata, "fasta"):
            seqlen = len(seq_record)
, p; X/ y/ W6 r# S* h4 t% k! r            self.genomeLength += seqlen
            for i in range(self.N):
5 p- ~5 w" c+ b$ a) P                # 生成断裂点
                breakpoint = self.generatebreakpoint(seqlen, self.averagefragmentlength)
, i* F3 i0 j% j5 k& M. z3 w                # 沿断裂点打断基因组
                self.breakgenome(seq_record, breakpoint, self.minfragmentlength, self.maxfragmentlength)
1 l' I& Z6 s* X+ c9 T1 I        # 模拟克隆时的随机丢失情况
        clonedfragmentList = self.clonefragment(self.fragmentList, self.cloneRetainprobability)
        # 模拟双端测序
        self.pairread(clonedfragmentList)
/ C8 Y" U6 L8 C) P1 v' R& ?# \        readsList_1 = [self.readsList for i in range(len(self.readsList)) if i % 2 == 0]
  s: ?( j, G$ l% W        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")
2 J/ z, f6 p2 F6 E4 X
    def resultsummary(self):
8 U" q- K0 @4 Y$ N4 T        print("基因组长度:" + str(self.genomeLength / 1000) + "kb")
        print("片段长度区间:" + str(self.minfragmentlength) + "-" + str(self.maxfragmentlength))
        print("N值:" + str(self.N))
1 k0 {3 o, W3 c% k; i        print("期望片段长度:" + str(self.averagefragmentlength))
6 K9 n4 F# O8 w5 d* f5 Z* t6 _        print("克隆保留率:" + str(self.cloneRetainprobability))
$ H' |5 e  l, k* r        print("片段数量:" + str(len(self.fragmentList)))
        print("reads长度:" + str(self.minreadslength) + "-" + str(self.maxreadslength))
        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)))
  L# E! N/ l- `8 b  d9 W2 m        print("覆盖率(1-e^-m):" + str(round(1 - exp(-m), 5)))
# -------------------------------------------主程序-------------------------------------------
# 模拟单端测序
3 L" R' z( D, l, Y2 s( n8 u9 TsequencingObj = Sequencing()
sequencingObj.singlereadsequencing("data/NC_025452.fasta", "result/virusSingleRead.fa")
sequencingObj.resultsummary()

# 模拟双端测序
4 s5 r: h+ V, F' a1 ]) RsequencingObj = Sequencing()
) B5 W* d% u6 e1 k, h) e# esequencingObj.pairreadsequencing("data/GCF_000146045.2_R64_genomic.fna", "result/yeastPairRead_1.fa", "result/yeastPairRead_2.fa")
4 G* ?2 {2 {" isequencingObj.resultsummary()
/ c* J/ K- i! z$ J6 ]# Pfrom Bio import SeqIO
from math import exp
. J: |: s- b0 S0 u2 r9 timport random

2 e6 z0 U$ u' F- e6 f$ @) \- zclass Sequencing:
    # N代表拷贝份数
- \$ n  J1 X" L+ I* }: I0 E    def __init__(self):
        self.fragmentList = []
        self.readsID = 1
( {. K% r( h0 M; S- Z        self.readsList = []
) _; q  B6 ~- g3 P9 Y        self.averagefragmentlength = 650
        self.minfragmentlength = 500
        self.maxfragmentlength = 800
+ C" f/ ]: ]+ O9 K$ c0 A        self.cloneRetainprobability = 1
$ z* B; @3 P- P! k2 r" H+ z/ A* I        self.minreadslength = 50
! [+ ~: R! t" i. U2 S        self.maxreadslength = 150
        self.N = 10
- s  f3 b8 E" R7 \- w+ p( q8 z$ Q1 H        self.genomeLength = 0
        self.allreadslength = 0
4 O0 x6 v9 V( C- s! {
    # 生成断裂点
    def generatebreakpoint(self, seqlen, averageLength):
/ h* u( K' A4 M        # 假设平均每500bp 产生一个断裂点(averageLength = 500),通过随机函数生成seqlen/500个随机断裂点(1到seqlen之间的随机整数)
        breakpoint = [random.randint(0, seqlen) for _ in range(int(seqlen / averageLength))]
        breakpoint.append(seqlen)
        breakpoint.append(0)
" f. Q" V" o7 ]; q& R, C        # 把随机断裂点从小到大排序
        breakpoint.sort()
8 [: e) v( b) W5 X( s: X        return breakpoint
) L! B2 q$ Y) X
9 S; f+ m! l' e& \0 F+ y    # 沿断裂点打断基因组,并删除不符合长度要求的序列片段,定义片段范围:200-1000bp
    def breakgenome(self, seq, breakpoint, minfragmentlength, maxfragmentlength):
9 ^* v* G3 r  a3 X        for i in range(len(breakpoint) - 1):
/ j. _* {' r6 Z& G4 S; ^" M4 {# Y            fragment = seq[breakpoint:breakpoint[i + 1]]
1 ?  O& [4 ]( t7 }            if maxfragmentlength > len(fragment) > minfragmentlength:
                self.fragmentList.append(fragment)
        return self.fragmentList
( C. d: m1 i0 L" y9 z
    # 模拟克隆时的随机丢失情况,random.random()生成0-1的保留概率
9 _0 l9 d% T3 P( u- ?    def clonefragment(self, fragmentList, cloneRetainprobability):
: I( {' I! C+ E2 b6 u: C        clonedfragmentList = []
        Lossprobability = [random.random() for _ in range(len(fragmentList))]
        for i in range(len(fragmentList)):
0 x# H$ S, t* `2 ~/ u& e            if Lossprobability <= cloneRetainprobability:
( V  m8 @- n" l1 f. N' r                clonedfragmentList.append(fragmentList)
8 Q) d2 ^2 h* \        return clonedfragmentList

    # 模拟单端测序,并修改reads的ID号
    def singleread(self, clonedfragmentList):
4 P; M5 H  d2 X  o# S        for fragment in clonedfragmentList:
            fragment.id = ""
            fragment.name = ""
9 R& y* w, M* q$ Y4 w! n            fragment.description = fragment.description[12:].split(",")[0]
            fragment.description = str(self.readsID) + "." + fragment.description
' Q: d; k8 I0 K1 ~7 C) L/ a+ }            self.readsID += 1
, t" H6 q! S# v& r7 F! B" Y            readslength = random.randint(self.minreadslength, self.maxreadslength)
            self.allreadslength += readslength
            self.readsList.append(fragment[:readslength])

: l4 H( F, A4 J  e    def singlereadsequencing(self, genomedata, sequencingResult):
        for seq_record in SeqIO.parse(genomedata, "fasta"):
1 \1 J$ D+ N  S; O            seqlen = len(seq_record)
7 ~! T( Y! z# z4 ]7 K            self.genomeLength += seqlen
! ~7 L1 U. W) C4 O0 Y/ v! f            for i in range(self.N):
8 E9 c1 B' W& J! f! y  N& `                # 生成断裂点
( {( _6 N; Q; }# O# t                breakpoint = self.generatebreakpoint(seqlen, self.averagefragmentlength)
+ Q! o  k3 G, k1 _6 Z                # 沿断裂点打断基因组
# z$ S3 S. G2 b                self.breakgenome(seq_record, breakpoint, self.minfragmentlength, self.maxfragmentlength)
        # 模拟克隆时的随机丢失情况
        clonedfragmentList = self.clonefragment(self.fragmentList, self.cloneRetainprobability)
        # 模拟单端测序
3 I. K# v; Y" \- Z7 X4 i8 \' h8 O        self.singleread(clonedfragmentList)
        SeqIO.write(self.readsList, sequencingResult, "fasta")
) x& }5 h" {/ V, O4 w7 J
    def pairread(self, clonedfragmentList):
3 q1 v; X3 [6 i( ]% B8 B        for fragment in clonedfragmentList:
            fragment.id = ""
5 `4 f% p  i* w3 m% l6 a1 [            fragment.name = ""
" E) }6 i: O0 M7 Q- r0 |5 l& v+ c            description = fragment.description[12:].split(",")[0]
* t+ [" i3 t) M* z' G3 ?2 E            fragment.description = str(self.readsID) + "." + description
/ }4 v) B2 t* v! ?! G5 O& j            readslength = random.randint(self.minreadslength, self.maxreadslength)
            self.allreadslength += readslength
# K5 q7 Y0 r. E. Y+ q9 x            self.readsList.append(fragment[:readslength])

            readslength = random.randint(self.minreadslength, self.maxreadslength)
6 z5 c3 g8 z, y; `+ N% V, S            self.allreadslength += readslength

9 `3 O, A* M5 M! G7 \9 P7 r8 L            fragmentcomplement = fragment.reverse_complement()
            fragmentcomplement.id = ""
+ M, O7 P9 P; ^' j, k            fragmentcomplement.name = ""
9 ^$ l) K# O" B9 k* a            fragmentcomplement.description = str(self.readsID) + "." + description
$ I+ q) y; ~5 r            self.readsList.append(fragmentcomplement[:readslength])
" t6 F7 @) m2 p* D. q" j; h4 k6 k
. }* r4 Y: c( H) x/ U6 P2 r5 P3 R            self.readsID += 1

, F& ], g. y. P+ Y2 n; ~    def pairreadsequencing(self,genomedata, sequencingResult_1, sequencingResult_2):
) q9 L: k2 n: P$ w        for seq_record in SeqIO.parse(genomedata, "fasta"):
1 L6 j6 {5 L. w& h: _9 @: L            seqlen = len(seq_record)
            self.genomeLength += seqlen
            for i in range(self.N):
                # 生成断裂点
- @* Q0 n6 I4 F) d! M, I; @                breakpoint = self.generatebreakpoint(seqlen, self.averagefragmentlength)
                # 沿断裂点打断基因组
                self.breakgenome(seq_record, breakpoint, self.minfragmentlength, self.maxfragmentlength)
        # 模拟克隆时的随机丢失情况
        clonedfragmentList = self.clonefragment(self.fragmentList, self.cloneRetainprobability)
+ h& ?2 D! x- \; M        # 模拟双端测序
        self.pairread(clonedfragmentList)
        readsList_1 = [self.readsList for i in range(len(self.readsList)) if i % 2 == 0]
( A8 @" C2 i! M' \; P; |5 P$ k0 D) i        readsList_2 = [self.readsList for i in range(len(self.readsList)) if i % 2 == 1]
        SeqIO.write(readsList_1, sequencingResult_1, "fasta")
# c0 }0 `1 u* M        SeqIO.write(readsList_2, sequencingResult_2, "fasta")
0 b/ n! F+ @  O8 C+ h1 n( X7 x- [
    def resultsummary(self):
        print("基因组长度:" + str(self.genomeLength / 1000) + "kb")
        print("片段长度区间:" + str(self.minfragmentlength) + "-" + str(self.maxfragmentlength))
6 F* P0 [) P; k: D$ ]5 _* p        print("N值:" + str(self.N))
        print("期望片段长度:" + str(self.averagefragmentlength))
' E  Q% D- X# ^  p( D! \        print("克隆保留率:" + str(self.cloneRetainprobability))
        print("片段数量:" + str(len(self.fragmentList)))
        print("reads长度:" + str(self.minreadslength) + "-" + str(self.maxreadslength))
! c( U1 r) Z' Z: O8 q        print("reads总数量:" + str(len(self.readsList)))
        print("reads总长度:" + str(self.allreadslength / 1000) + "kb")
+ h7 z* F, X0 q; n        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)))
# -------------------------------------------主程序-------------------------------------------
# 模拟单端测序
' e; h& n+ e) Y; {, x, SsequencingObj = Sequencing()
% _! `- n1 D# I, j) GsequencingObj.singlereadsequencing("data/NC_025452.fasta", "result/virusSingleRead.fa")
sequencingObj.resultsummary()

! A+ F/ T9 T8 s( a: W# 模拟双端测序
sequencingObj = Sequencing()
sequencingObj.pairreadsequencing("data/GCF_000146045.2_R64_genomic.fna", "result/yeastPairRead_1.fa", "result/yeastPairRead_2.fa")
+ J+ R$ t: @, ~- m% t2 w, gsequencingObj.resultsummary()
5 o1 ?* h% W5 T# N


7 f5 p$ c  J# Z8 x
8 N1 z2 ~: }, C8 V4 \8 Y& J3 m( M! U9 ~

3477959497

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