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

基于機器學習方法構建IDH野生型膠質母細胞瘤預測模型研究

2024-06-27 14:25:24許廣智張佳樂伊西才魏禮洲劉衛平
臨床神經外科雜志 2024年3期
關鍵詞:模型

許廣智 張佳樂 伊西才 魏禮洲 劉衛平

【摘要】 目的 建立異檸檬酸脫氫酶(IDH)野生型膠質母細胞瘤生存概率的列線圖模型及隨機生存森林模型。方法 回顧性分析2017年1月—2020年12月在空軍軍醫大學附屬西京醫院手術治療的127例IDH野生型膠質母細胞瘤患者臨床資料,進行預后因素分析并建立列線圖模型及隨機生存森林模型,通過C指數,校準曲線,決策曲線評價模型的區分度,校準度以及臨床凈獲益率。結果 使用Cox比例風險模型進行多因素分析發現,患者術前卡氏功能狀態評分(KPS)、是否接受同步放化療、年齡、O6-甲基鳥嘌呤甲基轉移酶(MGMT)蛋白表達,是獨立的預后因素(P<0.05)。通過Cox回歸模型建立列線圖預測模型;通過R軟件建立隨機生存森林模型,兩個模型均具有良好的區分度和校準度,隨機生存森林模型的臨床凈獲益優于列線圖模型。結論 建立的列線圖模型及隨機生存森林模型有助于臨床醫生判斷患者特定時間點的生存概率。

【關鍵詞】 IDH野生型膠質母細胞瘤;列線圖;隨機生存森林;預測模型;MGMT蛋白

【中圖分類號】 R739.41? 【文獻標志碼】 A? 【文章編號】 1672-7770(2024)03-0280-07

機器學習是近年來的熱點領域,被廣泛應用于生存分析中。最早的機器學習方法為生存樹的概念,類似于決策樹,生存樹也是通過樹節點選擇最佳分割,最大化兩個節點的生存差異,通常使用Log-rank檢驗作為分裂標準。隨著預測模型研究的進展,越來越多的預測模型被研發出來,在生存分析中最常用的預測模型是列線圖和隨機生存森林模型[1]。列線圖運用于生存分析中,原理是將Cox多因素分析建立的模型進行可視化繪制,通過對各因素分值的結合,確定個體發生某個臨床事件的概率。其優點為可以將連續變量以及其他影響因素進行整合,并可視化地呈現在列線圖中,便于臨床應用[2]。列線圖可以分析多種結局事件,如腫瘤的轉移、復發、死亡等,應用于生存分析中主要預測復發概率、生存概率等。列線圖雖然簡單實用,但由于其必須滿足比例風險假定,即假設某個事件發生率不隨時間變化而變化,并且由于技術的進步以及實際情況,部分事件發生率會隨著時間改變,會影響列線圖的準確性,因此基于Cox回歸分析的列線圖預測模型應用有一定限制[3]。隨機生存森林(random survival forest,RSF)是2008年首次提出的機器學習方法,是生存樹與隨機森林的結合。其特點為不要求數據滿足對數線性假定以及比例風險假定,并且其預測準確度不低于Cox比例風險模型,可用于高維數據的數據篩選。本研究通過納入2017年1月—2020年12月在空軍軍醫大學附屬西京醫院神經外科手術治療的127例異檸檬酸脫氫酶(isocitric dehydrogenase,IDH)野生型膠質母細胞瘤患者,構建生存概率的列線圖以及隨機生存森林模型,預測患者不同時間點的生存概率,并通過區分度、校準度、決策曲線評價兩種方法的適用性。

1 資料與方法

