基于电子化手绘螺旋测试的智能决策支持方法及系统

文档序号:29858015发布日期:2022-04-30 09:57阅读:134来源:国知局
基于电子化手绘螺旋测试的智能决策支持方法及系统

1.本发明涉及机器学习和人机交互技术领域,具体地说,本发明提出一种基于电子化手绘螺旋测试的智能决策支持方法及系统。


背景技术:

2.帕金森病作为当前世界第二大神经退行性疾病,在全球的患者人数达到七百万到一千多万之多。每一万个80岁以上的人中,就有大约190个患有帕金森病,并且帕金森病的患病率会随年龄增长而逐渐提高。目前临床主要由专科医生观察评估震颤、运动迟缓和强直等的主要症状,从而进行帕金森病的诊断。但是由于帕金森病在早期并未有明显表现,大部分帕金森病患者在疾病早期难以意识患有帕金森病,并且诊断易受到医生主观影响,导致误诊率较高,从而延误了早期的药物干预治疗。
3.目前临床有效诊断帕金森病的方式之一是观测患者绘制螺旋线的表现(如图1所示)。然而,传统的螺旋线绘制都是使用纸和笔,使得许多重要时间特征(例如绘制螺旋线时的平均速度、最大速度、速度的方差、平均加速度、最大加速度和加速度方差等)无法获得精准记录。其次,在2018年发表的研究报告《let’s draw:detecting and measuring parkinson’s disease on smartphones》中提到,受试者在重复进行单一测试时,手部的震颤有逐渐减轻地趋势,传统方式受到纸笔的限制,测试类别较为单一。另外,传统方式需要患者在特定人员的监督下完成,诊断受到时间和空间条件的约束,使用纸质材料记录的数据也不利于传播和保存。


技术实现要素:

4.针对上述问题,本发明提出了一种面向帕金森病辅诊的电子化手绘螺旋测试方法,该方法首先电子化显示螺旋线,并引导受试者完成两种螺旋线测试——静态螺旋测试sst和动态螺旋测试dst,同时采集相关数据,并基于采集到的数据提取特征,最后使用特征训练模型用于帕金森病辅助诊断。
5.具体而言,利用触摸屏显示引导用户描绘的螺旋参考线(分为静态和动态两种方式),患者利用手指绘制过程中记录轨迹坐标;并根据螺旋线的几何结构,提取表征帕金森病症状的有效特征,利用随机森林算法构建分类模型,实现高精准、强鲁棒的帕金森病分类结果作为确诊过程的重要辅助参考。
6.针对现有技术的不足,本发明提出一种基于电子化手绘螺旋测试的智能决策支持方法,其中包括:
7.步骤1、获取包括多条手绘螺旋轨迹的训练数据,且手绘螺旋轨迹具有其是否属于帕金森的类别标签;
8.步骤2、提取该训练数据中每条手绘螺旋轨迹的特征,以结合该类别标签,训练随机森林模型,将训练完成后的随机森林模型作为智能决策支持模型;
9.步骤3、将待分类的手绘螺旋轨迹的特征输入该智能决策支持模型,得到其所属的
类别,并将其作为待分类的手绘螺旋轨迹的决策支持结果。
10.所述的基于电子化手绘螺旋测试的智能决策支持方法,其中该手绘螺旋轨迹包括动态螺旋测试轨迹,该动态螺旋测试轨迹的提取过程为在触屏上生成闪烁的动态参考螺旋线,感测并存储受试者根据该动态参考螺旋线所绘制的坐标及各坐标对应的时间戳,作为该动态螺旋测试轨迹。
11.所述的基于电子化手绘螺旋测试的智能决策支持方法,其中该手绘螺旋轨迹包括静态螺旋测试轨迹,该静态螺旋测试轨迹的提取过程为在触屏上生成常亮的静态参考螺旋线,感测并存储受试者根据该静态参考螺旋线所绘制的坐标及各坐标对应的时间戳,作为该动态螺旋测试轨迹。
12.所述的基于电子化手绘螺旋测试的智能决策支持方法,其中该步骤2所提取的特征包括:平均速度,和/或最大速度,和/或速度方差,和/或平均加速度,和/或最大加速度,和/或加速度方差,和/或斜率变化率,和/或位移偏差均值,和/或位移偏差最大值,和/或位移偏差方差,和/或线圈间距均值,和/或线圈间距方差,和/或线圈间距最大值。
13.本发明还提出了一种基于电子化手绘螺旋测试的智能决策支持系统,其中包括:
14.初始模块,用于获取包括多条手绘螺旋轨迹的训练数据,且手绘螺旋轨迹具有其是否属于帕金森的类别标签;
15.训练模块,用于提取该训练数据中每条手绘螺旋轨迹的特征,以结合该类别标签,训练随机森林模型,将训练完成后的随机森林模型作为智能决策支持模型;
16.分类模块,用于将待分类的手绘螺旋轨迹的特征输入该智能决策支持模型,得到其所属的类别,并将其作为待分类的手绘螺旋轨迹的决策支持结果。
17.所述的基于电子化手绘螺旋测试的智能决策支持系统,其中该手绘螺旋轨迹包括动态螺旋测试轨迹,该动态螺旋测试轨迹的提取过程为在触屏上生成闪烁的动态参考螺旋线,感测并存储受试者根据该动态参考螺旋线所绘制的坐标及各坐标对应的时间戳,作为该动态螺旋测试轨迹。
18.所述的基于电子化手绘螺旋测试的智能决策支持系统,其中该手绘螺旋轨迹包括静态螺旋测试轨迹,该静态螺旋测试轨迹的提取过程为在触屏上生成常亮的静态参考螺旋线,感测并存储受试者根据该静态参考螺旋线所绘制的坐标及各坐标对应的时间戳,作为该动态螺旋测试轨迹。
19.所述的基于电子化手绘螺旋测试的智能决策支持系统,其中该训练模块所提取的特征包括:平均速度,和/或最大速度,和/或速度方差,和/或平均加速度,和/或最大加速度,和/或加速度方差,和/或斜率变化率,和/或位移偏差均值,和/或位移偏差最大值,和/或位移偏差方差,和/或线圈间距均值,和/或线圈间距方差,和/或线圈间距最大值。
20.本发明还提出了一种存储介质,用于存储执行所述任意一种基于电子化手绘螺旋测试的智能决策支持系方法的程序。
21.本发明还提出了一种客户端,用于上述任意一种基于电子化手绘螺旋测试的智能决策支持系系统。
22.由以上方案可知,本发明的优点在于:
23.本发明提出一种面向帕金森病辅诊的电子化手绘螺旋测试方法。该方法提供电子化螺旋线测试,实现快速、简单、低成本的帕金森辅助诊断,解决了传统帕金森诊断特征维
度低、测试形式单一、测试结果数据不易保存和传播、测试受到时间与空间条件约束等问题。从更多维特征中分析体现受试者手部震颤的表征,构建出高精度和强鲁棒性的识别模型。
附图说明
24.图1为现有技术的纸绘螺旋线示意图;
25.图2为本发明的流程图;
26.图3为基准阿基米德螺旋线。
具体实施方式
27.针对传统手绘螺旋测试形式单一、时间特征缺失、受时间与空间约束、数据不易保存等问题,本发明提出了一种面向帕金森病辅诊的电子化手绘螺旋测试方法。如图2所示,该方法主要包括数据采集、特征提取和模型训练三个部分。数据采集利用电子化方式,记录多维度数据、提供多样性测试提供(包括静态螺旋测试和动态螺旋测试);特征提取用于提取能够表征帕金森病症状的时间和空间特征;模型训练部分采用随机森林模型,利用提取的特征构建帕金森病辅助诊断模型。
28.1.数据采集
29.本发明中需要采集的数据分别来自两种测试,第一个是静态螺旋测试(static spiral test,sst),第二个是动态螺旋测试(dynamic spiral test,dst)。这两种测试使用完全一致的螺旋线,区别在于dst中使用的螺旋线模板每隔2秒钟才闪烁出现一次,受试者需要凭借记忆在螺旋线模板消失期间继续绘制螺旋线。
30.电子化手绘螺旋线的数据采集分为两个步骤:(1)首先在电子触摸屏上生成用于受试者临摹的参考螺旋线,由于不同设备的尺寸及分辨率不完全一致,为了保证在不同平台上采集到的螺旋线数据(例如速度、位移、间距等)的一致性,需要在不同设备端生成统一大小、形状、间距、圈数的参考螺旋线;(2)受试者以中心点为起始点,用食指根据参考线连续地进行螺旋线绘制,同时按照固定频率存储受试者所绘制轨迹的坐标及相应的时间戳。
31.2.特征提取
32.本发明利用坐标数据提取以下特征:平均速度、最大速度、速度方差、平均加速度、最大加速度、加速度方差、斜率变化率、位移偏差均值、位移偏差最大值、位移偏差方差、线圈间距均值、线圈间距方差、线圈间距最大值。
33.3.模型训练
34.使用提取到的特征,训练随机森林模型。使用模型精度和特征相关性验证特征的有效性,使用留一法得到的方差验证模型的鲁棒性。
35.为让本发明的上述特征和效果能阐述的更明确易懂,下文特举实施例,并配合说明书附图作详细说明如下。
36.1、数据采集
37.在发明内容一节中,我们提到采集到的数据分别来自dst和sst,由于两测试使用的螺旋线完全一致,对本小节中关于螺旋线的绘制、电子化显示,和数据采集的步骤这三个方面做统一阐述。
38.本发明使用阿基米德螺旋线作为受试者完成螺旋测试的临摹参考,应确保测试是在具有相同大小、形状、间距、圈数的参考螺旋线上完成的。因此,本节中,我们将围绕(1)阿基米德螺旋线地绘制、(2)在不同设备端显示、(3)数据采集步骤三个方面,分别阐述如何控制采集数据时其他变量的一致性。
39.阿基米德螺旋线的极坐标方程如式(1)所示:
40.r=a+b*θ
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(1)
41.其中,a为起始点与极坐标中心的距离,主要负责旋转整个螺旋线;b为控制螺旋线的螺距,b越大或变化越快时,螺旋线在相同角度下半径的r的增长越快,越稀疏;θ控制螺线圈的大小,θ越大则螺旋线的范围越大。
42.为了方便坐标表示和特征提取,结合式(1)将极坐标公式转化为直角坐标公式,见式(2)(3):
43.x=r*cosθ=(a+b*θ)*cosθ
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(2)
44.y=r*sinθ=(a+b*θ)*sinθ
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(3)
45.例如,本发明使用螺旋线坐标方程为式(4)(5):
46.x=(0+10.5*θ)*cosθ,θ∈[0,6π]
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(4)
[0047]
y=-(0+10.5*θ)*sinθ,θ∈[0,6π]
ꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(5)
[0048]
其中,a值为0,b值为10.5,θ的取值间隔为0.01π,范围在[0,6π]之间。
[0049]
例如,本发明中根据式(4)(5)绘制出的基准阿基米德螺旋线如图3所示,该螺旋线将作为受试者绘制时的模板:
[0050]
屏幕自适应螺旋线显示。本发明所使用的电子化螺旋测试可使用手机或者平板电脑等移动设备进行采集,针对不同设备在采集时会出现的分辨率不一致导致的图像大小、间距等发生改变的情况,我们便根据不同设备分辨率自适应生成统一大小的螺旋线图片。具体过程如下:
[0051]
a.获取屏幕分辨率。
[0052]
b.根据式(6)(7)计算屏幕实际可显示的长宽:由于设备的实际长宽与可显示的屏幕的长宽并不一致,所以需要计算出可显示的实际长宽。例如,本发明使用了两种设备,分别是mi pad 4,分辨率为1920*1200,和honor 30pro,分辨率为2340*1080。由于两设备大小不一,我们首先需要限制生成的螺旋线的大小(例如:2*2英寸)。为了使得两种设备上显示的螺旋线的物理特征(包括大小、形状、间距、圈数等)完全一致,我们使用式(8)(9)计算得到在不同的设备上显示相同的螺旋线时(例如:2*2英寸),需要占用的x和y方向上的分辨率的比例。使用该比例即可在不同设备端上显示出完全一致的螺旋线。
[0053]
width=resolution_width/xbpi
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(6)
[0054]
height=resolution_height/ybpi
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(7)
[0055]
其中,width表示实际可显示宽度,height表示实际可显示长度,resolution_width表示分辨宽,resolution_height表示分辨长,xbpi表示宽度(x)方向上每英寸的像素点个数,ybpi表示长度(y)方向上每英寸的像素点个数。
[0056]
x
ratio
=actual
width
/xbpi
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(8)
[0057]yratio
=actual
height
/ybpi
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(9)
[0058]
其中,x
ratio
和y
ratio
分别表示长(y)宽(x)方向上的所占比例,actual
width

