10X单细胞（10X空间转录组）基因表达的熵值分析

hello，大家好，今天我们来分享一个新的分析点，基因表达的熵值分析，这个分析主要是用来衡量不同样本同一基因或者通路的熵值变化，也是一个新的角度来看待差异基因。

SCEnt 是一个用于单细胞熵分析的包。 它可以计算细胞群内基因异质性或同质性的指标。 它可以使用这些指标对 scRNA-seq 数据执行特征选择。 文章在Measuring the Information Obtained from a Single-Cell Sequencing Experiment,我们主要来看看代码。

单细胞测序 (sc-Seq) 实验正在产生越来越大的数据集。 但是，大型数据集不一定包含大量信息。 在这里，引入了一个正式框架，用于评估从 sc-Seq 实验中获得的信息量，该框架可用于整个 sc-Seq 分析流程，包括质量控制、特征选择和聚类。 用一些简单的例子来说明这个框架，包括使用它来量化由建议的聚类解释的单细胞测序数据集中的信息量。 我们的信息理论框架提供了一种正式的方法来评估从 sc-Seq 实验中获得的数据的质量，并且对我们理解异质细胞群中基因表达模式的变异性具有广泛的意义 。（专业的词汇总是这个难以理解）

图片.png

• 图注：An information-theoretic view of sc-Seq data. Transcripts, or more generally counts, are assigned to cells after sequencing. If the population is pure, then the information unexplained by the hypothesis of homogeneity is zero (top left). In practice, the assignment process is stochastic, so will naturally result in a small non-zero information unexplained (bottom left). If the population is heterogeneous then transcripts are expressed preferentially in a subset of cells and the information unexplained is large (top right), reaching a maximum (at log2(N), where N is the number of cells) when only one cell expresses the gene (bottom right）

什么是熵值分析，Entropy Analysis。每个基因的表达可以被认为是一个概率分布。 在给定基因已表达的信息的情况下，表达计数可用于推断它在给定细胞中表达的概率。

每个基因的这些概率分布中的每一个都包含一定数量的信息，在数学中称为熵。 熵在历史上被描述为系统中编码的“surprise”的数量。 如果一个给定的基因只在一个细胞中表达，那么这个基因就被认为是低熵的，就好像那个基因被表达了，那么表达它的细胞是已知的，没有“surprise”。 但是，如果给定的基因在每个细胞中的表达均等，则该基因将具有高熵。 由于该基因的表达不会缩小它可能来自的细胞数量，这意味着答案将是一个“surprise”。

这意味着，基因在细胞群中表达得越均匀，表达分布的熵就越大。 因此，允许熵成为细胞群内同质性的度量。

每个概率分布都有一个熵值，有很多方法可以比较概率分布之间的熵值。 其中之一称为 Kullback-Liebler Divergence（或简称 KL Divergence）。 KL 散度测量如果使用一种分布来近似另一种分布将损失的熵量。 (关于KL散度已经分享了很多了，大家可以回看一些文章)。

在这种情况下，可以使用 KL Divergence 来查看如果使用均匀分布来表示基因的表达分布会丢失多少信息。 记住均匀分布将是完全同质基因的表达，KL 散度为我们提供了一个衡量从基因异质性中获得多少额外信息的方法。 本质上给出了细胞群中基因异质性的度量。

我们来看看具体怎么做的

Using SCEnt to Quantify Homogeneity and Heterogeneity

Here is some synthetic scRNA-seq data:

gene_counts
#>       gene1 gene2 gene3 gene4 gene5
#> cell1     0     5     2     3     0
#> cell2     0     5     0     3     0
#> cell3     0     3     2     3     0
#> cell4     0     2     1     3     0
#> cell5     1     0     3     3     5
#> cell6     2     0     0     3     0
#> cell7     3     0     1     3     0

gene_hom() 和gene_het() 函数可分别用于计算基因的同质性或异质性。 每一个都可以传递一个基因，它会返回一个值。

(gene1 <- gene_counts[,1])
#> cell1 cell2 cell3 cell4 cell5 cell6 cell7 
#>     0     0     0     0     1     2     3
gene_hom(gene1)
#> [1] 1.459148
gene_het(gene1)
#> [1] 1.348207

(gene2 <- gene_counts[,2])
#> cell1 cell2 cell3 cell4 cell5 cell6 cell7 
#>     5     5     3     2     0     0     0
gene_hom(gene2)
#> [1] 1.908613
gene_het(gene2)
#> [1] 0.8987422

(gene3 <- gene_counts[,3])
#> cell1 cell2 cell3 cell4 cell5 cell6 cell7 
#>     2     0     2     1     3     0     1
gene_hom(gene3)
#> [1] 2.19716
gene_het(gene3)
#> [1] 0.6101952

(gene4 <- gene_counts[,4])
#> cell1 cell2 cell3 cell4 cell5 cell6 cell7 
#>     3     3     3     3     3     3     3
gene_hom(gene4)
#> [1] 2.807355
gene_het(gene4)
#> [1] 0

(gene5 <- gene_counts[,5])
#> cell1 cell2 cell3 cell4 cell5 cell6 cell7 
#>     0     0     0     0     5     0     0
gene_hom(gene5)
#> [1] 0
gene_het(gene5)
#> [1] 2.807355

注意，gene4 是最大同质基因，gene5 是最大异质基因。 因此，他们的熵计算总是返回 0 或 2.807355 的极值，非零极值将等于 log2[N]，其中 N 是样本中的cells数。

可以传递整个基因表达矩阵，而不是单独提交每个基因。 但是，数据需要采用将基因表示为行而将细胞表示为列的格式。

gene_hom(t(gene_counts))
#>    gene1    gene2    gene3    gene4    gene5 
#> 1.459148 1.908613 2.197160 2.807355 0.000000

除了转置矩阵作为输入之外，为了方便起见，还有一个转置参数

gene_het(gene_counts, transpose = TRUE)
#>     gene1     gene2     gene3     gene4     gene5 
#> 1.3482070 0.8987422 0.6101952 0.0000000 2.8073549

gene_hom() 和gene_het() 函数还有其他参数，大家可以好好看看。 unit参数传递给熵计算以更改计算熵的单位，默认情况下设置为“log2”，以便所有输出都以位为单位。 如果基因表达已经是概率形式，则normalise参数就在那里。 通常不应使用此参数，除非在某些情况下可以加快计算速度。

Using SCEnt For Feature Selection

当想要使用基因表达进行预测时，群体内均匀表达的基因对任何分类都没有用。 相反，异质表达的基因更有可能用于预测。

函数scent_select() 将通过发现基因异质性并应用一些定义的阈值来对scRNA-seq 数据进行特征选择。 可选择的阈值有 bit_threshold、count_threshold 和 perc_threshold。 bit_threshold 采用一个数值，SCEnt 只会选择异质性值大于给定位值的基因。 count_threshold 取一个整数值，SCEnt 将返回只有前 n 个异质表达基因的数据，其中 n 是 count_threshold。 perc_threshold 的值介于 0 和 1 之间，SCEnt 将仅返回异质性大于由 perc_threshold 给出的种群异质性百分位数的基因。 与上面的熵计算不同的是，基因表达的矩阵应该以细胞为行，基因为列的格式。 同样，有一个内置的转置选项，以及要传递给熵计算的单位和归一化参数。

gene_counts
#>       gene1 gene2 gene3 gene4 gene5
#> cell1     0     5     2     3     0
#> cell2     0     5     0     3     0
#> cell3     0     3     2     3     0
#> cell4     0     2     1     3     0
#> cell5     1     0     3     3     5
#> cell6     2     0     0     3     0
#> cell7     3     0     1     3     0

根据提供的阈值，scent_select() 可以执行这三种形式的特征选择中的每一种。

scent_select(gene_counts, bit_threshold = 0.85)
#>       gene1 gene2 gene5
#> cell1     0     5     0
#> cell2     0     5     0
#> cell3     0     3     0
#> cell4     0     2     0
#> cell5     1     0     5
#> cell6     2     0     0
#> cell7     3     0     0
scent_select(gene_counts, count_threshold = 2)
#>       gene1 gene5
#> cell1     0     0
#> cell2     0     0
#> cell3     0     0
#> cell4     0     0
#> cell5     1     5
#> cell6     2     0
#> cell7     3     0
scent_select(gene_counts, perc_threshold = 0.25)
#>       gene1 gene2 gene3 gene5
#> cell1     0     5     2     0
#> cell2     0     5     0     0
#> cell3     0     3     2     0
#> cell4     0     2     1     0
#> cell5     1     0     3     5
#> cell6     2     0     0     0
#> cell7     3     0     1     0

尝试对多个值设置阈值会引发错误

scent_select(gene_counts, bit_threshold = 0.85, count_threshold = 2)
#> Error in scent_select(gene_counts, bit_threshold = 0.85, count_threshold = 2): 
#>  Only one threshold can be set at a time

这是为了避免操作顺序出现问题，而应该将scent_select() 传递回自身，以便应用阈值的顺序是明确的

gene_counts %>%
  scent_select(bit_threshold = 0.85) %>%
  scent_select(count_threshold = 2)
#>       gene1 gene5
#> cell1     0     0
#> cell2     0     0
#> cell3     0     0
#> cell4     0     0
#> cell5     1     5
#> cell6     2     0
#> cell7     3     0

其实熵值分析主要是来寻找变化大的基因，也就是同质化和异质化的基因，对我们理解样本的特征有很深的借鉴意义。

生活很好，有你更好