RNA-Seq探索基底型乳腺癌干细胞差异表达基因

李　婷, 李　璐, 王姝越, 戴晓峰, 白仲虎

(1. 江南大学 粮食发酵工艺与技术国家工程实验室， 无锡 214122； 2. 江南大学 生物工程学院， 无锡 214122； 3. 上海科技大学 生命科学与技术学院， 上海 201012)

乳腺癌是我国妇女中最常见的恶性肿瘤之一，其发病率高居榜首[1]。近几年，随着早期诊断技术以及治疗方式的发展，乳腺癌致死率已显著降低。乳腺癌是一种高度异质性疾病，至少包含5种截然不同的亚型，即腔A型(luminal A)、腔B型(luminal B)、HER2型、基底型(basal)和正常细胞样型(normal-like)[2]。目前，临床上针对luminal型和HER2型乳腺癌均已经有了较好的治疗方式，如激素受体阳性的luminal型常采用内分泌疗法并辅以他莫昔芬靶向治疗，该类乳腺癌患者预后较好；HER2型乳腺癌患者虽预后较差，但该类型患者对化疗及赫赛汀(Herceptin)靶向治疗均具有良好的反应；而基底型乳腺癌患者由于低表达雌激素受体(ER)、孕激素受体(PR)以及人类表皮生长因子受体2(HER2)，致使其对内分泌治疗及靶向治疗药物均不敏感，该类型患者恶性程度极高，常伴随复发、转移，5年存活率不到15%[3-4]。因此，针对基底型乳腺癌进行研究，找出其特定有效的靶标将对该类病人起到至关重要的作用。

肿瘤干细胞理论的提出为探索基底型乳腺癌的治疗方式带来了曙光。该理论认为乳腺癌中存在乳腺癌干细胞亚群，可自我再生并驱动肿瘤发生，是产生治疗抗性及复发的根源[5]。基底型乳腺癌许多特性与癌症干细胞十分相似[6-7]，许多研究表明针对肿瘤干细胞设计靶向创新型药物可能是治疗基底型乳腺癌的有效方法[8]。乳腺癌干细胞的识别方法中表面标志物CD44高表达CD24低表达(CD44+CD24-/low)分离法是最早被发现可以用来分选乳腺癌干细胞的方法，CD44+CD24-/low干细胞与恶性程度高度相关，转移部位乳腺癌表达CD44+CD24-/low干细胞比例显著高于原发灶[9]，且该亚群细胞具有间充质特性，与肿瘤的转移及化疗、放疗抗性密切相关[10]。因此对CD44+CD24-/low干细胞亚群进行深入研究，尤其是从整体转录水平研究其特征基因，将更加深入地了解该干细胞亚群的分子机制，找到干细胞相关靶点基因，为基底型乳腺癌患者提供靶向治疗的新策略，从而从根本上解决基底型肿瘤易复发、转移的难题。

本研究从乳腺癌干细胞方向出发，将流式分选的基底型乳腺癌细胞系HCC1937和SUM149PT的干细胞及非干细胞作为模型，通过高通量转录组测序(RNA-Seq)筛选出干细胞和非干细胞的差异表达基因(differentially expressed genes，DEGs)，并对DEGs进行生物信息学分析，为进一步探索基底型乳腺癌的靶向治疗提供了理论基础。

1　材料与方法

1.1　主要试剂和仪器

RPMI 1640(Hyclone)，Ham′s F-12、胎牛血清(FBS)、胰蛋白酶购于Gibco公司；D-PBS、HEPES(吉诺生物)，胰岛素(Insulin)购于无锡市第四人民医院；氢化可的松(Hydrocortinsone)购于TOCRIS；RNA提取试剂盒(TIANGEN)，反转录试剂盒、SYBR qRT-PCR试剂盒(TaKaRa)，引物设计由苏州金唯智测序合成；PE标记的鼠抗人CD24抗体、APC标记的鼠抗人CD44抗体(BD)，流式缓冲液(D-PBS+1%FBS)，流式细胞分选仪BD FACSAria II(BD)，qRT-PCR仪(Applied Biosystem)。