1.1 一般資料 共納入127例WHO Ⅳ級IDH野生型膠質母細胞瘤患者,其中男74例,女53例;年齡20~79歲,中位年齡55歲(四分位數間距為17);術前卡氏功能狀態評分(Karnofsky performance scale,KPS)評分≥70分的95例(74.8%)。納入標準:(1)所有患者首次入院,且由病理科專家診斷為WHO Ⅳ級IDH野生型膠質母細胞瘤;(2)術前無嚴重的心肝腎等系統異常;(3)所有患者知情同意,并簽署知情同意書;(4)所有患者均行腫瘤切除術。排除標準:(1)復發型或曾經確診為低級別膠質瘤的繼發性膠質瘤;(2)臨床資料不完整或病理資料不全;(3)妊娠或者哺乳期女性;(4)腫瘤病理標本不符合檢測要求,隨訪(聯系人)資料不全的患者;(5)術后有嚴重的手術并發癥或者嚴重的不良反應。本研究已通過西京醫院倫理委員會批準(批準號:XJLL-KY20222013)。

1.2 切除程度評價及隨訪方法 通過手術記錄,術后影像學資料來判定切除程度。通過電話隨訪的方式對患者進行出院后隨訪。總生存期定義為手術日至患者因IDH野生型膠質母細胞瘤死亡或隨訪截止時間。平均隨訪時間為18.2個月。

1.3 統計學方法 使用R4.2.1軟件及SPSS 23.0軟件對數據進行統計分析。使用R軟件的survival包建立Cox比例風險模型,多因素Cox比例風險模型采用單因素Cox分析中有意義的變量(P<0.05),并對多因素Cox比例風險模型進行比例風險(proportional hazards,PH)假定的檢驗,檢驗方法為使用R的Cox.zph函數,對Cox回歸模型的有效性進行Schoenfeld殘差的趨勢檢驗,檢驗納入的變量是否P>0.05并繪制可視化曲線。計算風險比(hazard ratio,HR)及其相應的95%置信區間(95% confidence interval,95%CI)。以P<0.05認為差異有統計學意義。

1.3.1 隨機生存森林的構建 使用R4.2.1進行統計分析。利用randomForestSRC包,以200棵樹為基礎構建隨機生存森林模型。基于使模型穩定的ntree納入Cox單因素分析有意義的變量。使用網格搜索法(grid search)計算各種組合的袋外錯誤率,將能夠達到最低袋外錯誤率的mtry和nodesize篩選出來,并構建隨機生存森林模型。再計算各變量的最小深度以及預測其重要性。

1.3.2 列線圖預測模型的構建 使用R4.2.1軟件的survival和rms包,利用內容二中構建的Cox多因素模型的β系數絕對值進行排序,然后將所有變量按照各自的β系數絕對值與β系數最大的系數絕對值相比,換算成相應分值,用rms包將各變量以其分值以平行線的形式繪制在同一個坐標系中,制作列線圖。在列線圖下方制作不同時間節點的總刻度所對應的預測生存率,可以從一定程度上預測個體在IDH野生型膠質母細胞瘤術后的生存概率。

1.3.3 預測模型的評價 使用1 000次Bootstrap重抽樣法,即自舉法進行內部驗證。通過一致性指數(concordance index,C-index)評價區分度,即將不同的患者區分開的能力,利用pec包中的cindex函數計算6個月、12個月、24個月的C指數(一致性指數)及校正C指數。rms包中的calibrate函數繪制校準曲線來評定模型的校準度。利用stdca.R程序繪制6個月、12個月、24個月的臨床決策曲線(decision curve analysis,DCA),其可以預測模型的臨床有效性,凈收益高,認為臨床效用更好,預測模型越遠離兩條極端線,其人群凈獲益更高。

2 結 果

2.1 臨床特征 腫瘤相關指標顯示,術中完全切除腫瘤有60例(47.21%),腫瘤未全切有67例(52.89%),術后接受同步放化療有72例(56.7%)、僅接受化療為9例(7.09%)、僅接受放療為3例(2.36%)、接受放化療治療但并未同步進行2例(1.57%)、沒有接受放化療41例(32.28%)。ATRX陽性患者共104例,占全部患者的81.9%;Ki-67陽性患者共118例,占全部患者的82.7%;O6-甲基鳥嘌呤甲基轉移酶(O6-methylguanine-DNA methyltransferase,MGMT)陽性患者共41例,占全部患者的32.3%;P53陽性患者共28例,占全部患者的22.0%。見表1。截至2022年8月5日因膠質瘤而死亡的為107例(84.3%),平均生存期為19.57個月,中位生存期為14.5(95%CI=12.7~16.3)個月。

