999精品在线视频,手机成人午夜在线视频,久久不卡国产精品无码,中日无码在线观看,成人av手机在线观看,日韩精品亚洲一区中文字幕,亚洲av无码人妻,四虎国产在线观看 ?

基于時間序列數據挖掘的地鐵車門亞健康狀態識別方法

2018-05-21 01:01:18支有冉許志興
計算機應用 2018年3期
關鍵詞:亞健康特征

薛 鈺,梅 雪,支有冉,許志興,史 翔

(1.南京工業大學 電氣工程與控制科學學院,南京 211816; 2.南京康尼機電股份有限公司,南京 210013)

近年來,城市人口的急劇增加,城市規模的增大,致使城市的交通系統面臨嚴峻的局勢。地鐵作為公共交通中重要的組成部分,是目前人們主要出行方式之一。列車車門由于其工作環境中經常受到擠壓和震動,所以它是整個車輛中故障頻發的部分[1]。車門打開和關閉過程中最重要的驅動部件就是電機,其工作狀態可以直接或間接反映出車門的工作狀態。最近這些年人們提出多種方法用于電機的故障診斷,例如基于神經網絡的方法[2]、基于信號處理的方法[3]和基于規則發現的診斷方法[4]等。車門亞健康狀態不同于故障,是介于車門正常運行和故障之間的中間狀態,對其進行識別存在重要意義。

由于車門電機參數為時間序列相關的數據,故可以利用時間序列數據挖掘算法和技術來對電機的運行狀態進行分析。目前,時間序列常用的表示方法有離散小波變換、分段線性表示和符號化方法。文獻[5]采用離散小波變換處理地震波時間序列,將轉化后小波的系數作為特征,實現對地震信號的跟蹤。文獻[6]等依據時間序列中斜率的變化,提出了基于一階濾波的時間序列分段線性表示方法(Piecewise Linear Representation method of time Series based on First-order Filtering, PLR_SFWF),將時間序列分段線性表示。文獻[7]挖掘字符串處理后時間序列中的頻繁模式,實現時間序列的聚類以及分類。雖然上述的算法對時間序列處理效果很好,但遺憾的是運算比較復雜,難以應用于地鐵車門電機數據處理中。本文針對地鐵門電機數據提出了一種新的亞健康識別算法。首先,采用多尺度滑動窗口的方法并結合拓展符號聚集近似(Extension of Symbolic Aggregate approXimation, ESAX)算法對地鐵門電機數據進行挖掘;然后,采用主成分分析(Principal Component Analysis, PCA)法降維并選擇較為敏感的參數作為亞健康判斷的特征量;最后,結合基礎特征利用分層模式識別模型對各類亞健康狀態進行識別,并以實測地鐵門電機數據為例驗證本文所提方法的有效性。

1 地鐵門電機數據特征提取

1.1 時間序列符號化

Lin等[8]提出符號化聚合近似(Symbolic Aggregate approXimation, SAX)算法,對于長度為n的時間序列X,若要將其用長度為w的字符串表示,需要進行下面四個步驟:

1)對數據X進行標準化處理。對于標準化后的時間序列,數據往往呈現高斯分布,計算公式如下:

(1)

3) 依據間斷點β={β1,β2,…,βα-1}將高斯空間等概率劃分,其中的區間個數為α。間斷點可以查表獲得,表1給出了α從3到5時的間斷點的值。

表1 高斯空間等概率劃分斷點(α為區間個數, β為斷點)Tab. 1 Breakpoints that divide a Gaussian distribution into equiprobable regions (α is the number of regions, β are the breakpoints)

1.2 改進時間序列符號化

SAX在應用中的良好效果使其成為近10年來被廣泛使用的一種符號化特征表示方法;然而,只計算均值往往會忽略時間序列中的極值點。圖1為車門在剛啟動階段時的電流曲線,通過SAX算法可以將其離散化為字符串DCBCDE(a=6,w=6)??梢园l現,圖1中很多重要的極值點信息被忽略了,這些極值點記錄車門運行過程中一些極端且不尋常的模式,這些模式往往是進行識別的關鍵。文獻[9]提出了拓展符號聚集近似(ESAX)算法并將其應用于金融數據處理,其思想是在計算均值的時候同時引入最大值和最小值,以此來保留時間序列數據中的極值點。

