临汾盆地黄土粒度分布的端元模型研究

2021-04-16赵格格田庆春杜五喜裴瑜鄂崇毅

海洋地质与第四纪地质 2021年2期

赵格格，田庆春，杜五喜，裴瑜，鄂崇毅

1.山西师范大学地理科学学院，临汾 041000 2.青海师范大学青藏高原资源与环境教育部重点实验室，西宁 810008 3.中国科学院青海盐湖研究所，西宁 810008

粒度作为反映古气候、搬运方式和物质来源的重要指标，被广泛地运用于黄土[1-3]、海洋[4]和湖泊[5]的古环境变化研究中。刘东生等[1]对中国黄土高原的典型马兰黄土样品的粒度分析表明，自西北向东南方向依次为砂黄土带、黄土带、黏化黄土带，粒度逐渐变细，黄土的粒度变化一般被解释为受控于风力强度的改变，被用作东亚冬季风变化的替代性指标。黄土-古土壤的沉积序列是在不同物源、多种搬运动力和复杂沉积环境下共同作用的结果，但是不同的粒级对环境变化的敏感程度是不同的，因此环境敏感粒级的提取对充分获得复杂的古环境变化信息来说是十分重要的。目前，提取敏感的环境粒度方法有端元模型[6-11]、粒级-标准偏差和主成分分析等，端元模型对黄土高原粒度数据进行反演应用已经非常广泛，且取得了大量成果[12-13]。

临汾盆地位于黄土高原东南缘，地处中国温带季风气候区，季风气候是这一地区最活跃的环境变迁因子，为获取地区高分辨率气候特征，一些研究者从粒度、磁化率、全铁以及有机质等指标对临汾盆地黄土-古土壤做了一定工作[14-15]，但关于沉积物环境敏感粒度组分的研究相对较少，缺乏使用其他数学方法来捕捉更高分辨率粒度组分的气候信息。不同区域颗粒组分指示的环境存在一定差异[16]，没有科学的分析很难对沉积环境有详细了解，而端元分析模型是识别多组分沉积物来源和查明搬运方式的有效工具，在古环境重建中得到广泛的应用。对此，本文选择临汾盆地丁村剖面黄土-古土壤序列为研究对象进行工作，通过对比参数化与非参数化端元分解方法的差异[17]，讨论两种方法对临汾盆地粒度数据的适用性后，选择参数化端元分析法提取环境敏感粒组的基础上，结合年代结果，揭示各端元的变化规律以及它们指示的古气候意义，并与深海氧同位素[18]为载体的气候代用指标进行对比分析，来揭示黄土高原东南缘临汾盆地的粒度分布特征及其环境意义。

1 区域概况

临汾盆地位于山西省南部，北依韩侯岭，南为峨眉台地，东至霍山，西抵罗云山，属于温带大陆性季风气候，年均气温11.7～13.2 ℃，年降水量453～547 mm，自然植被区带属暖温带落叶阔叶林带，山地与平原相接处，不同时代的黄土分布广泛，午城黄土和离石黄土主要分布于盆地边缘丘陵地带的冲沟当中，出露面积很小，马兰黄土分布广泛，覆盖在山麓丘陵、平原及高阶地之上，厚5～20 m。马兰黄土之上覆盖全新世黄土-古土壤，分布在盆地的南部和北部，蕴含着丰富的古气候演化信息。本研究选择襄汾县城南部汾河右岸丁村遗址附近的一黄土剖面（丁村剖面）为研究对象（如图1），经纬坐标为 35°49′34″N、111°25′32″E，海拔高度 406.5 m。

2 研究方法

2.1 样品采集

图1 丁村剖面位置图Fig.1 The site of the Dingcun profile

通过野外考察，我们选择丁村古人类遗址附近黄土-古土壤剖面为研究对象，命名为丁村剖面。剖面位于汾河三级阶地[19]，野外可清晰分辨出S1古土壤层，从剖面顶部以2 cm间隔进行样品的采集，一直采集到S1底部，共采集了310个样品，所采剖面总厚度约620 cm。同时使用长40 cm，直径4.5 cm的不锈钢管进行光释光年代样品的采集，以50～100 cm不等距间隔采集8组光释光（OSL）样品。

2.2 分析测试方法

2.2.1 粒度测试方法

粒度参数首先对样品进行前处理，前处理方法见参考文献[20]，然后使用Mastersizer 2000型激光粒度仪进行测试，测试范围0.02～2 000 μm，重复测量误差＜1%。