2.2 IDH野生型膠質母細胞瘤Cox比例風險模型分析結果 單因素Cox回歸分析顯示年齡越低、KPS指數≥70、MGMT表達陰性、接受同步放化療為總生存期(overall survival,OS)的保護因素(P<0.05)。通過R的cox.zph函數,對Cox比例風險模型進行Schoenfeld殘差趨勢檢驗。對上述Cox單因素模型有意義的變量,通過ggcoxzph函數進行可視化處理。可知年齡、KPS評分、MGMT蛋白表達、是否接受放化療這四個因素不因為時間的改變而改變發生率,因此符合PH假定可以納入Cox多因素分析中。

將上述單因素Cox分析中有意義的因素納入Cox多因素分析中。使用R語言的coxph函數,通過逐步回歸法,以死亡為結局變量、OS為時間變量進行分析。結果顯示患者年齡、術前KPS評分、是否接受同步放化療、MGMT蛋白表達,是總生存期的獨立影響因素(P<0.05)。見表2。

2.2.1 建立隨機生存森林預測模型 將年齡、是否接受同步放化療、KPS指數以及MGMT蛋白表達作為自變量,死亡作為結局事件,當ntree取200時,模型表現趨于穩定,計算時間適中。當mtry取2、nodesize取15時,袋外錯誤率達到最低。各變量對IDH野生型膠質母細胞瘤的重要性依次為是否接受放化療、年齡、KPS指數以及MGMT蛋白表達(圖1)。

2.2.2 隨機生存森林模型預測的前10個體的生存曲線 使用R語言的matplot函數,繪制數據中前10個體的生存曲線。隨機生存森林為不可視化的預測模型,對于個體患者的預測是通過將個體患者的年齡、KPS指數、是否接受放化療、MGMT蛋白表達輸入R中進行個體患者的生存曲線繪制(圖2)。

2.2.3 隨機生存森林模型評價 隨機生存森林在模型構建過程中可以形成帶外錯誤率,即對患者預后的錯分率,本模型袋外錯誤率為37.45%。隨機生存森林模型在6個月、12個月、24個月的 Bootstrap校正C指數分別為0.724、0.696、0.700。繪制隨機生存森林模型的Bootstrap法重抽樣1 000次驗證校準曲線,表明預測模型的校準度較好,即預測生存率和實際生存率接近(圖3)。

2.3 建立列線圖預測模型 根據Cox多因素回歸分析中的獨立預后因素年齡、KPS評分、是否接受放化療、MGMT蛋白表達構建列線圖。用這四個變量在Cox多因素回歸方程中的β系數,使用R語言的rms包構建列線圖(圖4)。例如1例48歲患者、接受了放化療、入院KPS評分小于70分、MGMT蛋白表達陽性,則對該患者的各項評分相加,47+65+70+48=230,則預測該患者6個月生存概率為73%、12個月生存概率為30%、24個月生存概率小于10%。

2.4 列線圖的內部驗證 通過Bootstrap重抽樣1 000次計算出三個點的置信區間作為三條豎線,繪制得到校準曲線。模型經過1 000次重抽樣內部驗證校準后,其校準曲線接近45°角,表明預測模型的校準度較好,即預測生存率和實際生存率接近(圖5)。列線圖模型在6個月、12個月、24個月的Bootstrap重抽樣校正C指數分別為0.701、0.647、0.653,表明模型具有良好的區分能力。

