通过高维细胞状态空间确定单个细胞轨迹对于理解从药物分化到病变细胞的细胞分化到表观遗传反应等生物学变化而言是一项艰巨的挑战。我们结合实验和理论来确定单个BRAF V600E的轨迹突变型黑色素瘤癌细胞处于药物纯净状态和药物耐受状态之间。尽管单细胞组学工具可以生成细胞状态概况的快照，但是通过该随机空间状态切换可能会使通过该空间的各个细胞轨迹的确定产生混淆。我们在药物治疗5天的各个时间点检测了一组信号传导，表型和代谢调节剂，以揭示细胞状态态势，并通过两条途径连接了纯药物性和耐药性状态。给定细胞采取的轨迹取决于谱系限制的转录因子的药物纯水平。每条轨迹都表现出独特的药物敏感性，从而更新了等基因细胞群体中适应性耐药性发展的范式。
从癌细胞中药物耐受状态的发展到干细胞分化的细胞过程可以描述为细胞状态变化。具体而言，某些癌细胞最初响应于靶向抑制剂，对这些致癌驱动动作1可经由非遗传机制演变成一个药物耐受性状态，或许耐药克隆的出现前述2，3，4，5。两种状态之间的癌细胞转移怎么能告知使用额外的药物分子设计的细节逮捕过渡6，7，8。以前的研究已经发现在信令，代谢，转录的耐药性机制的见解，和表观遗传水平5，9。然而，这些研究中的大多数要么仅在没有单细胞分辨率的情况下在整体水平上比较了药物耐受细胞和药物敏感性细胞，要么没有提供连接这两种状态的轨迹的详细时间分辨特征。我们假设在药物敏感性和药物耐受性状态之间，细胞可以有多个独立的通路。如果这是真的，那么寻找可以阻止不良细胞状态转变的药物组合的挑战就大大增加了。在这里，我们研究了一种高度可塑性的癌细胞系，该细胞系在用靶向抑制剂治疗后，在几天之内便从快速分裂的药物反应状态转变为药物耐受的慢循环状态。我们表明，单元格确实可以在两个状态之间采取多种类别的轨迹。
从功能的角度，单元状态变化通常在基因表达的伴随变化7，10，11，12，13，蛋白信号9，10，12，14，15，16，17，18，19，和细胞代谢20，21，22，23。高度多路复用单细胞的方法24，25，26，27可以绘制出与细胞状态改变相关联的单元状态景观提供强大的工具17，28，29，30，31。然而，即使对于同基因细胞系，捕获单个细胞在穿越这些景观时所走的轨迹也具有挑战性。这是因为多路复用单细胞组学方法仅在给定时刻提供占用的细胞状态空间的快照。在连续的时间点捕获的细胞之间的相似性测量可以通过风景暗示可能路径32，33，34，35。但是，单元格可能会从一个状态随机切换到另一种状态，因此单个单元格可能无法在状态之间采取平滑的轨迹。时间推移成像方法可以映射单个细胞的轨迹，但仅在两到三个分析物为每个小区等提供小区状态空间的有限视图36，37，38。因此，从动态的一系列细胞状态空间快照中提取细胞轨迹的能力将具有很高的价值。在这里，我们针对解决这一基本挑战的实验和理论方法相结合的报告。
我们利用源自患者的BRAF V600E突变型黑色素瘤癌细胞系39作为快速开发针对靶向抑制剂的药物耐受性的模型。下BRAF抑制，这些高度塑料细胞从药物响应性状态的药物耐受性状态迅速地转变10，16。我们使用集成的单细胞功能蛋白质组学和代谢测定法来表征这种转变，该测定法旨在广泛地采样与所选癌症特征和特定细胞状态相关的蛋白质和代谢产物。降维，信息理论分析和时间序列单细胞数据的可视化揭示了复杂的细胞状态空间格局，并暗示了在纯药物状态和药物耐受状态之间可能存在两条不同路径的可能性。进一步的实验测试这些路径是否构成独立的细胞轨迹。实际上，我们发现即使是同基因的肿瘤细胞也可以采取不同的，独立的轨迹来耐受药物。这两个轨迹与不同的信号传导和代谢网络相关，并且可以独立进行药物治疗。
我们通过测定所选蛋白质的面板，其特征个体黑色素瘤细胞药物适应，再加葡萄糖摄取，在BRAF V600E使用单细胞条码芯片（SCBC）期间BRAFi治疗的第5天突变M397细胞培养物10，17，26，40，41，42，43（图 1A）。在药物治疗的0、1、3和5天（D0对照，D1，D3和D5）之后，将单个细胞分离到SCBC内的纳升体积微腔中。将每个分离的细胞原位裂解以释放其细胞内含物。一个SCBC内的每个微型腔包含一个完整的条形码阵列，其中每个条形码元素是为特定的蛋白质捕获的抗体44或分子探针通过竞争测定设计用来测定特定代谢物42，43（图1A）。该板的设计是由BRAFi处理M397细胞的转录组学分析（补充图告知 1）和现有文献9，10，12，20，45。该小组广泛采样了癌症的各种功能和代谢标志以及细胞状态标志。
a单细胞整合蛋白质组和代谢分析实验设计。使用基于微流体的单细胞条形码（SCBC）技术，收集并处理BRAFi处理过程中来自不同时间点的细胞。每个细胞的特征在于六种不同类别的标记物的水平。b综合蛋白质组学和代谢分析数据集的热图表示。每行代表一个单独的单元格，每列（最后一列除外）代表一个单独的分析物，热图中的颜色代表所测量的分析物含量。最后一栏代表开始BRAFi治疗后的天数。在X轴上，标记分别对应于六个功能类别中的哪个类别进行着色。C小提琴图表示某些代表性标记在四个时间点上的分布。Y轴代表测得的标记物水平的自然对数。每个图都以所测量标记的功能类别的颜色为边界。
BRAFi天真（D0）M397细胞的单细胞分析揭示了基线时许多测定标记物的异质水平。参照图 1b，c和补充图 2，某些分析物在整个细胞群体中表现出高变异性。这些包括黑素细胞谱系转录因子MITF及其下游的黑素细胞状态标记MART1，代谢调节剂HIF1α和p-AMPKα以及增殖标记Ki67。Ki67的差异意味着该群体既包含快速循环细胞又包含缓慢循环细胞。相比之下，高葡萄糖摄取和乳酸脱氢酶（LDH）和PKM2代谢酶的表达在各个细胞之间相对一致。最初（在D1处）的药物治疗会抑制葡萄糖的摄取，并抑制大多数代谢调节剂和信号磷蛋白以及Ki67。这些癌症标志的抑制反映了最初的BRAF抑制后关键致癌信号通路的阻断。该药物还促进瞬时细胞分化，然后再去分化，这由D3中MART1表达的增加以及D5中的下调所证明。然而，在D1中，一小部分M397细胞亚群仍为Ki67-高，这意味着该细胞亚群中的药物反应较慢。在D3时，大多数分析物的方差急剧而短暂地增加，并随D5的收缩而缩小。此变化包括除p-LKB以外的所有代谢调节剂，除Slug以外的所有抗性状态标志物和调节剂，所有代谢酶和所有信号转导磷蛋白。根据先前的报告，D3处许多标记的波动幅度增加 在D1中，一小部分M397细胞亚群仍为Ki67高，这意味着该细胞亚群中的药物反应较慢。在D3时，大多数分析物的方差急剧而短暂地增加，并随D5的收缩而缩小。此变化包括除p-LKB以外的所有代谢调节剂，除Slug以外的所有抗性状态标志物和调节剂，所有代谢酶和所有信号转导磷蛋白。根据先前的报告，D3处许多标记的波动幅度增加 在D1中，一小部分M397细胞亚群仍为Ki67高，这意味着该细胞亚群中的药物反应较慢。在D3时，大多数分析物的方差急剧而短暂地增加，并随D5的收缩而缩小。此变化包括除p-LKB以外的所有代谢调节剂，除Slug以外的所有抗性状态标志物和调节剂，所有代谢酶和所有信号转导磷蛋白。根据先前的报告，D3处许多标记的波动幅度增加 和所有的信号磷蛋白。根据先前的报告，D3处许多标记的波动幅度增加 和所有的信号磷蛋白。根据先前的报告，D3处许多标记的波动幅度增加41，46，意味着近该时间点的一个或多个小区的状态变化。流式细胞仪分析也证实了这一点（补充图 3）。到D5，葡萄糖摄取增加回到D0水平，但方差增加。Ki67进一步降低，并且相对于D0而言，方差急剧减小。实际上，D5中的大多数细胞进入衰老状态，而凋亡细胞死亡的发生率却没有增加（补充图 4和 5）。）。另外，在第5天，上皮-间质转化相关转录因子Slug的变异和丰度增加，表明一些趋向于间质表型的细胞的出现。此外，与药物抗性相关的其他测定的蛋白质标记（AXL，N-钙黏着蛋白，NGFR和TNFR）的水平到D5时都更高。这些标记的变化也通过流式细胞仪分析得到了证实（补充图 6）。）。葡萄糖摄取的上调和许多抗性标志物表明细胞已通过D5启动了耐药性程序。因此，在单个分析物水平上观察时，单细胞整合蛋白质组学和代谢分析提供了在D1时初始药物反应，D3时药物诱导的细胞状态变化以及D5时出现药物耐受性的证据。细胞增殖（完全耐药）的增加已经显示在几周后发生。这些意见都与现有文献一致9，12。
同时可视化跨多个标记的时间相关，协调的变化需要算法，该算法可以降低数据集的高维性。我们应用了三种算法：FLOW-MAP 47，t-SNE 48和PHATE 49。所有的方法所提供的数据集（图的直观表示 2和补充图 7 - 9）。FLOW-MAP分析表明，黑色素瘤细胞主要根据药物暴露时间聚集（图 2）。，左上图），表示按时间顺序排列的细胞状态轨迹。大多数未经处理的M397细胞（在图表的左下方）的特征是，除了N-钙粘着蛋白，MITF，HIF1α，Ki67和MART1以外，所有被测分析物的水平均一（参见图2和补充图7的特定于分析物的 图）。 ）。这些非均匀表达的蛋白质大多数在整个D0细胞簇中表现出从左到右逐渐变化的差异，未处理细胞的一小部分亚群（D0簇的右侧）表现出Ki67，MITF和MART1的较低表达。这些特征指向一小组去分化的慢循环细胞。经过BRAFi处理，细胞群体最初分裂为占据FLOW-MAP的两个区域。在D1（绿色点）处，大多数细胞聚集在D0细胞的右上方，而一小部分亚群直接聚集在D0组的右侧。这种趋势在D3处继续，大多数细胞聚集在最大D1质量上方，而少数细胞聚集在小型D1组的右侧。通过D5（紫色），所有单元都聚集在图形的右侧。在第1天和第3天，细胞分叉意味着上下轨迹朝着药物耐受状态的可能性。t-SNE还表明了两类轨迹的可能性48点PHATE 49分析（补充图 8， 9）。因此，对单细胞数据集的计算分析表明，在BRAFi适应的早期阶段出现了分叉的药物反应。
每个点代表一个单独的单元格。每对单元之间的距离代表它们之间的总体多组相似性。距离足够近的单元对之间有一条边链接。主面板中的点的颜色（左上）代表相应单元格的BRAFi曝光时间（0、1、3或5天）。其他面板中的圆点颜色表示每个单元格中每个标记的丰度。MITF，MART1和Ki67水平面板中的虚线框显示了第0天细胞的一小部分亚群，这些亚细胞缓慢循环且黑素细胞表型较少。
为了进一步解剖与分叉药物响应轨迹相关的分子的变化的动力学，我们应用surprisal分析50，51，52，以我们的单细胞数据集。Surprisal分析是已被广泛地应用于理解大规模批量和单细胞组学数据集的热力学启发法46，50，52，53，54。这种方法是基于对系统的稳定状态的识别（正式来说最低自由能的状态），并且，从这个理论最小增加的自由能的任何约束（分析物模块）52，55。使用这种方法，我们确定了两个主要模块，每个模块代表一组显示跨细胞协调变化的分析物。基于这两个模块的所有20种分析物的预测表达与测得的单细胞数据集非常匹配（补充图 10和11），表明模块1和2概括了5个单元中所有细胞的所有分子标记的总体变化。天的药物治疗过程。
单元中模块的影响力得分（参考值52中定义的lambda值）表示与模块相关的分析物在该单元中富集或抑制的程度。通过将模块1和2的影响力得分颜色编码到FLOW-MAP图的每个节点上，可以使模块1和2可视化（图 3a）。我们发现，模块1的影响力分数沿上下两个路径从负（蓝色）逐渐增加到正（红色），中间时间点的符号变化明显（lambda1 = 0）（图 3a），左侧面板），表明沿过渡轨迹存在生物物理屏障。先前我们已经表明，这种符号变化可以暗示单元状态的转换以及不同单元状态之间的边界50。考虑到Ki67表达的负相关和NGFR / AXL表达与模块1得分的负相关（补充图 12），模块1得分变化的时间依赖性似乎反映了从药物反应性状态向慢反应状态的转变。在第1天和第3天之间处于循环，药物耐受状态。同样，当将模块2得分投影在FLOW-MAP上时，也显示出符号变化（lambda2 = 0），这表明存在一个生物物理屏障将上层和上层分开。较低的路径（图 3a，右图）。值得注意的是，黑素细胞表型转录因子MITF的表达及其下游蛋白MART1均与模块2得分呈负相关（图 3b和补充图 13）。），表明这两个路径的分离可能与细胞的黑素细胞谱系有关。通过额外的z分数归一化或删除前两个变化最大的标记Ki67和MART1（补充图14）， 也获得了相似的结果。总而言之，意外分析可以解析时间相关模块和路径特定模块。它还显示，随着细胞从单纯药物发展到耐受药物，它们占据了相当复杂的环境：包含由两个生物物理屏障隔开的四个不同的细胞状态（补充图 15）。
一个影响力的可视化得分从surprisal分析鉴定的两个调节模块。模块1是与时间相关的，而模块2具有特定于路径的模式。黑色虚线表示每个单元的相应模块得分接近零的区域。b各个标记水平与模块2得分之间的皮尔逊相关性。c，d模块2的Ki67和MITF表达水平的箱线图在第0天得分高和-低亚群。数据是中位数，标有第一和第三四分位数（框），以及顶部和底部四分位数（须）。每个实验是每组n = 16个生物学独立细胞的结果。Ë在第0天，通过qPCR在分选的MITF-High和MITF-Low细胞中测量Ki67相对表达。每个实验是每组n = 3个生物学独立样品的结果。f在未经处理的情况下，在第0天从分选的MITF-High和MITF-Low细胞中收集的时间加倍。每个实验是每组n = 3个生物学独立样品的结果。GBRAFi 5天期间MITF活性的单细胞延时显微镜分析。上图：BRAFi前后5天分类的GFP-High和GFP-Low细胞的延时图像。显示了来自三个生物学复制品的代表性图像。比例尺，100 µm。底部：MITF高（橙色）和MITF低（蓝色）单元的单单元MITF报告迹线。粗线代表平均响应。h模块2中的Slug，MITF，MART1和PFK相对表达水平分别为高和低亚群，在第5天从细胞中收集并从单细胞数据集中进行分析。每个实验都是每组n = 16个生物学独立样品的结果。一世通过qPCR在已用BRAFi处理5天的MITF-High和MITF-Low day-0细胞中进行qPCR测定，Slug，MITF，Mart1和PFK表达。每个实验是每组n = 3个生物学独立样品的结果。数据表示为平均值±SEM。源数据作为源数据文件提供。
意外分析为从单纯药物到耐药细胞状态的上下路径的存在提供了理论支持。但是，需要进行实验验证来确定单个细胞是否沿着一条路径或另一条路径遵循一条轨迹，或者细胞是否在路径之间随机切换。D0细胞数据上模块2的图暗示了生物学差异，这些差异甚至将未经处理的D0细胞也分成两个亚群（状态1和状态2）（补充图 15）。转录因子MITF及其直接下游靶标MART1的表达水平在区分这两个D0亚群的前四个标记中（补充图 16））。这一发现表明，用药物治疗的MITF-Low细胞可能遵循较低的路径，而MITF-High细胞则可能遵循较高的路径（补充图 17a）。因此，我们产生MITF绿色荧光蛋白（GFP）报道细胞系和分选GFP-高（MITF到高）和GFP-低（MITF-低）的亚群（补充图 17 B和18）。与我们的假设一致，与分类的MITF低的亚群相比，高MITF的高细胞显示Ki67和MITF的水平更高，并且加倍时间更短（图 3c–f）。该数据与已报道的黑素瘤表型从黑素细胞性高增殖状态转变为非黑素细胞性，更具侵入性状态的观察结果一致56。还证实了使用该报告系统可以分离D0细胞中的两个亚群。接下来，我们对MITF高和MITF低的亚群进行了时程实验，以通过流式细胞术分别分析了两个亚群中以下标记的表达，包括Ki67，MART1，p-ERK，NGFR，AXL和MITF。从高维空间到三维（3D）或二维（2D）空间的两个子群体的轨迹的可视化显示了两条轨迹的清晰分离（补充图 19））。这些数据表明，即使在同基因细胞系中，不同的亚群在BRAFi处理后的行为也可能不同。MITF高和MITF低的亚群可能代表在药物治疗后注定要遵循较高和较低路径的细胞。
为了量化在药物治疗期间分选的MITF高和MITF低的亚群之间随机相互转换的频率，我们在BRAFi治疗的5天时间内监测了大量单个细胞中的MITF活性。不出所料，MITF高细胞显示出比MITF低细胞更高的活性（由GFP报告基因定量）（图 3g），并且在两个轨迹之间没有明显的随机切换。
为了进一步确认分选的药物在服药5天后达到其各自的目的地状态，我们量化了在D5的上下路径之间差异表达的标记。对单细胞数据集的挖掘显示，在D5的两条路径（负值模块2和正值模块2）之间，差异表达了包括Slug，MITF，MART1和PFK在内的几个标记（图3h和补充图2b） 。 13和20a）。治疗5天后，通过分析这四个基因在分选的MITF-High和MITF-Low D0细胞中的表达（补充图 20b）），我们发现在处理5天后，它们在分选的MITF-Low细胞中的表达水平显着低于在MITF-High细胞中的表达水平（图 3i）。这些结果实验证明，在药物治疗后，MITF-High和MITF-Low细胞分别沿着上下路径呈现出不同的趋向于药物耐受性的轨迹（补充图 20a，左图）。
MITF被认为是内在药物耐受性的诱因57。为了研究MITF是否在药物反应中驱动分叉，我们生成了MITF稳定敲低的M397细胞系。在处理之前，敲低MITF会诱导细胞以特征性的低水平Ki67进入慢循环（补充图 21a，b），这表明MITF的下调将迫使这些细胞沿着较低的路径过渡。此外，当BRAFi处理5天后，MITF敲低的细胞显示出的显著较低水平SLUG，MITF，MART1和PFK相对于对照（补充图 21C），表明MITF沉默的细胞确实沿着较低的路径运动。因此，MITF被确定为区分我们确定的两种药物反应轨迹的重要分子驱动器。
临界点是不断变化的局势中的关键点，导致新的不可逆转的发展58。临界点分析已广泛用于理解物理系统中的状态转换。最近，更多的研究在生物系统的细胞状态转换的调查应用临界点分析46，59，60。在单元状态转换期间，将有一个临界点，在该临界点发生单元状态的关键变化。如果两个细胞状态被屏障隔开，则临界点可以理解为屏障的峰值，超过该峰值，细胞将不可逆转地过渡到新状态（补充图 22a））。确定这些临界点对于挖掘重要的监管机构至关重要，因为监管机构可以推动这一转变。下药这些稳压器可以提供用于停止这样的过渡（补充图的策略 22B）10，28，46，53，59，61。
对我们的单细胞数据集的意外分析表明，在D1-D3时间窗口（图3a，左图）中，上下路径均具有细胞状态转变（模块1分数的符号变化 ）的特征。为了确定沿两条路径中每条路径的引爆点，我们首先将来自所有时间点的单细胞数据聚类到FLOW-MAP上的14个不同的子簇中（图 4a）。此总体分析假设FLOW-MAP上的每个群集代表过渡处的中间状态。群集1、6、7、8、10、11和12与上方路径对齐，而群集2、3、9、13和14沿下方路径排列（图 4a）。由于有两条路径连接第0天和第5天的像元，因此我们希望有两个引爆点，每条引爆点一条。因此，假定每个群集在沿其自身路径距临界点的距离不同的位置包含单元。已有大量文献证明，当以单细胞分辨率进行分析时，沿着临界点过渡的临界点显示细胞之间的相关性降低，而基因之间的相关性随之增加60。此功能允许使用定量索引来预测高维状态空间中的关键转换。信令网络活动指数（SNAI）10和关键过渡指数（Ic）60都基于此类定量特征进行形式化的是两个已发布的索引，用于识别临界点附近的区域。使用这些索引，我们发现上路径的聚类7和下路径的聚类9在它们各自的路径内显示了这些索引的最大值（图 4b，c和补充图 23 – 25），表明聚类7 1和9最靠近沿两条路径的转折点。
将所有细胞聚类为四个时间点定义的亚群。左面板是FLOW-MAP，其中的单元格通过药物暴露时间进行了颜色编码。右面板是FLOW-MAP，其中的单元格颜色编码为从聚类分析中定义的14个子种群之一。b上部路径的临界点过渡分析。临界点索引SNAI在与上层路径相关联的每个子种群中计算，并用颜色编码到FLOW-MAP上。红色表示较高的SNAI值，而蓝色表示较低的SNAI值。在标记处显示的群集7在上方路径中显示最高SNAI值。C较低路径的临界点过渡分析。在与下部路径关联的每个子种群中计算临界点索引SNAI，并将其颜色编码到FLOW-MAP上。红色表示较高的SNAI值，而蓝色表示较低的SNAI值。簇9（在标记处显示）在较低路径中显示最高SNAI值。d标记-标记相关网络，是从聚类7单元内的SCBC数据中提取的。相关强度以每个边缘的颜色反映（橙色表示正相关，蓝色表示负相关）。Ë标记-标记相关网络，是从聚类9个单元中的SCBC数据中提取的。相关强度以每个边缘的颜色反映（橙色表示正相关，蓝色表示负相关）。f每个网络内每个节点的重要性评分，由节点度（网络内一个节点的连通性的量化）定义。颜色表示群集7或群集9网络中每个节点的节点度值。假定高分结点很重要，其中一些标有星号的结点将通过药物扰动进一步测试。
矿驱动单元状态转变的关键调节剂，我们接下来执行网络分析10，17，40中的簇的临界点7个9我们的假设是，监管网络内显示更高的连接应该是关键毂器，它维护网络。结果，给那些集线器调节器加药可以更有效地破坏网络，从而防止通过临界点过渡到耐药状态（补充图 22b）。这两个网络（针对群集7和群集9）的特征在于结构不同（图 4d，e），表示这些过渡受到不同方式的监管。为了确定网络的中心调节剂（可能是潜在的药物靶标），我们通过计算每个节点的节点度和中心得分来量化相关网络中每个分析物（节点）的连通性（请参见方法）。对于簇7（上部路径），我们发现两种转录因子和代谢酶（包括MITF，PFK，p-LKB，PKM2，LDH2和Slug）均通过两种评分指标显示出高水平的网络参与度（连接性）（图2）。 图4f和补充图 26）。对于簇9（下部路径），TNFR，N-钙黏着蛋白和p-NFκB-p65占主导地位（图 4f和补充图 26）。）。一个有趣的发现是，在聚类7中表现出高分的标记通常在聚类9中表现出低分，反之亦然，这表明两条路径的调控方式不同。重要的是，这些临界点分析结果在一系列聚类数范围内相对稳定（补充图 27）。
为了检查沿两条路径的过渡是否由不同的枢纽调节剂驱动，我们选择使用两种专门针对PKM2或核因子-κB（NFκB）的市售化合物干扰在聚类7和9中确定的各个枢纽节点。假设糖酵解酶PKM2和信号磷蛋白p-NFκB-p65的抑制将分别影响沿上下路径的转变（图 4f和补充图 26）。因此，我们将PKM2抑制剂（PKM2i）或NFκB抑制剂（NFκBi）与BRAFi结合使用来治疗分选的MITF高和MITF低细胞亚群。与我们的假设一致，MITF低的亚群对BRAFi +NFκBi组合更敏感（图。 5a），而MITF高亚群对BRAFi + PKM2i组合更敏感（图 5b）。通过相对于未经修饰的M397细胞，在MITF组合式细胞系上测试相同的药物组合，进一步证实了这一假设（图 5c，d）。因此，沿不同轨迹通过的细胞对PKM2和NFκB的抑制敏感性不同。
在服药前，针对MITF高和MITF低的亚群对 MITF-GFP报告细胞系进行了分类。然后将分选的细胞用BRAFi +NFκBi组合处理5天，然后收集细胞计数。绘制分类的MITF-High和MITF-Low细胞经过BRAFi +NFκBi联合治疗5天后的相对细胞存活率。将生存数据标准化为MITF-High样本。每个实验是每组n = 4个生物学独立样品的结果。b在给药前，将MITF-GFP报告细胞系针对MITF-高和MITF-低亚群进行分类。然后将分类的细胞用BRAFi + PKM2i组合处理5天，然后收集细胞计数。绘制分类的MITF-High和MITF-Low细胞经过BRAFi + PKM2i联合治疗5天后的相对细胞存活率。生存数据被标准化为MITF低样本。每个实验是每组n = 4个生物学独立样品的结果。CMIBRA击倒细胞和对照细胞用BRAFi +NFκBi组合处理5天，然后收集进行细胞计数。绘制分类的对照细胞和MITF-sh细胞经过BRAFi + NFKBi联合治疗5天后的相对细胞存活率。将存活数据相对于对照样品标准化。每个实验是每组n = 5个生物学独立样品的结果。d用BRAFi + PKM2i组合处理MITF组合式细胞和对照细胞5天，然后收集细胞计数。绘制分类的对照细胞和MITF-sh细胞经过BRAFi + PKM2i联合治疗5天后的相对细胞存活率。生存数据针对MITF-sh样本进行了标准化。每个实验都是n的结果 每组= 4个生物学独立的样本。对于a - d中的箱线图，数据为中位数，并标出第一和第三四分位数（框），并指示顶部和底部四分位数（须线）。ë与BRAFi，BRAFi +NFΚBi，BRAFi + PKM2i，和BRAFi +NFκBi+ PKM2i处理5天M397细胞共收集细胞数目计数。绘制了经过BRAFi，BRAFi +NFκBi，BRAFi + PKM2i或BRAFi + PKM2i +NFκBi治疗5天后的细胞相对细胞存活率。存活数据针对进行BRAFi单药治疗的细胞进行了标准化。每个实验是每组n = 4个生物学独立样品的结果。该P -值是由双尾非配对确定ŧ-检验，* P <0.05，** P <0.01，*** P <0.001。数据表示为平均值±SEM。源数据作为源数据文件提供。
考虑到两条轨迹的差分调节器相关性，我们进一步假设通过同时抑制PKM2和NFκB信号来共同阻断两条轨迹可能在防止向BRAFi耐受性过渡方面显示出累加效应。为了验证这一假设，我们使用三联药物组合（BRAFi + PKM2i +NFκBi）在体外处理M397细胞5天，并比较了单药疗法（仅BRAFi）和双药物组合（BRAFi + PKM2i）的最终细胞数和BRAFi +NFκBi）5天。与我们的预测一致，三药组合明显优于双药组合，后者优于单药疗法（图 5e）。）。此外，PKM2i或NFκBi单一疗法对M397细胞的生长抑制作用最小（补充图 28），这意味着这些药物可能通过选择性阻断BRAFi诱导的细胞状态向药物耐受状态的转变而起作用。这些结果表明，上，下路径是独立的，具有不同的调节剂，并且是可独立药物化的。
我们在这里探讨了是否可以通过所遍历的细胞状态空间景观的静态快照的动力学序列来确定连接细胞状态过渡的初始状态和最终状态的细胞轨迹。作为一个模型系统中，我们使用的高塑性，患者来源的M397 BRAF V600E突变的黑素瘤细胞系，其已被证明在治疗与BRAF抑制剂药物幼稚和耐药状态之间可逆地转变8，10，62。尽管单细胞组学工具已被证明对于解决单个给定时间点的组织细胞异质性具有巨大价值，但在此我们试图定量地将细胞异质性与细胞状态变化的动态异质性联系起来。
我们利用基于微流体的SCBC技术来表征药物反应前5天的细胞异质性。由于代谢活性和信号通路在早期药物反应期间均显示出功能变化，因此SCBC在此特别适用，因为它能够同时捕获单个细胞的代谢物和胞质蛋白（和磷蛋白）。但是，与单细胞RNA测序不同，单细胞蛋白质组学通常仅限于分析数十种功能蛋白和代谢产物。为了准确捕获BRAFi处理下M397细胞访问的细胞状态空间，我们首先利用转录组学分析和文献指导来定义20种分析物，包括表型标记，代谢活性，致癌信号和细胞增殖标记，所有这些在初始药物反应期间都会改变。使用此精选面板进行的单细胞分析可轻松解决在BRAFi处理的头几天中，细胞遍历的复杂细胞状态空间。当然，转向大量分析物肯定会提供更深入的表征63，64，65。
我们利用计算和理论方法29，32，33，48，并与其他细胞生物学实验集成，以从SCBC动力学序列快照中推断出单细胞轨迹。使用FLOW-MAP算法对数据集进行降维处理表明，细胞可能会通过两条路径（分别标记为上层和下层）通过细胞状态空间中的一条，从而将未加工过的细胞连接到药物耐受细胞。对相同数据的意外分析解决了时间相关的模块和路径相关的模块。路径依赖模块暗示沿一条路径传播的细胞被生物物理屏障与另一条路径隔开，该生物物理屏障似乎与转录因子MITF及其下游的黑素细胞标记MART1相关。这些分析进一步预测，特定细胞所采取的轨迹由药物治疗之前的MITF水平决定。这些预测已通过实验验证，这支持将计算可视化方法与理论生物物理方法相集成，以深入了解复杂的生物系统。这种方法应广泛适用于其他动态，复杂的生物系统，包括细胞分化，肿瘤发生等方面的研究。
在黑色素瘤的增殖和侵袭性的表型是公知的56，66。MITF，MART1，和Ki67已经被报告为区分这两种表型稳健标记56，66。我们发现，即使在我们的研究中使用的未经处理的等基因M397细胞系中，这两种不同的表型也可以共存。MITF高和MITF低的亚群在没有BRAFi治疗的情况下不仅表现出不同的加倍时间，而且在治疗后也遵循不同的药物反应轨迹。这些发现与黑色素瘤癌细胞表型从黑素细胞状态和高度增殖状态转变为非黑素细胞状态和更具侵入性状态的观察结果一致59。在那项研究中，增殖或侵袭性细胞系在培养中显示出固定的基因表达谱，但是当在体内移植时，每一类都产生了异质性肿瘤，其中含有具有两种表达谱的细胞。与在体外对固定基因表达谱的观察一致，我们在5天的治疗期间未观察到沿不同路径传播的细胞之间的显着相互转化。这些发现表明，这两种表型在体外BRAFi处理的几天中相对稳定。当然，我们的体外研究可能无法完全概括体内的黑色素瘤生物学，其中肿瘤的微环境可以产生强大的影响。对早期BRAFi治疗的黑色素瘤患者的临床样本的分析可能会进一步验证我们发现的临床意义。我们以前曾报道过，在BRAFi处理后，细胞开始从短暂分化开始，进入早期的药物耐受状态（在D3处MITF和MART1表达略有增加）。长时间的BRAF抑制（3周）会产生干细胞样耐药状态，其特征是细胞增殖急剧增加以及MITF和MART1表达丧失8，10。与这些发现一致，我们显示D5细胞（早期药物耐受状态）在BRAFi下显示出低水平的增殖标记Ki67。这些细胞大多数仍表达MART1和MITF，但已经开始启动抗性相关基因程序（葡萄糖摄取的上调和许多干细胞标志物，例如AXL和CDH2）。与批量研究的先前结果一致 67，我们发现短期药物治疗后，MITF高和MITF低的亚群都将变得更加黑素细胞化。这些表明，向较黑素细胞状态的瞬时分化可能是黑素瘤细胞尽管初始细胞状态，但响应BRAF抑制而利用的一般早期作用机制。此外，我们还发现，通过敲低MITF转录因子可以很容易地诱导向MITF-低侵入性表型的转变。这表明复杂的细胞状态态势很可能由很少的主调节剂调节。它还强调了MITF作为调节黑素瘤表型确定的分子驱动因素的重要性。这些发现大大增加了我们对黑色素瘤表型调控的了解，并且通过单细胞分析得到了独特的揭示。
我们的单细胞分析表明，未经处理的细胞同时含有MITF低和MITF高细胞亚群，它们倾向于采取不同的途径发展耐药性。因此，单元的初始状态可能会确定该单元可以采取的路径。甚至在同基因细胞系中两种不同的药物反应轨迹的共存也可以解释所谓的混合反应，这在黑素瘤患者的治疗期间通常会观察到。这样的替代逃逸路径也可以解释为什么黑素瘤对BRAFi靶向治疗如此难治。有趣的是，对于两种途径中的每一种，通过临界点分析和网络分析都鉴定出不同的药物敏感性：发现上部途径易于抑制糖酵解酶PKM2，而下部路径对NFκb-p65抑制敏感。这些不同的药物敏感性结果也与以前有关黑色素瘤浸润表型的大量研究一致：据报道，MITF低，浸润（或间充质）黑色素瘤细胞更依赖于NFκB信号传导。10，我们研究的单细胞分辨率揭示了该观察结果背后的确切分子和细胞动力学。PKM2和NFκB通路的共同抑制显示出抑制肿瘤生长的优异作用。然而，这两种基因都是正常细胞中必不可少的调节剂，它们的抑制作用可能对非恶性组织产生毒性。值得注意的是，MITF的表达水平已显示出与BRAFi敏感性 68相关。因此，将BRAFi的剂量从当前的细胞抑制水平增加到细胞毒性水平可以消除MITF-High亚群及其各自的途径。总之，已解决的异质药物反应轨迹更新了目前对耐药性发展的理解，并可以提供一种有效的方法来鉴定有效的治疗组合。
本研究中使用的患者来源的黑色素瘤细胞系M397先前是根据UCLA IRB批准号11–003254 39生成的。细胞在37℃，5％CO培养2在RPMI 1640与升-谷氨酰胺（Life Technologies），补充了10％的胎牛血清（Omega）和0.2％的抗生素（来自Lonza的MycoZapTM Plus-CL）。定期使用GenePrint®10 System（Promega）对细胞系进行身份验证，以使其早日通过。通过OncoMap 3或Iontrone检查目标基因中突变的存在，并通过PCR和Sanger测序确认。将全部来自Selleck Chemicals LLC的BRAF抑制剂（vemurafenib），PKM2i（化合物3 K）和NFκBi（JSH-23）以指定的浓度溶于二甲基亚砜（DMSO）中，然后再施加到细胞培养基上。将M397细胞以60％融合度铺在10-cm组织培养板上，并用3 µM BRAF抑制剂处理指定的天数。
通过在使用前立即使DNA-抗体或DNA-探针共轭混合物溶液流动，将每个微腔室内的DNA微阵列转换为抗体或纳米探针微阵列。在各个时间点进行胰蛋白酶消化之前，我们用磷酸盐缓冲盐水（PBS）洗涤了死细胞。收集的细胞用Gluc-Bio处理，然后随机装入SCBC内的微腔室进行分析。每个微腔室都有一个测定成分和一个单独的裂解缓冲液容器，并在细胞加载后照相。然后将SCBC在冰上冷却以裂解细胞。在室温下2小时的蛋白质和代谢物捕获期后，冲洗微腔并将阵列上捕获的蛋白质或代谢物转换为荧光读数，并通过Genepix扫描仪（Molecular Devices）进行数字化。
通过自定义的MATLAB代码，提取给定条形码中所有条的平均荧光信号，并将其与该阵列的显微照片进行匹配，以准备一个表格，其中包含微室地址，细胞数以及每种测定的测得的荧光水平蛋白质或代谢产物。收集来自具有单个单元格的微腔室的SCBC读数，以形成m × n矩阵表，其中每一行（m）代表特定的微腔室地址，每一列（n）代表特定分析物的丰度。该矩阵表用于进一步分析。在第0天，第1天，第3天和第5天分别分析了156个，185、162和171个单细胞。
所有FLOW-MAP可视化都是使用GitHub（https://github.com/zunderlab/FLOWMAP/）上提供的FLOWMAPR R包（版本1.2.0）创建的。生成的图具有seed.X = 1，并且没有聚类或下采样。使用FLOWMAPR软件包中描述的蓝色调色板或使用喷射彩虹调色板在Gephi（https://gephi.org/）中生成最终图形。可根据要求提供用于生成精确的FLOW-MAP图的代码。
为了确定细胞状态的边界，我们认为使用突变分析将细胞内许多标志物的变化反卷积为仅几个模块的变化。每个模块代表一组标记，这些标记从一个单元到另一个单元共同变化。这种分析可以极大地简化变更的复杂性，并将其范围缩小到仅几个模块，然后可以将其进一步分解以进行详细的生物学发现。计算上，如先前所述52进行了意外分析。简而言之，在单元c处测得的分析物i的水平ln X i（c）表示为稳态项ln X i 0（Ç）和几个约束（模块）λĴ（Ç）× G ^ IJ表示从稳定状态偏差。每个偏差项是依赖性细胞重量（影响分数）的约束的产物λ Ĵ（Ç）和分析物的该约束（模块）的细胞不依赖性贡献ģ IJ。为了执行意外分析，我们计算矩阵ln X i（c）的奇异值分解。这个因素这个矩阵中，其确定两组参数中需要的surprisal分析的方式：所述拉格朗日乘数（λ Ĵ给定时间点的所有约束（模块），以及所有时间和所有约束i的所有分析物i的G ij（与时间无关）分析物模式。在图 3中，具有最高正模块2得分的前10％的单元被定义为Module2-High单元，而具有最低负10％单元被定义为Module2-Low单元。
在具有硬件自动对焦（ZDC2）的Olympus IX8倒置荧光显微镜和保持37°C，5％CO 2的培养环境的环境室中拍摄电影。使用自动采集软件（METAMORPH，分子设备）每15分钟从多级位置采集差分干涉对比（DIC）和GFP图像。
使用开源软件ilastik（1.3.2版）分割核染色图像，以获取分割的核体。五帧（193帧中的五帧）用作每个位置图像分割的训练集。使用ilastik 1.3.2的像素分类功能将所有193帧的像素分割为“背景”和“单元”。然后，将分段的h5文件与原始核染色电影一起用于ilastik“手动跟踪工作流程”，以从18个电影位置获得110条单细胞核痕迹。根据细胞跟踪结果，使用自定义Python代码从相应的GFP图像中提取所有107条单细胞痕迹的GFP荧光数据（减去背景）。然后，按前50帧的平均GFP水平降序排列107条单细胞GFP迹线。其中， 3克。
聚类之前，按时间点（即，第0天，第1天，第3天和第5天）分隔所有单单元格数据。然后应用Rclusterpp聚类，将细胞聚集成14个亚群。使用“肘方法”确定簇数，所述“肘方法”是基于平方和之内的总和度量69。使用所有默认设置（https://github.com/nolanlab/Rclusterpp），使用Rclusterpp R软件包（0.2.5版）生成了Rclusterpp群集。所有聚类算法均使用在以下标记上聚类的细胞进行：Ki67，Mart1，HIF1a，LDH，AMPKA，p-ERK1，PFK，p-ACAC，Slug和p-LKB。可根据要求提供用于集群的代码。
临界点分析是在两个单独的运行中实施的。每次运行仅在来自两条路径之一的集群上进行。利用两个定量指标来预测高维状态空间中的关键跃迁：IC 60和SNAI 10。这两个指标都根据临界点的数学特征进行了形式化：标记物与标记物的相关性增加以及细胞异质性的增加。SNAI值定义为蛋白质-蛋白质相关性决定因素的倒数。Ic值定义为所有蛋白质对蛋白质相关系数对的平均值与所有细胞对细胞相关系数对的平均值之比，并进行计算。可根据要求提供用于计算单个小区群集的SNAI / Ic索引的代码。
使用Hmisc R程序包（版本4.2-0，可从https://cran.r-project.org/web/packages/Hmisc/index.html获得）在14个群集中的每个群集上计算成对相关矩阵。计算Spearman的相关性。从Hmisc程序包输出的相关性生成成对的相关值矩阵。Bonferroni校正了p使用-值通过统计显着性过滤相关网络，并使用自定义MATLAB代码绘制相关网络。使用igraph R包（版本1.2.4.1）计算每个相关网络中每个标记的中心得分和节点度。每个标记的两个分数均从0缩放到1，以进行并排比较，并绘制图形以显示计算相关方法之间中心行为的标记间差异。可根据要求提供用于执行相关网络分析的代码。
根据制造商的规程，使用RNeasy Mini Kit或RNeasy plus Micro Kit（Qiagen）从细胞中提取RNA。使用iScript cDNA合成试剂盒（Bio-Rad）从提取的总RNA中合成第一链cDNA。使用特异性引物通过基于SYBR Green的实时定量逆转录PCR分析人Slug，MITF，MART1和PFK转录物的表达（补充表 1）。将数据标准化为RPL19的表达，并表示为倍数变化。
靶向MITF编码序列和对照shRNA的短发夹RNA（shRNA）购自Santa Cruz。通过基于慢病毒的载体及其包装载体psPAX2和pMD2.G的瞬时转染，在HEK-293T细胞中产生了编码对照shRNA和MITF shRNA的慢病毒。收集病毒，在48小时后通过0.45 µm注射器过滤器过滤，然后在900 g和30°C下用补充了10 µg / mL聚乙烯的病毒上清液自旋感染M397细胞90分钟。从转导后3天开始，使用嘌呤霉素选择转导的细胞。
将人酪氨酸酶启动子从pLightSwitch Prom S700747（SwitchGear Genomics，Carlsbad，CA）亚克隆到慢病毒载体主链的BamH1和EcoRI位点，从而驱动Zsgreen基因的表达。如上所述产生了慢病毒颗粒，以稳定地转导M397细胞。通过单细胞分选和扩增进一步产生了克隆细胞系。然后通过BD FACSAria融合细胞分选仪将细胞分类为GFP高和GFP低种群，以进行进一步处理和分析。
在黄色荧光蛋白（YFP）和DIC通道中的EVOS FL自动成像系统（Fisher Scientific）上以10倍（奥林巴斯，10X FL PH，0.3 NA）采集图像。MITF高井和MITF低井之间的光或激光强度，曝光和增益设置为相同。
然后将所有细胞用BD Bioscience的Fix-Perm缓冲液固定，然后对MART1，NGFR，AXL，p-ERK和Ki67的细胞内染料偶联抗体染色。使用购自Thermo Fisher的Attune NxT流式细胞仪进行流式细胞仪分析，并使用FlowJo软件分析数据。为了可视化细胞状态转换轨迹，分别通过意外分析和t-SNE将六维流式细胞术数据投影到3D和2D空间上。
洗涤细胞并从培养板上用胰蛋白酶消化，然后在500× g和4°C 下离心 5分钟以沉淀细胞。然后将细胞沉淀重悬于含有1％BS牛血清白蛋白的PBS中，然后进行荧光激活的细胞分选。门控策略显示在补充图 29中。
使用衰老β-半乳糖苷酶染色试剂盒（Cell Signaling，9860）根据制造商的方案定量衰老细胞的百分比。简而言之，将板上的细胞用PBS冲洗，然后用固定液固定15分钟。固定后，将板用PBS冲洗两次，然后在染色溶液中于37°C孵育过夜。在相差显微镜下检查板。
通过处理在相应条件下培养的指示细胞系进行细胞凋亡测定。在流式细胞仪分析之前，在室温下用膜联蛋白V-异硫氰酸荧光素和碘化丙啶对细胞染色15分钟。使用未染色的对照确定门。
将M397细胞以最佳融合度铺在装有新鲜培养基的六孔板上。每两天补充一次介质（含药物或DMSO）。染色时，将4％的多聚甲醛施加到菌落上以固定细胞，并使用0.05％的结晶紫溶液对菌落染色。
|
|
