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

基于深度學習的密閉式豬舍內溫濕度預測模型

2020-10-29 12:16:32謝秋菊蘇中濱
農業機械學報 2020年10期
關鍵詞:模型

謝秋菊 鄭 萍 包 軍 蘇中濱

(1.東北農業大學電氣與信息學院, 哈爾濱 150030;2.東北農業大學農業農村部生豬養殖設施工程重點實驗室, 哈爾濱 150030;3.東北農業大學動物科學技術學院, 哈爾濱 150030)

0 引言

隨著商業利益和消費需求的不斷增長,環境可控式的密閉豬舍養殖方式發展迅速。在密閉式豬舍的有限空間中環境品質尤為重要,它是除遺傳基因之外制約養豬生產發展的關鍵因素[1]。在眾多因素中,溫度是維持豬體溫恒定、影響豬的健康水平和繁殖能力的重要外部環境因素之一,直接影響豬體的熱平衡[2]。豬體通過產熱和散熱平衡來維持恒定的體溫,進行正常的生命活動[3]。豬舍空氣濕度影響豬的體熱調節,高溫高濕環境嚴重影響豬的日增質量,同時造成細菌滋生,誘發疾病[4]。因此,營造適宜的舍內熱濕環境可以有效地提高豬群的健康水平、繁殖能力,降低飼料消耗,提高生長效率。研究低延遲、精準化的豬舍環境調控方法是實現豬舍內環境精準控制的關鍵。

豬舍內熱濕環境受豬舍建筑樣式、舍內飼養豬的數量和質量、通風模式[5]、豬的活動量、生產方式、舍外天氣情況等多種因素相互耦合的影響。國內外許多學者針對畜禽舍內溫濕度控制進行了大量研究[6-17]。

這些研究大多基于傳感器網絡技術實時監測豬舍中的環境狀況,并做出判斷,實現舍內環境調控。傳感器只能對當前舍內環境狀況進行監測,無法提前對舍內的環境因子變化趨勢作出預判,因而只依靠傳感器監測數據來決定控制系統的運行狀態。然而,由于豬舍內的熱濕環境隨著舍內的氣流而逐漸變化,因此導致環境調控設備運行后對環境調控達到預期效果會存在一定的滯后性,這種滯后性因受豬舍空間結構、建筑樣式、舍內通風、加熱等控制設備運行效率的影響而不同。試驗測試表明,在開啟通風或加熱等設備20 min后舍內熱濕環境可以得到明顯改善。因此,需要對豬舍內環境因子變化提前作出預測,從而提前決定環境控制設備的運行狀態,實現舍內環境優化控制、降低控制效果的滯后。目前,隨著傳感器及大數據技術的發展,豬舍環境監測的數據量不斷增加,傳統的數據驅動預測方法,如線性[18]和非線性的方法[19-20]、CFD方法[8-9]、神經網絡方法[21-23]等,已經無法實現大數據量的輸入時序序列數據的精確預測。

豬舍內溫濕度監測數據是由豬舍環境監測傳感器在一定時間間隔內獲取的大量時序數據。長短時記憶(Long short-term memory,LSTM)網絡是由GERS等[24]提出的深度學習的典型代表方法之一。LSTM 通過對時序數據歷史特征的提取實現數據未來變化的預測,它適合于多個領域時序數據處理及預測應用[25-28]。然而,根據傳感器監測獲取的大量豬舍內溫濕度時序數據,應用LSTM預測模型的研究尚未見報道。

本文基于深度學習方法,結合傳感器監測的大量時序序列歷史數據,綜合考慮各環境因素對豬舍內溫濕度變化的影響,建立具有多層結構的長短時記憶LSTM網絡的溫濕度變化預測模型,以實現當前時間豬舍內溫濕度的變化預測。

1 影響因素相關性分析

