三角帆蚌育珠蚌群体生长性状和转录组差异分析

张爱菊，张根芳，顾志敏，周志明

( 1.浙江省淡水水产研究所，中国水产科学研究院东海水产研究所浙江研究中心，农业农村部淡水渔业健康养殖重点实验室，浙江省淡水水产遗传育种重点实验室，浙江 湖州 313001；2.金华职业技术学院，浙江 金华 321000 )

三角帆蚌(Hyriopsiscumingii)广泛分布于长江、淮河等流域，是我国特有的河蚌资源，又是淡水珍珠养殖的优良品种。然而, 由于现有的供片蚌和育珠蚌质量参差不齐以及养殖群体种质退化等问题日益突出，使得淡水珍珠产业可持续发展受到制约。因此，开展三角帆蚌的遗传改良和品种选育对淡水珍珠产业有着极为重要的意义。迄今为止，国内关于三角帆蚌不同地理种群间形态学及其育珠性状的比较研究等方面已有较多报道[1-4]，均表明三角帆蚌不同地理种群的形态及育珠相关性状差异明显。然而，不同三角帆蚌育珠蚌选育群体形态及育珠相关性状的差异分析尚未有报道。

转录组学是从RNA水平上研究基因表达情况的学科[5]，它将基因组学研究带入了高速发展的崭新时代[6]。转录组学可用于分析不同组织或生理状态下某些基因表达水平的差异，从而发掘与特定生理功能相关的未知基因，具有信噪比高、分辨率高、应用范围广、价格低廉等优势。目前，用于第二代高通量测序的技术主要有3种：Roche公司的454测序、Illumina公司的Solexa测序以及ABI公司的SOLiD测序[7]。目前，软体动物中大多数研究主要集中于利用Roche-454测序或Illumina测序技术研究海洋物种的转录组学研究[8-10]，而在淡水育珠贝类中，仅见利用Roche-454测序技术研究紫色和白色珍珠层三角帆蚌组织中转录组表达差异[11]和利用Illumina测序技术研究三角帆蚌组织总RNA转录组[12]的报道。

为探究3个已经过2代选育无显著形态学差异的三角帆蚌育珠蚌选育群体——鹰蚌、超长蚌和纵纹蚌群体幼蚌之间的差异，笔者运用聚类分析、主成分分析两种多元统计方法比较分析这3个群体幼蚌的生长性状指标，并利用Illumina测序技术比较鹰蚌和超长蚌群体肝胰腺组织转录组表达的差异情况，以期筛选出具有明显生长优势、适合作为插核或插片的育珠蚌群体，为三角帆蚌优良育珠品种的选育和育珠实践提供科学指导。

1 材料与方法

1.1 试验材料

以鹰蚌、超长蚌、纵纹蚌3个三角帆蚌幼蚌选育群体为研究对象，于2014年10月在同一养殖池塘中随机抽取这3个选育群体各40只，暂养于试验池塘，2 d后进行测量。取样时，暂养池内水温23 ℃、pH 7.69、透明度30 cm、亚硝酸氮质量浓度4.98 mg/L、氨氮质量浓度1.38 mg/L。2015年3月再次于同一养殖池塘取样，随后暂养，方法同上所述。此时暂养池内水温10 ℃、pH 7.64、透明度30 cm、亚硝酸氮质量浓度3.05 mg/L、氨氮质量浓度1.16 mg/L。两次取样的幼蚌均为2014年春季同批繁殖，且均取自浙江金华威旺养殖新技术有限公司试验基地。3个群体均从良种场的自然养殖群体中通过家系选育的方法筛选而出，且已经过2代选育。

1.2 生长性状差异分析

1.2.1 数据测量

