白斑狗魚雌、雄基因組混池重測序研究

張俊杰，趙瑞陽，蔣 麗，胡 瓊，李 菁，唐 露，李勝忠

(1.新疆農業大學動物科學學院，烏魯木齊 830052；2.中國水產科學研究院，北京 100141)

白斑狗魚(EsoxLucius)屬鮭形目(Salmoniformes)狗魚科(Esocidae)狗魚屬(Esox)，自然分布于北美洲及歐亞大陸北部淡水流域，在我國僅見于新疆北部的額爾齊斯河流域[1,2]。白斑狗魚作為自然水體中的主要頂級捕食者，其生長迅速，體型較大，具有重要經濟價值[3]。自人工繁殖取得成功以來，白斑狗魚已經成為新疆重要的特色養殖魚類，并逐漸向全國推廣，取得了顯著的經濟效益。

研究表明，白斑狗魚雌魚的生長速度明顯快于雄魚，而且性成熟時間更晚，壽命更長，體型更大[4-7]。Luczynski等[8]通過性逆轉與雌核發育結合，確定了白斑狗魚屬于XX/XY(雌性同型)的染色體決定類型，但近年來對白斑狗魚的研究主要集中在生物學特性、繁殖生物學、生態學等方面，其性別鑒定只能在性成熟后通過觀察和擠壓生殖孔等方法進行[9-13]。基因組重測序技術能夠在個體或群體樣本的全基因組水平上篩選出相關性狀的基因序列差異，已在多種魚類中獲得性別特異分子標記[14-17]。本實驗通過對雌、雄白斑狗魚進行混池重測序(BSA-seq)，探究白斑狗魚雌、雄性別的基因組差異，以期篩選出雌、雄特異分子標記，實現對白斑狗魚的性別控制技術及全雌魚單性生產。

1 材料與方法

1.1 實驗材料

實驗所用白斑狗魚購自烏魯木齊北園春市場，體重500～800 g，解剖后通過性腺觀察進行生理性別鑒定，然后剪取鰭條，在無水乙醇中保存，用于后續基因組DNA的提取。

1.2 實驗試劑

實驗所用試劑主要有：組織裂解液[1 mL Tris-HCl(1 mol/L，pH 8.0)，1 mL EDTA(0.5 mol/L，pH 8.0)，20 mL NaCl(0.5 mol/L)，10 mL SDS(10%)]、乙酸銨(7.5 mol/L)和PCR反應預混液(天根生化科技北京有限公司)。所用引物均由生工生物工程(上海)股份有限公司合成。

1.3 實驗方法

1.3.1 白斑狗魚DNA制備

選取鑒定性別后的雌、雄白斑狗魚各20尾，使用醋酸銨法提取基因組DNA，經全波長酶標儀(InfiniteM200，瑞士TECAN公司)測定濃度和A260 nm/A280 nm值后，送至上海凌恩生物科技有限公司進行測序。用同樣方法提取和測定雌、雄各24尾白斑狗魚的DNA，并稀釋至50 ng/μL，用于性別特異標記的PCR驗證。

1.3.2 混池重測序與數據質控

