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

基于超聲射頻信號的支持向量機雙參量B線識別方法*

2023-09-20 06:49:36張皓宇馬泉龍鐘徽
應用聲學 2023年5期
關鍵詞:分類信號效果

張皓宇 馬泉龍 張 蕾 鐘徽

(西安交通大學生命學院 生物醫學信息工程教育部重點實驗室 西安 710049)

0 引言

隨著現代超聲醫學的不斷發展進步,曾經被視作超聲盲區的肺臟超聲檢查也由于其無輻射、可床旁診斷等優點成為臨床上檢測肺臟疾病、監測重癥病人生理指標的重要手段。肺超聲中有一種被稱為B 線的重要超聲偽影,多用于檢測、診斷肺部疾病,在超聲影像上表現為從胸膜線產生并延伸至屏幕底部的離散、激光樣垂直高回聲偽影[1],并與肺滑動同步移動。1998 年,Lichtenstein等[2]提出B 線定義,并用肺臟超聲圖像中的B 線來區分肺水腫和慢性阻塞性肺病。2004 年,Jambrik等[3]通過分析121名患者的胸部X光片與肺超聲之間的相關性,證實了B線作為檢測血管外肺水指標的可行性。2010 年,Jambrik等[4]通過肺損傷動物模型實驗驗證了B 線數量與肺干濕比之間的顯著相關性,提出了一種簡單、半定量、無創的評估肺超聲B 線的肺積水指數。肺超聲影像中B 線數量的增加與肺水腫、肺纖維化、肺炎等肺病變之間顯著的相關性已經成為臨床上對患者進行床旁無創、快速診斷的重要依據,醫生可以按照肺超聲檢查國際推薦標準[5],快速評估患者的病情,并做出相應的診斷。但在相同的評估標準之下,不同的醫生對同一患者的評估、甚至同一醫生在不同時間對同一患者的評估都可能會產生一定偏差,這與評估者的臨床經驗、評估時的生理狀態等因素有關。因此,能夠輔助醫生對肺超聲進行評估的算法可以大大提高醫生診斷的穩定性和準確性。

2016年,Moshavegh等[6]提出了一種自動檢測和表征B 線的方法,該方法首先檢測圖像中分割胸膜線,再向下對B線進行檢測。Corradi等[7]通過比較肺超聲圖像灰度值變化的平均值和標準差,提出了一種檢測機械通氣心臟手術患者肺水腫的計算機輔助定量方法。2020年,Cristiana等[8]提出了一種基于深度卷積神經網絡的自動評估肺超聲B 線的深度學習算法,并將評估結果與超聲專家評估結果進行比較,表現出較好的敏感性和特異性。

現有的B 線識別方法大多基于超聲圖像,但超聲圖像是由原始超聲回波射頻(Radio frequency,RF)信號經過預處理、平滑、對比度增強等處理得到的便于醫生肉眼觀察的圖像,在由RF 信號轉變為圖像的過程中,真實信號丟失或弱化了一部分的組織信息以換取超聲圖像更有利于觀察的性能,故本文提出的利用原始RF 信號對肺超聲結果進行評估可以提高信息的利用率,避免不同超聲儀器參數設置的影響,可以得到普適性更高的評估算法。與基于卷積神經網絡的自動評估B 線深度學習算法相比,本文采用的非線性支持向量機(Support vector machine,SVM)復雜度更小,提取特征參量的運算量也更小,有望應用于臨床超聲檢查中,實現實時的B線識別功能。

1 理論和方法

1.1 特征參數選取

本文選用的超聲RF 數據每幀包含256 條掃描線,每條掃描線有2336個數據點,采樣率為40 MHz。由于B 線的特征是垂直于成像探頭的激光樣線,與掃描線方向相同。故將每一條掃描線上的全部RF數據點看作一個識別對象,分別對每個識別對象提取信息熵、排列熵、峰度、偏度和能量5 種特征參量,并評估這些參量B線識別的效果。

(1) 信息熵

