一种持续性注意水平的敏感性特征指标的选取方法与流程

文档序号:12777614阅读:530来源:国知局
一种持续性注意水平的敏感性特征指标的选取方法与流程

本发明涉及应激检测技术,特别涉及一种持续性注意水平的敏感性特征指标的选取方法。



背景技术:

与常规列车操作方式相比,在动车中的采用了列车自动控制系统(atc),使得动车司机的作业方式由传统的列车操纵型转变为对动车实时运行控制信息、车辆状态信息、运行环境信息的高负荷、长时间地持续监控。因此,动车司机的持续性注意水平成为影响司机作业可靠性的关键因素。选取高敏感于持续性注意水平的脑电特征参数,为动车司机的持续性注意水平进行测评及识别提供量化指标,对于保障动车运行安全是十分必要与迫切的。

目前,国内外就铁路司机的研究主要针对于疲劳,而对持续性注意水平的研究则较为鲜见。在现有技术中,研究人员已经采用人眼持续闭合时间比率(perclos)和人眼平均闭合速度(aecs)作为测评指标对司机疲劳进行研究,研究结果表明:持续闭合时间比率越大,人眼平均闭合速度越低,则司机的疲劳程度越大,并在此基础上开发了司机疲劳监测系统。研究人员还基于类似的3项眼动指标,提出一种基于贝叶斯网络信息融合的机车司机驾驶疲劳监测系统。研究人员还用机车司机驾驶状态监测视频文件,通过对机车司机的姿态(规定手势的完成状态)来判定司机的疲劳水平。研究人员还利用脑电信号(eeg)和心电信号(eog)研究了司机睡眠剥夺与持续性注意间的关系,结果表明:睡眠时间的减少,疲劳程度的增加,将导致司机持续性注意的大幅下。研究人员的研究也再次证明了脑电信号作为直接反应大脑活动的神经生理信号,与司机的当前精神状态具有高度相关性,当大脑处于较高唤醒水平时,脑电信息熵越大,相反当大脑处于较低唤醒水平时,脑电信息熵越小。

脑电信号作为一种神经电信号与大脑活动状态相关,能够有效反映驾驶员精神状态。因此,亟待研究和提出一种持续性注意水平的敏感性特征指标的选取方法,从而为动车司机持续性注意水平的量化测评及后期持续性注意水平识别模型构建提供依据。



技术实现要素:

有鉴于此,本发明提供一种持续性注意水平的敏感性特征指标的选取方法,从而实现了敏感于动车司机的持续性注意水平脑电特征指标的选取。

本发明的技术方案具体是这样实现的:

一种持续性注意水平的敏感性特征指标的选取方法,该方法包括:

通过多个电极采集脑电信号,记录与所述脑电信号对应的反应时间参数,并采集主观疲劳测评数据和行为绩效测评数据;

根据采集到的脑电信号计算得到脑电熵值参数;

根据所述脑电熵值参数,采用kruskal-wallis检验方法与relief算法,分别选取敏感于不同持续性注意水平的指标集合;

将得到的两个指标集合的交集g作为持续性注意水平的敏感性特征指标。

较佳的,所述根据采集到的脑电信号计算得到脑电熵值参数包括:

将采集到的脑电信号划分成多个时间长度为k分钟的分析单元,对每个分析单元的脑电信号以0~35hz的带宽进行整体滤波处理;

对每个分析单元的脑电信号,以预设步长为时间窗从左到右逐段滑动,并以预设的时间窗重叠率将一个分析单元的脑电信号分割成多个时间窗信号;

将每一个时间窗信号内乘于等长度的汉明窗,得到中间变量f[n];

对f[n]进行快速傅里叶变换,得到脑电信号在频域的幅值分布f(k);

从幅值分布f(k)中分别提取预设频段的平均幅值;对于每一个分析单元的脑电信号,当时间窗从前至后逐段滑动时,得各个预设频段各自的幅值序列;

对于各频段的幅值序列,去除正负3倍标准差外的异常数据后,得到长度都为l的各频段幅值序列集合;

根据各频段的幅值序列集合,计算得到各频段的样本熵、近似熵和香农熵,作为分析单元的脑电信号的脑电熵值参数。

较佳的,所述k的值为1;所述预设步长为2秒;所述时间窗重叠率为50%。

