全长转录组测序技术解析不同转移性能肝癌细胞系的转录本表达谱与结构变异

杨淑欣，彭民武，罗娟娟，徐景祥，万绍贵

（赣南医学院基础医学院分子病理中心，江西 赣州 341000）

肝癌是一类发病率高、死亡率高的恶性肿瘤，也是全球范围内最常见的癌症之一，而肝细胞癌是原发性肝癌中最主要的亚型［1］。中国是全球肝癌发病率最高的国家［2］，同时肝癌在我国常见肿瘤中发病率位居第二［3］。尽管肝细胞癌的治疗取得了巨大进展，但由于术后复发和转移率高，5年生存率仅为15%~30%［4］。常见的肝内和肝外转移是导致肝细胞癌患者临床预后不良的主要原因［5］。因此，揭示肝癌转移的分子机制对改善其临床治疗具有重要意义。

可变剪切是一种关键的转录后基因表达调控方式，有助于扩大蛋白质种类复杂性和调控mRNA代谢［6-7］，研究报道mRNA的可变剪切在多种肿瘤中能够改变蛋白质组的多样性，异常的可变剪切参与肿瘤细胞增殖和转移的调控过程，在肿瘤发生、发展、药物治疗及耐药等方面发挥重要作用［8-10］。ARHGEF2基因的可变剪切转录本V1和V3变体与肝癌的转移和预后密切相关，提示可变剪切调控可能是肝癌转移的重要分子机制［11］。

高通量测序技术的发展促进了关于肿瘤基因变异的相关研究，该技术为研究者在基因表达差异、可变剪切和新转录本鉴定等研究提供极大的便利［12-13］。尽管高通量测序以前所未有的准确度和深度提供了巨大的测序能力，然而，传统的二代测序技术由于短读长的特点，其测序数据需要通过生物信息学算法拼接来分析可变剪接，因此无法完全准确地反映全长剪接转录本的序列和定量表达情况［14-15］，这一缺点将限制其在肿瘤中对可变剪接的大规模研究和转化应用。最新的三代测序技术克服了二代测序技术的不足，比如纳米孔测序技术，具有长读长的优势，平均测序长度能够达到10 kb［16］。长读长测序技术能够更加全面准确地识别新转录本、可变剪切和基因融合等结构变异，并进行转录本定量分析［17-18］。

本研究利用纳米孔三代测序技术对两株不同转移潜能的肝癌细胞系MHCC97H和MHCC97L进行全长转录组测序分析，旨在从转录本水平上系统地揭示肝癌转移可能存在的分子机制。

1 材料与方法

1.1 细胞系本实验使用的高转移潜能肝癌细胞系（MHCC97H）及低转移潜能肝癌细胞系（MHCC97 L）由赣南医学院血管生物学研究平台提供。

1.2 全长转录组文库构建及上机测序利用RNA simple总RNA提取试剂盒（天根生化科技北京有限公司）分别提取MHCC97H及MHCC97L细胞系总RNA，每株细胞做3个重复，共6个样本。然后在Qsep-100 Advance毛细管电泳仪上使用RNA卡夹对RNA样本进行完整度质控检测。利用Maxima H Minus Reverse Transcriptase试剂盒（Thermo Fisher公司，美国）进行反转录，并使用cDNA-PCR测序试剂盒（SQK-PCB109，Oxford Nanopore Technologies公司）制备cDNA测序文库，具体步骤按照试剂盒说明书进行，简要过程如下：逆转录得到的cDNA加Switch Oligo，再合成互补链，然后经DNA末端修复加A尾，AMPure XP磁珠纯化并添加测序接头。将上述cDNA测序文库上样到Nanopore测序芯片（FLO-MIN106D，R9.4.1）中，在Nanopore MinION平台上通过MinKNOW2.2软件操作进行全长转录组测序。

1.3 数据过滤与生信分析Nanopore MinION测序下的原始数据（Raw data）格式为包含所有原始测序信号的FAST5格式，通过MinKNOW 2.2软件包中的Guppy软件对原始数据进行电流信号到碱基序列信息的转换，并将数据转换为FASTQ格式。委托北京百迈客生物科技有限公司对上述FASTQ格式数据进行过滤与生信分析，进一步过滤短片段和低质量的原始数据，以及去除接头序列得到有效数据（Clean data），过滤条件为所测序列长度＞500 bp，测序质量Q score＞7。基于Oxford Nanopore Technologies单分子实时测序技术的全长转录组测序无须打断RNA片段，反转录得到全长cDNA。该平台的超长读取包含了单条完整转录本序列信息，后期分析无需组装，所测即所得。根据所测转录本与参考基因组的进行比对分析，将比对到参考基因组gtf文件未注释区域的序列定义为新基因。真核生物的基因往往具有多个转录本，将由可变剪切等导致的结构不同的转录本定义为新转录本。

