高產聚蘋果酸黑色素短梗霉CGMCC18996全基因組組裝注釋及關鍵蛋白分析

王舸楠，李佳謙，李雨桐，陳世偉，王淑賢，趙廷彬，賈士儒，喬長晟,,

（1.天津科技大學生物工程學院，天津 300457；2.工業發酵微生物教育部重點實驗室暨天津市工業微生物重點實驗室，天津 300457；3.天津市微生物代謝與發酵過程控制技術工程中心，天津 300457；4 .天津慧智百川生物工程有限公司，天津 300457）

聚蘋果酸（polymalic acid，PMLA）是以蘋果酸為唯一單體的均聚高分子聚合物，屬于聚酯類聚合物，具有高生物相容性、高水溶性、生物可吸收性、化學可衍生性、可降解性和無免疫原性等多種優良性能[1]。PMLA分子的側鏈帶有可修飾性較強的羧基，一些分子能夠通過官能團反應引入其聚合鏈上，同時PMLA水解產生的蘋果酸是理想的調酸劑，這使PMLA在生物醫藥、食品和生物材料領域具有潛在的應用前景[2-6]。

研究表明產黑色素短梗霉（Aureobasidium melanogenum）是一種具有較強產PMLA能力的類酵母真菌[4]，但目前對其PMLA生物合成研究并不透徹，對通路中關鍵酶的研究并不清晰。產量還有很大的提升空間[7]。同時，該菌基因組信息依舊較少，限制了該菌株的開發利用。因此，對該菌種進行基因組測序及組裝可為改造菌種，提高產量提供理論依據。

目前，主流的基因測序分為兩種，單分子測序（如Pacbio三代測序平臺）和高通量測序（如二代Illumina平臺），其中單分子測序具有讀長長組裝結果好的特點，但其測序成本與錯誤率較高；高通量測序采用雙端測序策略，結果兼備質量高、價格低的特點，但其測序長度較短，組裝會產生較多的片段（contigs）。因此采用三代測序結果組裝再通過二代測序結果進行糾錯（基因組polish），能夠很大程度降低contigs數目同時保持基因組準確率，在完成基因組組裝后，可通過轉錄組測序結果對基因組進行結構注釋，以同時獲得高質量基因組注釋結果[8]。

本研究通過PacBio Sequel II及Illumina NovaSeq 6000測序平臺對高產PMLA產黑色素短梗霉基因組進行測序，通過不同組裝軟件對測序的下機文件進行組裝及優化；結合轉錄組數據對組裝結果進行基因結構注釋。之后對基因組注釋結果進行不同數據庫的功能注釋，同時分析PMLA合成關鍵蛋白。通過這種方法，期望得到適合于基因組分析以及后續分子生物學實驗的高質量基因組，為產黑色素短梗霉的開發利用提供一定生物信息學參考，并同時為其他類似物種的基因組組裝提供思路。

1 材料與方法

1.1 材料與試劑

菌株黑色素短梗霉，已保存至中國普通微生物菌種保藏管理中心（保藏號：CGMCC18996）；基因組提取試劑盒 湖南艾科瑞生物工程有限公司；乙腈（色譜純）北京鼎國昌盛生物技術有限公司；磷酸二氫鉀、磷酸氫二鉀、磷酸（均為色譜純）天津市大茂化學試劑廠。

1.2 儀器與設備

UC7超薄切片機 德國Leica公司；LC-10高效液相色譜儀 日本島津公司；ZWYR-D2403恒溫培養振蕩器上海智城分析儀器制造有限公司；LRH-250A生化培養箱 泰宏醫療器械有限公司；NanoDrop One超微量分光光度計 美國賽默飛世爾科技公司。

1.3 方法

1.3.1 基因組提取與測序

本實驗的基因組提取、建庫以及測序委托諾禾致源科技股份有限公司（北京）完成，實驗流程均按照諾禾致源公司測序方法進行（https://www.novogene.com/novo/sd_cpjs_6.html，https://www.novogene.com/novo/ed_cpjs_2.html）。

1.3.2 數據處理及生物信息學分析

使用曙光云計算服務器（曙光信息產業股份有限公司）進行生物信息學分析。最終得到相關的測序下機文件以及組裝的基因組文件已上傳至國家生物信息中心（項目號：PRJCA011444）。

1.3.2.1 原始下機文件獲取及質控

通過Illumina NovaSeq 6000平臺獲取高通量測序產生的fastq基因組測序文件，通過PacBio Sequel II測序平臺獲取單分子測序產生的bam基因組測序文件，同時驗證高通量測序和單分子測序文件的md5值與平臺提供是否一致。