信息熵是反映信號混亂程度的參數[9],信號越混亂,信息熵越大。肺中富含氣體,超聲波在氣相表面發生全反射,無法探測氣相后的組織信息,在超聲圖像上產生A 線和黑色區域,黑色區域在RF 信號中表現為0 左右的低幅值數據,由于混響效應[10]而產生的B 線在圖像上表現為明亮的彗星尾征,在RF 信號中表現為幅值較高的回波數據。本研究認為與正常肺氣相區域相比,B 線區域RF 數據的信息熵在宏觀上應具有更大的值,選取信息熵作為特征參量有助于B線識別算法的實現。信息熵的計算公式如下:

其中,H為信息熵,X為RF 信號時間序列幅值,f(Xi)是信號幅值為Xi的概率密度,單位為bit。

(2) 排列熵

排列熵是一種檢測動力學突變和時間序列隨機性的方法,可以反映出時間序列微小的細節變化,能夠定量評估信號序列中含有的隨機噪聲,越規律的時間序列,排列熵越小[11]。本研究認為胸膜線下,病變導致的血管外肺水增多使得肺水腫區域富含氣體與液體,而超聲波在該區域的多次反射形成了B 線,故B 線區域與其他區域相比,其RF 數據在微觀上應具有一定的規律性與重復模式,選擇排列熵作為特征參量有助于區分識別B 線。排列熵的計算步驟如下:

對一個含有2336 個數據點的元素X進行空間重構,得到矩陣Y為

式(2)中,m為嵌入維數,t為延遲時間,k=N -(m-1)t,矩陣Y的每一行都是一個重構分量,共有k個重構分量。

(a) 將每一個重構分量按照升序重新排列,得到向量中各元素位置的列索引構成一組符號序列:

(b)M維相空間映射不同的符號序列總共有m!種,計算每一種符號序列出現的次數除以總數m!,作為該符號序列出現的概率。得到概率{P1,P2,···,PK},按公式(4) 計算得到元素X的排列熵:

式(4)中,Hpe為元素X的排列熵,k為重構分量總數Pj為第j個符號序列出現的概率。

峰度與偏度是研究信號包絡特征以及分級定征的常用統計量,本研究提取各個元素的峰度與偏度并進一步評估其在B線識別方面的效果。

(3) 峰度

峰度又稱峰態系數,表征概率密度分布曲線在平均值處峰值高低的特征數。峰度反映了RF 信號包絡峰部的尖度,是和正態分布相比較而言的統計量,如果峰度大于3,峰的形狀比較尖,比正態分布更陡峭,反之亦然。RF數據元素X包絡的峰度按公式(5)計算:

式(5)中,Kurt(X)為RF數據元素X包絡的峰度,μ為RF數據元素X包絡的均值。

(4) 偏度

偏度是統計數據分布偏斜方向和程度的度量,反映了RF 信號包絡幅度值非對稱性程度,偏度小于0 說明RF 信號包絡中出現了少量值較小的數據點,偏度大于0說明RF信號包絡中出現了少量值較大的數據點。RF 數據元素X包絡的偏度按公式(6)計算:

式(6) 中,Skew(X)為RF 數據元素X包絡的偏度,μ為RF 數據元素X包絡的均值,σ為RF數據X包絡的方差,k2、k3分別表示二階和三階中心矩。

(5) 能量

信號能量是表征信號時域特性的一個重要參數,從超聲圖像可以直觀地看出與非B線區域相比,B 線區域更加明亮,灰度值更高,RF 數據同樣表現出相似的特點,即B 線區域幅值更大,信號能量更高,故本研究選擇信號能量作為B 線識別的特征參量,并進一步研究其B線識別性能。RF 信號元素X的能量按公式(7)計算:

式(7)中,E(X)為RF 信號元素X的能量,N為元素X中數據點的數量,Xj是元素X中第j個數據點的幅值。

1.2 基于非線性SVM的分類識別