圖1 通過SAX字符化車門電機電流曲線 Fig. 1 Current data of motor represented by SAX

圖2為使用ESAX算法對電流曲線進行離散化。當α=6,w=6時,曲線離散字符化結果為BDFECBABBBCCCDDDEE。這種改進的字符化算法能夠保留圖中的極值點信息,缺點是使得原來長度為6的字符串變為長度為18的字符串。

圖2 通過ESAX字符化車門電機電流曲線 Fig. 2 Current curve of motor represented by ESAX

為了使ESAX算法進一步適用于車門電機數據,對其作出改進,結合滑動窗口重疊分割法來離散時間序列[10]?;瑒哟翱诘脑砣鐖D3所示,其將時間序列分成若干等寬的短序列,序列個數為N-s+1(為了避免序列末端數據剩余,取步長r=1)。

圖3 滑動窗原理 Fig. 3 Principle diagram of sliding window

1.3 基于多尺度滑動窗口的特征提取

地鐵車門所受的阻力主要來源于機械阻力、空氣阻力、密封條反彈力和密封條摩擦力四個方面,而這些阻力受到溫度和濕度等環境的影響。車門在開關門過程中的開始階段和結束階段阻力變化較大,從而使得這部分電機的轉速和轉矩信號變化明顯。圖4為車門開關門時電機的轉矩曲線,依據曲線的變化,將開關門過程大致分為3個階段,即啟動段、勻速段和減速段。開關門過程中啟動段和減速段蘊含著更多的車門狀態信息,對這兩個階段的數據挖掘有著更加重要的意義。在本文中,采用一種多尺度的滑動窗口的方法對開關門過程中的轉速、轉角和電流信號進行挖掘,在勻速階段采用較大的滑動窗口來將時間序列離散成字符,在啟動段和減速段采用比較小的滑動窗口進行離散。在盡可能多地挖掘時間序列數據特征的同時降低數據的維數。

通過符號化特征表示之后,可以識別車門電機參數形態特征,接著通過選取合適的距離度量的函數,可以找出數據之間的差別。歐氏距離[11]和最小距離[8]是目前最常見的兩種計算離散字符序列距離的方法。本文通過對最小距離方法的改進,計算時間序列的相似性。對于字符化后的時間序列AEsax={a1,a2,…,an}和BEsax={b1,b2,…,bn},計算公式如下:

(2)

其中:c為壓縮率,c=3n/w,n為滑動窗口長度,w為字符化后字符的個數;ai和bi來自字符表V={V1,V2,…,Vn},且

(3)

其中k為間斷點β={β1,β2,…,βn}之間的距離。

本文基于多尺度滑動窗口的符號化特征提取算法如下:

1) 初始化相關數據。字符種類數α=6,字符個數w=5,電機的轉角數據滑動窗口長度s=120;轉速和電流數據在啟動段和減速段滑動窗口長度s=30,勻速段滑動窗口長度s=120。

2)標準化車門正常狀態下的模板曲線,利用滑動窗口對其中的轉角、轉速和電流曲線依次分段,并通過對子序列均值、最大值和最小值的計算,將時間序列離散成ESAX字符串。

3)對亞健康和其對應正常數據的電機轉速、轉角和轉矩曲線預處理,包括問題數據的去除和數據長度的處理,然后采用步驟2)的方法進行離散,獲得離散后的字符串序列。

4)利用式(2)分別計算亞健康狀態下轉速、轉角和轉矩和模板之間的距離,即d(ESAX,ESAXT),將其組合獲得距離表示的特征集。

步驟2)中模板曲線對車門亞健康狀態的識別有著重要意義。由于地鐵門受其工作環境和人為因素的影響,使得其正常狀態之間的曲線也存在差異,比如溫度的降低或者氣壓的增加,車門開關過程中受到的阻力會變大,其電機電流將會增加。因此有必要考慮各亞健康對應的模板曲線,以減小實驗過程中門的調整所帶來的誤差,提高特征的可靠性。本文將車門調整之后正常運行一段時間數據的均值作為模板曲線。

