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

基于SVM的下肢運動表面肌電信號的特征提取與辨識分析

2015-11-07 03:12:06趙慧龍張向剛秦開宇電子科技大學航空航天學院人機耦合技術實驗室成都611731
載人航天 2015年2期
關鍵詞:動作

張 羿,趙慧龍,張向剛,秦開宇(電子科技大學航空航天學院人機耦合技術實驗室,成都611731)

基于SVM的下肢運動表面肌電信號的特征提取與辨識分析

張羿,趙慧龍,張向剛*,秦開宇
(電子科技大學航空航天學院人機耦合技術實驗室,成都611731)

為了對蹲、站、行走支撐與行走擺動四類人體典型運動行為進行有效的分類辨識,三位健康且未經過專業訓練的受試者被邀請參加運動實驗,實驗對人體下肢股內側肌的表面肌電信號(Electromyography,EMG)進行實時采集和記錄。通過時域、頻域、時頻域方法,對特征值進行提取,發現下肢股內側肌的表面肌電信號在下蹲、起立、行走支撐期和行走擺動期四類運動方式下的動態特征具有明顯的差異性。基于上述結論,介紹了一種基于支持向量機(Support Vector Machine,SVM)的誤識別樣本二次分類方法,對上述典型運動類型進行了辨識分析。與傳統單次樣本識別結果相比,本文所發展的基于SVM特征值分析的誤識別樣本二次分類方法能較好地提高識別效果,樣本辨識結果顯示時域和時頻域的識別效果最好,時頻域方法抗外界干擾能力較好。

下肢;表面肌電信號;特征提取;支持向量機;時頻域

1 引言

近些年,宇航員地外生存的健康保障課題已被認定為最具挑戰的最高優先級技術[1]。從近地走向深空,是載人航天的內在驅動力。我國載人航天在完成“三步走”戰略目標之后,探索的腳步將逐步邁向地球高軌、月球、近地小行星、火星等目的地[2]。然而,載人航天實踐表明,中長期空間站駐留,將會導致宇航員心血管功能障礙、骨質丟失、肌肉萎縮、免疫功能下降、空間運動病等生理效應[3]。人體肌肉神經系統的反應特征能有效地評估上述生理效應的變化規律,安全運動鍛煉能提高人體健康水平和身體素質。因此,提供一種安全、可靠且有效的個性化運動方案可以幫助宇航員提高心血管系統健康,預防功能退化性疾病,如骨質疏松、肌肉萎縮等[4]。對表面肌電信號(Surface Electromyogram,SEMG)的研究能夠幫助分析肌肉的勞損程度,預測肢體運動信息,評估宇航員在地外生活中日常運動對肌肉活動能力的影響。表面肌電信號是通過電極引導、放大而從肌肉表面記錄下來的神經肌肉系統活動時的非平穩一維時間序列生物電信號,它能夠反映肌肉強度和人的運動特征[5]。人的運動意圖通常通過神經興奮刺激肌肉細胞的收縮與舒張實現[1,2,6],由于不同的肢體運動中肌肉收縮的模式不同,導致相應的表面肌電信號的特征也具有差異,在人的意圖控制下產生的表面肌電信號能較好地反映出肢體運動或動作特征,評估人體運動意圖[1,3]。

傳統的基于表面肌電信號的人體運動特征分析方法往往僅針對時域或頻域信息對人體不同的運動狀態進行特征分析與提取[7-9]。這些方法都將肌電信號視為平穩或者短時平穩信號進行處理,然而實際上肌電信號應為非平穩信號[1,5,6]。為了有效地分析SEMG的非平穩性,小波變換被應用到了SEMG的特征提取[10,11]。它是將時域和頻域分析相結合,對表面肌電信號的時間和頻率所包含的信息研究。由于運動刺激產生的表面肌電信號的非平穩隨機性、時變和易干擾等特性,小波變換也很難對由運動產生的SEMG特征提取進行一般性概括[6,12]。