1.3.1 转录本表达量分析转录组测序可以模拟成一个随机抽样的过程，为了让片段数目能真实地反映转录本表达水平，需要对样品中Mapped Reads的数量进行归一化。采用CPM（counts per million）［19］作为衡量转录本或基因表达水平的指标，CPM计算公式如下：CPM=reads mapped to transcript/total reads aligned in sample×1 000 000（“reads mapped to transcript”表示比对到某一转录本上的reads数，“total reads aligned in sample”表示比对到参考转录组的片段总数）。使用edgeR R package（3.8.6）进行两株细胞间的差异表达分析。使用Benjamini和Hochberg软件计算PValue和FDR（False discovery rate），将Pvalue＜0.01和差异倍数≥1.5的基因认定为差异表达基因。

1.3.2 差异表达转录本KEGG注释对差异表达基因的信号通路注释分析有助于进一步解读基因的功能。KEGG（Kyoto Encyclopedia of Genes and Genomes，http：//www.genome.jp/kegg/）是系统分析基因功能、基因组信息的数据库，它有助于研究者把基因及表达信息作为一个整体网络进行研究。本研究利用KOBAS［20］软件进行差异表达转录本的KEGG信号通路注释。

1.3.3 可变剪切分析通过Astalavista软件［21］获取每个样品存在的可变剪切类型，主要的基因可变剪切类型主要包括外显子跳跃（Exon skipping，ES）、内含子保留（Intron retained，IR）、可变5ʹ剪切位点（Alternative 5ʹsplice-site，A5S）、可 变3ʹ剪切 位 点（Alternative 3ʹsplice-site，A3S）、外 显 子 互 斥（Mutually exclusive exons，MEE），从Astalavista软件分析结果中，对转录本发生上述5种可变剪切事件情况进行统计。

1.3.4 融合基因分析融合基因是指将两个或多个基因的编码区首尾相连，置于同一套调控序列（包括启动子、增强子、核糖体结合序列及终止子等）控制之下，构成的嵌合基因。融合基因的表达产物为融合蛋白。使用Tofu（版本：13.0.0；参数：default）比对及寻找融合转录本，检测融合转录本的分析原理有：（1）比对到2个或多个位点；（2）每个位点必须比对上至少5%的转录本长度，最小比对长度为1 bp；（3）所有位点比对到的总长度必须占转录本总长度的95%以上；（4）两位点间距离必须达到10 kb以上。

2 结 果

2.1 全长转录组数据质控及分析通过Nanopore MinION平台测序，在两株细胞系（各3个重复样本）的6个文库中共得到5 801 060条有效reads，N50和平均读长分别为985 bp和815 bp，最大读长达到134 443 bp。过滤核糖体RNA后的reads数为4 284 539条，全长序列为3 796 820条，全长序列占总mRNA reads数的88.6%，详细信息见表1。

表1 全长序列数据统计表

2.2 差异基因与转录本表达分析本研究中两株细胞共注释到9 807个基因，包括643个已有参考基因组中未注释的新基因。获得28 532条转录本，其中已知转录本为26 061条，及与已有转录本具有结构差异的转录本2 471条。在MHCC97H vs.MHCC97L的比较分析中发现293个差异表达基因，包括9个新注释的基因。这些差异表达基因中，有116个基因在MHCC97H中高表达，剩余177个基因则在MHCC97H中低表达（图1A）。同时在转录本水平比较分析发现两株细胞共有74条差异转录本，其中10条为新发现的转录本。这些差异表达的转录本中有22条在MHCC97H中高表达，剩余52条则在MHCC97H中低表达（图1B），其中ITM2A的两个转录本仅在MHCC97H细胞中表达，而IL24的四个转录本仅在MHCC97L细胞中表达（图2）。

图1 MHCC97H和MHCC97L细胞的（A）差异表达基因热图和（B）差异表达转录本热图

图2 部分差异表达转录本在MHCC97H和MHCC97L细胞中的表达差异

差异表达转录本的KEGG通路注释发现，大部分差异表达转录本被注释在Human Diseases这个类别的肿瘤相关信号通路中（图3），其中有7条转录本注释到癌症蛋白多糖（Proteoglycans in cancer）通路，分别有5条转录本注释到结直肠癌（Colorectal can⁃cer）、癌症信号通路（pathways in cancer）及胶质瘤（Glioma）等通路。

图3 MHCC97H和MHCC97L细胞中差异表达转录本的KEGG数据库注释结果

