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

放射多組學協(xié)同學習預測鼻咽癌自適應放療觸發(fā)機制

2024-04-09 01:41:36邱成羽李兵林世杰盛嘉寶滕信智張將程煜婷張馨勻周塔葛紅張遠鵬蔡璟
智能系統(tǒng)學報 2024年1期
關鍵詞:劑量特征模型

邱成羽,李兵,林世杰,盛嘉寶,滕信智,張將,程煜婷,張馨勻,周塔,4,葛紅,張遠鵬,4,蔡璟,4

(1.香港理工大學 健康科技與資訊學系, 香港 999077; 2.南通大學 醫(yī)學信息學系, 江蘇 南通 226019; 3.鄭州大學附屬腫瘤醫(yī)院, 河南 鄭州 450008; 4.香港理工大學深圳研究院, 廣東 深圳 518057; 5.香港理工大學 生物醫(yī)學工程學系, 香港 999077)

鼻咽癌(nasopharyngeal carcinoma, NPC)是指發(fā)生于鼻咽腔頂部和側(cè)壁的惡性腫瘤。根據(jù)國際癌癥研究機構的數(shù)據(jù)[1],2020年,在世界范圍內(nèi)估計有1 930萬新的癌癥病例,其中大概有13.3萬例鼻咽癌新發(fā)病例,死亡人數(shù)約為8萬人,表現(xiàn)出較高的死亡率。在中國,鼻咽癌的發(fā)病人數(shù)明顯高于世界平均水平(1.9萬/10萬)[2-3],主要發(fā)生在香港和南方五省(廣東、廣西、湖南、福建和江西),是當?shù)仡^頸部惡性腫瘤的首位。

隨著人工智能技術的發(fā)展,影像組學被廣泛應用于癌癥的臨床分期、治療方法選擇、療效評估以及預后分析等領域,并且取得了一定的研究成果。在鼻咽癌領域,Liu等[4]從鼻咽癌患者的CT圖像提取影像組學數(shù)據(jù)并與患者的唾液量相結合,通過預測患者的預期唾液量來構建預測患者未來急性口干燥癥的可能性和程度的模型,該模型能在早期準確地預測患者的唾液量,并提前預防口干癥狀;金哲等[5]從1 393例鼻咽癌患者的MRI圖像中提取影像組學特征并結合患者的臨床指標構建鼻咽癌遠端轉(zhuǎn)移預測模型,在預測鼻咽癌遠端轉(zhuǎn)移風險方面具有較高效能。

但影像組學數(shù)據(jù)也有自己的局限性[6-9],由于各個醫(yī)院之間所采用的機器不同、協(xié)議不同、參數(shù)不同,以及負責感興趣區(qū)域(volume of interest,VOI)標注的臨床醫(yī)師的不同,使得影像組學特征的再現(xiàn)性受到影響,進而影響模型的準確性和可重復性。因此,影像組學數(shù)據(jù)融合其他組學或臨床指標的方法逐漸開始受到重視。將影像組學與其他組學相結合的研究方式,可以充分利用來自不同組學特征進行互補,也可以避免單組學特征冗余度過高的問題,以提高模型的預測精度和泛用性。例如, Zhou等[10]結合影像組學和基于EQD2的劑量組學,建立了基于CT的放射性肺炎(padiation poisoning, RP)預測模型,發(fā)現(xiàn)劑量組學結合影像組學的訓練效果要優(yōu)于傳統(tǒng)的劑量-體積直方圖(dose-volume histogram, DVH)、影像組學+DVH,且具有統(tǒng)計學意義;Cai等[11]將水平計算機斷層掃描(CT)融合劑量學,所建立的融合劑量學模型均優(yōu)于單一的CT/劑量模型;Cui等[12]在深度學習中整合了PET(positron emission tomography)影像組學、細胞因子和miRNA的多組學信息,構建了預測正常組織并發(fā)癥概率/腫瘤控制概率的模型,要優(yōu)于傳統(tǒng)的預測模型。

