Analyses of differentially expressed genes of Eucalyptus urophylla in response to single and dual factor stresses of water and nutrient
-
摘要:目的
比较尾叶桉Eucalyptus urophylla对水肥单、双因素胁迫的分子响应机制,为桉树抗逆育种提供思路和分子基础。
方法以无性系ZQUA44为试验材料,设3个水肥胁迫处理和1个对照(CK)。水肥双因素胁迫:20%~40%田间持水量,基肥[钙镁磷肥(CMPF) 250 g];水分单因素胁迫:20%~40%田间持水量,基肥(CMPF 250 g+复合肥150 g);养分单因素胁迫:60%~80%田间持水量,基肥(CMPF 250 g);CK:60%~80%田间持水量,基肥(CMPF 250 g+复合肥150 g)。水分单因素胁迫和对照在种植2个月后施尿素100 g,8月份施复合肥150 g,第2年春天施复合肥100 g。对受胁迫处理的尾叶桉进行生长指标测定,并取叶片进行转录组测序,对测序结果进行拼接,利用生物信息学手段获得差异表达基因,对差异基因进行GO功能注释和KEGG通路分析。最后随机选择4个具有不同表达模式的差异表达基因进行qRT-PCR,验证转录组测序数据的可信度。
结果3种胁迫处理均在一定程度上抑制尾叶桉的生长发育,胁迫处理下尾叶桉各项生长指标多显著低于对照。水肥双因素胁迫、水分单因素胁迫和养分单因素胁迫处理共得到5 547个差异表达基因,分别有2 585、1 472和1 490个差异表达基因,各有1 195、222和665个基因表达上调和1 390、1 250和825个基因表达下调;3种胁迫处理分别获得155、75和108个基因编码转录因子。水肥单、双因素处理的差异表达基因均显著富集在苯丙烷的生物合成路径中。qRT-PCR结果与转录组测序结果基本一致,说明测序结果可信。
结论尾叶桉在3种胁迫处理下的转录变化具有明显的普遍性和特异性,综合胁迫对植物生长的影响比单独胁迫的影响更大,不同胁迫条件下植物可能具有不同的响应方式。
Abstract:ObjectiveTo compare the molecular response mechanisms of Eucalyptus urophylla under single and dual factor stresses of water and nutrient, and provide insight and molecular basis for resistance breeding of Eucalyptus.
MethodThe experiment selected E. urophylla clone ZQUA44 as matetial, and set three stress treatments and one control (CK). Water and nutrient deficiency stress: 20%−40% field moisture capacity, applying 250 g calcium magnesium phosphate fertilizer (CMPF) as base fertilizer; Water deficiency stress: 20%−40% field moisture capacity, applying 250 g CMPF and 150 g compound fertilizer as base fertilizer; Nutrient deficiency stress: 60%−80% field moisture capacity, applying 250 g CMPF as base fertilizer; CK: 60%−80% field moisture capacity, applying 250 g CMPF and 150 g compound fertilizer as base fertilizer. The urea of 100 g was applied at two months after planting, 150 g compound fertilizer in August, and 100 g compound fertilizer in the next spring in water deficiency stress and control treatments. High-throughput transcriptome sequencing was performed to E. urophylla leaves in different stress treatments, and the growth indexes of plant were determined. After assembling, differentially expressed genes were obtained using bioinformatics methods, and GO functional annotations and KEGG pathway analyses were performed. Finally, four differentially expressed genes with different expression patterns were randomly selected for qRT-PCR to validate the reliability of transcripcome sequencing data.
ResultAll three stress treatments partly suppressed the growth and development ofE. urophylla. Most of the growth indexes of E. urophylla in three stress treatments were significantly lower than those in control. A total of 5 547 differentially expressed genes were obtained in water and nutrient deficiency treatment, water deficiency treatment and nutrient deficiency treatment. The 2 585, 1 472 and 1 490 differentially expressed genes were obtained, including 1 195, 222 and 665 up-regulated genes and 1 390, 1 250 and 825 down-regulated genes, respectively. The 155, 75, 108 gene encoding transcription factors were respectively obtained in three stress treatments. The differentially expressed genes in three stress treatments were all significantly enriched in phenylpropane biosynthetic pathway. The qRT-PCR results were basically consistent with transcriptome sequencing results, indicating that the sequencing results were credible.
ConclusionThe transcriptional changes of E. urophylla in three stress treatments have obvious universality and specificity, and the effects of dual factor stress on plant growth are more than that of single factor stress. Plant may have different responce modes in different stress conditions.
-
林木遗传改良通常需在田间试验基础上,估算参试样本(种源、家系、单株或无性系)目标性状的遗传参数,并运用一定的评价方法筛选出具有优异目标性状的种源、家系或无性系[1]。目前,林木遗传评估基本采用混合线性模型,即根据观测表型值和设计的关联矩阵,利用最佳线性无偏估计(Best linear unbiased estimation, BLUE)方法获得固定效应估计值,利用最佳线性无偏预测(Best linear unbiased prediction, BLUP)方法获得随机效应预测值。而育种值(Breeding value)即属随机效应预测值,且BLUP值依赖于混合线性模型的关联设计矩阵。然而,对林木田间试验而言,个体间存在竞争,且基因型与环境互作较为显著,即空间效应明显。如竞争效应和空间效应没有被考虑到混合线性模型中,遗传评估则会产生偏差,其中空间正相关易导致加性遗传方差偏大,中等竞争效应则抑制偏差[2]。因而,在林木遗传评估模型中,应纳入竞争效应和空间效应。
空间效应(或称环境异质性)可以是空间连续的,其表现为土壤和小气候效应的相似模式,也可以是不连续的,表现为育林措施和测量方法的不同效应,此外还可以是随机的,其表现为微环境的异质性[3]。针对空间效应的空间分析(Spatial analysis)已被广泛用于农作物品种田间试验研究。与传统试验设计相比,空间分析可减少试验误差并提高农作物品种评估的准确性。空间分析方法通常有趋势面分析法和空间自相关分析法等。其中,趋势面分析法是采用空间坐标的多项式函数来模拟环境的变异[4],而空间自相关分析法则通过行和列的空间自相关分析剔除地点内的空间相关误差来提高遗传试验分析的精确性[5]。
在林木研究上,Costa等[6]利用空间自相关分析法分析了辐射松Pinus radiata等3个树种的子代和无性系试验林,发现空间自相关分析可更有效拟合试验林的空间变异,显著降低了区组和小区的方差,且极大提高了遗传评估的准确度与遗传增益。另外,Dutkowski等[7-8]在更大数据量的树种与试验林上进一步采用空间分析,也发现空间分析可显著提高种源、家系、亲本和无性系遗传效应估算的准确性。但林木是多年生植物,除空间差异外,往往还存在个体竞争[9]。虽然Cappa等[10]采用遗传干预模型(Genotypic interference model, GIM)评估了林木遗传加性效应和竞争加性效应,但是没有考虑空间效应,因此模型比较受限。值得注意的是,Stringer等[11]提出空间效应与竞争效应的联合模型进行糖蔗品种试验的遗传分析,模型采用等根二阶自回归法(Equal-roots second-order autoregressive, EAR2)拟合环境的空间效应,而遗传的竞争效应则采用处理干预模型(Treatment interference model, TIM)。此外,Cappa等[12]提出拟合空间效应与竞争效应的林木单株混合模型,并采用贝叶斯法估计遗传参数。上述模型均在一定程度上提高了林木遗传评估的准确性,但模型涵盖的遗传竞争效应均为单一值,而实际上林木个体周边的遗传竞争效应可能不同,因此有必要提出新的单株混合模型,同时拟合不同的遗传竞争效应和异质的空间效应。基于此,本研究根据现代线性统计理论,构建了异质空间效应和不同竞争效应的联合模型,并通过数据模拟的方式,结合实际测定数据进行验证,同时运用REML 法估算遗传参数,对不同模型进行比较,进而探讨异质空间效应和不同竞争效应的分析模型在林木遗传测定中的有效性,以期为林木精确的遗传评估和后续可靠的良种选育提供理论参考和技术支撑。
1. 方法
1.1 统计模型
单株混合模型为:y=Xb+Zgug+e。式中: y为个体性状表型值构成的向量; b为固定效应构成的向量; u为随机遗传效应构成的向量,
${\boldsymbol u}_{\rm{g}}' =\left( {{u}_{{{\rm{g}}_{\rm{d}}}}',{u}_{{{\rm{g}}_{{\rm{en}}}}}',} \right. $ $\left. {u}_{{{\rm{g}}_{{\rm{wn}}}}}', {u}_{{{\rm{g}}_{{\rm{sn}}}}}',{u}_{{{\rm{g}}_{{\rm{nn}}}}}' \right)$ ,${u}_{{{\rm{g}}_{\rm{d}}}}'\text{、} {u}_{{{\rm{g}}_{{\rm{en}}}}}'\text{、} {u}_{{{\rm{g}}_{{\rm{wn}}}}}'\text{、} {u}_{{{\rm{g}}_{{\rm{sn}}}}}'\text{和 } {u}_{{{\rm{g}}_{{\rm{nn}}}}}'$ 分别为直接遗传效应、东向近邻效应、西向近邻效应、南向近邻效应和北向近邻效应;e为剩余误差向量;X为固定效应的关联矩阵;Zg为随机效应的关联矩阵。遗传效应的方差–协方差矩阵为:
$$ {\rm{var}}\left( {{\boldsymbol{u}_{\rm{g}}}} \right) = \left[ \begin{array}{l} \begin{array}{*{20}{c}} {\delta _{{{\rm{g}}_{\rm{d}}}}^2} & {{\delta _{{{\rm{g}}_{{\rm{den}}}}}}} & {{\delta _{{{\rm{g}}_{{\rm{dwn}}}}}}}\\ {{\delta _{{{\rm{g}}_{{\rm{den}}}}}}} & {\delta _{{{\rm{g}}_{{\rm{en}}}}}^2} & {{\delta _{{{\rm{g}}_{{\rm{ewn}}}}}}}\\ {{\delta _{{{\rm{g}}_{{\rm{dwn}}}}}}} & {{\delta _{{{\rm{g}}_{{\rm{ewn}}}}}}} & {\delta _{{{\rm{g}}_{{\rm{wn}}}}}^2} \end{array}\begin{array}{*{20}{c}} {{\delta _{{{\rm{g}}_{{\rm{dsn}}}}}}} & {{\delta _{{{\rm{g}}_{{\rm{dnn}}}}}}}\\ {{\delta _{{{\rm{g}}_{{\rm{esn}}}}}}} & {{\delta _{{{\rm{g}}_{{\rm{enn}}}}}}}\\ {{\delta _{{{\rm{g}}_{{\rm{wsn}}}}}}} & {{\delta _{{{\rm{g}}_{{\rm{wnn}}}}}}} \end{array}\\ \begin{array}{*{20}{c}} {{\delta _{{{\rm{g}}_{{\rm{dsn}}}}}}} & {{\delta _{{{\rm{g}}_{{\rm{esn}}}}}}} & {{\delta _{{{\rm{g}}_{{\rm{wsn}}}}}}}\\ {{\delta _{{{\rm{g}}_{{\rm{dnn}}}}}}} & {{\delta _{{{\rm{g}}_{{\rm{enn}}}}}}} & {{\delta _{{{\rm{g}}_{{\rm{wnn}}}}}}} \end{array}\begin{array}{*{20}{c}} {\delta _{{{\rm{g}}_{{\rm{sn}}}}}^2} & {{\delta _{{{\rm{g}}_{{\rm{snn}}}}}}}\\ {{\delta _{{{\rm{g}}_{{\rm{snn}}}}}}} & {\delta _{{{\rm{g}}_{{\rm{nn}}}}}^2} \end{array} \end{array} \right] \otimes {\boldsymbol{I}_m}, $$ 其中,对角线上依次为直接遗传效应、东向近邻效应、西向近邻效应、南向近邻效应和北向近邻效应的方差分量,非对角线为上述两两效应之间的协方差分量。var代表方差,δ2 和 δ 分别代表方差和协方差,Im代表单位矩阵。
为减少计算量和促进模型收敛,可通过因子分析法[1]将遗传效应的方差–协方差矩阵改写为如下矩阵:
$$ {\rm{var}}\left( \boldsymbol{u} \right) = \left( {\boldsymbol \varGamma \times \boldsymbol\varGamma ' + \boldsymbol\varPsi } \right) \otimes \boldsymbol{A}, $$ 其中,Γ 是因子载荷矩阵,Ψ 是特殊方差矩阵,A是加性遗传相关矩阵。
在常规的空间分析[7]中,剩余误差(R)可分解为空间自相关误差方差(
$\delta _\xi ^2$ )和随机测量误差方差($\delta _\eta ^2$ ),因此剩余误差的方差协方差矩阵为:$$ \boldsymbol{R} = \delta _\xi ^2\left( {\sum\nolimits_c {\rm{(}}{\rho _{_c}}{\rm{)}} \otimes \sum\nolimits_r {\rm{(}}{\rho _{_r}}{\rm{)}}} \right) + \delta _\eta ^2\boldsymbol{I}, $$ 其中,
$\sum\nolimits_c{\rm{(}}{\rho _{_c}}{\rm{)}}$ 、$\sum\nolimits_r{\rm{(}}{\rho _{_r}}{\rm{)}}$ 分别是列、行自回归的相关矩阵,ρc、ρr分别是列、行自相关值,I是单位矩阵。$\mathop \sum \nolimits {\rm{(}}\rho {\rm{)}}$ 自回归的相关矩阵,具体如下:$$ \sum {(\rho )} = \left| {\begin{array}{*{20}{c}} 1 & \rho & {{\rho ^2}} & \cdots & {{\rho ^{n - 1}}}\\ \rho & 1 & \rho & \cdots & \vdots \\ {{\rho ^2}} & \rho & 1 & \cdots & \vdots \\ \vdots & \vdots & \vdots & \ddots & \vdots \\ {{\rho ^{n - 1}}} & \cdots & \cdots & \cdots & 1 \end{array}} \right|\text{。} $$ 1.2 数据模拟方法
采用R软件[1]和R程序包breedR[13]进行数据模拟。为简化数据模拟,参考下述松树数据,以该松树胸径均值(13)和估算的加性方差(
$\delta _a^2$ 为2.46)、剩余残差($\delta _e^2$ 为3.00)作为模拟初始值,4次重复,5株行式小区,10个父本,10个母本,共30行、40列。为简化不同模型的比较,设置加性遗传效应与竞争遗传效应相关值[简称竞争相关值(Cr)]为–0.6, 空间行、列自相关值分别为0.6、–0.6,比较模型包括随机完全区组设计模型(RCBM)、空间模型(SM)、空间与测量误差模型(SUM)和空间与竞争模型(SCM)。为探究遗传竞争效应和环境空间效应的影响,设置竞争相关值(Cr= –0.3, –0.6, –0.9)以及空间行相关值( ρr= 0.3, 0.6, 0.9)、列自相关值(ρc= –0.3, –0.6, –0.9),通过正交设计 [1],产生9组数据。1.3 实测数据
数据来源于文献[1]中的例6.4.4.2,该数据有46个全同胞家系,采用完全随机区组设计,设7个区组,5株行式小区,试验林共35行40列,测定性状是10年生胸径。
1.4 模型分析方法
所有模型分析均采用ASReml软件V4.0[14],除总体均值作为固定因子外,所有其他试验因子均作为随机因子。空间相关的残差采用行列AR1自相关结构,遗传竞争效应的方差采用单因子XFA1结构。拟合模型的原始数据图和残差图采用R包AAfun[15]生成,空间变异图variogram由ASReml运行生成,所有估计的参数由ASReml运行获得。
单株狭义遗传力(
$h_i^2$ )的计算公式如下:$$ h_i^2 = \delta _a^2/\left( {\delta _a^2 + \delta _e^2} \right), $$ 式中,
$\delta _a^2$ 为加性方差,${\delta _e^2}$ 为随机误差方差(剩余残差)。AIC值的计算公式如下:
$$ {\rm AIC} = - 2 \times ({\rm ln}L - P)\text{,} $$ 式中,AIC值是赤池信息量准则,lnL是模型迭代收敛似然值,P是模型中随机参数的数量。
2. 结果与分析
2.1 不同模型比较
对模拟数据设置随机区组设计模型(RCBM)、空间模型(SM)、空间与测量误差模型(SUM)和空间与竞争模型(SCM)4种模型,模型的拟合结果(表1)显示:模型RCBM的参数最少(4个),收敛似然值(lnL)最小(–1 927),AIC值最大(3 863),剩余残差(7.56)也最大,所估算的单株狭义遗传力也最小(0.24±0.08);模型SM的参数为6个,收敛似然值(lnL)稍微增加,为–1 904,AIC值下降至3 820,空间自相关误差为7.04,加性方差为3.44,比RCBM提高了45%,但估计的行列相关值远低于模拟设置值,行自相关值为0.25,列自相关值为–0.07;模型SUM的参数中,随机测量误差方差为5.72,空间自相关误差方差为1.83,行列相关值与模拟设置值一致,加性方差为2.43,估计的单株狭义遗传力为0.30±0.10;模型SCM 参数最多(10个),似然值(lnL)最大(–1 853),但AIC值(3 734)和随机测量误差方差(3.13)均最小,估计的单株狭义遗传力最高,为0.43±0.13。此外,模型SCM估计的4个近邻效应方差在0.34~0.74,并且均显著。根据似然值lnL最大和AIC值最小原则,可知SCM为最佳模型。
表 1 模拟数据不同模型估计参数的比较1)Table 1. Comparison of estimated parameters from different models for simulated data模型 参数/
个lnL AIC $\delta _{{\rm{rep}}}^2$ $\delta _{{\rm{plot}}}^2$ $\delta _e^2$ $\delta _{\rm{\eta }}^2$ $\delta _{\rm{\xi }}^2$ ρr ρc $\delta _a^2$ $\delta _{{{\rm{g}}_{{\rm{en}}}}}^2$ $\delta _{{{\rm{g}}_{{\rm{wn}}}}}^2$ $\delta _{{{\rm{g}}_{{\rm{sn}}}}}^2$ $\delta _{{{\rm{g}}_{{\rm{nn}}}}}^2$ $h_i^2$ RCBM 4 –1 927 3 863 0.060.85 0.101.03 7.5613.52 2.372.67 0.24±0.08 SM 6 –1 904 3 820 0.080.81 0.070.80 7.048.92 0.258.92 –0.075.87 3.442.52 … SUM 7 –1 897 3 808 0.070.89 0.080.84 5.729.02 1.834.21 0.615.84 –0.708.00 2.432.70 0.30±0.10 SCM 10 –1 853 3 734 3.134.74 2.204.75 0.535.24 –0.709.63 2.372.81 0.592.22 0.341.60 0.682.39 0.742.54 0.43±0.13 1) lnL为迭代收敛似然值,AIC 为赤池信息量准则, $\delta _{{\rm{rep}}}^2$为区组方差, $\delta _{{\rm{plot}}}^2$为小区方差, $\delta _e^2$为剩余残差, $\delta _{\rm{\eta }}^2$、 $\delta _{\rm{\xi }}^2$分别为随机测量误差方差和空间自相关误差方差,ρr、ρc 分别为行、列自相关值, $\delta _a^2$为加性方差, $\delta _{{{\rm{g}}_{{\rm{en}}}}}^2$、 $\delta _{{{\rm{g}}_{{\rm{wn}}}}}^2$、 $\delta _{{{\rm{g}}_{{\rm{sn}}}}}^2$、 $\delta _{{{\rm{g}}_{{\rm{nn}}}}}^2$分别为东向、西向、南向和北向近邻遗传竞争方差, $h_i^2$为单株狭义遗传力;表中的上标值是 t 检验统计量,为估计值与标准误的比值;…表示单株遗传力因没有随机测量误差而不能计算。 残差图显示,RCBM的残差图(图1C)与原始数据图(图1A)非常相似,表明RCBM模型未能很好地拟合空间异质性。虽然SM模型比RCBM有所改进,但也只拟合了残差小的空间部分(图1D)。SUM和SCM模型,图形比较相似,有明显的行列空间趋势,但SCM模型更能消除部分残差高的区域(图1E、1F)。从空间变异图可以看出,SM模型的图形比较平坦(图2B),表明捕获的空间变异很小,而SUM和SCM模型的图形比较相似(图2C、2D),在行方向上存在明显的空间变异趋势,在列方向的突起尖峰表示列相关为负值,存在环境的竞争效应,表明SUM和SCM模型捕获空间变异的效果不错。此外,残差值与近邻表型均值的散点图(图2A)显示两者存在明显的负相关(r= –0.30),表明存在遗传竞争效应。上述结果与 表1的结果相一致。
2.2 SCM模型拟合结果的稳定性分析
为了解不同初始参数值条件下,SCM模型拟合的效果,根据竞争相关值、行列自相关值设计了9组数据,结果(表2)表明,所有数据拟合的加性方差为2.25~2.52,随机测量误差方差为2.90~3.96,空间自相关误差方差为1.31~2.43,除了列相关值(pc= –0.3)拟合得一般外,其他均与初始值接近。各近邻效应方差多数在0.45~0.80,基本上都显著;此外,估计的单株狭义遗传力差异不大,均在0.40左右。可见,SCM 模型对于空间效应和竞争效应的拟合效果稳定性较好。
表 2 不同参数初始值的SCM模型估计结果1)Table 2. Estimated results of SCM with different initial values for parameters参数初始值 ρr ρc $\delta _{\rm{\eta }}^2$ $\delta _{\rm{\xi }}^2$ $\delta _a^2$ $\delta _{{{\rm{g}}_{{\rm{en}}}}}^2$ $\delta _{{{\rm{g}}_{{\rm{wn}}}}}^2$ $\delta _{{{\rm{g}}_{{\rm{sn}}}}}^2$ $\delta _{{{\rm{g}}_{{\rm{nn}}}}}^2$ $h_i^2$ Cr–0.3/ρr0.3/ρc–0.3 0.362.24 –0.675.45 3.965.20 1.312.77 2.472.72 0.742.43 0.241.32 0.562.01 0.652.43 0.38±0.12 Cr–0.3/ρr0.6/ρc–0.9 0.586.00 –0.9133.47 3.285.61 2.395.05 2.452.75 0.652.36 0.331.71 0.492.32 0.762.58 0.43±0.13 Cr–0.3/ρr0.9/ρc–0.6 0.9125.41 –0.738.56 3.245.37 1.803.95 2.522.77 0.912.60 0.472.03 0.522.35 0.832.83 0.44±0.13 Cr–0.6/ρr0.3/ρc–0.9 0.403.35 –0.9131.91 3.356.09 2.375.08 2.252.89 0.532.28 0.411.89 0.602.63 0.622.44 0.40±0.12 Cr–0.6/ρr0.6/ρc–0.6 0.535.24 –0.709.63 3.134.74 2.204.75 2.362.81 0.592.22 0.341.60 0.682.39 0.742.54 0.43±0.13 Cr–0.6/ρr0.9/ρc–0.3 0.9440.34 0.010.04 3.214.77 2.213.69 2.273.29 0.452.12 0.451.88 0.802.61 0.512.18 0.41±0.11 Cr–0.9/ρr0.3/ρc–0.6 0.383.04 –0.687.71 3.965.20 1.312.77 2.472.72 0.742.43 0.241.32 0.562.01 0.652.43 0.44±0.08 Cr–0.9/ρr0.6/ρc–0.6 0.545.36 –0.668.38 2.904.56 2.414.51 2.353.24 0.452.10 0.431.87 0.802.58 0.592.33 0.43±0.13 Cr–0.9/ρr0.9/ρc–0.9 0.9120.45 –0.9224.47 3.246.82 2.433.05 2.383.64 0.512.41 0.572.48 0.652.77 0.532.51 0.42±0.10 1) Cr为竞争相关值,ρr、ρc分别为行、列自相关值, $\delta _{\rm{\eta }}^2$、 $\delta _{\rm{\xi }}^2$分别为随机测量误差方差和空间自相关误差方差, $\delta _a^2$为加性方差, $\delta _{{{\rm{g}}_{{\rm{en}}}}}^2$、 $\delta _{{{\rm{g}}_{{\rm{wn}}}}}^2$、 $\delta _{{{\rm{g}}_{{\rm{sn}}}}}^2$、 $\delta _{{{\rm{g}}_{{\rm{nn}}}}}^2$分别为东向、西向、南向和北向近邻遗传竞争方差, $h_i^2$为单株狭义遗传力;表中的上标值是t检验统计量,为估计值与标准误的比值。 2.3 真实数据的模型拟合结果比较
真实数据拟合的结果(表3)显示,常规的RCBM模型,估计的区组方差为0.19,加性方差为2.46,剩余残差为2.96;仅有空间行列自相关的模型SM,加性方差被高估为7.58;加上区组因子后的空间模型SM1,这种偏差被修正,加性方差为2.74,空间自相关误差方差为2.84,但行、列自相关值仅为0.11、0.13;“空间+测量”误差模型SUM,有较好的拟合结果,加性方差为2.21,随机测量误差方差为2.80,空间自相关误差方差为0.80;空间与竞争模型SCM为最优模型,加性方差为2.25,随机测量误差方差和空间自相关误差方差均比SUM有所降低,分别为2.41和0.72,所估算的近邻竞争方差在0.01~0.13,但没达到显著水平。从似然值lnL和AIC值可知,SCM模型仍是最优。
表 3 真实数据不同模型估计参数的比较1)Table 3. Comparison of estimated parameters from different models for real data模型 参数/
个lnL AIC $\delta _{{\rm{rep}}}^2$ $\delta _{{\rm{plot}}}^2$ $\delta _e^2$ $\delta _{\rm{\eta }}^2$ $\delta _{\rm{\xi }}^2$ ρr ρc $\delta _a^2$ $\delta _{{{\rm{g}}_{{\rm{en}}}}}^2$ $\delta _{{{\rm{g}}_{{\rm{wn}}}}}^2$ $\delta _{{{\rm{g}}_{{\rm{sn}}}}}^2$ $\delta _{{{\rm{g}}_{{\rm{nn}}}}}^2$ $h_i^2$ RCBM 4 –1 605 3 219 0.191.54 0.00 2.968.04 2.463.82 0.45±0.09 SM 4 –1 578 3 165 0.816.84 0.9839.29 0.9842.75 7.5824.29 … SM1 6 –1 600 3 212 0.191.49 0.00 2.846.84 0.118.92 0.132.69 2.743.69 … SUM 7 –1 569 3 153 0.010.14 0.00 2.808.37 0.801.26 0.9738.63 0.9839.76 2.213.81 0.44±0.10 SCM 10 –1 559 3 138 2.416.43 0.722.96 0.9120.33 0.9326.38 2.253.75 0.130.75 0.091.13 0.010.32 0.051.00 0.48±0.10 1) lnL为迭代收敛似然值,AIC为赤池信息量准则, $\delta _{{\rm{rep}}}^2$为区组方差, $\delta _{{\rm{plot}}}^2$为小区方差, $\delta _e^2$为剩余残差, $\delta _{\rm{\eta }}^2$、 $\delta _{\rm{\xi }}^2$分别为测量随机误差方差和空间自相关误差方差,ρr、ρc分别为行、列自相关值, $\delta _a^2$为加性方差, $\delta _{{{\rm{g}}_{{\rm{en}}}}}^2$、 $\delta _{{{\rm{g}}_{{\rm{wn}}}}}^2$、 $\delta _{{{\rm{g}}_{{\rm{sn}}}}}^2$、 $\delta _{{{\rm{g}}_{{\rm{nn}}}}}^2$分别为东向、西向、南向和北向近邻遗传竞争方差, $h_i^2$为单株狭义遗传力;表中的上标值是t检验统计量,为估计值与标准误的比值;…表示单株遗传力因没有随机测量误差而不能计算。 3. 讨论与结论
本研究利用R软件及其程序包breedR模拟数据,结合实测数据,采用XFA1结构拟合加性效应、近邻竞争效应和AR1结构拟合空间效应,建立了RCBM、SM、SUM和SCM 4种模型并进行比较,发现SCM是最优模型,该模型可大大降低随机误差方差,并估算到四周近邻竞争方差,设置参数的不同初始值,SCM估计的参数结果均较为稳定。
对林木研究而言,试验林面积一般较大,采用常规试验设计难以避免环境异质性,传统的分析模型因忽略环境的异质性而导致遗传参数和育种值估算的偏差,因此Dutkowski 等[8]提出林木遗传试验数据应尽可能采用空间分析。对于环境的空间变异,一般使用3次B样条函数(Cubic B spline)或行列自回归拟合。考虑到残差结构的灵活性,本文采用行列自回归法进行空间分析。从模拟数据的似然值lnL和AIC值结果可以看出,空间分析模型比常规分析模型拟合的结果更好,但是没有随机误差方差的空间模型,随机误差被纳入空间相关误差,从而导致空间相关误差方差被高估,因此,也造成空间变异趋势被抑制,使得行、列自相关值被低估,这从残差图和空间变异图也可得到佐证。然而,仅有行、列自相关的空间模型会导致加性方差被严重高估[3],本文也得到一致的结果。Magnussen[2] 指出,环境异质性一般存在于试验林早期,林木个体周围相似的微环境会导致个体间存在正相关,但随着林木生长,个体间的竞争日趋严重,将会导致个体间存在负相关。因此,对于多年生的林木试验林,个体间存在竞争效应是必然的,这是造成个体致死率的重要因素[9]。换言之,林木性状的表型中,会同时存在空间效应和竞争效应,忽略任何一个因素,都可能会造成遗传评估的偏差,并进一步影响后续的选择工作。对于遗传竞争效应,拟合的方法有GIM[10]、TIM[11]和CSM[12],但上述研究中遗传竞争效应均为单一值,忽略了林木个体周边的遗传竞争效应可能不同。本文提出的以加性效应及其四周近邻竞争效应的方差–协方差矩阵作为G结构,替代常规的加性G结构,从而达到拟合近邻竞争效应值不唯一的目的。考虑到降低运算量和模型收敛,采用因子分析法来拟合G结构。模拟数据的拟合结果表明,新的SCM模型,拟合方差虽然没有增加,但随机误差方差大大降低了,分别由7.56(RCBM)、5.72(SUM)降低到3.13(SCM),分别降低了58.6%、45.3%,这与同时考虑竞争效应和空间效应的研究模型结果[12]相一致。此外,拟合数据的SCM模型单株狭义遗传力保持在0.40左右,高于RCBM(0.24)和SUM模型(0.30)。而针对参数不同初始值的模拟数据,SCM模型的拟合结果都比较稳定,表明其可用于林木的遗传分析。最后,利用松树的实测数据验证SCM模型,估算结果与模拟结果比较一致,笔者认为,新单株混合模型即SCM模型可用于林木遗传分析,在提高遗传评估的准确性和可靠性方面具有优势。
-
图 1 4种不同水肥处理中尾叶桉无性系生长指标的变化
T1:水肥双因素胁迫,T2:水分单因素胁迫,T3:养分单因素胁迫;各小图相同指标柱子上不同小写字母表示不同处理间差异显著(P < 0.05,Duncan’s法)
Figure 1. Growth index changes of Eucalyptus urophylla clones in four different water and fertilizer treatments
T1: Water and nutrient deficiency stress, T2: Water deficiency stress, T3: Nutrient deficiency stress; Different lowercase letters on columns of the same index in each figure indicate significant differences among different treatments (P < 0.05, Duncan’s method)
图 3 3种胁迫处理中尾叶桉差异表达基因韦恩图
T1:水肥双因素胁迫,T2:水分单因素胁迫,T3:养分单因素胁迫;括号中的百分数表示该基因在所有上调或下调基因中占的比例
Figure 3. Venn diagram of differentially expressed genes of Eucalyptus urophylla in three stress treatments
T1: water and nutrient deficiency stress, T2: water deficiency stress, T3: Nutrient deficiency stress; The percentage in parenthesis indicates the proportion of the gene in all up-regulated or down-regulated genes
表 1 qRT-PCR验证的引物信息
Table 1 Primers for qRT-PCR validation
基因
Gene引物序列(5′→3′)
Primer sequenceEucgr.F02337.1 F: 5′-CTGTTCGGTGAATTGCCACG-3′ R: 5′-AGGGATCTGCTGCGAGAATG-3′ Eucgr.F00180.1 F: 5′-CTCGTCGCTCGAATGGAAGA-3′ R: 5′-CTTCTTGGAGTGAGGCTCCG-3′ Eucgr.K02129.1 F: 5′-GGGAACGACTGACAGCTGAA-3′ R: 5′-CCGGATGTATGCTCGGAGAC-3′ Eucgr.B03214.1 F: 5′-CCAACACCACTTCGACCTCTC-3′ R: 5′-TGACCTTGAAGGAGAGGGACT-3′ 表 2 3种胁迫处理中尾叶桉叶片转录组测序质量1)
Table 2 Transcriptome sequencing quality of Eucalyptus urophylla leaves in three stress treatments
处理
Treatment重复
Replicate原始序列数
No. of raw reads去杂序列数
No. of clean reads去杂碱基数/(×109)
No. of clean bases错误率/%
Error rateQ20/% Q30/% GC含量/%
GC contentCK 1 102 875 922 98 982 508 14.85 0.01 97.84 94.29 47.41 2 91 836 470 88 460 334 13.07 0.02 96.61 91.37 47.30 3 88 160 956 85 873 298 12.82 0.02 95.22 88.34 45.19 T1 1 118 961 564 112 278 350 16.84 0.02 95.26 88.29 47.78 2 99 694 016 96 296 560 14.44 0.01 98.23 95.22 47.15 3 96 089 830 93 089 066 13.96 0.02 97.23 92.76 50.93 T2 1 93 283 460 90 948 242 13.64 0.01 98.39 95.68 45.45 2 91 547 658 88 602 608 13.27 0.02 95.06 87.93 45.32 3 92 936 506 90 084 620 13.51 0.01 98.15 95.05 46.51 T3 1 94 953 336 91 803 584 13.55 0.02 96.68 91.55 47.04 2 90 612 280 87 845 130 13.11 0.02 95.14 88.20 48.26 3 97 341 022 94 914 876 14.24 0.01 97.47 93.33 44.54 1) T1:水肥双因素胁迫,T2:水分单因素胁迫,T3:养分单因素胁迫;Q20:Phred数值大于20的碱基占总碱基的百分比,Q30:Phred数值大于30的碱基占总碱基的百分比;GC含量:去杂序列中G与C在4种碱基中占的百分比
1) T1: Water and nutrient deficiency stress, T2: Water deficiency stress, T3: Nutrient deficiency stress; Q20: Percentage of bases with phred greater than 20 in total bases; Q30: Percentage of bases with phred greater than 30 in total bases; GC content: Percentage of G and C in four bases in clean reads -
[1] ZHOU Q, GUO J J, HE C T, et al. Comparative transcriptome analysis between low- and high-cadmium-accumulating genotypes of pakchoi (Brassica chinensis L.) in response to cadmium stress[J]. Environ Sci Technol, 2016, 50(12): 6485-6494. doi: 10.1021/acs.est.5b06326
[2] JOHNSON S M, LIM F L, FINKLER A, et al. Transcriptomic analysis of Sorghum bicolor responding to combined heat and drought stress[J]. BMC Genomics, 2014: 15. doi: 10.1186/1471-2164-15-456.
[3] AHUJA I, DE VOS R C H, BONES A M, et al. Plant molecular stress responses face climate change[J]. Trends Plant Sci, 2010, 15(12): 664-674. doi: 10.1016/j.tplants.2010.08.002
[4] MITTLER R. Abiotic stress, the field environment and stress combination[J]. Trends Plant Sci, 2006, 11(1): 15-19. doi: 10.1016/j.tplants.2005.11.002
[5] RIZHSKY L, LIANG H J, SHUMAN J, et al. When defense pathways collide: The response of Arabidopsis to a combination of drought and heat stress[J]. Plant Physiol, 2004, 134(4): 1683-1696. doi: 10.1104/pp.103.033431
[6] RAMPINO P, MITA G, FASANO P, et al. Novel durum wheat genes up-regulated in response to a combination of heat and drought stress[J]. Plant Physiol Biochem, 2012, 56: 72-78. doi: 10.1016/j.plaphy.2012.04.006
[7] RIZHSKY L, LIANG H J, MITTLER R. The combined effect of drought stress and heat shock on gene expression in tobacco[J]. Plant Physiol, 2002, 130(3): 1143-1151.
[8] MESEKA S, MENKIR A, BOSSEY B, et al. Performance assessment of drought tolerant maize hybrids under combined drought and heat stress[J]. Agronomy-Basel, 2018, 8(12). doi: 10.3390/agronomy8120274
[9] MACHADO S, PAULSEN G M. Combined effects of drought and high temperature on water relations of wheat and sorghum[J]. Plant Soil, 2001, 233(2): 179-187
[10] 张胜. 侧柏对干旱与自然低温胁迫响应的分子机制研究[D]. 杨凌: 西北农林科技大学, 2017. [11] 潘文, 张卫华, 杨会肖, 等. 尾叶桉苯丙氨酸解氨酶基因的克隆、表达与单核苷酸多态性分析[J]. 华南农业大学学报, 2018, 39(2): 89-94. doi: 10.7671/j.issn.1001-411X.2018.02.014 [12] 赛佳. 根冠比[EB/OL]. [2019-12-20]. https://baike.sogou.com/v8081475.htm?fromTitle. [13] 曹维. 苗期玉米叶片响应非生物胁迫的转录组分析[D]. 扬州: 扬州大学, 2017. [14] 杨献光, 梁卫红, 齐志广, 等. 植物非生物胁迫应答的分子机制[J]. 麦类作物学报, 2006, 26(6): 158-161. doi: 10.3969/j.issn.1009-1041.2006.06.034 [15] 周峰. 干旱胁迫对植物生长及其生理的影响[J]. 北京农业, 2014. [16] 白禹, 卯丹, 王继玥. 黄酮类物质的生物合成及对干旱胁迫响应的研究进展[J]. 安徽农业科学, 2018, 46(16): 24-26. doi: 10.3969/j.issn.0517-6611.2018.16.008 [17] 王冰, 程宪国. 干旱、高盐及低温胁迫下植物生理及转录因子的应答调控[J]. 植物营养与肥料学报, 2017, 23(6): 1565-1574. doi: 10.11674/zwyf.17312 [18] ABE H, URAO T, ITO T, et al. Arabidopsis AtMYC2 (bHLH) and AtMYB2 (MYB) function as transcriptional activators in abscisic acid signaling[J]. Plant Cell, 2003, 15(1): 63-78. doi: 10.1105/tpc.006130
[19] SEO P J, XIANG F, QIAO M, et al. The MYB96 transcription factor mediates abscisic acid signaling during drought stress response in Arabidopsis[J]. Plant Physiol, 2009, 151(1): 275-289. doi: 10.1104/pp.109.144220
[20] 关淑艳, 焦鹏, 蒋振忠, 等. MYB转录因子在植物非生物胁迫中的研究进展[J]. 吉林农业大学学报, 2019, 41(3): 253-260. [21] 戚莹雪, 王蕾, 王尧尧, 等. 植物WRKY类转录因子在非生物胁迫下的功能与作用机制[J]. 分子植物育种, 2019, 17(18): 5973-5979. [22] MA T, LI M, ZHAO A, et al. LcWRKY5: An unknown function gene from sheepgrass improves drought tolerance in transgenic Arabidopsis[J]. Plant Cell Rep, 2014, 33(9): 1507-1518. doi: 10.1007/s00299-014-1634-3
[23] CHEN Y F, LI L Q, XU Q, et al. The WRKY6 transcription factor modulates PHOSPHATE1 expression in response to low Pi stress in Arabidopsis[J]. Plant Cell, 2009, 21(11): 3554-3566. doi: 10.1105/tpc.108.064980
[24] 于燕杰, 张大兵, 袁政. WOX蛋白家族调控干细胞发育分子机制的研究进展[J]. 植物学报, 2016, 51(4): 565-574. doi: 10.11983/CBB15137 [25] HIRAKAWA Y, KONDO Y, FUKUDA H. TDIF peptide signaling regulates vascular stem cell proliferation via the WOX4 homeobox gene in Arabidopsis[J]. Plant Cell, 2010, 22(8): 2618-2629. doi: 10.1105/tpc.110.076083
[26] JI J, STRABLE J, SHIMIZU R, et al. WOX4 promotes procambial development[J]. Plant Physiol, 2010, 152(3): 1346-1356. doi: 10.1104/pp.109.149641
-
期刊类型引用(2)
1. 李一华,徐贤柱. 丝状真菌培养及其产物的研究进展. 生物化工. 2022(02): 148-151 . 百度学术
2. 贺新生,张锐杰,王亚洲,贺亚宁,赵栎棋. 中国炭角菌属Xylaria物种多样性. 食用菌. 2022(04): 14-23 . 百度学术
其他类型引用(0)