较佳的,所述幅值分布f(k)为:

其中,wn=cos(2π/n)-jsin(2π/n)。

较佳的,所述预设频段为θ、α和β频段;

所述θ频段为4~8hz,α频段为8~13hz,β频段为13~30hz。

较佳的,所述根据各频段的幅值序列集合,计算得到各频段的样本熵、近似熵和香农熵包括:

步骤31,对于每一个频段,将该频段的幅值序列由一维向量构造成m维向量ym(i);

步骤32,计算任意两向量对应元素之间的绝对差d,将其中最大的绝对差dmax作为两向量间的距离;

步骤33,根据预设阈值r,对每个向量ym(i)统计dmax≤r×sd的数目,并计算得到该数目与距离总数(l-m+1)的比值其中,sd为幅值序列的标准差;

步骤34,将取对数,根据如下的公式计算得到其均值bm(r):

步骤35,将该频段的幅值序列由原来的一维向量构造成m+1维向量,根据上述步骤31~35得到bm+1(r);

步骤36,将bm(r)和bm+1(r)代入如下的公式中,计算得到该频段的样本熵sampen和近似熵apen:

sampen(m,r)=-ln[bm(r)/bm+1(r)];

apen(m,r)=bm(r)-bm+1(r);

步骤37,对原有幅值序列按升序排列,其中每个元素所出现的概率分布为{p1,p2,...,pj},代入如下的公式中,得到该频段的香农熵h(p):

步骤38,对各个频段均按照上述步骤31~37得到各个频段的样本熵、近似熵与香农熵,将所得到的样本熵、近似熵和香农熵作为分析单元的脑电信号的脑电熵值参数。

较佳的,所述采用kruskal-wallis检验方法选取敏感于不同持续性注意水平的指标集合包括:

将所述脑电熵值参数组成样本组集合其中表示第k组样本中第t个参数,将样本组集合中的各样本组进行混合并按升序排列;

求取各变量的秩,构造统计量h:

其中,n为样本组数;q为参数样本的总量;rk为第k组的样本中的秩总和;nk为第k组的样本量;

在给定的显著性水平α下,计算统计量h的p值;若p>α,则接受原假设,表明各组脑电熵值参数样本无显著差异;若p<α,则拒绝原假设,表明各组脑电熵值参数样本之间存在显著差异;

选取差异性最为显著的h项脑电熵值参数作为持续性注意水平敏感性指标,即指标集合z={z1,z2,...,zh}。

较佳的,所述采用relief算法选取敏感于不同持续性注意水平的指标集合包括:

将所述脑电熵值参数组成样本集合其中任意两个样本xi与xj在脑电熵值参数r(1≤r≤9q)上的差为:

其中,xmax与xmin为脑电熵值参数r(1≤r≤9q)在样本集合x中的最大值与最小值;

从所述样本集合x中随机抽取一个样本xk,对与xk同类及异类的样本,均根据d值按升序排序,并分别选取前z个样本,用m(k)表示与xk同类的z个样本组成的集合,用m(c)表示与xk异类且属于c类的z个样本组成的集合;

使用如下的公式更新脑电熵值参数r的权重wr:

其中,wr(i)表示第i(1<i≤n)次迭代脑电熵值参数r的权重;p(·)为所属类型的样本个数在样本总数所占的比率;n为迭代次数;

经过n次迭代之后,求取各个脑电熵值参数相应权重的均值,选取其中权重最大的h项脑电熵值参数作为持续性注意水平敏感性指标,即指标集合z'={z1',z'2,...,z'h}。

如上可见,在本发明中的持续性注意水平的敏感性特征指标的选取方法中,由于先通过多个电极采集脑电信号,记录与所述脑电信号对应的反应时间参数,并采集主观疲劳测评数据和行为绩效测评数据;根据采集到的脑电信号计算得到脑电熵值参数,然后采用kruskal-wallis检验方法与relief算法,分别选取敏感于不同持续性注意水平的指标集合,最终持续性注意水平的敏感性特征指标,从而实现了敏感于动车司机的持续性注意水平脑电特征指标的选取,为动车司机持续性注意水平的量化测评及后期持续性注意水平识别模型构建提供了依据。

附图说明

图1为本发明实施例中的持续性注意水平的敏感性特征指标的选取方法的流程示意图。

图2为本发明实施例中的脑电熵值参数的示意图。