2.5 決策曲線 繪制6個月、12個月、24個月的隨機生存森林模型與列線圖模型的決策曲線用以評估兩個模型的臨床凈獲益。在三個時間節點中,隨機生存森林模型與列線圖模型的曲線與兩條極端曲線較遠,因此使用兩種預測模型的臨床凈獲益率較高(圖6)。在12個月、24個月的DCA曲線中,隨機生存森林的曲線相對列線圖模型,距離兩條極端曲線更遠,因此在這兩個時間節點,隨機生存森林模型的臨床凈獲益優于列線圖模型。

3 討 論

IDH突變的WHO Ⅳ級成人彌漫性膠質瘤中位生存期為39.4個月,遠高于既往報道的WHO Ⅳ級成人彌漫性膠質瘤的平均生存時間14.6個月[4]。因此2021年新版的中樞神經系統分類將成人彌漫性膠質瘤分成IDH突變型以及IDH野生型兩類,本研究著重分析IDH野生型膠質母細胞瘤[5]。IDH野生型膠質母細胞瘤作為一種新的分類方式,其患者人口學特征、臨床特征以及免疫組化標志物的情況尚待總結分析[6]。

單因素及多因素預后分析表明,在IDH野生型膠質母細胞瘤中,術前KPS評分、年齡、是否接受同步放化療均為獨立影響因素,而MGMT蛋白表達,是分子標志物中的獨立影響因素(P<0.05)。這個結果與WHO Ⅳ級成人彌漫性膠質瘤的Cox回歸分析結果相似[7]。年齡是IDH野生型膠質母細胞瘤的一個很重要的預后因素。由于老年患者伴隨的循環系統疾病、衰老的內臟功能以及能量儲備等原因,導致老年患者并非標準治療方案的最佳適用者,因此部分老年患者無法接受手術或者放化療治療,這也與老年患者預后差有關。是否接受了同步放化療治療是IDH野生型膠質母細胞瘤一個重要預后因素。有研究表明,在年齡大于65歲的高齡人群中,聯合放化療相對于單獨放療生存收益更高[8]。手術切除程度也是一個很重要的預后指標,可以減輕占位效應、獲得病理標本,并且可以很好地減少體內腫瘤細胞數量。在WHO Ⅳ級成人彌漫性膠質瘤中,單純的腫瘤切除術可以使生存期增加約6個月[9]。KPS評分是一個與預后關系密切的指標,其評分越高,患者健康情況越好、對治療的耐受性更高、更能接受完整的標準化治療。MGMT蛋白是一種修復蛋白,當MGMT蛋白不足時會導致DNA修復過程受到影響,導致未修復的DNA受到損傷。既往研究表明,MGMT蛋白的低表達是WHO Ⅳ級成人彌漫性膠質瘤的獨立保護因素[1012],同時,本研究結果證實MGMT蛋白低表達在WHO Ⅳ級IDH野生型膠質母細胞瘤中仍然是獨立的保護因素。

機器學習作為一種新興的熱門學科,在醫學統計領域發揮了重要的作用,可以使用不同的算法提高模型的預測能力。在傳統的生存分析研究中,Kaplan-Meier法以及Cox比例風險模型是最常見的方法。隨著對預測模型的研究深入,將預測模型以及生存分析相結合,提出了各種不同的預測模型,包括支持向量機、決策樹、神經網絡等。其中列線圖與隨機生存森林是最常見的兩種方法。

列線圖可以將多種預后相關因素整合,很便捷地將研究結果可視化,在膠質瘤預后分析中有很重要的作用[13]。其通過構建Cox多因素比例風險模型,根據各個因素的回歸系數大小來評價其對結局的貢獻程度,通過對各個因素賦值,再計算其相加的總和,來計算出結局事件的發生概率,即生存概率[6]。列線圖在預測得到的生存概率與實際生存概率之間有一定的差異,但其在可視化個體預測上具有很重要的意義,因此常被用于臨床實踐中。在WHO Ⅳ級成人彌漫性膠質瘤的預測模型研究中,確定與預后相關的變量是非常重要的過程[14];通常要考慮變量的臨床與統計意義,其中最常見的臨床特征包括年齡、KPS評分、接受的治療、手術切除程度等,而其中最常見的分子標志物為IDH突變狀態以及MGMT啟動子甲基化狀態。IDH野生型膠質母細胞瘤作為一種預后很差的腫瘤,對于新診斷的患者及其家屬而言,了解其預后不同時間點的生存概率,有助于患者家屬及醫生共同選擇更合適的治療方案。列線圖由于其便捷性,可以簡單計算出6個月、12個月、24個月的生存概率,對患者及家屬直觀了解預后有很大幫助。對于患者而言,了解生存概率有助于緩解對于疾病的恐懼,也有助于制定相應的應對策略。

