今天我们来将一个数据分布的情况，接触10X单细胞或者10X空间转录组不久之后，就认识到了一些数据的特征，比如经常提到的dropout/missing data/zero-inflated，当时深以为然，也没有过多思考，直觉上感觉是对的，但是随着分析的深入，很多分析点和数据分布的特点格格不入，所以最近看了几篇文献，纠正自己看法的同时，也给大家一些警示。

首先第一篇文献，

v2-65b2c9b8c6686c44c23081f0e0548330_720w.jpg

这篇文献主要讲了以下几个事情，

• 如果single cell UMI data 先分好细胞类型之后（注意文章中用的数据都是UMI data），数据中的 zero-proportion 是和普通的泊松分布吻合的非常好的。不需要Zero-inflation，这个其实说明了zero并不特殊，就是简单的随机过程的结果。

图片.png

原文Figure1 C： 左边是Umap可视化的图，可以看到有三种细胞类型。右边三幅横坐标是基因表达值，纵坐标是该基因表达值为 0 的细胞比例。实线是泊松分布的理论估计值，每个点都是一个 gene，黑色代表三个细胞类型混在一起的结果，右边是分开之后的结果。可以看到在分完细胞类型之后这个zero proportion吻合的相当好。

• 文章另外说明了如果先做 pre-processing 会引入一些噪声

  
  
  图片.png

原文 Figure 2A,B ：左边是raw data，左三是SF normalized 之后，反而有了可区分的noise。且从右图可以看到其实是引入了测序深度这个noise。

• 所以在上一点的基础上利用 zero-proportion 直接来做了一个迭代分 cluster 的方法，而不是先做pre-processing。每次迭代先选取最体现细胞异质性（也就是上图中最偏离期望值的基因），根据这些gene来划分细胞（也就是拟合混合泊松），然后迭代。比较有意思的事情是传统的方法挑选出来的 gene mean 普遍比较小（毕竟是用CV挑的，mean小自然就大了），而且因为mean比较小其实variance也比较小。

  
  
  图片.png

原文Figure 3C: 可以看到和传统的比，作者认为挑选出来的gene mean的分布感觉是比较合理的，因为那些gene mean小的不那么准确。

很简单的统计模型得到有意思的结论！

第二篇文献

图片.png

这篇文章个人觉得视作是一个 review 会更合适，没有提出新的工具和方法，但是提供了一个（不那么新的）思路，将对表达值的建模分为真实表达模型和measurement模型，整理总结了在这个角度下对许多之前的工作的理解。 其中第二段对 “dropout”, “missing data”, “imputation”, “zero-inflation” 的阐述很好，尤其是叫大家别用 dropout 这个概念

背景

基于二代的测序方法概括的说是想办法把生物系统内的分子（DNA/RNA）扩增了，然后通过测序得到对这些分子在生物系统内含量的定量化（其中值得一提的是使用 UMI 可以大大减少扩增过程带来的影响），因此自然的得到的定量化就是关于真实含量的一个函数。

概念澄清

首先作者号召大家不要在使用 dropout 这个词了，因为其含义过于混乱

最早 dropout 是使用在 PCR 时对某些等位基因无法成功扩增从而导致数据中缺失或者错误。而在单细胞测序中，通常被认为是某些 gene 的表达因为技术原因没有在某些细胞中测到。但是这个概念在很多地方有不同的意思（如作者在文中所述，颇为混乱）

Since then, use of the term “dropout” has spread rapidly, but its meaning varies among papers and presentations (and sometimes even within papers and presentations!).

Dropout 有时候指 zero value，有时候指表达但是没测到的，有时候指对所有表达 值的影响（不只是0，因为其他值都被低估了），所以作者认为不用 Dropout 这个词，你好我好大家好。

missing data 和传统数据科学中这个概念的矛盾

有些文章（特别是 imputation 的一些文章）中也会将 0 称为 missing data，这显然是和统计中的 missing data 有所矛盾，一般统计中的 missing data 是指没有observed 的数据，比如不知道样品来自的性别，那这个性别变量就是 missing data。

作者举了一个小例子，比如你需要在马路上统计通过的车辆数量，即使没有看到任何车辆，这个 0 仍然是有意义的。

zero-inflated 也许是没有必要的。