在多組學研究中,處理多組學數(shù)據(jù)傳統(tǒng)方法之一是將來自不同組學的特征向量連接成特征新向量,然后在連接后的向量上應用單視角算法進行處理。而這樣簡單的特征拼接處理容易帶來以下問題:1)忽略了不同組學之間的信息互補關系;2)特征數(shù)量會隨著組學數(shù)量的增加而增加,造成維度災難。最終導致過擬合,模型精度較差等問題[13]。多組學問題應當基于多視角學習的原則設計[14]:1)假設視角之間的分歧是分類誤差的上限,力求最大限度地保證每個視角的預測一致性;2)假設每個視角都包含其他視角沒有的信息,從每個視角中提取差異,同時保留共享信息。因此,采用合適的多組學融合方案是多組學研究的關鍵。對此,本研究針對影像組學、劑量組學和輪廓組學,設計了一種多組學協(xié)同學習算法(multi-omics collaborative learning, MOCL)。MOCL通過一致性約束挖掘不同組學特征之間的互補模式,另外還通過香農(nóng)熵自適應學習不同組學特征的權重。通過在真實臨床影像數(shù)據(jù)上的實驗結果表明,MOCL對鼻咽癌遠端轉(zhuǎn)移有較好的預測作用。另外,和引入的對比算法的比較結果表明,MOCL在多組學協(xié)同學習上具有一定的優(yōu)勢。

1 放射多組學協(xié)同學習

1.1 問題提出

鼻咽癌患者在治療過程中,從其CT、MRI等常規(guī)放射影像學資料中,可以獲取的組學特征包括影像組學、劑量組學以及輪廓組學。已有研究表明,影像組學、劑量組學以及輪廓組學特征均和患者的治療預后密切相關[15-16],只使用單一組學構建模型往往會出現(xiàn)模式信息有限的問題,如何有效地挖掘不同組學特征所包含的互補信息對于準確預測患者預后來說顯得非常重要。傳統(tǒng)的學習策略將各組學進行拼接,人為地割裂了各個組學之間潛在的聯(lián)系,只是簡單地將單獨學習的結果進行集成,很難達到較為滿意的學習效果。因此,基于多組學協(xié)同學習的問題描述如下。

假設{X,Y}為一個多組學訓練集,其中,每個組學特征空間可以表示為Xk∈RN×d,則多組學協(xié)同學習的問題可以形式化描述為

式中:第1項為各個組學特征空間上的損失之和;第2項為各組學之間的一致性約束,用于保障同一對象在不同組學特征空間預測的一致性;第3項為泛化項,用于保障模型的泛化能力。

1.2 多組學特征空間損失

在使用機器學習進行線性分類的過程中,對于每個樣本,可以使用一個轉(zhuǎn)置矩陣將其映射到一個低維空間中,即模型對該樣本的預測值f(X)。但需要對模型的預測值f(X)與真實值Y的差異程度的進行估算。兩者相差值越小,說明模型的預測精度就越好。在描述多組學的特征空間損失時,具體定義如下:

式中:K是數(shù)據(jù)集中所包含的組學的數(shù)量,ωk是第k個組學所占的權重,Xk為第k個組學的特征向量,Ak代表第k個組學特征的轉(zhuǎn)置矩陣。通過最小化式(2),能保證訓練過程中得到的標簽和真實結果誤差盡可能最小。

1.3 一致性約束

為了確保各組學特征空間之間預測結果的一致性,即認為同一個鼻咽癌患者,在不同視圖下,其放療預后的結果應該盡量保持一致,設置一個先驗知識,具體定義如下:

其中是針對每個組學特征空間單獨預先訓練得到的一個轉(zhuǎn)換矩陣。通過最小化式(3),能保證從任意組學特征空間進行預測時,都能與其余組學特征空間得到的預測結果的差異最小,保證同一樣本在不同組學特征空間預測的結果盡量一致。