2.2.2 OSL测试方法

OSL样品在青海师范大学光释光年代学室测定，使用Risø TL/OSL-DA-20-C/D型热/光释光仪进行测试（具体步骤见参考文献[21]），选择180～250 μm粒级组分的粗颗粒钾长石为测年物质，采用富钾长石pIRIR法[22-23]，对钾长石进行多次IR激发的同时增加样品的测量温度可以将没有异常衰减的IRSL信号分离出来，每一个测量温度下产生的IRSL信号可以获得一个Do值和年龄。

2.2.3 端元分析法

端元分析可以把粒度数据中的一些混合信息分解成很多端元，这些端元是由特定的动力机制分选出来的，通过对粒度数据进行不同的处理方式而分为两种方法：（1）参数组分分离是通过对单个样品进行拟合，设定最佳参数，以达到分离出样品亚组分的目的。Paterson等[17]研究出了参数化的端元分析法，认为通过特定的动力方式所沉积的黄土具有单峰形态，与曲线拟合方法的结果较为一致。在程序中设置了多种拟合方法，输入全样本数据并设置拟合函数值后，即可得到最终的端元拟合结果，即为最优结果。（2）非参数组分分离是通过对整体样品进行拟合，目的是分离出几种具有相同频率分布的粒度组分，Weltje[11]等开发了一种综合考虑多个数据并根据可能的动力来源进行分解的非参数方法，该方法应用比较广泛，但分解结果可能具有多个峰值形状，与曲线拟合方法的结果有些不同，对动力来源的解释存在困难。对丁村剖面黄土的粒度数据，应用沉积物粒度端元模型进行分析，Weltje等[7,10-11]认为，沉积物粒度组分是在不同的搬运或改造动力方式共同作用下由不同物质来源的颗粒物沉积形成，因此提出了沉积物端元组分的概念，是根据一组样品的沉积学粒度特征，从统计学意义上对各不同端元组分进行鉴别，Paterson等[17]在此基础上基于MATLAB软件进行了改良，使用参数组分分离方法对单个样品进行拟合。

3 年代结果

年代学结果见表1，由于OSL只能得到地层的年代框架，因此以5个准确的OSL数据为控制点，选择前人提出的粒度-年龄模型[24]，带入中值粒径得到了丁村剖面年代框架（图2），570 cm处（S1底界）年代约为128.5 ka，与前人建立的黄土地层年代标尺S1/L2界限年龄为128.8 kaBP较为一致[25]，说明我们的年代-深度模型是正确的。剖面顶部年代约20 ka（0～20 ka的地层可能受到侵蚀），说明本剖面包含与深海氧同位素一致的MIS5、MIS4和MIS3这3个阶段。通过野外观察与测年数据，将剖面自上而下分为4层：0～60 cm，L1L1黄土层；60～190 cm，L1S1弱古土壤层；190～420 cm，L1L2黄土层；420～570 cm，S1古土壤层，为浅红棕色黏土。

4 结果分析

4.1 临汾黄土粒度分类结果

谢帕德三角分类法[26-27]被广泛地接受，并应用到黄土-古土壤的粒度分类研究当中[28]。沉积物分类图可以显示剖面样品的粒度组成情况，可以非常明显地看出剖面沉积物粒径属于哪一类别。故使用谢帕德三角分类法对临汾盆地黄土-古土壤进行分类，结果如图3。丁村剖面黄土-古土壤粒径主要分布在粉砂和黏土质粉砂中，分布较为集中。

4.2 端元模型结果分析

4.2.1 端元数量的确定

本文在MATLAB环境下，运用Paterson等人[17]编写的Analysize-masters程序，导入丁村剖面的全样粒度数据，在假设端元数为1～10的基础上，进行参数化和非参数化分析，然后通过对线性相关和角度离差进行综合分析。参数化方法是在Weibull分布函数中增加了一个附加位置参数的Gen.Weibull分布函数来进行分解，是为了在形状上有很大的灵活性并且能很好地控制偏度。关于选择最优端元数量，主要考虑3个指标[29]：线性相关度，是为了说明原数据集与端元的总体相关度，相关度越高拟合程度越好，当线性相关度大于0.9时，满足拟合要求，线性拟合也很良好；角度离差，是原样粒度曲线与端元在进行形状拟合时产生的偏差，偏差值越小则端元曲线在形状拟合的误差越小，当偏差值小于5°时拟合结果越理想；端元相关度，是分解得到的每个端元之间的相关度，相关度越高各个端元之间的分布有较高的重合，说明其独立性差、拟合过度。