圖4 電機電流曲線分段 Fig. 4 Curve segmentation of motor current

2 分層亞健康狀態識別模型

在實際的地鐵車門亞健康狀態診斷中,亞健康類別與亞健康的狀態并不存在完全的一一對應的關系,部分亞健康例如電機組松動、下擋銷干涉對電機的運行過程中所造成的阻力干擾比較小,對曲線的變化不如其他亞健康狀態明顯。將這些靠近正常狀態的亞健康狀態歸為Ⅱ類亞健康,而將遠離正常狀態的亞健康歸為Ⅰ類亞健康。本文采用一種由粗到精細分層分類的算法[12],對2類亞健康逐層進行區分,算法流程如圖5所示。

分層亞健康狀態識別算法含有2層識別過程:針對正常數據和亞健康數據,首先采用字符化處理對數據進行挖掘,獲得距離表示的特征集,接著使用主成分分析(Principal Component Analysis, PCA)進行降維,去除冗余信息并選擇貢獻率較高的特征用于第一層分類;第一層的訓練識別的過程主要針對Ⅰ類亞健康,這類亞健康往往擁有著很高的區分度。在第二層識別中,通過融合數據的基礎特征和距離表示的特征集,進一步區分出Ⅱ類亞健康數據和正常數據。本文的分層亞健康識別模型采用鑒別能力和推廣特性較好的支持向量機(Support Vector Machine, SVM)作為分類器。

圖5 分層亞健康識別模型 Fig. 5 Hierarchical sub-health state identification model

SVM的本質是采用核函數將特征非線性映射到高維空間,并在高維空間尋找最大間隔分類面。SVM常見的核函數包括線性核、高斯核和Sigmoid核,本文選用高斯核。SVM在進行分類時通常會引入松弛變量ξ和懲罰系數C,求解軟邊界二次規劃問題。其中C>0為一個自定義的懲罰因子,它決定對錯分樣本的懲罰程度,用來控制樣本偏差與機器泛化能力之間的平衡。C越大,懲罰就越大,暗示著對離群點的重視程度也越高,其取值與具體應用有關。當各亞健康特征的區分特性較好時,可以適當提高懲罰因子,以提高SVM分類能力。當亞健康特征區分度較差時,必須適當減小懲罰因子,保證分類器有著比較好的泛化性能。因此,在第一層分類模型中選擇較大的C來區別出Ⅰ類亞健康;第二層分類模型中,減小C的取值,以得到較為準確的識別結果。

在每一層的分類識別過程中,需要同時對多種亞健康狀態進行分類。單個支持向量機只能解決二分類問題,對于多分類問題,本文使用一對一的方法來實現多分類的功能。通過在各個亞健康之間構造決策函數,對t個類別共需構造t(t-1)/2個決策超平面。當對未知亞健康數據進行分類時,采用投票法,即得票最多的類為樣本所屬的類別。

傳統特征[13]通過總結時間序列中一個單一的值來描述時間序列,本文稱這些特征為基礎特征。這些特征計算簡單并且快速,實時性效果比較好,因而在實際中被廣泛地使用。常用的基礎特征類型包括三種:1)簡單的靜態值,例如時間序列均值、最大最小值、方差、標準差、斜率;2)頻率主導特征,例如時間序列的傅里葉變換以及穿越橫軸的次數和頻率;3)不同序列之間的相關特征,例如協方差、相關系數、歐氏距離等。本文采用的基礎特征主要包括兩部分:一是運動曲線靜態值特征,通過對三種曲線的分析,啟動段提取的特征包括最大轉速以及其對應的轉角、加速度和最大電流。勻速段提取的特征包括最大轉速、平均轉速和電流的有效值。減速段提取的特征包括最大加速度、最大電流、截止轉速、電流有效值和開關到位堵轉電流以及時間。二是曲線與模板曲線的相關特征,包括協方差和相關系數的計算。其中相關系數反映當前曲線與模板曲線的相似程度,協方差反映當前曲線與模板曲線的總體誤差。

3 實驗結果分析

3.1 數據來源及分析

