APP下载

国外引进亚麻种质资源遗传多样性分析

2019-12-31郭栋良江海霞叶佳丽颜清城谢丽琼

新疆农业科学 2019年11期
关键词:蒴果油用亚麻

郭栋良, 江海霞, 张 喻,叶佳丽,王 玥,颜清城,谢丽琼

(新疆大学生命科学与技术学院,乌鲁木齐830046)

0 引 言

【研究意义】亚麻(LinumusitatissimumL.)是一种古老的植物纤维和油料作物,至今为止仍然是一个全球重要的经济作物[1]。韧皮纤维优质[2],亚麻籽粒含有丰富不饱和脂肪酸,是重要的医疗保健食品[3]。丰富遗传多样性的种质资源是作物育种的基础,国内亚麻育种长期使用骨干亲本使得遗传多样性较低,品种选育难以取得较大的突破[4]。从国外引进种质可以拓宽亚麻种质基因库,为育种提供丰富的遗传亲本。对引进亚麻种质资源的农艺性状进行详细的分析,对于高效利用亚麻种质资源是十分必要的。【前人研究进展】 目前,形态学标记和DNA分子标记已被广泛用于评估作物种质资源的遗传多样性[5-11]。对亚麻遗传多样性已有较多研究,邓欣等[12]用25个RAPD引物将10个亚麻种质分为3个类群,纤用亚麻和油用亚麻的亲缘关系较远; 王玉富等[13]将150份国外引进亚麻种质资源依据农艺性状分为4类,为引进亚麻资源的利用提供依据。李明等[14]用7对引物对85份亚麻品种进行AFLP分析,发现纤用亚麻的遗传背景狭窄,其遗传多样性远低于油用亚麻;王利明等[15]对国外引进的256份油用亚麻主要农艺性状进行主成分和系统聚类分析,筛选出3份农艺性状优良,利用价值较高的材料。邓欣等[4]对15个数量性状主成分分析,基于PC1和PC2 把535份亚麻分为纤用和油用两个基因库。2015年Habibollahi等[16]用SSR标记分析了23个亚麻种质的遗传多样性,发现材料间的遗传距离和地理距离之间有显著的相关性。李秋芝等[17]对国内外300份亚麻种质资源农艺性状进行鉴定,筛选出早熟、高产、抗倒伏的优良种质。【本研究切入点】尽管分子标记已广泛应用于亚麻遗传多样性研究,但采用的分子标记主要以DNA标记和数量性状为主,而遗传多样性是基因型多样性与环境的共同作用[18],表型性状能够快速对引进种质资源进行鉴定评价。研究对国外引进亚麻种质表型性状进行多样性分析。【拟解决的关键问题】选用144份来源于不同国家的亚麻种质,利用24个农艺性状对种质资源进行遗传多样性、相关性、主成分和聚类分析。为亚麻育种奠定基础。

1 材料与方法

1.1 材 料

供试144份亚麻种质来源产于26个国家,大部分为印度、美国、澳大利亚和巴基斯坦等国家,所有供试亚麻种质均从加拿大、美国亚麻种质资源库中引入。选用144份来源于不同国家的亚麻种质,利用24个农艺性状对种质资源进行遗传多样性、相关性、主成分和聚类分析。表1

表1 144份亚麻种质来源组成Table1 Origin of 144 flax accessions

来源Origin品种数Numbers来源Origin品种数Numbers来源Origin品种数Numbers印度 India16波兰 Poland4伊朗 Iran2美国 The United States14德国 Germany3意大利 Italy2澳大利亚 Australia11亚美尼亚 Armenia3捷克共和国 Czech Republic1巴基斯坦 Pakistan9阿富汗 Afghanistan2尼泊尔 Nepal1法国 France9埃及 Egypt2瑞典 Sweden1俄罗斯 Russia8爱尔兰 Ireland2危地马拉 Guatemala1加拿大 Canada7比利时 Belgium2希腊 Greece1埃塞俄比亚 Ethiopia6罗马尼亚 Romania2以色列 Israel1匈牙利 Hungary6土耳其 Turkey2未知 Unknown26

1.2 方 法

