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

基于TCGA在宮頸癌腫瘤微環境中挖掘免疫療效相關的預后基因

2022-03-08 14:28:49吳南昌王靜盧東林韋云寶盧海娉何雪芬陳升才
右江醫學 2022年1期

吳南昌 王靜 盧東林 韋云寶 盧海娉 何雪芬 陳升才

【摘要】 目的 結合生物信息學方法,使用腫瘤基因組圖譜(TCGA)數據庫中宮頸癌表達譜數據,通過ESTIMATE算法分析腫瘤微環境挖掘有價值的預測免疫治療療效相關的預后基因。

方法 從TCGA數據庫下載宮頸癌數據,運用R軟件ESTIMATE包將宮頸癌基因表達譜數據分為高低免疫/高低基質評分組,通過limma包計算各組差異基因,最終兩組的交集基因將用于生存分析、功能富集分析及蛋白互作網絡構建。使用R Survival軟件包分析高低免疫/高低基質評分組與預后的關系,同時使用Wilcox檢驗分析高低免疫/高低基質評分組與臨床特征之間的關系。

結果 高免疫評分組與預后呈正相關(P=0.035);高基質評分組中位總生存期延長,但無統計學意義(P=0.391)。高免疫評分組、高基質評分組與M分期呈正相關(P<0.05)。設定閾值為|log2(FC)|>2且FDR<0.05,共鑒定出485個交集差異基因(上調475個,下調10個)。差異基因生存分析顯示157個基因與預后有關。基因本體分析(GO)和京都基因與基因組百科全書途徑(KEGG)共富集788個基因本體術語和36條通路(FDR<0.05)。差異基因構建PPI網絡并利用Cytoscape提取3個關鍵核心網絡,共包括66個基因,其中12個也與預后相關。

結論 基于TCGA數據庫分析腫瘤微環境獲得的免疫預后相關基因可能是潛在的免疫治療療效預測的生物標志物,值得深入研究。

【關鍵詞】 宮頸癌;腫瘤微環境;免疫基質評分;ESTIMATE算法;癌癥基因組圖譜

中圖分類號:R737.33?? 文獻標志碼:A?? DOI:10.3969/j.issn.1003-1383.2022.01.002

Exploration on prognosis genes related to immunization efficacy based on TCGA in cervical cancer tumor microenvironment

WU Nanchang, WANG Jing, LU Donglin, WEI Yunbao, LU Haiping, HE Xuefen, CHEN Shengcai

(Department of Gynecology, Affiliated Hospital of Youjiang Medical University for Nationalities, Baise 533000, Guangxi, China)

【Abstract】 Objective To combine the bioinformatics method and use the Cancer Genome Atlas (TCGA) database cervical cancer expression profile data to analyze the tumor microenvironment by using the ESTIMATE algorithm to explore valuable prognostic genes related to the prediction of immunotherapy efficacy.

Methods Cervical cancer data were downloaded from the TCGA database. The cervical cancer gene expression profile data were divided into high and low immune/high and low matrix score groups by R software ESTIMATE package. The differential genes were calculated by limma package. And the final two sets of intersection genes would be used in survival analysis, functional enrichment analysis, and protein interaction network construction. R Survival package was used to analyze the relationship between the high and low immune/high and low matrix score groups and the prognosis, and the relationship between the high and low immune/high and low matrix score groups and clinical features was analyzed by Wilcox test.

Results The high immune score group was positively correlated with prognosis (P=0.035); the median overall survival was prolonged in the high matrix score group, but difference was not statistically significant (P=0.391). The high immune score group and the high matrix score group were positively correlated with M stage (P<0.05). The threshold was set to |log2(FC)|>2 and FDR<0.05. A total of 485 intersection difference genes were identified (475 up-regulation and 10 down-regulation). Differential gene survival analysis showed that 157 genes were associated with prognosis.The gene ontology (GO) analysis and the Kyoto Encyclopedia of Genes and Genome Pathways (KEGG) enriched 788 gene ontology terms and 36 pathways (FDR<0.05). The differential gene constructed PPI network and used Cytoscape to extract 3 key core networks, including 66 genes, 12 of which were also associated with prognosis.

Conclusion Immune prognosis-related genes obtained from tumor microenvironment based on TCGA database may be a potential biomarker for predicting the efficacy of immunotherapy, which is worthy of further study.

【Key words】 cervical cancer; tumor microenvironment; immune matrix score; ESTIMATE algorithm; TCGA