實驗所依賴的地鐵門臺架為塞拉門結構,通過對其進行調整來模擬各種亞健康狀態。數據采集由內置的采集設備對電機開關門過程中的各項數據實時采集,并通過無線傳輸的方式(由外置監測設備轉發)將數據信息傳到數據中心服務器。電機數據采集設備(編碼器)由采樣電路、霍爾傳感器及無線傳輸模塊等組成,采集的數據量包括轉角、轉速和電流,數據采集流程如圖6所示。

圖6 數據采集流程 Fig. 6 Data acquisition process

主要針對實際運行過程中地鐵門可能發生的各種常見亞健康進行分析,通過對門結構進行微調整來模擬這些常見的亞健康狀態,本實驗包括的亞健康類型主要有V型異常、電機組件松動、對中尺寸變化小、對中尺寸變化大、上滑道外移、下擋銷橫向干涉、下擋銷縱向干涉和壓輪過壓這8種。由于每次的實驗過程中地鐵門的調整都會對門的正常狀態造成影響,所以每一組亞健康都對應一組正常數據。

在實驗中,選取8種亞健康狀態下的電機數據,每種亞健康狀態選取約40組數據(每組數據又分別包含開門數據和關門數據),每組正常數據也約40組,分別與亞健康數據一一對應。圖7為其中的一種亞健康狀態(對中尺寸變化小)以及其對應的正常狀態曲線。從圖7中可以看出,亞健康曲線與正常曲線相似程度比較高,其中轉角的正常曲線和亞健康曲線基本重合;轉速曲線和電流曲線在勻速段基本重合,在啟動段和減速段差異較為明顯。其他各組亞健康情況基本類似,所以選取了一組亞健康進行說明。

圖7 車門電機數據 Fig. 7 Data of door motor

3.2 距離表示的特征集參數分析

利用多尺度滑動窗口的特征提取算法對8種亞健康以及對應的正常數據進行分析,每組數據分別可以獲得600多個開門特征及關門特征。最終獲得的整體的特征已經達到1 200個之多,構成的特征維數很高,各組特征之間的相關性較大,且存在大量冗長信息(主要存在于開關門過程中的均速段)。直接將這些數據放入分類器進行模式識別,不僅給分類的過程造成較大的負擔,同時還影響分類的準確率。主成分分析(PCA)[14]通過映射的方法將特征變換為維數較少的新特征,在保留主要特征的同時極大地降低特征的維數。

針對開關門整體特征參數,通過求解協方差矩陣的特征值和對應的特征向量,確定各主要成分。圖8列出了前10個主成分的貢獻率隨主成分變化趨勢圖。

由圖8分析得到,從第5個主成分開始,貢獻率的變化趨勢已趨于平穩并逐漸接近于0。累計經過主成分變換后各特征值的貢獻率,前4個特征值累計貢獻率達到73%,相對的前20個累計貢獻率為90%。由此可見,前4個主成分可以很好地表示原本數據的大部分特征。為了進一步反映各主成分對車門亞健康狀態的敏感程度,圖9給出了前4個主成分在不同的亞健康狀態下的分布。由圖9可知,第一個特征擁有著最高的貢獻率,能準確地區分對中尺寸變化大、對中尺寸變化小、上滑道外移和正常這4個狀態,對其他幾種狀態存在著混疊現象;其他3個特征僅對部分的亞健康狀態有著很好的區分度,對其余狀態存在著明顯的混疊以及波動較大。結合這些特征的分布情況,可以發現區分度最好的亞健康狀態為對中尺寸變化小、對中尺寸變化大、壓輪過壓和正常狀態,其次為上滑道外移、V型異常和下擋銷縱向干涉,而電機組松動和下擋銷橫向干涉的區分度較差。綜上所述,將電機組松動和下擋銷橫向干涉歸為Ⅱ類亞健康,其他6種歸為Ⅰ類亞健康;在分層識別算法,首先對6種Ⅰ類亞健康進行區分,接著結合基礎特征,對正常數據和剩下的2種Ⅱ類亞健康狀態作進一步識別。

圖8 貢獻率隨特征變化趨勢 Fig. 8 Trend of contribution rate with principal components changing