表1 丁村剖面OSL年代结果Table 1 OSL dating results of Dingcun profile

图2 丁村剖面年龄模型Fig.2 Age model of Dingcun profile

图3 谢帕德三角分类Fig.3 Shephard’s classification

然后通过综合分析两种方法的线性相关度、角度离差和端元相关度确定了最优的端元数量，其拟合特性见表2、图4。非参数方法的线性相关度小于参数化方法；非参数化方法的端元标准偏差均值却远大于参数化方法；非参数化方法的端元相关度大于参数化方法。因此，我们最终确定使用参数化分解方法，选取3个端元对丁村剖面粒度数据进行端元反演。

4.2.2 端元组分分析

通过端元模型对临汾盆地丁村剖面黄土粒度进行分析，反演得出了3个端元粒度频率分布曲线图（图5），将3个端元粒度组分分别命名为EM1、EM2和EM3，分布都接近正态分布，分选由差变好。EM1在2.51 μm的粒径时含量占比最多，曲线分布宽平，范围更广，峰值较低；EM2在12.61 μm的粒径时含量占比最多，曲线分布范围较窄，较为陡峭；EM3在35.56 μm的粒径时含量占比最多，曲线分布范围最窄，更集中在粗粒端，曲线窄而尖，峰值更高且陡峭，变化幅度大。从3个端元含量变化图中可以看出（图5），EM1的含量变化在3个端元中相对较为平缓，冷暖波动不明显，整体曲线没有明显的变化趋势；EM2含量变化也有比较大的波动，在0～60 cm出现低值，60～190 cm为高值，190～420 cm为低值，420～570 cm为高值，EM3含量变化具有较大的波动。

5 讨论

5.1 端元的环境意义

EM1的众数粒径为2.51 μm，与中值粒径的相关系数仅为0.090，几乎没有相关性（图6），一些学者研究认为黄土沉积后又经过了风化成壤作用，使用端元分析，从沉积物中分解出在风化与成壤作用下形成的平缓端元，因此众数粒径为2.51 μm的EM1是通过高空长距离搬运而来的远源细粒物，且这些细粒物是一直存在于大气中的本底粉尘，同时黄土沉积后，因为一系列的成壤风化过程也会对几微米以下的粒径产生影响[30-31]；与中国黄土众数粒径2～6 μm的细颗粒组分具有一致性[32]，故推测EM1可能是黄土在夏季风作用下的风化和成壤过程中形成的。

EM2的众数粒径为12.61 μm，与中值粒径的相关性较强（R2=0.492）（图6），根据大气动力学的研究结果，粒径范围2～16μm的细粉砂在风力的作用下可以被搬运到距离地面几千米的高空，它们长期漂浮在高空大气中，在风力减小时降落到全球任一地方沉积[33]；而且Vriend等[34]对黄土高原地区多个剖面进行研究，对黄土粒度运用端元分析后，也提出了细颗粒端元是被高空远距离搬运的，称之为远源端元；王兆夺等[13]对长城塬剖面的研究提出了细粒度组分，认为物源较远，是在搬运动力减弱的情况下沉积，因此它可能代表了西风带控制的由高空气流搬运的远源粉尘。基于上述研究，推测EM2可能代表的是高空西风搬运的远源沉积物。

表2 丁村剖面黄土粒度参数化与非参数化端元拟合特性Table 2 Fittings of parametric and nonparametric EMA for the Dingcun profile

图4 丁村剖面黄土粒度非参数化与参数化端元分解结果与相关判别指标Fig.4 Frequency distribution curves of each end-member in the Dingcun profile

图5 EM1、EM2、EM3频率分布与含量变化曲线Fig.5 Frequency distribution and content variations curves of EM1、EM2、EM3

图6 EM1、EM2、EM3与中值粒径的相关性分析Fig.6 Correlations of EM1、EM2、EM3 with median grain size