即使 UMI data, zero proportion 很高，但是是在泊松分布或者负二项分布之内的。

Modeling scRNA-seq data

本文关键的核心就是将 measurement 和 actual expression 分开（当然事实上很多文章就是这样做的）。也就是下式，认为真实表达服从一个分布，然后基于这个表达产生了一个measurement的分布。而且作者通过一系列的empirical 实验说明 measurement这部分基本可以被泊松分布给建模。

图片.png

首先作者对之前的文章总结，之前的许多文章的假设如下表所示（其实大部分文章就是这样写的，本来就是分开的）。

v2-fe669d41ac8848c9072f221cc958c312_720w.jpg

文章还对 multiple genes 的一些模型做了类似的总结，这边略去细节，感兴趣的可以看原文。其实思路还是一样的，只不过会有个 link function （因为一般model时，都会将 expression matrix 投到低维上去，所以就需要这样一个 link function 来做高维与低维之间的映射，这部分文章很多，以 VAE 自编码器来理解比较容易，也许之后可以写一写）

v2-b78980b587bd257a2fb96259b1f7ae07_720w.jpg

结果部分大家可以读一读，这里就不过多分享了。

然后是第三篇文献

图片.png

这篇文章其实就说了一件事，droplet scRNA-seq 方法得到的数据中 0 的比例可以很好的被负二项（Negative Binomial）分布拟合，没有 zero-inflated。

结论

其实非常简单，首先根据数据拟合负二项分布，也就是得到负二项分布的参数（泊松过程同理）。这边使用了 CEL-seq 的数据，在 notebook 中可以看到其他方法的结果，为什么这样算，可以稍微了解一下负二项分布，都是些基本的性质。

Load CEL-seq data

# Calculate the mean
mean_array = np.array(data.mean(axis=0)).ravel()

# Calculate the varience 
var_array = np.array(data.var(axis=0)).ravel()

# Calculate the dispersion parameter
r = np.square(mean_array)/((var_array-mean_array))

然后根据参数就可以计算出负二项分布的理论 0 的比例

# Calculate the zero proportion 
zero_propotion = ((data == 0).sum(0))/data.shape[0]

# Calculate the zero proportion predicted by Negative Binomial distribution 
NB_predict = np.power(r/(r+mean_array),r)

# Calculate the zero proportion predicted by Poisson distribution 
Poisson_predict = np.exp(-mean_array)

把理论值和真实数据中的结果作比较。

plt.figure(figsize=(12,5))
plt.subplot(1, 2, 1)
plt.scatter(np.log(mean_array), zero_propotion, s=10, label ='Observed')
plt.scatter(np.log(mean_array), Poisson_predict, s=10, label ='Poisson distribution')
plt.xlabel('Log(mean expression)')
plt.ylabel('Zero proportion')
plt.title('Poisson Distribution')
plt.legend()
plt.subplot(1, 2, 2)
plt.scatter(np.log(mean_array), zero_propotion, s=10, label ='Observed')
plt.scatter(np.log(mean_array), NB_predict, s=10, label ='Negative Binomial distribution')
plt.xlabel('Log(mean expression)')
plt.ylabel('Zero proportion')
plt.title('Negative Binomial Distribution')
plt.legend()
plt.show()

v2-cf79ca51bdadfd3c0de17f9fd7d30827_720w.jpg

再看一下散点图就会发现负二项吻合的极好

v2-6575b54962fa28c406af92c152316f54_720w.jpg

一些思考

这个负二项拟合的很好可以告诉我们什么呢？

还记得前文我们提到这个表达值可以拆成两部分，一部分是 expression model，一部分是 measurement model， measurement model 部分是一个 泊松我觉得从测序的原理上来说非常有道理，那么如果最后 observation model 是个 负二项说明什么呢？回忆一下之前的那张表。

v2-fe669d41ac8848c9072f221cc958c312_720w.jpg

可以看到其实可以说明大概 expression model 是个 Gamma 分布，这也是很有道理的。

另外前面分享也就是 "Demystifying" drop-outs" in single cell UMI data." bioRxiv (2020) 中说到，细胞的异质性也会影响 zero proportion，分完细胞类型之后泊松就拟合的很好了。

感觉读了以后，对自己之前的认识有了很多的改变，看来学习确实是终身的事情

生活很好，有你更好