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

基于雙通道卷積神經網絡的航班延誤預測模型

2018-08-27 10:58:18吳仁彪李佳怡屈景怡
計算機應用 2018年7期
關鍵詞:特征提取模型

吳仁彪,李佳怡,屈景怡

(中國民航大學 天津市智能信號與圖像處理重點實驗室,天津 300300)(*通信作者電子郵箱qujingyicauc@163.com)

0 引言

隨著航空業的迅速發展,航班延誤問題日益嚴重。由航班延誤引發的旅客與機場、航空公司之間的矛盾已經成為一種嚴重的社會問題,引起極大反響。航班延誤預測可通過已知信息進行到港延誤等級預測,使相關部門提前作好應對措施,提升旅客滿意水平。同時,通過預測可提前知曉航班到達序列,可進一步改進航班調度和機位指派,對航班延誤的治理具有重要意義[1-2]。

國內外相關學者針對航班延誤預測問題已展開大量研究[3-6],已有方法中多采用決策樹、貝葉斯網絡、神經網絡等機器學習算法[7-10]。程華等[8]研究基于決策樹方法的到港延誤預測問題,通過實驗發現使用C4.5決策樹進行延誤等級預測可得到80%的正確率;曹衛東等[9]考慮同一架飛機的連續航班情況,利用貝葉斯網絡進行航班延誤波及分析,可反映連續航班延誤原因分布及波及分布;Khanmohammadi等[10]采用二進制編碼航班特征后輸入神經網絡,在實現延誤預測的同時可獲得延誤原因。以上模型較少考慮氣象信息對航班延誤的影響,而且傳統機器學習算法仍是針對小數據集的信息獲取和挖掘,處理能力有限。為使正確率進一步提高則需加入專家知識或通過聚類算法等方式擴展信息[11],因此,有學者提出將大數據技術應用于航班延誤預測[12-13],可使模型自動進行特征提取,無需人為干預,并可以提升數據處理能力。

Kim等[13]率先采用大數據方法預測機場延誤狀態,使用長短時記憶法(Long Short Term Memory, LSTM)考慮時間維度特征,并通過實驗證明了模型的可擴展性;但全連接神經網絡進行特征提取最多只能訓練五層網絡模型,沒有克服深度加深時模型無法收斂的缺點。深度卷積神經網絡(Convolutional Neural Network, CNN)通過權值共享,可減少訓練參數,訓練更深的模型,從而最大限度從海量數據中提取信息。近幾年,在圖像分類與識別、語音識別等領域已取得突破性進展[14-16]。其中,Highway網絡[17]和殘差網絡(Residual Network, ResNet)[18]通過對卷積神經網絡傳統結構的改進,進一步增加了網絡的深度,可取得更好的特征提取效果。

本文分數據預處理、特征提取和分類預測三部分構建航班延誤預測模型。在數據預處理部分考慮到氣象因素對航班延誤的重要影響,進行了航班數據和氣象數據的融合。在特征提取部分考慮到傳統機器學習算法處理大數據的困難,以及網絡深度對模型性能的影響,參考ResNet網絡結構,提出了一種雙通道卷積神經網絡(Dual-Channel Convolutional Neural Network, DCNN)結構。該網絡結構可克服梯度消失的缺陷,訓練深度網絡,增強數據處理性能,從而進一步提升模型預測能力。

1 雙通道卷積神經網絡

1.1 傳統卷積神經網絡

傳統卷積神經網絡主要由卷積層和池化層構成,如要訓練深層網絡則需進行如圖1(a)所示的卷積層堆疊。每個卷積層后接激活函數進行非線性變換,如式(1)所示:

A=f(W?X)

(1)

其中:A表示輸出特征矩陣;W表示需要學習的卷積核;X表示卷積層的輸入特征矩陣;?表示Hadamard乘積運算,即矩陣對應元素相乘求和,是一種相關運算;f代表使用激活函數進行非線性變換。本文選用線性整流單元(Rectified Linear Unit,ReLU)[19]作為激活函數,函數定義如下:

f(x)=max(x,0)

(2)

依據式(1)~(2),圖1(a)的非線性映射函數可表示為:

Y=f(W(L)?(…f(W(1)?X)…))

(3)

其中L表示卷積層堆疊層數。

1.2 雙通道卷積神經網絡

雙通道卷積神經網絡結構如圖1(b)所示,由多個Block堆疊而成,每個Block中包含直通通道和卷積通道。

直通通道是對數據的簡單線性映射變換,以保證數據的無損傳輸,其函數定義為:

Y=X

(4)

參考ResNet網絡,DCNN同樣進行了特征矩陣的累加;但不同于ResNet網絡中的一次累加[18],DCNN每個Block進行了兩次累加,即每個Block中包含兩個卷積模塊,且分別和直通通道進行一次累加。其目的在于融合不同網絡深度提取的特征,可同時考慮底層至高層的數據特征。出于保持兩個通道輸出數據完整性的考慮,每次累加前后不接激活函數。則卷積通道的輸出Yn1、Yn2可由式(5)進行計算:

(5)

其中:Yn1和Yn2分別代表兩個卷積模塊的輸出特征矩陣,W(1)、W(2)、W(3)和W(4)分別代表4個卷積層的卷積權值矩陣。更一般的情況,每個Block可包含m(m≥2)個卷積模塊,每個卷積模塊可包含l(l≥2)層卷積層,故卷積模塊的非線性映射函數可表示為:

H(X,{W(l)})=W(l)?(…f(W(2)?f(W(1)?X))…)

(6)

為保證信息在網絡中傳輸的暢通性,引入卷積衰減因子對卷積通道輸出的特征矩陣進行稀疏性限制。卷積衰減因子的作用類似于Highway網絡中的控制門[17],但Highway網絡中控制門的參數需要進行大量計算,且當網絡深度加深時參數學習變得困難。本文模型直接引入卷積衰減因子,無需增加格外計算,同時也起到了控制門的作用。出于以直通通道為主、卷積通道為輔的原則,將直通通道輸出權值設置為1;同時出于以淺層特征為主、深度特征用于微調的目的,另卷積衰減因子0≤λ2≤λ1≤1,增大淺層網絡貢獻率,維持網絡的穩定性。

圖1 網絡結構

由圖1(b)可看出,第n(n≥2)個Block的映射函數可表示為:

Yn=λ1Yn1+λ2Yn2+Yn-1

(7)

其中λ1和λ2分別表示兩個卷積模塊的卷積衰減因子。

每個Block的輸出為下一個Block的輸入,若以第n個Block的輸出作為整個網絡的結果,則由級聯的Block構成的DCNN的非線性映射函數可表示為:

Y=λ1Yn1+λ2Yn2+λ1Y(n-1)1+λ2Y(n-1)2+…+

λ1Y11+λ2Y12+X

(8)

由式(8)可看出,DCNN可等效于對不同深度網絡的“投票”,通過自身的學習選擇不同深度網絡的特征組合,且賦予相應的權值,控制不同輸出特征的比例。

1.3 誤差項計算

在卷積神經網絡訓練時,首先進行前向傳播計算各層輸出特征矩陣并計算損失函數J,然后依據BP(Back Propagation)鏈式法則反向計算各層誤差項δ(l),并根據誤差項計算權值梯度,計算公式如下:

(9)

最后,按照式(10)和式(11)更新網絡權值參數,反復訓練直至更新后的權值參數使網絡收斂:

(10)

W(t+1):=W(t)+V(t+1)

(11)

其中:t表示迭代次數;μ為動量因子(momentum),決定了歷史權重修正量對本次權重修正的貢獻大小;V(t)是動量項(momentum variable);η為學習率(learning rate);λ為權重衰減系數(weight decay)。

在計算各層誤差項時,傳統卷積神經網絡按式(12)反向計算各隱藏層誤差項:

(12)