密閉豬舍是一個由磚、保溫板等為圍護結構的封閉空間。豬舍內環境溫濕度變化受豬舍建筑樣式、圍護結構,舍內飼養豬的數量和質量,舍內加熱及通風系統,舍外天氣狀況等的影響。建立豬舍內溫濕度變化預測需要對影響因素進行分析,確定影響因素,為建立優化的預測模型奠定基礎。本文在實際監測的數據基礎上,考慮除建筑樣式和圍護結構相關的影響因素,利用皮爾遜(Pearson) 相關分析方法 (P<0.05),確定豬舍內溫濕度變化的相關影響因素,如表1所示。

表1 影響因素相關性系數Tab.1 Correlations coefficients between factors

從表1可以看出,豬舍內溫度(Pig living space temperature, PLS T)與豬舍外溫度(Outdoor temperature, Out T)、糞坑溫度(Pit temperature, Pit T)和豬的活動量(Activity, Act)呈現較大的正相關性,相關系數分別為0.489、0.518和0.395;同時,豬舍內溫度與豬舍外相對濕度(Pig living space relative humidity, PLS RH)和豬的總質量(Weight, Wt)呈現負相關性,相關系數分別為-0.315和-0.436。因此表明,豬舍內溫度受豬舍外溫度、舍內豬的總質量、糞坑溫度、舍外相對濕度影響較大;此外,豬舍內溫度與豬舍內相對濕度和通風量(Ventilation, Vent)也具有一定的相關性,相關系數分別為0.284和0.147。

豬舍內相對濕度與豬的活動量呈現正相關性,相關系數為0.378,表明當豬由于進食或排泄等的活動量增大時可導致舍內相對濕度增大;同時,豬舍內相對濕度與豬舍內溫度和糞坑溫度呈現正相關性,相關系數分別為0.284和0.293;與豬的總質量和豬舍外相對濕度(Outdoor relative humidity, Out RH)呈現負相關性,相關系數分別為-0.253和-0.178。

因此,本文選取豬舍內相對濕度、豬舍外相對濕度、豬舍外溫度、豬的活動量、豬的總質量、通風率、糞坑溫度作為豬舍內溫度預測的外部影響因素;選取豬舍內溫度、豬舍外溫度、豬舍外相對濕度、豬的活動量、豬的總質量、通風率、糞坑溫度作為豬舍內相對濕度預測的外部影響因素。

2 基于LSTM的豬舍內溫濕度預測模型

2.1 模型結構

預測模型由時序數據輸入層、隱含層和輸出層組成。模型結構如圖1所示。

圖1 LSTM預測模型結構圖Fig.1 Structure diagram of LSTM prediction model

輸入層負責接收采集到的豬舍內溫度、相對濕度歷史數據,舍外溫度、相對濕度,糞坑溫度,通風率,豬的活動量及總質量數據,并且將其進行標準化處理。標準化處理公式為

(1)

式中z′i——標準化后的數據

zi——輸入的時序數據

n——時序序列數據的長度

標準化后的豬舍內溫度和相對濕度歷史數據轉換成時序序列x1和x2。豬舍內溫度和相對濕度的外部影響因素被標準化處理后參與完全連接層的運算。輸入層節點數量為2t,t為豬舍內溫度和相對濕度歷史時序序列長度。

隱含層包含LSTM、完全連接和回歸3個子層。LSTM子層負責學習時間序列數據每步之間的長期依賴關系,以執行有助于改善長序列梯度流的交互。LSTM子層由元胞狀態c、隱含狀態h、輸入門i、遺忘門f、候選門g和輸出門o組成,如圖2所示[29-30]。

圖2 LSTM子層組成結構圖Fig.2 Consist structure diagram of LSTM sublayer

元胞狀態c包含從前一步學習來的信息,并且在每一步的學習中,LSTM子層將從元胞狀態c中加入或刪除信息。每一步的隱含狀態h都包含上步LSTM層的輸出信息。為了使信息能夠更好地更新,LSTM層增加了輸入門i、遺忘門f、元胞候選門g和輸出門o。輸入門i決定元胞的更新程度,遺忘門f決定元胞狀態的重置水平,候選門g負責向元胞中加入一些新的信息,輸出門o決定將有多少元胞信息加入隱含狀態。