2.3 可变剪切转录本分析鉴定在MHCC97H和MHCC97L细胞中共鉴定到1 008次可变剪切事件，其中MHCC97H细胞中鉴定到619次可变剪切事件，包括83次3’端可变剪切（Alternative 3'splice site，A3S）、87次5’端可变剪切（Alternative 5'splice site，A5S）、403次外显子跳跃（Exon skipping，ES）、28次内含子保留（Intron retention，IR）及18次外显子互斥（Mutually exclusive exon，MEE）；MHCC97L细胞中鉴定到870次可变剪切事件，包括120次3’端可变剪切（Alternative 3'splice site，A3S）、107次5’端可变剪切（Alternative 5'splice site，A5S）、572次外显子跳跃（Exon skipping，ES）、42次内含子保留（Intron retention，IR）及29次外显子互斥（Mutually exclusive exon，MEE）。两种细胞中各种可变剪切事件的占比基本一致（图4），可变剪切事件的具体信息详见表2。

表2 部分基因的可变剪切信息概要

图4 MHCC97H和MHCC97L细胞可变剪切分析结果

2.4 融合基因分析在MHCC97H和MHCC97L细胞系中，共鉴定出5条融合转录本（表3），这5条融合转录本均为新发现的融合突变类型，之前未见文献报道。其中包括肝癌中常见异常基因CTNNB1的融合转录本及两个新转录本ONT.5884与ONT5894形成的融合基因。

表3 融合转录本统计结果

3 讨论

纳米孔测序技术因具有超长读长的特点，使得其在全长转录组测序方面表现出独特优势。本研究通过纳米孔三代测序技术对具有不同转移潜能的肝细胞癌细胞系MHCC97H和MHCC97L进行全长转录组分析，共得到了5 801 060条Clean reads，N50和平均读长分别为985 bp和815 bp，最大读长达到134 443 bp，体现出纳米孔测序技术在鉴定肝癌细胞全长转录本方面的优势。通过生物信息学分析共注释到9 807个基因，包括643个新基因，获得28 532条转录本，其中已知转录本26 061条，本研究中鉴定出的新转录本2 471条。在这两种不同转移潜能的肝癌细胞系中存在着74条差异转录本，包括10条新发现的转录本。在结构变异分析中发现了1 008次可变剪切转录事件，其中占比最多的是外显子跳跃，这个结果与CHEN等利用Pacibio三代测序技术在肝细胞癌样本和MIHA细胞中的研究结果一致［11］。本研究结果为在转录本水平的基因表达调控参与肝细胞癌转移分子机制提供了初步的研究思路与线索。

由于二代测序技术短读长的限制，无法对转录本进行结构精确定量和差异表达分析。同一个基因转录而来的前体mRNA（pre-mRNA）通过可变剪切可形成不同的剪接异构体，最终形成不同的蛋白质而发挥不同的功能［22］。三代测序技术为深入研究转录本结构提供了强大的工具［23］，基于二代短读长测序数据只能进行基因表达量的计算和差异表达分析，但基于三代长读长测序数据不仅能够同时进行基因和转录本表达量的计算和差异表达分析，还能对基因和转录本的结构进行精确的可变剪切分析。可变剪切作为一种转录后的调控机制，在肿瘤的发生发展中经常出现异常［24］，比如促进肿瘤进展和转移的特异性可变剪切体［25-26］。与正常肝脏相比，肝细胞癌肿瘤组织中存在高度差异性可变剪切，其中许多可变剪切差异与肝细胞癌患者的存活率密切相关［27-28］。CHEN等通过Pacbio三代长读长测序技术对肝细胞癌患者样本和MIHA细胞系进行全长转录组分析鉴定出了肝癌细胞特有的isoform，研究结果体现了三代测序技术在鉴定可变剪切事件上的优势［11］。本研究通过纳米孔测序，在全长转录组水平全面系统地刻画了MHCC97H及MHCC97L细胞的可变剪切事件，为后续进一步研究可变剪切转录本的功能及肝癌转移分子机制提供了重要线索和数据信息。

本研究通过全长转录组测序更加精确地比较不同转录本在样本间的差异。MHCC97H与MHCC97L细胞是从MHCC97人肝癌细胞株再次分离培养筛选得到的具有高低转移潜能差异的两株细胞［29-30］。由于这两株细胞均来源于MHCC97，虽然具有不同的转移潜能，但总体差异较小，这可能是本研究仅筛选到74个差异转录本的原因。在这些差异表达转录本中，ITM2A基因的两个转录本仅在高转移肝癌细胞株MHCC97H中表达，而IL24基因的四个转录本仅在低转移肝癌细胞株MHCC97L中表达。有研究报道ITM2A表达与肝癌的预后密切相关，可能作为肝癌临床预后的标志物［31］，IL24基因能够抑制肺癌细胞的迁移和侵袭［32］。

本研究利用纳米孔三代测序技术对两种不同转移潜能的肝癌细胞系进行全长转录组测序，在转录本水平比较了两种细胞的表达差异，同时揭示了其可变剪切等基因结构变异。研究结果进一步体现出三代长读长测序技术在鉴定可变剪切事件上的优势，可变剪切转录本可能作为肝细胞癌治疗的新型潜在分子靶标，同时为后续进一步揭示肝癌转移的分子机制提供了新的思路和线索。