形态学测量参数为壳长、壳宽、壳高等11个可量性状，采用游标卡尺等测量工具(精确至0.01 mm)，测量方法参照文献[1]的方法，测量部位与参数见图1。测量性状有：壳长(AB)、壳高(OH)、壳宽(WI)、全高(FG)，以及壳顶至前端长(OA)、壳顶至后端长(OB)、壳顶至铰合部前缘长(OC)、壳顶至铰合部后缘长(OD)、壳顶至壳顶上方突起长(OE)、帆状顶点至前端长(FA)、帆状顶点至后端长(FB)。质量参数采用电子天平等工具测量(精确至0.1 g)，为湿质量、壳质量、内脏团质量等3个参数。

图1 形态学测量位点

1.2.2 数据计算

为消除蚌体规格大小对参数值的影响，将每只蚌的所有形态学参数分别除以壳长(AB)予以校正，得出10个形态学比例性状。同时，参考文献[4，13]的方法计算体质量指数(IBW)、壳质量指数(ISW)、肥满度(CF)和壳宽指数(IST)。按下式计算每个选育群体的绝对生长率和相对生长率：

ρ1=(L2-L1)/t

ρ2=(L2-L1)/L1/t×100%

式中，ρ1为绝对生长率(cm/月或g/月)，ρ2为相对生长率(%/月)，L1为第一次采样时壳长、壳宽、壳高或体质量(cm或g)，L2为第二次采样时壳长、壳宽、壳高或体质量(cm或g)，t为采样间隔月数(月)。

1.2.3 数据处理

应用Excel 2007与SPSS 16.0软件进行数据处理。数据表示为平均值±标准差，以P<0.05为差异显著，P<0.01为差异极显著。对10个形态学比例性状参数与4个指数参数进行统计学分析。采用的多元统计方法主要为聚类分析、主成分分析。其中，聚类方法为欧式距离的最短系统聚类法[1,14]。

1.3 转录组差异表达分析

1.3.1 组织样本获取

形态学测量参数和质量参数测量完成后，挑选2015年选取的批次中壳体完整、活力好的鹰蚌和超长蚌2个群体幼蚌各9只，清除其体表附着生物及污物，活体解剖，迅速取肝胰腺组织，用液氮速冻后分别放在对应的群体冻存管中，存于-80 ℃冰箱中备用。

1.3.2 总RNA的提取

试验用枪头在DEPC-H2O中过夜后高压灭菌，置于100 ℃烘箱中烘干备用。采用TRIzol总RNA提取试剂(Invitrogen, USA) 分别提取2个选育群体肝胰腺样本的总RNA。

1.3.3 转录组测序

参照文献[12]的方法，由上海美吉生物医药科技有限公司提供测序服务。

1.3.4 测序数据处理及分析

Illumina HiseqTM2000测序得到的原始图像数据处理参照文献[12]的方法进行质控后，使用TopHat软件分别将鹰蚌组和超长蚌组质控后得到的clean reads与文献[12]的de novo组装结果比对和组装。组装好的结果按照文献[12]的方法进行功能注释与功能分类。

根据每个基因/转录本在样本中的FPKM值，分析差异表达分析。对筛选出的差异表达基因，使用EdgeR软件对所有基因在各组样本中的表达进行差异显著性分析，筛选标准为：FDR(False Discovery Rate)<0.05，并且︱logFC︱≥1。使用Goatools软件进行GO富集分析，使用KOBAS软件进行KEGG富集分析,满足条件FDR<0.05所对应的GO条目和KEGG代谢通路即为显著富集的GO条目代谢通路和KEGG代谢通路。

2 结 果

2.1 3个选育群体间生长性状比较

2.1.1 性状参数比较

方差分析结果显示，壳顶至铰合部后缘长/壳长、壳顶至壳顶上方突起长/壳长在3个群体间均无显著性差异(P>0.05)；鹰蚌群体的壳宽/壳长为(0.23±0.003)，显著高于超长蚌的(0.20±0.018)和纵纹蚌的(0.21±0.003)(P<0.05)，超长蚌与纵纹蚌壳宽/壳长之间差异不显著(P>0.05)；壳高/壳长、全高/壳长、壳顶至前端长/壳长、壳顶至后端长/壳长、壳顶至铰合部前缘长/壳长、帆状顶点至前端长/壳长在3个群体间差异显著(P<0.05)，均呈现纵纹蚌>鹰蚌>超长蚌的趋势；鹰蚌群体帆状顶点至后端长/壳长显著高于纵纹蚌和超长蚌(P<0.05)。综合分析可知， 鹰蚌群体偏宽型，纵纹蚌次之；纵纹蚌群体偏高型，鹰蚌次之；超长蚌偏长型，鹰蚌次之。