每一步元胞狀態和隱含狀態計算式為

ct=ft°ct-1+it°gt

(2)

ht=ot°σ(ct)

(3)

式中 °——向量元素相乘符號

it、ft、gt、ot——輸入門、遺忘門、候選門、輸出門

σ——狀態激活函數

在每一個時間步t中, 輸入門it、遺忘門ft、輸出門ot和候選門gt可描述為

it=σg(wixt+Riht-1+bi)

(4)

ft=σg(wfxt+Rfht-1+bf)

(5)

gt=σc(wgxt+Rght-1+bg)

(6)

ot=σg(woxt+Roht-1+bo)

(7)

式中σc——狀態激活函數,采用雙曲正切函數tanh來計算狀態激活函數

σg——門狀態函數

激活函數由公式σ(x)=(1+e-x)-1計算得到;LSTM層中的權重是由具有可學習性的輸入矩陣W、回歸權重矩陣R以及每個分量偏差矩陣b組成。這些矩陣表示為

完全連接子層中的所有節點都與LSTM子層中的節點以及外部影響數據相連接,以獲得上一層和外部影響因素對輸入變量x1和x2學習所得的所有特征。完全連接層的每一步都是獨立執行的,輸入變量與權重矩陣W相乘后,再加上偏置矩陣向量b。此層輸出變量數目等于輸入變量的數目2,即溫度和相對濕度2個輸入變量。

回歸子層與完全連接子層相連,用于計算序列到序列網絡的均方根誤差損失(Half-mean-squared-error loss)。每步的未歸一化的損失預測響應計算式為

(8)

式中L——計算序列到序列網絡的均方根誤差損失

S——序列的長度

tkj、ykj——第k步的序列值、網絡預測值

當網絡訓練時,觀察值的平均損失由小批量網絡訓練時計算而得。

輸出層將隱含層輸出量O1和O2反歸一化處理成預測值豬舍內溫度y1和豬舍內相對濕度y2,反標準化公式采用式(1)對觀察值進行標準化處理時的變換所得。

2.2 模型評價

本文通過統計分析方法中的決定系數R2、均方根誤差eRMSE和平均絕對誤差百分比eMAPE來評價LSTM預測模型的性能。

2.3 模型參數設置

在LSTM子層中,初始網絡參數設置為:隱含層神經元數為200,狀態激活函數為雙曲正切函數tanh,門激活函數為sigmod,采用自適應矩估計Adam優化器進行訓練,最大步數設置為250,初始學習速率為0.005,時間間隔(Time lag)為1。

為了獲得最優的LSTM預測模型,需要對預測模型中的LSTM子層和完全連接層中的神經元數、網絡的學習速率2個主要參數進行對比及組合優化。設LSTM子層和完全連接層中的神經元個數取值集合為{200,400,800,1 200,1 600,2 000},網絡的學習速率取值集合為{0.001, 0.003, 0.005, 0.008, 0.010}。

2.4 模型數據集選取

時序序列數據集分為訓練數據集和驗證數據集。訓練數據集訓練LSTM模型自學習及調整連接參數等,驗證數據集對模型相關性能參數進行驗證。本文中所用數據集包含從2015年8月28日至9月30日(代表夏季數據集)和2015年12月24日至2016年1月13日(代表冬季數據集)共計1 320組監測數據。將這些數據分成4個數據集,在夏季數據集中前94%的監測數據(從2015年8月28日至9月28日)用作模型訓練,剩余6%的監測數據(從2015年9月29日至30日)用于模型驗證;在冬季數據集中前90.4%的監測數據(從2015年12月24日至2016年1月11日)用于模型訓練,剩余9.6%的監測數據(從2016年1月12日至13日)用于模型驗證。

