单分子全长转录本异构体测序揭示了心肌细胞中与疾病相关的 RNA 异构体
抽象的
选择性剪接产生不同的 RNA 异构体,控制真核生物的表型复杂性。它的故障是许多疾病的基础,包括癌症和心血管疾病。在全基因组范围内对 RNA 同种型进行比较分析一直很困难。在这里,我们建立了一个实验和计算管道,该管道执行从头转录本注释并准确量化来自 cDNA 序列的转录本异构体,全长异构体检测准确度为 97.6%。我们生成了一个可搜索的、定量的人类转录组注释,其中包含 31,025 个已知和 5,740 个新的转录亚型 ( http://steinmetzlab.embl.de/iBrowser/ )。通过在 RNA 结合基序蛋白 20 ( RBM20) 与侵袭性扩张型心肌病 (DCM) 相关的突变,我们确定了 107 个心脏基因中的 121 个差异表达的转录亚型。我们的方法能够对复杂的转录结构进行定量剖析,而不仅仅是识别单个外显子的包含或排除,如发现由 RBM20 突变错误剪接的IMMT同种型为例。因此,我们实现了独立于现有转录亚型注释的直接差异表达测试的途径,提供更直接的生物学解释和更高分辨率的转录组比较。
介绍
几乎所有的人类多外显子基因都被交替剪接1 , 2 , 3 , 4,允许单个基因产生多个 RNA 异构体,从而产生不同的蛋白质异构体5 , 6,从而驱动真核生物1 的表型复杂性。许多疾病,包括癌症、神经系统疾病和心血管疾病,都与选择性剪接失调有关,例如心脏特异性选择性剪接调节因子RBM20 的突变,导致扩张型心肌病 (DCM) 7 , 8 , 9. 研究选择性剪接的常规方法是短读长 RNA-seq 1 , 4 , 10 , 11 , 12 , 13,它对 RNA 片段进行测序,因此只能检测整个基因或单个外显子的表达变化。尽管 Kallisto 14或 Salmon 15等软件工具允许使用短读长 RNA-seq 进行转录水平定量,这需要参考注释或从头转录组组装。后者是发现新同种型所需的,从短读数据中尤其具有挑战性。因此,短读长 RNA-seq 对于定量检测单个 RNA 异构体并不理想。包括Pacific Biosciences(PacBio)和Oxford Nanopore Technologies(ONT)在内的长读长测序技术的出现,为更全面地分析可变剪接提供了一种新工具。PacBio 的一项初始全基因组研究使用了 476,000 个读数16,这使得从头转录物识别成为可能,但不足以进行准确的量化。因此,后续的用PacBio或ONT研究集中于使用富集方法量化所选择的基因,缺少基因组范围的型材17,18,19,20。最近的研究已经使用长读测序小鼠量化转录物21,22和乳腺癌细胞23,然而样品之间没有定量差异进行了评估。此外,长读长测序刚刚开始在单细胞转录组学中用于研究转录异质性24 , 25. 应用长读长 RNA-seq 来确定转录本的全局表达变化是具有挑战性的,因为:(1) 低测序深度通常不允许足够的基因组覆盖进行准确定量;(2) 高测序错误率需要新的计算方法来实现转录本量化、分类和可视化;(3) 由于测序伪影,例如逆转录酶 (RT) 的不希望的模板转换和寡核苷酸 (dT) 21 的脱靶引发,已识别的未注释转录本的错误发现率很高混淆下游分析。因此,迄今为止,还不可能进行定量、全基因组差异表达分析来大规模检测同种型变化。这使得很难确定许多涉及剪接的疾病的分子机制。在这里,我们证明了单分子全长 RNA 测序能够识别疾病相关的转录亚型。
结果和讨论使用 ONT 测序生成长读长 RNA-seq 数据集
我们使用 cDNA ONT 测序为人诱导多能干细胞衍生的心肌细胞 (iPSC-CM) 生成了一个大型长读长数据集,并开发了一个工作流程,该工作流程可准确定量地测量和比较基因组范围内的全长剪接异构体,而不依赖于现有注释。来自具有和不具有 DCM 相关RBM20突变的人类 iPSC-CM(每个具有两个独立克隆)的 RNA 以及掺入对照(ERCC 26和亮片27),被转化为全长 cDNA 并使用 ONT MinION 技术进行测序(方法部分、补充图 1和补充表 1)。我们总共生成了 2100 万个高质量读数(平均 qscore ≥6,方法部分,补充图 2),并量化了总共 11,707 个基因的 36,765 个转录亚型,覆盖了蛋白质编码的 53.8%(10,682/19,847)基因在人类基因组(图 1,补充图 3,和补充表格 1 - 3)。相比之下,我们具有约 12 亿对末端读数的互补短读数据识别了 16,726 个蛋白质编码基因(读数 >1)。
图 1:全长人类 iPSC-CM 转录组的复杂图景。
以 19 号染色体区域 (chr19:57531255-57906445) 为例,对人类 iPSC-CM 中全长剪接亚型的全基因组测量。基因位点表示为水平黑线,分别在基因组轴上方或下方的 + 或 – 链上带有注释基因。折叠 GENCODE 全面和基本的注释以深紫色和紫色轨道呈现。已知的转录异构体显示为绿色轨迹,以前未识别的转录异构体显示为红色。已知转录本异构体未通过转录本过滤标准以浅绿色显示。每个识别出的转录本的阅读计数以圆圈中的数字给出。基于新的剪接事件(用灰色框和文本符号表示的位置),新的转录本异构体被分类为新的外显子组合(NC)、新的外显子(NE)、新的内含子保留 (IR)、新的外显子跳跃 (ES) 或新的替代剪接位点 (ASS)。不显示通过成绩单阅读的可疑小说。
全尺寸图片开发 FulQuant 以定义全长转录组当前使用长读长测序定量同种型的方法的主要问题是测序伪影和导致高错误转录发现的高测序错误。为了应对这些挑战,我们建立了全长转录本量化 (FulQuant) 的计算方法,该方法集成了 ONT cDNA 序列的准确全长读取识别、转录本量化和可视化。这种新方法允许从头转录本注释,并采用严格而复杂的标准过滤读数、比对和转录本,以排除 RT 和测序错误中的伪影,并根据其剪接位点将读数分组到转录本中,生成一组高度可信的具有丰度估计的转录本(方法部分和补充图 4)。为了评估 FulQuant 识别同种型的准确性,我们确定了其在合成测序掺入序列(亮片,掺入占总 RNA 的 3%)上的性能,其浓度范围约为 10 6倍,涵盖了观察到的基因表达的动态范围跨越人类转录组27. 这些spike-ins 代表72 个人工基因位点,具有多个编码156 种替代异构体的外显子,以模拟不同序列、结构、丰度和长度的人类剪接异构体。我们鉴定了 156 种亮片同种型中的 82 种,包括 50% 最集中的同种型中的 90% (70),占所有亮片 RNA 分子的 99.95%。我们的同种型鉴定真阳性率,定义为鉴定的参考转录物与所有报告的转录物的比率,基于所有亮片同种型为 97.6%,最大的亮片同种型为~7 kb(图 2a)。总之,人类 GENCODE 注释中所有剪接转录本的 98.1% 低于这个 7 kb 阈值(图 2b)。对于短于 1728 bp 的同种型,真阳性率为 100%,其长度大于所有人类剪接转录本的 70%。FulQuant 的识别准确度高于现有的长读取管道 FLAIR 28,当使用默认设置(补充文本)对我们的亮片数据集进行基准测试时,它产生了 30% 的真阳性率。FulQuant 在转录本识别方面的改进主要是由于去除了合成对照中出现的 RT 和测序中的人工制品造成的假阳性转录本。使用 FulQuant,亮片的量化与输入浓度呈线性相关(平均 Pearson 相关系数 0.85,补充图 5)。此外,98.8% 的转录本边界位于带注释的转录本的 20 bp 内(补充图 6),表明可以可靠地捕获全长转录本,而不依赖于注释。同样,当我们将 FulQuant 应用于我们的人类 iPSC-CM 数据集时,我们观察到技术和生物学重复之间的高度相关性(补充图 7),以及全长转录本的准确边界定义(补充表 2和3)。FulQuant 的量化在生物复制之间产生了比 FLAIR 更好的一致性(补充图 7),这可能是由于假阳性转录本的数量减少了。
图 2:FulQuant 方法准确识别 iPSC-CM 中的转录异构体。
a FulQuant 方法的长度相关精度(绿色)和召回率(蓝色),使用全部或最丰富的 50% 合成亮片对照进行估计。b GENCODE 中转录本、我们鉴定的同种型和亮片对照的长度分布。a、b与相同的 x 轴垂直对齐,以说明亮片对人类转录本估计的准确性的相关性。c基于 GENCODE 综合注释和本研究中鉴定的转录本(分为已知、新颖和全部),每个基因的同种型数量分布。d包含重复序列的已知外显子的比例与已识别的新外显子的比例相比。电子按新剪接事件类型划分的新转录本异构体的百分比:新外显子 (NE)、新内含子保留 (IR)、新替代剪接位点 (ASS)、新外显子跳跃 (ES) 和已知外显子的未注释组合 (NC) . 模棱两可的成绩单是那些有对齐问题的成绩单。f每个新转录物类别中预测的非 NMD 产品的比率。g , h分别具有已知外显子和新外显子未注释组合的基因示例。综合 GENCODE 注释中的已知同种型以灰色轨迹显示。未显示具有 <5 个读数的同种型。源数据作为源数据文件提供。
全尺寸图片人类 iPSC-CM 全长转录组的全局分析我们的数据揭示了已知和新型转录本亚型的全基因组复杂可变剪接模式。通过将我们的转录异构体与人类 GENCODE 注释 (V24) 进行比较,我们确定了 36,765 个转录异构体,其中 5740 个(15.6%)是新的。例如,我们鉴定了锌指基因ZNF211 的八个已知和四个新转录本(图 1)。与综合注释集 CHESS (v2.2) 29 相比,3028 (8.2%) 是新颖的(补充文本和图 8)。我们的全基因组全长同种型数据集可通过可搜索浏览器获得,该浏览器是我们开发并免费提供的 ( http://steinmetzlab.embl.de/iBrowser/ )。
我们鉴定的转录本异构体的长度从 177 bp 到 10,718 bp,中位数为 1768 bp(图 2b和补充图 9)。长度分布类似于 GENCODE 全长蛋白质编码和 lincRNA 转录本(补充图 10),但是我们没有捕获大的转录本,例如TTN同种型(大小约为 100 kb,心肌的主要结构成分)可能是由于 cDNA 合成的 RT 限制。确定了每个基因的两个可变剪接转录本同种型的中位数(图 2c)。值得注意的是,通过将我们的方法仅应用于一种细胞类型,即 iPSC-CM,我们已经鉴定了 24.1% 的蛋白质编码同种型(n = 80,901) 涵盖 GENCODE 综合注释中的 9946 个基因,该注释累积了来自各种人类细胞系、组织和群体的数据。这反映了我们的长读长方法的综合性质和深度,并表明人类基因组中所有剪接异构体的数量和复杂性可能被低估了。
异常剪接可能会产生新的转录本,可用于诊断和药物开发。值得注意的是,在我们的数据集中检测到的新外显子富含由 RepeatMasker 定义的重复序列(图 2d),其中 21% 的新外显子与 SINE/Alu 序列重叠。除了更高程度的内含子保留(IR,4.1%)外,大多数新型转录本异构体还包含替代剪接位点(ASS,29.2%)、新型外显子跳跃(ES,12.7%)、新型外显子(NE,8.6%)、和未注释的已知外显子组合(NC,42.9%,图 2e,g,h)。为了评估我们管道的性能,我们使用 RT-PCR 的独立片段分析成功验证了 33 个新剪接事件(外显子跳跃、内含子保留和新外显子)中的 26 个(方法部分,补充图)。 11和12以及补充表 6)。在蛋白质编码潜力方面,我们确定的新型转录本异构体平均每四个基因添加 1 个潜在的蛋白质编码异构体(图 2f和补充文本)。支持我们的数据的分辨率,我们还确定从相邻基因(补充图45个潜在剪接多顺反子转录物 13被从线粒体基因组(补充图转录)和未剪接转录多顺反子 14)30,31。
人类 iPSC-CM 中协调外显子使用事件的分析外显子的协调表达是许多调节机制的基础23。它的分析也提供了重要的见解选择性剪接的策略32,33。特别是在心脏中,TPM1、2和 3等基因中的互斥外显子对已被用作模型系统来研究剪接调节34。我们的长读长数据首次系统地分析了心脏细胞中整个转录本亚型的相邻和远距离外显子之间的相互依赖性。我们在 722 个基因中鉴定了 2793 个显着共调控的外显子对(调整后的P < 0.001,相互排斥和包含关联,方法部分和补充表 4)。例如,我们在心脏特异性基因TPM1 中检测到有据可查的外显子对关联,该基因包含相邻的互斥外显子 2 和 3,它们代表平滑肌和骨骼肌35之间的肌肉类型转换。此外,我们在TPM1转录本亚型中鉴定了 6 个进一步的共关联事件,包括 4 个远端对,例如外显子 2 和 11 之间的相互包含对和外显子 7 和 11 之间的互斥对(图 3a)。另一个突出的例子是MYL7(图 3b),一个心脏功能的关键基因,包含 18 个相互包含 (6) 和排斥 (12) 的外显子对。该基因中所有测试的外显子对 (25) 中超过 70% 显着相关,并且大多数 (61%) 位于远端,表明该基因的剪接调控复杂。总的来说,我们观察到一个或两个外显子内的相邻共关联比远端共关联更常见(补充图 15)。由于新生 RNA 的剪接以极快的速度共转录发生(“50% 的剪接在 3' 剪接位点合成后约1.4 秒内完成”)36,我们推断这在理论上有利于相邻的新合成的快速相关剪接外显子。
图 3:人类 iPSC-CM 中的外显子相关事件。
TPM1 ( a ) 和MLY7 ( b ) 显示为包含多个相邻和远端相互包含或排斥的替代外显子的基因的例子。替代外显子对之间的共关联强度(p值的log10 )在互斥(红色,-)和包含对(蓝色,+)的三角形图中说明,而灰色表示未测试的外显子对。紫色和灰色注释分别用于折叠和展开的 GENCODE 综合注释。蓝色读数与图1 的对数刻度相同 暗度与读取计数相关。仅显示了用于我们的共关联分析的外显子。为清楚起见,未显示低丰度的转录本异构体。源数据作为源数据文件提供。
全尺寸图片鉴定RBM20突变体中错误拼接的全长转录本亚型异常剪接与多种疾病有关,但确定致病事件具有挑战性。我们使用全长转录本测序的方法为研究不同状态之间的差异异构体表达提供了机会,这将使健康和疾病之间的直接比较成为可能。要制定一个量化的方式来确定样品类型之间的差异表达的转录亚型,我们应用我们的管道研究全长剪接异构体是如何通过在心脏特异性剪接突变的调节影响RBM20 37,38,39。20元突变(例如 R634Q 在其富含精氨酸/丝氨酸的结构域内)导致 DCM 的侵袭性形式,并且短读长测序表明RBM20靶外显子在大鼠37、小鼠40、猪41和人类心脏细胞中被错误剪接38 , 42. 然而,在这些目标中的许多目标中,尚不清楚这些差异剪接的外显子属于哪种异常剪接的转录亚型,以及是否可能会或可能不会产生故障蛋白质,或者可能只是导致正确剪接蛋白质的缺失。除了-6 M读取为R634Q突变体,我们还产生9M读取的P633L突变体,一个新的致病突变(脯氨酸到亮氨酸改变在氨基酸位置633),我们最近发现42,43. 我们使用 R634Q 和野生型数据进行了从头转录本注释。虽然这种方法不能识别仅在 P633L 突变体中看到的新同种型,但它使我们能够定量比较所有突变体中在 R634Q 和野生型中检测到的同种型的表达水平。为了将携带突变体 R634Q 和 P633L 的 iPSC-CM 的全长转录组与野生型RBM20进行比较,我们使用 DESeq2 44对量化的转录本进行了差异异构体表达分析(补充表 3)。我们在 107 个差异表达的基因中鉴定了 121 个转录亚型(调整后的P < 0.001,图 4a,补充表 5,和方法部分)。两种突变体都显示出相似的表达水平(相关系数 = 0.77),并且在使用差异表达分析相互比较时没有显着差异(补充图 16)。通过使用 Kallisto(方法部分、补充文本和补充图17和18)对短读长数据进行相同的分析运行,并使用 Kallisto 进行转录同种型定量,我们对通过长读长测序鉴定的超过 80% 的候选进行了概括 。对这 107 个基因的基因本体分析显示心脏功能富集,如横纹肌发育、心脏收缩调节、离子转运和基于肌动蛋白丝的过程(调整后的P < 0.05),支持了心脏功能中的RBM20(补充文字和补充图 19)。
图 4:RBM20突变体 R643Q 和 P633L 中差异表达的全长转录本异构体。
WT 和突变体中平均全长同种型表达水平的散点图。每个点代表一个异构体,红点突出显示差异表达的异构体(调整后的P < 0.001)。对于异构体IMMT基因标记。b含有差异表达同种型的基因。顶部,每个同种型(三角形)的差异表达分析的显着性水平,显着性为红色。底部,每个同种型(点)的表达水平的倍数变化,红色显着。红色的基因名称是基因差异表达仅占其同种型的比例,在总基因表达水平上没有显着变化。仅显示具有三个以上表达同种型的基因。ç IMMT基因与 GENCODE 注释,短读数据和本研究的长读数据中鉴定的转录亚型。从上到下的轨迹:IMMT基因座周围 GENCODE 转录本的折叠视图、IMMT GENCODE 转录本的扩展视图、突变体和 WT 的短读覆盖、短读数据中推断的突变体和 WT 之间的 delta-PSI,以及已知(绿色)和在长读数据中识别的新(红色)转录本。GENCODE 转录本类型用红色(蛋白质编码)、蓝色(处理过的)和金色(哈瓦那合并)进行颜色编码。编码区域由 GENCODE 注释轨迹中的填充颜色突出显示。d RBM20中IMMT同种型表达的并排比较R634Q、P633L 突变体和 WT。两种新的红色转录本同工型IMMT-SL1和IMMT-SL2仅在 WT 样品中表达。已知的 GENCODE 转录本为绿色,原始长读为蓝色。错误拼接的外显子 6 用星号表示。IMMT-205和IMMT-208由于读取计数低(所有样本中为 1)而被省略。源数据作为源数据文件提供。
全尺寸图片值得注意的是,并非所有具有差异表达转录本同种型的 107 个基因在总基因表达水平上都有显着变化,如使用标准差异基因表达分析评估的那样(方法部分和图 4b))。在 107 个基因中,27 个表达单一同种型,其余 80 个具有多种同种型的基因中,69 个基因显示总基因表达差异。在这 69 个基因中,60 个基因在主要同种型表达上存在差异(> 基因总转录本的 50%),而 9 个基因仅差异表达次要同种型。具有多种同种型的 11 个基因的总基因表达没有变化,但差异表达了一种或多种同种型。在评估转录本异构体的使用时,这 11 个基因中有 8 个显示突变体和野生型之间异构体比率的显着差异(Fisher 精确检验,调整后的P值 < 0.005,方法部分和补充图 20),清楚地表明剪接水平的复杂调控,独立于启动子或同种型稳定性差异,因为总基因表达不受影响(补充表 5)。这种复杂的同种型调控表明,在基因水平而非同种型水平上测量差异 RNA 表达可能无法检测到特定同种型的显着表达变化,不需要事先注释的方法是有利的。具有特定转录本同种型差异表达(而非总基因表达)的基因的两个突出例子是TPM1 37和IMMT 38,之前报告说两者都是从大鼠和/或人类心脏转录组的短读测序中错误拼接的。我们的长读长数据清楚地表明,并非所有同种型在突变体中都受到同等影响(图 4b)。对于IMMT,我们对短读数据外显子包含(PSI)分析的百分比确定仅外显子6与野生型和突变体之间的差异剪接的,如在以前的研究中报道的37,38,41,42(图 4C和补充文本)。然而,在全长转录亚型误剪接的作用不是由以前的短读分析捕获的37,38。
具体来说,目前尚不清楚有多少和哪些转录本受到RBM20突变的影响。对于IMMT,基于参考 GENCODE 注释,假设转录本IMMT-205是唯一受影响的同种型是合乎逻辑的,因为这是唯一没有外显子 6的IMMT注释转录本。基于短读数据的转录本量化差异表达分析确实将IMMT-205报告为IMMT(方法部分)的唯一重要影响。这是一种误导,因为我们基于长读长的方法可以清楚地识别两种新的IMMT转录亚型,IMMT-SL1和IMMT-SL2,与它们相应的已知转录本IMMT-207和IMMT-204相比,缺少外显子 6 。这两种同工型在野生型中表达,在 R634Q 中不存在,仅在 P633L 中略有表达,与三种未差异表达的同工型形成对比(图 4a、d)。在所有样品中仅鉴定了IMMT-205 的一个读数。基于域分析,我们发现外显子 6 的缺失可能导致蛋白质具有新功能,因为该外显子编码的 32 个氨基酸的缺失破坏了蛋白质中一个现有的卷曲螺旋结构域和一个固有的无序区域(补充图. 21 )。因此,新型IMMT的发现 isoforms 说明了无需事先注释的全长转录组识别和定量的效用。
讨论
总之,我们提出了人类细胞全长转录本亚型的全基因组分析及其在关键心脏基因RBM20中的剪接突变背景下的差异分析。虽然在剪接缺陷的RBM20突变37、38、40、41、42 中描述了不同的外显子使用,但我们发现了选择转录本上复杂的异构体失调,这可能反映了启动子活性、异构体稳定性和单个基因的可变剪接变化的复杂组合。在 RBM20 突变体中。我们的 FulQuant 工作流程允许直接识别和量化转录本异构体,无需事先注释。以两本小说为例IMMT转录本异构体是交替剪接的、无注释的全长异构体测序,对于识别以前遗漏的许多新的、信息丰富的转录本异构体的存在至关重要。这些观察结果,结合人类转录本异构体景观的整体、复杂的复杂性(如图 1 所示)和在基因组浏览器中),证明了更广泛地采用长读长测序进行转录组研究的迫切需要,以提供先前确定的可变剪接外显子的背景。了解全长同种型对于理解表达基因的功能产物至关重要。我们直接量化和比较全长转录本同工型的工作流程,包括各类新型转录本,可以很容易地扩展到其他细胞类型和组织,以及已提议异常剪接在其中发挥作用的许多疾病。
方法人类 iPSC 中的基因组编辑
我们的出版物42详细介绍了基因组编辑和 iPSC-CM 生成的方法。简而言之,我们从斯坦福心血管研究所生物库获得了人类 iPSC。将退火(T4 连接缓冲液,NEB)和磷酸化(T4 PNK,NEB)引导 RNA(gRNA,使用http://crispr.mit.edu 上的工具设计)寡核苷酸引入 pSpCas9(BB)-2A 的 BbsI 位点-GFP 质粒并转化到 STBL3大肠杆菌细胞中并通过 Sanger 测序验证。将细胞以低密度铺在基质胶包被的 6 孔板中,在 Essential 8 (E8) 培养基中接种一天后,将培养基更换为含有 Rock 抑制剂(Tocris Cat. No. 1254)的 E8。对于每个孔,CRISPR/Cas9 载体(1 μ,pSpCas9(BB)-2A-GFP)和 4 μg 单链 DNA 供体(补充信息)通过 Lipofectamine 3000 转染引入细胞。转染后 36-48 小时使用 FACSAria IIu(DB Biosciences)流式细胞仪与 100-μm 喷嘴分离 GFP + 细胞,用 Rock 抑制剂维持在 E8 中,用于前 3 天,初始密度为 2–3 × 10 3细胞/孔,随后转移到常规 E8,直到集落大小达到 ~0.5 mm。单独的 iPSC 克隆被分离并重新接种到 E8 中 24 孔板的孔中,并含有 Rock 抑制剂。在目标基因组区域的 PCR 扩增后通过测序确认纯合编辑(补充信息)。R634Q 电池可应要求提供。
iPSC-CM 生成通过调节 WNT 信号,iPSCs 分化为单层心肌细胞。iPSCs 以低密度接种在 Matrigel 涂层板上,使其在 4 天后(分化的第 1 天)达到 70-80% 汇合。用 3 ml RPMI 1640 (Life Technologies 11875-093) 和 1X B27® Minus 胰岛素 (Life Technologies 0050129SA) 诱导分化,并补充 6 μM CHIR (TOCRIS 4953)。在第 6 天,将培养基更换为 3 ml 含有 1X B27® Minus Insulin 的 RPMI 1640,并补充有 5 μM IWR。从第 8 天起,将细胞保存在含有 1X B27® 无血清补充剂(Life Technologies 17504-044)的 RPMI 1640(Life Technologies 11875-093)中。从第 12-15 天开始,用含 1X B27 的 RPMI 1640 无葡萄糖(Life Technologies 11879-020)处理细胞,然后在含 1X B27 的 RPMI 1640 中恢复 2 天,并随后以 3Mi 细胞/孔的密度重新种植。然后将细胞在含有 1X B27 的 RPMI 1640 中再保持 4 周,然后收集用于实验。
逆转录使用 ZymoResearch RNA Clean and Concentrator-5 清洗和处理 RNA。总共 4.5 μL 混合物,含有 5 ng 纯化的 iPSC-CM RNA,含 2.8% ERCC 和 3% Sequins A 版,1 μL 10 μM oligodT(/5SpC3/A*A*G*CAGTGGTATCAACGCAGAGGTACTTTTTTTTTTTTTTVTTTTTTTTTTTTTTTTTTTTTT1TTL dNTP、0.1 μL Takara 重组核糖核酸酶抑制剂 (40U/μL) 在 72°C 下孵育 3 分钟,在 4°C 下孵育 10 分钟,然后在 25°C 下孵育 1 分钟。另一种 5.5 μL 混合物,包含 2 μL 5x SuperScript II 缓冲液、2 μL 5 M 甜菜碱、0.5 μL 100 mM DTT、0.5 μL SuperScript II (200 U/μL)、0.25 μL Takara 重组核糖核酸酶抑制剂、0.1 μL 100 μM TSO、0.09 μL 无核酸酶水和 0.06 μL 1 M MgCl 2加入 10 μL 的最终反应体积。最后的 10 μL RT 反应在 42 °C 下孵育 90 分钟,50 °C 下 2 分钟和 42 °C 下 2 分钟4进行 10 个循环。然后通过在 70°C 下孵育 15 分钟并保持在 4°C 来停止反应。然后用 1 μL 1:10 稀释的 NEB RecJf(目录号 M0264S)在 37°C 下处理 30 分钟,在 65°C 下处理 20 分钟。
全长 cDNA 的 PCR 扩增将 25 μL 2x NEB Q5 HotStart 主混合物、0.25 μL 10 μM PCR 引物(ONT_Index1_ISPCR 和 ONT_Index2_ISPCR)5和 14.5 μL 无核酸酶水加入 RT 反应中。PCR 反应如下:98 ℃ 30 s,98 ℃ 10 s,67 ℃ 15 s,72 ℃ 6 min,72 ℃ 2 min,20 个循环。 4℃。PCR 反应使用 0.75 体积的 Ampure XP 珠进行纯化,并用 25 μL 无核酸酶的水洗脱。使用高灵敏度 DNA 芯片在安捷伦生物分析仪上分析全长 cDNA。
牛津纳米孔测序按照制造商的描述,总共将 1 μg cDNA 与 ONT 测序接头连接,并在以下步骤中进行了修改: (1) 末端修复和 dA 拖尾;(2) ONT测序接头的连接;(3) 添加测序系链。使用 MinION 设备上的 R9.4 和 R9.5 流通池进行测序。数据用 ONT 的 MinKNOW (1.4.2) 软件记录 48 小时。碱基识别是使用 ONT 的 Albacore 软件(版本 2.1.3)与选项“-disable_filtering -kit SQK-LSK108”和相应的流动池版本(补充表 1)进行的。
FulQuant 用于全长转录本异构体定量除非另有说明,否则分析是在 R (3.5.3) 中或使用自定义 Bash 脚本进行的。我们评估了对齐识别 (PID) 百分比对平均读取质量(由 ONT 定义的 mean_qscore)的依赖性,并确定了平均读取质量得分截止值 6,其中超过 80% 的读取具有 PID > 80%(补充图 2))。在去除低质量读数后,我们在读数的两端(200 bp 窗口)修剪了测序接头 ISPCR(AAGCAGTGGTATCAACGCAGAGTAC)和 polyA 序列(≥10 个连续 A/T,允许一个不匹配)。一端具有接头序列而另一端具有polyA 序列的读数被认为是全长的。我们使用 minimap2(版本:2.9-r751-dirty)和参数:'-K500m -secondary=no -a -x splice -splice-flank=yes' 将读数与人类基因组 (GRCh38) 对齐。从下游分析中删除满足以下标准的比对:(i) 由补充比对组成,(ii) 两端的软剪切大小大于 30 bp,(iii) 剪接点附近的插入大小(由雪茄 N 识别)> 5 bp 和 (iv) 未拼接对齐(单例)。
转录组注释是对 WT 和 R634Q 的数据进行的。我们使用先前描述的聚类方法确定了内含子-外显子和外显子-内含子连接的标签剪接位点45使用参数'HFWINSIZE = 5 DISTHRES = 8 PTHRESHOLD = 3 SUMHOLD = 5'。共享相同标签剪接位点的读数被折叠成一致的转录本异构体。我们包括末端外显子,但不考虑 5' 和 3' 端的转录变异,因为读取两端的错误率很高。事实上,对于末端外显子,我们精确地注释了它们的剪接位点,而不是它们的转录起始(TSS)和终止(TES)位点。我们使用一种转录本同种型中所有读数的中位数 5' 和 3' 位置作为其 TSS 和 TES 位点。以下规则用于去除假阳性转录本:(1) 读取计数 >3,(2) 全长读取 >40%,(3) 转录终止位点 >10 bp 来自已知基因组 polyA 序列,(4) 链信息可用。此外,我们要求在每个基因座中的所有转录本:(1) 有超过 3% 的最丰富转录本的读取计数,(2) 不是错配的结果,以及 (3) 不是任何已识别转录本的 5' 截断产物,除非计数高于 25 % 较长的转录本。用于执行上述步骤的代码附在 带有详细注释的补充文件。对于最终的注释,我们还删除了只有 2 个外显子或可能是已知转录本的截断产物的转录本。
我们将我们的注释与 GENCODE 人类综合注释 (V24) 进行了比较,以将转录本分类为已知的和新颖的。映射到已知 GENCODE 注释的过滤转录本被拯救。我们进一步将新转录本注释为新外显子、新组合、外显子跳跃、内含子保留、替代剪接位点和多个事件。可以在补充信息和代码中找到每个新转录本类别的详细定义 。
外显子连接分析我们首先将所有转录本中存在的外显子定义为组成型。对于基因内的所有外显子对,我们计算了以下情况下的读取数:(1) 两个外显子都存在,(2) 只存在外显子 1,(3) 只存在外显子 2,以及 (4) 两个外显子都存在不存在。使用这个计数矩阵,我们在使用BH 方法调整p值之后使用 Fisher 精确检验测试了它们的共关联。远端对被定义为至少相隔一个组成型外显子的外显子对。
预测无意义介导的衰变事件对于包含在 GENCODE 注释中定义的已知起始密码子的所有新转录本,我们计算了转录本上第一个终止密码子和最后一个剪接点之间的距离。如果这个距离 >55 bp 46,我们将转录本定义为 NMD 产物。
短读数据分析使用 STAR v2.5.1b 将读数映射到 GRCh38 后,使用 featureCounts v1.6.0 确定基因表达水平。对于每个外显子,我们分别将包含和排他的读数计算为包括感兴趣的外显子的读数和包括上游和下游外显子但不包括感兴趣的外显子的读数。计算了表示包含读取与包含和排除读取总和之间的比率的拼接 (PSI) 百分比。使用 Kallisto (0.46.2) 与 GENCODE 或我们的注释作为索引的输入注释进行转录量化。
基因和转录本异构体水平差异表达分析通过计算共享该转录本相同剪接位点的读数数量来实现转录本量化。通过总结该基因座中的转录物计数来实现基因定量。我们将 CRISPR 生成的 R643Q 和 P633L 突变样本与未经 CRISPR 编辑处理的样本 (WT_NC) 和未经 CRISPR 编辑处理的样本 (WT) 进行了比较,其中使用 DESeq2 (1.26.0) 对基因和转录本进行标准差异表达分析异构体。由于测序深度有限,我们只考虑了所有样本中总read数超过10的基因和转录本异构体来测试长读长数据。使用 BH 调整的p定义重要候选者突变体和 WT 之间的比较值截止值为 0.001。对于短读长数据,使用 DESeq2 的相同测试程序应用于使用 Kallisto 获得的转录本和基因量化,无需任何过滤。GO富集分析使用PANTHER通过http://geneontology.org的网络界面进行,表达的基因作为参考列表以控制背景。为了评估突变型和野生型之间转录本使用的差异,我们总结了生物学重复并对具有一种或多种差异表达同种型的 80 个候选基因进行了 Fisher 精确检验。使用 BH 程序调整P 值。
使用片段分析验证新的剪接事件我们使用片段分析来验证 33 个新的剪接事件(补充表 6和补充图 10)。设计针对新剪接事件的寡核苷酸引物,使得具有和不具有剪接事件的扩增子将在凝胶上显示可检测的大小差异。所有引物序列都可以在补充表6 中找到 。我们进行了标准的 RT-PCR 并使用标准的 Bioanalyzer 分析分析了片段大小。只有包含具有预期大小的条带的配置文件才被视为阳性验证。对于 7 个不成功的案例,3 个案例没有产生任何 PCR 产物,4 个案例包含与预测大小相差至少 20 bp 的片段。
报告摘要有关研究设计的更多信息,请参见与本文链接的 自然研究报告摘要。
数据可用性
本研究中生成的长读长数据已存放在 ArrayExpress 中,登记号为E-MTAB-7334。本研究中使用的短读取数据已存放在序列读取存档 (SRA) 中,登录号为PRJNA579336 42。我们的全基因组全长同种型数据集可在http://steinmetzlab.embl.de/iBrowser/ 上作为自定义基因组浏览器使用。 本文提供了源数据。
代码可用性
GitHub https://github.com/czhu/FulQuant (doi: 10.5281/zenodo.4960275)上提供了执行转录过滤基本步骤的 R 代码。
参考1.
王,ET 等。人体组织转录组中的替代异构体调节。自然 456 , 470–476 (2008)。
ADS CAS PubMed PubMed Central Article Google Scholar
2.Pan, Q., Shai, O., Lee, LJ, Frey, BJ & Blencowe, BJ 通过高通量测序深入调查人类转录组中可变剪接的复杂性。纳特。基因 40 , 1413–1415 (2008)。