本文首先對人體下肢下蹲、起立、行走擺動和行走支撐四個典型基本動作的肌電(EMG)信號進行了實驗調查。為了有效降低跨個體差異性(Inter subject variability),每位受試者對每個動作試驗進行3次試驗,同時每次試驗也對該動作重復20次。而考慮到同個體差異性(Intra subject variability),三位健康且未經過專業運動訓練的受試者(年齡范圍在20~30歲)被邀請參加試驗,在試驗前1小時僅被允許攝入適量的食物。三位受試者分別進行了下蹲、起立、行走擺動和行走支撐四類動作的實驗,單通道肌電采集系統(如圖1所示)被采用來獲取運動時股內側肌的表面肌電信息。其次,為了揭示上述四類基本動作的可分類性程度,時域、頻域和時-頻域三類特征提取方法被采用,支持向量機(Support Vector Machine,SVM)分別對上述四個動作的時域、頻域和時-頻域特征進行分類驗證與對比分析,同時也提出了一種基于支持向量機(SVM)特征值分析的誤識別樣本二次分類方法(Double-round Feature Recognition and Filtering,DFRF),較好地改善了特征的可辨識率(包括時域、頻域和時-頻域特征)。

2 實驗介紹

本實驗在電子科技大學航空航天學院人機耦合實驗室完成。為觀察不同受試者在下肢運動時的股內側肌活動反應,三名受試者被邀請,每位受試者設定代號為A~C,受試者基本資料如表1所示。

表1 受試者基本資料Table 1 Data of the sub jects

實驗所采用的SEMG傳感器包括信號的檢測、信號放大、信號濾波、A/D轉換,通過串口輸入到PC機上。如圖1所示:

測試對象為正常人群,且該被測試對象身體健康、體態正常、身體各關節運動狀況良好,無任何關節運動障礙疾病。同時,為了減小同一對象和主體間的差異性,在進行SEMG采集前24 h內,被測試者須避免任何形式的劇烈運動,并保證下肢肌肉的舒適感。實驗要采集下肢的下蹲、起立、行走擺動期和行走支撐期四種動作的數據,如圖2和圖3所示,其中圖2中的A—B為下蹲,B—C為起立,圖3中的A—C為擺動期、C—E為支撐期。由于四種動作與膝關節密切相關,控制膝蓋彎曲等動作,因此本文針對股內側肌的肌肉運動反應特征進行采集與記錄[13]。

圖1 肌電信號傳感器Fig.1 Sensor of the SEMG signal

與其他生物電信號一樣,SEMG本身非常微弱,信噪比較低,易受周圍噪聲信號的干擾,盡管在電路設計中已充分考慮,但在實際SEMG采集過程中仍需采用一定的措施加以抑制:信號采集之前,用濕巾擦洗貼片位置并晾干;在粘貼電極片時,電極片的正負兩端沿肌纖維的長軸方向縱向放置,并將電極片粘貼在肌腹處,用手輕輕按壓,保證電極片與皮膚接觸良好并牢牢固定;表面肌電電極片均為一次性粘貼,不能重復使用[13]。

數據采集過程中,受試者做有規律的運動,動作周期大約為3 s,每次采集大約1 min,對每個受試者的每個動作每次采集20次,每次采集完之后測試者要休息10min,然后再接著采集下一組,每個動作每個人采集3組,然后從中挑選出采集效果比較好的數據進行分析。

原始數據進行10階中位值濾波以去除采集過程中由于通信錯誤等原因產生的奇異點,然后進行20~500 Hz帶通濾波以去除低頻和高頻的噪音并加絕對值,最后用50階滑動平均濾波器反復濾波10次以得到平滑的時域信號,其實驗結果如圖4所示(其中灰色區域為所有樣本的范圍,即所有樣本都在灰色區域中,黑色加粗線為所有樣本的平均值)。