图3为本发明实施例中的主观行为与绩效得分的示意图。

图4为本发明实施例中p值分布脑地形图和前10项敏感性指标示意图。

图5为本发明实施例中电极位置分布图和前10项敏感性指标的示意图。

具体实施方式

为使本发明的目的、技术方案及优点更加清楚明白,以下参照附图并举实施例,对本发明进一步详细说明。

本实施例提供了一种持续性注意水平的敏感性特征指标的选取方法。

图1为本发明实施例中的持续性注意水平的敏感性特征指标的选取方法的流程示意图。如图1所示,本发明实施例中的持续性注意水平的敏感性特征指标的选取方法主要包括如下所述的步骤:

步骤11,通过多个电极采集脑电信号,记录与所述脑电信号对应的反应时间参数,并采集主观疲劳测评数据和行为绩效测评数据。

在本发明的技术方案中,为了测试驾驶员对突发时间的简单反应时间,需要对被测试人员(例如,志愿者)进行相应的测试,从而获得相应的脑电信号、与脑电信号对应的反应时间参数以及主观疲劳及行为绩效测评数据。

例如,较佳的,在本发明的较佳实施例中,可以先召集多个被测试人员,例如,在本发明的一个具体实验中,选取了22名中国铁路总公司93、94期动车司机班学员(年龄为36±2岁,驾龄为6年以上,累积安全里程大于10万公里)作为被测试人员。所选被测试人员无任何精神或心理疾病史,且视力或矫正视力正常,实验前24h禁止被测试人员饮用含有酒精或咖啡因的功能性饮料。

在本发明的一个具体实验中,采用crh1e动车模拟器,该模拟器采用单通道大屏前向视景系统,其屏幕分辨率为2048×768pix。模拟器声音环境由7.1的数字音频发生系构成,该系统可高保真模拟仿真动车运行时的背景音环境。机车操作台包括atp和pis系统、列车运行控制监视装置(lkj2000)、cir无线电通讯和后视视频摄像系统以及带按钮,开关和指示灯的控制面板。

实验线路选取北京南至上海虹桥高铁线路,线路全长为1463km,最高限速为220km/h,要求司机以低于限速5km/h的速度(标准速度)保持动车运行。

对随机探测信号的反应时间是衡量持续性注意水平的有效指标,为实时探测连续驾驶作业中司机的持续性注意水平,实验采用在前方呈现随机红点的方式来模拟列车运行中随机事件的发生。刺激红点在5个可能位置按照随机时间间隔(例如,240±10s)出现。要求被测试人员在红点出现后立即点击反应键,系统自动记录被测试人员反应时间。整个实验持续4小时,全程无休息。

为保证在实验开始前被测试人员具有较高的持续性注意水平,实验均安排在早上8:00进行。为了消除练习效应的影响,实验开始前预先进行15min的模拟驾驶练习,然后开始4小时的正式实验。在实验进行至46min及实验结束后,要求司机填写kss(karolinskasleepinessscale)量表,即采用kss量表对被测试人员的疲劳程度进行测量,得到主观疲劳测评数据。

较佳的,在本发明的具体实施例中,所述kss量表为9分量表,用于评估司机的疲劳状态:1-极度清醒;2-非常清醒;3-清醒;4-有些清醒;5-既不清醒也不困倦;6-开始出现困倦的征兆;7-困倦,容易控制;8-困倦,困难但可以控制;9-昏昏欲睡,无法控制。主观疲劳量表可以分别于第46min及实验结束后测试2次。

较佳的,在本发明的具体实施例中,所述驾驶行为的行为绩效测评数据包括:当前限速与列车当前运行速度(采样频率10hz),以及对随机探测事件(红点)的反应时间及正确率。

较佳的,在本发明的具体实施例中,脑电信号(eeg)可以由32导脑电仪连续采集,采用10~20导联头皮电极系统,记录水平与垂直眼电。选取fcz电极作为参考电极,脑电采样率设置为1000hz,采集频率带宽设置为0.5-100hz,保持所有电极阻抗小于5kω。

另外,在本发明的技术方案中,还可以对持续性注意水平进行划分。

