scRNA-seq能够在单个细胞分辨率下研究细胞异质性对扰动的响应。然而，由于技术限制，扩大高通量筛选（HTSs，highthroughput screens）来测量许多药物的细胞反应仍然是一个挑战。因此，目前依然需要借助常规的bulk-RNA传递信息。chemCPA是一种新的编码器-解码器结构，用于研究未知药物的扰动效应。作者将该模型与迁移学习的架构相结合，演示了如何在现有的bulkRNA-HTS数据集上进行训练可以提高泛化性能。更好的泛化减少了对单细胞分辨率的昂贵需求，最终加速药物发现。

来自：Predicting Cellular Responses to Novel Drug Perturbations at a Single-Cell Resolution

背景概述

单细胞测序允许同时分析数百万个细胞，增加了探索细胞异质性的分辨率。利用scRNA-seq和高通量筛选(HTSs)，人们现在可以在细胞分辨率上研究不同扰动（即药物剂量组合，drug-dosage combinations）对转录组的影响。与传统的HTs不同，scRNAseq HTs可以识别基因表达和细胞异质性的细微变化，这是制药和药物发现的基石。

因此，需要计算方法来解决现有实验方法的有限探索能力，并发现有希望的候选药物。合适的方法应该预测对未知组合扰动的反应。就医学影响而言，对未知扰动的预测可能是最理想的。同时，它要求模型在多个不同的细胞环境中正确捕获复杂的化学相互作用。这种泛化能力还不能仅从单细胞HTs中学习，因为它们被认为不能覆盖所需的化学相互作用的广度。在这项工作中，作者利用跨数据集的信息来缓解这个问题。

相关工作

在过去的几年里，深度学习(DL)已经成为分析和解释scRNA-seq数据的重要工具。特别是表征学习，不仅可用于识别细胞异质性和整合，还能将query映射到reference，而且还可用于模拟单细胞扰动响应。

与线性模型不同，DL适用于捕获非线性细胞类型特异性反应，并且易于扩展到全基因组测量。最近，Lotfollahi等人引入了CPA方法对scRNA-seq数据进行扰动建模。CPA不能推广到看不见的化合物，这阻碍了其应用于尚未通过scRNA-seq数据测量的药物的虚拟筛选，而这是有效药物发现所必需的。

另一方面，对于bulkRNA数据，已经提出了几种方法来预测de novo化学物质的基因表达谱。至关重要的是，LINCS项目引入的L1000数据集极大地促进了基于表型的化合物筛选的进展。然而，目前尚不清楚如何将这些方法转化到具有少量化合物信息的单细胞数据集，并且可用于不同的基因集合。

化学组合扰动自编码器

考虑一个数据集 D = { ( x i , y i ) } i = 1 N = { ( x i , ( d i , s i , c i ) ) } i = 1 N D=\left\{(x_{i},y_{i})\right\}_{i=1}^{N}=\left\{(x_{i},(d_{i},s_{i},c_{i}))\right\}_{i=1}^{N} D={(xi​,yi​)}i=1N​={(xi​,(di​,si​,ci​))}i=1N​其中$x_i\in R^n$描述了$n$维基因表达，$y_i$是属性集。对于scRNA-seq扰动数据，我们通常考虑药物和剂量属性（drug and dosage attributes），$d_i\in${drugs in $D$}，$s_i\in R$。细胞$i$的细胞系记为$c_i$。注意，这组属性依赖于可用数据，并且可以扩展为协变量，如患者或物种。

预测反事实组合的一种可能方法是将细胞的基因表达$x_i$从其属性$y_i$不变地编码为latent $z_i$，称为basal state。$z_i$可以与属性$z_{d_i}$和$z_{c_i}$组合，编码任意属性组合$y_i^{\prime }\ne y_i$，并解码回对应这组新选择属性的基因表达状态${\hat{x}}_i$。

为此，作者将化学组合扰动自编码器chemCPA分为三个部分：1.基因表达编码器和解码器，2.属性嵌入器，3.对抗性分类器，参见图1的说明。