分析结果显示，壳质量指数、壳宽指数、体质量指数、肥满度在3个群体间差异显著(P<0.05)，鹰蚌群体的壳质量指数(21.26±0.402)和体质量指数(68.09±0.583)均低于超长蚌的(28.59±0.895，74.76±0.991)和纵纹蚌的(22.19±0.383，68.27±0.638)群体的(P>0.05)，而其壳宽指数(0.13±0.001)显著高于超长蚌的(0.12±0.001)和纵纹蚌的(0.12±0.001)(P<0.05)，且其肥满度(36.08±0.481)高于超长蚌的(33.76±0.577)和纵纹蚌的(32.15±0.527)群体(P>0.05)。故推测鹰蚌群体具有相对较肥满的内脏团，且在壳宽生长上具有明显优势，可能为更适合内脏囊插核的三角帆蚌养殖群体。

2.1.2 聚类分析

对所有样本14个参数进行聚类分析得出(图2)，鹰蚌群体与纵纹蚌群体首先聚在一起，表示鹰蚌与纵纹蚌的距离较短，在形态上可能更为接近。

图2 3个三角帆蚌选育群体聚类分析

2.1.3 主成分分析

主成分分析结果显示，只有前4个特征根大于1，因此SPSS只提取了前4个主成分。第一主成分的总方差贡献率为28.2%，第二主成分为16.74%，第三主成分为16.68%，第四主成分为8.61%，贡献率累积达到70.24%(表1)。4个主成分的累计贡献率较高，基本可以反映不同选育群体之间的生长差异。由表1可见，主成分1在体质量指数、壳质量指数和壳高/壳长上有较大的载荷，主要从体质量和壳高上反映形态和育珠性能；主成分2在壳宽指数、壳宽/壳长上有较大载荷，从壳宽方面反映形态和育珠性能。主成分3在帆状顶点至前端长/壳长上有较大载荷，表现为三角帆蚌纵向的性状方面；主成分4在全高/壳长、壳顶至后端长/壳长上有较大载荷，也主要表现为三角帆蚌纵向的性状方面。

表1 4个主成分的载荷值和方差贡献率

注：*表示载荷值>0.700.

3个育珠蚌选育群体第1、2主成分的散布情况见图3。由图3可见，3个群体均有较多的重叠区域，尤其是鹰蚌群体和纵纹蚌群体，进一步验证其形态更为相近。

图3 3个三角帆蚌选育群体第1、2主成分散布

2.1.4 生长速率比较

超长蚌群体和鹰蚌群体的壳长绝对生长率显著高于纵纹蚌群体(P<0.05)，但鹰蚌群体和超长蚌群体之间差异不显著(P>0.05)；鹰蚌群体的壳宽绝对生长率显著高于其他2个群体(P<0.05)，而其壳高绝对生长率显著高于纵纹蚌群体(P<0.05)，但与超长蚌群体差异不显著(P>0.05)；鹰蚌群体和超长蚌群体的体质量绝对生长率极显著高于纵纹蚌群体(P<0.01)(图4)。同时，3个选育群体的壳长、壳宽、壳高和体质量的相对生长率呈现与对应参数的绝对生长率相同的趋势(图5)。因此，综合分析认为，鹰蚌群体在壳宽和壳高生长上具有明显优势，超长蚌群体在壳长生长方面具有明显优势。

2.2 转录组差异分析

2.2.1 测序产量统计及质量评价

测序共计获得61.92×106条原始序列，碱基数达到62.5×108bp，近94.20%的序列测序质量值均在Q20以上。质控后获得58.09×106条高质量序列，碱基数56.79×108bp，近98.48%的序列测序质量值均在Q20以上。