试验材料于2016年10月播种于云南大理州经济作物研究所麻类示范种植中心试验地(101°52'E,25°44'N;海拔1 118.4 m,)。小区面积为2 m2,区长1 m,宽2 m,5行区条播,有效播种粒数为2 000粒/m2。区间道0.4 m,组间道1 m,试验区四周设1.5 m宽保护行,按常规生产方式进行田间管理。农艺性状调查参照《亚麻种质资源描述规范和数据标准》[19]进行。选取30株统计调查株高、工艺长度、工艺长度/株高(株型)、茎粗、分枝数、分枝角度、蒴果数、种子长度、种子宽度、种子长度/种子宽度(种子形状)、千粒重、子叶色、子叶形状、叶色、萼片色、花冠形状、花瓣色、柱头色、花药色、花丝色、果实形状、分枝习性、抗倒伏性和形态一致性,共计24个农艺性状。

1.3 数据处理

运用Excle计算最大值、最小值、平均数、标准差、变异系数和遗传多样性指数,Shannon-weave多样性指数计算公式为:H=-ΣPi·InPi,其中Pi为某一性状第i级出现的频率。质量性状按表2分级直接计算各级频率,数量性状需要进行质量化处理,从第1级[Xi<(X-2s)]到第10级[Xi>(X+2s)],每0.5 s为1级[20]。运用Graph Pad Prism5软件对亚麻种质资源主要农艺性状进行相关性分析。主成分分析采用SPSS 22软件,聚类分析采用MeV4.9.0(Multi Experiment Viewer,http://www.tm4.org/)作图。

2 结果与分析

2.1 亚麻农艺性状描述统计

研究表明,13个质量性状的遗传多样性指数变幅在0.37~1.20,平均为0.83。其中抗倒伏性、花冠形状和柱头色的遗传多样性指数较大,这些性状在每一级中的分布比较均匀;抗倒伏材料占所有材料的11.81%。花药色、子叶形状和果实形状的遗传多样性指数较小,表明这几个性状分布比较集中。花药色为蓝色的材料占所有材料的91.67%。花色以蓝色为主,从浅蓝色到深蓝色连续过渡,占所有材料的93.74%,白花占5.56%。引进亚麻种质之间存在较大差异。表2

表2 亚麻种质资源质量性状分级标准及遗传多样性
Table 2 Grading standard and genetic diversity of qualitative traits in flax germplasm resources

性状Traits分级Grading scale遗传多样性指数H'频率分布 Distribution of frequency1234567子叶色CC1:浅绿,2:黄绿,3:绿,4:深绿色0.91 64.5812.5022.220.69---子叶形状CS11:椭圆,2:长椭圆0.58 73.6126.39-----叶色LC1:浅绿,2:绿,3:深绿0.97 21.5358.3320.14----萼片色SC11:浅绿,2:绿,3:深绿0.88 45.8347.926.25----花冠形状CS21:圆锥形,2:漏斗形,3:五角星形,4:碟形,5:轮形1.16 18.0655.56015.9710.42--花瓣色PC1:白,2:粉,3:红,4:黄,5,浅蓝,6,深蓝,7:紫0.82 5.5600066.6727.080.69柱头色SC21:白色,2:粉红色,3:蓝色,4:浅紫色,5:紫色1.03 38.89049.312.789.03--花药色AC1:微黄,2:橘黄,3:浅灰,4:蓝0.37 3.473.471.3991.67---花丝色FC1:白,2:蓝,3:紫0.64 66.6733.330----果实形状CS31:扁圆形,2:球形,3:卵形0.60 9.0381.949.03----分支习性BH1:紧凑型,2:中间型,3:松散型0.95 59.7219.4420.83----抗倒伏性LR1:0级,2:1级,3:2级,4:3级1.20 11.8143.0636.119.03---形态一致性NM1:一致,2:连续变异,3:非连续变异0.73 50.6948.610.69----

CC: Cotyledon , CS1: Cotyledon shape, LC: Leaf , SC1: Sepal , CS2: Corolla shape, PC: Petaline , SC2: Stigma , AC: Anther , FC: Filament , CS3: Capsule shape, BH: Branch habit, LR: Lodging resistance, UM: in morphologic, the same as below

研究表明,11个数量性状遗传多样性指数变幅在1.97~2.09,平均为2.03,其中工艺长度/株高的遗传多样性指数最高,蒴果数的遗传多样性指数最低。亚麻种质资源主要农艺性状变异系数呈现较大差异,纤维产量相关性状如株高和工艺长度变异系数分别为19.88%和27.58%,其中株高最高材料达128.13 cm, 株高最低材料仅为39.00 cm, 这些材料分别适用于纤用和油用亚麻育种。种子产量相关性状如株高、分枝数、蒴果数、种子尺寸和千粒重等性状中蒴果数的变异系数最大为40.74%,变幅为12.80~106.67个;种子长度和种子宽度的变异系数较小,分别为10.64%和8.78%,变幅为2.80~5.19 mm和1.73~2.72 mm;种子宽度/种子长度的变异系数最小为6.22%,变幅为0.46~0.62,说明种子形状变异较小。千粒重变异系数为27.72%,变幅为1.72~10.51 g;其中多分枝、大粒种子等性状适用于国内油用亚麻改良。表3

2.2 相关性

研究表明,株高与工艺长度、工艺长度/株高、茎粗、分枝角度、种子宽度/长度呈显著正相关,与分枝数、种子长度、种子宽度呈显著负相关。蒴果数与茎粗、分枝数、分枝角度、种子长度、种子宽度、千粒重均呈现显著正相关,与工艺长度、工艺长度/株高、种子宽度/长度呈现极显著负相关,而与株高未呈现显著相关性。千粒重与分枝数、蒴果数、种子长度、种子宽度呈现显著正相关,与株高、工艺长度、工艺长度/株高、种子宽度/长度呈现显著负相关。工艺长度/株高与株高、工艺长度、分枝角度、种子宽度/长度呈现显著正相关,与分枝数、蒴果数、种子长度、千粒重呈现显著负相关。种子宽度/长度与株高、工艺长度、工艺长度/株高呈现显著正相关,与分枝数、蒴果数、种子长度、千粒重呈现显著负相关。表4

表3 亚麻种质资源数量性状变异参数
Table 3 The variation parameter of quantitative traits in flax germplasm resources

性状Traits最大值Max.最小值Min.极差Range平均数Mean标准差SD变异系数(%)CV遗传多样性指数H'株高PH128.1339.0089.1378.6115.6319.882.03工艺长度TL105.5022.2283.2851.5814.2327.582.01工艺长度/株高TL/PH0.840.440.400.650.0811.922.09茎粗SD3.941.422.522.500.4216.832.02分枝数BN10.004.205.806.981.1115.962.06分枝角度BA39.2515.9023.3523.383.7015.832.02蒴果数CN106.6712.8093.8742.3117.2340.741.97种子长度SL5.192.802.394.150.4410.642.07种子宽度SW2.721.730.992.180.198.782.05种子宽度/长度SW/SL0.620.460.160.530.036.222.03千粒重1 000-GW10.511.728.795.701.5827.722.03平均Mean-----18.372.03

PH:Plant height,TL:Technical length,TL/PH:Technical length/plant height, SD:Stem diameter,BN: number,BA:Branch angle,CN:Capsule number,SL:Seed length,SW:Seed width,SW/SL:Seed width/Seed length, 1,000-GW:1,000-grain weight, the same as below

表4 亚麻种质资源主要农艺性状相关性
Table 4 Correlation of main agronomic traits in flax germplasm resources

TraitPHTLTL/PHSDBNBACNSLSWSW/SLTL0.915 8∗∗∗ TL/PH0.396 2∗∗∗0.722 4∗∗∗SD0.549 3∗∗∗0.392 1∗∗∗-0.015 3BN-0.210 6∗-0.243 9∗∗-0.168 4∗0.494 8∗∗∗BA0.197 3∗0.232 3∗∗0.189 8∗0.399 2∗∗∗0.264 9∗∗CN-0.092 1-0.313 8∗∗∗-0.538 9∗∗∗0.613 3∗∗∗0.654 8∗∗∗0.265 1∗∗SL-0.407 7∗∗∗-0.397 1∗∗∗-0.212 1∗0.104 10.418 7∗∗∗0.042 20.331 0∗∗∗SW-0.269 3∗∗-0.244 6∗∗-0.086 10.075 90.347 7∗∗∗-0.004 90.227 1∗∗0.827 7∗∗∗SW/SL0.336 2∗∗∗0.350 4∗∗∗0.246 7∗∗-0.078 1-0.247 5∗∗-0.082 4-0.255 6∗∗-0.573 2∗∗∗-0.019 81 000-GW-0.397 1∗∗∗-0.395 5∗∗∗-0.210 1∗0.102 30.431 9∗∗∗-0.037 50.345 0∗∗∗0.911 9∗∗∗0.868 8∗∗∗-0.361 0∗∗∗

注:*、**和***分别表示显著水平为 0.05、 0.01 和 0.001

Note:“*”, “**” and “***” indicate significance atP<0.05,P<0.01 andP<0.001, respectively

2.3 主成分

研究表明,在PCA分析中,前9个成分解释了73.57%的表型变异,其中,前2个主成分解释了大约32.31%的表型变异。PC1代表“油用亚麻特征性状因子”如株高(r=-0.594 1),工艺长度(r=-0.662 8),分枝数(r=0.631 3),蒴果数(r=0.605 8),种子长度(r=0.859 3),种子宽度(r=0.686 1),种子宽度/长度(r=-0.531 6),千粒重(r=0.842 7),占总方差变异的20.25%。PC2代表“纤用亚麻特征性状因子”如株高(r=0.658 1),工艺长度(r=0.658 9),茎粗(r=0.797 6)和分枝角度(r=0.601 0),占总方差变异的12.06%。PC3-9主要解释花瓣色,花药色,抗倒伏性,形态一致性,果实形状等,分别占9.17%、7.26%、6.08%、5.39%、4.75%、4.55%和4.06%的总方差变异。表5

表5 亚麻种质资源主要农艺性状的主成分
Table 5 Principal components analysis of main agronomic traits in flax germplasm resources

性状Traits主成分Principal component123456789株高PH-0.594 10.658 1-0.073 1-0.019 30.157 3-0.019 00.138 20.106 6-0.128 5工艺长度TL-0.662 80.658 90.137 6-0.103 80.094 20.111 00.032 20.105 1-0.000 8工艺长度/株高TL/PH-0.486 00.409 30.429 5-0.204 4-0.071 40.297 5-0.132 70.068 80.239 7茎粗SD0.137 10.797 6-0.303 90.197 50.280 9-0.163 90.090 9-0.023 7-0.122 5分枝数BN0.631 30.384 7-0.215 00.157 60.047 0-0.185 20.070 6-0.074 20.220 7分枝角度BA0.076 30.601 0-0.173 20.017 8-0.047 6-0.009 9-0.344 8-0.059 30.121 8蒴果数CN0.605 80.299 7-0.482 30.304 70.197 9-0.284 50.043 7-0.003 2-0.021 5种子长度SL0.859 30.087 90.345 1-0.087 50.139 70.148 90.067 7-0.034 4-0.104 8种子宽度SW0.686 10.147 30.480 9-0.079 4-0.015 7-0.020 20.371 00.183 2-0.012 9种子宽度/长度SW/SL-0.531 60.048 10.073 80.040 6-0.259 3-0.284 50.414 90.345 20.164 5千粒重1 000-GW0.842 70.090 20.392 4-0.024 70.051 60.067 90.212 1-0.019 5-0.041 9子叶色CC0.070 70.040 3-0.001 4-0.396 20.435 20.097 1-0.014 40.158 6-0.419 0子叶形状CS10.506 70.100 70.390 6-0.069 90.008 80.086 5-0.285 00.206 60.057 2叶色LC-0.102 20.385 30.418 4-0.200 3-0.081 0-0.156 80.294 3-0.366 6-0.002 9萼片色SC1-0.105 8-0.050 30.321 5-0.040 70.422 0-0.139 90.014 5-0.485 90.457 7花冠形状CS20.421 10.370 50.028 4-0.016 0-0.403 90.081 8-0.325 60.045 6-0.222 3花瓣色PC-0.130 0-0.140 50.465 30.639 50.210 4-0.018 2-0.284 0-0.075 3-0.147 1柱头色SC2-0.223 60.010 60.337 60.441 6-0.030 6-0.036 90.245 10.256 2-0.163 8花药色AC-0.219 6-0.000 60.320 20.737 5-0.032 2-0.110 2-0.072 5-0.019 1-0.048 0花丝色ASC-0.354 3-0.246 10.108 2-0.087 00.635 8-0.048 1-0.071 00.140 2-0.139 4果实形状CS30.220 00.086 50.106 8-0.047 70.291 0-0.137 7-0.261 30.550 40.446 7分支习性BH0.092 5-0.312 3-0.350 50.106 70.336 20.289 50.275 10.053 80.101 4抗倒伏性LR0.244 4-0.008 1-0.260 40.301 4-0.085 20.549 00.186 30.171 70.250 5形态一致性UM-0.150 60.313 3-0.037 40.231 10.121 00.674 40.077 1-0.184 9-0.002 1特征值Eigen value4.861 12.893 42.199 81.743 41.458 41.294 41.141 11.092 20.973 5贡献率(%)Contributive rate20.254 712.056 09.165 67.264 46.076 55.393 44.754 44.550 74.056 2累计贡献率(%)Cumulative contributive rate20.254 732.310 641.476 348.740 654.817 160.210 564.964 969.515 673.571 8

研究表明,144份亚麻种质主要划分为两个群体。PC1从负到正,株高、工艺长度逐渐减小,分枝数、蒴果数逐渐增加。现代栽培亚麻主要分为油用亚麻和纤用亚麻,前者表现为株高较矮而分枝蒴果数多等特征,后者株高工艺长度长而少分枝、蒴果等特征。分析各类群株高、分枝数等重要农艺性状差异,发现Group I、Group II两大类群间各农艺性状均呈现极显著差异。与Group II相比,GroupI表现为株高、工艺长度高,分枝数、蒴果数低等特征。Group II表现为株高、工艺长度低,分枝数、蒴果数高特征。Group I明显表现为纤用亚麻特征,Group II表现为油用亚麻特征。图1,图2

图1 基于PC1和PC2的二维主成分
Fig. 1 Two-dimensional principal component analysis based on PC1 and PC2

A:种子宽度/长度性状极端材料表型分布; B:千粒重性状极端材料表型分布; C:株高与蒴果数曲线回归分析; D-G:亚麻不同类群农艺性状差异比较。D:株高; E:工艺长度; F:分枝数;G:蒴果数

A: The phenotypic distribution of extreme material of seed width / length traits; B:The phenotypic distribution of extreme material of 1000 grain weight; C: Regression analysis on the curve of plant height and capsule number; D-G: Comparison of agronomic characters of different flax groups. D: Plant height; E: Technical length; F: Branch number,The number of branches; G: Capsule number

图2 亚麻种质资源各类群农艺性状表现特征
Fig. 2 Characteristics of agronomic traits for every group of flax germplasm resources

2.4 聚类分析

基于24个农艺性状标准化数据将供试亚麻材料聚为 2 大类群 ,形态学标记先将纤用亚麻和油用亚麻区分开。Group I 表现为株高、工艺长度较高,而分枝数、蒴果数、千粒重较低等纤用亚麻特征,Group II表现为株高、工艺长度较矮等油用亚麻特征,而兼用亚麻主要集中在油用亚麻中。Group II类群又分为3个亚类群,Group II-1包含36份材料,表现为分枝数和蒴果数中等而种子长度、种子宽度和千粒重最大等特征。Group II-2包含40份材料,主要表现为株高和工艺长度最矮且分枝数和蒴果数最低等特征。 Group II-3包含33份材料,主要表现为茎粗、分枝数、蒴果数较高等特征。兼用材料主要集中在Group II-3中,表现出株高、工艺长度、分枝数、蒴果数均较高的特征。

Group I中材料可作为纤用亚麻种植或作为培育高产纤用亚麻的亲本。Group II-1中大粒种质可用来改良油用亚麻籽粒产量。多数矮杆材料集中在Group II-2中,例如80、41和14号,株高仅为40 cm左右,这些材料可用来改良油用材料的株高。Group II-3中部分兼用亚麻材料株高、工艺长度、分枝数、蒴果数均表现优良,可作为优质兼用亚麻进行利用。图3

图3 基于24个农艺性状的聚类分析
Fig. 3 Cluster analysis based on 24 agronomic traits

3 讨 论

具有丰富遗传多样性的种质资源是作物遗传改良的基础[21],研究引进的144份亚麻种质丰富了我国亚麻种质基因库,为亚麻育种提供重要的基因资源。纤用亚麻更关注株高因为株高直接与纤维产量相关[22],增加株高往往提高倒伏的威胁[23]。研究抗倒伏材料占所有材料的11.81%,油用亚麻产量形成往往受多性状调节,包括“高度因子”、“分枝因子”和“种子大小因子”等。研究引进种质包括多个矮杆材料,株高仅为40 cm左右,可用于油用亚麻矮杆育种。国内栽培油用亚麻种子千粒重约6 g左右,研究引进种质千粒重最大达10.51 g,可用于改良油用亚麻种子大小。然而,种子大小和种子数量之间存在制衡[24],使得协调“分枝因子”和“种子大小因子”进一步提高油用亚麻种子产量变得复杂。

纤用亚麻和油用亚麻在农艺性状上存在显著差异[25],与油用亚麻相比,纤维亚麻一般表现为株高较高,分枝数和蒴果数较少,种子产量较低等特征[26]。各农艺性状相关性分析中,株高、工艺长度等 “高度因子”与分枝数、蒴果数、种子长度、种子宽度、千粒重等“分枝因子”、“种子大小因子”呈现显著负相关,其中只有株高与蒴果数未呈现显著相关性,这与王利民[27]和张丽丽等[28]的研究结果一致,其主要原因是兼用亚麻株高、工艺长度、分枝数、蒴果数均较高而对株高和蒴果数相关性造成干扰,致使株高与蒴果数未能表现出显著关联。研究第一次采用工艺长度/株高来描述分枝位高低,采用种子宽度/长度描述种子形状。工艺长度/株高与株高、工艺长度等呈现极显著正相关,且与工艺长度相比,工艺长度/株高与蒴果数负相关程度最高,工艺长度/株高可以作为衡量亚麻用途的重要性状。种子宽度/长度与株高等呈现极显著正相关,与千粒重、种子长度等呈现极显著负相关,大粒种子更偏向于长粒,小粒种子更偏向于圆粒。

主成分和聚类分析常用于作物农艺性状评价[9, 29, 30],研究PCA分析中,前两个主成分解释大约32.31%的表型变异。PC1代表“油用亚麻特征性状因子”,PC2代表“纤用亚麻特征性状因子”,以PC1和PC2绘制二维坐标图,144份亚麻种质主要划分为两个群体。Group I表现为株高、工艺长度高,分枝数、蒴果数低等纤用亚麻特征。Group II表现为株高、工艺长度低,分枝数、蒴果数高等油用亚麻特征。利用24个农艺性状将144份亚麻材料聚为 2 大类群,形态学标记最先将纤用亚麻和油用亚麻区分开,表明亚麻驯化过程中产量相关性状受到主要选择。Group II-3中包含部分株高、工艺长度、分枝数、蒴果数均表现优良的兼用亚麻,这些材料或许受到纤用或油用亚麻基因流的影响。基于农艺性状的主成分和聚类分析均将亚麻种质聚为纤用亚麻和油用亚麻两大类群,这与前人的研究结果一致,认为亚麻种质根据品种用途分为不同的群体[14, 31]。

育种者的选择对亚麻遗传多样性产生重要影响。亚麻驯化过程中育种者紧紧围绕亚麻产量相关性状进行选择,纤用亚麻和油用亚麻产量相关性状间相互制衡,双向驯化导致纤用亚麻和油用亚麻在农艺性状上存在显著差异,主要农艺性状选择的差异造就纤用和油用亚麻基因组中的差异。

4 结 论

4.1 国外引进亚麻种质资源存在着较为丰富的多样性,质量性状中抗倒伏性的多样性指数最高为1.20,数量性状中蒴果数、千粒重、工艺长度的变异系数均较高,分别为40.74%,27.72% 和27.58%。这些产量相关性状在亚麻育种中有较大的选择空间。

4.2 主成分分析和聚类分析将144份亚麻材料聚为纤用和油用亚麻两个群体,表明驯化过程中纤用和油用亚麻差向选择形成显著的表型差异。

4.3 各农艺性状相关性分析中,"高度因子"与"分枝因子"、"种子大小因子"呈现显著负相关,表明"油用亚麻特征性状"和"纤用亚麻特征性状"相互制衡,体现种子产量和纤维产量权衡中的复杂性。

猜你喜欢

蒴果油用亚麻
光伏支架对油用牡丹生长特性影响的研究
蓝亚麻蒴果离散元模型仿真参数标定
油用牡丹助农增收
油用牡丹助农增收
与碑石对饮(外一首)
与碑石对饮(外一首)
亚麻·亚麻
吃亚麻籽有方法
新申亚麻&孙林:亚麻时尚的民俗融合
“逐梦森林”:新申亚麻&孙林时装发布会