隨機森林是于2001年提出的基于決策樹的有監督學習方法[15], 隨機生存森林于2008年被提出,是基于隨機森林的算法[16]。隨機生存森林是隨機森林的擴展,是完全非參數模型,能評價變量間的復雜影響,無需要限制性假設,并能計算出變量的重要性。多數情況下隨機生存森林與經典Cox比例風險模型相比性能更優良。隨機生存森林是傳統的二元決策樹的拓展。隨機生存森林通過自舉法(Bootstrap)對樣本和變量進行抽樣生成大量的決策樹。對每個樣本(觀測對象)來說,所有決策樹依次對其進行預測。隨機生存森林可以處理大量輸入變量,并且可以評估變量的重要性。其造模時使用無偏估計,模型泛化能力強,當數據缺失較多時,仍可以維持一定的精度[17]。

傳統的區分能力和校準能力不能判斷使用一個模型做臨床決策是否有益,或者不同模型之間哪一個會更有臨床意義,尤其是一個模型區分度高,一個模型準確度高的情況下[18]。決策曲線可以判斷一個模型是否值得使用,當遇到假陽性或者假陰性的情況時,既然兩種情況都不能避免,那應該找出一個得到凈收益的方法,通過計算風險閾值的真陽性減去假陽性得出不同風險閾值概率的凈收益[19]。DCA可以預測模型的臨床有效性,凈收益高,被認為其臨床效用更好[20]。

本研究通過Cox比例風險模型的基礎上繪制列線圖,用來判斷具有不同年齡、KPS評分、接受放化療方式以及MGMT蛋白表達的患者,其6個月、12個月、24個月的生存概率;同時,通過R軟件構建隨機生存森林模型,用來判斷不同時間點患者的生存概率。通過Bootstrap法重抽樣1 000次驗證預測模型的區分度與校準度。隨機生存森林模型與列線圖預測模型相比,具有更好的臨床凈獲益率,并且隨機生存森林還具有因素篩選及因素重要性排序的功能。但列線圖具有操作簡便、不需要將變量輸入軟件的特點,因此也有很高的實用性。

與其他IDH野生型膠質母細胞瘤的預測模型研究相比,本研究主要納入指標為預后相關臨床因素及免疫組化檢測的蛋白表達結果。本研究發現ATRX蛋白表達、Ki-67蛋白表達、P53蛋白表達這3個在成人彌漫性膠質瘤中常用的分子標志物,在IDH野生型膠質母細胞瘤中與預后不相關。本研究還證實了MGMT蛋白是IDH野生型膠質母細胞瘤的獨立預后因素,并將其納入預測模型中。MGMT蛋白表達檢測迅速并且價格低廉,因此在臨床運用中價值很高,包含MGMT蛋白表達的預測模型更切合實際臨床工作。

本研究的局限性為單中心研究,并且IDH野生型膠質母細胞瘤是一種發病率較低的惡性腫瘤[21],因此樣本量有一定局限,這與兩個模型的C指數較低,隨機生存森林模型袋外錯誤率較高有關,并且本研究沒有探究MGMT蛋白表達與MGMT啟動子甲基化的關系,后續研究會考慮納入更多因素,納入臨床多中心相關數據使模型更加完善。本研究納入的臨床特征是常見與預后相關的特征,后續研究會納入更多臨床特征,以研究臨床特征與蛋白表達之間相關性。