2.2.2 de novo 序列组装结果

利用Trinity软件进行de novo拼接，将58.09×106条高质量序列拼接成92 347条非冗余unigene。unigene长度为351～28 861 bp，平均长度1150.61 bp。

图4 3个三角帆蚌选育群体的绝对生长率变化Tukey Test分析，柱状图上同个参数中不同小写字母的表示群体组间差异显著(P<0.05)，不同大写字母的表示群体组间差异极显著(P<0.01)，下同.

图5 3个三角帆蚌选育群体的相对生长率变化

2.2.3 unigene功能注释

利用GO数据库，可以将基因通过GO功能分类，共有11 174条unigenes获得了GO信息，归为生物学过程，分子功能和细胞组件3个分支。细胞组件分支中，归类到细胞组分和细胞条目的unigenes最多，均为5181条，分别占22.47%。其次为细胞器和细胞器组分，分别为3670条(15.92 %)和2168条(9.40%)；分子功能分支中，结合功能的unigenes最多，5772条(45.27%)，其次为催化活性(5262条，41.27%)和转运活性(494条，3.87%)；生物学过程分支中，归为细胞生理过程的unigenes最多(17.78%)，代谢过程次之，占14.82%。

2.2.4 差异表达基因的筛选

鉴于鹰蚌群体和超长蚌群体明显的壳体生长优势，主要对这两个群体的转录组差异表达进行阐述。以超长蚌每个基因在样本中的FPKM值为对照，共筛选出差异表达unigenes 63 098条，其中显著性差异表达基因388个。与超长蚌群体相比，鹰蚌群体中显著上调的基因有244个，显著下调的基因有144个，分析显示，大部分基因与代谢、免疫功能等生物学过程相关。

2.2.5 GO功能富集分析

GO功能富集分析可以给出两样本间差异表达基因所参与的主要生物学过程。通过GO功能富集分析，共有3571个差异表达基因富集上了140个GO号，其中显著富集的GO号80个，归类到生物学过程类别的GO号41个，分子功能类别的GO号20个，归类到细胞组分类别的GO号19个。分析结果显示，显著富集的GO条目主要参与各类物质代谢过程的调节、生物合成过程的调节以及基因表达的调节等生物学过程。

2.2.6 KEGG通路富集分析

KEGG代谢通路富集分析可以给出两样本间差异表达基因所参与的最主要代谢通路。本次研究结果显示，2个选育群体中的差异表达基因共富集上了241条通路，其中注释到癌症通路中差异表达基因最多，共有46个；其次为泛蛋白介导的蛋白酶解通路、HTLV-I感染通路和弓形虫(Toxoplasmagondii)病通路，差异表达基因分别为45、40、38个。KEGG 通路富集分析结果显示，共有2条显著富集的通路，为半胱氨酸和蛋氨酸代谢通路和柠檬酸循环通路，注释到前者通路中差异表达基因为 12个，而注释到后者通路中差异表达基因为13个。进一步分析后得出，这2条显著富集的通路分别属于代谢通路中的氨基酸代谢和糖代谢通路。

3 讨 论

3.1 3个群体幼蚌生长性状差异比较分析

钱荣华等[1]利用多元分析方法对我国5大淡水湖的三角帆蚌进行分析，从形态学上探讨了地理种群间的差异；董志国等[2]运用3种多元分析方法对3个地理种群杂交与自交后的9个群体的10个形态比例性状进行了比较研究；闻海波等[3-4]先后对长江和淮河两个水系的三角帆蚌亲本形态特征及3个地理种群三角帆蚌形态等性状进行了比较。这些研究均表明，三角帆蚌不同地理种群的形态及育珠相关性状差异明显。

