韓躍，汪春

（黑龍江八一農(nóng)墾大學(xué)工程學(xué)院，大慶163319）

東北地區(qū)屬東亞季風(fēng)氣候區(qū)，地處中高緯歐亞大陸的東岸，是中國(guó)的老工業(yè)基地和最大的商品糧食基地。由于全球氣候變暖，東北地區(qū)夏季氣溫不僅僅表現(xiàn)在平均氣溫的升高，氣候變率也在增大，自1986年后，東北地區(qū)夏季氣溫升高顯著[1-3]。尤其近兩年來(lái)，東北地區(qū)夏季7月份平均氣溫達(dá)24℃以上，最高氣溫可達(dá)30～34℃。而作物生長(zhǎng)的適宜溫度在28℃以下，當(dāng)外界氣溫達(dá)到25℃時(shí)，溫室自然通風(fēng)情況下的氣溫一般都在32℃以上，遠(yuǎn)遠(yuǎn)超過(guò)了作物生長(zhǎng)的適宜溫度[4]。

Venlo型溫室，是一種源于荷蘭的小屋面雙坡面玻璃溫室[5]，我國(guó)為發(fā)展設(shè)施農(nóng)業(yè)建設(shè)，已在許多地區(qū)引進(jìn)應(yīng)用。近年來(lái)東北地區(qū)夏季氣溫升高，如何經(jīng)濟(jì)有效的降低溫室內(nèi)的溫度是保證作物安全越夏生長(zhǎng)并提高溫室生產(chǎn)效率的重要環(huán)節(jié)。

為了更好地提高東北夏季溫室的生產(chǎn)效率，通過(guò)測(cè)定溫室內(nèi)空氣溫度及其影響因素?cái)?shù)據(jù)，對(duì)Venlo型玻璃溫室溫度建立了時(shí)間序列模型，可對(duì)溫室溫度進(jìn)行預(yù)測(cè)，可為夏季東北Venlo型玻璃溫室生產(chǎn)提供理論依據(jù)，提高農(nóng)業(yè)信息化水平[6]。

1 材料與方法

1.1 實(shí)驗(yàn)材料

試驗(yàn)以黑龍江八一農(nóng)墾大學(xué)的Venlo型玻璃溫室為研究對(duì)象，運(yùn)用JTJW-2建筑熱工溫度與熱流巡回檢測(cè)儀，測(cè)溫傳感器為銅-康銅熱電偶。Venlo型玻璃溫室南北走向，東西長(zhǎng)24m，雙屋脊跨度為8.4m，頂高5m。Venlo型玻璃溫室具備交錯(cuò)式頂開(kāi)窗，雙跨具備6扇頂窗，每扇小窗長(zhǎng)3 m，寬1 m，頂窗關(guān)閉時(shí)與地面夾角約為30度，溫室頂部覆蓋材料為4mm厚浮法玻璃，四周覆蓋材料為6mm厚浮法玻璃。Venlo型玻璃溫室南面?zhèn)葔ε鋫溥B續(xù)型鋁合金側(cè)開(kāi)外翻窗，長(zhǎng)8m，寬1m，室內(nèi)配備長(zhǎng)8m，高1.5 m的濕簾機(jī)構(gòu)，東墻配備兩個(gè)風(fēng)機(jī)，長(zhǎng)寬均為1.4 m。系統(tǒng)布置了18個(gè)均勻的溫度測(cè)點(diǎn)，東西走向劃分了3個(gè)平面，一個(gè)面上各安置6個(gè)測(cè)點(diǎn)，其中縱向測(cè)點(diǎn)分別距地1.2m和2.4m。大棚結(jié)構(gòu)及室溫測(cè)量的布點(diǎn)位置如圖1。試驗(yàn)時(shí)間為2011年7月9日-2011年8月24日，主要采集溫室內(nèi)的空氣溫度，及氣象站采集的室外溫度，風(fēng)速，濕度，太陽(yáng)輻射，記錄時(shí)間間隔為1 h。期間無(wú)遮陽(yáng)并保持自然通風(fēng)，即頂窗和側(cè)窗保持開(kāi)啟，風(fēng)機(jī)每天8：00～17：00開(kāi)啟，并根據(jù)天氣情況，手動(dòng)開(kāi)關(guān)濕簾降溫系統(tǒng)，運(yùn)用LM8000四合一風(fēng)速儀檢測(cè)風(fēng)機(jī)風(fēng)速及通過(guò)濕簾后的濕簾表面空氣溫度。

圖1 溫室結(jié)構(gòu)及布點(diǎn)位置圖Fig.1 Locationmap of the greenhouse structure and the distribution point

1.2 建模方法