本文使用Matlab對LSTM模型進行訓練及驗證。在LSTM模型中共用到8個變量,分別作為輸入變量和外部影響變量進行模型特征提取,實現下一時刻豬舍內溫度和相對濕度的預測。其中,豬舍內溫度、相對濕度歷史數據作為模型的輸入變量;豬舍內相對濕度、豬舍外溫度、豬舍外相對濕度、豬舍糞坑溫度、舍內通風率、豬的活動量、豬的總質量作為豬舍內溫度預測的外部影響變量;豬舍內溫度、豬舍外相對溫度、豬舍外濕度、豬舍糞坑溫度、舍內通風率、豬的活動量、豬的總質量作為豬舍內相對濕度預測的外部影響變量。

3 預測結果實例驗證

3.1 數據來源及數據監測

本文LSTM預測模型中用到的訓練及驗證數據來源于普渡大學動物研究教育中心豬舍[31]開展的研究項目(the USDA National Institute of Food and Agriculture Hatch project)。該豬舍位于美國印第安那州西拉法葉市,地處大陸性濕潤氣候,冬季寒冷,1月氣溫最低,平均溫度為-7.2℃,相對濕度為69%;夏季炎熱,7月溫度最高,平均溫度為29.4℃,平均相對濕度為71%。

該豬舍為密閉式機械通風、鋼木混結構。豬舍共有12個南北朝向的小間,舍外新鮮空氣通過豬舍東西兩側進風口送入舍內,然后通過棚頂和過道入口送入各個小間。單個小間尺寸為11 m×6.1 m×2.7 m(長×寬×高),內部建有雙列鋼質圍欄,中間為過道,每列有6個豬欄,可以容納60頭育肥豬。小間內是水泥漏縫地板,地板下方有1.8 m深的儲糞池。豬舍外部實景及結構如圖3所示[32]。

圖3 豬舍外部實景和結構圖Fig.3 Photo and structure diagram of pig building

每個小間分別安裝有天然氣加熱器(Guardian 60型,L.B. White Co.,美國)和2個直徑分別為356 mm(V4E35型,Multifan,美國)和508 mm(V4E50型, Multifan,美國)的墻上恒速風機,用于舍內加熱和通風。舍內溫度和糞坑溫度數據分別通過安裝在距離豬舍地面1.5 m和儲糞池的排風口處的T型熱電耦監測得到。豬舍內相對濕度通過安裝在豬欄中間距離地面1.5 m處的相對濕度傳感器(HX92BC型,Omega,美國)測量得到。小間內墻上風機所受壓力通過氣壓測量儀(260型,Setra Systems,Inc.,美國)測得。

豬舍的通風率包括舍內墻上風機的通風率和糞坑風機通風率。在豬舍下方的儲糞池外墻上部安裝有2個直徑為250 mm的變速風機(P4E30型,Multifan,美國),為每間豬舍提供最小量通風,糞池風機的通風率由風速儀直接測量得到[32];墻上風機的通風率通過持續測量墻上風機的運行時間、小間內靜壓力以及后期風機測量模型[33]間接計算得到。

在每個小間墻上安裝有紅外活動量傳感器(Visonic SRN 2000 Detector型,Visonic Inc.,美國)將探測到的紅外光線信號轉換成電壓信號,用于豬的活動量(以電壓度量)監測[34]。

每個小間中豬的數量是由人工每個星期進行計數得到。豬的質量由豬的初始質量和生長模型計算而得到,豬的生長模型曲線表示為

W0~29=0.008 592 3D2+0.107 91D+5.35

(9)

W30~150=0.002 339 7D2+0.562 14D-3.31

(10)

式中W0~29、W30~150——入欄0~29 d和30~150 d的豬質量

D——入欄的時間,d

由生長曲線計算得到的豬的質量與實際質量的擬合度為0.999 8,說明由生長模型計算得到的豬的質量完全可以反映豬的質量變化情況。這些數據根據需要,插值處理成時間間隔為1 h的數據。

豬舍外部溫度與相對濕度等天氣數據來源于距離豬舍5 km的普渡大學園藝研究與教育中心,以及安裝在距離豬舍南側70 m的Davis 6152 Wireless Vantage Pro2型氣象站(Davis Instruments,美國)監測得到。