同時為了從不同的角度來探索一致或互補的信息,本文提出了一種基于信息論的多視圖加權機制,即自適應權重。具體來說,引入了一個基于“香農(nóng)熵”的視圖加權懲罰項來自適應學習每個組學的權重系數(shù),對于包含模式信息較少的特征空間,為其分配較小的權重,反之,則為其分配較大的權重其定義如下:

通過最小化式(4),能使得包含更多有效特征的組學在訓練中被給予更大的權重。

1.4 泛化性能提升

為了減少在迭代過程中可能出現(xiàn)的過擬合現(xiàn)象,基于流形學習的方法,采用緊致度圖來盡可能避免過擬合問題。在緊致度圖中,同一標簽的兩個不同的樣本對應的兩個節(jié)點由一條無向邊連接,其權重定義如下:

其中σ是熱核參數(shù)。如果Wij能夠保證具有相同標簽的樣本彼此應該相互靠近,那么兩個具有相同標簽的樣本fi和fj就應該被賦予更大的權重Wij[17]。選擇映射的合理標準是最小化以下目標函數(shù):

其中:fi和fj表示xi和xj進行變換后的結果,通過最小化該式,能進一步確保fi和fj在空間中的相對接近。因此,可以通過式(6),保持樣本在變換空間中的緊致性。再引入拉普拉斯矩陣之后,式(6)可以進一步變化為

式中:Lk是拉普拉斯矩陣,定義為L=Z-W;Z是對角矩陣,它的對角項可以由計算。

1.5 目標函數(shù)和優(yōu)化

基于式(2)~(4)、(7)所設計的組件,多組學協(xié)同學習算法的目標函數(shù)可以表示為

其中,η、λ、β作為平衡參數(shù)且>0,以調(diào)配目標函數(shù)中不同部分的影響。

定理1針對J(ωk,Ak),當ωk固定時,求得最佳的AK計算如下:

證明給定任意ωk,式(8)中的目標函數(shù)可以重新表達為

證畢。

定理2針對J(ωk,Ak),當Ak固定時,可以求得最佳ωk:

證明給定任意的Ak,式(8)中的目標函數(shù)可以被重新表達為

通過將J(ωk)對ωk的偏導數(shù)設置為0,可以得到

證畢。

1.6 MOCL算法描述

根據(jù)Ak和ωk的 計算原則,具體算法計算流程如下。

算法1MOCL訓練

輸入多組學訓練樣本{Xk,Y},正則化參數(shù)η、λ、β,迭代次數(shù)t

輸出Ak和ωk

1) 計算影像組學、輪廓組學和劑量組學的轉(zhuǎn)置矩陣;

2) 為每一個組學特征空間計算Lm;

3) 以1/K為每個組學特征空間賦予初始權重ωk重復;

4) 為每個組學特征空間計算Ak;

5) 為每個組學特征空間計算ωk直到

6) 返回Ak和ωk。

算法2MOCL測試

輸入測試集,Ak和ωk

輸出

2 實驗分析

2.1 數(shù)據(jù)來源和預處理

2.1.1 患者數(shù)據(jù)

本研究對2012—2015年期間在香港伊利沙伯醫(yī)院接受放療的311名NPC患者進行了回顧性分析。研究對象包括符合以下條件的患者:1)初次診斷為經(jīng)組織活檢證實的原發(fā)性鼻咽癌,無遠處轉(zhuǎn)移和合并其他類型腫瘤;2)接受了治愈性的同步放化療(concurrent chemoradiation,CCRT)或CCRT加輔助化療(AC);3)接受了螺旋調(diào)強放療治療。被排除的患者包括:1)在接受CCRT治療之前接受過新輔助化療的患者;2)接受放療但未同時接受化療的患者;3)未注射造影劑進行計劃對比增強CT(contrast enhancement CT,CECT)或計劃對比增強T1-w(CET1-w)磁共振成像的患者;4)沒有完整的臨床/影像數(shù)據(jù)的患者。個體患者以在放療后是否需要進行自適應復查(adaptive review)來進行分類,其中需要進行自適應復查的標記為1,不需要的標記為0。

2.1.2 感興趣區(qū)域

