王建茹，李曉輝

0 引 言

動(dòng)脈粥樣硬化(atherosclerosis，AS)是一種由免疫系統(tǒng)介導(dǎo)的以動(dòng)脈壁炎癥和斑塊形成為特征的慢性動(dòng)脈疾病，是多種心腦血管疾病的共同病理基礎(chǔ)[1]。頸動(dòng)脈粥樣硬化(carotid atherosclerosis，CAS)是AS的一種臨床類型，因其解剖位置表淺，探查方便，常被作為全身AS的“窗口”，來(lái)判定AS的整體情況[2]。CAS發(fā)病機(jī)制與炎癥免疫反應(yīng)、脂質(zhì)代謝及非編碼RNA等關(guān)系密切，但其機(jī)制仍未完全闡釋清楚。眾所周知，巨噬細(xì)胞作為免疫細(xì)胞中重要的成員之一，在AS病理過(guò)程中的致病作用一直是研究的熱點(diǎn)和重點(diǎn)[3]。因此，深入探索巨噬細(xì)胞在CAS中的病理機(jī)制，挖掘關(guān)鍵基因，對(duì)于更好地理解CAS的發(fā)病機(jī)制以及防治靶點(diǎn)具有積極的作用。

單細(xì)胞RNA測(cè)序(single-cell RNAsequencing，scRNA-Seq)作為新一代的高通量測(cè)序技術(shù)，可進(jìn)行單細(xì)胞的轉(zhuǎn)錄組分析，識(shí)別新的或罕見(jiàn)的細(xì)胞亞型，為深入探析不同細(xì)胞類型中單細(xì)胞的行為、機(jī)制及其與機(jī)體的關(guān)系提供了新的方法和思路，在心血管疾病中方面產(chǎn)生了革命性的影響[4-8]。近年來(lái)，已有多項(xiàng)研究利用scRNA-Seq技術(shù)在深化理解AS病理機(jī)制，促進(jìn)潛在藥物靶點(diǎn)研發(fā)和臨床診斷標(biāo)志物的篩選等方面做出了重要貢獻(xiàn)[9-15]。

現(xiàn)階段，高通量技術(shù)聯(lián)合生物信息學(xué)的方法已被廣泛用于CAS的研究中，為其病理機(jī)制、發(fā)展轉(zhuǎn)歸和診療等提供了新見(jiàn)解。本研究結(jié)合scRNA-Seq技術(shù)和生物信息學(xué)分析探索了人CAS組織中巨噬細(xì)胞的特征基因，以期更好的理解CAS的病理過(guò)程及巨噬細(xì)胞介導(dǎo)CAS的潛在作用機(jī)制，為其診斷、治療提供新的參考。

1 材料與方法

1.1 數(shù)據(jù)來(lái)源及預(yù)處理從GEO數(shù)據(jù)庫(kù)(http://www.ncbi. nlm.nih.gov/geo/)中下載CAS的scRNA-Seq數(shù)據(jù)集GSE159677。該數(shù)據(jù)集的樣本來(lái)自3例CAS患者，包括3個(gè)動(dòng)脈粥樣硬化斑塊(scRNA-Seq-AS組)和3個(gè)匹配的頸動(dòng)脈近端相鄰組織(scRNA-Seq-對(duì)照組)。利用hdf5r包和Seurat包分別對(duì)每個(gè)樣本的數(shù)據(jù)進(jìn)行讀取，然后再分別對(duì)兩組的樣本數(shù)據(jù)進(jìn)行合并后保存。

同時(shí)，下載基因芯片數(shù)據(jù)集GSE43292及其平臺(tái)文件GPL6244-17930。該數(shù)據(jù)集的組織來(lái)自32例CAS患者，包括32個(gè)動(dòng)脈粥樣硬化斑塊組織(AS組)和32個(gè)匹配的遠(yuǎn)端完整的頸動(dòng)脈組織(對(duì)照組)。依據(jù)平臺(tái)文件對(duì)表達(dá)數(shù)據(jù)進(jìn)行注釋，將探針矩陣轉(zhuǎn)化為基因矩陣，剔除基因表達(dá)量低的基因，基因?qū)?yīng)多個(gè)探針則利用avereps函數(shù)取均值。然后，再利用normalizeBetweenArrays函數(shù)對(duì)數(shù)據(jù)進(jìn)行矯正，并保存以進(jìn)行后續(xù)的特征基因驗(yàn)證和免疫細(xì)胞浸潤(rùn)分析。