根據豬舍所處地理位置、氣候特點以及監測數據的完整性等,本文分別選取2個朝南向的小間2和小間4,從2015年8月28日至9月30日和2015年12月24日至2016年1月13日2個階段的采集間隔為1 min的舍內、外環境數據,將這些數據處理為時間間隔為1 h的數據進行模型訓練與驗證,實際監測數據變化,如表2和圖4所示。

表2 實際監測數據變化范圍Tab.2 Variations of collected data

圖4 豬舍內外環境監測數據Fig.4 Indoor and outdoor environmental data of pig building

3.2 模型參數性能對比及選取

3.2.1不同學習速率

從不同網絡的學習速率取值集合{0.001,0.003,0.005,0.008,0.010}中對模型進行測試,以得到高時效、訓練損失較小的優化模型。不同學習速率下網絡損失如表3所示,當學習速率取值增大時,網絡訓練損失迅速下降。雖然在學習速率0.010時網絡訓練損失值最小,但是,網絡損失變化率最小,在夏季和冬季分別為0.500 0和0.517 2;然而在學習速率為0.008時的預測模型在夏季和冬季的網絡訓練損失變化率分別為0.843 3和0.701 0,基本達到最大值,表明在學習速率大于0.008后,預測模型的網絡訓練損失下降不再明顯。因此綜合考慮網絡訓練損失和訓練時間效率,本文選取0.008作為LSTM預測模型的網絡訓練學習速率。

表3 不同學習速率下的網絡訓練損失Tab.3 Training losses at different learning rates

3.2.2不同節點數

網絡訓練的學習速率為0.008時,對LSTM子層和完全連接子層的節點數進行測試,以獲得較優的預測性能。如圖5a、5b所示,舍內溫度和相對濕度夏季預測模型的eRMSE和eMAPE雖然在節點數為2 000和1 600時達到最小值,然而當節點數大于400時,eRMSE和eMAPE下降較緩慢,并且相對濕度預測的eRMSE變化曲線有逐漸上升的趨勢;冬季溫度和相對濕度預測模型都在節點數為400時,eRMSE和eMAPE達最小值。如圖5c所示,夏季舍內溫度和相對濕度預測模型的決定系數R2分別在節點數為2 000和1 600時達到最大值,然而當節點大于400時,溫度預測模型的R2增加較為緩慢;冬季溫度和相對濕度預測模型都在節點數為400時,R2達最大值。因此,根據以上溫度和相對濕度預測模型的eRMSE、eMAPE和R2的性能變化分析,同時結合模型節點數與運行效率考慮,本文選取節點數400作為預測模型LSTM子層和完全連接層神經元的數目,優化模型的預測性能。

圖5 不同節點數時LSTM預測模型的性能變化Fig.5 Performances of LSTM mode with different neurons

3.3 預測結果驗證

選取2015年9月29—30日和2016年1月12—13日4 d監測的豬舍內環境數據以及舍外的天氣數據作為豬舍內夏季、冬季溫度和相對濕度驗證數據對預測模型進行驗證。

3.3.1溫度預測驗證

如圖6所示,豬舍內溫度預測值和實際測量值變化趨勢非常一致,能夠很好地反映舍內溫度變化情況。

圖6 舍內溫度預測值與實測值對比Fig.6 Comparisons of indoor temperature predictions and observations

如圖6a所示,夏季豬舍內溫度高峰出現在14:00—15:00,溫度低谷出現在01:00—02:00;豬內實測溫度變化范圍為21.5~26.6℃,預測溫度變化范圍為21.6~25.3℃;預測值與實測值的最大誤差為1.9℃,最小誤差為0℃,平均誤差0.6℃;預測值與實測值的平均絕對百分比誤差eMAPE為2.504%,決定系數R2為0.703。

如圖6b所示,冬季豬舍內溫度高峰出現在22:00—01:00,溫度低峰出現在04:00—05:00;溫度實際測量值變化范圍為17.6~22.1℃,預測溫度變化范圍為17.6~22.1℃;預測值與實際測量值的最大誤差為1.3℃,最小誤差為0,平均誤差為0.3℃;預測值與實測值的平均絕對百分比誤差eMAPE為1.392%,決定系數R2為0.938。