1.3.2.2 基因組組裝及評估

使用SPAdes-3.15.5[9]組裝軟件對二代測序clean reads進行組裝，為了得到最優N50及L50值，將k-mer大小設定為77、87、97、107、117及127進行組裝。并通過QUAST[10]軟件對組裝結果進行評估。三代測序通過Flye[11]、Canu[12]、NextDenovo（https://github.com/Nextomics/NextDenovo）以及MECAT2[13]組裝軟件進行組裝，通過quickmerge[14]將最優結構與其他組裝結果進行融合，之后對得到的融合后結果進行去重復處理進一步降低contigs數目。

1.3.2.3 基因組注釋

使用BRAKER2[15]軟件進行基因組注釋，首先通過HISAT2將轉錄組序列比對到基因組上，生成bam文件，再將生成的bam文件與基因組文件作為輸入參數使用BRAKER2進行注釋。將獲得的氨基酸序列通過InterProScan[16]與eggNOG-mapper[17]進行注釋得到UniProtKB AC/ID（基因名），將這些基因名進行京都基因與基因組百科全書（Kyoto Encyclopedia of Genes and Genomes，KEGG）（https://www.kegg.jp/kegg/）、基因本體論（Gene Ontology，GO）[18]及直系同源集（Clusters of Orthologous Groups，COG）[19]數據庫注釋，并通過antiSMASH[20]對次級代謝產物合成基因簇進行預測。

1.3.3 菌體形態

透射電鏡拍攝方法如下，將電子顯微鏡固定劑（2.5%戊二醛溶于磷酸緩沖液）添加到經離心分離的菌株中，重懸菌體，使用 0.1 mol/L磷酸漂洗液漂洗3 次，1%鋨酸固定液固定1 h；之后使用0.1 mol/L磷酸漂洗液漂洗3 次，丙酮脫水3 次。使用純EMBed 812樹脂（90529-77-4）包埋，移入65 ℃烘箱中聚合48 h；使用超薄切片機切至60～80 nm薄片，用2%乙酸鈾飽和醇溶液和檸檬酸鉛染色15 min。最后，使用JEM1200電鏡進行觀察拍攝，得到放大4000 倍的產黑色素短梗霉CGMCC18996圖像。

1.3.4 蛋白質結構預測

蛋白質結構預測使用Alphafold2-2.2.0[21]，將基因組注釋的氨基酸序列作為輸入參數，運行模式參數選擇monomer，其余均選用默認選項進行蛋白質結構預測，最終在生成文件中選擇最優結果（即rank0結果）。之后將蛋白預測結果與蛋白質數據庫（Protein Data Bank，PDB）[22]中已上傳的晶體結構進行比對（磷酸烯醇式丙酮酸羧化酶激酶（phosphoenolpyruvate carboxykinase，PCKA）：1YLH[23]；蘋果酸合成酶（malate synthase，MASY）：3CUZ[24]）。

2 結果與分析

2.1 通過二代測序數據進行基因組預組裝

基于Illumina NovaSeq 6000測序平臺對產黑色素短梗霉基因組進行測序，共得到35.85 Gb×2的雙端測序結果。且該測序平臺得到的raw reads具有較高的質量評分（Q＞35）。之后，選取質控后clean reads進行后續的組裝實驗。

不同k-mer長度對產黑色素短梗霉clean reads進行組裝結果見表1，組裝結果顯示，產黑色素短梗霉基因組大小約為44 Mb。其中，在k-mer長度為127時得到N50值為908694，L50值為14，基因組覆蓋度為121.2×。基于該結果可以判斷該菌基因組長度理論值在44 Mb左右；且隨著k-mer長度的增加，N50值增大，L50值減小；若繼續增大k-mer值會進一步優化組裝參數。但由于SPAdes軟件進行組裝的k-mer值最大為127，這可能是考慮到繼續加大k-mer值對測序深度要求較高，從而提高測序成本。

2.2 三代測序組裝結果及基因組結構注釋

基于PacBio平臺的單分子測序共產生46 Gb大小的bam基因組測序文件，轉換為fasta文件后，根據二代組裝結果設置基因組大小為44 Mb進行組裝。各組裝軟件組裝結果如表2所示，其中，選用Canu組裝的最優結果通過quickmerge軟件與其他組裝結果進行融合，并結合二代測序文件進行基因組糾正（polish），在刪除重復contigs后的最終組裝結果見表3。

表2 三代組裝結果Table 2 Results of third-generation assembly

表3 產黑色素短梗霉基因組裝結果Table 3 Results of genome assembly of A.melanogenum