溫室環(huán)境控制系統(tǒng)以等時(shí)間間隔采集到的溫度及其他各變量數(shù)據(jù)，是按照時(shí)間次序排列的一系列數(shù)據(jù)，這些數(shù)據(jù)記錄序列是多組獨(dú)立的時(shí)間序列[7]。因此對(duì)Venlo型玻璃溫室溫度的預(yù)測(cè)模型建立，采用時(shí)間序列分析方法。具體步驟如下[8-9]:（1）對(duì)通過(guò)實(shí)驗(yàn)測(cè)得的溫室溫度及其影響因素的原始數(shù)據(jù)序列，進(jìn)行線性差值法處理異常、缺失的數(shù)據(jù)。（2）溫室無(wú)遮陽(yáng)下三種降溫狀態(tài)，即自然通風(fēng)狀態(tài)、自然通風(fēng)+風(fēng)機(jī)排風(fēng)狀態(tài)、自然通風(fēng)+風(fēng)機(jī)排風(fēng)+濕簾降溫狀態(tài)，分別對(duì)各狀態(tài)下的溫室溫度、HOBO氣象站數(shù)據(jù)、自然通風(fēng)率、風(fēng)機(jī)開(kāi)啟時(shí)的總通風(fēng)率、濕簾開(kāi)啟時(shí)經(jīng)過(guò)濕簾表面的空氣溫度對(duì)序列進(jìn)行ADF檢驗(yàn)（單位根檢驗(yàn)），檢驗(yàn)是否為平穩(wěn)序列。如果不是，則再進(jìn)行一次一階差分，然后再進(jìn)行平穩(wěn)性檢驗(yàn)；如果是平穩(wěn)序列，則根據(jù)自相關(guān)系數(shù)估計(jì)值和偏相關(guān)系數(shù)估計(jì)值選擇模型類(lèi)型。（3）對(duì)差分后的平穩(wěn)時(shí)間序列進(jìn)行協(xié)整檢驗(yàn)并進(jìn)行模型預(yù)估。（4）對(duì)估計(jì)的預(yù)測(cè)模型進(jìn)行誤差修正。（5）進(jìn)行模型診斷檢驗(yàn)和預(yù)測(cè)精度評(píng)價(jià)。

2 溫室溫度預(yù)測(cè)模型構(gòu)建

2.1 異常與缺失數(shù)據(jù)的補(bǔ)足

異常點(diǎn)（outliner）是指一個(gè)時(shí)間序列觀測(cè)中，遠(yuǎn)離序列一般水平的極端大值和極端小值，有時(shí)稱(chēng)為離群點(diǎn)或野值。將異常點(diǎn)去掉，便形成了缺省值，而由于數(shù)據(jù)采集時(shí)，也會(huì)由于儀器故障，操作失誤等原因造成缺省值[10]。異常點(diǎn)檢驗(yàn)與缺省值補(bǔ)足方法如下：

假定正常數(shù)據(jù)的序列值｛Xt｝是平滑的，而異常點(diǎn)是突變的，若

2.2 平穩(wěn)性檢驗(yàn)

室外溫度TWw（℃）、室外太陽(yáng)輻射強(qiáng)度RWw（w/m2）、自然通風(fēng)速率V1（m3/s2）、風(fēng)機(jī)開(kāi)啟時(shí)的總通風(fēng)速率V2（m3/s2）、室外空氣相對(duì)濕度SWw（%）以及通過(guò)濕簾降溫的濕簾表面空氣溫度Ts（℃），這些是影響室內(nèi)空氣溫度T（℃）的主要因子，其他因素均包含在隨機(jī)誤差項(xiàng)δ中引入模型，Venlo型玻璃溫室無(wú)遮陽(yáng)下三種降溫狀態(tài)，即自然通風(fēng)狀態(tài)、自然通風(fēng)+風(fēng)機(jī)排風(fēng)狀態(tài)、自然通風(fēng)+風(fēng)機(jī)排風(fēng)+濕簾降溫狀態(tài)，分別可以表示為：

其中Venlo型玻璃溫室自然通風(fēng)狀態(tài)下的通風(fēng)速率[11]（V1，m3/s2）為：

式中，h為進(jìn)風(fēng)口中心和出風(fēng)口中心的垂直距離（m），TW為室外空氣溫度（℃）；ΔT為室內(nèi)外溫度差（ΔT=TWW-T），TW為室內(nèi)溫度（℃）;uW為室外風(fēng)速；g為重力加速度（m/s2）；Cd，Cu分別為流量系數(shù)和綜合風(fēng)效系數(shù)，分別為0.64，0.09；St（m2）為頂窗開(kāi)啟面積（Sh）和側(cè)窗開(kāi)啟面積（Sl）之和，φ因子表示為：

式中b=Sh/Sl=2.25，經(jīng)計(jì)算φ≈1，溫室內(nèi)自然通風(fēng)速率表示為：