本研究涉及8個不同的器官結構的VOI,包括原發(fā)性NPC腫瘤(GTVnp)和轉(zhuǎn)移淋巴結(GTVn)的總腫瘤體積、同側(cè)腮腺(IpsiPG)、對側(cè)腮腺(ContraPG)、腦干(BS)、脊髓(SC)、淋巴結計劃靶體積的高劑量和低劑量區(qū)(分別為處方劑量為70-Gy的PTVn_high_dose和處方劑量為60-Gy的PTVn_low_dose)。圖1結出了本研究中涉及的每個VOI的位置。

其中GTVnp是與計劃MRI圖像對準后,在軸向CT層面上手動勾畫出來的,GTVn是由一名豐富經(jīng)驗且獲得認證的頭頸放射腫瘤學家按照國際頸部水平CT描繪指南在CECT圖像上繪制的。為了區(qū)分患者的同側(cè)腮腺和對策腮腺,本文按照以下策略進行:1)確定腮腺(PG)體積上特定體素點與GTVnp表面上所有體素點之間的最小幾何距離;2)選擇另一體素點,重復(1)過程,直到確定整個PG上所有體素點相對于GTVnp的最小距離;3)計算所有最小距離的中值,并以此確定各個患者的PG與GTVnp的總體接近度。中值較小的認為是同側(cè)腮腺,較大的認為是對側(cè)腮腺。所有分割均使用Varian ARIA和Eclipse治療計劃系統(tǒng)v13(加利福尼亞州帕洛阿爾托的Varian醫(yī)療系統(tǒng)公司)進行。

2.1.3 影像組學特征

在實驗中,由GTVnp、GTVn、IpsiParotid和ContraParotid 4種來自不同器官結構的VOI參與影像組學的特征提取的計算。使用Python中的SimpleITK(v1.2.4)和Pyradiomics(v2.2.0)對原始醫(yī)學圖像進行提取。一般來說,影像組學特征可以分為3類:一階統(tǒng)計特征、形態(tài)特征、二階及高階紋理特征。一階統(tǒng)計特征反應所測體素的對稱性、均勻性以及局部強度分布變化。包括中值、平均值、最小值、最大值、標準差、偏度、峰度等。形態(tài)特征定量描述感興趣區(qū)的幾何特性,如腫瘤的表面積、體積、表面積和體積比、球形度、緊湊性和三維直徑等,這些特征可以描述腫瘤三維的大小和形態(tài)信息。二階及高階紋理特征能夠反映圖像體素灰度之間的空間排列關系,可以分為灰度共生矩陣(gray-level covarionce matrix, GLCM)、灰度游程長度矩陣(gray-level run-length matrix,GLRLM)和灰度尺寸區(qū)域矩陣(gray-level size zone matrix, GLSZM)。在對CECT、CET1-w和T2-w磁振圖像進行影像組學特征計算的過程中,以高斯-拉普拉斯濾波器(核大小:1、3、6 mm)和小波濾波器(HHH、HLL、LHL、LLH、LHH、HLH、HHL、LLL)進行濾波。在本研究中,將4個VOI的形態(tài)學特征全部分離出來,為每個研究的VOI提供了總共6 348個影像組學特征。同時將4個VOI(GTVnp、GTVn、IpsiPG和ContraPG)的伸長、平整度、最小軸長、長軸長、小軸長、最大2D直徑列、最大2D直徑行、最大2D直徑切片、最大3D直徑網(wǎng)格體積、球度、表面積、表面體積比、體素體積14個形態(tài)學特征組合起來的56個特征作為一組。影像組學和形態(tài)學特征的詳細定義可在Pyradiomics文檔中找到。

2.1.4 劑量組學特征