圖9 前4個主成分在各亞健康狀態下的分布 Fig. 9 Distribution of first four principal components in different sub-health states

3.3 亞健康狀態的識別

3.3.1 基于距離表示特征的亞健康狀態識別

為了進一步驗證距離表示特征集的有效性,采用幾個常用分類器對其進行分類。選取3.1節中的8種正常數據和亞健康數據,其中一半用作訓練集,另一半用作測試集,得到總的訓練集數據330個,測試集數據326個,實驗分為3次,分別依據車門開門特征、關門特征以及整體特征進行識別。實驗過程中選擇降維后的前20個主成分進行訓練識別,分類器分別為1NN(k-Nearest Neighbor,其中k=1)[15]、人工神經網絡(Artificial Neural Network, ANN)[2]和支持向量機(Support Vector Machine, SVM)[12]。其中:SVM選擇高斯核,懲罰系數C為10;ANN使用多層感知器拓撲結構,節點輸出函數采用Sigmoid函數。

表2中的識別率為測試集正確識別的個數占測試集總個數的比例,基于距離的特征集通過對電機數據的字符化處理和序列中極值點的挖掘,其識別的準確率要高于傳統特征集。從各特征的識別準確率來看,開關門整體特征要比單個特征擁有更好的識別結果。各分類器識別的準確率相差不大,SVM識別的準確率比其他分類器略高。由于Ⅱ類亞健康與正常數據比較接近,錯誤識別主要發生在某些正常數據被錯分成Ⅱ類亞健康。

表2 不同特征的識別結果 %Tab. 2 Results of classification based on different feature %

3.3.2 基于融合特征的分層亞健康狀態識別

采用分層亞健康識別模型對同樣的數據進行識別。在第1層訓練中,將Ⅱ類亞健康和正常數據看作一類,依據基于距離特征集,對其他6種Ⅰ類亞健康進行識別;在第2層的識別中,融合基于距離的特征集和基礎特征集,對正常數據和Ⅱ類亞健康進一步判斷。各層分類器均采用SVM作為分類器,核函數使用高斯核,第1層中懲罰系數C為10,第2層中懲罰系數C為1。測試集數據每層的識別結果如表3所示。

表3 分層亞健康狀態識別結果 %Tab. 3 Results of hierarchical sub-health state identification %

表3中第1層的識別率表示第1層中測試集正確識別的個數占測試集總個數的比例;第2層對Ⅱ類亞健康進一步判斷,識別率表示兩層中正確識別的個數之和占測試集總數的比例。分析表3可知,與3.3.1節的識別結果相比,分層亞健康識別模型的開門、關門和整體特征的識別率都有所提高。在第1層中,將區分度較差的亞健康數據和正常數據歸為一類,先區分較為明顯的Ⅰ類亞健康,故三種特征的識別率都比較高,其中整體特征在這層的識別結果全部正確;第2層中,亞健康差別較小,錯誤明顯增多,識別率下降較大。最終本文方法的識別率可達到99%,測試樣本中僅有3個下檔銷橫向干涉被識別錯誤。由此可見,采用分層識別的模型能夠對亞健康數據做到更加準確的區分。

4 結語

本文提出了一種基于時間序列數據挖掘的軌道車輛門亞健康狀態識別的方法,通過對實測地鐵門電機數據的分析得到以下結論:

1)采用多尺度滑動窗口的方法并結合ESAX字符化算法對地鐵門電機數據進行字符化,有效地捕捉了時間序列中的極值點信息,同時能夠多層次地挖掘出車門運行過程中不同階段的信息。

2)將亞健康數據與模板曲線之間的距離作為特征集,減小實驗中車門調整所帶來的誤差;同時使用主成分分析對特征降維并分析了前4種主成分在不同的亞健康狀態下的分布,將實測的8種車門亞健康劃分為Ⅰ類和Ⅱ類。

3)通過車門8種亞健康狀態結果分析,表明分層亞健康識別模型比單個特征集有更好的分類效果,能夠準確識別各種亞健康狀態。

由于本文算法字符化后的字符串的長度較長,存在著效率較低的缺點,接下來的工作就是進一步提高算法實現的效率,以及對更多類型的亞健康狀態進行分析。