當(dāng)室外風(fēng)速大于2m·s-1時(shí)，溫室內(nèi)自然通風(fēng)速率可表示為：

試驗(yàn)期間，風(fēng)機(jī)排風(fēng)時(shí)間均在8：00～17：00間，所以期間總通風(fēng)率為風(fēng)機(jī)排風(fēng)速率與自然通風(fēng)速率之和，其中風(fēng)機(jī)排風(fēng)速率約為每單位秒下的風(fēng)機(jī)面積與風(fēng)機(jī)風(fēng)速的乘積，經(jīng)測(cè)試，風(fēng)機(jī)風(fēng)速約為6 m·s-1，所以計(jì)算得，風(fēng)機(jī)排風(fēng)速率為V2=9.23m3/s2。

式（3）中各變量若均為平穩(wěn)序列，則根據(jù)自相關(guān)系數(shù)估計(jì)值和偏相關(guān)系數(shù)估計(jì)值選擇模型類(lèi)型。如果各變量是不平穩(wěn)序列，對(duì)非平穩(wěn)的時(shí)間序列建模會(huì)產(chǎn)生偽回歸現(xiàn)象[12]，因此需要對(duì)各變量進(jìn)行平穩(wěn)性檢驗(yàn)。

如果，經(jīng)檢驗(yàn)非平穩(wěn)的時(shí)間序列Zt的自回歸模型中含有單位根，其自回歸模型為：

式（8）中λ是序列Zt的重要構(gòu)成部分，εt是具有方差為σε2和零均值的獨(dú)立同分布的序列，Φi為Zt的自回歸系數(shù)。方程（4）可簡(jiǎn)化為3個(gè)模型:①λ具有常數(shù)項(xiàng)C（model“C，NT”）；②λ具有常數(shù)項(xiàng)C，趨勢(shì)項(xiàng)T（model“C，T”）；③Zt不含λ（model“NC，NT”）。

表1 單位根檢驗(yàn)的臨界Table1 Criticality of the unit root test

運(yùn)用eviews進(jìn)行ADF單位根檢驗(yàn)結(jié)果如表1，表2。經(jīng)結(jié)果顯示：各變量檢驗(yàn)統(tǒng)計(jì)量值均有小于相應(yīng)臨界值，從而拒絕原假設(shè)，均不含單位根，即均是平穩(wěn)的時(shí)間序列。

表2 單位根檢驗(yàn)結(jié)果（t-統(tǒng)計(jì)量）Table2 Results of the unit root test（t-statistic）

2.3 選擇模型類(lèi)型

2.3.1 自相關(guān)檢驗(yàn)和模型識(shí)別

溫室自然通風(fēng)狀態(tài)下，根據(jù)方程回歸殘差的相關(guān)圖和Q統(tǒng)計(jì)量，判斷出ARMA模型的自回歸過(guò)程可能是2階，移動(dòng)平均過(guò)程MA是低階的。因此，估計(jì)下列兩種模型形式：ARMA（2，1）和ARMA（2，2）。

表3 方程回歸殘差的相關(guān)圖和Q統(tǒng)計(jì)量Table3 Correlogram of equation of regression residuals and the Q statistic

2.3.2 ARMA模型估計(jì)

識(shí)別ARMA模型的形式后，可以使用EViews估計(jì)方程參數(shù)。模型整體上是顯著的，且擬合效果良好。模型ARMA（2，1）的AIC準(zhǔn)則和SC準(zhǔn)則都小于ARMA（2，2）。因此認(rèn)為模型ARMA（2，1）較好。

表4 ARMA（2，1）模型估計(jì)結(jié)果Table4 Results of ARMA（2，1）model estimation

利用滯后多項(xiàng)式寫(xiě)出自然通風(fēng)狀態(tài)下的溫室室溫模型ARMA（2，1）的估計(jì)結(jié)果：

（1-0.817 8L+0.181 4L2）T=10.562+0.537 7Tw+0.012 2Rw+0.02Sw+0.004 6V1+（1-0.104L）ε （9）

2.3.3 模型診斷檢驗(yàn)

運(yùn)用eviews軟件包對(duì)所估計(jì)的ARMA（2，1）模型生成的殘差進(jìn)行自相關(guān)檢驗(yàn)。殘差序列的樣本自相關(guān)函數(shù)都在95%的置信區(qū)域以內(nèi)，從滯后1階至12階的自相關(guān)函數(shù)相應(yīng)的概率值P都大于0.9，遠(yuǎn)大于檢驗(yàn)水平0.05，因此模型ARMA（2，1）估計(jì)結(jié)果的殘差序列不存在自相關(guān)。

同自然通風(fēng)狀態(tài)下的溫室溫度預(yù)測(cè)模型方法和步驟相同，建立自然通風(fēng)+風(fēng)機(jī)通風(fēng)狀態(tài)下溫室室溫預(yù)測(cè)模型為：