結合轉錄組測序文件進行基因組結構注釋共注釋出15684 個基因，并找到基因編碼區與氨基酸預測區，因結合轉錄組測序文件進行結構注釋，這些基因可能包含有可變剪切和重復注釋的結構，會增加注釋出的基因數目。因此，進行功能注釋和基因名稱注釋后需刪除重復的基因名；最終，共獲得6202 個基因注釋結果。

2.3 GO、KEGG、COG以及antiSMASH次級代謝注釋結果

對注釋出的6202 個基因進行GO、KEGG與COG數據庫注釋，結果如圖1所示。其中，COG注釋結果顯示大部分基因與碳水化合物轉運及代謝、氨基酸轉運代謝、轉錄后修飾、RNA加工及修飾有關。KEGG注釋結果顯示大部分基因所處代謝通路與核糖體、過氧化物體、RNA轉運有關。GO注釋結果顯示大部分基因與RNA、過氧化物體以及線粒體有關。最終，antiSMASH次級代謝物預測結果共發現4 個非核糖體肽合成酶（nonribosomal peptide synthetase，NRPS）基因簇、5 個一類聚酮合酶（polyketide synthetase，pks）基因簇、3 個β-內酯類合成基因簇以及7 個萜類合成基因簇，其中1 個一類pks基因簇和黑色素合成有關（相關性100%），1 個一類pks基因簇與黑麥酮酸類化合物合成有關（相關性18%）。

圖1 COG（a）、KEGG（b）、GO（c）數據庫注釋結果Fig.1 Results of annotation in the COG (a),KEGG (b),and GO (c) database

2.4 產黑色素短梗霉菌體透射電鏡結果

如圖2所示，在低產組中菌株具有較大細胞核（N）以及周圍存在的具有雙層膜結構線粒體（M），在高產組中出現了類似乙醛酸循環體的圓形結構，提示黑色素短梗霉中可能存在有乙醛酸循環途徑。

圖2 產黑色素短梗霉PMLA低產組（a、c）及PMLA高產組（b、d）透射電鏡結果Fig.2 TEM results of low-and high-PMLA producing A.melanogenum

2.5 PMLA合成相關基因結構預測

在對基因結構與基因名注釋后，找到PCKA、MASY的蛋白質序列，其晶體結構預測結果如圖3、4所示，其中，PCKA、MASY與PDB上傳晶體結構（PCKA：1YLH[23]，MASY：3CUZ[24]）進行比對，結果顯示預測的蛋白結構與PDB數據庫中序列比對結果基本一致，氨基酸比對結果中序列與其晶體結構中小分子配體的結合位點也具有較高的一致性。同時蛋白的保守作用位點如PCKA中86R、140V、146G、287G、288D、289D，MASY中440C、275C、276G、277R、278W在比對結果中一致。

圖3 Alphafold對PCKA預測結果（a）及與PDB中已上傳晶體結構的比對（b）和序列對比結果（c）Fig.3 Results of PCKA prediction by Alphafold software (a),crystal structure of PCKA versus that in PDB database (b),and amino acid alignment results (c)

圖4 Alphafold對MASY預測結果（a）及與PDB中已上傳晶體結構的比對（b）和序列對比結果（c）Fig.4 Results of MASY prediction by Alphafold software (a),crystal structure of MASY versus that in PDB database (b),and amino acid alignment results (c)

3 討論

結合基因組及轉錄組測序結果可以組裝并注釋出質量較高的基因組。本研究通過三代測序（Pacbio sequel II平臺）、二代測序（Illumina NovaSeq 6000）平臺對產黑色素短梗霉基因組進行測序，因三代組裝需要預估基因組大小，因此，首先通過二代數據進行基因組預組裝共得到44 Mbp基因組大小。之后考察了不同組裝軟件對產黑色素短梗霉基因組的組裝效果，在一般的默認選項下，Canu軟件得到了較好的組裝效果。對該組裝結果通過二代數據修正并去重后得到包含26 個contigs、N50為2204220、GC值為50.09%、大小為42 Mb的較高質量基因組組裝結果。通過轉錄組數據對該組裝結果進行結構注釋，共找到6202 個基因。產黑色短梗霉屬于出芽短梗霉的亞種，其在不同環境中也會呈現出酵母狀與菌絲體狀的不同形態[25]。將基因組組裝與注釋結果與酵母和絲狀真菌的模式生物基因組進行比較；其中釀酒酵母（Saccharomyces cerevisiae）基因組大小為12.15 Mb，編碼6016 個蛋白，GC含量為38.15%[26]；構巢曲霉基因組大小為30.30 Mb，編碼10008 個蛋白，GC含量為50.10%[27]。因此，產黑色素短梗霉基因組組成更偏向構巢曲霉，且本實驗室先前通過無參轉錄組注釋發現很多與構巢曲霉同源的基因[28]。因構巢曲霉是絲狀真菌的模式生物，該結果提示產黑色素短梗霉可能同樣適用構巢曲霉的分子轉化方法。

