一种单细胞转录组RNA污染去除方法、介质和设备

相关申请的交叉引用

本申请是基于申请号为2022113640956，申请日为：2022年11月02日，发明名称为：一种单细胞转录组RNA污染鉴定方法、介质和设备的分案申请。

技术领域

本发明涉及生物数据分析方法，尤其涉及一种单细胞转录组RNA污染去除方法、介质和设备。

背景技术

基于微流控技术的单细胞转录组测序能够在单个实验中对数万个细胞的基因表达进行量化。其主要是基于序列标签来识别单细胞，其核心技术是给每个细胞添加一个独特的序列标签，在测序时把携带相同标签的核酸序列视为来自同一个细胞。10X Genomics单细胞转录组测序平台为目前应用广泛的一种技术，该平台利用微流控、油滴包裹和barcode标签等技术来实现高通量的细胞分选与捕获，能够一次性分离、并标记500至数万个单细胞，测序后可获得每个细胞的转录组信息，具有细胞通量高、建库成本低、捕获周期短等优势。

典型的单细胞转录组测序实验流程如下，首先制备细胞悬浮液，在相应的平台仪器上利用微流控芯片将细胞悬浮液与磁珠混合，并用油滴包裹。每个微珠带有独特的核苷酸序列，即barcode标签，可以对单个细胞进行标记。而每条barcode标签还连有一段由核苷酸序列组成的分子标识符(unique molecular identifier，UMI)，每个UMI可以标记一个mRNA转录本。经过逆转录、PCR扩增、生成文库和测序，在测序数据中，根据barcode标签和UMI标记，即可确定结果中每条序列是否来自于同一个细胞、同一个mRNA，这一方法可以减少PCR对不同分子的偏好影响。通过匹配和计数barcode和UMI，汇总基因表达信息在计数矩阵中，由此获得单个细胞的转录组表达谱。

基于液滴的单细胞转录组数据分析的前提是假设所有获得的RNA序列都是细胞内源性的。然而，在实际实验过程中，难免混入其他非细胞RNA，这些RNA可能来源于样本中已经破裂甚至破碎的细胞，或者来源于试剂溶液中，也可能来源于实验人员等环境因素。实验过程中，非细胞RNA与内源性细胞RNA会被一起包裹在液滴中，可能与细胞RNA一起被标记、扩增和计数。这导致不同细胞群之间转录本的交叉污染，构成了单细胞转录组数据的背景污染，使单细胞转录组数据的生物学解释偏离真实情况。目前还未有系统判定是否存在RNA污染的方法。

发明内容

为了解决背景技术中提到的至少一个技术问题，本发明的目的在于提供一种单细胞转录组RNA污染鉴定方法、介质和设备，分析单细胞转录组数据中是否存在背景RNA污染。

为实现上述目的，本发明提供如下技术方案：

一种单细胞转录组RNA污染鉴定方法，包括以下步骤：

S1，根据基因表达量从大到小对barcode进行排序；

S2，求得基因表达量变化幅度由小变大的A类拐点和基因表达量变化幅度由大变小的B类拐点；

S3，取规定区间内对应基因表达量最大的A类拐点为起点拐点，对应基因表达量最小的B类拐点为终点拐点；将起点拐点到终点拐点之间的所有barcode设为背景barcode；

S4，提取背景barcode的表达谱，统计不同基因在背景barcode中的表达比例；将表达比例大于污染阈值GP的基因设定为污染基因；

S5，若存在污染基因，则判定单细胞转录组数据存在RNA污染。

进一步的，所述A类拐点的求解方法如下：

SA21，以barcode的排名为X轴，基因表达量为Y轴绘制散点图；

SA22，在所述散点图上，间隔规定距离取最近的一个点，求相邻两个点之间的斜率；

SA23，当斜率的变化趋势为由小到大，且在该趋势持续过程中，斜率首次小于设定的斜率阈值时，将对应的点设为A类拐点。

进一步的，所述B类拐点的求解方法如下：

SB21，以barcode的排名为X轴，基因表达量为Y轴绘制散点图；

SB22，在所述散点图上，间隔规定距离取最近的一个点，求相邻两个点之间的斜率；

SB23，当斜率的变化趋势为由大到小，且在该趋势持续过程中，斜率首次小于设定的斜率阈值时，将对应的点设为B类拐点。

进一步的，在S21绘制散点图之前，对barcode的排名和基因表达量先进行对数处理。

进一步的，所述表达比例的计算方法如下：

P＝N

其中，P为基因的表达比例；N

进一步的，所述规定区间的确定方法如下：设定基因表达量阈值G1和G2，G1＞G2；所述规定区间为基因表达量大于G2且小于G1。

进一步的，所述基因表达量阈值G1＝800，G2＝80。

进一步的，所述污染阈值GP为70％。

一种计算机存储介质，其上存储有计算机程序，该程序被处理器执行时实现如上所述的单细胞转录组RNA污染鉴定方法。

一种终端设备，包括存储器、处理器以及存储在所述存储器上并可在所述处理器上运行的计算机程序，所述处理器执行所述计算机程序时实现如上所述的单细胞转录组RNA污染鉴定方法。

与现有技术相比，本发明的有益效果是：

本发明通过基因表达量的变化幅度求得两类拐点，并进一步确定背景barcode，进而提取背景barcode的表达谱，统计不同基因在背景barcode中的表达比例，若存在某一基因的表达比例大于污染阈值，则判定单细胞转录组数据存在RNA污染，从而提高后续分析的可靠程度。