由于ReLU激活函數的導數值為1,權值矩陣W中元素的絕對值普遍小于1,具有較強的稀疏性。因而在式(12)中,δ(1)?δ(L),且L越大,傳至第一層時所乘小于1的因子越多,δ(1)越小。此時,按式(9)計算權值梯度值越小。當網絡不斷加深,就容易出現梯度消失現象。

雙通道卷積神經網絡由于加入直通通道,可通過對誤差項反向傳播時的調控,弱化梯度消失現象。取m=2,l=2,其誤差項計算過程如下:

(13)

其中:δ(1)、δ(2)、δ(3)和δ(4)代表各卷積層誤差項,δ(n)和δ(n+1)代表相鄰兩個Block的誤差項,δ(n2)和δ(n1)分別代表第n個Block中兩個卷積模塊的誤差項。

若令式(12)中L=4,即CNN和DCNN具有相同的深度,對比式(12)、(13)可知:

δn≥δn+1≥δ(L)?δ(1)

(14)

由此可見,相比CNN,DCNN反向傳播時誤差項衰減變緩,從而保證了梯度的穩定性,使深度網絡的權值參數得到更加充分的訓練和學習。

2 模型描述

航班延誤預測模型如圖2所示,主要由數據預處理、特征提取和分類預測三部分組成。特征提取使用深層DCNN完成,以訓練更深層網絡,提煉更多有效信息。現對數據預處理和分類預測部分進行描述。

圖2 模型描述

2.1 數據預處理

為方便后文敘述,現對數據定義如下。

定義1 航班數據F中包含〈C,Sd,Sa,Ad,Aa,tsd,tad,tsa,taa〉。其中:C代表航空器型號;Sd代表起飛機場所在州;Sa代表目的機場所在州;Ad代表起飛機場;Aa代表目的機場;tsd代表計劃起飛時間;tad代表實際起飛時間;tsa代表計劃到達時間;taa代表實際到達時間。

定義2 氣象數據O中包含〈As,ts,T,H,Wd,Ws,P,S,V〉。其中:As代表觀測站所在機場;ts代表觀測時間;T代表溫度;H代表濕度;Wd代表風向;Ws代表風速;P代表氣壓;S代表天空狀況;V代表能見度。

定義3 航班到港延誤時間Dt定義為實際到達時間和計劃到達時間的差值,另f為單條航班數據,則:

Dt=f.taa-f.tsa

(15)

數據預處理過程主要完成數據編碼、數據融合和特征矩陣的構造。由于航班數據F和氣象數據O分別具有不同的量綱和標識方式,因此對名義變量進行編碼。編碼采用頻次統計的方式,每一變量中出現頻次最高的變量值編碼為0,其次為1,依此類推。這種編碼方式可保證特征矩陣的稀疏性,增強魯棒性。

由于本文只關注延誤數據,因此對數據集中的取消和改航數據進行過濾處理,同時依據Sd和Sa進行各州州時轉換,保證全部時間的時區一致性。為探討目的機場氣象數據對航班延誤的影響,通過Spark并行化平臺[20]實現航班數據和氣象數據的融合,輸出總的航班延誤數據。

數據融合流程如下:

步驟1 輸入航班數據F和氣象數據O,使兩份數據均以彈性數據集的形式均勻分布于集群的各個節點上。

步驟2 建立Map1階段。將單條航班數據f從F中提取,設定目的機場Aa和計劃到達時間tsa為關聯主鍵join_key1,輸出數據〈join_key1,f〉。

步驟3 建立Map2階段。將單條氣象數據o從O中提取,設定觀測站所在機場As和觀測時間ts為關聯主鍵join_key2,輸出數據〈join_key2,o〉。

步驟4 進行數據的重新緩存。將各個彈性分布式數據集中的〈join_key1,f〉和〈join_key2,o〉進行重新分區,具有相同關聯主鍵的數據緩存于相同節點的數據集中。

步驟5 建立Reduce1階段。將緩存于同一分區的f和o進行合并構建〈f,o〉,即相同主鍵數據進行銜接,輸出融合后數據w。

步驟6 建立Reduce2階段。將輸出的w從各節點中獲取,構建融合后的輸出矩陣W。