本次试验的三角帆蚌是自洞庭湖、鄱阳湖等天然水域中引进的亲本，经过两代选育形成的3个不同家系的群体，在形态上既相似又有一定程度的差异，成熟鹰蚌个体腹缘具有鹰嘴状凹刻，成熟超长蚌个体的壳长/壳高明显大于鹰蚌和纵纹蚌，而纵纹蚌的成熟个体贝壳上则具有显著的纵纹。但这3个群体的幼蚌则无类似明显的形态学特征。为此，采用多元统计分析方法，对其生长性状参数进行了比较分析。聚类分析结果中，鹰蚌群体与纵纹蚌群体距离最短，形态更为接近，而超长蚌群体相对这2个群体形态差异较大。主成分散布图进一步验证了这一结果。通常，形态特征是由遗传因子与环境因子共同作用的结果，本试验的3个三角帆蚌均采自同一地区的同一池塘，基本可以排除环境因子因素。究竟是什么遗传因素导致这些群体之间形态的变异，以及这些群体之间的形态变异究竟多大程度地反应在遗传物质与等位基因的变异上，还有待于进一步的深入研究。

3.2 贝壳形态、生长性状与选育

珍珠质量除受养殖环境因子影响外[15]，还与多种因素相关，而其中育珠蚌亲本选择是一个重要因素。Wada[16]认为，育珠蚌亲本的选择应主要以珍珠层颜色、贝壳形态及生长等指标为主。白志毅等[17]认为，三角帆蚌选育应以壳质量(或体质量)和壳宽为主要选育目标，李梦军等[18]群体选育结果也证明了这一观点。由于本研究的目的是筛选出具有明显生长优势、适合作为插核或插片的育珠蚌群体，不考虑珍珠颜色的选育，故认为三角帆蚌亲本选育应以贝壳形态及生长等指标为主。鉴于鹰蚌和超长蚌群体幼蚌较快的生长速率，主要考虑贝壳形态方面的因素。本研究对3个三角帆蚌选育群体幼蚌的贝壳形态测定与分析结果表明，鹰蚌幼蚌的壳宽指数和肥满度均高于超长蚌和纵纹蚌，因此鹰蚌群体的育珠性能优于其他2个群体，推断其作为育珠蚌更适合内脏囊插核，将有利于大颗粒优质有核珍珠的培育。同时，从贝壳形态上看，超长蚌为长扁型，在壳长生长方面具有明显优势，作为育珠蚌更适合外套膜插片，可以增加插入的细胞小片数量，有利于无核珍珠的培育。

3.3 高通量测序技术与选育

近年来，运用高通量测序技术辅助水产经济物种的遗传选育工作在我国越来越受到重视。目前，对于三角帆蚌的转录水平研究主要采用cDNA文库/EST测序技术[19]、cDNA文库/抑制消减杂交技术[20]、cDNA文库/Roche测序技术[11]、Illumina/Solexa测序技术[12]，但这些研究主要集中于三角帆蚌整体水平或者珍珠颜色选育等方面，而关于有核珍珠/无核珍珠育珠蚌选育方面的研究略有匮乏。

目前，软体动物外套膜转录组高通量测序的研究较多，且主要用于生物矿化分析[8,11]。然而，鉴于本研究的目的以及鹰蚌和超长蚌各自的壳体生长优势，同时考虑肝胰腺在软体动物众多生命大分子的合成及分解代谢等方面的重要功能[21]，很大程度上影响其生长发育等因素，本研究进一步应用高通量测序技术对同龄鹰蚌和超长蚌2个选育群体的幼蚌肝胰腺转录组进行深度测序。结果发现，筛选出差异表达unigenes 63 098条，其中显著性差异表达基因388个，显著富集的GO条目主要参与各类物质代谢过程的调节、生物合成过程的调节以及基因表达的调节等生物学过程，显著富集的KEGG通路为半胱氨酸和蛋氨酸代谢通路和柠檬酸循环通路，均属于代谢通路。故推测这2条代谢通路及其中的差异表达基因可能是造成鹰蚌和超长蚌2个选育群体生长性能差距的主要遗传因素。

以上测序结果和分析为三角帆蚌的遗传育种提供了丰富的分子生物学信息，有利于今后的家系选育工作。然而，这些结果仅是在转录组水平生物信息学分析上的一种推测，今后仍需进一步进行验证。