圖4 實驗數據結果(左上:行走擺動,右上:行走支撐,左下:下蹲,右下:起立)Fig.4 Experiment results(upper left:sw ing phase,upper right:support phase,lower left:stand-squat transition,lower right:squatstand transition)

頻域信號如圖5所示,上述動作被選取任意20個樣本,通過20~500 Hz的帶通濾波器之后,計算出其功率譜密度,然后對20個樣本求平均。

圖5 頻域肌電信號圖Fig.5 Frequency-domain feature of EMG signals

根據上述實驗結果,我們發現上述四類動作在時域中均存在較大的差別性(如圖4所示),且下蹲、起立和步態這三個動作之間的幅值有比較大的差別,而步態中的擺動期和支撐期在面積上有一定差別;上述動作的功率譜在幅度上也存在差異,下蹲、起立和擺動支撐之間有較大差異,但是下蹲和起立之間差異較小,且支撐和擺動之間差異也較小(如圖5所示)。

表2 時域特征均值及標準差Table 2 Average and standard deviation of time-domain features

3 特征提取及分析

本文針對實驗觀察結果分別進行了時域、頻域和時頻域三種方法的比較分析。

時域分析方法采用了平均絕對值肌電值(MAV)、均方根(RMS)、絕對值積分(iEMG)和過零點數(ZC)四種提取方法[7]。表2介紹了上述方法對動作特征的統計結果。由表2可知,在時域中下蹲、起立與支撐和擺動三者之間的MAV、RMS和iEMG差異比較大,說明這三個特征在這三者之間的類間距離很大,具有相當的可分性,而MAV的標準差小說明聚類性好;支撐擺動之間的iEMG和ZC差異比較大,且標準差小,說明它們之間的類間距很大且聚類性很好,具有很好的可分性。

對于頻域分析方法,本文采用了常用的平均功率頻率(MNF)、中值頻率(MDF),作為表面肌電信號的頻域特征值指標[14]。表3揭示了上述動作樣本的MNF和MDF的標準值與偏差。由表3可知,在頻域中下蹲、起立、支撐和擺動四者之間的MDF差異比較大,說明這個特征在這四者之間的類間距離很大,具有可分類性,而擺動和支撐的MDF標準差小說明擺動和支撐的聚類性也很好。

表3 頻域特征均值及標準差Table 3 Average and standard deviation of frequencydomain features

最后,本文采用了基于小波變換的時頻域特征提取方法,對第2節所獲得的實驗數據進行了分析。圖6展示了本文所采用的5層小波分解方式。給定一個長度為N的信號S,離散小波分解(DWT)最多可以把信號分解成log2N個頻率級。第一步分解開始于信號S,分解后分解系數由兩部分組成:低頻系數向量D1和高頻系數向量H1,向量D1是由信號S與低通分解濾波器經過卷積運算得到的,向量H1是由信號S與高通分解濾波器經過卷積運算得到的。在下一步分解中,用同樣的方法把低頻系數D1分成兩部分,即把上面的S用D1代替,分解后返回尺度2的低頻系數D2和高頻系數H2,再往下分解依此類推。最后得到一個低頻系數列D5和五個高頻系數列H5,H4,H3,H2,H1。基于上述分解,本文選用db3小波函數進行實驗,以下是對上面6個系數列提取奇異值特征的結果。表4記錄了四種動作樣本的以上特征的平均值和標準差。通過觀察上表可知下蹲、起立、擺動和支撐的H4的差異比較大且標準差相對不大,所以H4的奇異值可分性比較好。

表4 小波變換后各層系數的奇異值的均值及標準差Table 4 The average and standard deviation for singular values of hierarchical coefficien ts after w avelet transform

4 基于SVM的分類

設給定的訓練集{(x1,y1),(x2,y2),…,(xm,ym)},其中xi∈Rn,y∈{-1,1},即xi為一個樣本矢量,y表示分類的標記,也稱為屬性矢量。如果訓練集的所有向量均能由某超平面正確劃分,并且距離最近的異類向量之間的距離最大,則該超平面為最優平面,其中距離超平面最近的異類向量被稱為支持向量。該超平面記為式(1)[15]:

其中,w是為最優分類面的權系數向量,b為常數,X為輸入樣本向量。我們定義一個樣本點到超平面的間隔為:

現在對w和b歸一化,即用w/‖w‖和b/‖w‖分別代替w和b,則間隔可以寫作式(2):

其中,σg表示幾何間隔。幾何間隔與樣本誤分次數之間的關系可表示為式(3):

其中R=max‖xi‖,σ是樣本集合到分類面的間隔。由以上可得知構造最優化超平面的問題轉化為求min‖w‖,但是必須得加約束,否則幾何間隔有可能為無窮大。因此,約束條件可被表示為:

引入Lagrange乘子,得到函數:

其中ai為非負Lagrange乘子。

根據最優化理論,存在唯一的最優解,最優解在該函數的極點處,如式(4)對w和b求偏導數可得式(4):

上述是線性可分的情況,如果所給的樣本是線性不可分的,線性分類器的求解程序會無限循環下去。這必然使得它的應用范圍大大縮小,但是它的很多優點我們不愿放棄,這就引入了核函數,它將輸入空間映射到更高維的空間中。

常用的核函數包括線性核,多項式核和徑向基核,分別為式(5)~(7)所示:

其中γ,r和σ均為核參數。

5 結果

5.1時域結果分析

圖7為四種動作在時域中各自的特征值分布情況:其中藍色—下蹲,紅色—起立,綠色—擺動,黑色—支撐(經過統一尺度的變換)。

圖7時域特征值分布圖Fig.7 Distribution of time-domain features

圖8為用SVM進行識別后的結果,紅點為提前標定的類別即為實際的類別,而黑點為用SVM識別的結果,每個樣本的紅點和黑點互相重疊的說明識別的結果是正確的,不重疊則說明是不正確的。

圖8 時域識別結果圖Fig.8 Recognition resu lts of tim e-domain features

由圖8可知,下蹲動作中有兩個樣本識別錯誤,被錯劃分成了起立動作,而擺動樣本有一個被錯劃分成了支撐動作。四種動作總的識別率為94.23%。此結果驗證了下蹲、起立、行走擺動和行走支撐在時域中的可分類性。

5.2頻域結果分析

圖9為四種動作在頻域中各自的特征值分布情況:其中藍色—下蹲,紅色—起立,綠色—擺動,黑色—支撐,橫軸表示特征值類型,縱軸表示特征值的大小(經過統一尺度的變換)。

圖9 頻域特征值分布圖Fig.9 Distribution of frequency-domain features

由圖9和表3可以觀察到,其頻域特征不具有可分類性。圖10為用SVM進行識別后的結果,其中紅點為提前標定的類別即為實際的類別,而黑點為用SVM識別的結果,每個樣本的紅點和黑點互相重疊的說明識別的結果是正確的,不重疊則說明是錯誤的。

由圖10可知,起立動作中有一個樣本識別錯誤,被錯劃分成了下蹲動作,而擺動樣本有二個被分別錯劃分成了支撐動作和起立動作,但是支撐中有很多被錯劃分成了擺動。四種動作總的識別率為78.80%。結果驗證了下蹲、起立、行走擺動和行走支撐在時域中不具有可分類性。

圖10頻域識別結果圖Fig.10 Recognition results of frequencydomain features

圖11為頻域特征值大小與頻域特征的關系。

圖11 頻域特征值大小Fig.11 The am p litudes of frequency-domain features

其揭示了擺動和支撐之間的識別率很低,然而支撐(綠色)和擺動(黑色)的斜率差較大,因此可通過斜率區分支撐與擺動狀態。

圖12為加上斜率(slope)的特征值分布結果:

圖12斜坡特征值分布圖Fig.12 Distribution of slope features

圖13為加上斜率后的識別結果。