2.2 分類預測

將航班延誤時間Dt按表1進行延誤等級劃分,從而將預測問題轉化為分類問題。根據文獻[21]的要求,將延誤時間為15 min以內的航班設定為不延誤航班,其他等級參考文獻[8]進行劃分。

分類過程中使用全局平均池化進行數據降維,使用全連接層將最后輸出的特征矩陣轉換為(1×1)×q維,q為分類類別數。最后,采用softmax進行分類。每一個類別j的輸出概率值為p(y=j|x),則函數為:

(16)

(17)

表1 航班延誤等級劃分

3 模型優化

3.1 批歸一化

批歸一化(Batch Normalization, BN)技術[22]一定程度上可提高網絡的訓練速度。其基本思想是在利用每一批樣本對網絡進行訓練時,分別對各層網絡的輸出數據進行歸一化處理。

假設某一層網絡輸入數據的集合為{x1,x2,…,xk},該批樣本的數量為k,則BN技術首先對該批次數據按式(18)進行歸一化處理:

(18)

其中:E[x]表示該批次輸入的樣本均值;Var[x]代表該批次輸入的樣本方差;ε為一個加到方差上的數值常量,一般取一個接近于零的正數,用來保證重構變換響應數值的穩定性。

歸一化后的數據符合均值為0、方差為1的分布,大部分落在激活函數的非飽和區域,防止了梯度消失,加速訓練;但進行歸一化后不可避免地破壞了上一層網絡的原始表達,因此引入重構參數γ和β對批歸一化后的數據按重構函數(19)進行變換,得到最終該層輸出數據{y1,y2,…,yk}:

(19)

利用BN技術對網絡內部數據進行歸一化處理,等效于在網絡當中引入一個BN網絡層。按文獻[22]的實驗結果,BN層加在激活函數之前效果更佳,因此在本文模型中,將BN層直接加在卷積通道的每個卷積操作之后。同時為保證直通通道的暢通性,不在直通通道上加入該操作。

3.2 Padding

Padding技術[23]通過在特征矩陣外圍補零,控制卷積層輸入、輸出特征矩陣維度。Padding可增強對邊緣數據的特征提取,同時對航班延誤數據進行維度控制,使特征提取部分的深層網絡得以訓練。

設輸入特征矩陣(input)、卷積核(convolution kernel)、Padding、步長(stride)、輸出特征矩陣(output)在不同軸上長度相等,分別為i、c、p、s、out,則在是否采用Padding時,各長度關系如式(20)所示:

out=

(20)

由于航班延誤數據具有先驗信息較少的特點,經數據預處理后維度較小,在不采用Padding技術時,特征提取部分維度會迅速下降,無法訓練深層網絡,因此在本模型中采用p=1的Padding技術,同時和c=3、s=1的卷積核相結合,可保證特征矩陣維度不變。

4 實驗

實驗運行的軟件環境是在64位Ubuntu14.04操作系統下搭建的Spark開源框架[20]和Caffe[24]深度學習開源框架,本文對官方版本的Caffe進行了部分修改,并使用GPU加速計算。實驗計算機配置:Inter Xeon處理器,16 GB內存,GeForce GTX TITAN X顯卡。

本文模型所用數據為美國交通運輸統計局(Bureau of Transportation Statistics, BTS)提供的航班準點數據(Airline On-Time Performance data, AOTP)和美國國家氣象數據中心(National Climatic Data Center, NCDC)提供的本地質量控制氣象數據(Quality Controlled Local Climatological Data, QCLCD)。其中AOTP包含從1987年至今的航班數據,QCLCD包含美國1 600個氣象觀測站的氣象數據,氣象數據每小時進行一次更新。實驗中使用兩數據源中2015—2016年數據,輸入數據源包含航班數據103 405 986條,44維特征;氣象數據11 235 340條,48維特征。融合后數據為103 402 300條,編碼為217維特征。最終經矩陣變換輸入網絡中的數據維度為103 405 986×15×15。按表1進行等級劃分時,類別1到5的比例為46∶6∶2∶2∶1,其中不延誤的比例很大,為偏分類問題。為使數據盡量均衡,在網絡訓練迭代時對類別1的數據進行隨機抽樣,抽取與類別2相同個數,使輸入樣本比例變為6∶6∶2∶2∶1。