1.2　细胞培养及乳腺癌干细胞分选

1.2.1细胞培养

细胞系SUM149PT来源于Asterand，HCC1937来源于ATCC (American Type Tissue Culture Collection)，SUM149PT细胞系用含有5% FBS、10 mmol/L HEPES、5 μg/mL Insulin和1 μg/mL Hydrocortinsone的Ham′s F-12培养基，HCC1937用含10% FBS的RPMI 1640培养基，于37℃、5%CO2培养箱中培养。

1.2.2乳腺癌干细胞分选及纯度鉴定

收集对数期细胞，D-PBS洗涤1次，用流式缓冲液重悬细胞并计数，调节细胞浓度至107个/mL。设置阴性对照组、实验组，对照组加入100 μL细胞悬液，实验组加入500 μL细胞悬液以及CD24-PE、CD44-APC各100 μL，混匀后室温避光孵育，每隔10 min混匀细胞1次，孵育20 min后流式缓冲液洗涤细胞1次，加入流式缓冲液重悬，流式细胞仪上机分选，得到各细胞系的CD44+CD24-/low干细胞亚群及non-CD44+CD24-/low非干细胞亚群，每个细胞系分选3次，分选后的细胞再次用流式细胞仪检测各亚群细胞纯度，保证各亚群细胞纯度都在95%以上，纯度合格后离心收集分选后的细胞，Trizol重悬，液氮速冻后转移至-80℃保存用于RNA-Seq。

1.3　转录组测序

RNA-Seq由上海欧易生物医学科技有限公司所做，简要步骤如下：分别提取分选的HCC1937、SUM149PT干细胞与非干细胞总RNA，使用DNase消化DNA后富集mRNA，用打断试剂将mRNA打断成短片段并以此为模板合成cDNA，使用试剂盒纯化双链cDNA再进行末端修复、加A尾，连接接头后选择合适的片段大小进行PCR扩增，构建好的文库用Agilent 2100 Bioanalyzer质检合格后，使用Illumina HiSeqTM 2500测序仪进行测序，上样量为50 ng。

1.4　测序质量评估