1.2scRNA-Seq數(shù)據(jù)的質(zhì)控和降維聚類利用R語(yǔ)言讀取GSE159677中的scRNA-Seq-AS組和scRNA-Seq-對(duì)照組樣本數(shù)據(jù)，并對(duì)重復(fù)基因取均值。利用Seurat包對(duì)樣本數(shù)據(jù)進(jìn)行過(guò)濾。然后，使用NormalizeData函數(shù)對(duì)數(shù)據(jù)進(jìn)行標(biāo)準(zhǔn)化，再利用FindVariableFeatures函數(shù)計(jì)算細(xì)胞間基因的變異系數(shù)(標(biāo)準(zhǔn)差)。使用RunPCA函數(shù)對(duì)scRNA-Seq-AS組的樣本數(shù)據(jù)進(jìn)行主成分分析( principal component analysis，PCA)，選取P< 0.05的主成分，使用RunTSNE函數(shù)進(jìn)行t分布隨機(jī)鄰接嵌入(t-distributed stochastic neighbor embedding，t-SNE )聚類。

1.3細(xì)胞類型注釋對(duì)scRNA-Seq-AS組的樣本數(shù)據(jù)降維聚類后，使用Seurat包中的FindAllMarker函數(shù)找出每個(gè)細(xì)胞亞群(Cluster)相對(duì)其它Cluster差異表達(dá)的基因，然后利用SingleR包對(duì)不同的細(xì)胞亞群進(jìn)行注釋。

1.4差異基因篩選利用FindAllMarker函數(shù)，以|log2(fold change, FC)|>1和矯正后的P<0.05為閾值，對(duì)scRNA-Seq-AS組和scRNA-Seq-對(duì)照組間的DEGs，以及scRNA-Seq-AS組不同Cluster間的DEGs進(jìn)行篩選[16-17]。然后，將組間的DEGs和scRNA-Seq-AS組巨噬細(xì)胞的特異性DEGs取交集(即共同DEGs)用于后續(xù)分析。

1.5基因功能富集分析利用clusterProfiler包對(duì)共同DEGs進(jìn)行GO和KEGG通路富集分析，并利用ggplot2、GOplot包進(jìn)行可視化。