因為傳統(tǒng)的劑量-體積直方圖(DVH)不包含輻照器官內(nèi)空間劑量分布的信息,所以本文采用RT劑量數(shù)據(jù)對8種不同器官結構的VOI進行劑量學特征計算。相比之下,劑量組學能夠表征所研究的8個VOI內(nèi)局部輻射劑量分布的空間格局,在癌癥預后和治療反應的各種預測模型中得到了廣泛的研究[18-19]。本研究采用Gabryovic等[20]的方法計算了8種VOI的DVH曲線點的劑量學特征,示例包括但不限于最大劑量、最小劑量、平均劑量、至少接受一定劑量水平的VOI的體積以及一定體積的VOI接受的最小劑量。此外,本文提取了每個研究VOI內(nèi)的空間劑量分布,以全面描述沉積劑量的異質(zhì)性,如沿3個成像軸(x、y和z方向)的劑量梯度。這些特征的定義在Buettner等[21]之前的出版物中描述過。同時,本文將每個研究VOI內(nèi)的三維劑量分布都轉(zhuǎn)換為3D圖像,以便后續(xù)使用Pyradiomics包計算與影像組學特征相似的劑量組學特征;所提取的特征包括一階劑量統(tǒng)計量、灰度依賴矩陣(gray-level dependence matrix, GLDM)、GLCM、GLRLM、GLSZM和相鄰灰度差分矩陣(neighboring gray tone difference matrix, NGTDM)等。本研究共從8個VOI中提取了1 608個劑量組學特征。

2.1.5 輪廓組學特征

在本研究中,是否需要實施自適應放療(adaptive radiotherapy, ART)是由頭頸部不同內(nèi)部器官幾何關系的變化決定的。故我們提取了能描述4對器官結構的VOI (GTVnp與IpsiPG、GTVnp與ContraPG、GTVnp與SC、PTVn_low_dose與SC)之間復雜幾何關系的特征。這些特征是從RT等高線數(shù)據(jù)中提取的,作為“輪廓組學”特征。對于每個VOI對,從距離描述符重疊體積直方圖(overlapping volume nistogram, OVH)中計算一系列輪廓特征;例如,在治療計劃階段,SC與PTVn_low_dose之間的最大和最小距離分別作為零和滿體積時OVH上的距離。在本研究中,OVH的計算使用了Wu等[22]設計的算法。此外還使用投影重疊體積(projected overlap volume, POV),即一個VOI與另一個VOI在特定投影角度的平行投影重疊,來進一步揭示VOI對潛在的輪廓特征。本研究共從4對VOI中提取了132個輪廓特征。表1總結了計算所研究的4種組學特征所涉及的VOIs的來源(“—”表示無相應組學特征)。

2.2 特征選擇

在最小化模型過度擬合的風險時,特征選擇(feature selection, FS)能減少建模階段的冗余和不相關性較高的特征,是機器學習中必不可少的步驟。在這項研究中,仿照了Li等[23]的方法選取了6種常用的無監(jiān)督FS算法和4種有監(jiān)督FS算法產(chǎn)生24組FS組合。進行多次隨機分割作為訓練集和驗證集,以每種FS組合的可分辨性和再現(xiàn)性的乘積作為評價指標,且各組學獨立進行特征選擇過程,互相無影響。最佳的特征選擇方法依次是費舍爾得分Fisher Score和Pearson Score、T_score和Pearson Score、T_score和Lap_score,組成311×35、311×4和311×8的數(shù)據(jù)集。

2.3 實驗設置

為了綜合評估MOCL算法對于多組學數(shù)據(jù)的預測性能,選取了3種傳統(tǒng)的單視角機器學習算法SVM、C4_5決策樹、Adaboost元算法和Zhang等[24]提出的DICS、Yang等[25]提出的WeightReg兩種多視角算法進行比較。針對單視角算法,采用傳統(tǒng)的特征拼接方法,即影像組學、劑量組學、輪廓組學的3組特征向量拼接成一組新特征向量,進行歸一化后直接使用單視角算法對拼接后的向量進行處理。本實驗的評價指標設置為接收者操作特征曲線下的面積(area under curve, AUC)和F1_score。