附图说明

图1为本发明一实施例的RNA污染鉴定方法流程图。

图2为本发明一实施例的散点图。

图3为本发明一实施例的拐点示意图。

图4为本发明一实施例的RNA污染去除方法流程图。

具体实施方式

下面对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

实施例一：

请参阅图1，本实施例提供一种单细胞转录组RNA污染鉴定方法，包括以下步骤：

计算分析基因表达量(UMI counts数，简称C

S1，根据基因表达量C

S2，求得基因表达量C

具体的，所述A类拐点的求解方法如下：

SA21，以barcode的排名为X轴，基因表达量为Y轴绘制散点图；为了放大基因表达量C

SA22，在所述散点图上，沿X轴，间隔0.1取最近的一个点，间隔长度可以自定义设置，还可以是0.2或0.3。

求相邻两个点之间的斜率k

表1：基因表达量排名及斜率

SA23，当斜率k

由上表1可得知，排名5007的barcode对应的斜率在连续增大的情况下首次大于-1，因此认定为A类拐点。

同理，所述B类拐点的求解方法如下：

SB21，以barcode的排名为X轴，基因表达量为Y轴绘制散点图；

SB22，在所述散点图上，间隔规定距离取最近的一个点，求相邻两个点之间的斜率；

SB23，当斜率的变化趋势为由大到小，且在该趋势持续过程中，斜率首次小于设定的斜率阈值时，将对应的点设为B类拐点。

同理，由上表1可得知，排名1996的barcode对应的斜率在连续降低的情况下首次小于-1，因此认定为B类拐点。

如图3所示，通过步骤S2，可以在散点图上求得若干A类拐点和B类拐点，对拐点的数量进行统计以便后续步骤使用。

S3，设定基因表达量阈值G1和G2，G1＞G2；基因表达量C

S4，提取背景barcode的表达谱，统计不同基因在背景barcode中的表达比例P。表达了某基因(假设为基因A)的barcode数量用N

P＝N

其中，P为基因的表达比例；N

值得一提的是，表达比例P的计算还可采用基因表达量计算。而此步骤中采用barcode数量而不是基因表达量计算表达比例P的原因是，污染本身应具有普遍性，即较多背景barcode中能检测到污染基因的表达，若用基因表达量计算比例P，容易受到个别高表达背景barcode(如细胞碎片barcode)的影响，且与普遍性这一特性存在异质性。总之，使用barcode数量计算比例P的结果更稳定，更符合污染的定义中的普遍性特点，更接近真实情况。

S5，若存在污染基因，则判定单细胞转录组数据存在RNA污染，反之则判定无RNA污染。

实施例二：

本实施例二在实施例一判定单细胞转录组数据存在RNA污染的基础上，还进一步去除细胞中的RNA污染。请参照图4，通过如下步骤实现：

T1，根据污染基因在背景barcode中的表达谱，计算每个污染基因的基因表达量的平均值u1和标准差a1；如下表2所示为一实施例的污染基因A在背景barcode中的表达谱：

表2：污染基因A在背景barcode中的表达谱

可以求得以上5个污染基因表达量以下数据：

平均值u1＝5.8；

标准差a1＝2.99。

T2，输入单细胞转录组数据，首先使用10X Genomics官方的cellranger软件筛选出细胞barcode，再使用Seurat软件对细胞barcode进行分群。针对每一细胞群，分别计算污染基因在该群细胞中基因表达量的平均值u2和标准差a2。

将细胞barcode分群后进行再进行计算的原因为，单细胞转录组数据中的每个barcode可能代表不同的细胞类型，且不同细胞类型的barcode数量不同，若未进行分群处理，污染基因的分布难以形成可计算的模型，去污染的误差较大。进行分群处理后再去除细胞污染，可以最大程度上保留细胞的原始表达量，尽可能准确地去除污染部分。

例如：分群之后，得到一个细胞群包括5个细胞barcode，其表达谱针对上述污染基因A的表达谱如下表3所示：

表3：污染基因A在某细胞群中的表达谱

可以求得以上5个污染基因表达量以下数据：

平均值u2＝41.8；

标准差a2＝13.7。

T3，根据公式Cp＝(C

表4：求得的污染表达量Cp以及对应的实际表达量C

注：以上污染表达量Cp以及对应的实际表达量C

对每个污染基因和每个细胞都进行如上处理后，则获得去除污染后的单细胞转录组数据。

实施例三：

一种计算机存储介质，其上存储有计算机程序，该程序被处理器执行时实现如实施例一所述的单细胞转录组RNA污染鉴定方法和/或实施例二所述的单细胞转录组RNA污染去除方法。

实施例四：

一种终端设备，包括存储器、处理器以及存储在所述存储器上并可在所述处理器上运行的计算机程序，所述处理器执行所述计算机程序时实现如上所述的单细胞转录组RNA污染鉴定方法和/或实施例二所述的单细胞转录组RNA污染去除方法。

对于本领域技术人员而言，显然本发明不限于上述示范性实施例的细节，而且在不背离本发明的精神或基本特征的情况下，能够以其他的具体形式实现本发明。因此，无论从哪一点来看，均应将实施例看作是示范性的，而且是非限制性的，本发明的范围由所附权利要求而不是上述说明限定，因此旨在将落在权利要求的等同要件的含义和范围内的所有变化囊括在本发明内。