航班延誤預測模型采用文獻[25]提出的權值初始化策略。使用隨機梯度下降法對模型進行訓練,批處理數量(min-batch size)為512,動量因子為0.9,權重衰減系數為0.000 1,訓練時總的迭代次數為5 000,初始學習率為0.1。DCNN中選用尺寸為3×3、步長為1的卷積核進行卷積,選用尺寸為3×3、步長為2的濾波器進行池化。為探索網絡深度對特征提取的影響分別使用16,24,32,40,48和56層模型結構,網絡層數配置如表2所示。在實驗對比時,DCNN與傳統卷積神經網絡均采用softmax分類方法,以保證實驗條件的一致性。

表2 網絡層數配置

4.1 卷積衰減因子影響分析

從式(13)可以看出卷積衰減因子的取值影響了反向傳播時梯度的大小,從而控制網絡收斂。若卷積衰減因子設置較大則會導致梯度爆炸現象,違背了設置雙通道卷積神經網絡的初衷;若設置過小則會導致返回的梯度值減小,使網絡收斂變慢,增加了模型的訓練時間。實驗中,采用如表3所示的卷積衰減因子進行訓練,以探索設置的最佳準則。其中,NAN表示網絡未收斂。

表3 卷積衰減因子對正確率影響

從表3可以看出:在淺層網絡時,設置λ1=0.2,λ2=0.1預測效果較好;而網絡加深時,同樣的設置正確率卻下降。原則上,越深的網絡結構應該對應越小的卷積衰減因子,但設置過小會使網絡在有限的迭代次數內沒有收斂,實驗結果無法達到最好。在深層網絡中,設置為λ1=0.1,λ2=0.05,達到最優,因此本文后續實驗將卷積衰減因子設置為λ1=0.1,λ2=0.05。

4.2 氣象數據影響分析

文獻[21]顯示,我國上海、深圳、南京等重要城市航班延誤因素中氣象因素占比高達40%以上。從實際情況分析,風速和風向的變化會引起風切變等極端氣象情況,使飛機無法起降;降雨、降雪等極端惡劣天氣也都伴隨著風向、風速和濕度的改變;風暴天氣則會大幅影響能見度,因此,本文在模型的數據預處理過程進行了航班數據和氣象數據的融合。為驗證數據融合的重要性,首先取60 min作為判斷是否延誤的閾值,作二分類預測。以刪除某個氣象特征時,延誤召回率較不刪除時的下降值表征此氣象特征對延誤預測的影響程度。分別去掉某個氣象特征后,召回率下降值如圖3所示。

從實驗結果可以看出,風向、風速和能見度對延誤預測影響較大,濕度和氣壓其次,天空狀況對預測影響最小。DCNN具有很強的大數據處理能力,因此本著盡可能提高正確率的思路,不對特征作降維處理。

對航班數據和融合后數據分別進行延誤預測。表4給出了航班數據和融合后數據的預測準確率。融合后數據預測準確率在各層均增加,增加范圍基本保持在1個百分點,由此可見氣象數據是航班延誤預測時重要的先驗信息。

表4 數據融合前后正確率對比

4.3 正確率

為驗證直通通道對模型“暢通性”的影響,對使用傳統卷積神經網絡和使用直通通道改進后的卷積神經網絡進行實驗對比。使用傳統卷積神經網絡時,為保證可對比性,除最后一層卷積層,其他各層也相應增加ReLU激活函數進行非線性變換并相應加入BN層和Padding技術進行優化。

由表5可知,傳統卷積神經網絡在32層時效果達到最佳,正確率為0.912;但從40層開始,傳統卷積神經網絡正確率呈下降趨勢,到56層時,正確率已降至0.790。而雙通道卷積神經網絡隨深度的增加預測正確率逐漸提升,當網絡達到56層時,其預測正確率已達到0.921。