參考文獻(References)

[1] REN J B, LONG J, QIN Y, et al. Fault criticality evaluation of metro door based on WLSM and FWGM [C]// EITRT 2013: Proceedings of the 2013 International Conference on Electrical and Information Technologies for Rail Transportation, LNEE 288. Berlin: Springer, 2014, Ⅱ: 293-300.

[2] 朱興統,熊建斌.基于PCA和BP神經網絡的故障診斷仿真系統[J].自動化與儀器儀表,2015(12):47-48. (ZHU X T, XIONG J B. The fault diagnosis system based on PCA and BP neural network [J]. Automation & Instrumentation, 2015(12): 47-48.)

[3] 李海林,郭崇慧,楊麗彬. 基于時間序列數據挖掘的故障檢測方法[J].數據采集與處理,2016,31(4):782-790. (LI H L, GUO C H, YANG L B. Fault detection algorithm based on time series data mining [J]. Journal of Data Acquisition and Processing, 2016, 31(4): 782-790.)

[4] 胡為,高雷,傅莉.基于最優階次HMM的電機故障診斷方法研究[J].儀器儀表學報,2013,34(3):524-530. (HU W, GAO L, FU L. Research on motor fault detection method based on optimal order hidden Markov model [J]. Chinese Journal of Scientific Instrument, 2013, 34(3): 524-530.)

[5] GRUBB H J, WALDEN A T. Characterizing seismic time series using the discrete wavelet transform [J]. Geophysical Prospecting, 2010, 45(2): 183-205.

[6] 林意,王智博.基于一階濾波的時間序列分段線性表示方法[J].計算機工程,2016,42(9):151-157.(LIN Y, WANG Z B. Time series piecewise linear representation method based on first-order filtering [J]. Computer Engineering, 2016, 42(9): 151-157.)

[7] LIN J, KHADE R, LI Y. Rotation-invariant similarity in time series using bag-of-patterns representation [J]. Journal of Intelligent Information Systems, 2012, 39(2): 287-315.

[8] LIN J, KEOGH E, WEI L, et al. Experiencing SAX: a novel symbolic representation of time series [J]. Data Mining and Knowledge Discovery, 2007, 15(2): 107-144.

[9] LKHAGVAK B, SUZUKI Y, KAWAGOE K. Extended SAX: extension of symbolic aggregate approximation for financial time series data representation [EB/OL]. [2017- 03- 16]. http://www.ieice.org/iss/de/DEWS/DEWS2006/doc/4A-i8.pdf.

[10] 田再克,李洪儒,孫健,等.基于改進MF-DFA和SSM-FCM的液壓泵退化狀態識別方法[J].儀器儀表學報,2016,37(8):1851-1860. (TIAN Z K, LI H R, SUN J, et al. Degradation state identification method of hydraulic pump based on improved MF-DFA and SSM-FCM [J]. Chinese Journal of Scientific Instrument, 2016, 37(8): 1851-1860.)

[11] LAI C-P, CHUNG P-C, TSENG V S. A novel two-level clustering method for time series data analysis [J]. Expert Systems with Applications, 2010, 37(9): 6319-6326.

[12] 陳立江,毛峽,MITSURU I,等.基于Fisher準則與SVM的分層語音情感識別[J].模式識別與人工智能,2012,25(4):604-609. (CHEN L J, MAO X, MITSURU I, et al. Multi-level speech emotion recognition based on Fisher criterion and SVM [J]. Pattern Recognition and Artificial Intelligence, 2012, 25(4): 604-609.)

[13] SIIRTOLA P, KOSKIMAKI H, HUIKARI V, et al. Improving the classification accuracy of streaming data using SAX similarity features [J]. Pattern Recognition Letters, 2011, 32(13): 1659-1668.

[14] 律方成,金虎,王子建,等.基于主成分分析和多分類相關向量機的GIS局部放電模式識別[J].電工技術學報,2015,30(6):225-231. (LYU F C, JIN H, WANG Z J, et al. GIS partial discharge pattern recognition based on principal component analysis and milticlass relevance vector machine [J]. Transactions of China Electrotechnical Society, 2015, 30(6): 225-231.)

