浙江省結核分枝桿菌耐藥突變特征及其與基因型相關性分析

吳坤陽 陸燁瑋 張明五 朱業蕾 李相辰 潘軍航 王曉萌 王蔚 江敏敏 彭小軍 王煒欣 高俊順 柳正衛

結核病的耐藥問題依然嚴峻，尤其是耐多藥結核病(multidrug-resistant tuberculosis，MDR-TB)對結核病防控構成嚴重威脅。據世界衛生組織統計，2020年71%的病原學確診肺結核患者接受了利福平耐藥檢測，發現了13.2萬例耐多藥或利福平耐藥結核病患者(multidrug/rifampicin-resistant tuberculosis, MDR/RR-TB)和2.6萬例準廣泛耐藥或廣泛耐藥結核病患者(pre-extensively or extensively drug-resistant tuberculosis, pre-XDR/XDR-TB)[1]。

2013—2014年浙江省開展的第四次結核病耐藥監測顯示，浙江省結核病初治耐藥率和耐多藥率均低于全國水平，但仍高于全球平均水平[2]。全基因組測序(whole-genome sequencing，WGS)技術目前已廣泛應用于結核分枝桿菌(Mycobacteriumtuberculosis，MTB)的研究當中[3-5]，該技術既可以克服傳統表型藥物敏感性試驗(簡稱“藥敏試驗”)及其他分子藥敏試驗檢測的局限性，也可以用于耐藥性分析和基因型鑒定[6]。故筆者運用WGS技術對2018—2019年浙江省第五次結核病耐藥監測時收集的MTB菌株耐藥相關基因突變及其與基因型的相關性進行分析，以評估WGS藥敏試驗與表型藥敏試驗檢測結果的一致性，并了解浙江省目前結核病耐藥現狀，為本地區耐藥結核病的早診斷和早治療提供科學依據，探索WGS替代傳統表型藥敏試驗進行耐藥監測工作的可能性。

材料和方法

一、菌株來源

2018—2019年，浙江省參照世界衛生組織《結核病藥物耐藥性監測指南》的方法，在浙江省第四次耐藥監測的原有30個耐藥監測點開展了第五次結核病耐藥監測[2]。按監測指南要求，對納入調查的996例涂陽肺結核患者痰液標本進行細菌分離培養，共獲得996株分枝桿菌菌株，剔除47株污染菌、33株未生長菌株、81株非結核分枝桿菌，最終835株MTB菌株用于研究分析。MTB標準株H37Rv作為監測菌株由中國疾病預防控制中心國家結核病參比實驗室提供。

二、研究方法

1.菌株復蘇培養：將耐藥監測分離保存的菌株凍存管自然解凍后，吸取0.1～0.15 ml菌液均勻接種于改良羅氏培養基上，每株菌株接種于2支培養管，置于(36±1) ℃恒溫培養箱中孵育復蘇，具體操作參考《結核病實驗室檢驗規程》[7]。

2.表型藥敏試驗：將復蘇后的菌株進行比例法固體藥敏試驗。取培養基上生長1～2周的菌落，制備成10-2mg/ml和10-4mg/ml菌懸液，各取0.01 ml分別均勻接種至對照培養基和藥敏培養基表面。接種完成后，將培養基置于(36±1) ℃恒溫培養箱中培養4周后報告結果。若受試菌株耐藥百分比>1%，則判斷為耐藥。本研究由于只有異煙肼檢測到的耐藥菌株數量較多，且異煙肼可以較為準確地評估WGS耐藥性分析方法的性能，故本研究只使用異煙肼藥敏羅氏培養基檢測菌株的異煙肼表型藥敏，具體操作和耐藥百分比計算參照《結核病實驗室檢驗規程》[7]進行。

3.WGS：轉種的菌落初生長至2周時，刮取培養基斜面2環菌放于含有300 μl去離子水的1.5 ml離心管中，80 ℃滅活30 min。委托杭州廣科安德生物科技有限公司對MTB分離菌株進行測序，利用高通量二代Illumina測序技術平臺，構建DNA文庫，完成PE150雙端測序。