表5 不同網絡結構正確率對比

為進一步探索模型可訓練深度,對72,104和136層網絡進行訓練,訓練結果如表6所示。由結果可知,在深層網絡模型中,傳統卷積神經網絡性能進一步下降,而雙通道卷積神經網絡性能一直保持穩定。

表6 深層網絡正確率對比

4.4 損失值

由表5可知不同網絡結構的差異性在深度為40層時開始凸顯,現著重討論40層以后兩種網絡結構的損失值變化。圖4所示為損失值隨迭代次數變化示意圖。對比可知本文模型訓練時損失值可減小至0.3以下,而傳統模型一直維持在0.5以上,說明DCNN相比傳統的卷積神經網絡能夠更好地對深層網絡進行訓練。

圖4 不同網絡結構的損失值隨迭代次數變化示意圖

4.5 梯度值

采用DCNN結構進行特征提取的航班延誤預測模型之所以能夠得到更高的準確率,得益于正向信息的無損傳輸和反向梯度的穩定性,從而使式(10)和式(11)更好地進行卷積權值更新。

圖5所示為40、48和56層網絡結構第一層對應梯度值(簡稱梯度)隨迭代次數變化示意圖。圖5(a)為傳統卷積神經網絡梯度變化,可以明顯看出隨深度的增加,梯度銳減且梯度隨迭代次數變化極不穩定,最小時可達到10-14左右。這說明傳統卷積神經網絡模型在權值更新時存在梯度消失問題,難以進行充分的特征提取。圖5(b)中相應梯度基本穩定在0.01至0.001之間,由此可見本文的DCNN模型可克服梯度消失缺陷,進行更加完備的信息表達。

4.6 算法對比

文獻[8]和文獻[13]均與本文一樣采用分類預測的方式進行延誤等級判斷。其中:文獻[8]采用國內31 613條航班數據建模,使用C4.5決策樹算法預測,并與樸素貝葉斯(Native Bayes)算法進行對比;文獻[13]采用與本模型相同數據源的單個機場數據進行預測,利用傳統神經網絡與LSTM組合的方式構建淺層模型。表7所示為本文深層模型與這兩種模型中所提算法的正確率對比。

C4.5決策樹和樸素貝葉斯算法都是傳統機器學習算法的典型代表,兩種算法均具有較好的分類性能且計算簡單,但當數據量增長時兩個算法都會遇到計算瓶頸,無法充分利用大數據優勢,因此它們做的是在小樣本情況下的航班延誤預測,正確率較低。文獻[13]中,作者同樣采用了加深深度的方式提高正確率,這使得該文中基于LSTM算法建立的模型與基于C4.5決策樹和樸素貝葉斯算法建立的模型相比,預測正確率得到大幅提升;但是由于沒有克服深度網絡中梯度消失的缺點,基于LSTM的模型在加深深度時難以收斂,最終只是采用了淺層的網絡結構。

圖5 不同網絡結構的梯度值隨迭代次數變化示意圖

Tab. 7 Correct rate comparison of different algorithms

對比可知,與傳統算法相比,深度卷積神經網絡可進行基于大數據的信息挖掘,因此在處理數據量與維度上具有較強優勢。同時與淺層的LSTM算法相比,深度雙通道卷積神經網絡的正確率更高。實驗表明采用本文模型進行基于大數據的航班延誤等級預測可有效提高正確率。

5 結語

本文創新性地提出使用數據融合和深度雙通道卷積神經網絡構建航班延誤預測模型。通過將延誤時間進行等級劃分,將預測問題轉化為分類問題,分數據預處理、特征提取和分類預測三步完成航班延誤等級預測,并采用批歸一化和Padding技術進行模型優化。經實驗驗證,結論如下:

1)融合氣象數據后,CNN和DCNN進行特征提取后的模型性能都有所提升,證明氣象數據對于航班延誤預測十分重要。

2)當特征提取部分的網絡加深時,使用CNN的模型性能變差,而利用DCNN的模型仍保持穩定,證明傳統卷積神經網絡難以訓練深度模型,而改進后的網絡結構可以保證深度網絡的訓練。