• 图1：chemCPA架构，模型包含3部分：1.基因表达编码器和解码器，2.属性嵌入器，3.对抗性分类器。分子编码器$G$可以是任意的graph-based和language-based模型，前提是生成固定长度的embedding $h_{drugs}$。MLPs $S$和$M$被训练用于映射embedding到扰动latent space。这里，$z_{d_i}$被添加到basal state和协变量embedding $z_{c_i}$中。在这项工作中，后者对应于细胞系。basal state $z_i=E_{\theta }(x_i)$通过对抗性分类器$A_{\varphi }^j$训练为不变量，解码器$D_{\psi }$产生高斯似然$N(x_i|{\mu }_i,{\sigma }_i)$。

基因编码解码器

模型基于结合对抗训练的编码解码器架构。编码器$E_{\theta }:R^n\to R^l$是一个MLP将测量的基因表达$x_i\in R^n$映射到$z_i\in R^l$。通过对抗分类器，$z_i$被训练成不包含任何关于属性$y_i$的信息。这让我们可以控制latent space，我们在latent space中用选择的附加属性embedding更新$z_i$并获得$z_i^{\prime }$。

解码器$D_{\psi }:R^l\to R^{2n}$是MLP，取$z_i^{\prime }$作为输入并计算基因表达的分布$P$，如果输入的$x_i$是原始数据则$P$遵循NB分布，如果输入的$x_i$是预处理数据则$P$遵循高斯分布。假设是高斯分布，则得到参数化输出$\mu =D_{\psi }^{\mu }(z^{\prime }),{\sigma }^2=D_{\psi }^{{\sigma }^2}(z^{\prime })$。虽然chemCPA支持这两种设置，但作者在实验中观察到高斯分布具有更好的收敛性，因此重建损失为：$$\begin{gathered} L_{rec}(\theta ,\psi )=N(x_i|{\mu }_i,{\sigma }_i)=\frac{1}{2}[In(D_{\psi }^{{\sigma }^2}(z_i^{\prime }))+\frac{(D_{\psi }^{\mu }(z_i^{\prime })-x_i{)}^2}{D_{\psi }^{{\sigma }^2}(z_i^{\prime })}] \\ {z}_i^{\prime }=E_{\theta }(x_i)+z_{attribute} \end{gathered}$$接下来，作者提供了关于如何有意义地解释潜在空间的方法以及如何编码药物和细胞系属性的方法。

属性embedding和加性latent space

假设潜在空间中的扰动响应具有加性结构：$$z_i^{\prime }=z_i+z_{c_i}+{\hat{s}}_i{z}_{d_i}$$其中，$z_{c_i}$和$z_{d_i}$对应细胞系和药物属性，${\hat{s}}_i$为剂量的编码。

这种线性的选择使模型对生物学家等用户来说是可解释的，因为它允许单独分析和去除组分。加性结构的另一个优点是它的排列不变性，它允许添加新的协变量。

由于药物和细胞系属性的不同，作者在潜空间中分别编码药物和细胞系属性。对于细胞系，使用与Lotfollahi等人(2021)相同的方法，其中为每个细胞系$c$编码了$l$维潜在表示$z_c$。对于药物，作者提出了一个新的嵌入网络$P_{\phi }$。

扰动网络

网络$P_{\phi }$映射分子的表示，由$G,M,S$组成，见图1(2)。分子编码器$G$将分子映射到固定长度的向量$h_{d_i}\in R^m$。后续步骤中，扰动编码器$M:R^m\to R^l$以分子向量作为输入，生成药物扰动表示$z_{d_i}$。剂量$S:R^{m+1}\to R$使用$h_{d_i}$，并将其与剂量$s_i$一起映射到缩放后的剂量表示${\hat{s}}_i$。作者选择$S$映射回标量值$\hat{s}$，因为这使我们能够以一种简单的方式计算药物反应曲线。此外，这种编码方式与$z_{d_i}$编码药物的一般效应（与剂量无关）的场景相匹配。合在一起，得到：$${\hat{s}}_i\times z_{d_i}=P_{\phi }(g_i,s_i)=S(h_{d_i},s_i)\times M(h_{d_i})$$分子编码器$G$可以是任何将分子表示映射到固定大小嵌入的编码网络。由于scRNA-seq HTs中可用的药物数量有限，作者建议依靠预训练的编码模型并在训练期间冻结$G$。同时，作者发现RDKit特征的表现良好，并报告了使用RDKit作为分子编码器$G$的chemCPA的结果。