與圖10結果相比,二次誤樣本的頻域識別結果(如圖10和圖13所示)有了很大改善,且僅一個樣本被錯話分成了擺動,另有兩個被錯劃分成了起立。四種動作總的識別率為88.4%,相比加策略之前有了很大的提高。

圖13 二次誤樣本的頻域識別結果圖Fig.13 The recognition results of frequency-domain features

5.3時頻域結果分析

圖14為四種動作在頻域中各自的特征值分布情況:(經過統一尺度的變換)。

圖14 時頻域特征值分布圖Fig.14 Distribution of time-frequency domain features

由圖14時域特征值分布可以更直觀的看到,以上時域特征具有比較好的可分性,它更加直觀地體現了我們從表4中得到的結論。

圖15為用SVM進行識別后的結果,其中紅色的點為提前標定的類別即為實際的類別,而藍色的點為用SVM識別的結果,每個樣本的紅色點和黑色點互相重疊的說明識別的結果是正確的,不重疊則說明是不正確的。

由圖15可知,擺動和支撐之間的識別率比較低,而下蹲和起立的識別率還是比較高的,四種動作總的識別率為88.8%。

圖16揭示了時頻域特征值大小與頻域特征的關系。

圖16結果分析可知,擺動和支撐之間的識別率很低,但是仔細觀察圖16中的綠線和黑線可以發現,H3和H4之間的斜率相差比較大,所以我們可以用圖16中H3和H4的斜率來識別支撐和擺動這兩個動作。圖17為加上斜率(slope)的特征值分布結果:

圖17加斜率后特征值分布圖Fig.17 Feature distribution w ith the slope

圖18為加上斜率后的識別結果。通過比較第一次識別結果(如圖15所示)和本文提出的算法識別結果(如圖18所示),可驗證二次誤樣本的加斜率辨識算法能有效提高四種動作的時頻域特征的識別率,在測試樣本中的識別率達96.1%。

圖18 加斜率后的識別結果Fig.18 The recognition resultsw ith the slope

6 結論

本文利用時域、頻域和時頻域三種特征提取方法,對人體下肢EMG反應特征進行了提取與分析,比較了三種方法對蹲、站、行走支撐與行走擺動四類人體典型運動行為可分類性,也提出并驗證了一種基于支持向量機(SVM)特征值分析的誤識別樣本二次分類方法(DFRF)。與第一輪誤識別結果相比,該方法能有明顯的改善。驗證結果也表明,時域和時頻域方法在離線分析中正確率最高,為了考慮外界干擾等因素時頻域方法具有較好的可靠性。

[1]強靜,宋堯.國外航天生物醫學發展戰略[J].載人航天信息,2013,3:1-5. Qiang Jing,Song Rao.Foreign space biomedical development strategy[J].Manned Space Information,,2013,3:1-5.(in Chinese)

[2]白延強,吳大蔚.長期載人航天中的醫學挑戰與對策[J].航天醫學與醫學工程,2008,21(3):210-214. Bai Yanqiang,Wu Dawei.Long term manned spaceflight medical challenges and Countermeasures[J].Space Medicine &Medical Engineering,2008,21(3):210-214.(in Chinese)

[3]果琳麗,王平,朱恩涌.載人月球基地工程[M].北京:中國宇航出版社.2013:502-510. Guo Linli,Wang Ping,Zhu Enyong.Manned lunar base engineering[M].Beijing:China Aerospace Press,2013:502-510.(in Chinese)

[4]Xu Z,Tian Y,Li Z,et al.Proceedings of the International Conference on Information Engineering and Applications(IEA)2012[M].London:Springer,2013:495-503.

[5]Pal PR,Goswami G,Mohanty N P.Computer Networks and Intelligent Computing[M].Springer Berlin Heidelberg,2011:502-510.

[6]Ageberg E,Bennell K L,Hunt M A,et al.Validity and inter-rater reliability ofmedio-lateral kneemotion observed during a single-limbmini squat[J].BMCMusculoskeletal Disorders,2010,11(1):265.