本研究通過Illumina NovaSeq 6000平臺對產黑色素短梗霉進行了大約35 GB下機文件大小的高通量測序，該測序量的測序深度較大，理論測序深度為1000×。通過二代組裝軟件進行組裝后，最終組裝出N50值為908694的基因組，該結果僅達到三代測序組裝結果N50的41.4%，因此，考慮到組裝結果以及建庫與測序成本[29]，三代測序進行組裝輔以二代測序進行基因組修正的測序方法更具有性價比。

通過COG、KEGG與GO數據庫注釋結果表明，產黑色素短梗霉中基因表達主要位于核糖體、線粒體以及過氧化物體等細胞器中，其中，有研究表明線粒體中的三羧酸循環、乙醛酸循環體中的乙醛酸循環以及細胞質中的還原性三羧酸途徑與PMLA的生物合成有關[7,30-31]。功能注釋結果中出現了大量與過氧化物體有關的基因，而能夠產生蘋果酸的乙醛酸循環體也屬于一類過氧化物體。因此，該結果說明在產黑色素短梗霉中可能存在乙醛酸循環體。實驗室先前研究發現乙醛酸體中的MASY可通過乙醛酸途徑生成蘋果酸，相比于線粒體中存在的蘋果酸/天冬氨酸穿梭體系，乙醛酸循環體中蘋果酸可能可以直接通過單層乙醛酸體膜進入細胞質中而被非核糖體肽合成酶（non-ribosomal peptide synthase，NRPS）聚合成PMLA[28]。透射電鏡結果也顯示，在加入CaCO3的高產PMLA組中，菌體內部出現了較多的圓形細胞器，提示菌體中存在乙醛酸循環體的可能。antiSMASH注釋結果表明在產黑色素短梗霉中存在4 個NRPS基因簇，該結果也同樣反映出NRPS蛋白聚合蘋果酸形成PMLA的可能性。同時，注釋出的與黑色素合成的pks基因簇可解釋該菌種在生長過程中逐漸變黑的現象。

實驗室先前研究表明[28]，PCKA、MASY和NRPS基因在PMLA高產組中發生了較大幅度的上調。因此，通過對基因組注釋得到的PCKA和MASY的編碼區進行了蛋白結構預測，并將預測結果與PDB上已有的晶體結構進行比對。比對結果表明，蛋白預測結果與已有晶體結構基本吻合，且保守結構域相似，提示預測蛋白與其他菌株中已存在蛋白可能具有相似的功能。PCKA是糖異生途徑的限速酶，一般情況下催化草酰乙酸轉化為磷酸烯醇式丙酮酸，但有研究表明PCKA在真菌細胞質中可反向催化磷酸烯醇式丙酮酸生成草酰乙酸[23]。而細胞質中的草酰乙酸會進一步轉化為PMLA的前體物質蘋果酸。MASY可催化乙醛酸循環中的第二部反應，將乙醛酸轉化為蘋果酸[24]，而乙醛酸體的單層膜結構可能使蘋果酸被動運輸進入細胞質，進而被細胞質中存在的NRPS蛋白聚合成為PMLA[28]。本實驗沒有對NRPS蛋白預測的原因是因為其氨基酸序列較長（5000 aa左右），軟件無法預測，其次是在PDB上還沒有找到近似的蛋白晶體結構，比對較為困難。

本研究對產黑色素短梗霉進行二代、三代基因組測序并組裝得到最優組裝結果，后通過轉錄組數據對基因組進行結構注釋，將得到的結構注釋結果進行不同數據庫的功能注釋。同時對產黑色素短梗霉進行了透射電鏡拍攝，電鏡結果提示菌株中存在乙醛酸循環體結構，該結構可能與PMLA合成有關，最終，實驗通過對PMLA合成相關的PCKA和MASY蛋白進行結構預測并與PDB上已有電鏡結果進行比對。發現產黑色素短梗霉中這兩種蛋白與PDB上已有晶體結構基本一致，且這兩種蛋白功能最終都與蘋果酸代謝有關。本實驗結果可為產黑色短梗霉菌株的PMLA代謝提供一定的參考，相關的測序下機文件以及組裝的基因組文件已上傳至國家生物信息中心（PRJCA011444），為后續的菌株開發利用提供基礎。