（1-0.528 7L）T=19.593 8+0.340 9Tw+0.005 3Rw+0.004 1Sw-0.002 2V2+（1-0.093 2L）ε （10）

自然通風(fēng)+風(fēng)機(jī)通風(fēng)+濕簾降溫狀態(tài)下的溫室室溫預(yù)測(cè)模型如下:

（1-0.719 5L）T=10.935 6+0.585 5Tw+0.004 2Rw+0.006 4S-0.004 5V2+0.069 2T+（1-0.147L）ε （11）

3 模型預(yù)測(cè)與檢驗(yàn)

根據(jù)自變量RW、TW、SW、V、TS的數(shù)據(jù)，預(yù)測(cè)三種狀態(tài)下7 d的因變量溫室室溫T數(shù)據(jù)。取7月9日10：00-7月16日10：00的數(shù)據(jù)為已知條件進(jìn)行預(yù)測(cè)數(shù)值T′與實(shí)際值T比較，部分?jǐn)?shù)據(jù)如表5顯示，溫室溫度預(yù)測(cè)值的最大絕對(duì)誤差、最大相對(duì)誤差、平均絕對(duì)誤差和平均相對(duì)誤差分別為0.9℃、4.4%、0.1℃和0.1%，三種模型均預(yù)測(cè)效果比較理想。

綜上所述，預(yù)測(cè)模型通過(guò)了理論驗(yàn)證和實(shí)際數(shù)據(jù)檢驗(yàn)，可以用于溫室溫度的預(yù)測(cè)。

表5 預(yù)測(cè)值與實(shí)際值誤差表Table5 The error table of the predicted and actual value

4 結(jié)論

4.1 異常點(diǎn)檢驗(yàn)和插值法對(duì)實(shí)驗(yàn)測(cè)量的數(shù)據(jù)修正融合，實(shí)現(xiàn)了溫室空氣溫度的精確測(cè)量。

4.2 通過(guò)對(duì)模型各變量進(jìn)行ADF單位根檢驗(yàn)，表明各變量均為平穩(wěn)時(shí)間序列?；谏鲜鎏卣鞑捎米韵嚓P(guān)分析，然后進(jìn)行模型預(yù)估，并通過(guò)引入殘差自相關(guān)檢驗(yàn)確保模型的預(yù)測(cè)精度。

4.3 上述結(jié)果表明所建室內(nèi)空氣溫度時(shí)間序列模型具有較高的預(yù)測(cè)精度，研究結(jié)果為全開(kāi)型玻璃溫室夏季的溫度控制提供了良好的理論依據(jù)。

[1]李菲，李輯，管兆勇.我國(guó)東北夏季氣溫年代際變化特征及與太平洋海溫異常關(guān)系的研究[J].氣象與環(huán)境學(xué)報(bào)，2010，3（26）：19-26.

[2]李蘭，王盤(pán)興，陳長(zhǎng)勝，等.東北夏季（6—8月）氣溫異常的時(shí)空特征分析[J].南京氣象學(xué)院學(xué)報(bào)，2005，6（28）：801-806.

[3]李輯，龔強(qiáng).東北地區(qū)夏季氣溫變化特征分析[J].氣象與環(huán)境學(xué)報(bào)，2006，22（1）：6-10.

[4]J.J.編著.溫室管理[M].鄭光華，譯.北京：科學(xué)出版社,1984.

[5]李文榮.論設(shè)施農(nóng)業(yè)的創(chuàng)新與發(fā)展[J].農(nóng)機(jī)化研究，2007（8）：183-186.

[6]王熙，莊衛(wèi)東，汪春.黑龍江墾區(qū)農(nóng)業(yè)信息化進(jìn)展及展望[J].黑龍江八一農(nóng)墾大學(xué)學(xué)報(bào)，2010，22（5）：18-21.

[7]Peter J Brockwel,l Richard A Davis.Time series：theory andmethods[M].New York：Springer-Verlag,1991.

[8]潘迪特SM，吳憲民.時(shí)間序列及系統(tǒng)分析與應(yīng)用[M].北京：機(jī)械工業(yè)出版社，1988.

[9]BoxG E P.時(shí)間序列分析：預(yù)測(cè)與控制[M].顧嵐,譯.北京：中國(guó)統(tǒng)計(jì)出版社，1997.

[10]肖枝洪，郭明月.時(shí)間序列分析與SAS應(yīng)用[M].武漢：武漢大學(xué)出版社，2009.

[11]H Temrati，TBoulard，A Bekkaoui．Natural ventilation and microclimatic performance of a large scale Banana Green-house[J].Journal of Agricultural Engineering Research.2001，80（3）：261-271．

[12]高鐵梅．計(jì)量經(jīng)濟(jì)分析方法與建模[M].北京：清華大學(xué)，2004.