引用本文： 向潤, 謝少華, 廖瓊, 李強, 邵偉康, 李娟. 基于全外顯子組測序的胸腺癌基因組分析和生物標志物探索的回顧性隊列研究. 中國胸心血管外科臨床雜志, 2024, 31(2): 288-303. doi: 10.7507/1007-4848.202310011 復制
版權信息： ?四川大學華西醫院華西期刊社《中國胸心血管外科臨床雜志》版權所有，未經授權不得轉載、改編
胸腺是位于胸骨柄后方的前縱隔上部的次級淋巴腺。胸腺上皮腫瘤（thymic epithelial tumors，TET）是由胸腺上皮細胞分化而來的各種腫瘤，可分為胸腺瘤和胸腺癌（thymic carcinoma，TC）。TC是一種罕見且具有高度侵蝕性的惡性腫瘤[1]，整體預后較差。TC在 Ⅰ～Ⅱ期、Ⅲ期、Ⅳ期的5年生存率分別為88.2%、51.7%和37.6%[2-3]，而胸腺瘤在Ⅰ期、Ⅱ期、Ⅲ期、ⅣA期和ⅣB期的5年生存率分別為100.0%、98.4%、88.7%、70.6%和 52.8%[4]。
TC可與胸腺瘤的病理特征有顯著差異，它們具有惡性組織學特征以及不同的免疫組織化學和遺傳學特征[5-6]。然而，TC與轉移至胸腺且組織學表現相似的原發性肺惡性腫瘤之間較難區分，大大增加了診斷難度[7-8]。在分子技術層面，目前的研究[9]已經證實在所有成人惡性腫瘤中，胸腺腫瘤的腫瘤突變負荷（tumor mutation burden，TMB）最低，并且胸腺瘤和TC的分子畸變模式不同。但由于疾病罕見且數量有限，目前對TC分子特征的報道仍然較少。
在TC的治療和診斷中，可供參考的治療方法很少，根治性手術是TC的最佳治療選擇，鉑類化療是無法手術患者的標準治療方法，但晚期患者的反應期通常很短。鉑類化療失敗后可用的治療方法很少[10]。這是因為對這類腫瘤的生物學認識不足。因此，有必要探索其他新治療方法或其他生物標志物，幫助延長患者生存期，改善預后。
最近，二代測序（next-generation sequencing，NGS）技術作為一種快速和具有成本-效益的手段應用于描述個體患者基因組中的突變模式，揭示與癌癥進展有關的突變癌癥基因[7]。很大程度上，NGS了解疾病的分子基礎是一種有效可行的手段，既往研究也逐漸將NGS技術應用于胸腺上皮腫瘤（TET)的診斷和治療。例如，Shimada等[11]利用NGS-DNA技術靶向測序53個基因和聚合酶鏈式反應（PCR）技術分析TET的基因組圖譜，確定RAS突變可能是TET治療的候選基因靶點。Szpechcinski 等[12]使用NGS-DNA技術靶向檢測了15個實體瘤常見突變基因，共納入34個TC和19個胸腺瘤樣本，發現TC表現出比胸腺瘤更大的遺傳失調，并展示了兩種腫瘤之間的基因組異質性。Tsukaguchi等[13]通過NGS鑒定攜帶PIK3CA突變的黏液性TC患者，使用lenvatinib治療失敗，從而推斷PIK3CA可能與lenvatinib耐藥相關。但是TC作為一種罕見腫瘤，其分子驅動因素在很大程度上仍然未知。目前相關研究對象大多為TET或胸腺瘤，僅單獨針對TC的研究隊列極少。其次，既往研究基因檢測方式大多為靶向NGS檢測部分基因的突變信息，總體納入研究的基因數較少，在全面生物標志物的研究方面存在欠缺，目前亟需大隊列TC的全面基因組研究以獲取更全面的分子突變信息。
綜上，目前有必要對TC患者進行更全面的基因測序研究。因此，我們利用全外顯子組測序（whole exome sequencing，WES）技術檢測24例TC患者基因的體細胞突變，試圖通過闡明TC的基因突變譜，更深入地了解這種罕見腫瘤的分子遺傳信息，確定參與TC發展的致癌基因和對TC患者的個體化靶向治療，并提供更精準的治療方案。
1 資料與方法
1.1 樣本采集
本研究回顧性納入2015年1月—2021年2月四川省腫瘤醫院胸外科、四川省癌癥中心收治的24例TC患者。進行NGS檢測的樣本為甲醛固定后石蠟包埋（formalin fixed paraffin embedded，FFPE）的腫瘤組織，患者匹配的外周血白細胞作為胚系DNA對照，并于2021年6月開始WES檢測。
1.2 全外顯子建庫測序
用AllPrep DNA/RNA mini Kit（Qiagen 80204）試劑盒從血液淋巴細胞中提取基因組DNA（gDNA），gDNA使用KAPA Hyper Prep kit（Kapa Biosystems）按照說明書進行文庫構建。用MagMaxTM Cell-Free DNA Isolation Kit（Thermo fisher A29319）試劑盒提取ctDNA。基于平鋪方式對全外顯子基因組進行設計，gDNA文庫使用安捷倫定制的全外顯子捕獲探針對目標區域進行了富集。將測序樣品上機illumina Hiseq X10平臺進行高通量測序，測序策略為PairEnd150。腫瘤組織測序數據量35 g，要求測序深度≥200 倍，白細胞對照測序數據量15 g，要求測序深度≥100 倍。
1.3 序列比對
使用fastp[14]對原始數據進行質量控制。高質量的雙端讀數通過BWA-MEM[15]與人類參考基因組hg19進行對比。分別使用Samtools （
1.4 單核苷酸替換變異檢測
基于配對的排列文件（腫瘤和匹配的胚系），使用MuTect2 （
利用maftools[19]中的somaticInteraction函數對各基因中體細胞突變進行相關性分析。首先構建一個線性模型，然后再對相關性水平進行估計，最后篩選相關性最高的前25個基因進行展示。
1.5 腫瘤突變負荷計算
為了確定TMB值，我們對NGS檢測到的體細胞非同義SNVs的數量進行量化，并使用可靠的算法將該值推算至整個外顯子組。然后，我們在得到腫瘤樣本的絕對突變數與正常樣本的突變點后，用以下公式計算腫瘤樣本的TMB：絕對突變數×1 000 000/外顯子堿基數。簡言之，TMB以每Mb的突變數衡量[20]。
1.6 微衛星不穩定狀態分析
微衛星不穩定（microsatellite instability，MSI）是根據微衛星長度的變化可將其分為穩定型（MSS：長度沒有變化）及不穩定型（MSI：長度發生變化）。用MSIsensor軟件進行微衛星狀態的檢測，軟件算法原理為：首先，對于在腫瘤和正常樣本中測序深度都≥20倍的微衛星位點，統計其每種重復長度的reads數目分布情況；然后使用χ2檢驗對微衛星位點上的分布進行統計檢驗，若差異有統計學意義，則認為該位點不穩定；最后，統計不穩定位點的比例（MSI評分）并按照MSI評分>10判定為MSI狀態（微衛星不穩定），MSI評分≤10判定為MSS狀態（微衛星穩定）[21]。
1.7 突變等位基因腫瘤異質性計算
變異等位基因頻率（variant allele frequency，VAF）是每個位置的交替等位基因觀測值與讀數深度的比率。為了包括VAF在0.05～1.00的所有體細胞變異，我們修改了突變等位基因腫瘤異質性（mutant-allele tumor heterogeneity，MATH）得分，計算方法是100×中位絕對偏差/VAF的中位數[22]。
1.8 人類白細胞抗原分型和腫瘤新抗原負荷計算
對于每個樣本，使用Optitype推斷出4位數的人類白細胞抗原（human leukocyte antigen，HLA）類型，輸入文件是正常胚系樣本的bam文件[23]。腫瘤新抗原負荷（tumor neoantigen burden，TNB）是指在一個基因組區域中每堿基包含的新抗原數量[24]。pVACseq過濾結果中肽的數量記錄為新抗原的數量，TNB是每兆堿基中新抗原數量，絕對TNB是每兆中新抗原對應的突變數量。
1.9 腫瘤純度、倍性和同源重組缺陷分析
用Seqenza軟件包對腫瘤樣本的純度和倍性進行推斷。其基本原理如下：Sequenza基于基因組測序reads的分割數據構建概率模型。通過計算腫瘤和正常配對樣本的平均深度比和B等位基因頻率（在種系雜合位置測量的兩個等位基因部分中較小者）。在此模型下給出純度和倍性的估計值[25]。
通過scarHRD軟件包確定WES數據中同源重組缺陷（homologous recombination deficiency，HRD）水平，包括端粒等位基因不平衡（telomeric allelic imbalance，TAI）、雜合性缺失（loss of heterozygosity，LOH）、大片段遷移（large-scale state transitions，LST）的數量。
1.10 富集分析和統計學分析
富集分析使用R（v. 3.6.1）中的clusterProfiler（v. 3.18.1）[26]包進行統計分析和繪圖。基因組百科全書（KEGG）富集分析采用Ingenuity Pathway Analysis進行分析。統計分析軟件采用R軟件。連續變量組間差異檢驗采用Wilcoxon秩和或Kruskal-Wallis檢驗。分類變量組間差異檢驗采用χ2檢驗用。樣本量較小的組間差異檢驗用Student’s t檢驗。總生存期（OS）為從隨機化開始到死亡（任何原因）的時間。生存分析使用Cox比例風險模型估計風險比（HR）及其相應的95%置信區間（CI）。若未明確說明， P≤0.05為差異有統計學意義。
1.11 倫理審查
本研究嚴格按照《赫爾辛基宣言》的要求進行，經四川省腫瘤醫院醫學倫理委員會審查批準（倫理審批號：SCCHEC-02-2021-036）。所有患者均在入組前簽署書面知情同意書。
2 結果
最終納入24例TC患者，其中男16例、女8例，中位年齡55（42～74）歲。
2.1 胸腺癌患者的基因突變圖譜
我們分析了24例患者的SNV數據（圖1a）。在所有樣本中TTN基因的突變頻率最高（42%），其次是HSPG2（29%）、OBSCN（29%）和TSSC1（25%）。通過對所有SNV突變變體的分析，發現錯義突變是最常見的突變分類（圖1b）。單核苷酸多態性（SNP）是變異類型中最常見的類型。在SNV大類中，C-T的發生頻率最高。樣本K004002T的變異最多，為3 607個，每個樣本發生變異次數的中位數為72.5。突變頻率最高的前10個基因是TTN、OBSCN、HSPG2、DNAH5、MEGF8、KIF21B、CAD、PLEC、KIAA0100、DYNC1H1。另外，所有患者的顛換、轉換突變比例中，C-T占所有突變類型的比例最高（圖1c），同時轉換的比例遠高于顛換。最后，對所有樣本的各基因中體細胞突變進行相關性分析，并對相關性最高的前25個基因進行展示（圖1d）。UBR4、KALRN與其他基因發生共突變事件的比例較高。
圖1
胸腺癌患者基因突變情況總覽
a：TC患者的單核苷酸突變的整體突變情況；b：突變基因的類型分析；c：所有突變類型的比例、顛換及轉換的比例；d：體細胞突變的相互作用關系
對所有樣本突變分析結果（圖2a）顯示，發生染色體（CNV）突變頻率最高的是ZNF276基因（54%，loss），其次是BEND3（50%，loss）、DHODH（50%，loss）和VAC14（50%，loss）。大多數發生CNV的類型均為loss, 僅有NBPF10（33%）的突變類型為gain。同時對染色體的位置與CNV的擴增或缺失情況進行評估），2q32.2位置的G-score最大，2q11.1的G-score最小（圖2b）。
圖2
胸腺癌患者染色體拷貝數變異總覽和富集分析
a：胸腺癌患者的CNV整體突變情況；b：染色體的位置與CNV的擴增或缺失情況進行評估；c：KEGG通路富集分析；A：心律失常性右室心肌病；B：ECM-受體的相互作用；C：病灶粘連；D：擴張型心肌病；E：ABC轉運蛋白；F：胰島素分泌；G：肥厚型心肌病；H：磷脂酰肌醇信號系統；I：谷氨酸能突觸；J：肌醇磷酸代謝；K：催產素信號通路；L：胃酸分泌；M：醛固酮的合成和分泌；N：通路為軸突導向；O：cGMP?PKG信號；P：GABA能突觸；Q：肌動蛋白細胞骨架的調控；R：心肌細胞的腎上腺素能信號傳導；S：溶酶體；T：蛋白質的消化和吸收；CNV：染色體拷貝數變異；KEGG：基因組百科全書
KEGG通路分析顯示，發生SNV突變的基因富集到的通路。根據矯正后 P 值對結果進行排序，最顯著的參與通路是心律失常性右室心肌病、細胞外基質（ECM）-受體的相互作用、病灶粘連和擴張型心肌病等（圖2c）。
2.2 TMB、MATH和MSI
MSI分析結果顯示，6例患者為MSI，占患者總數的25%，其余患者MSI分數均相對較低，為MSS型。評估MSI狀態有助于確定臨床預后，并可能預測腫瘤對化學治療的反應。TMB結果顯示，大多數患者的TMB水平較低，只有3例患者的TMB值較高，中位TMB為2.045，平均TMB為9.86。晚期克隆擴增、空間隔離和不完全選擇性掃描導致遺傳上不同的細胞群都會出現腫瘤內異質性，MATH評分可以反映這種腫瘤內異質性。絕大多數患者的MATH值很高，只有少數患者的MATH值<20（表1）。
表1
胸腺癌隊列MSI、TMB、MATH的結果
2.3 HLA分型和TNB分析
HLA分析結果表示，大多數樣本都有雜合的HLA等位基因。所有樣本中共檢測到14個HLA-A位點、24個HLA-B位點和15個HLA-C位點。HLA-A*11:01、B*46:01和C*01:02是TC的易感基因（表2）。
表2
胸腺癌隊列HLA分析結果
TNB結果顯示絕大多數患者的TNB值很低，只有1例患者的TNB值較高（表3）。
表3
胸腺癌隊列新抗原負荷結果
2.4 同源重組缺陷、腫瘤純度和倍性
腫瘤組織樣本包含癌細胞和正常細胞的混合物。因此，腫瘤數據分析必須考慮標本的腫瘤純度和倍性，以便于在其他分析中剔除這些因素帶來的影響。大多數樣本的倍性為2.0，純度大于0.8（表4）。
表4
胸腺癌隊列的腫瘤純度和倍性結果
HRD-LOH是同源重組缺陷相關的LOH，全部為雜合性缺失，其中樣本K004002T和K034958T在HRD-LOH、LST、TAI的數量均為最多，其余樣本的個數均較少（表5）。
表5
胸腺癌隊列HRD-LOH的結果
2.5 分組分析和預后分析
為了更好地探索臨床特征與基因突變的關系，將患者分別按照性別、年齡、是否吸煙、病理分型、Masaoka分期、治療方案、是否生存、腫瘤直徑、是否復發分成9個組，分析每組樣本的數量。將MATH、HRD、CNV、純度、倍性等生物標志物進行組間統計學分析（表6）。
表6
本研究的分組方式及對應的樣本數統計
CNV突變負荷結果在多個分組內均有顯著差異。3組吸煙者的CNV負荷（P=0.026）和CNV-loss負荷（P=0.018）均高于不吸煙者（圖3a）。6組手術患者和放化療患者的CNV負荷（P=0.0079）、CNV-gain負荷（P=0.04）、CNV-loss負荷（P=0.0078）均顯著高于新輔助放化療患者（圖3b）。7組生存患者的CNV負荷（P=0.036）和CNV-loss負荷（P=0.018）顯著低于死亡患者（圖3c），死亡患者的CNV-gain負荷（P=0.44）有高于生存患者的趨勢。9組是否復發（圖3d）與7組是否死亡有相似的結論（P=0.0092，P=0.0033，P=0.16）。
圖3
胸腺癌患者拷貝數變異負荷在不同分組的差異
a：CNV負荷、CNV gain負荷和CNV loss負荷分別在3組是否吸煙組； 6組不用治療方案組（b）、7組是否生存組（c）和9組是否復發組（d）的統計分析；CNV：染色體拷貝數變異
MATH 第6組不同治療方案組差異有統計學意義（P=0.041，圖4a），新輔助放化療組的MATH值顯著低于手術組。第8組腫瘤直徑50～100 mm患者的MATH值顯著低于腫瘤直徑<50 mm患者和>100 mm患者（P=0.033）。腫瘤的倍性（圖4b）在1組性別方面差異有統計學意義，女性樣本的倍性顯著高于男性樣本（P=0.027），但是不同性別樣本的腫瘤純度差異無統計學意義（P=0.98）。HRD結果在7組是否生存和9組是否復發的差異有統計學意義。7組死亡患者的HRD-LOH數量（P=0.033）、HRD總數（P=0.024）、LST數量（P=0.026）顯著高于生存患者（圖4c）。9組復發患者的HRD-LOH數量（P=0.046）、HRD總數（P=0.016）、LST數量（P=0.017）顯著高于非復發患者（圖4d）。但是兩者TAI差異均無統計學意義（P=0.92）。
圖4
胸腺癌患者多個生物標志物在組間的差異
a：MATH分別在6組不同治療方案患者和8組腫瘤大小的統計分析；b：腫瘤倍性和細胞性分別在1組性別的統計分析中；c：HRD-LOH、HRD總數、LST數量和TAI數量分別在7組是否生存的統計分析；d：HRD-LOH、HRD總數、LST數量和TAI數量分別在9組是否復發的統計分析；e：在4組不同病理患者間的預后分析；f：在9組是否復發的預后分析；MATH：突變等位基因腫瘤異質性；HRD：同源重組缺陷；LOH：雜合性缺失；TAI：端粒等位基因不平衡； LST：大片段遷移
通過Cox回歸進行單因素預后分析，分別計算所有臨床特征和患者OS之間的關系。其中，TC的不同病理分型（圖4e）和是否復發（圖4f）是與OS顯著相關的因素。胸腺神經內分泌癌的OS最短，胸腺鱗狀細胞癌的OS相對最長（P=0.0075）。未復發患者的OS顯著長于復發患者（P<0.001）。
2.6 TCGA數據庫對比
通過本隊列TC數據與TCGA數據庫中胸腺瘤數據的對比，得到兩癌種突變圖譜的差異和相似點。胸腺瘤和TC數據的SNV突變圖譜表明（圖5a），突變頻率最高的前4個基因為GTF2I（41%）、TTN（11%）、HRAS（8%）和HSPG2（6%）。僅在胸腺瘤隊列中，GTF2I基因的突變頻率最高（49%），其次是HRAS（8%）、MUC16（5%）和TTN（5%），其中GTF2I僅在胸腺瘤數據中檢測到。TC樣本的SNV平均突變數高于胸腺瘤樣本。
圖5
胸腺癌患者多數據集的突變圖譜對比
a：本隊列TC數據與TCGA胸腺瘤數據的SNV突變情況對比；b：本隊列TC數據與TCGA胸腺瘤數據的CNV突變情況對比；SNV：單核苷酸突變；CNV：染色體拷貝數變異
胸腺瘤和TC數據的CNV突變圖譜表明，突變頻率最高的前4個基因為RIPK1（24%）、AKAP12（23%）、BEND3（22%）和NUP43（22%）。其中，HLA-E 和DUS2基因僅在TC數據中被檢測到突變。其次在TC數據中僅有NBPF10基因是CNV-gain型的突變，但在胸腺瘤數據中，發生CNV-gain突變的基因數量遠高于TC隊列。這表明胸腺瘤和TC的分子畸變模式不同（圖5b）。
3 討論
我們的研究結果顯示，TTN是突變頻率最高的基因（42%），其次是HSPG2（29%）和OBSCN（29%），這早期研究[27-28]結果相似。TTN是目前最大的蛋白，由364個外顯子組成，位于染色體2q31上[29]。TTN突變在人類癌癥中起重要作用。研究[30]證實，TTN基因突變是甲狀腺癌的獨立危險因素［HR=4.558，95%CI（1.808，11.494），P=0.001］。TTN可與lncRNA TTN-AS1結合，在乳腺癌、肺癌、消化系統腫瘤等許多癌癥中具有促癌作用[31]。更重要的是，在結直腸癌中，TTN/OBSCN“雙重打擊”是一種“免疫熱”亞型，具有潛在更好的免疫治療效果，是良好預后的預測因子[32-33]。在本研究中，TTN和OBSCN的突變頻率極高，可能是TC的危險因素，并對患者免疫治療起重要作用。
與大多數成人腫瘤類型相比，TET的TMB相對較低，且TC的TMB水平遠高于胸腺瘤[34-36]。既往研究中，胸腺瘤的平均TMB為0.663[37]。本隊列中TC平均TMB為9.86，顯著高于文獻[34-36]報道的胸腺瘤的平均TMB，同時本研究中TC隊列的突變數遠遠高于胸腺瘤隊列（圖5a）。臨床研究[38]表明，高TMB患者有免疫治療的機會，且對免疫治療的反應明顯優于低TMB患者。這是因為高TMB患者可能含有更多激活免疫細胞的新抗原，因此在使用免疫檢查點抑制劑（ICIs）治療后具有更長的生存期和更高的反應率[39]。2017年，新英格蘭醫學雜志發表了27種腫瘤的中位TMB與ICIs治療的客觀緩解率（ORR）之間的線性圖，發現TMB與ORR存在顯著的正相關性（P<0.001），相關系數為0.74，表明在這27種腫瘤中TMB能夠預測ICIs的療效[40]。2020年6月，帕博利珠單抗第二個不限癌種的新適應證獲得美國食品藥品監督管理局（FDA）批準，用于治療具有高組織腫瘤突變負荷（tTMB-H）無法切除或轉移性實體瘤患者。根據TC具有較高TMB的特點，未來或許可以從免疫治療中獲益。綜上，TC患者進行TMB檢測是非常必要的。但是并非所有的突變都是具有免疫原性的，后續對TC免疫細胞的組成和功能的深入了解對有限管理TC的進展和免疫反應至關重要[35]。
MSI-H在泛癌種中占比為0%～16.5%。在子宮內膜癌中MSI占比高達16.5%，但是在黑色素瘤（0%）、肺癌（0.36%）和很多血液腫瘤中MSI-H是罕見的[41]。2017年5月，美國FDA批準帕博利珠單抗（pembrolizumab）用于無論腫瘤來源，微衛星不穩定性高（MSI-H）的實體腫瘤患者。曾有文獻[42]報道1例MSI型TC患者接受阿維單抗（Avelumab）和阿昔替尼（Axitinib）聯合治療后，獲得持久的部分反應。在我們的研究中，24例患者中有6例是MSI型，占總人數的25%，與其他癌種相比TC的MSI型患者占比相對較高，因此從免疫治療中獲益的可能性很高，同時也表明TC患者可能需要進行MSI檢測來指導用藥。
MATH值 是一個腫瘤異質性的定量指標。以前的研究[43-44]表明顯示，高MATH與較低的總生存率有關。 不幸的是，在本研究中，絕大多數患者的MATH值較高，只有少數患者的MATH值<20。這表明TC預后相對較差。
HRD檢測是通過檢測腫瘤細胞基因組中與DNA修復機制有關的基因變異情況，來評估細胞中HRD的程度，從而判斷腫瘤的基因組穩定性[45-46]。目前已經在很多癌種中證實，HRD與腫瘤患者的預后密切相關，HRD陽性腫瘤可能預示患者預后更差。LOH現象在癌癥中普遍存在，已有研究[47]表明LOH導致抑癌基因失活，從而影響癌癥的發生和進展。該結論在本研究中也得到證實，復發患者HRD數量顯著高于未復發患者（圖4d），且復發患者OS顯著短于未復發患者（圖4f）。同理，死亡患者的HRD數量顯著高于生存患者（圖4c）。
盡管本研究初步取得了一些重要的發現，但樣本量相對較小，一定程度上限制了我們對結果的解釋和及其可推廣性。因此，我們計劃在后續研究中擴大樣本量，以加強和驗證我們的研究結果，更全面地闡述研究問題，增加結論的可信度。
總之，本研究通過NGS WES子技術檢測了TC的體細胞突變，并利用生物信息學方法進一步分析了多維度的生物標志物。這項研究提供了迄今為止第一個關于中國TC突變情況的全面分析。本隊列中，在SNV突變層面，TTN基因的突變頻率最高（42%），其次是HSPG2（29%）、OBSCN（29%）和TSSC1（25%）。發生CNV突變頻率最高的是ZNF276基因（54%，Loss），其次是BEND3（50%，Loss）、DHODH（50%，Loss）和VAC14（50%，Loss）。通過研究免疫治療的可能性發現，TC患者的TMB相比胸腺瘤較高，且MSI型患者的占比較高，表明TC患者有機會接受免疫治療，并有較長的有效生存期。這些發現有助于更好地了解TC的分子發病機理，并可能為治療TC提供指導。我們將在未來的研究中收集更多的臨床樣本并驗證這些結果。
利益沖突：邵偉康受雇于無錫臻和生物科技有限公司；其他作者聲明無潛在利益沖突。
作者貢獻：向潤負責文獻檢索，論文初稿撰寫；謝少華、廖瓊負責論文初稿審閱、對文章知識性內容作批判性審閱；邵偉康負責數據分析，初稿撰寫和校對；李娟、李強負責論文設計指導，參與制定研究思路，全文審校。
胸腺是位于胸骨柄后方的前縱隔上部的次級淋巴腺。胸腺上皮腫瘤（thymic epithelial tumors，TET）是由胸腺上皮細胞分化而來的各種腫瘤，可分為胸腺瘤和胸腺癌（thymic carcinoma，TC）。TC是一種罕見且具有高度侵蝕性的惡性腫瘤[1]，整體預后較差。TC在 Ⅰ～Ⅱ期、Ⅲ期、Ⅳ期的5年生存率分別為88.2%、51.7%和37.6%[2-3]，而胸腺瘤在Ⅰ期、Ⅱ期、Ⅲ期、ⅣA期和ⅣB期的5年生存率分別為100.0%、98.4%、88.7%、70.6%和 52.8%[4]。
TC可與胸腺瘤的病理特征有顯著差異，它們具有惡性組織學特征以及不同的免疫組織化學和遺傳學特征[5-6]。然而，TC與轉移至胸腺且組織學表現相似的原發性肺惡性腫瘤之間較難區分，大大增加了診斷難度[7-8]。在分子技術層面，目前的研究[9]已經證實在所有成人惡性腫瘤中，胸腺腫瘤的腫瘤突變負荷（tumor mutation burden，TMB）最低，并且胸腺瘤和TC的分子畸變模式不同。但由于疾病罕見且數量有限，目前對TC分子特征的報道仍然較少。
在TC的治療和診斷中，可供參考的治療方法很少，根治性手術是TC的最佳治療選擇，鉑類化療是無法手術患者的標準治療方法，但晚期患者的反應期通常很短。鉑類化療失敗后可用的治療方法很少[10]。這是因為對這類腫瘤的生物學認識不足。因此，有必要探索其他新治療方法或其他生物標志物，幫助延長患者生存期，改善預后。
最近，二代測序（next-generation sequencing，NGS）技術作為一種快速和具有成本-效益的手段應用于描述個體患者基因組中的突變模式，揭示與癌癥進展有關的突變癌癥基因[7]。很大程度上，NGS了解疾病的分子基礎是一種有效可行的手段，既往研究也逐漸將NGS技術應用于胸腺上皮腫瘤（TET)的診斷和治療。例如，Shimada等[11]利用NGS-DNA技術靶向測序53個基因和聚合酶鏈式反應（PCR）技術分析TET的基因組圖譜，確定RAS突變可能是TET治療的候選基因靶點。Szpechcinski 等[12]使用NGS-DNA技術靶向檢測了15個實體瘤常見突變基因，共納入34個TC和19個胸腺瘤樣本，發現TC表現出比胸腺瘤更大的遺傳失調，并展示了兩種腫瘤之間的基因組異質性。Tsukaguchi等[13]通過NGS鑒定攜帶PIK3CA突變的黏液性TC患者，使用lenvatinib治療失敗，從而推斷PIK3CA可能與lenvatinib耐藥相關。但是TC作為一種罕見腫瘤，其分子驅動因素在很大程度上仍然未知。目前相關研究對象大多為TET或胸腺瘤，僅單獨針對TC的研究隊列極少。其次，既往研究基因檢測方式大多為靶向NGS檢測部分基因的突變信息，總體納入研究的基因數較少，在全面生物標志物的研究方面存在欠缺，目前亟需大隊列TC的全面基因組研究以獲取更全面的分子突變信息。
綜上，目前有必要對TC患者進行更全面的基因測序研究。因此，我們利用全外顯子組測序（whole exome sequencing，WES）技術檢測24例TC患者基因的體細胞突變，試圖通過闡明TC的基因突變譜，更深入地了解這種罕見腫瘤的分子遺傳信息，確定參與TC發展的致癌基因和對TC患者的個體化靶向治療，并提供更精準的治療方案。
1 資料與方法
1.1 樣本采集
本研究回顧性納入2015年1月—2021年2月四川省腫瘤醫院胸外科、四川省癌癥中心收治的24例TC患者。進行NGS檢測的樣本為甲醛固定后石蠟包埋（formalin fixed paraffin embedded，FFPE）的腫瘤組織，患者匹配的外周血白細胞作為胚系DNA對照，并于2021年6月開始WES檢測。
1.2 全外顯子建庫測序
用AllPrep DNA/RNA mini Kit（Qiagen 80204）試劑盒從血液淋巴細胞中提取基因組DNA（gDNA），gDNA使用KAPA Hyper Prep kit（Kapa Biosystems）按照說明書進行文庫構建。用MagMaxTM Cell-Free DNA Isolation Kit（Thermo fisher A29319）試劑盒提取ctDNA。基于平鋪方式對全外顯子基因組進行設計，gDNA文庫使用安捷倫定制的全外顯子捕獲探針對目標區域進行了富集。將測序樣品上機illumina Hiseq X10平臺進行高通量測序，測序策略為PairEnd150。腫瘤組織測序數據量35 g，要求測序深度≥200 倍，白細胞對照測序數據量15 g，要求測序深度≥100 倍。
1.3 序列比對
使用fastp[14]對原始數據進行質量控制。高質量的雙端讀數通過BWA-MEM[15]與人類參考基因組hg19進行對比。分別使用Samtools （
1.4 單核苷酸替換變異檢測
基于配對的排列文件（腫瘤和匹配的胚系），使用MuTect2 （
利用maftools[19]中的somaticInteraction函數對各基因中體細胞突變進行相關性分析。首先構建一個線性模型，然后再對相關性水平進行估計，最后篩選相關性最高的前25個基因進行展示。
1.5 腫瘤突變負荷計算
為了確定TMB值，我們對NGS檢測到的體細胞非同義SNVs的數量進行量化，并使用可靠的算法將該值推算至整個外顯子組。然后，我們在得到腫瘤樣本的絕對突變數與正常樣本的突變點后，用以下公式計算腫瘤樣本的TMB：絕對突變數×1 000 000/外顯子堿基數。簡言之，TMB以每Mb的突變數衡量[20]。
1.6 微衛星不穩定狀態分析
微衛星不穩定（microsatellite instability，MSI）是根據微衛星長度的變化可將其分為穩定型（MSS：長度沒有變化）及不穩定型（MSI：長度發生變化）。用MSIsensor軟件進行微衛星狀態的檢測，軟件算法原理為：首先，對于在腫瘤和正常樣本中測序深度都≥20倍的微衛星位點，統計其每種重復長度的reads數目分布情況；然后使用χ2檢驗對微衛星位點上的分布進行統計檢驗，若差異有統計學意義，則認為該位點不穩定；最后，統計不穩定位點的比例（MSI評分）并按照MSI評分>10判定為MSI狀態（微衛星不穩定），MSI評分≤10判定為MSS狀態（微衛星穩定）[21]。
1.7 突變等位基因腫瘤異質性計算
變異等位基因頻率（variant allele frequency，VAF）是每個位置的交替等位基因觀測值與讀數深度的比率。為了包括VAF在0.05～1.00的所有體細胞變異，我們修改了突變等位基因腫瘤異質性（mutant-allele tumor heterogeneity，MATH）得分，計算方法是100×中位絕對偏差/VAF的中位數[22]。
1.8 人類白細胞抗原分型和腫瘤新抗原負荷計算
對于每個樣本，使用Optitype推斷出4位數的人類白細胞抗原（human leukocyte antigen，HLA）類型，輸入文件是正常胚系樣本的bam文件[23]。腫瘤新抗原負荷（tumor neoantigen burden，TNB）是指在一個基因組區域中每堿基包含的新抗原數量[24]。pVACseq過濾結果中肽的數量記錄為新抗原的數量，TNB是每兆堿基中新抗原數量，絕對TNB是每兆中新抗原對應的突變數量。
1.9 腫瘤純度、倍性和同源重組缺陷分析
用Seqenza軟件包對腫瘤樣本的純度和倍性進行推斷。其基本原理如下：Sequenza基于基因組測序reads的分割數據構建概率模型。通過計算腫瘤和正常配對樣本的平均深度比和B等位基因頻率（在種系雜合位置測量的兩個等位基因部分中較小者）。在此模型下給出純度和倍性的估計值[25]。
通過scarHRD軟件包確定WES數據中同源重組缺陷（homologous recombination deficiency，HRD）水平，包括端粒等位基因不平衡（telomeric allelic imbalance，TAI）、雜合性缺失（loss of heterozygosity，LOH）、大片段遷移（large-scale state transitions，LST）的數量。
1.10 富集分析和統計學分析
富集分析使用R（v. 3.6.1）中的clusterProfiler（v. 3.18.1）[26]包進行統計分析和繪圖。基因組百科全書（KEGG）富集分析采用Ingenuity Pathway Analysis進行分析。統計分析軟件采用R軟件。連續變量組間差異檢驗采用Wilcoxon秩和或Kruskal-Wallis檢驗。分類變量組間差異檢驗采用χ2檢驗用。樣本量較小的組間差異檢驗用Student’s t檢驗。總生存期（OS）為從隨機化開始到死亡（任何原因）的時間。生存分析使用Cox比例風險模型估計風險比（HR）及其相應的95%置信區間（CI）。若未明確說明， P≤0.05為差異有統計學意義。
1.11 倫理審查
本研究嚴格按照《赫爾辛基宣言》的要求進行，經四川省腫瘤醫院醫學倫理委員會審查批準（倫理審批號：SCCHEC-02-2021-036）。所有患者均在入組前簽署書面知情同意書。
2 結果
最終納入24例TC患者，其中男16例、女8例，中位年齡55（42～74）歲。
2.1 胸腺癌患者的基因突變圖譜
我們分析了24例患者的SNV數據（圖1a）。在所有樣本中TTN基因的突變頻率最高（42%），其次是HSPG2（29%）、OBSCN（29%）和TSSC1（25%）。通過對所有SNV突變變體的分析，發現錯義突變是最常見的突變分類（圖1b）。單核苷酸多態性（SNP）是變異類型中最常見的類型。在SNV大類中，C-T的發生頻率最高。樣本K004002T的變異最多，為3 607個，每個樣本發生變異次數的中位數為72.5。突變頻率最高的前10個基因是TTN、OBSCN、HSPG2、DNAH5、MEGF8、KIF21B、CAD、PLEC、KIAA0100、DYNC1H1。另外，所有患者的顛換、轉換突變比例中，C-T占所有突變類型的比例最高（圖1c），同時轉換的比例遠高于顛換。最后，對所有樣本的各基因中體細胞突變進行相關性分析，并對相關性最高的前25個基因進行展示（圖1d）。UBR4、KALRN與其他基因發生共突變事件的比例較高。
圖1
胸腺癌患者基因突變情況總覽
a：TC患者的單核苷酸突變的整體突變情況；b：突變基因的類型分析；c：所有突變類型的比例、顛換及轉換的比例；d：體細胞突變的相互作用關系
對所有樣本突變分析結果（圖2a）顯示，發生染色體（CNV）突變頻率最高的是ZNF276基因（54%，loss），其次是BEND3（50%，loss）、DHODH（50%，loss）和VAC14（50%，loss）。大多數發生CNV的類型均為loss, 僅有NBPF10（33%）的突變類型為gain。同時對染色體的位置與CNV的擴增或缺失情況進行評估），2q32.2位置的G-score最大，2q11.1的G-score最小（圖2b）。
圖2
胸腺癌患者染色體拷貝數變異總覽和富集分析
a：胸腺癌患者的CNV整體突變情況；b：染色體的位置與CNV的擴增或缺失情況進行評估；c：KEGG通路富集分析；A：心律失常性右室心肌病；B：ECM-受體的相互作用；C：病灶粘連；D：擴張型心肌病；E：ABC轉運蛋白；F：胰島素分泌；G：肥厚型心肌病；H：磷脂酰肌醇信號系統；I：谷氨酸能突觸；J：肌醇磷酸代謝；K：催產素信號通路；L：胃酸分泌；M：醛固酮的合成和分泌；N：通路為軸突導向；O：cGMP?PKG信號；P：GABA能突觸；Q：肌動蛋白細胞骨架的調控；R：心肌細胞的腎上腺素能信號傳導；S：溶酶體；T：蛋白質的消化和吸收；CNV：染色體拷貝數變異；KEGG：基因組百科全書
KEGG通路分析顯示，發生SNV突變的基因富集到的通路。根據矯正后 P 值對結果進行排序，最顯著的參與通路是心律失常性右室心肌病、細胞外基質（ECM）-受體的相互作用、病灶粘連和擴張型心肌病等（圖2c）。
2.2 TMB、MATH和MSI
MSI分析結果顯示，6例患者為MSI，占患者總數的25%，其余患者MSI分數均相對較低，為MSS型。評估MSI狀態有助于確定臨床預后，并可能預測腫瘤對化學治療的反應。TMB結果顯示，大多數患者的TMB水平較低，只有3例患者的TMB值較高，中位TMB為2.045，平均TMB為9.86。晚期克隆擴增、空間隔離和不完全選擇性掃描導致遺傳上不同的細胞群都會出現腫瘤內異質性，MATH評分可以反映這種腫瘤內異質性。絕大多數患者的MATH值很高，只有少數患者的MATH值<20（表1）。
表1
胸腺癌隊列MSI、TMB、MATH的結果
2.3 HLA分型和TNB分析
HLA分析結果表示，大多數樣本都有雜合的HLA等位基因。所有樣本中共檢測到14個HLA-A位點、24個HLA-B位點和15個HLA-C位點。HLA-A*11:01、B*46:01和C*01:02是TC的易感基因（表2）。
表2
胸腺癌隊列HLA分析結果
TNB結果顯示絕大多數患者的TNB值很低，只有1例患者的TNB值較高（表3）。
表3
胸腺癌隊列新抗原負荷結果
2.4 同源重組缺陷、腫瘤純度和倍性
腫瘤組織樣本包含癌細胞和正常細胞的混合物。因此，腫瘤數據分析必須考慮標本的腫瘤純度和倍性，以便于在其他分析中剔除這些因素帶來的影響。大多數樣本的倍性為2.0，純度大于0.8（表4）。
表4
胸腺癌隊列的腫瘤純度和倍性結果
HRD-LOH是同源重組缺陷相關的LOH，全部為雜合性缺失，其中樣本K004002T和K034958T在HRD-LOH、LST、TAI的數量均為最多，其余樣本的個數均較少（表5）。
表5
胸腺癌隊列HRD-LOH的結果
2.5 分組分析和預后分析
為了更好地探索臨床特征與基因突變的關系，將患者分別按照性別、年齡、是否吸煙、病理分型、Masaoka分期、治療方案、是否生存、腫瘤直徑、是否復發分成9個組，分析每組樣本的數量。將MATH、HRD、CNV、純度、倍性等生物標志物進行組間統計學分析（表6）。
表6
本研究的分組方式及對應的樣本數統計
CNV突變負荷結果在多個分組內均有顯著差異。3組吸煙者的CNV負荷（P=0.026）和CNV-loss負荷（P=0.018）均高于不吸煙者（圖3a）。6組手術患者和放化療患者的CNV負荷（P=0.0079）、CNV-gain負荷（P=0.04）、CNV-loss負荷（P=0.0078）均顯著高于新輔助放化療患者（圖3b）。7組生存患者的CNV負荷（P=0.036）和CNV-loss負荷（P=0.018）顯著低于死亡患者（圖3c），死亡患者的CNV-gain負荷（P=0.44）有高于生存患者的趨勢。9組是否復發（圖3d）與7組是否死亡有相似的結論（P=0.0092，P=0.0033，P=0.16）。
圖3
胸腺癌患者拷貝數變異負荷在不同分組的差異
a：CNV負荷、CNV gain負荷和CNV loss負荷分別在3組是否吸煙組； 6組不用治療方案組（b）、7組是否生存組（c）和9組是否復發組（d）的統計分析；CNV：染色體拷貝數變異
MATH 第6組不同治療方案組差異有統計學意義（P=0.041，圖4a），新輔助放化療組的MATH值顯著低于手術組。第8組腫瘤直徑50～100 mm患者的MATH值顯著低于腫瘤直徑<50 mm患者和>100 mm患者（P=0.033）。腫瘤的倍性（圖4b）在1組性別方面差異有統計學意義，女性樣本的倍性顯著高于男性樣本（P=0.027），但是不同性別樣本的腫瘤純度差異無統計學意義（P=0.98）。HRD結果在7組是否生存和9組是否復發的差異有統計學意義。7組死亡患者的HRD-LOH數量（P=0.033）、HRD總數（P=0.024）、LST數量（P=0.026）顯著高于生存患者（圖4c）。9組復發患者的HRD-LOH數量（P=0.046）、HRD總數（P=0.016）、LST數量（P=0.017）顯著高于非復發患者（圖4d）。但是兩者TAI差異均無統計學意義（P=0.92）。
圖4
胸腺癌患者多個生物標志物在組間的差異
a：MATH分別在6組不同治療方案患者和8組腫瘤大小的統計分析；b：腫瘤倍性和細胞性分別在1組性別的統計分析中；c：HRD-LOH、HRD總數、LST數量和TAI數量分別在7組是否生存的統計分析；d：HRD-LOH、HRD總數、LST數量和TAI數量分別在9組是否復發的統計分析；e：在4組不同病理患者間的預后分析；f：在9組是否復發的預后分析；MATH：突變等位基因腫瘤異質性；HRD：同源重組缺陷；LOH：雜合性缺失；TAI：端粒等位基因不平衡； LST：大片段遷移
通過Cox回歸進行單因素預后分析，分別計算所有臨床特征和患者OS之間的關系。其中，TC的不同病理分型（圖4e）和是否復發（圖4f）是與OS顯著相關的因素。胸腺神經內分泌癌的OS最短，胸腺鱗狀細胞癌的OS相對最長（P=0.0075）。未復發患者的OS顯著長于復發患者（P<0.001）。
2.6 TCGA數據庫對比
通過本隊列TC數據與TCGA數據庫中胸腺瘤數據的對比，得到兩癌種突變圖譜的差異和相似點。胸腺瘤和TC數據的SNV突變圖譜表明（圖5a），突變頻率最高的前4個基因為GTF2I（41%）、TTN（11%）、HRAS（8%）和HSPG2（6%）。僅在胸腺瘤隊列中，GTF2I基因的突變頻率最高（49%），其次是HRAS（8%）、MUC16（5%）和TTN（5%），其中GTF2I僅在胸腺瘤數據中檢測到。TC樣本的SNV平均突變數高于胸腺瘤樣本。
圖5
胸腺癌患者多數據集的突變圖譜對比
a：本隊列TC數據與TCGA胸腺瘤數據的SNV突變情況對比；b：本隊列TC數據與TCGA胸腺瘤數據的CNV突變情況對比；SNV：單核苷酸突變；CNV：染色體拷貝數變異
胸腺瘤和TC數據的CNV突變圖譜表明，突變頻率最高的前4個基因為RIPK1（24%）、AKAP12（23%）、BEND3（22%）和NUP43（22%）。其中，HLA-E 和DUS2基因僅在TC數據中被檢測到突變。其次在TC數據中僅有NBPF10基因是CNV-gain型的突變，但在胸腺瘤數據中，發生CNV-gain突變的基因數量遠高于TC隊列。這表明胸腺瘤和TC的分子畸變模式不同（圖5b）。
3 討論
我們的研究結果顯示，TTN是突變頻率最高的基因（42%），其次是HSPG2（29%）和OBSCN（29%），這早期研究[27-28]結果相似。TTN是目前最大的蛋白，由364個外顯子組成，位于染色體2q31上[29]。TTN突變在人類癌癥中起重要作用。研究[30]證實，TTN基因突變是甲狀腺癌的獨立危險因素［HR=4.558，95%CI（1.808，11.494），P=0.001］。TTN可與lncRNA TTN-AS1結合，在乳腺癌、肺癌、消化系統腫瘤等許多癌癥中具有促癌作用[31]。更重要的是，在結直腸癌中，TTN/OBSCN“雙重打擊”是一種“免疫熱”亞型，具有潛在更好的免疫治療效果，是良好預后的預測因子[32-33]。在本研究中，TTN和OBSCN的突變頻率極高，可能是TC的危險因素，并對患者免疫治療起重要作用。
與大多數成人腫瘤類型相比，TET的TMB相對較低，且TC的TMB水平遠高于胸腺瘤[34-36]。既往研究中，胸腺瘤的平均TMB為0.663[37]。本隊列中TC平均TMB為9.86，顯著高于文獻[34-36]報道的胸腺瘤的平均TMB，同時本研究中TC隊列的突變數遠遠高于胸腺瘤隊列（圖5a）。臨床研究[38]表明，高TMB患者有免疫治療的機會，且對免疫治療的反應明顯優于低TMB患者。這是因為高TMB患者可能含有更多激活免疫細胞的新抗原，因此在使用免疫檢查點抑制劑（ICIs）治療后具有更長的生存期和更高的反應率[39]。2017年，新英格蘭醫學雜志發表了27種腫瘤的中位TMB與ICIs治療的客觀緩解率（ORR）之間的線性圖，發現TMB與ORR存在顯著的正相關性（P<0.001），相關系數為0.74，表明在這27種腫瘤中TMB能夠預測ICIs的療效[40]。2020年6月，帕博利珠單抗第二個不限癌種的新適應證獲得美國食品藥品監督管理局（FDA）批準，用于治療具有高組織腫瘤突變負荷（tTMB-H）無法切除或轉移性實體瘤患者。根據TC具有較高TMB的特點，未來或許可以從免疫治療中獲益。綜上，TC患者進行TMB檢測是非常必要的。但是并非所有的突變都是具有免疫原性的，后續對TC免疫細胞的組成和功能的深入了解對有限管理TC的進展和免疫反應至關重要[35]。
MSI-H在泛癌種中占比為0%～16.5%。在子宮內膜癌中MSI占比高達16.5%，但是在黑色素瘤（0%）、肺癌（0.36%）和很多血液腫瘤中MSI-H是罕見的[41]。2017年5月，美國FDA批準帕博利珠單抗（pembrolizumab）用于無論腫瘤來源，微衛星不穩定性高（MSI-H）的實體腫瘤患者。曾有文獻[42]報道1例MSI型TC患者接受阿維單抗（Avelumab）和阿昔替尼（Axitinib）聯合治療后，獲得持久的部分反應。在我們的研究中，24例患者中有6例是MSI型，占總人數的25%，與其他癌種相比TC的MSI型患者占比相對較高，因此從免疫治療中獲益的可能性很高，同時也表明TC患者可能需要進行MSI檢測來指導用藥。
MATH值 是一個腫瘤異質性的定量指標。以前的研究[43-44]表明顯示，高MATH與較低的總生存率有關。 不幸的是，在本研究中，絕大多數患者的MATH值較高，只有少數患者的MATH值<20。這表明TC預后相對較差。
HRD檢測是通過檢測腫瘤細胞基因組中與DNA修復機制有關的基因變異情況，來評估細胞中HRD的程度，從而判斷腫瘤的基因組穩定性[45-46]。目前已經在很多癌種中證實，HRD與腫瘤患者的預后密切相關，HRD陽性腫瘤可能預示患者預后更差。LOH現象在癌癥中普遍存在，已有研究[47]表明LOH導致抑癌基因失活，從而影響癌癥的發生和進展。該結論在本研究中也得到證實，復發患者HRD數量顯著高于未復發患者（圖4d），且復發患者OS顯著短于未復發患者（圖4f）。同理，死亡患者的HRD數量顯著高于生存患者（圖4c）。
盡管本研究初步取得了一些重要的發現，但樣本量相對較小，一定程度上限制了我們對結果的解釋和及其可推廣性。因此，我們計劃在后續研究中擴大樣本量，以加強和驗證我們的研究結果，更全面地闡述研究問題，增加結論的可信度。
總之，本研究通過NGS WES子技術檢測了TC的體細胞突變，并利用生物信息學方法進一步分析了多維度的生物標志物。這項研究提供了迄今為止第一個關于中國TC突變情況的全面分析。本隊列中，在SNV突變層面，TTN基因的突變頻率最高（42%），其次是HSPG2（29%）、OBSCN（29%）和TSSC1（25%）。發生CNV突變頻率最高的是ZNF276基因（54%，Loss），其次是BEND3（50%，Loss）、DHODH（50%，Loss）和VAC14（50%，Loss）。通過研究免疫治療的可能性發現，TC患者的TMB相比胸腺瘤較高，且MSI型患者的占比較高，表明TC患者有機會接受免疫治療，并有較長的有效生存期。這些發現有助于更好地了解TC的分子發病機理，并可能為治療TC提供指導。我們將在未來的研究中收集更多的臨床樣本并驗證這些結果。
利益沖突：邵偉康受雇于無錫臻和生物科技有限公司；其他作者聲明無潛在利益沖突。
作者貢獻：向潤負責文獻檢索，論文初稿撰寫；謝少華、廖瓊負責論文初稿審閱、對文章知識性內容作批判性審閱；邵偉康負責數據分析，初稿撰寫和校對；李娟、李強負責論文設計指導，參與制定研究思路，全文審校。