EM3的众数粒径为35.56 μm，与中值粒径之间表现出非常明显的相关性且呈正相关（图6），相关系数很高（R2=0.866），且EM3的粒径范围为3.34～60.45 μm，而中值粒径的范围为 9.44～21.85 μm，EM3的粒径变化幅度大于中值粒径。Pye[35]总结出20～70 μm的中粗粉砂和细砂在尘暴中做短距离的悬移运动，一般在大气中的运动范围只有几百米的距离，待风力减小时降落到地面上，大于30 μm的颗粒物质在较强风暴条件下最多迁移至距离源区30 km的地方，极端风暴条件下也只是迁移至距离源区几百千米的地方，而且搬运的方式主要是短时间的近地面悬移；孙东怀在对黄土高原粒度组分的研究中[31]，通过计算发现粗粒组分的平均粒径在21～54 μm范围内，认为该范围内的黄土由冬季风沉积，且搬运距离近，为低空搬运；在Prins等[36]对黄土高原多个剖面的研究中，把粗粒组分划分为两个端元，一个为40 μm左右的组分，该组分的黄土可能是短距离的悬移搬运，认为可以反映冬季风的变化情况，另一个为大于60 μm的组分，该组分可能为短距离跃移，也可以反映冬季风的变化情况。因此基于上述原因，推测EM3可以反映冬季风的变化情况，可能代表的是较强的冬季风所搬运的近源沉积物以及冬季短时间尘暴所带来的物质。

5.2 古环境和古气候变化

根据Gen.Weibull函数分布的参数化端元分析法，并将临汾盆地各环境指标与深海氧同位素记录[18]对比（如图7），对末次间冰期以来气候演化特征进行分析：

128.8 ka之前（570 cm以下），磁化率和频率磁化率在该阶段为低值，EM1和EM2细颗粒组分的含量减少，中值粒径和EM3粗颗粒组分的含量增加，说明该区冬季风增强，沉积物主要由冬季风带来，对应于MIS6晚期。

128.8～75 ka（570～420 cm），对应于 MIS5，磁化率和频率磁化率在该阶段急剧增高，在整个曲线出现最高值，EM1和EM2细颗粒组分的含量增加，中值粒径减小和EM3粗颗粒组分的含量减少，表明该阶段东亚冬季风较弱而夏季风较强，磁化率与EM3组分的变化与深海氧同位素的变化较为相似，表明临汾盆地与全球气候表现一致，处于温暖湿润时期。深海氧同位素在该阶段分为a、b、c、d、e五个明显的冷暖变化亚阶段，呈现暖-冷-暖-冷-暖的气候变化趋势，本区气候在该阶段虽然也有冷暖变化，但是没有表现出特别明显的冷暖波动，可能与临汾盆地距离冬季风的源区较远有关，沉积速率相对较低，加上后期强烈的成壤作用[37]和淋溶作用[38]，S1古土壤大多由1层或3层古土壤叠置而成[39]，可能影响其气候记录。

75～57 ka（420～190cm)，对应于 MIS4，磁化率和频率磁化率在该阶段急剧降低，在整个曲线出现最低值，EM1和EM2细颗粒组分的含量减少，中值粒径和EM3粗颗粒组分的含量增加，表明该阶段冬季风强劲，气候转冷，深海氧同位素在该阶段指示全球冰量增加，全球处于寒冷的时期，但内部有短暂回暖，可能受北大西洋千年尺度气候影响[40]，但具体的对应关系还有待进一步研究。

57～25 ka(190～60 cm)，对应于 MIS3，磁化率和频率磁化率在该阶段缓慢增加，EM1和EM2细颗粒组分的含量增加，中值粒径和EM3粗颗粒组分含量减少，表明该阶段气候较温暖，且与三阶段全球弱暖期一致，出现两次暖期，一次冷期，呈现暖-冷-暖气候变化趋势，但是该区域的冷暖变化趋势比深海氧同位素大，可能与东亚冬季风强度快速变化有关[41]。

图7 丁村剖面磁化率、频率磁化率、中值粒径、EM1、EM2和EM3与深海氧同位素记录[18]的对比Fig.7 Comparison of the magnetic susceptibility、frequency magnetic susceptibility、median grain size、EM1、EM2 and EM3 of the Dingcun profile with SPECMAP δ18O record[18]

25～20 ka(60 cm以上)，EM1和EM2细颗粒组分的含量减少，中值粒径和EM3粗颗粒组分的含量增加，磁化率和频率磁化率值降低，与深海氧同位素指示的全球冰量增加相一致，说明气候开始变冷，对应于MIS2早期。

临汾盆地黄土沉积物各指标波动特征与深海氧同位素一致，说明临汾盆地气候演化受全球冰量控制，但在每个阶段的内部波动上又表现出一定的区域特征，对我们研究黄土高原东南缘临汾盆地的粒度分布特征及其环境具有一定意义。