3)深度DCNN算法與傳統算法相比具有更強的數據處理能力,使數據量和維度大幅提升,最終可得到更高的正確率。

本文的預測結果與傳統卷積神經網絡相比有較大提升,有利于進一步改善空管、航空公司等部門的決策,提升民航綜合服務能力;但受數據來源限制,本模型在構建時未使用國內航班數據,目前正與國內空管部門爭取合作,希望下一步使用國內數據進行建模。后續還將結合混淆矩陣分析預測結果,使實驗數據更加完善。

猜你喜歡
特征提取模型
一半模型
特征提取和最小二乘支持向量機的水下目標識別
重要模型『一線三等角』
重尾非線性自回歸模型自加權M-估計的漸近分布
基于Gazebo仿真環境的ORB特征提取與比對的研究
電子制作(2019年15期)2019-08-27 01:12:00
基于Daubechies(dbN)的飛行器音頻特征提取
電子制作(2018年19期)2018-11-14 02:37:08
Bagging RCSP腦電特征提取算法
3D打印中的模型分割與打包
FLUKA幾何模型到CAD幾何模型轉換方法初步研究
基于MED和循環域解調的多故障特征提取
主站蜘蛛池模板: 亚洲AⅤ波多系列中文字幕| 99资源在线| 国产福利一区二区在线观看| 亚洲欧洲日本在线| 亚洲va欧美va国产综合下载| 国产精品19p| AV在线天堂进入| 国产乱码精品一区二区三区中文| 免费人成网站在线观看欧美| 在线看片免费人成视久网下载| 91免费观看视频| 波多野结衣一二三| 另类欧美日韩| 国产午夜人做人免费视频中文| 国产一级毛片在线| 日本a级免费| 国产乱人视频免费观看| 国产杨幂丝袜av在线播放| 99热国产这里只有精品9九| 伊人久久久久久久久久| 999国产精品永久免费视频精品久久| 无码一区二区波多野结衣播放搜索| 性视频一区| 99福利视频导航| 亚洲人成日本在线观看| 中国一级特黄大片在线观看| 国产99精品久久| www.99在线观看| 国产欧美网站| 国产在线专区| 强奷白丝美女在线观看 | 国产精品99久久久| 日韩国产精品无码一区二区三区 | 99热这里都是国产精品| 啪啪永久免费av| 一级毛片a女人刺激视频免费| 亚洲性视频网站| 日本a∨在线观看| 91年精品国产福利线观看久久| 在线五月婷婷| 欧美精品v欧洲精品| 人妻中文字幕无码久久一区| 波多野结衣国产精品| 久久久91人妻无码精品蜜桃HD| 六月婷婷精品视频在线观看| 91区国产福利在线观看午夜| 国产一区二区三区在线无码| 欧美69视频在线| 精品撒尿视频一区二区三区| 国产成人在线无码免费视频| 色成人综合| 呦视频在线一区二区三区| 黄色片中文字幕| AV熟女乱| 亚洲日韩每日更新| 中文字幕日韩视频欧美一区| 久久精品aⅴ无码中文字幕| 大陆精大陆国产国语精品1024| 最新午夜男女福利片视频| 亚洲国产清纯| 亚洲国产欧美目韩成人综合| 自拍偷拍欧美日韩| 亚洲精品va| 久久久久久国产精品mv| 亚洲国产日韩欧美在线| 久久久久无码精品| 一边摸一边做爽的视频17国产| 亚洲国产精品一区二区第一页免 | 国产一区二区人大臿蕉香蕉| 国产亚洲视频免费播放| 国产成人精品免费视频大全五级| 女人爽到高潮免费视频大全| 国产免费黄| 国产成人永久免费视频| hezyo加勒比一区二区三区| 国产午夜一级淫片| 特级毛片8级毛片免费观看| 无码网站免费观看| 视频二区中文无码| 国产后式a一视频| 日本午夜影院| 9cao视频精品|