例如,较佳的,在本发明的具体实施例中,可以选取上述实验中前后两个时间段作为高、低持续性注意水平的划分。该划分方式是人因工程学中的常用方法,在此不再赘述。在本发明的技术方案中,还可采用选取前后不同时段的行为指标与主观疲劳测评对驾驶持续性注意水平进行差异性验证。因此,在本发明的较佳实施例中,可以在上述实验中选取16-46min,210-240min两个时间段来对应动车司机不同的持续性注意水平。对随机探测信号的反应时间是衡量持续性注意水平的有效指标,因此,本发明可以通过上述前后两个时段司机对随机信号的反应时间的差异性对比分析,来验证该持续性注意水平划分的合理性。

步骤12,根据采集到的脑电信号计算得到脑电熵值参数。

在本发明的技术方案中,可以使用多种方式根据采集到的脑电信号计算得到脑电熵值参数。以下将以其中的一种具体实现方式为例,对本发明的技术方案进行详细的介绍。

例如,较佳的,在本发明的具体实施例中,所述根据采集到的脑电信号计算得到脑电熵值参数包括:

步骤121,将采集到的脑电信号划分成多个时间长度为k分钟的分析单元,对每个分析单元的脑电信号以0~35hz的带宽进行整体滤波处理,以去除工频电及部分肌电等伪迹干扰。

例如,如果采集到的脑电信号时间长度(简称时长)为(m*k)分钟,则可以将该采集到的脑电信号划分成m个分析单元,每个分析单元的脑电信号的时长为k分钟。

在本发明的技术方案中,所述k的取值可以根据实际应用情况的需要而预先设置。例如,较佳的,在本发明的具体实施例中,所述k的值为1。

步骤122,对每个分析单元的脑电信号,以预设步长为时间窗从左到右逐段滑动,并以预设的时间窗重叠率将一个分析单元的脑电信号分割成多个时间窗信号。

在本发明的技术方案中,所述预设步长的取值可以根据实际应用情况的需要而预先设置。例如,较佳的,在本发明的具体实施例中,所述预设步长为2秒。

在本发明的技术方案中,所述时间窗重叠率的取值也可以根据实际应用情况的需要而预先设置。例如,较佳的,在本发明的具体实施例中,所述时间窗重叠率为50%。

因此,例如,当分析单元的时长为1分钟、预设步长为2000毫秒、时间窗重叠率为50%时,则可以将一个分析单元的脑电信号分割成60个时间窗信号,每个时间窗信号的脑电信号的时长为2000毫秒,且各个时间窗信号之间的重叠率为50%。

步骤123,将每一个时间窗信号内乘于等长度的汉明窗,得到中间变量f[n]。

在本发明的技术方案中,为了消除旁瓣效应对快速傅里叶变换(fft)的影响,可以将每一个时间窗信号都内乘于等长度的汉明窗,得到中间变量f[n]。

步骤124,对f[n]进行快速傅里叶变换,得到脑电信号在频域的幅值分布f(k)。

例如,较佳的,在本发明的具体实施例中,所述幅值分布f(k)为:

式中,wn=cos(2π/n)-jsin(2π/n)。

步骤125,从幅值分布f(k)中分别提取预设频段的平均幅值;对于每一个分析单元的脑电信号,当时间窗从前至后逐段滑动时,得各个预设频段各自的幅值序列。

在本发明的技术方案中,所述预设频段的取值可以根据实际应用情况的需要而预先设置。例如,较佳的,在本发明的具体实施例中,所述预设频段为θ、α和β频段;所述θ频段为4~8hz,α频段为8~13hz,β频段为13~30hz。

因此,例如,在本发明的较佳实施例中,当一个分析单元的脑电信号为时长为1分钟的脑电信号、各个时间窗信号之间的重叠率为50%时,可以得到上述θ、α和β等3个频段各自的幅值序列。

步骤126,对于各频段的幅值序列,去除正负3倍标准差外的异常数据后,得到长度都为l的各频段(例如,θ、α和β频段)幅值序列集合(例如,θ、α和β频段的幅值序列集合{uθ(j),uα(j),uβ(j);1≤j≤l})。

步骤127,根据各频段的幅值序列集合,计算得到各频段的样本熵、近似熵和香农熵,作为分析单元的脑电信号的脑电熵值参数。

在本发明的技术方案中,可以先根据上述各频段的幅值序列集合,计算得到各频段的样本熵、近似熵和香农熵,然后再将所述样本熵、近似熵和香农熵作为分析单元的脑电信号的脑电熵值参数。