利益沖突:所有作者均聲明不存在利益沖突。

[參 ?考 ??文 ??獻]

[1] Mijderwijk HJ,Nieboer D,Incekara F,et al.Development and external validation of a clinical prediction model for survival in patients with IDH wild-type glioblastoma[J].J Neurosurg,2022:110.

[2] 羅治文,陳曉,張業繁,等.機器學習算法和Cox列線圖在肝細胞癌術后生存預測中的應用價值[J].中華消化外科雜志,2020,19(2):166178.

Luo ZW,Chen X,Zhang YF,et al.Application value of machine learning algorithms and Cox nomogram in the survival prediction of hepatocellular carcinoma after resection[J].Chin J Dig Surg,2020,19(2):166178.

[3] 陳金鳳.生存分析在隨訪研究中的應用[J].實用老年醫學,2021,35(9):896899.

Chen JF.Application of survival analysis in follow-up study[J].Pract Geriatr,2021,35(9):896899.

[4] Wong QH,Li KK,Wang WW,et al.Molecular landscape of IDH-mutant primary astrocytoma Grade IV/glioblastomas[J].Mod Pathol,2021,34(7):12451260.

[5] Louis DN,Perry A,Wesseling P,et al.The 2021 WHO classification of tumors of the central nervous system:a summary[J].Neuro Oncol,2021,23(8):12311251.

[6] Alzial G,Renoult O,Paris F,et al.Wild-type isocitrate dehydrogenase under the spotlight in glioblastoma[J].Oncogene,2022,41(5):613621.

[7] Gülten G,Yaln N,Baltalarl B,et al.The importance of IDH1,ATRX and WT-1 mutations in glioblastoma[J].Pol J Pathol,2020,71(2):127137.

[8] Perry JR,Laperriere N,OCallaghan CJ,et al.Short-course radiation plus temozolomide in elderly patients with glioblastoma[J].N Engl J Med,2017,376(11):10271037.

[9] Rammeloo E,Schouten JW,Krikour K,et al.Preoperative assessment of eloquence in neurosurgery:a systematic review[J].J Neurooncol,2023,165(3):413430.

[10]Butler M,Pongor L,Su YT,et al.MGMT status as a clinical biomarker in glioblastoma[J].Trends Cancer,2020,6(5):380391.

[11]許廣智,張佳樂,伊西才,等.IDH野生型膠質母細胞瘤患者預后影響因素分析[J].臨床神經外科雜志,2022,19(2):130134.

Xu GZ,Zhang JL,Yi XC,et al.Prognostic factors of IDH wild-type glioblastoma patients[J].J Clin Neurosurg,2022,19(2):130134.

[12]Castresana JS,Meléndez B.Glioblastoma biology,genetics and possible therapies[J].Cells,2023,12(16):2063.

[13]Jekel L,Brim WR,von Reppert M,et al.Machine learning applications for differentiation of glioma from brain metastasis-a systematic review[J].Cancers,2022,14(6):1369.

[14]Zheng H,Yan TN,Han YS,et al.Nomograms for prognostic risk assessment in glioblastoma multiforme:Applications and limitations[J].Clin Genet,2022,102(5):359368.

[15]Sylman JL,Mitrugno A,Atallah M,et al.The predictive value of inflammation-related peripheral blood measurements in cancer staging and prognosis[J].Front Oncol,2018,8:78.

[16]Kim TG,Park W,Kim H,et al.Baseline neutrophil-lymphocyte ratio and platelet-lymphocyte ratio in rectal cancer patients following neoadjuvant chemoradiotherapy[J].Tumori,2019,105(5):434440.

[17]李淼.隨機生存森林在不同維度肺癌患者預后預測中的應用[D].太原:山西醫科大學,2021.

Li M.Application ofrandom survival forest in prognosisprediction of lung cancer patients with different dimensions[D].Tai Yuan:SHANXI MEDICAL UNIVERSITY,2021.