4.測序菌株質控及篩選：首先，使用fastp軟件(https://github.com/OpenGene/fastp)去除原始測序數據中的接頭序列，并進行低質量堿基的過濾處理[8]。接下來，將處理后得到的高質量序列數據輸入Kraken軟件(http://ccb.jhu.edu/software/kraken/)進行物種鑒定，去除鑒定結果為其他物種或者鑒定到MTB比例過低(<80%)的污染樣本。最后，使用bwa軟件(http://bio-bwa.sourceforge.net/)將剩余樣本的測序數據與H37Rv標準株基因組(NC_000962.2)序列進行比對。篩選出平均測序深度>20X、基因組覆蓋度>95%的滿足分析要求的808株樣本測序數據進行后續的WGS數據分析。

5.WGS數據分析：將808株測序樣本數據輸入至TB-Profiler軟件的本地版本(https://github.com/jodyphelan/TBProfiler)，以鑒定突變位點基因型，并檢測其對14種抗結核藥物(包括異煙肼、利福平、乙胺丁醇、吡嗪酰胺、鏈霉素、氟喹諾酮類、阿米卡星、卡那霉素、卷曲霉素、乙硫異煙胺、對氨基水楊酸、利奈唑胺、貝達喹啉、氯法齊明)的耐藥譜[9-10]。使用TB-Profiler軟件檢測各菌株相較于H37Rv標準株基因組(NC_000962.2)的單核苷酸多態性(single nucleotide polymorphism, SNP)，將檢測到的SNP與已知的基因型特異性SNP比對獲得菌株的基因型及其亞型，再與已知的耐藥基因突變數據庫比對獲得菌株的耐藥突變譜。以是否發生耐藥相關突變作為判斷藥物敏感性的標準。

三、統計學處理

利用R 4.0.5軟件進行數據的統計分析。以α=0.05為檢驗水準，通過χ2檢驗或矯正χ2檢驗分析不同藥物主要突變類型的發生頻率在不同基因型間(北京基因型和非北京基因型)的差異是否有統計學意義。為進一步了解WGS與表型藥敏試驗結果的一致性，本研究對WGS藥敏試驗結果中耐藥菌株數較多的異煙肼結果進行了表型藥敏試驗檢測結果的比對分析，以敏感度、特異度、陽性預測值、陰性預測值、一致率和Kappa值為指標計算WGS技術的耐藥預測效能。其中，Kappa值≥0.75為一致性很好，0.40≤Kappa值<0.75為一致性較好，Kappa值<0.40為一致性較差。平均測序深度為平均每份樣本的有效測序數據中覆蓋到該位點的測序reads數，平均測序深度越深表明該突變位點的檢測結果越可靠，一般要求>20；平均突變頻率(%)為檢測到該突變的reads數與覆蓋到該位點的總reads數的比例，鑒于MTB的基因組是單倍型，其固定下來的耐藥突變頻率一般>90%，故較低的突變頻率可反映異質性耐藥的情況及耐藥突變在基因組上進行固定的過程。

結 果

一、MTB菌株耐藥突變特征

808株MTB菌株的WGS藥敏試驗結果顯示，153株發生了除利奈唑胺、貝達喹啉、氯法齊明以外的11種抗結核藥物耐藥相關基因突變，總突變率為18.9%；其中，MDR-MTB菌株占2.1%(17/808)，pre-XDR-MTB菌株占0.9%(7/808)。11種發生耐藥相關基因突變的藥物中，對異煙肼的耐藥率最高，達到9.3%；其次是鏈霉素和氟喹諾酮類，分別為6.8%和5.1%；而對其他藥物的耐藥率均在5%以下。具體見表1、2。

表1 808株結核分枝桿菌對11種抗結核藥物的耐藥基因突變情況

續表1

表2 808株結核分枝桿菌對11種抗結核藥物的耐藥基因突變類型構成情況

二、MTB菌株對一線抗結核藥物的耐藥相關基因突變特征(表1，2)

1.異煙肼耐藥相關基因突變特征：75株異煙肼耐藥菌株中，不同耐藥相關基因的突變頻率分布為katG(78.7%， 59/75)、inhA(22.7%，17/75)、ahpC啟動子區(6.7%，5/75)。耐藥基因中共產生21種突變類型，其中最主要的突變類型為katG-315-S/T；此外，發現6株(8.0%，6/75)異煙肼耐藥菌株存在≥2個位點聯合突變，以katG聯合inhA或ahpC啟動子區突變居多。

2.利福平耐藥相關基因突變特征：28株利福平耐藥菌株中，全部在利福平耐藥相關rpoB基因上發生耐藥突變，其中2株同時在rpoC基因上發生耐藥突變。在rpoB基因上共發現9種突變類型，其中92.9%(26/28)的突變發生在rpoB基因耐藥決定區(rifampicin-resistance determining region，RRDR)，最主要的突變類型為rpoB-450-S/L(57.1%)。另外發現有5株(17.9%，5/28)利福平耐藥菌株存在≥2個位點的聯合突變；其中rpoB基因聯合rpoC基因發生耐藥突變的2株菌株的突變類型分別是rpoB-450-S/L+rpoC-527-L/V和rpoB-170-V/F+rpoC-527-L/V。

3.乙胺丁醇耐藥相關基因突變特征：16株乙胺丁醇耐藥菌株中，發現embB(93.8%，15/16)和embA啟動子區(12.5%，2/16)發生耐藥突變。其中embB基因306、354、423、406、504和1024位點分別占56.3%(9/16)、12.5%(2/16)、6.3%(1/16)、6.3%(1/16)、6.3%(1/16)和6.3%(1/16)，embA啟動子區中-16和-12位點分別占6.3%(1/16)和6.3%(1/16)。耐藥基因中共發生9種突變類型，其中最主要的突變類型為embB-306-M/V(43.8%)，且其中有1株(6.3%)存在位點聯合突變，為embB基因和embA啟動子區聯合突變，突變類型是embB-406-G/A+embA啟動子-12-C/T。

4.吡嗪酰胺耐藥相關基因突變特征：12株吡嗪酰胺耐藥菌株中，耐藥突變均發生在pncA基因上，且突變位點分散，幾乎遍布整個編碼區域。突變類型共有12種，未發現主要突變類型。

5.鏈霉素耐藥相關基因突變特征：55株鏈霉素耐藥菌株中，耐藥突變發生在rpsL、rrs和gid基因，共產生11種突變類型，主要突變類型為rpsL-43-K/R(65.5%)，此外發現1株(1.8%)存在rrs基因和gid基因位點聯合突變。

三、對二線抗結核藥物產生的耐藥相關基因突變特征(表1，2)

1.氟喹諾酮類耐藥相關基因突變特征：41株耐藥菌株中，耐藥突變主要在gyrA和gyrB基因上，以gyrA基因為主(90.2%，37/41)，主要突變類型為gyrA-94-D/G(36.6%)和gyrA-90-A/V(29.3%)。

2.阿米卡星、卡那霉素、卷曲霉素耐藥相關基因突變特征：三者耐藥菌株數分別為3、4和3株，共計4株菌株。其中有3株對3種藥物同時耐藥，另有1株僅對卡那霉素耐藥，分別在rrs基因和eis啟動子區發現了耐藥突變，突變類型以rrs-1402-C/A和rrs-1484-G/T為主；其中有2株(2/4)發生rrs-1402-C/A和rrs-1484-G/T聯合突變，1株(1/4)卡那霉素單獨耐藥菌株上發現eis啟動子區耐藥突變eis啟動子-10-G/A。

3.乙硫異煙胺耐藥相關基因突變特征：20株乙硫異煙胺耐藥菌株中，在inhA和ethA基因上發現主要的突變類型為inhA-15-C/T。

4.對氨基水楊酸耐藥相關基因突變特征：8株對氨基水楊酸耐藥菌株中，耐藥相關基因為thyA和folC，發現主要突變類型為thyA-75-H/N和folC-150-S/G。

四、遺傳背景與耐藥突變的相關性分析

808株MTB菌株中，586株(72.5%)為北京基因型菌株，222株(27.5%)為非北京基因型菌株。受限于本研究中耐藥菌株數較少的情況，僅對突變率較高的4種耐藥基因katG、rpsL、gyrA和rpoB的突變情況及突變率較高的4種耐藥基因位點katG-315、rpsL-43、gyrA-94、rpoB-450在不同基因型中的分布進行相關性分析。結果顯示，鏈霉素耐藥相關rpsL基因突變僅在北京基因型中發現，與在非北京基因型中的分布差異有統計學意義(表3)。同時，katG-315位點和rpsL-43位點在北京基因型中的突變率均明顯高于非北京基因型(表4)。

表3 808株結核分枝桿菌4種耐藥基因突變在不同基因型中的分布情況

表4 808株結核分枝桿菌4種耐藥基因突變位點在不同基因型中的分布情況

五、WGS藥敏試驗與表型藥敏試驗檢測結果一致性分析

對808株菌株進行異煙肼表型藥敏試驗檢測，共獲得750株菌株的表型藥敏結果，檢出56株(7.5%)異煙肼耐藥菌株。以表型藥敏試驗結果為參照標準，WGS對異煙肼的耐藥性檢測敏感度、特異度、陽性預測值、陰性預測值、一致率和Kappa值分別為87.5%(49/56)、98.0%(680/694)、77.8%(49/63)、99.0%(680/687)、97.2%(729/750)和0.808(表5)。

表5 全基因組測序以表型藥敏試驗結果為參照標準對異煙肼耐藥性的檢測效能

在21株WGS與表型藥敏試驗結果不一致的菌株中，katG-315-S/T和ahpC啟動子-48-G/A是導致結果不一致的主要位點，分別有5株和3株，并且發現2株存在位點聯合突變，分別為inhA-15-C/T+katG-110-A/V和ahpC啟動子-48-G/A+katG-607-E/K。同時，在結果不一致的菌株中，有7株菌株是WGS未能檢出的野生型耐藥菌株(表6)。

表6 21株兩種藥敏試驗方法檢測結果不一致菌株的基因突變類型分析

討 論

WGS技術預測MTB菌株對于不同抗結核藥物耐藥性的性能已得到了相關研究的驗證[11-12]，同時也獲得了世界衛生組織的認可和推薦[13]。該方法可以克服目前臨床分子檢測方法存在的不足，如GeneXpert MTB/RIF和GenoType MTBDRplus等方法只能檢測到部分抗結核藥物及其中的部分耐藥基因突變位點[14-15]，對指導結核病臨床用藥和控制耐藥結核病意義重大。浙江省第四次結核病耐藥監測結果表明本省耐藥問題依然嚴峻[2]，故本研究針對浙江省第五次結核病耐藥監測的MTB菌株，使用WGS技術揭示了目前浙江省主要的耐藥突變類型及其與菌株基因型之間的相關性，獲得了本地區比較全面的耐藥突變特征。

本研究中發現的利福平、乙胺丁醇、氟喹諾酮類藥物的主要突變類型和突變頻率與已有研究一致[16-19]，但異煙肼耐藥菌株出現了更高比例的inhA基因耐藥突變(22.7%)，明顯高于深圳市(8.6%)和全國范圍(15.2%)的研究[18,20]，提示本地區可能存在較高比例的異煙肼耐藥患者，抗結核藥物治療難度加大，這可能與inhA基因突變常發生異煙肼低水平耐藥[21]，使得標準劑量的異煙肼對inhA基因突變的患者治療無效，只能使用高劑量并與其他抗結核藥物聯合治療才可以取得一定的療效有關[22]。因此，建議浙江省可以通過WGS技術檢測inhA基因耐藥突變，有針對性地使用高劑量異煙肼療法以提高本地區耐藥結核病的治愈率。另外，本研究還發現6株(8.0%)異煙肼耐藥菌株存在katG聯合inhA或ahpC啟動子區突變的情況，提示患者體內可能存在異質性耐藥[23]，也體現出WGS技術在異質性耐藥檢測上相較于其他分子檢測方法的優勢[24-25]。其中，有1株異煙肼耐藥菌株為ahpC+inhA+katG的聯合突變形式，值得后續進一步研究分析。相關文獻報道顯示，吡嗪酰胺主要的耐藥相關基因為pncA，其突變類型較多[26]，主要突變類型為pncA啟動子-11-T/C[18-19]，與本研究發現的12株吡嗪酰胺耐藥基因突變均在pncA基因上、突變位點分散且幾乎遍布了整個編碼區域的結果一致，但本研究并未發現pncA啟動子-11-T/C的突變類型。由于受到樣本量的限制，本研究發現的對其他二線抗結核藥物耐藥菌株較少，無法更深入地挖掘其耐藥突變特征。總而言之，國內不同地區MTB耐藥譜中主要突變類型基本一致，但突變頻率可能存在一定的差異，這可能受到研究樣本量、菌株基因型及不同地區結核病臨床治療偏好等因素影響。而本研究的結果進一步說明，充分且全面地了解本地區MTB菌株耐藥突變特征，對評估分子檢測方法的適用性及制定臨床診療方案具有重要作用。

不同地區的流行菌株存在遺傳背景差異，可能會導致同一藥品在不同地區產生耐藥突變譜的差異[27]。本研究分析了突變率較高的耐藥基因及突變位點與基因型的相關性，發現586株(72.5%)為北京基因型菌株，222株(27.5%)為非北京基因型菌株，所占比例與全國其他地區相似，如北京[16]和上海[12]；而且，鏈霉素耐藥相關rpsL基因突變僅在北京基因型菌株中出現；同時，katG-315位點和rpsL-43位點在北京基因型中的突變率也明顯高于非北京基因型；但在樣本中未發現突變率較高的利福平耐藥相關基因rpoB和氟喹諾酮類耐藥相關基因gyrA與基因型之間的明顯相關性，這可能與樣本量不足有關，需要擴大樣本量以進一步證實北京基因型與上述耐藥相關基因的相關性。

WGS耐藥性分析方法是否可以在結核病耐藥監測工作中發揮更加重要的作用，還需要進行系統性評估。本研究對異煙肼的WGS分子藥敏結果與表型藥敏試驗檢測結果進行了交叉比對，并結合異煙肼相關耐藥突變特征進行分析，展示了一套系統性的評估流程，為WGS耐藥性分析方法的普及和推廣進行了探索性工作。本研究中，WGS對異煙肼耐藥檢測的敏感度(87.5%)、特異度(98.0%)和一致率(97.2%)均較高，Kappa值達到了0.808，表明兩種方法對異煙肼的耐藥性檢測結果一致性較高，可重復性強，提示在耐藥監測工作中WGS具有取代表型藥敏試驗的潛力。但是，也發現WGS對部分耐藥突變位點的預測性能較弱，如ahpC啟動子-48-G/A、inhA-15-C/T和katG-315-S/T這3個位點在本研究中的不一致率分別為3/4、16.7%和11.6%，這些性能較弱的突變位點可導致假陽性結果的產生，降低了陽性預測值(77.8%)，對整體分析性能存在著制約作用。ahpC啟動子-48-G/A在本研究中較高的不一致率提示該突變位點可能在浙江地區存在較低的預測性能，需要收集更多的耐藥菌株進一步驗證。

本研究也存在一定的局限性。雖然收集到808株MTB菌株的測序數據，但從部分抗結核藥物檢測到耐藥相關基因突變的菌株并不多，無法進行進一步分析。另外，受到不同藥物耐藥菌株數量的限制，本研究僅對耐藥率較高的異煙肼進行了WGS耐藥性分析的性能評估，但認為該評估方法可應用于其他抗結核藥物并值得推廣。

綜上所述，本研究對浙江省第五次結核病耐藥監測數據進行WGS技術檢測，揭示了MTB菌株對一線及二線抗結核藥物的耐藥相關基因突變譜及其與基因型的相關性，在分子層面展示了本地區的耐藥情況，為耐藥結核病分子檢測方法的優化提升以及臨床用藥的精準指導提供了數據和理論支持，同時對WGS技術在耐藥監測工作中普及和推廣進行了探索性的研究。

利益沖突所有作者均聲明不存在利益沖突

作者貢獻吳坤陽和陸燁瑋：醞釀設計實驗，實施具體研究，收集分析數據，撰寫論文以及審閱修改；張明五、朱業蕾和李相辰：參與設計完善實驗步驟，幫助實施研究，分析解釋部分數據，撰寫論文部分內容；潘軍航、王曉萌、王蔚、江敏敏、彭小軍、王煒欣和高俊順：指導部分研究、幫助采集分析數據；柳正衛：審閱論文，指導研究，為研究提供經費、行政、技術和材料支持