[15] 陳法法,湯寶平,蘇祖強,等.基于等距映射與加權KNN的旋轉機械故障診斷[J].儀器儀表學報,2013,34(1):215-220. (CHEN F F, TANG B P, SU Z Q, et al. Rotation machinery fault diagnosis based on isometric mapping and weightedKNN [J]. Chinese Journal of Scientific Instrument, 2013, 34(1): 215-220.)

XUEYu, born in 1992, M. S. candidate. His research interests include data mining, pattern recognition.

MEIXue, born in 1975, Ph. D., associate professor. Her research interests include image processing, pattern recognition.

ZHIYouran, born in 1984, Ph. D., associate professor. His research interests include fire science, pattern recognition.

XUZhixin, born in 1970, Ph. D. His research interests include fault diagnosis, data mining.

SHIXiang, born in 1956, M. S. His research interests include fault diagnosis, data mining.

猜你喜歡
亞健康特征
抓住特征巧觀察
新型冠狀病毒及其流行病學特征認識
如何表達“特征”
不忠誠的四個特征
當代陜西(2019年10期)2019-06-03 10:12:04
甩掉亞健康 我們在行動
少先隊活動(2018年5期)2018-12-29 12:12:58
跟蹤導練(一)
抓住特征巧觀察
亞健康科技工作者
豬群亞健康
《內經》治未病與亞健康的干預
主站蜘蛛池模板: 自拍偷拍欧美日韩| 伊人久久大香线蕉aⅴ色| 欧美日韩在线成人| 婷婷五月在线| 国产激爽大片在线播放| 婷婷激情五月网| 亚洲午夜福利精品无码不卡| 超级碰免费视频91| 人妻熟妇日韩AV在线播放| 国产激情无码一区二区三区免费| 国产国产人免费视频成18| 国产精品一老牛影视频| 亚洲日韩Av中文字幕无码| 国产精品女在线观看| 欧美日本在线播放| 色婷婷在线播放| 国产香蕉97碰碰视频VA碰碰看| 啪啪免费视频一区二区| 国内丰满少妇猛烈精品播| 熟妇人妻无乱码中文字幕真矢织江| 国产精品无码影视久久久久久久| 日韩欧美国产精品| 精品国产电影久久九九| 久久国产精品嫖妓| 久久久久人妻精品一区三寸蜜桃| 青青操国产视频| 91av国产在线| 国产女人水多毛片18| 欧美一区日韩一区中文字幕页| 波多野结衣一区二区三区AV| 人妻无码中文字幕第一区| 午夜电影在线观看国产1区| 一级毛片中文字幕| 女人爽到高潮免费视频大全| 99在线国产| aaa国产一级毛片| 日本91在线| 免费人欧美成又黄又爽的视频| 激情五月婷婷综合网| 无套av在线| 婷婷亚洲天堂| 真实国产精品vr专区| 国产精品亚洲欧美日韩久久| 日韩专区欧美| 亚洲男人在线| 97se亚洲综合在线天天| 国产主播一区二区三区| 在线免费观看AV| 午夜日本永久乱码免费播放片| 国产成人一级| 亚洲精品成人7777在线观看| 97人人模人人爽人人喊小说| 亚洲乱码精品久久久久..| 亚洲制服中文字幕一区二区| 亚洲国产欧美中日韩成人综合视频| 色综合天天娱乐综合网| 欧美激情第一区| 久久国产热| 免费精品一区二区h| 亚洲天堂网站在线| 国产制服丝袜91在线| 久久夜色精品| 国产综合在线观看视频| 成人在线欧美| 五月天香蕉视频国产亚| 久久久久久高潮白浆| 波多野结衣一二三| 黄网站欧美内射| 日韩在线观看网站| 亚洲国产日韩一区| 国产精品第页| 欧美精品成人一区二区在线观看| 丁香婷婷激情网| 国产成人免费| 国产精品久久久精品三级| 国产1区2区在线观看| 婷婷成人综合| 日韩精品一区二区三区swag| 免费xxxxx在线观看网站| 国产97色在线| 在线另类稀缺国产呦| 亚洲日韩精品伊甸|