1.6蛋白互作網(wǎng)絡(luò)的構(gòu)建及特征hub基因的篩選將共同DEGs上傳至STRING數(shù)據(jù)庫(kù)(https://string-db.org/)，選擇“Multiple proteins”模式，互作得分設(shè)置為中等置信度≥0.400，進(jìn)行PPI分析，并將用Cytoscape3.7.2軟件構(gòu)建PPI網(wǎng)絡(luò)。利用CytoHubba插件的MCC、DMNC、MNC、Degree、EPC算法來(lái)篩選PPI網(wǎng)絡(luò)中的巨噬細(xì)胞特征hub基因(即5種算法中Top10基因的交集基因)[17-18]。

1.7巨噬細(xì)胞特征hub基因的驗(yàn)證從預(yù)處理后的GSE43292數(shù)據(jù)集中提取1.6中篩選出來(lái)的各樣本巨噬細(xì)胞特征hub基因的表達(dá)量，繪制箱體圖對(duì)特征基因的組間表達(dá)差異情況進(jìn)行可視化，并利用pROC包繪制特征基因的ROC曲線，計(jì)算AUC值。

1.8免疫細(xì)胞浸潤(rùn)分析CIBERSORT和ssGSEA是廣泛應(yīng)用的量化免疫細(xì)胞浸潤(rùn)狀態(tài)的方法[19]。利用這兩種方法對(duì)預(yù)處理后GSE43292的表達(dá)譜矩陣進(jìn)行模擬分析，進(jìn)而獲得所有樣本的免疫細(xì)胞的浸潤(rùn)模式，利用Wilcoxon檢驗(yàn)以P< 0.05為閾值篩選兩組間的差異免疫細(xì)胞；然后，提取AS樣本中巨噬細(xì)胞含量的相對(duì)比例。

1.9特征hub基因與AS中巨噬細(xì)胞浸潤(rùn)狀態(tài)的相關(guān)性分析采用Pearson法對(duì)GSE43292數(shù)據(jù)集中特征hub基因的表達(dá)量與1.8中獲取的AS樣本巨噬細(xì)胞含量的相對(duì)比例進(jìn)行相關(guān)性分析，并利用cor.test函數(shù)計(jì)算相關(guān)性系數(shù)。相關(guān)性系數(shù)大于0為正相關(guān)，相關(guān)性系數(shù)小于0為負(fù)相關(guān)，相關(guān)系數(shù)的絕對(duì)值代表強(qiáng)、弱或無(wú)相關(guān)，以P≤0.05為具有統(tǒng)計(jì)學(xué)意義。

2 結(jié) 果

2.1 scRNA-Seq數(shù)據(jù)質(zhì)控結(jié)果預(yù)處理GSE159677數(shù)據(jù)集中的樣本后，scRNA-Seq-對(duì)照組和scRNA-Seq-AS組樣本中共獲得2756個(gè)細(xì)胞，scRNA-Seq-AS組樣本中獲得2095個(gè)細(xì)胞。利用Seurat包計(jì)算所有樣本中細(xì)胞的基因數(shù)(nFeature)、轉(zhuǎn)錄本測(cè)序Count數(shù)(nCount)、線粒體基因百分比(percent.mt)，過(guò)濾掉每個(gè)細(xì)胞測(cè)到的基因數(shù)目小于50，線粒體基因組占比大于5%的細(xì)胞。對(duì)篩選細(xì)胞間標(biāo)準(zhǔn)化方差較大的1500個(gè)基因進(jìn)行后續(xù)分析，這些基因代表了細(xì)胞間的異質(zhì)性。見(jiàn)圖1。

a、b、c、d、：分別為AS樣本中細(xì)胞的基因數(shù)、轉(zhuǎn)錄本測(cè)序Count數(shù)、所有細(xì)胞的線粒體比例、細(xì)胞間標(biāo)準(zhǔn)方差較大的前1500個(gè)變異基因；e、f、g、h：分別為scRNA-Seq-對(duì)照樣本和scRNA-Seq-AS樣本中細(xì)胞的基因數(shù)、轉(zhuǎn)錄本測(cè)序Count數(shù)、所有細(xì)胞的線粒體比例、細(xì)胞間標(biāo)準(zhǔn)方差較大的前1500個(gè)變異基因

2.2降維聚類及注釋結(jié)果利用PCA法對(duì)過(guò)濾后的樣本進(jìn)行了線性降維，發(fā)現(xiàn)有20個(gè)主成分的P<0.05。然后，選取這20個(gè)主成分進(jìn)行t-SNE聚類分析，結(jié)果發(fā)現(xiàn)1927個(gè)細(xì)胞被分為了12個(gè)亞群，不同細(xì)胞亞群被標(biāo)記為了相應(yīng)的顏色；SingleR包注釋后12個(gè)細(xì)胞亞群被分為了6個(gè)亞群，包括內(nèi)皮細(xì)胞(3個(gè)細(xì)胞亞群)、巨噬細(xì)胞(3個(gè)細(xì)胞亞群)、單核細(xì)胞(3個(gè)細(xì)胞亞群)、T細(xì)胞(1個(gè)細(xì)胞亞群)、平滑肌細(xì)胞(1個(gè)細(xì)胞亞群)、組織干細(xì)胞(1個(gè)細(xì)胞亞群)。見(jiàn)圖2。

a：每個(gè)主成分的P值；b：各細(xì)胞亞群的聚類分布情況；c：注釋后各細(xì)胞亞群分布的t-SNE圖

2.3共同DEGs的篩選及其富集分析差異分析的結(jié)果顯示，巨噬細(xì)胞亞群和其它細(xì)胞亞群之間的有260個(gè)DEGs，其中上調(diào)88個(gè)，下調(diào)172個(gè)；組間有109個(gè)DEGs，其中上調(diào)35個(gè)，下調(diào)74個(gè)；兩者取交集后共獲得47個(gè)DEGs，其中上調(diào)27個(gè)，下調(diào)20個(gè)。利用clusterProfiler包、org.Hs.eg.db包等將篩選出的共同DEGs進(jìn)行功能富集分析和可視化，見(jiàn)圖3、圖4。GO功能富集分析結(jié)果顯示，依據(jù)P<0.05確定了434個(gè)GO條目；生物學(xué)過(guò)程(biological process，BP)條目為355條，主要涉及中性粒細(xì)胞激活、中性粒細(xì)胞脫顆粒、膠原蛋白分解代謝過(guò)程等；細(xì)胞組分條目(cellular component，CC)為43條，主要涉及含膠原蛋白的細(xì)胞外基質(zhì)、溶酶體腔等；分子功能(molecular function，MF)條目為36條，主要涉及鐵結(jié)合、鐵氧化酶活性、氧化還原酶活性等。KEGG通路富集分析結(jié)果顯示，依據(jù)P<0.05共映射出了6條信號(hào)通路，包括溶酶體、鐵死亡、凋亡、PPAR信號(hào)通路、自噬、補(bǔ)體與凝血級(jí)聯(lián)反應(yīng)。

圖 3 共同DEGs的GO富集分析

2.4特征hub基因的篩選及驗(yàn)證經(jīng)STRING數(shù)據(jù)庫(kù)對(duì)共同DEGs進(jìn)行PPI分析后，利用Cytoscape軟件構(gòu)建了PPI網(wǎng)絡(luò)。該網(wǎng)絡(luò)包括84條邊和38個(gè)節(jié)點(diǎn)，其中上調(diào)基因23個(gè)，下調(diào)基因15個(gè)。對(duì)cytoHubba插件中MNC、DMNC、MCC、EPC、Degree算法排序前10的基因取交集共獲得5個(gè)交集基因，即凝聚素(CLU)、組織蛋白酶D(CTSD)、組織蛋白酶B(CTSB)、組織蛋白酶Z(CTSZ)、組織蛋白酶L(CTSL)，其中CLU為下調(diào)基因，其余為上調(diào)基因。如圖5所示，與其他細(xì)胞亞群相比，CLU在巨噬細(xì)胞的表達(dá)和分布最小，而CTSD、CTSB、CTSZ、CTSL則相反。隨后，該研究利用芯片數(shù)據(jù)集GSE43292驗(yàn)證了5個(gè)特征hub基因在AS病理過(guò)程中的表達(dá)情況，結(jié)果與 scRNA-Seq數(shù)據(jù)分析的結(jié)果一致，見(jiàn)圖6；同時(shí)，5個(gè)特征hub基因表現(xiàn)出了良好的診斷AS的效能。見(jiàn)圖7。

圖 4 共同DEGs的KEGG通路富集分析

a：各細(xì)胞亞群中5個(gè)特征hub基因表達(dá)的氣泡圖；b、c、d、e、f：分別為CLU、CTSD、CTSB、CTSL、CTSZ在每個(gè)細(xì)胞分布的t-SNE圖

*P<0.001

2.5免疫細(xì)胞浸潤(rùn)分析CIBERSORT結(jié)果所示，對(duì)照和AS組織樣本之間存在8種差異的免疫細(xì)胞，其中M0巨噬細(xì)胞在AS中的比例顯著增加。ssGSEA結(jié)果所示，對(duì)照和AS組織樣本之間存在27種差異的免疫細(xì)胞，其中巨噬細(xì)胞在AS中的比例顯著增加。兩種算法的結(jié)果均提示，巨噬細(xì)胞在AS組織中的浸潤(rùn)比例較對(duì)照樣本高。見(jiàn)圖8。

圖 7 特征hub基因在GSE43292數(shù)據(jù)集中ROC分析的結(jié)果

a：CIBERSORT算法分析的組織樣本中22種免疫細(xì)胞浸潤(rùn)差異分析的小提琴圖；b：ssGSEA算法分析的組織樣本中29種免疫細(xì)胞浸潤(rùn)差異分析的小提琴圖

2.6特征hub基因與巨噬細(xì)胞的相關(guān)性分析CIBERSORT算法結(jié)果顯示，巨噬細(xì)胞含量的相對(duì)比例與CLU(r=-0.39,P=0.029)呈負(fù)相關(guān)性；與CTSD(r=0.83,P<0.001)、CTSB(r=0.76,P<0.001)、CTSL(r=0.85,P<0.001)、CTSZ(r=0.82,P<0.001)均表現(xiàn)出了明顯的強(qiáng)正相關(guān)性。ssGSEA算法結(jié)果顯示，巨噬細(xì)胞含量的相對(duì)比例與CLU(r=-0.51,P=0.002)呈一定的負(fù)相關(guān)性；與CTSD(r=0.87,P<0.001)、CTSB(r=0.91,P<0.001)、CTSL(r=0.86,P<0.001)、CTSZ(r=0.88,P<0.001)均表現(xiàn)出了明顯的強(qiáng)正相關(guān)性。

3 討 論

CAS作為AS最常見(jiàn)的臨床類型，可導(dǎo)致腦梗、心梗等不良事件[20-21]。研究顯示，巨噬細(xì)胞作為AS中一種豐富的功能異質(zhì)性免疫細(xì)胞，在CAS的病理過(guò)程中也同樣發(fā)揮著重要作用[20]。巨噬細(xì)胞通過(guò)調(diào)控脂質(zhì)代謝、炎癥反應(yīng)、LC3相關(guān)吞噬作用、免疫代謝、胞葬等多種途徑在AS中發(fā)揮重要作用，被認(rèn)為是治療AS的一個(gè)有潛力的作用靶點(diǎn)[1,3,22-24]。因此，深入挖掘巨噬細(xì)胞在CAS中潛在的分子機(jī)制和臨床診斷標(biāo)志物，對(duì)于有效安全的管理CAS具有重要的意義。

本研究對(duì)數(shù)據(jù)集GSE159677進(jìn)行分析，發(fā)現(xiàn)CAS樣本中巨噬細(xì)胞亞群有47個(gè)DEGs。富集分析顯示，這些DEGs可通過(guò)調(diào)控中性粒細(xì)胞的活化及相關(guān)免疫應(yīng)答、膠原蛋白分解及代謝、鐵結(jié)合及鐵氧化酶活性、細(xì)胞凋亡、鐵死亡、PPAR通路等生物學(xué)功能和通路，來(lái)介導(dǎo)CAS的病理過(guò)程。既往研究報(bào)道，中性粒細(xì)胞明膠酶相關(guān)脂質(zhì)沉積蛋白高表達(dá)于活化的中性粒細(xì)胞，可誘導(dǎo)巨噬細(xì)胞分泌促炎介質(zhì)加劇癥狀性CAS局部及全身性的炎癥反應(yīng)[25]。在CAS斑塊組織中，浸潤(rùn)的巨噬細(xì)胞通過(guò)促進(jìn)膠原蛋白分解、鐵死亡等，加速斑塊的不穩(wěn)定性[26-27]。

本文對(duì)共同DEGs進(jìn)行PPI分析，篩選出5個(gè)巨噬細(xì)胞特征hub基因即CLU、CTSD、CTSB、CTSZ、CTSL。同時(shí)，利用芯片數(shù)據(jù)集對(duì)hub基因進(jìn)行分析，結(jié)果發(fā)現(xiàn)hub基因在芯片數(shù)據(jù)和scRNA-Seq數(shù)據(jù)中的差異表達(dá)趨勢(shì)一致，并在診斷CAS方面表現(xiàn)出了良好的診斷效能。最后，進(jìn)行免疫細(xì)胞浸潤(rùn)分析，發(fā)現(xiàn)巨噬細(xì)胞在CAS組織中處于高浸潤(rùn)狀態(tài)。相關(guān)性分析顯示，在CAS環(huán)境中CLU的表達(dá)水平與巨噬細(xì)胞的相對(duì)含量存在一定的負(fù)相關(guān)性；而CTSD、CTSB、CTSL、CTSZ分別與巨噬細(xì)胞具有明顯的強(qiáng)正相關(guān)性。CLU又稱載脂蛋白J，通過(guò)調(diào)節(jié)巨噬細(xì)胞參與CAS病程變化的機(jī)制可能包含兩個(gè)方面，一是抑制巨噬細(xì)胞浸潤(rùn)和促炎性M1巨噬細(xì)胞的極化[28]；二是促進(jìn)巨噬細(xì)胞源泡沫細(xì)胞的膽固醇逆流[29]。巨噬細(xì)胞分泌的CTSB參與了CAS斑塊的進(jìn)展和破裂[30]。CTSL的表達(dá)隨著CAS的嚴(yán)重程度而增加，并與斑塊不穩(wěn)定和帽狀結(jié)構(gòu)破裂有關(guān)，而巨噬細(xì)胞是易損和破裂的頸動(dòng)脈斑塊中表達(dá)CTSL的主導(dǎo)細(xì)胞[31]。也就是說(shuō)，CTSL可能通過(guò)調(diào)控巨噬細(xì)胞凋亡，參與壞死核的形成和不穩(wěn)定斑塊的發(fā)展。研究顯示，CTSD可以促進(jìn)巨噬細(xì)胞泡沫化進(jìn)而加重AS病變[32]?，F(xiàn)階段，尚未有CTSD、CTSZ調(diào)控巨噬細(xì)胞的功能介導(dǎo)CAS病理過(guò)程的報(bào)道。綜上所述，這5個(gè)hub基因的高表達(dá)在CAS中發(fā)揮著重要作用，將有助于闡釋CAS可能的分子機(jī)制，并為CAS的診斷和防治提供一定的思路。

本研究尚有一定的局限性：①由于GEO等公共數(shù)據(jù)庫(kù)里人CAS樣本較少，未能對(duì)多數(shù)據(jù)集、多批次的數(shù)據(jù)進(jìn)行整合分析，導(dǎo)致納入分析的樣本有限。②本研究雖然利用外部基因芯片數(shù)據(jù)集對(duì)挖掘出的5個(gè)巨噬細(xì)胞特征基因進(jìn)行了驗(yàn)證，但結(jié)果仍需在后續(xù)的基礎(chǔ)和臨床研究中進(jìn)行佐證。

總之，本研究利用生物信息分析和scRNA-Seq數(shù)據(jù)相結(jié)合的方法從免疫學(xué)的角度探索了CAS可能的病理機(jī)制，挖掘出5個(gè)巨噬細(xì)胞特征基因(即CLU、CTSD、CTSB、CTSL、CTSZ)。本研究的核心意義不僅有助于更好地了解CAS的病理機(jī)制，為后續(xù)以巨噬細(xì)胞為切入點(diǎn)來(lái)探索CAS的免疫學(xué)機(jī)制提供了方向，還為CAS潛在臨床診斷標(biāo)志物的開發(fā)提供了一定的借鑒。