宮頸癌是威脅我國女性健康最常見的婦科惡性腫瘤之一。全球每年有超過50萬新發的宮頸癌患者,約30萬例死于宮頸癌,其發病率及病死率均居女性惡性腫瘤第4位[1]。手術為主、放化療為輔的傳統治療對早期宮頸癌治療效果不錯,但對于晚期、轉移和復發宮頸癌5年生存率低于20%[2~3]。近年來,免疫治療是腫瘤治療的新手段,其機制是通過重新激活自身的抗腫瘤免疫系統對腫瘤細胞強且持續性的殺傷,顯著提高生存時間。免疫治療目前研究最全面的是免疫檢查點抑制劑:其代表藥為程序性死亡蛋白1(PD-1)抑制劑(帕姆單抗),在多種癌癥中被證明有效,但總體客觀有效率僅為20%~30%[4]。目前用于指導免疫治療的分子靶標主要局限于程序性死亡蛋白配體1(PD-L1)表達水平、高度微衛星不穩定(MSI-H)[5]、錯配修復系統缺陷(dMMR)[6]、腫瘤突變負荷(TMB)[7~8],雖然取得一些進展,但仍存在很大爭議,因此,通過全面分析腫瘤微環境篩選更合理的分子標志物以指導免疫治療顯得尤為重要。腫瘤的發生發展是致癌因素持續存在刺激基因突變與腫瘤細胞經免疫逃逸后形成微環境的共同結局[9]。腫瘤微環境是指腫瘤細胞與腫瘤局部浸潤的免疫細胞、基質細胞等共同構成的局部病理環境。既往研究提示,腫瘤微環境中高免疫、基質細胞浸潤水平往往伴隨良好預后。最近,YOSHIHARA等[10]開發了一種通過分析腫瘤內部免疫和基質細胞的特定基因表達特征來計算免疫和基質評分進而預測腫瘤純度的算法。該算法在前列腺癌[11]、乳腺癌[12]和結腸癌[13]中得到應用并被證實有效,然而在宮頸癌領域尚未見報道。本研究利用腫瘤基因組圖譜(TCGA)數據庫中宮頸癌隊列和ESTIMATE算法衍生的免疫評分[10],探索了與宮頸癌腫瘤微環境免疫預后相關的基因,為挖掘有價值的預測免疫治療療效生物標志物提供新思路。

1 材料與方法

1.1 基因表達數據集