SVM是一種監督學習方式的二分類模型,其決策邊界是對學習樣本求解的最大邊距超平面,是一種具有稀疏性和穩健性的分類器,可以對線性可分的樣本進行較準確的二元分類。對于線性不可分的學習樣本,非線性SVM 利用核函數,通過非線性映射算法將二位線性不可分的樣本映射到高維特征空間中,使得樣本點再高維空間線性可分,采用與線性SVM 相同的間隔最大化原則學習得到非線性SVM模型。

對于n維特征空間,目標超平面為

其中,w=(w1,w2,···,wn)為法向量;b為位移項。求解最大超平面:

利用拉格朗日乘子法,結合特征空間多維特征,得到其對偶問題:

其中,κ(xi,xj)稱為核函數,用于降低升維后特征空間的計算復雜度。本文主要使用的核函數為高斯核函數:

本研究從肺超聲RF 數據提取的多個特征參數組成的二維原始空間均為線性不可分,故采用非線性SVM進行識別分類。

2 結果與分析

2.1 實驗數據的獲取

本實驗選用西安交通大學醫學院實驗動物中心提供的實驗6只白兔,平均體質量2.3 kg。實驗器材和藥品包括兔臺、兔子固定箱、5 mL 注射器、留置針、0.9%氯化鈉溶液、分析純油酸、異氟烷。實驗儀器采用Sonix Touch 超聲診斷系統,線陣超聲探頭,頻率4~10 MHz。美國MATRX動物呼吸麻醉機,氧氣流量400 mL/min,異氟烷麻醉揮發罐檔位2~4。

實驗前對實驗兔禁食8 H,并稱重,將實驗兔俯臥位放入兔子固定箱,佩戴呼吸麻醉機面罩,調整氧氣流量至400 mL/min,麻醉揮發罐檔位調至3 檔,待兔子麻醉后撤下兔子固定箱,仰臥位固定于兔臺,胸部備皮。由左右兩側第四、第五肋間對肺進行超聲檢查,留取正常無B線、具有明顯胸膜線和A線的肺部超聲圖(圖1(a))和RF 數據。按0.1 mL/kg 經耳緣靜脈注射分析純油酸,制備急性肺損傷模型,每隔1 min 進行一次肺超聲檢查,記錄肺超聲圖像由存在少量B 線(圖1(b))到存在彌漫性B線(圖1(c))過程的RF數據和對應的超聲圖像共100組。

圖1 肺超聲圖Fig.1 Images of lung ultrasound

2.2 單參數分析

本文共選取100 幀RF 原始數據,依據B 線從胸膜線產生并延伸至屏幕底部、與肺滑行同步運動、呈現離散、激光樣垂直高回聲的偽影表現的特征,對所有識別對象標定B 線或非B 線標簽。分別統計B 線區域與非B 線區域各個特征參量的均值和標準差并使用獨立樣本t 檢驗判斷其是否有顯著差異性,如表1 所示。最后將各個參數隨機平均分為10份,輪流取其中9份作為訓練集,剩下1份作為測試集輸入貝葉斯分類器,得到單參數分類的平均靈敏度、特異性、準確率和曲線下面積(Area under curve,AUC)如表2所示。

表1 特征參數統計特征Table 1 Statistical characteristics of characteristic parameters

表2 特征參數單參數貝葉斯分類Table 2 Bayesian classification results of single characteristic parameters

由表1 可見,對本文提取的5 個特征參數分別進行獨立樣本t 檢驗后,結果顯示B 線與非B線的上述所有特征參數之間均具有顯著差異性(P <0.05)。可以發現,B 線排列熵均值(1.974)顯著低于非B 線排列熵均值(2.454);B 線信息熵均值(8.589)顯著高于非B 線信息熵均值(8.134)。這是因為信息熵與排列熵分別從兩個不同方面反映了信號的復雜程度,信息熵從信號幅值出現的概率密度評估信號在宏觀上的混亂程度,排列熵卻并不關注信號具體幅值,而是分析信號重構分量內數據升序排列順序出現的概率,排列熵側重于被分析信號內相似片段的重復性,而對相似片段內信號幅值較小的波動變化以及異常點并不敏感,相比于信息熵,排列熵弱化了B 線信號衰減等其他因素導致的信號幅值變化的影響,突出了B 線信號在微觀上的重復規律性,這使得表1 中數據顯示與非B 線區域相比,信息熵顯示B 線區域信號混亂程度更高而排列熵卻顯示其更具有規律性。