actual
height
分别表示实际生成的螺旋线的长于宽,xbpi表示宽度(x)方向上每英寸的像素点个数,ybpi表示长度(y)方向上每英寸的像素点个数。
[0059]
数据采集。获得螺旋线之后,引导受试者从螺旋线的中心开始,使用食指沿顺时针方向绘测螺旋线。测试时的绘制时间戳和相应坐标,将按特定频率(如60hz)被记录,最终得到螺旋线数据集的组成包括时间(ti)和横纵轴坐标(xi,yi):
[0060]
s={i=0,

,t|ti,xi,yi}
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(10)
[0061]
其中,s为时间戳和坐标数据集,i为索引,t为总用时,ti为时间节点,(xi,yi)为i时刻的坐标。
[0062]
2、特征提取
[0063]
针对获取受试者绘制螺旋线数据,围绕手部震颤程度、手部发力能力和手部稳定性等方面进行特征提取。
[0064]
起始点偏移值:起始点偏移值表示受试者绘制螺旋线时的起始点位置距离螺旋线模板起始点的位移。若该值过大则可能表示患者手部抖动、震颤较为剧烈导致绘制食指无法有效对准螺旋线模板起始点。计算方法如下:
[0065][0066]
其中,dist表示距离,(xi,yi)为受试者绘制的螺旋线起始点坐标,(x0,y0)为螺旋线模板起始点坐标。
[0067]
螺旋线总体位移偏差的总和、均值、最大值、和方差:圈总偏移表示受试者绘制的螺旋线和基准螺旋线之间总体位移偏差,可以体现受试者在绘制螺旋线时跟随基准螺旋线的能力,螺旋线位移偏差由受试者所绘测的螺旋线的各点与其基准螺旋线上相对应的角度相同的最近点之间的欧几里得距离相加之和得到。计算步骤为:
[0068]
遍历坐标数据集,根据式(12)计算当前坐标(xi,yi)到基准螺旋线的最近距离。
[0069][0070]
其中,dist表示距离,n为总坐标个数,i为索引,(xi,yi)为当前坐标,(x
i’,y
i’)为具有相同角度的基准点坐标。
[0071]
累加各个坐标到基准螺旋线的距离即可得到总体位移偏差。
[0072]
根据式(13)计算平均位移。
[0073][0074]
其中,avg
dist
表示平均位移,dist为总位移,n为总坐标个数。
[0075]
根据式(14)计算位移方差。
[0076][0077]
其中,std
dist
表示位移方差值,n为总坐标个数,i为索引,disti为总位移,avg
dist
为位移平均值。
[0078]
螺旋线的线圈间隔距离的平均值、最大值和方差:线圈间距表示受试者绘制的螺
旋线之间的距离。通过分析螺旋线间距离可以观测螺旋线的平滑程度从而体现受试者绘制螺旋线时的手部稳定性和跟随螺旋线的能力。计算方法如下:
[0079]
遍历坐标数据集(x,y)={(xi,yi)|i=1,

,t},根据公式(11)(15)计算各点角度θi与半径ri。
[0080][0081]
其中,angle表示角度,(xn,yn)为当前坐标,(x0,y0)为起始点坐标。
[0082]
将具有相同角度值的半径集合得到:
[0083]
去除临近值:遍历字典中所有角度找出每一个角度对应的所有半径值r={r|r∈r},由于绘画时速度不一致,相同角度上可能对应多个临近的半径值,找出临近的半径值,并得到平均值。
[0084]
合成缺失值:遍历字典中所有角度遍历当前角度对应的所有半径值r={r|r∈r},若缺失半径值,遍历当前角度θ相邻的i={i|1,..,5}个角度所对应半径值,找到对应的半径缺失值。
[0085]
计算半径差值得到各线圈间距和其中的最大值。
[0086]
根据公式(13)计算平均间距。
[0087]
根据公式(14)计算间距方差。
[0088]
绘测平均速度、最大速度、速度方差,绘测速度表示的是受试者在描绘螺旋线时的快慢程度。若速度较慢可能表示其绘测时手部状态不稳定,需要放慢速度以控制手部发力来更好地跟随模板。计算方法如下:
[0089]
遍历坐标数据集,使用式(12)得到每五个点之间的总位移。
[0090]
由于采样频率一致,五个点之间所用时间相同,位移的大小即可用来衡量速度的快慢。原理见下式(16)(17)。
[0091]
velocity=dist/time
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
(16)
[0092][0093]
式(17)中,velocity表示速度,dist表示位移,time表示时间。
[0094]
式(18)中,v1,v2表示两速度,dist1,dist2表示两位移,t1,t2表示量时间。
[0095]
绘测加速度、最大加速度、加速度方差:绘测加速度表示的是受试者在描绘螺旋线时的加速度大小。加速度的变化幅度可以反映出受试者绘制螺旋时的手部稳定程度。计算方法如下:
[0096]
根据(4)中得到的位移改变量来衡量加速度。
[0097]
平均加速度、最大加速度和加速度方差均可用使用式(12)计算得到的位移改变量计算得出。
[0098]
编号特征描述1.maxdist最大位移偏差
2.stdevdist位移偏差方差3.avgdist位移偏差平均值4.stdevvelo_dstdst速度方差5.maxvelo_dstdst速度最大值6.stdevacc_dstdst加速度方差7.minacc最小加速度8.minacc_dstdst最小加速度9.stdevacc加速度方差10.stdevgap_dstdst间距方差11.stdevgap间距方差12.avggap间距均值13.maxgap_dstdst最大间距14.adjust_dstdst起点偏移
[0099]
表1提取的特征:在特征一列中带有_dst后标的为dst中提取的特征,无后标的为sst中提取的特征
[0100]
3、模型训练
[0101]
获得以上特征后,使用基尼指数,构建随机森林模型。
[0102]
实验验证:
[0103]
本节在一个基准数据集上,对比两种测试分别合成特征和结合两种测试构建模型时精度的差别,以下我们将在实验数据集、实验结果等方面进行阐述。
[0104]
数据集:
[0105]
实验验证所用数据集是来自伊斯坦布尔大学医学院公开的帕金森患者手绘螺旋图,该数据集包含两组对象,帕金森患者与正常对照。共有62名帕金森患者和15名健康人士参与螺旋测试中。记录有的数据包括横纵坐标、时间戳等属性。
[0106]
特征提取:根据坐标和时间戳这些数据,根据第2节中的特征提取方法,可以得到以下几个特征,包括最大位移偏差、位移偏差方差、位移偏差平均值、dst速度方差、dst速度最大值、dst加速度方差、最小加速度、dst最小加速度、加速度方差、dst间距方差、间距方差、间距均值、dst最大间距、dst起点偏移,见表1。
[0107]
实验结果:
[0108]
为了检验所提取特征的有效性,我们1)使用了最大互信息系数(mic)和皮尔森系数来检验所提取特征与结果的相关性,计算结果见表2。其中,皮尔森系数只能检测特征与标签的线性相关性;而mic系数具有较好的普适性、公平性和对称性。所谓普适性,是指在样本量足够大(包含了样本的大部分信息)时,能够捕获各种各样的有趣的关联,而不限定于特定的函数类型(如线性函数、指数函数或周期函数),或者说能均衡覆盖所有的函数关系。2)使用了多种模型验证特征的有效性。针对两种不同的测试(dst和sst)所提取的特征,我们将其分为了三类分别进行测试,包括单独使用sst特征、dst特征和二者结合的特征,使用多种模型精度来验证合成特征的有效性。以下表格(3)(4)(5)中的数据展示了使用不同特征集构建的多种模型的准确率、精度、回归率和受试者工作特征(即roc)曲线下面积。同时,我们使用了交叉验证方法再次验证特征的有效性和模型的鲁棒性。
[0109][0110][0111]
表2特征-标签相关性表格
[0112][0113]
表3sst特征模型实验结果
[0114][0115]
表4dst特征模型实验结果
[0116][0117]
表5合并特征模型实验结果
[0118]
在特征相关性验证方面,我们使用了最大信息系数(mic)和皮尔森系数来验证特征的相关性。在表2中我们可以看到,所提取特征的mic值都高于0.6,由于mic值的高鲁棒性,可以有效反映出所提取特征与标签的相关性。虽然,某些特征的皮尔森系数的绝对值小于0.5,其原因可能是皮尔森系数只能够检测与标签具有线性相关性的特征。
[0119]
在模型训练中,我们将特征分为三组(包括sst特征、dst特征,和两测试的结合特征),分别使用三种模型检验特征的有效性。在模型测试阶段,三种模型都使用了交叉验证法对模型的精度及泛化能力进行检验。实验结果显示在使用合并特征(sst+dst)进行预测的时候,三种模型的精度都是最高的,说明了动态螺旋测试在实验中的重要性,因此电子化螺旋测试具有比传统诊断更好的辅助诊断表现。同时,交叉验证的精度证实了特征的有效性,方差值证实了模型的高鲁棒性。
[0120]
以下为与上述方法实施例对应的系统实施例,本实施方式可与上述实施方式互相配合实施。上述实施方式中提到的相关技术细节在本实施方式中依然有效,为了减少重复,这里不再赘述。相应地,本实施方式中提到的相关技术细节也可应用在上述实施方式中。
[0121]
本发明还提出了一种基于电子化手绘螺旋测试的智能决策支持系统,其中包括:
[0122]
初始模块,用于获取包括多条手绘螺旋轨迹的训练数据,且手绘螺旋轨迹具有其是否属于帕金森的类别标签;
[0123]
训练模块,用于提取该训练数据中每条手绘螺旋轨迹的特征,以结合该类别标签,训练随机森林模型,将训练完成后的随机森林模型作为智能决策支持模型;
[0124]
分类模块,用于将待分类的手绘螺旋轨迹的特征输入该智能决策支持模型,得到其所属的类别,并将其作为待分类的手绘螺旋轨迹的决策支持结果。
[0125]
所述的基于电子化手绘螺旋测试的智能决策支持系统,其中该手绘螺旋轨迹包括动态螺旋测试轨迹,该动态螺旋测试轨迹的提取过程为在触屏上生成闪烁的动态参考螺旋线,感测并存储受试者根据该动态参考螺旋线所绘制的坐标及各坐标对应的时间戳,作为该动态螺旋测试轨迹。
[0126]
所述的基于电子化手绘螺旋测试的智能决策支持系统,其中该手绘螺旋轨迹包括静态螺旋测试轨迹,该静态螺旋测试轨迹的提取过程为在触屏上生成常亮的静态参考螺旋线,感测并存储受试者根据该静态参考螺旋线所绘制的坐标及各坐标对应的时间戳,作为该动态螺旋测试轨迹。
[0127]
所述的基于电子化手绘螺旋测试的智能决策支持系统,其中该训练模块所提取的特征包括:平均速度,和/或最大速度,和/或速度方差,和/或平均加速度,和/或最大加速度,和/或加速度方差,和/或斜率变化率,和/或位移偏差均值,和/或位移偏差最大值,和/或位移偏差方差,和/或线圈间距均值,和/或线圈间距方差,和/或线圈间距最大值。
[0128]
本发明还提出了一种存储介质,用于存储执行所述任意一种基于电子化手绘螺旋测试的智能决策支持系方法的程序。
[0129]
本发明还提出了一种客户端,用于上述任意一种基于电子化手绘螺旋测试的智能决策支持系统。
当前第1页1 2 
网友询问留言 已有0条留言
  • 还没有人留言评论。精彩留言会获得点赞!
1