為了進行公正的比較,每一種算法都進行了相對公平的參數(shù)設定和評價指標。所有的算法都是在搭載了win10操作系統(tǒng)的機器上運行的。具體流程如圖2所示,其中支持向量機(support vector machines, SVM)算法和MOCL算法都額外使用10折交叉驗證尋找最佳參數(shù);由于Adaboost元算法是通過迭代的方式不斷優(yōu)化分類器的性能,故取20個分類器中最優(yōu)分類器的實驗精度作為10次隨機劃分的結果并求均值。表2給出了實驗過程中,各算法參數(shù)的設定值。

表2 各算法參數(shù)設定值Table 2 Parameter settings for each algorithm

圖2 算法流程Fig.2 Algorithmic process

2.4 實驗結果和分析

首先,給出MOCL及其對比算法的分類準確率,如表3所示。可以看到和傳統(tǒng)的單視角機器學習算法相比,采用自適應權重和一致性約束的MOCL算法能取得較為優(yōu)秀的實驗結果。其中SVM和Adaboost元算法也都有較好的結果,可能是因為所使用的數(shù)據(jù)中,影像組學的特征占比較大,對于拼接后的特征向量,劑量組學和輪廓組學特征帶來的影響比較小。同時相對于兩種針對多視角的算法,MOCL算法明顯取得了更加優(yōu)秀的實驗結果。DISC采用聯(lián)合非負矩陣分解將多視圖數(shù)據(jù)分解為公共部分和視圖特定部分,而這種方法的結果對于子空間維數(shù)相當敏感,在面對不同的設定時結果差異很大;WeightReg算法設計了一種基于判別回歸的框架將多視圖數(shù)據(jù)映射到統(tǒng)一的低維判別子空間,并引入了一組可學習的權重參數(shù)保留原始視圖的相關和互補信息在投影子空間中,這種設計在處理不同視角特征數(shù)目不等的真實數(shù)據(jù)集合時可能存在一定的不足。

表3 MOCL及其對比算法結果平均值(標準差)Table 3 Mean of results for MOCL and its comparison algorithims (standard deviation)

2.5 消融分析

為了確定采用多組學協(xié)作學習算法構建模型和使用單一組學特征構建模型確實存在差異,將MOCL算法進行消融,即去除其中有關自適應權重和一致性約束的計算,將其改造為針對單一組學特征空間的算法。對3個組學特征進行單獨10次隨機劃分后、訓練和驗證,并以10次驗證結果的均值作為最終結果,結果如表4所示。可以看出,綜合考慮AUC和F1_score兩項評價指標,采取多組學協(xié)作學習的方式構建模型相比單獨考慮單一組學特征空間,能取得更好的預測精度。

表4 MOCL與其消融算法結果對比平均值(標準差)Table 4 Mean of results for MOCL and its ablation algorithms (standard deviation)

2.6 統(tǒng)計分析

為了進一步確定每一種對比算法和MOCL算法是否存在統(tǒng)計學上的差別,依次將對比算法與MOCL進行了T檢驗,設定α為0.05,具體結果如表5所示。可以看出,MOCL算法與SVM算法在AUC不存在顯著的差異,能達到類似的性能,但從F1的角度評估與SVM還有一定差距。而與其余算法相比,MOCL改良算法至少在一項評價指標上顯現(xiàn)出顯著差異,證明無論是與單視角算法對比還是已有的多視角算法對比,MOCL都有較為優(yōu)秀的分類性能和穩(wěn)定性。MOCL消融算法和MOCL的顯著差異也更一步體現(xiàn),使用多組學協(xié)作學習的方式相較于應用單一組學的方式能取得更好的模型性能。

表5 MOCL與其對比算法統(tǒng)計學比較Table 5 Statistical comparison of MOCL with its comparison algorithm

3 結束語