宮頸癌患者的3級基因表達譜數據來源于TCGA數據門戶網站(https://tcga-data.nci.nih.gov/tcga/)。另外,我們還下載了包括年齡、組織學類型、生存時間和分期等臨床特征信息。下載的數據使用ESTIMATE算法來計算免疫評分和基質評分[10]。本研究嚴格遵循TCGA批準的出版指南,因此不需要倫理委員會批準。

1.2 差異表達基因分析

所有宮頸癌樣本根據中位值分為高、低免疫評分組和高、低基質評分組后分別進行差異表達計算。將Log2|FC|>2且FDR<0.05設定為選擇閾值,通過limmar包篩選顯著差異基因,并使用R軟件的pheatmap包生成熱圖。使用Venn包選擇出高、低免疫評分組和高、低基質評分組交集的基因。

1.3 生存分析

為表明免疫/基質評分以及差異基因與患者預后之間的關系,根據Kaplan-Meier曲線法在R中運行Survival軟件包繪制Kaplan-Meier圖。

1.4 功能分析

為了解差異基因的潛在功能,免疫基質組交叉DEG的功能富集分析采用基因本體分析(GO)和京都基因與基因組百科全書途徑(KEGG)。FDR<0.05作為截止值。

1.5 免疫、基質評分與臨床特征的關系

使用Wilcox檢驗分析免疫、基質評分與臨床特征之間的關系,包括T期[(T3 + T4) vs (Tis + T1 + T2)]、N期(N1 vs N0)、M期(M1 vs M0)、FIGO臨床期(ⅠA~ⅡA2 vs ⅡB~ⅣB)、年齡、淋巴血管侵犯(YES vs NO)和腫瘤組織學分級(G1~2 vs G3~4)。P<0.05為差異有統計學意義。

1.6 PPI網絡的構建

進一步尋找基因之間的相互作用,利用STRING數據庫[14]構建蛋白質-蛋白質相互作用(PPI)網絡,并借助Cytoscape[15]平臺重建,只有包含超過10個節點的單個網絡用于進一步分析。計算網絡的每個節點的連接程度,然后使用MCODE插件基于拓撲找到簇以定位密集連接的區域。

2 結? 果

2.1 差異表達的免疫相關基因分析

將TCGA數據庫307例宮頸癌樣本的表達譜數據分成高、低免疫評分組和高、低基質評分組,并比較選擇差異表達的免疫相關基因。高免疫評分組有735個上調基因,412個下調基因;高基質評分組有1156個上調基因和20個下調基因(圖1A、1B)。此外,為提高結果的可靠程度,韋恩圖顯示,免疫與基質變化方向一致的交集差異基因共485個(475個上調,10個下調)。見圖1C、1D。

2.2 免疫基質評分與臨床特征的關系

當設定P<0.05時,結果提示免疫、基質評分與M分期有統計學意義。低免疫、低基質腫瘤細胞與轉移相關,高免疫、高基質腫瘤細胞不易發生轉移(圖2A、2B)。

2.3 免疫和基質評分與宮頸癌預后的關系

為了解免疫/基質評分是否與患者預后相關,根據307例宮頸癌患者的免疫/基質評分分別納入高分組和低分組并構建Kaplan-Meier生存曲線。結果提示,高免疫評分組病例的中位總生存期(overall survival,OS)高于低免疫評分組,免疫評分與OS呈正相關(P=0.035);高基質評分組較低基質評分組中位生存期更長,但差異無統計學意義(P=0.391)。見圖2C、2D。

2.4 交叉基因生存分析

485個交叉基因在通過R的Survival進行Kaplan-Meier曲線分析,結果157個基因與總生存率相關(P<0.05),其中參與構建核心網絡的有12個(P<0.05,圖3)被認為是潛在的免疫相關的預后基因。

2.5 功能富集分析

GO通路和KEGG共富集788個GO術語和36條通路(FDR<0.05,如補充S1、S2)。結果顯示前5條KEGG通路為:細胞因子-細胞因子-受體相互作用、趨化因子信號通路、造血細胞系、原發性免疫缺陷、細胞黏附分子(圖4A)。前五個GO術語為:T細胞激活、淋巴細胞活化的調節、白細胞-細胞黏附、T細胞活化的調節、白細胞-細胞黏附的調節(圖4B)。

2.6 PPI核心網絡構建

將485個基因導入String并結合Cytoscape構建PPI網絡,揭示不同基因之間可能存在功能上的互作,包括239個節點和1267個邊緣。包含15個節點以上的三個模塊核心網絡被選出進一步分析(圖5A~C)。模塊A(圖5A)由33個節點及528個邊緣構成,C3AR1、GNG8、GNGT2和FPR2與其他基因緊密關聯。模塊B(圖5B)具有18個節點和141個邊緣,ITGB2、FCER1G和GPR84具有更高的連通度值。模塊C有15個節點和70個邊緣(圖5C),重要節點包括CD3E、CD3D和CD3G。3個模塊共包含66個基因,其中12個與預后密切相關,提示可以作為免疫相關的預后基因。

3 討? 論

除貝伐單抗聯合化療一線治療外,缺乏有效的二線治療方案是晚期、復發及轉移性宮頸癌患者死亡的主要原因。目前,新興的免疫治療在多種實體瘤中取得良好療效,尤其是免疫檢查點(PD-1/PD-L1)抑制劑,包括宮頸癌。KEYNOTE 028(第一階段研究)[16]和KEYNOTE 158(第二階段研究)[17]分別以10 mg/(kg·2W)和200 mg/3W劑量的帕姆單抗評估了復發性、轉移性宮頸癌,結果提示帕姆單抗具有持久的抗腫瘤活性和可控制的安全性,FDA于2018年加速批準帕姆單抗用于化療期間或之后疾病進展且PD-L1陽性的晚期宮頸癌患者。隨后NCCN婦瘤指南指出MSI-H/dMMR/PD-L1陽性的晚期/復發和耐藥宮頸癌患者二線用藥推薦使用帕姆單抗。雖然免疫治療在個別晚期宮頸癌收到奇效,但總體有效率僅為10%,尋找更加合適的標志物篩選受益人群是急需解決的難題。因此,我們試圖通過研究宮頸癌腫瘤微環境鑒定可以預測宮頸癌免疫治療療效的分子靶標。本研究下載TCGA數據庫宮頸癌307例樣本的轉錄數據,經ESTIMATE算法計算,人為將它們分為高、低免疫評分組和高、低基質評分組進行差異表達分析。我們的研究結果再次證明了腫瘤微環境免疫及基質占比高低與患者預后相關,并且最終鑒定了485個基因與宮頸癌腫瘤微環境免疫顯著相關。

近年,在直腸癌中建立的免疫評分研究提示,免疫評分高,生存時間長,預后較好,可作為預后的獨立影響因素[18]。本研究分析每個樣本免疫和基質的評分分別與個體預后及臨床病理特征的關系,發現腫瘤微環境高免疫、高基質評分與良好預后及腫瘤病理M分期呈正相關,然而基質細胞與預后無統計學意義。該結果與既往結論基本一致,進一步說明腫瘤微環境對腫瘤的惡性增殖、免疫逃逸及逃避凋亡等有促進作用。

將樣本分為高、低免疫評分和高、低基質評分的差異基因交集結果作為最終的免疫相關基因,通過GO、KEGG分析,結果顯示大部分基因參與腫瘤微環境的眾多關鍵環節(圖4),如GO中的T細胞激活、淋巴細胞活化的調節、T細胞活化的調節及KEGG中細胞因子-細胞因子-受體相互作用、趨化因子信號通路、細胞黏附分子。這與先前的報道一致,即免疫細胞與細胞外基質細胞的分子功能在構建宮頸癌腫瘤微環境中有重要作用[19]。

為了進一步確定預后免疫相關基因,我們對485個差異基因進行總體生存分析,其中168個與宮頸癌患者的預后有關。此外,我們選取了包含節點大于15個基因的3個PPI蛋白互作網絡進一步分析,發現網絡中12個核心基因同時具備預后特征,因此,我們對它們特別感興趣。6個核心基因已經被報道(CCL5、CXCL9、CD3D、CD3E、CD8A、IL2RB)參與腫瘤微環境的構建,甚至影響免疫治療的效果,這表明我們基于TCGA數據庫分析腫瘤微環境相關免疫基因具有重大意義。其余6個基因在腫瘤微環境領域尚未見報道,包括CD3G、CLEC4C、GNG8、GNGT2、P2RX1、P2RY13,可作為宮頸癌免疫治療潛在的療效評價生物標志物。研究表明在b16-f10腫瘤中靶向大自噬/自噬基因becn1/beclin1激活mapk8/jnk-jun/c-jun信號通路并釋放CCL5細胞因子誘導NK細胞浸潤抑制腫瘤生長[20]。ALLEN等[21]在卵巢癌中發現Th17細胞可通過IL-17與干擾素γ的協同作用,刺激CXCL9和CXCL10趨化因子產生,使腫瘤浸潤效應T細胞向腫瘤微環境中募集,進而提高抗腫瘤免疫能力。SHI等[22]通過分析4個肌肉浸潤性膀胱癌(MIBC)公共數據集,發現CD3D/CD4比值在MIBC中是一個穩定的獨立預后因子,高CD3D/CD4比值有較好的生存率,CD4表達明顯高于CD3D提示細胞具有免疫抑制作用,可用于免疫治療療效評估的潛在指標。CHEN等[23]基于腫瘤基因組圖譜數據庫對軟組織肉瘤(STS)分析,發現IL-33和ST2的表達水平與CD3E、CD8A和趨化因子的表達水平呈正相關,提示IL-33/ST2軸可能在1型極化CD8+T細胞的免疫應答中起重要作用,提示這些分子可作為免疫治療的預后標志物。ZHU等[24]分析了92例ICC組織中PD-L1表達與CD8+T細胞浸潤的關系,發現腫瘤IFN-γ與PD-L1和CD8A基因表達顯著相關,此外,用CD8+T細胞分泌的重組IFN-γ刺激HCCC9810和RBE細胞可在體外快速誘導PD-L1的上調,表明高表達CD8A可能是選擇抗PD-L1/PD-1治療的潛在標志物。人類IL-2rβ鏈(IL2RB)基因缺陷是致命的免疫失調的原因,CAMPBELL等[25]研究采用野生型IL2RB干細胞移植改善了1例患者T淋巴細胞的IL-2反應性,為使用IL-2療法治療免疫性疾病和癌癥提供了新觀點。

總之,本研究基于TCGA數據庫的宮頸癌表達譜數據通過ESTIMATE算法得到高低免疫、基質評分組,將兩者的交集差異基因進行生存分析、GO及KEGG富集,同時利用String構建PPI網絡,最終提取了一份宮頸癌腫瘤微環境免疫相關的預后基因。腫瘤微環境由多種因素相互作用,極其復雜,目前對于免疫治療療效預測和評估的生物標志物有限,準確篩選獲益人群仍是難題,我們的研究通過分析宮頸癌TME中主要成分免疫、基質細胞免疫評分進而獲得免疫相關的基因列表,它們可能是先前被忽略的潛在有力的免疫治療療效預測及預后評估的生物標志物,值得進一步全面深入研究。

參 考 文 獻

[1] ?BRAY F,FERLAY J,SOERJOMATARAM I,et al.Global cancer statistics 2018:GLOBOCAN estimates of incidence and mortality worldwide for 36 cancers in 185 countries[J].CA:A Cancer J Clin,2018,68(6):394-424.

[2] ?PFAENDLER K S,TEWARI K S.Changing paradigms in the systemic treatment of advanced cervical cancer[J].Am J Obstet Gynecol,2016,214(1):22-30.

[3] ?TEWARI K S,SILL M W,PENSON R T,et al.Bevacizumab for advanced cervical cancer:final overall survival and adverse event analysis of a randomised,controlled,open-label,phase 3 trial (Gynecologic Oncology Group 240)[J].Lancet,2017,390(10103):1654-1663.

[4] ?IWAI Y,HAMANISHI J,CHAMOTO K,et al.Cancer immunotherapies targeting the PD-1 signaling pathway[J].J Biomed Sci,2017,24(1):26.

[5] ?LE D T,URAM J N,WANG H,et al.PD-1 blockade in tumors with mismatch-repair deficiency[J].N Engl J Med,2015,372(26):2509-2520.

[6] ?LE D T,DURHAM J N,SMITH K N,et al.Mismatch repair deficiency predicts response of solid tumors to PD-1 blockade[J].Science,2017,357(6349):409-413.

[7] ?GOODMAN A M,KATO S,BAZHENOVA L,et al.Tumor mutational burden as an independent predictor of response to immunotherapy in diverse cancers[J].Mol Cancer Ther,2017,16(11):2598-2608.

[8] ?YARCHOAN M,HOPKINS A,JAFFEE E M.Tumor mutational burden and response rate to PD-1 inhibition[J].N Engl J Med,2017,377(25):2500-2501.

[9] ?FIDLER I J.The pathogenesis of cancer metastasis:the 'seed and soil' hypothesis revisited[J].Nat Rev Cancer,2003,3(6):453-458.

[10] ?YOSHIHARA K,SHAHMORADGOLI M,MARTNEZ E,et al.Inferring tumour purity and stromal and immune cell admixture from expression data[J].Nat Commun,2013,4:2612.

[11] ?SHAH N, WANG P, WONGVIPAT J, et al.Regulation of the glucocorticoid receptor via a BET-dependent enhancer drives antiandrogen resistance in prostate cancer[J].Elife, 2017,6:e27861.

[12] ?PRIEDIGKEIT N,WATTERS R J,LUCAS P C,et al.Exome-capture RNA-sequencing of decade-old breast cancers and matched decalcified bone metastases identifies clinically actionable targets[J].bioRxiv,2017.

[13] ?ALONSO M H,AUSS S,LOPEZ-DORIGA A,et al.Comprehensive analysis of copy number aberrations in microsatellite stable colon cancer in view of stromal component[J].Br J Cancer,2017,117(3):421-431.

[14] ?SZKLARCZYK D,FRANCESCHINI A,WYDER S,et al.STRING v10:protein-protein interaction networks,integrated over the tree of life[J].Nucleic Acids Res,2015,43(database issue):D447-D452.

[15] ?SHANNON P,MARKIEL A,OZIER O,et al.Cytoscape:a software environment for integrated models of biomolecular interaction networks[J].Genome Res,2003,13(11):2498-2504.

[16] ?OTT P A,BANG Y J,PIHA-PAUL S A,et al.T-cell-inflamed gene-expression profile,programmed death ligand 1 expression,and tumor mutational burden predict efficacy in patients treated with pembrolizumab across 20 cancers:KEYNOTE-028[J].J Clin Oncol,2019,37(4):318-327.

[17] ?CHUNG H C,ROS W,DELORD J P,et al.Efficacy and safety of pembrolizumab in previously treated advanced cervical cancer:results from the phase II KEYNOTE-158 study[J].J Clin Oncol,2019,37(17):1470-1478.

[18] ?MLECNIK B,BINDEA G,ANGELL H K,et al.Integrative analyses of colorectal cancer show immunoscore is a stronger predictor of patient survival than microsatellite instability[J].Immunity,2016,44(3):698-711.

[19] ?DENKO N,SCHINDLER C,KOONG A,et al.Epigenetic regulation of gene expression in cervical cancer cells by the tumor microenvironment[J].Clin Cancer Res,2000,6(2):480-487.

[20] ?NOMAN M Z,BERCHEM G,JANJI B.Targeting autophagy blocks melanoma growth by bringing natural killer cells to the tumor battlefield[J].Autophagy,2018,14(4):730-732.

[21] ?ALLEN F,BOBANGA I D,RAUHE P,et al.CCL3 augments tumor rejection and enhances CD8+ T cell infiltration through NK and CD103+ dendritic cell recruitment via IFNγ[J].Oncoimmunology,2018,7(3):e1393598.

[22] ?SHI M J,MENG X Y,WU Q J,et al.High CD3D/CD4 ratio predicts better survival in muscle-invasive bladder cancer[J].Cancer Manag Res,2019,11:2987-2995.

[23] ?CHEN H,CHEN Y,LIU H,et al.Integrated expression profiles analysis reveals correlations between the IL-33/ST2 axis and CD8+T cells,regulatory T cells,and myeloid-derived suppressor cells in soft tissue sarcoma[J].Front Immunol,2018,9:1179.

[24] ?ZHU Y,WANG X Y,ZHANG Y,et al.Programmed death ligand 1 expression in human intrahepatic cholangiocarcinoma and its association with prognosis and CD8+ T-cell immune responses[J].Cancer Manag Res,2018,10:4113-4123.

[25] ?CAMPBELL T M,BRYCESON Y T.IL2RB maintains immune harmony[J].J Exp Med,2019,216(6):1231-1233.

(收稿日期:2021-08-03 修回日期:2021-11-16)

(編輯:潘明志)

主站蜘蛛池模板: 亚洲嫩模喷白浆| 亚洲国产欧美国产综合久久| 国产大全韩国亚洲一区二区三区| 国产极品嫩模在线观看91| 狠狠操夜夜爽| aaa国产一级毛片| 911亚洲精品| www亚洲天堂| 国产人人乐人人爱| 日本午夜精品一本在线观看| 国产亚洲精品在天天在线麻豆| 亚洲福利视频网址| 欧美视频在线不卡| h网址在线观看| 日本黄网在线观看| 精品无码日韩国产不卡av| 少妇被粗大的猛烈进出免费视频| 91精品国产福利| 人妻一区二区三区无码精品一区| 欧美色视频日本| 欧美福利在线观看| 91色综合综合热五月激情| 国产综合亚洲欧洲区精品无码| 国产精品亚洲专区一区| 2021国产精品自产拍在线观看| 国产精品亚洲一区二区三区在线观看 | 国产精品乱偷免费视频| 92午夜福利影院一区二区三区| 国产毛片基地| 欧美无专区| 久青草免费在线视频| 亚洲视频a| 国产欧美日韩一区二区视频在线| 在线观看av永久| 国产精品不卡片视频免费观看| 亚洲欧美激情另类| 国产人人射| 午夜高清国产拍精品| 夜夜操国产| 国产丰满大乳无码免费播放 | 人妻一本久道久久综合久久鬼色| 国产啪在线| 一本大道东京热无码av| 国产综合精品一区二区| 免费中文字幕在在线不卡| 中文字幕佐山爱一区二区免费| 毛片网站观看| 99手机在线视频| 国产在线自乱拍播放| 国产96在线 | 99re在线视频观看| 国产精品尹人在线观看| 国产精品爽爽va在线无码观看| 亚洲香蕉在线| 欧美a√在线| 亚洲人成网站色7777| 四虎精品黑人视频| 国产又粗又爽视频| 乱人伦视频中文字幕在线| 欧美成在线视频| 亚洲AⅤ波多系列中文字幕| 色婷婷亚洲综合五月| 欧美成人看片一区二区三区| 乱人伦中文视频在线观看免费| 亚洲无码A视频在线| 欧美一级高清免费a| 国产微拍一区二区三区四区| 三上悠亚一区二区| 谁有在线观看日韩亚洲最新视频| 亚洲v日韩v欧美在线观看| 91蜜芽尤物福利在线观看| 狠狠色丁香婷婷| 日韩人妻精品一区| 欧美精品亚洲精品日韩专| 中国美女**毛片录像在线 | 国产在线欧美| 午夜爽爽视频| 久爱午夜精品免费视频| 日韩AV无码免费一二三区| 亚洲一区二区三区国产精华液| 91久久精品日日躁夜夜躁欧美| 香蕉久久国产超碰青草|