通过该设计，我们可以为药物和细胞系选择一组新的属性，并计算新的latent状态为：$z_i^{\prime }=z_i+z_{attribute}$。基于扰动网络，chemCPA可以预测实验未观察到的分子的药物扰动。

学习不变表示的目的：数据集中的基因表示本质都是受到协变量控制的，解耦后才能和协变量结合产生有意义的预测结果。因此，接下来，将描述如何“从属性信息中剥离$z_i$”以获得basal状态表示。

对抗分类器用于学习不变的basal状态

生成basal也被称为解开$z_i,z_{d_i},z_{c_i}$的纠缠，使用对抗分类器$A_{\varphi }^{drug}$和$A_{\varphi }^{cov}$。对抗网络$A_{\varphi }^j:R^l\to R^{N_j}$以$z_i$为输入，目标是预测细胞$i$应用的药物及其细胞系$c_i$。在训练这些分类器以提高分类性能的同时，作者还在编码器$E_{\theta }$的训练目标中添加了reversed sign分类损失。因此，编码器试图产生一个不包含有关属性信息的潜在表示$z_i$。注意，这种basal、药物和协变量信息的明确分离，称之为解纠缠，是使问题易于处理的近似。同时，这种分离有助于将扰动效应归因于特定来源，例如药物或细胞系，这与生物学应用和下游分析有关。

使用交叉熵损失：$$\begin{gathered} L_{class}^{drugs}=CE(A_{\varphi }^{drug}(z_i),d_i) \\ {L}_{class}^{cov}=CE(A_{\varphi }^{cov}(z_i),c_i) \end{gathered}$$作者在对抗性分类器的损失函数中添加了一个零中心梯度惩罚，以最小化：$$L_{pen}^j=\frac{1}{k}\sum _k||{\partial }_{z_i}{A}_{\varphi }^j(z_i{)}_k|{|}_2^2$$

当应用于生成对抗网络时，这种梯度惩罚被证明使判别器对噪声更具鲁棒性，并实现局部收敛。在训练期间，作者在以下相互竞争的目标之间交替更新步骤：$$\begin{gathered} L_{AE}(\theta ,\psi ,\phi |\varphi )=L_{rec}(\theta ,\psi ,\phi )-{\lambda }_{dis}\sum _j{L}_{class}^j(\theta |\varphi ) \\ {L}_{Adv}(\varphi ,\theta )=\sum _j{L}_{class}^j(\varphi |\theta )+{\lambda }_{pen}{L}_{pen}^j(\varphi ) \end{gathered}$$其中，${\lambda }_{dis}$平衡重建和编码器的不变表示，以产生解纠缠的$z_i$。

数据集和迁移学习

作者分别使用sci-Plex3和L1000数据集进行单细胞数据的主要评估和bulk实验的预训练。

数据集：L1000数据包含了大约130万份针对978种不同基因的bulkRNA的观察结果。它包括对近2万种不同药物的测量，其中一些是FDA批准的，而另一些是合成化合物，对任何疾病都没有证实的效果。与scRNA-seq数据相比，L1000允许探索更多样化的分子空间，这使其成为预训练的理想选择。

sci-Plex3包含649,340个细胞的7561个药物敏感基因的测量值。在三种人类癌细胞系（A549、MCF7和K562）上，研究了188种药物在4种不同剂量（10 nM、100 nM、1μM和10μM)下的单化合物扰动。请注意，所有细胞系和大约150种化合物与L1000数据重叠。此外，Srivatsan等人(2020)为所有化合物分配了19种不同的作用模式，也称为pathway，描述了细胞暴露于药物分子所导致的变化。

迁移学习：当我们用高斯似然损失训练chemCPA时，首先对数据集进行归一化，然后log变换。根据实验，进一步减少了单细胞数据中包含的基因数量。对于扩展了基因集的情况（因为通常考虑2000个基因），可以在$E_{\theta }$前加一个非线性层将微调数据集的2000基因压缩到预训练上的977基因。在微调过程中训练所有的层，包括新添加的层。