表2 中單參數的貝葉斯分類器分類結果顯示信息熵(靈敏度72.547%,特異性91.337%,準確率85.074%,AUC=0.91)和排列熵(靈敏度89.953%,特異性78.303%,準確率82.186%,AUC=0.90)的識別效果最好,峰度、偏度和能量3 個參數識別特異性較高,靈敏度較低且與信息熵和排列熵相比差距較大。表明排列熵與信息熵對于判斷RF 信號是否為B線的能力較上述其他特征參數更高。貝葉斯分類靈敏度較低的3個特征參數(峰度、偏度、能量)其標準差較大,B 線與非B 線相應特征參數數值重合區間內識別對象數量較多,導致大多數B 線數據被分類為非B 線,從而使得分類模型表現出靈敏度非常低而特異性非常高的特點。

2.3 基于非線性SVM的雙參數B線識別

由2.2 節可知,基于B 線和非B 線RF 數據提取的特征參數雖然有顯著差異性,但其數值均有一定程度上的重合,而這樣線性不可分的性質對貝葉斯分類器有較差的影響。為了獲得更好的分類效果并且尋找上述參數中識別效果最好的雙參數組合,選擇單參數貝葉斯分類中識別效果前三的參數(信息熵、排列熵和峰度)組成3 組雙特征參數組合,以及識別效果最差的兩個參數(偏度和能量)組成一組雙特征參數組合,將上述4 個特征參數組合輸入非線性SVM 進行分類,本研究使用的SVM 基于MATLAB 中LibSVM 工具箱,采用高斯徑向核函數。

不同參數組合如下:組合1 為排列熵加信息熵;組合2 為信息熵加偏度;組合3 為排列熵加偏度;組合4 為峰度加能量。對不同的參數組合采用十折交叉驗證的方法測試其分類效果:將數據集隨機平均分為10份,輪流將其中9 份作為訓練集,剩余1 份作為測試集,進行10 次模型訓練和B 線識別,計算不同參數組合的受試者工作特征AUC,平均靈敏度、平均特異性和平均準確率如表3所示。

表3 不同輸入組合十折交叉驗證結果統計Table 3 Ten fold cross validation results for different input combinations

由表3 可知,4 個組合輸入SVM 的識別效果最好的是信息熵加排列熵組合,將該組合作為輸入的B 線識別平均靈敏度為90.521%,平均特異性為98.106%,平均準確率為96.328%,AUC等于0.95,均為4 種組合中的最高值,可以較好的實現對B 線的識別。識別效果最差的是偏度加能量組合,將該組合作為輸入的B 線識別平均靈敏度為39.359%,平均特異性為83.228%,平均準確率為68.605%,AUC等于0.70。信息熵加峰度雙特征參數組合對B 線的平均識別靈敏度(靈敏度=65.071%)較信息熵加排列熵(靈敏度=90.521%)以及排列熵加峰度(靈敏度=84.138)較低,但由于其平均特異性(特異性=94.205%)較排列熵加峰度組合高,故組合二、組合三的識別平均準確率相近,分別為84.494%和83.083%。

與表2中單參數貝葉斯分類器識別效果對比可知,特征參數兩兩組成雙參數特征參數組合輸入非線性SVM 得到的B 線識別效果要優于該組合中單參數識別效果較差的特征參數的識別效果。當組成雙特征參數組合的兩個特征參數的單參數識別效果均較好時,比如信息熵加排列熵組合,其雙參數SVM 對B 線的識別效果顯著優于其單獨的識別效果。同時,當組成雙特征參數組合的兩個特征參數單參數識別效果均較差時,比如偏度加能量組合,其雙參數SVM 的識別效果雖然依然較差,但較其各自單參數識別效果有顯著的改善。