[7]楊鵬,陳玲玲,郭欣,等.基于支持向量機的表面肌電信號識別[C]//2006中國控制與決策學術年會論文集. 2006. Yang Peng,Chen Lingling,Guo Xin,et al.Support Vector Machine based surface EMG signal recognition[C].Shenyang:China control and decision Symposium,2006.(in Chinese)

[8]Rekhi N S,Arora A S,Singh S,et al.Multi-Class SVM Classification of Surface EMG Signal for Upper Limb Function[C]//Bioinformatics and Biomedical Engineering,2009. ICBBE 2009.3rd International Conference on.IEEE,2009:1-4.

[9]杜樹新,吳鐵軍.模式識別中的支持向量機方法[J].浙江大學學報:工學版,2003,37(5):521-527. Du Shuxin,Wu Tiejun.Support vector machine method in pattern recognition[J].Hangzhou:Journal of Zhejiang University:Engineering Science Edition,2003,37(5):521-527.(in Chinese)

[10]Xu Z,Tian Y,Li Z,et al.Proceedings of the International Conference on Information Engineering and Applications(IEA)2012[M].Springer London,2013:495-503.

[11]蔡立羽,王志中,張海虹.小波變換在表面肌電信號分類中的應用[J].生物醫學工程學雜志,2000,17(3):281-284. Cai Liyu,Wang Zhizhong,Zhang Haihong.Application of wavelet transform in the classification of sEMG[J],2000,25(1):105-108.(in Chinese)

[12]胡曉,王志中,任小梅.小波變換和非線性分析在表面肌電信號中的應用及進展[J].北京生物醫學工程,2006,25(1):105-108. Hu Xiao,Wang Zhizhong,Ren Xiaomei.Application and progress of sEMG signal by wavelet transformand nonlinear analysis[J].Beijing:Beijing Biomedical Engineering,2006,25(1):105-108.(in Chinese)

[13]賈雪琴,王旭,李景宏,等.基于小波變換和K-L展開的單通道表面肌電信號識別[J].東北大學學報:自然科學版,2006,27(8):859-862. Jia Xueqin,Wang Xu,Li Jinghong,etal.Single channel surface EMG signal recognition based on wavelet transform and K-L expansion[J].Shenyang:Journal of Northeastern University:Natural Science Edition,2006,27(8):859-862.(in Chinese)

[14]邱青菊.表面肌電信號的特征提取與模式分類研究[D].上海:上海交通大學,2009. Qiu Qingjun.Research on feature extraction and pattern classificationof surface EMG signal[D].Shanghai:Shanghai Jiao Tong University,2009.(in Chinese)

[15]王士允.基于表面肌電信號的膝關節康復機器人控制技術研究[D].南京:南京理工大學,2013. Wang Shiyun.Research on control technology of knee rehabilitation robot based on surface EMG signal[D].Nanjing:Nanjing University of Science and Technology,2013.(in Chinese)

Support Vector M achine based Feature Extraction and
Recognition of EMG Signals during Lower Extrem ity M ovement

ZHANG Yi,ZHAO Huilong,ZHANG Xianggang,QIN Kaiyu
(The Human-machine Coupling Technology Laboratory,School of Aeronautics and Astronautics,University of Electronic Science and Technology of China,Chengdu 611731,China)

To effectively identify the dynamic behaviors of squat,stand,support,and swing phases ofwalking from EMG signals in the lower extremity,three healthy and untrained subjectswere invited to participate in the exercise experiment.During experiments,the EMG signals of vastusmedialis in terms of squat,stand,support,and swing phases ofwalking were investigated and the features of the time-domain,frequency-domain,and time-frequency domain were extracted respectively.The experimental data showed that those dynamic characteristics were distinctively different.Based on the observed experimental results,the authors proposed a Double-round Feature Recognition and Filtering(DFRF)approach for the error recognition samples by using Support Vector Machine-based Regression(SVM).Compared to the previous recognitionmethods,it is advantageous that the error recognition rate can be well improved since those error recognition samples evaluated by the firstround feature recognition are further filtered following the second-round feature recognition.The proposedmethod also showed a promising result,indicating that the EMG featuresextracted by time-domain and time-frequency domain have better recognition performance than those by frequency-domain.