较佳的,在本发明的具体实施例中,所述根据各频段的幅值序列集合,计算得到各频段的样本熵、近似熵和香农熵包括:

步骤31,对于每一个频段,将该频段的幅值序列由一维向量构造成m维向量ym(i);

例如,在本发明的较佳实施例中,当当前频段为θ频段时,可以通过如下所述的公式将该频段的幅值序列由一维向量构造成m维向量:

ym(i)=[uθ(i),uθ(i+1),...,uθ(i+m-1)](2)

其中,{uθ(j),1≤j≤l}为θ频段的幅值序列集合,l为θ频段的幅值序列集合的长度,ym(i)为m维向量,i=1,2,...,l-m。

其中,m维向量ym(1),ym(2),...,ym(l-m),代表着从第i个点开始连续的m个uθ值。

对于其它频段的幅值序列,可以以此类推,在此不再赘述。

步骤32,计算任意两向量对应元素之间的绝对差d,将其中最大的绝对差dmax作为两向量间的距离。

步骤33,根据预设阈值r,对每个向量ym(i)统计dmax≤r×sd的数目,并计算得到该数目与距离总数(l-m+1)的比值其中,sd为幅值序列的标准差。

步骤34,将取对数,根据如下所述的公式计算得到其均值bm(r):

步骤35,将该频段的幅值序列由原来的一维向量构造成m+1维向量,根据上述步骤31~35得到bm+1(r)。

步骤36,将bm(r)和bm+1(r)代入如下所述的公式中,计算得到该频段的样本熵sampen和近似熵apen:

sampen(m,r)=-ln[bm(r)/bm+1(r)](4)

apen(m,r)=bm(r)-bm+1(r)(5)

步骤37,对原有幅值序列按升序排列,其中每个元素所出现的概率分布为{p1,p2,...,pj},代入如下所述的公式中,得到该频段的香农熵h(p):

步骤38,同理,对各个频段均按照上述步骤31~37得到各个频段的样本熵、近似熵与香农熵,将所得到的样本熵、近似熵和香农熵作为分析单元的脑电信号的脑电熵值参数。

通过上述的步骤31~38,即可得到分析单元的脑电信号的脑电熵值参数。

对于每个分析单元,都可以通过如上所述的步骤121~127,得到各个分析单元的脑电信号的脑电熵值参数。

例如,在本发明的较佳实施例中,对于通过q个电极采集到的脑电信号,通过上述的步骤121~127,最终可得到9q项脑电熵值参数(即持续性注意熵值参数),记为xi(i=1,2,...,9q),如图2所示。

步骤13,根据所述脑电熵值参数,采用kruskal-wallis检验方法与relief算法,分别选取敏感于不同持续性注意水平的指标集合。

在本发明的技术方案中,所述kruskal-wallis检验方法与relief算法为现有技术中的常用方法。

具体来说,在本发明的较佳实施例中,采用kruskal-wallis检验方法选取敏感于不同持续性注意水平的指标集合包括:

将所述脑电熵值参数组成样本组集合其中表示第k组样本中第t个参数,将样本组集合中的各样本组进行混合并按升序排列;

求取各变量的秩,构造统计量h:

其中,n为样本组数;q为参数样本的总量;rk为第k组的样本中的秩总和;nk为第k组的样本量;

在给定的显著性水平α下,计算统计量h的p值;若p>α,则接受原假设,表明各组脑电熵值参数样本无显著差异;若p<α,则拒绝原假设,表明各组脑电熵值参数样本之间存在显著差异;

经过上述计算之后,即可选取差异性最为显著的h项脑电熵值参数作为持续性注意水平敏感性指标,即指标集合z={z1,z2,...,zh}。

具体来说,在本发明的较佳实施例中,采用relief算法选取敏感于不同持续性注意水平的指标集合包括:

将所述脑电熵值参数组成样本集合其中任意两个样本xi与xj在脑电熵值参数r(1≤r≤9q)上的差为:

其中,xmax与xmin为脑电熵值参数r(1≤r≤9q)在样本集合x中的最大值与最小值;

从所述样本集合x中随机抽取一个样本xk,对与xk同类及异类的样本,均根据d值按升序排序,并分别选取前z个样本,用m(k)表示与xk同类的z个样本组成的集合,用m(c)表示与xk异类且属于c类的z个样本组成的集合;