本文認為高斯核函數將特征空間新向量表示為原始空間映射函數的內積,構建新的高維線性可分特征空間,并建立一個線性超平面對樣本進行分類,新的特征空間包含了特征參數各自對于分類的相關性,使參數組合具有了與其各自對分類的相關性相比更高的相關性,提高了分類準確率。

2.4 后處理

由于本文基于RF 數據提取的B 線區域與非B線區域的信息熵與排列熵是線性不可分的,采用了高斯核的非線性SVM 進行分類,模型中引入了懲罰因子,因此在更高維度的特征空間中進行分類時,分類模型對于邊界間的樣本是有一定容忍度的,即允許一定的錯誤分類的情況出現,這是因為對于懲罰因子過大時,識別模型對錯分情況容忍度很低,會出現過擬合;懲罰因子過小時,可能導致無法達到分類效果,出現欠擬合現象。與之對應的,本文采用信息熵、排列熵雙參數非線性SVM 識別后的結果顯示,一張超聲圖像中識別錯誤的情況大致分為兩類:一是將非B 線區域中的若干掃描線識別為B線;二是將B線區域中的若干掃描線識別為非B線。圖2 展示了SVM識別效果圖,其中識別為B線的掃描線用藍色半進行了標注。

圖2 SVM 識別圖(識別為B 線的線束用藍色在圖中標注)Fig.2 SVM identification images (B-Line regions were marked in blue)

分析識別錯誤的識別對象,發現出現識別錯誤的線束大多為圖像中位置分散的單獨掃描線或少量連續掃描線,而B 線的寬度總是≥3 條掃描線的寬度,同一條B 線包含的掃描線在圖像中的位置是連續的,并且B 線從胸膜線產生,一直延續到屏幕底端。根據以上性質,本文在SVM識別的基礎上提出了一種后處理的算法,可以有效地提高B 線識別效果。

對于SVM 模型的識別結果,后處理算法首先遍歷其識別為B 線的掃描線位置,找到被識別為B線,但與其他B線間隔>5 (每根掃描線的寬度為1)且寬度≤2 的掃描線區域,將其識別結果更正為非B線。再次遍歷識別結果,找到所有間隔≤2的相鄰B 線,判斷其間隔中非B 線在屏幕底部位置是否依然具有較強的能量,若是,則認為此掃描線滿足B線由胸膜線一直延伸到屏幕底端的性質,將其識別結果更正為B線,若否,則不做更改。

圖3 展示了后處理算法提高B線識別效果的對比圖。可以看出,后處理算法很好地找到B 線區域中和非B 線區域中識別錯誤的掃描線,并有效地更正了識別結果,得到了更加準確的B線位置與范圍。表4 展示了將排列熵加信息熵輸入SVM 后所有識別結果輸入后處理算法得到的B 線識別平均靈敏度、特異性和準確率。可以看出,后處理算法對B線識別的靈敏度有顯著提升,達到95.23%,準確率也有提高(準確率96.88%),特異性有較小下降(特異性97.22%),綜合3項指標可知后處理算法可以顯著有效地提升B 線識別能力,以較小的運算量改善SVM識別效果。

表4 后處理前后B 線識別效果對比Table 4 Comparison of B-Line recognition effect before and after post-processing

圖3 原始B 超圖像、SVM 識別B 線效果圖和后處理效果圖對比(圖中藍色標注的線束為識別為B 線的線束)Fig.3 Comparison of original B mode images,SVM recognition images and post-processing images

3 討論