由于影像組學特征的再現(xiàn)性容易受到多種因素影響,使用多組學融合進行研究的方法越來越收到重視。而如何挖掘不同組學之間的互補信息是多組學研究的關鍵。針對此,本文設計了一種基于放射多組學的協(xié)同學習算法用于預測鼻咽癌患者是否需要進行ART,并在來自香港伊利沙伯醫(yī)院311名患者的真實數(shù)據(jù)上完成驗證。相比已有的一些針對多視角的算法,能夠更好地適應多組學數(shù)據(jù),有較好的實驗精度。但相比SVM、Adaboost兩種經(jīng)典的機器學習算法,沒有顯示出統(tǒng)計學上的明顯差異。可能是因為所使用的數(shù)據(jù)集中,劑量組學和輪廓組學特征數(shù)量較小,在計算過程中權重影響較低,在接下來的研究中會著重于劑量組學和輪廓組學的特征提取和選擇,盡可能擴大其特征數(shù)量。總的來說,本算法在預測鼻咽癌患者是否需要進行ART時,可以為臨床決策提供較為可靠的參考意見,更早發(fā)現(xiàn)有可能發(fā)生轉(zhuǎn)移的病人,提高病人的存活率。

猜你喜歡
劑量特征模型
一半模型
結合劑量,談輻射
·更正·
全科護理(2022年10期)2022-12-26 21:19:15
90Sr-90Y敷貼治療的EBT3膠片劑量驗證方法
重要模型『一線三等角』
重尾非線性自回歸模型自加權M-估計的漸近分布
如何表達“特征”
不忠誠的四個特征
當代陜西(2019年10期)2019-06-03 10:12:04
抓住特征巧觀察
3D打印中的模型分割與打包
主站蜘蛛池模板: www.99精品视频在线播放| 毛片网站在线播放| 91精品国产福利| 国产欧美日韩专区发布| 亚洲精品动漫在线观看| 中文字幕乱码二三区免费| 午夜少妇精品视频小电影| AV在线麻免费观看网站| 欧美中出一区二区| 在线欧美一区| 人人爱天天做夜夜爽| 国产精品区视频中文字幕| 婷婷综合在线观看丁香| 欧美一级大片在线观看| 国产正在播放| av免费在线观看美女叉开腿| 国产导航在线| 精品福利国产| 99久久成人国产精品免费| 国产真实乱了在线播放| 日韩午夜福利在线观看| 日韩av高清无码一区二区三区| 伊人成人在线| 国产精品青青| …亚洲 欧洲 另类 春色| 58av国产精品| 91成人免费观看在线观看| 亚洲一级毛片在线播放| 亚洲综合香蕉| 亚洲精品第五页| 视频二区亚洲精品| 亚洲一区二区三区国产精品| 日韩精品久久无码中文字幕色欲| 激情国产精品一区| 日本一区中文字幕最新在线| 国产99精品久久| 亚洲国产日韩视频观看| 国产精品成人一区二区不卡 | 国产欧美视频在线观看| 久久精品国产国语对白| 永久免费av网站可以直接看的 | 天堂中文在线资源| 特级毛片8级毛片免费观看| 亚洲侵犯无码网址在线观看| 色综合五月婷婷| 亚洲日本中文字幕天堂网| 18禁黄无遮挡免费动漫网站| 伊人激情综合网| 97精品久久久大香线焦| 丰满人妻中出白浆| 欧美日一级片| 亚洲人人视频| 欧美va亚洲va香蕉在线| 中文字幕中文字字幕码一二区| 国产日韩欧美一区二区三区在线| 国产精品永久久久久| 久久精品国产999大香线焦| 91色爱欧美精品www| 五月婷婷伊人网| 三区在线视频| AV不卡国产在线观看| 国产成人精品无码一区二| 精品福利视频导航| 亚洲 日韩 激情 无码 中出| 国产精品第一区在线观看| 在线观看国产黄色| 98超碰在线观看| 国产在线一区视频| 精品亚洲麻豆1区2区3区| 三级欧美在线| 在线网站18禁| 拍国产真实乱人偷精品| 国产成人午夜福利免费无码r| 久久久久久久久亚洲精品| 97在线免费| 日韩欧美在线观看| 欧美国产精品不卡在线观看| 国产精品无码AV中文| 日本在线亚洲| 丰满人妻中出白浆| 国产肉感大码AV无码| 无码中字出轨中文人妻中文中|