將20尾雌性和20尾雄性白斑狗魚的DNA分別混合，構建雌、雄2個基因組DNA樣品池，將其序列分別隨機打斷，電泳回收長度為300～500 bp的DNA片段，兩端加上接頭后進行cluster制備，構建重測序文庫，在Illumina Hiseq 4000平臺進行測序(上海凌恩生物科技有限公司)。將重測序獲得的原始測序數據(reads)，通過Trimmomatic(https://github.com/timflutre/trimmomatic)軟件進行過濾和質量檢測，得到待分析測序數據(clean reads)。

1.3.3 性別標記篩選關聯分析與差異片段篩選

使用bwa(https://github.com/lh3/bwa)和samtools(https://github.com/samtools/samtools)將重測序獲得的clean reads與NCBI中的白斑狗魚參考基因組序列(NCBI登錄號：GCA_000721915.3)進行比對和排序。將clean reads和參考基因組使用bwa進行比對，再用samtools將比對后的基因組文件進行格式轉換和排序，接著使用GATK(https://github.com/broadinstitute/gatk)軟件包對上一步的結果文件進行SNP和Indel檢測，通過MarkDuplicates進行去重以及HaplotypeCaller分析，得到SNP和Indel的數據集，將比對質量值大于20的SNP位點作為最終SNP位點，并對得到的每個SNP和Indel位點進行注釋。使用SNP-Index算法進行標記與性別之間的關聯分析，篩選性別特異SNP位點，并對性別相關侯選區域進行選擇，在這些候選區域中進行候選基因的注釋，隨后在STRING(https://string-db.org)網站對這些基因進行功能互作分析。

將比對到參考基因組上的雌、雄基因組序列進行相互比對，找出差異序列作為性別特異候選序列。同時將不能比對到白斑狗魚參考基因組上的雌、雄基因組序列使用SOAPdenovo(http://soap.genomics.org.cn)進行拼接，然后將拼接的雌、雄序列進行相互比對，排除共有序列后，剩余序列也作為該性別特異候選序列，最后逐一將這些特異序列使用BLAST進行比對。

1.3.4 性別特異候選序列的驗證

隨機挑取部分雌、雄候選特異序列，使用primer 5.0設計PCR引物，以雌、雄個體基因組DNA為模板，用PCR擴增進行兩個階段的驗證。首先用4尾雌魚和4尾雄魚的DNA進行驗證驗，再根據驗證結果，用24尾雌魚和24尾雄魚的DNA進行驗證。PCR反應體系如下：DNA模板1 μL，PCR反應預混液6.25 μL，正反引物各0.5 μL(10 μmol/L)，加ddH2O使體積達到12.5 μL。PCR反應程序為：95 ℃預變性5 min；95 ℃變性30 s，57 ℃退火30 s，72 ℃延伸1 min，35個循環；72 ℃再延伸5 min，最后4 ℃保存。用1.5%瓊脂糖凝膠電泳檢測PCR擴增產物。

2 結果

2.1 混池重測序數據質量評估

雌、雄樣本混池后進行重測序得到的原始測序讀長數目分別為304 730 496條和305 121 064條；堿基數分別為45 709 574 400 bp和45 768 159 600 bp；堿基測序錯誤率(Q30)分別為90.66%和89.39%；測序深度50 X。對reads進行過濾后得到高質量的雌、雄混池數據，clean reads數目分別為282 400 688條和277 391 000條；堿基數分別為41 394 508 708 bp和40 501 950 971 bp；Q30分別提升至94.44%和93.98%，統計結果如表1所示。

表1 白斑狗魚混池重測序數據質量統計Tab.1 Quality statistics of mixed pool resequencing data of E.lucius

2.2 性別特異標記篩選

通過表2可以看到，雌性池SNP位點總計2 076 647個，純合SNP位點11 167個，占0.54%；雌性池Indel位點總計105 810個，缺失突變有54 234個，比插入突變多2 658個，雌性特異候選序列有524條。雄性池SNP位點總計2 076 942個，純合SNP有12 573個，占0.61%；雄性池Indel位點共計647 466個，缺失突變有337 940個，比插入突變多28 414個，雄性特異候選序列500條。篩選得到了10個比對質量值大于20的性別特異SNP位點，其中轉換位點7個，顛換位點3個。同時雌性SNP位點的堿基信息和雄性SNP位點堿基信息具有顯著差異，具體情況如表3所示。

表2 白斑狗魚雌、雄混池SNP和Indel位點以及性別特異候選序列Tab.2 SNP and Indel loci and sex specific candidate sequences in the female and male mixed pool of E.lucius

表3 重測序獲得的白斑狗魚性別特異SNP位點Tab.3 Sex specific SNP loci of E.lucius obtained by resequencing

2.3 性別特異候選基因分析

使用SNP-Index算法獲得Δ(SNP-Index)的分布和局部加權回歸(locally weighted regression，LOESS)擬合曲線(圖1)。經bootstrap方法篩選得到26個差異顯著(P<0.01)的Δ(SNP-Index)區間作為性別特異候選區域，SNP數目介于89～946。26個候選區域中存在不同數目的基因，基因總數為101，詳細信息見表4。經STRING網站分析發現有27個基因存在部分相互作用，這些基因分別在生物體結構發育、細胞增殖和能量代謝等方面發揮作用，其中Pdhx和Prkaa1基因調控主要進行的信號通路(圖2)。

2.4 性別特異片段的PCR實驗驗證

將雌、雄混池數據比對到參考基因組后篩選得到雌性特異序列35條、雄性特異序列14條，未比對到參考基因組數據(unmapping reads)經組裝后篩選得到雌性特異序列489條、雄性特異序列486條。在1 024條性別特異片段中隨機選取了276條可設計引物的片段進行PCR引物設計，其中雌性引物116對，雄性160對。通過對4尾雌魚和4尾雄魚的驗證后，篩選出20對條帶清晰且具有性別偏向性的引物(圖3)。然后又通過對24尾雌魚和24尾雄魚的驗證，最終得到1對引物(編號為2043138)(表5)的擴增產物出現性別差異，雌性陽性比率為4/24，雄性陽性比率為16/24，具有顯著的雄性偏向性(圖4)。將白斑狗魚參考基因組和該序列進行BLAST比對，結果表明該序列不存在于白斑狗魚參考基因組，組裝出來的序列長度為773 bp(圖5)，其中PCR擴增產物長度為384 bp(圖4)。

圖1 白斑狗魚染色體SNP標記的Δ(SNP-Index)分布和LOESS擬合曲線Fig.1 Δ(SNP-Index)distribution and LOESS fitting curve of SNP loci on the chromosome of E.lucius

表4 白斑狗魚染色體性別特異候選區域信息統計Tab.4 Data statistics of sex special chromosome candidate region of E.lucius

續表

注：“染色體區域”中信息依次表示為“染色體編號：起始堿基數-結束堿基數”

圖2 白斑狗魚27個性別特異候選基因的相互作用Fig.2 Interaction of 27 sex special candidate genes of E.lucius藍綠色線-相互關系來源于數據庫；紫色線-相互關系來源于實驗測定；綠色線-預測為基因鄰接關系；紅色線-預測為基因融合；藍色線-預測為基因同現；黃色線-文本挖掘；黑色線-共表達；灰色線-同源蛋白質

圖3 3個性別特異片段在4尾雌性和4尾雄性白斑狗魚的PCR驗證結果Fig.3 PCR result for 3 sex special fragments in 4 female and 4 male E.lucius性別特異片段編號從左向右依次為2019708、2022414、2043138

圖4 性別特異片段編號2043138在24尾雌性和24尾雄性白斑狗魚的PCR擴增結果Fig.4 PCR result for sex special fragment 2043138 in 24 female and 24 male E.lucius

圖5 白斑狗魚性別特異片段編號2043138完整序列Fig.5 Total sequence for sex special fragment 2043138 of E.lucius選中黑色部分為白斑狗魚性別特異片段編號2043138的PCR產物區域。

表5 白斑狗魚性別特異片段(編號2043138)引物設計Tab.5 Primer for sex special fragment 2043138 of E.lucius

3 討論

魚類性別相關研究不僅可以加強魚類性別發育機制的基礎研究，對于水產養殖業來說，進行魚類性別控制也具有重要的應用價值[18]。但魚類是低等脊椎動物，性染色體分化程度較低，多數魚類性染色體間差異非常小，通過核型分析往往很難鑒定魚類性別。由于魚類性別決定方式非常復雜，而且環境因素也常作為重要影響因子，和遺傳因子共同作用影響性別的形成過程，甚至出現生理性別與遺傳性別不一致的情況[19]。通過篩選魚類基因組中性別特有DNA片段或位點，開發出性別分子標記進行魚類性別鑒定成為一種簡便快速的方法。隨著基因組測序技術的發展，開發出混合群體分離分析(又稱混池重測序，Bulked Segregant Analysis，BSA)技術，為魚類性別研究提供了一種有效方法[20]。林曉煜等[14]結合混池重測序和個體重測序，篩選出一個大黃魚(Larimichthyscrocea)雄性特異SNP位點，Lin等[16]通過混池重測序和個體重測序，在大黃魚基因組鑒定出雄性特異的15 bp缺失片段。

本研究利用BSA方法得到雌、雄混池數據，篩選得到26個性別特異候選區域，包含101個性別相關的基因，有27個基因之間存在不同程度的相互作用。其中Pdhx(pyruvate dehydrogenase complex component X，丙酮酸脫氫酶復合物X)和Prkaa1(AMP-activated protein kinase catalytic subunit alpha-1，AMP激活蛋白激酶α1基因)與其他基因互作聯系較多，Pdhx基因[21]催化丙酮酸不可逆脫羧作用生成乙酰輔酶a，同時形成NADH，而且丙酮酸的還原氧化還與增強細胞增殖能力有關。Prkaa1基因[22]參與機體的能量代謝過程，在能量平衡的調節過程中起重要作用，且Prkaa1多態性與肌內脂肪含量顯著相關。而在白斑狗魚的雌、雄個體之間存在性別二型性，個體生長速率與能量代謝和脂肪含量之間可能存在性別差異，所以推測Pdhx基因和Prkaa1基因通過調控細胞增殖和能量代謝的過程與個體性別發育存在一定的關聯。

經過雌、雄混池序列分析，以及PCR擴增的實驗驗證，獲得1條源于unmapping reads具有雄性偏向性的序列，全長773 bp。基因組比對產生unmapping reads的主要原因是由于基因組結構變異、堿基數目大于組裝軟件所允許的錯配數目和gap數目、測序誤差、樣本污染，同時非模式生物的基因組還可能受到核參考基因組質量差、不完全共生體、細胞器基因組的影響，由此unmapping reads仍可能有新的信息發現，例如它們可能提供有關病原體、共生體或參考基因組中缺失的序列/基因的信息[23-25]。將該序列經BLAST后未能獲得與其他物種數據庫中相匹配的序列，也未能在Pan等[26]上傳的白斑狗魚參考基因組數據庫中搜索得到該序列甚至是相似序列，僅和鯽(Carassiusauratus)、恐龍魚(Erpetoichthyscalabaricus)、金線鲃(Sinocyclocheilusrhinocerous)的uncharacterized protein K02 A2.6-like具有61.11%、63.24%、65.75%的相似度，和南亞野鯪(Labeorohita)的voltage-dependent L-type calcium channel subunit alpha-1D protein具有63.01%的相似度。同時該序列在雌、雄各15尾白斑狗魚個體重測序數據中經BLAST比對后，結果顯示該序列似乎并不存在于某一特定染色體，而是對應于不同個體的不同染色體(結果另文發表)，其生物學功能還需要進一步研究。

本研究通過對白斑狗魚雌、雄混池基因組的測序和篩選分析，對unmapping reads組裝后篩選獲得了1條雄性偏向性序列，可為以后進一步對白斑狗魚的性別發育機制研究提供一定理論基礎。