因此,無論在冬季還是在夏季,溫度預測值與實測值的誤差都較小,能夠較準確地反映豬舍內溫度變化情況,本文所建立的溫度預測模型,可以作為豬舍內溫度優化調控策略的參考,實現舍內風機或加熱器等溫度控制設備的預先啟停控制,在一定程度上彌補環境控制效果滯后問題。

3.3.2相對濕度模擬驗證

豬舍內相對濕度預測值與實際測量值對比如圖7所示,豬舍內相對濕度預測值與實測值的整體變化趨勢一致。

圖7 舍內相對濕度預測值與實測值對比Fig.7 Comparisons of indoor relative humidity predictions and observations

如圖7a所示,由于在9月29日06:00起完全開啟直徑為508 mm的風機進行通風,通風量達到9 800 m3/min以上,此時刻的舍外溫度為20.6℃,舍內相對濕度在10:00出現低峰值28.2%,相對濕度變化較大;夏季舍內相對濕度實際測量值變化范圍為28.2%~50.1%,平均值為46.8%;相對濕度預測值變化范圍為34%~50.7%,平均值47.1%;相對濕度預測值與實測值的最大誤差為7%,最小誤差為0,平均誤差為1.6%;相對濕度的預測值與實測值eMAPE、eRMSE、R2分別為3.769%、2.346%、0.747。

如圖7b所示,冬季豬舍內相對濕度變化呈現正弦周期曲線變化趨勢,高峰值出現在13:00—15:00,低峰值出現在02:00—06:00;相對濕度實測值變化范圍為27.7%~55.2%,相對濕度預測值變化范圍為35.4%~51.7%;相對濕度預測值與實測值的最大誤差為13.5%,最小誤差為0,平均誤差為3.1%;相對濕度的預測值與實測值的eMAPE、eRMSE、R2分別為7.281%、4.245%、0.542。

根據豬舍環境控制要求及豬舍環境管理經驗,在豬舍各環境因素中舍內相對濕度重要性稍次于舍內溫度,允許在一定范圍內舍內相對濕度稍微大一些的變化。雖然,冬季相對濕度預測值與實測值最大誤差為13.5%,但是仍然滿足舍內相對濕度控制要求。因此,本研究建立的舍內相對濕度預測模型能夠反映舍內實際相對濕度環境變化情況。

3.3.3結果分析

本文建立的LSTM溫濕度預測模型是結合傳感器監測所得豬舍內外環境歷史時序數據和相關影響因素時序數據,對建立的預測模型進行訓練,獲取最優化的網絡參數;而文獻[20]是通過能量和質量平衡方程,建立豬舍內動態變化的熱濕小氣候模型,結合傳感器監測數據,利用多元非線性回歸的方法和一些經驗值來確定模型中的相關參數,會有一定的局限性。

文獻[20]中,溫度模擬值與實測值的最大誤差為2.4℃,平均誤差為0.6℃,平均決定系數R2為0.812;相對濕度模擬值與實測值的最大誤差為13.3%,最小平均誤差為3.35%,平均決定系數R2為0.905。

本文建立的預測模型中,溫度預測值與實測值最大誤差為1.9℃,平均誤差為0.5℃,平均決定系數R2、平均eRMSE和平均eMAPE分別為0.821、0.605℃和1.984%;相對濕度預測值與實測值的最大誤差為13.5%,平均誤差為2.3%,平均決定系數R2、平均eRMSE和平均eMAPE分別為0.645、3.297%和5.525%。

因此,綜合考慮各項參數,本文建立的溫濕度預測模型性能較優。

4 結論

(1)通過實測數據對LSTM溫濕度預測模型參數進行訓練優化,確定網絡學習速率為0.008,完全連接子層神經元數目為400時,模型可以達到較優的性能。舍內溫度預測的平均eRMSE和eMAPE分別為0.605℃和1.984%,舍內相對濕度預測的平均eRMSE和eMAPE分別為3.297%和5.525%。