本研究基于超聲RF 數據進行肺超聲特殊征象B 線的識別,選取了RF 數據的5 種特征參數,包括信息熵、排列熵、包絡峰度、偏度和RF 信號能量。為了確定對于B 線分類效果最好的特征參量以及最好的B 線識別效果,對各個參數進行了獨立樣本t 檢驗,采用單參數貝葉斯分類器與雙參數非線性SVM,并選擇了4 種不同的特征參數組合輸入SVM,對100幀RF 數據進行了分類識別實驗,進行了各特征參數的統計學分析,比較了各個參數單參數分類的效果、不同參數組合輸入非線性SVM 的識別效果以及單參數與雙參數的識別效果,最后在信息熵與排列熵雙參數SVM 識別的基礎上提出并采用了一種后處理算法,并對后處理算法的效果進行了比較分析。

結果顯示B 線與非B 線的所有特征參數之間均具有顯著差異,但不同特征參數對于B 線識別的相關性差別較大,根據單參數貝葉斯分類結果,標準差較小的排列熵和信息熵相關性較高,而標準差較大的峰度、偏度和信號能量相關度較低,說明對于兩種類別RF 信號的特征參數,具有顯著差異性且標準差較小的特征參數更有利于B線的識別。

與單參數貝葉斯分類相比,雙參數SVM 分類的識別靈敏度、特異性、準確率以及AUC 均有顯著的提高,且單參數分類效果好的特征參數作為雙特征參數輸入組合的一部分時,分類準確率也高于那些含有單參數分類效果較差的特征參數的參數組合,說明多參數輸入分類器時會提供各自對于分類識別的能力,并且參數組合不僅僅包含單參數分類能力,還會由于參數間隱含的相關關系,使得參數組合與分類的相關度更高,進一步提高識別準確率。但本研究并未比較雙參數組合與更高維度的特征參數組合輸入SVM 時的分類效果,后期研究會采用更高維度的特征參數組合,探究更高維特征空間對于B線識別的能力。

4 種輸入非線性SVM 的特征參數組合中,信息熵加排列熵的分類效果最好,平均靈敏度為90.521%,平均特異性為98.106%,平均準確率達到96.328%,AUC達到0.95,均大于其他3 種輸入特征參數的組合。進一步采用后處理算法對信息熵加排列熵的分類結果進行處理,顯著改善B線識別效果,平均靈敏度達到95.23%,平均特異性為97.22%,平均準確率為96.88%,實現了對B線的較準確識別。

本文提出的這種基于超聲RF 數據SVM 雙參數(信息熵和排列熵)B線識別方法相較于Cristiana等[8]提出的基于卷積神經網絡的深度學習算法對計算機算力的要求更低,計算速度更快,有利于對B線的實時快速識別。Cristiana等[8]分析識別效果時計算的靈敏度、特異性等指標是以某超聲圖像是否存在B 線為標準,比較模型識別與專家判斷結果是否一致為依據,若模型識別與專家判斷某超聲圖像中存在B 線則認為識別正確,而本文是將超聲圖像的每條掃描線包含的所有RF 數據作為一個最小的識別對象,分析識別效果時計算的靈敏度、特異性與準確率是以每一條掃描線是否屬于B 線為標準,比較模型識別與專家判斷結果是否一致,識別結果不僅包含被識別的圖像中是否存在B 線的信息,并且指明了B 線的具體位置與范圍,可以為定量研究B線與肺病變嚴重程度的關系提供重要依據。

本研究選用的超聲RF 數據來自油酸致肺損傷動物模型實驗,實驗中出現B 線的具體病理可能并不統一,具有一定的局限性。在后續研究中會與相關醫院超聲科合作,獲得可供研究的具體病例的RF 數據,將B 線識別與產生B 線的具體病理特征結合起來,以提高B線識別的可靠性與特異性。

4 結論