[18]Zhang Z,Jin ZP,Liu DY,et al.A nomogram predicts individual prognosis in patients with newly diagnosed glioblastoma by integrating the extent of resection of non-enhancing tumors[J].Front Oncol,2020,10:598965.

[19]Van Calster B,Wynants L,Verbeek JFM,et al.Reporting and interpreting decision curve analysis:a guide for investigators[J].Eur Urol,2018,74(6):796804.

[20]Vickers AJ,van Calster B,Steyerberg EW.A simple,step-by-step guide to interpreting decision curve analysis[J].Diagn Progn Res,2019,3:18.

[21]Ostrom QT,Cioffi G,Gittleman H,et al.CBTRUS statistical report:primary brain and other central nervous system tumors diagnosed in the United States in 2012—2016[J].Neuro Oncol,2019,21(Suppl 5):v1v100.

猜你喜歡
模型
一半模型
一種去中心化的域名服務本地化模型
適用于BDS-3 PPP的隨機模型
提煉模型 突破難點
函數模型及應用
p150Glued在帕金森病模型中的表達及分布
函數模型及應用
重要模型『一線三等角』
重尾非線性自回歸模型自加權M-估計的漸近分布
3D打印中的模型分割與打包
主站蜘蛛池模板: 五月婷婷综合网| 久久男人视频| 色综合日本| 亚洲有码在线播放| 亚洲中文字幕手机在线第一页| 国内精品伊人久久久久7777人| 园内精品自拍视频在线播放| 婷婷丁香在线观看| 亚洲人成电影在线播放| 日本不卡在线视频| 内射人妻无码色AV天堂| 欧美色图第一页| 2018日日摸夜夜添狠狠躁| 精品人妻一区二区三区蜜桃AⅤ| 色噜噜久久| 无码av免费不卡在线观看| 国产在线一二三区| www精品久久| 婷婷色一二三区波多野衣| 91精品国产91欠久久久久| 欧美在线精品一区二区三区| 成人字幕网视频在线观看| 亚洲国产精品日韩欧美一区| 欧美日韩第二页| 免费国产好深啊好涨好硬视频| 日韩成人在线网站| 欧美精品伊人久久| 欧美色视频网站| 一级毛片在线直接观看| 欧美成人日韩| 日韩大片免费观看视频播放| 久久夜夜视频| 亚洲精品无码在线播放网站| 婷婷亚洲视频| 国产精品无码一二三视频| 国产菊爆视频在线观看| 亚洲精品无码在线播放网站| 国产欧美在线观看一区| 一本大道东京热无码av| 免费一级大毛片a一观看不卡| 国产一级裸网站| 毛片网站观看| 亚洲人成亚洲精品| 亚洲专区一区二区在线观看| 8090午夜无码专区| 伊人久久综在合线亚洲91| 国产美女丝袜高潮| 久久中文字幕不卡一二区| 亚洲av色吊丝无码| 婷婷六月综合| 91精品国产自产在线老师啪l| 黄色国产在线| 国产波多野结衣中文在线播放| 国产美女91呻吟求| 久久成人免费| 欧美中出一区二区| 中文字幕日韩丝袜一区| 国产麻豆aⅴ精品无码| 免费va国产在线观看| 亚洲品质国产精品无码| 久久精品亚洲专区| 亚洲国产日韩在线观看| 九九视频在线免费观看| 国产美女免费| 一级毛片a女人刺激视频免费| 亚洲国产欧美国产综合久久| 久久不卡精品| 小说 亚洲 无码 精品| 国产欧美视频在线观看| 综合网久久| 婷婷六月激情综合一区| 免费一级毛片| 亚洲综合色吧| 国产一区二区三区免费观看| 凹凸国产熟女精品视频| 激情六月丁香婷婷四房播| 亚欧美国产综合| 日日碰狠狠添天天爽| 国产主播一区二区三区| 精品91在线| 无码粉嫩虎白一线天在线观看| 日本免费高清一区|