(2)無論在夏季還是在冬季,采用建立的預測模型舍內溫度與相對濕度預測值與實測值變化趨勢一致。預測溫度與實測溫度的平均決定系數R2為0.821,舍內溫度預測值與實測值最大誤差為1.9℃,平均誤差為0.5℃;預測相對濕度與實測相對濕度的平均決定系數R2為0.645,舍內相對濕度預測值與實測值的最大誤差為13.5%,平均誤差為2.3%。

(3)提出的基于深度學習的溫濕度預測模型較精確地反映了豬舍內溫度與相對濕度的變化情況,將預測模型的預測輸出加入到控制策略中,可以預測下一小時的豬舍環境狀況、提前對環境控制設備發出控制指令,控制設備提前運行可緩解控制效果滯后的問題,為豬舍內環境優化控制策略提供可靠的參考。

猜你喜歡
模型
一半模型
一種去中心化的域名服務本地化模型
適用于BDS-3 PPP的隨機模型
提煉模型 突破難點
函數模型及應用
p150Glued在帕金森病模型中的表達及分布
函數模型及應用
重要模型『一線三等角』
重尾非線性自回歸模型自加權M-估計的漸近分布
3D打印中的模型分割與打包
主站蜘蛛池模板: 久久久久久国产精品mv| 午夜a视频| 久久久精品无码一区二区三区| 欧美在线综合视频| 天天躁日日躁狠狠躁中文字幕| 一本大道无码日韩精品影视| 欧美三级自拍| 91色国产在线| 国产jizzjizz视频| 久久不卡精品| 亚洲欧美日韩成人在线| 久久99热66这里只有精品一| 亚洲女人在线| 欧美在线精品怡红院| 久久综合激情网| 乱系列中文字幕在线视频| 久久综合九色综合97网| 免费A∨中文乱码专区| 亚洲Aⅴ无码专区在线观看q| 国产第一页免费浮力影院| 91无码视频在线观看| 免费在线国产一区二区三区精品| 日本人妻丰满熟妇区| 91亚洲精选| 国产精品无码久久久久AV| 国产大片喷水在线在线视频| 亚洲高清在线天堂精品| 精品国产一区二区三区在线观看 | 国产在线拍偷自揄观看视频网站| 欧美一级高清免费a| 久久动漫精品| 亚洲一级毛片在线观播放| 亚洲日韩精品伊甸| 成人福利免费在线观看| 国产一区二区三区免费| 国产超薄肉色丝袜网站| 国产菊爆视频在线观看| 欧美日韩成人在线观看| 精品视频在线观看你懂的一区| 欧美日韩专区| 伊人激情久久综合中文字幕| 国产精女同一区二区三区久| 亚洲美女一级毛片| 国产成人综合日韩精品无码不卡| 麻豆国产在线观看一区二区| 成人久久18免费网站| 伊人激情综合| 国产欧美日韩视频一区二区三区| 亚洲精品视频免费看| 九九九精品成人免费视频7| 久久免费精品琪琪| 97成人在线观看| 一级黄色欧美| 欧美日韩国产综合视频在线观看| 亚洲日韩精品综合在线一区二区| 国产亚卅精品无码| 男女男免费视频网站国产| 99热这里只有精品在线观看| 国产精品黄色片| 亚洲香蕉伊综合在人在线| 国产真实乱人视频| 鲁鲁鲁爽爽爽在线视频观看| 色综合综合网| 无码精品国产VA在线观看DVD| 一级香蕉视频在线观看| 亚洲欧美日韩成人高清在线一区| 欧美啪啪网| 亚洲区第一页| 亚洲中文久久精品无玛| 久久99国产精品成人欧美| 九色国产在线| 2020国产精品视频| 日本三级黄在线观看| 日韩精品成人网页视频在线 | 中文成人无码国产亚洲| 在线网站18禁| 日韩无码视频播放| 日韩黄色在线| 日韩美女福利视频| 婷婷综合色| 91网在线| 婷婷亚洲天堂|