lower extremity;EMG;feature extraction;support vectormachine;time-frequency domain

V445,TP391

A

1674-5825(2015)02-0179-08

2014-12-01;

2015-03-13

中央高校基本科研業務費項目(ZYGX2013J123)

張羿(1984-),男,博士,講師,研究方向為生理信號處理、生理參數的監護、建模與控制。E-mail:yi.zhang@uestc.edu.cn

張向剛(1973-),男,博士,副教授,研究方向為人機系統與機器人。E-mail:csxgzhang@uestc.edu.cn

猜你喜歡
動作
動作不可少(下)
巧借動作寫友愛
下一個動作
動作描寫要具體
畫動作
讓動作“活”起來
動作描寫不可少
非同一般的吃飯動作
動作喜劇電影周
電影故事(2015年30期)2015-02-27 09:03:12
神奇的手
主站蜘蛛池模板: 色天天综合| 国产理论一区| 国产精品刺激对白在线| 中文字幕无码中文字幕有码在线| 亚洲AV成人一区国产精品| 国产毛片不卡| 91亚洲免费视频| 欧美啪啪一区| 欧美一级一级做性视频| 女同久久精品国产99国| 久久特级毛片| 成人免费午夜视频| 九九久久精品免费观看| 无码国产伊人| 精品国产成人a在线观看| 成人精品亚洲| 国产经典三级在线| 免费一级毛片不卡在线播放| 亚洲天堂久久新| 国产亚洲欧美日韩在线一区| 老司国产精品视频91| 国产成人麻豆精品| 精品无码一区二区在线观看| 成人国产精品网站在线看| 岛国精品一区免费视频在线观看| 99热这里只有免费国产精品| 777国产精品永久免费观看| 99re视频在线| 免费无码在线观看| 亚洲精品无码日韩国产不卡| 台湾AV国片精品女同性| 欧美国产视频| 久久semm亚洲国产| 国产麻豆va精品视频| 免费亚洲成人| 视频二区国产精品职场同事| 伊人久久大香线蕉综合影视| 国产网站免费观看| yy6080理论大片一级久久| 色一情一乱一伦一区二区三区小说| 中文精品久久久久国产网址| 亚洲AⅤ波多系列中文字幕 | 999福利激情视频| 婷婷成人综合| 成人国内精品久久久久影院| 人禽伦免费交视频网页播放| 精品成人免费自拍视频| 日韩a级片视频| 天天色综合4| 在线日本国产成人免费的| 狠狠色婷婷丁香综合久久韩国| 欧美日韩在线第一页| 狠狠综合久久| 中文字幕在线视频免费| 亚洲系列无码专区偷窥无码| 在线毛片网站| 在线中文字幕日韩| 国产麻豆福利av在线播放| 另类重口100页在线播放| 国产精品污污在线观看网站| 亚欧成人无码AV在线播放| 亚洲va精品中文字幕| 国产成人精品男人的天堂下载| 精品超清无码视频在线观看| 国产精品极品美女自在线网站| 欧美成人免费一区在线播放| 日韩第一页在线| 日韩在线播放欧美字幕| 中文成人在线视频| 欧美日韩国产成人在线观看| 国产欧美日韩另类| 亚洲成人精品在线| 九九九国产| 18禁黄无遮挡网站| 99久久国产精品无码| 国产一二三区视频| 一本大道香蕉高清久久| 国产精品三级av及在线观看| 国产日本视频91| 波多野结衣中文字幕久久| 国产18在线| 日韩在线欧美在线|