本文提出了一種基于超聲RF 信號的肺臟超聲特殊征象B 線的識別方法,提取肺超聲RF 數據的信息熵、排列熵、峰度、偏度和能量作為特征參數,利用單參數貝葉斯分類器以及非線性SVM,獲得并比較了不同參數以及不同參數組合對于B 線識別的效果,并將識別效果最好的信息熵加排列熵參數組合的識別結果輸入后處理算法,最終得到B 線識別平均靈敏度為95.23%,平均特異性為97.22%,平均準確率為96.88%,能夠以較高的靈敏度、特異性和準確率實現對肺超聲B 線的識別。與基于圖像的B線識別方法相比,本文提出的基于RF 數據的B 線識別方法受到超聲儀器參數調整的影響更小,可以適用于多種不同超聲儀器,普適性較強,應用場景較廣。同時,基于RF數據提取的特征參量包含更多組織結構信息,對B線的識別更加可靠有效,可以為臨床經驗有限的醫生提供準確的輔助診斷信息,減緩肺超聲診斷的學習曲線,提高診斷準確率。

猜你喜歡
分類信號效果
按摩效果確有理論依據
分類算一算
信號
鴨綠江(2021年35期)2021-04-19 12:24:18
完形填空二則
分類討論求坐標
迅速制造慢門虛化效果
基于FPGA的多功能信號發生器的設計
電子制作(2018年11期)2018-08-04 03:25:42
數據分析中的分類討論
抓住“瞬間性”效果
中華詩詞(2018年11期)2018-03-26 06:41:34
教你一招:數的分類
主站蜘蛛池模板: 91精品专区国产盗摄| 国产精品所毛片视频| 国产精彩视频在线观看| h网址在线观看| 精品无码国产自产野外拍在线| 色天堂无毒不卡| 爆乳熟妇一区二区三区| 97se亚洲综合在线天天| 亚洲欧洲日产无码AV| 欧美日韩福利| 亚洲日本中文综合在线| 99久久国产综合精品2023| 日韩av在线直播| 亚洲高清在线播放| 亚洲欧洲日产国码无码av喷潮| 综合人妻久久一区二区精品| 在线视频一区二区三区不卡| 免费看美女毛片| 国产国产人在线成免费视频狼人色| 女高中生自慰污污网站| 伊人激情久久综合中文字幕| 性激烈欧美三级在线播放| 亚洲国产精品一区二区高清无码久久| 国产成人精品一区二区秒拍1o| 久久婷婷五月综合97色| 国产在线精彩视频二区| 精品自窥自偷在线看| 黄色在线网| 又爽又大又黄a级毛片在线视频| 尤物亚洲最大AV无码网站| 老司机精品99在线播放| 国产永久在线视频| 97久久精品人人| 色综合中文| 国产一级毛片在线| 亚洲精品不卡午夜精品| 91最新精品视频发布页| 71pao成人国产永久免费视频 | 亚洲丝袜第一页| 婷婷六月综合| 欧美亚洲网| 欧美成人一级| 久久香蕉欧美精品| 亚洲乱码精品久久久久..| 亚洲中文字幕无码爆乳| 这里只有精品国产| 日本黄色a视频| 日韩免费无码人妻系列| 人人爱天天做夜夜爽| 国内毛片视频| www.亚洲一区| 国产精品九九视频| 色综合天天娱乐综合网| 亚洲人精品亚洲人成在线| 国产激情无码一区二区APP| 制服丝袜亚洲| 欧美成人精品一级在线观看| 成人免费黄色小视频| 伊人久久综在合线亚洲2019| 精品国产自在现线看久久| 国产第一页亚洲| 国产精品美女免费视频大全| 亚洲国产精品无码久久一线| 国内熟女少妇一线天| 青青草原偷拍视频| 色丁丁毛片在线观看| 日韩一区二区三免费高清| 亚洲国产av无码综合原创国产| 国产真实乱人视频| 国产视频大全| 中文字幕波多野不卡一区| 免费a在线观看播放| 无码aⅴ精品一区二区三区| 综合色天天| 巨熟乳波霸若妻中文观看免费| 国产综合日韩另类一区二区| 国产乱子伦无码精品小说| 制服丝袜 91视频| 久久精品亚洲中文字幕乱码| 国产女人在线观看| 丁香婷婷激情网| 国产精品毛片一区视频播|