使用如下所述的公式更新脑电熵值参数r的权重wr:

其中,wr(i)表示第i(1<i≤n)次迭代脑电熵值参数r的权重;p(·)为所属类型的样本个数在样本总数所占的比率;n为迭代次数;

经过n次迭代之后,求取各个脑电熵值参数相应权重的均值,选取其中权重最大的h项脑电熵值参数作为持续性注意水平敏感性指标,即指标集合z'={z1',z'2,...,z'h}。

步骤14,将得到的两个指标集合的交集g作为持续性注意水平的敏感性特征指标:

g={(z∩z')=(z1,z2,...,ze)}(10)

通过上述的步骤11~14,即可通过选取得到持续性注意水平的敏感性特征指标。

另外,在本发明的技术方案中,可以通过实验验证持续性注意水平的划分的合理性。

例如,在本发明的较佳实施例中,为了验证上述步骤11中的较佳实施例中选取16-46min,210-240min2个时间段对应动车司机不同持续性注意水平假设的合理性,可以对上述两时段的主观疲劳测评数据、3项行为数据(反应时、速度偏差与有效检测率)进行计算及配对样本t检验,其结果如图3所示。

结果表明,主观疲劳kss得分在第2时段显著高于第1时段(2.13vs7.23,t(21)=-6.86,p<.01);在速度偏差方面,第2时段的偏差显著高于第1时段(3.22vs5.13,t(21)=-5.21,p<.01);在对随机探测信号反应方面,第2时段的反应时显著高于第2时段(565.43vs602.58,t(21)=-3.11,p<.05);第2时段的有效检测率显著低于第1时段(90.62vs78.86,t(21)=2.52,p<.05)。

对随机信号的反应时间和有效检测率是评估持续性注意能力的有效指标,与第1时段(16-46min)相比,动车司机在第2时段(210-240min)次任务中对红点随机信号检测能力显著下降(反应时增加,有效检测率下降),这表明本实验中以第1时段(16-46min)作为高持续性注意水平、第2时段(210-240min)作为低持续性注意水平的划分是合理的。

另外,在本发明的技术方案中,还可以

对于实验所采集的脑电信号,在使用步骤12中的处理方法对其进行快速的傅里叶变换,然后对θ(4~8hz)、α(8~13hz)和β(13~30hz)3个频段进行样本熵、近似熵与香农熵处理,共可得到288组特征指标。采用上述的kruskal-wallis检验方法,得出在32个电极上两阶段脑电数据间的p值,如图4(a)所示;由此选出对持续性注意水平最敏感的10个指标,如图4(b)所示。

根据上述的relief算法,可得到各个脑电熵指标的权重均值,选取权重值最大的10个熵指标作为持续性注意水平的敏感性指标。其电极位置分布如图5(a)所示,前10项敏感性指标如图5(b)所示。

以kruskal-wallis检验和relief算法筛选所得到的敏感性指标交集作为持续性注意的特征指标,即位于fp1、f7电极处β波段的香农熵h(p)以及位于fz电极α波段的样本熵sampen。

综上可知,在本发明中的持续性注意水平的敏感性特征指标的选取方法中,由于先通过多个电极采集脑电信号,记录与所述脑电信号对应的反应时间参数,并采集主观疲劳测评数据和行为绩效测评数据;根据采集到的脑电信号计算得到脑电熵值参数,然后采用kruskal-wallis检验方法与relief算法,分别选取敏感于不同持续性注意水平的指标集合,最终持续性注意水平的敏感性特征指标,从而实现了敏感于动车司机的持续性注意水平脑电特征指标的选取,为动车司机持续性注意水平的量化测评及后期持续性注意水平识别模型构建提供了依据;并且通过前后两时段的主观疲劳测评数据与对随机信号的反应时差异性对比分析,验证了长时间连续动车驾驶作业会引发动车司机精神疲劳及持续性注意水平下降。

以上所述仅为本发明的较佳实施例而已,并不用以限制本发明,凡在本发明的精神和原则之内,所做的任何修改、等同替换、改进等,均应包含在本发明保护的范围之内。

当前第1页1 2 
网友询问留言 已有0条留言
  • 还没有人留言评论。精彩留言会获得点赞!
1