测序得到的原始数据(Raw reads)使用FastQC(http://www.bioinformatics.babraham.ac.uk/projects/fastqc/)进行质量评估，并通过软件NGS QC Toolkit v2.3.3(http://59.163.192.90:8080/ngsqctoolkit/)对原始数据进行质控(QC)，过滤得到Clean reads，用tophat/bowtie2将clean reads与所选的参考序列(GRCh38)进行比对。统计reads在参考序列上的分布情况及覆盖度，判断比对结果是否通过第2次质控。

1.5　差异表达基因的筛选

比对结果通过QC后，对基因进行表达分析，基因表达量使用FPKM(Fragments Per kb per Million reads)方法计算，公式如下(以基因A为例)[11]：

以校正的P值(false discovery rate, FDR)≤0.05、差异倍数的对数绝对值(︱log2FoldChange︱，︱FC︱)≥1为阈值，筛选干细胞和非干细胞的DEGs。

1.6　差异表达基因的GO分析和KEGG富集

筛选出DEGs后，我们对DEGs进行GO和KEGG富集分析，并用超几何分布检验(hypergeometric test)方法计算每个GO条目和每个通路中DEGs富集的显著性，P≤0.05认为DEGs显著富集。

1.7　差异表达基因的qRT-PCR验证

提取对数期细胞总RNA并反转录成cDNA，取1 μg cDNA为模板进行PCR(方法参照文献[3])检测目的基因的相对表达量。以GAPDH为内参基因，采用2-△△Ct方法计算基因表达丰度。各基因qRT-PCR引物序列如下：GAPDH引物(上游：5′-CCCACTCCTCCACCTTTGAC-3′，下游：5′-ATGAGGTCCACCACCCTGTT-3′)，TGFBR2引物(上游：5′-AGGGCAAGTGCGTAAATGGAG-3′，下游：5′-TGGGCAGATTTCACGGCTG-3′)，IL6引物(上游：5′-ACTCACCTCTTCAGAACGAATTG-3′，下游：5′-CCATCTTTGGAAGGTTCAGGTTG-3′)和MMP2引物(上游：5′-TACAGGATCATTGGCTACACACC-3′，下游：5′-GGTCACATCGCTCCAGACT-3′)。

2　结果与分析

2.1　乳腺癌细胞系干细胞的分选与纯度鉴定

为了在全转录水平探索干细胞的分子机理，挖掘其关键DEGs，我们根据此前文献报道，选取干细胞比例在50%左右的基底型乳腺癌细胞系SUM149PT和HCC1937进行流式分选[12]，分别收集各细胞系的CD44+CD24-/low干细胞和non-CD44+CD24-/low非干细胞亚群(图1-A、B)，每个细胞系分选3次。将分选的细胞进行二次流式检测分选的细胞亚群纯度，结果表明分选的HCC1937、SUM149PT干细胞亚群与非干细胞亚群纯度都在95%以上(图1-C、D、E和F)，满足下一步试验需求。

图1 流式分选乳腺癌细胞系的干细胞与非干细胞亚群

A： SUM149PT细胞系的分选，红色框线内为干细胞亚群，蓝色框线内为非干细胞亚群； B： HCC1937细胞的分选； C： 分选后的SUM149PT干细胞纯度检测； D： 分选后SUM149PT非干细胞纯度检测； E： 分选后HCC1937干细胞纯度检测； F： 分选后的HCC1937非干细胞纯度检测

2.2　测序质量评估

2.2.1测序数据预处理

分选的HCC1937以及SUM149PT干细胞和非干细胞分别提取总RNA，对各个细胞亚群进行RNA-Seq，每组3个重复。测序后得到SUM149PT干细胞、SUM149PT非干细胞、HCC1937干细胞、HCC1937非干细胞的原始片段(Raw reads)分别为63 736 099、63 966 929、63 364 586和64 750 000条，Raw reads经质控过滤得到的Clean reads分别为62 485 489、62 680 296、61 810 252和63 261 312条，分别占总Reads的98.04%、97.99%、97.55%和97.70%，其中过滤后的Q30比率分别为95.76%、95.80%、95.47%和95.58%(表1)，表明测序质量较好可以用于下一步分析。

2.2.2参考序列比对分析

质控后的Clean reads与参考序列(GRCh38)进行比对，结果表明(如表2)，能定位到基因组上的测序序列数量(Total mapped reads)占总reads数的88%以上，唯一比对到基因组序列(Uniquely mapped reads)占总reads的80%以上，这表明过滤的数据不存在污染并且选取的参考基因组合适。

2.2.3测序随机性评估及reads覆盖度分析

将reads在基因上的位置标准化到相对位置，通过统计基因的不同位置比对上的reads数来检测RNA-Seq测序过程中mRNA序列打断的随机性。结果如图2，不同样本的reads在基因5′端到3′端分布相对均匀，打断随机性好。此外，我们计算了基因覆盖度(每个基因被reads覆盖的百分比)，结果表明覆盖度在90%～100%的基因(枚红色)占70%以上(图3)。

表1 不同样本过滤得到的测序产量

1，2，3: 各细胞系干细胞和非干细胞亚群转录组测序的3个重复；P2：干细胞亚群；P3：非干细胞亚群。下同

表2 不同样本的比对分析

2.3　差异表达基因的筛选

在去除FPKM<0.5的低表达基因后，我们以FDR≤0.05、︱FC︱≥1为标准，筛选HCC1937干细胞和非干细胞、SUM149PT干细胞和非干细胞的DEGs，结果如图4-A，HCC1937中干细胞和非干细胞的DEGs有413个(165个基因在干细胞中上调，248个基因下调)，SUM149PT中干细胞和非干细胞的DEGs有832个(135个基因在干细胞中上调，697个基因下调)。通过VENNY分析，我们得到HCC1937和SUM149PT细胞系中共有的上调DEGs 39个(图4-B)，︱FC︱在1.0～3.5倍，下调DEGs 95个(图4-C)，︱FC︱在1.0～7.4倍。

图2 不同样品的测序随机性分布Fig 2 Randomness assessment of the samples

A-C：SUM149PT干细胞测序随机性评估；D-F：SUM149PT非干细胞测序随机性分布；G-I：HCC1937干细胞测序随机性分布；J-L：HCC1937非干细胞测序随机性分布。每行代表样品的3个重复

2.4　差异表达基因GO功能分析

我们将筛选得到的134个DEGs进行GO功能富集分析，得到的1841个GO功能注释中有1415个为生物过程(76.86%，biological processes, BP)，其中富集显著的BP有血管生成、细胞迁移和间叶细胞增殖的正调控等，这说明干细胞DEGs很可能主要是通过调节生物学过程来促进自身的自我更新及较强的转移侵袭特性；此外富集到细胞组分(cellular_component，CC)的GO条目占8.33%，最显著富集的有胞外区、细胞外间隙、顶端质膜和黏附功能等；分子功能(molecular functions， MF)的GO条目占15.26%，主要表现在细胞因子活性、细胞黏附分子结合和钙离子结合等(表3)。

2.5　差异表达基因通路富集分析

在大多数生物学过程中基因常通过相互作用参与其中，为了进一步探索这些干细胞中的DEGs生物学功能，我们通过KEGG数据库去研究DEGs主要参与的信号通路。134个DEGs共得到77条显著富集的信号通路(P≤0.05)，DEGs显著富集的通路包括肿瘤相关通路、蛋白多糖与肿瘤通路、PI3K-Akt信号通路，调节肿瘤干细胞多能性通路和黏着斑通路等(图5)。其中最富集的通路以及富集基因最多的通路都是肿瘤相关通路，14个富集到该通路中的DEGs中，有8个基因上调。

图3 不同样品的基因覆盖度Fig 3 Gene coverage distribution in different cell subset

A-C：SUM149PT干细胞亚群的基因覆盖度；D-F：SUM149PT非干细胞亚群的基因覆盖度；G-I：HCC1937干细胞亚群的基因覆盖度；J-L：HCC1937非干细胞亚群的基因覆盖度。每行代表样品的3个重复

图4 干细胞和非干细胞的DEGs

A： HCC1937和SUM149PT细胞系中干细胞和非干细胞的DEGs，Venny分析共有的DEGs； B: 共有的DEGs在干细胞中上调的数目； C： 共有的DEGs在干细胞中下调的数目

图5 DEGs富集最显著的10条KEGG通路

2.6　差异表达基因的qRT-PCR验证

为验证转录组测序筛选的DEGs准确性，我们随机选取了3个DEGs，采用qRT-PCR检测了它们在乳腺癌细胞系SUM149PT和HCC1937干细胞和非干细胞间的mRNA表达水平差异。结果(如图6)表明，qRT-PCR检测的TGRBR2、IL6和MMP2 3个DEGs的变化趋势和转录组测序结果一致。

图6 qRT-PCR验证转录组结果

A： SUM149PT干细胞和非干细胞的DEGs验证； B： HCC1937干细胞和非干细胞的DEGs验证

表3 DEGs在各个生物学过程中富集最显著的10个GO条目(A. 生物学过程，B. 细胞组分，C. 分子功能)

3　讨论与结论

基底型乳腺癌多发于年轻和绝经前女性，肿瘤状态常为晚期，该类患者对内分泌疗法不敏感，且尚无较好的靶向治疗策略[4]。随着肿瘤干细胞研究的深入，越来越多的证据表明针对此类亚群细胞进行研究找出其潜在靶点可能是治疗基底型乳腺癌的新曙光。RNA-Seq的发展为我们全面深入地研究肿瘤干细胞的分子机制带来了希望。在本研究中，我们采用干细胞标志物CD44+CD24-/low流式分选了基底型乳腺癌细胞系的干细胞和非干细胞亚群，通过转录组测序分析了转录水平上干细胞和非干细胞之间DEGs，用以全面深入地探索乳腺癌干细胞的基因调控机制。

肿瘤干细胞是指存在于肿瘤当中未分化的一部分亚群，多处于休眠状态，具有自我更新及多向分化潜能，是导致肿瘤复发、转移以及抗药性的根源所在[5]。在众多的乳腺癌干细胞分离方法中，我们发现采用CD44+CD24-/low标志物进行分选的方法研究最为广泛，AL-Hajj等首次发现该方法可以分离乳腺癌干细胞，这一亚群细胞只需200个便可在免疫缺陷型小鼠中成瘤[13]，且该亚群细胞具有基底样型特征，可提高肿瘤的侵袭能力并带来化疗抗性[10]。因此从CD44+CD24-/low亚群出发，靶向杀伤这一亚群细胞可能是治疗基底型乳腺癌最有效的治疗手段。在对转录组分析的结果中发现，在HCC1937和SUM149PT两个细胞系共有的DEGs中，与上皮-间充质(EMT)过程相关的IL6、MMP2和FOXC2在干细胞中显著高表达。IL6是白介素家族的重要一员，有报道称IL6可促使乳腺细胞发生上皮-间充质(EMT)转化，促使肿瘤产生药物抗性[14]；MMP2是基质金属蛋白酶家族成员，它通过酶解作用激活TGF-β促使EMT进而促使肿瘤转移[15]；差异最显著的基因FOXC2在EMT过程中被激活[16]，也有报道称shRNA敲除FOXC2后可显著降低乳腺癌的转移[17]。而我们的转录组结果中IL6、MMP2和FOXC2在干细胞中均显著高表达，这表明这些基因很可能是干细胞的关键靶点，可以用来指导肿瘤的靶向治疗。我们在这些DEGs中还发现，GATA3，FOXA1和ESR1转录因子在干细胞中显著下调，这些基因都曾被报道过是管腔型乳腺癌的重要分型转录因子[18]。结合这些报道，可以发现从干细胞入手挖掘到的基因不仅可以用来指导靶向治疗，还可以用来帮助识别乳腺癌亚型分类。

为了进一步研究DEGs所涉及的生物学功能及参与的关键信号通路，我们对134个DEGs进行了GO功能和KEGG通路富集。在GO分析中，我们发现大多数DEGs都富集在肿瘤恶性程度密切相关的GO条目中，如血管生成、细胞迁移、增殖的正调控，细胞黏附分子结合等。在KEGG富集通路中，值得关注的是最富集的通路与通路中富集DEGs最多地出现了重叠，即肿瘤相关通路，该通路富集的14个DEGs中有8个基因(AXIN2、F2R、FGFR1、IL6、LAMC3、PRKCA、MMP2、TGFBR2)在干细胞中显著上调。在最富集的前10条通路中也发现了PI3K-Akt信号通路，该通路被报道在肿瘤干细胞中起着重要作用，AKT的激活会促使乳腺癌细胞转移和侵袭[8]。

本研究从乳腺癌干细胞入手，通过RNA-Seq初步探索了干细胞和非干细胞DEGs以及这些基因所富集的关键GO条目和KEGG信号通路，为深入研究乳腺癌干细胞的分子机制以及靶向治疗基底型乳腺癌提供了理论基础，接下来仍需对关键的干细胞DEGs开展更进一步的研究。

[1]CHEN W, ZHENG R, BAADE P D, et al. Cancer statistics in China, 2015 [J]. CA: a Cancer Journal for Clinicians, 2016, 66(2):115-132.

[2]CHEANG M C, VODUC D, BAJDIK C, et al. Basal-like breast cancer defined by five biomarkers has superior prognostic value than triple-negative phenotype [J]. Clinical Cancer Research, 2008, 14(5):1368-1376.

[3]白仲虎, 李 璐, 戴晓峰. FOXA1:乳腺癌分子分型的重要转录因子 [J].生物学杂志, 2017, 34(1):5-10.

[4]RAKHA E A, REISFILHO J S, ELLIS I O. Basal-like breast cancer: a critical review [J]. Journal of Clinical Oncology Official Journal of the American Society of Clinical Oncology, 2008, 26(15):2568-2581.

[5]PAJONK F, VLASHI E, MCBRIDE W H. Radiation resistance of cancer stem cells: the 4 R′s of radiobiology revisited [J]. Stem Cells, 2010, 28(4):639-648.

[6]SMALLEY M, ASHWORTH A. Stem cells and breast cancer: a field in transit [J]. Nature Reviews Cancer, 2003, 3(11):832-844.

[7]TURNER N C, REISFILHO J S. Basal-like breast cancer and the BRCA1 phenotype [J]. Oncogene, 2006, 25(43):5846.

[8]GINESTIER C, LIU S, DIEBEL M E, et al. CXCR1 blockade selectively targets human breast cancer stem cells in vitro and in xenografts [J]. Journal of Clinical Investigation, 2010, 120(2):485-497.

[9]BALIC M, LIN H, YOUNG L, et al. Most early disseminated cancer cells detected in bone marrow of breast cancer patients have a putative breast cancer stem cell phenotype [J]. Clinical Cancer Research, 2006, 12(19):5615-5621.

[10]SHERIDAN C, KISHIMOTO H, FUCHS R K, et al. CD44+/CD24-breast cancer cells exhibit enhanced invasive properties: an early step necessary for metastasis [J]. Breast Cancer Research, 2006, 8(5):R59.

[11]TRAPNELL C, WILLIAMS B A, PERTEA G, et al. Transcript assembly and quantification by RNA-Seq reveals unannotated transcripts and isoform switching during cell differentiation [J]. Nature Biotechnology, 2010, 28(5):511-515.

[12]MAROTTA L L, ALMENDRO V, MARUSYK A, et al. The JAK2/STAT3 signaling pathway is required for growth of CD44+CD24 stem cell like breast cancer cells in human tumors [J]. Journal of Clinical Investigation, 2011, 121(7):2723-2735.

[13]AL-HAJJ M, WICHA M S, BENITO-HERNANDEZ A, et al. Prospective identification of tumorigenic breast cancer cells [J]. Proceedings of the National Academy of Sciences, 2003, 100(7): 3983-3988.

[14]SULLIVAN N J, SASSER A K, AXEL A E, et al. Interleukin-6 induces an epithelial|[ndash]|mesenchymal transition phenotype in human breast cancer cells [J]. Oncogene, 2009, 28(33):2940.

[15]GIALELI C, THEOCHARIS A D, KARAMANOS N K. Roles of matrix metalloproteinases in cancer progression and their pharmacological targeting [J]. Febs Journal, 2011, 278(1):16-27.

[16]BATTULA V L, EVANS K W, HOLLIER B G, et al. Epithelial-mesenchymal transition-derived cells exhibit multi-lineage differentiation potential similar to mesenchymal stem cells [J]. Stem Cells, 2010, 28(8):1435-1445.

[17]MANI S A, YANG J, BROOKS M, et al. Mesenchyme Forkhead 1(FOXC2) plays a key role in metastasis and is associated with aggressive basal-like breast cancers [J]. Proceedings of the National Academy of Sciences of the United States of America, 2007, 104(24):10069-10074.

[18]韩明利, 吴诚义. 乳腺癌干细胞与基因分型 [J]. 中国普外基础与临床杂志, 2010, 17(3):302-305.