1、 學號： 畢業設計（論文）（2015屆）題 目 基于時間序列在糧食產量中的方法研究 學 生 學 院 專業班級 校內指導教師 專業技術職務 校外指導老師 專業技術職務 基于時間序列在糧食產量中的方法研究摘要：糧食是我們生產和生活中的基本消費品，我國民生國計的首要大事就是解決我國的糧食產量問題。本文介紹了幾種時間序列的建模方法，來分析預測出我國的糧食產量問題。并且通過分析我國2000到2014年的糧食生產總量數據的特點，建立了自回歸移動平均模型ARIMA(p,d,q)模型。最后，通過使用Eviews6.0操作軟件成功的計算完成了關于我國糧食產量的預測問題。得出的結果如下，在未來的幾年內，我國糧食生

2、產在不受到自然災害等因素影響的前提下，將會緩慢增長。通過分析，顯示農業科技技術和重大的自然災害對我國糧食產量的影響較為嚴重，為了保證糧食產量的增長要發展好農業技術和做好重大自然災害的預防措施。關鍵詞：時間序列； 糧食產量問題； ARIMA模型。IIStudy on the method of time series based on grain outputAbstract：Food is the most basic consumer goods for human survival, the problem of the grain output of a country relates

3、to the national economy and the national economy and the countrys national economy and the nations livelihood. To make food forecast, this paper presents several modeling methods for time series. Through the analysis of the characteristics of the total grain production data in 2000-2014, the model o

4、f auto regressive moving average ARIMA(p,d,q)is established. Finally, the forecast of grain output of our country is calculated by Eviews6.0 software. The results show that the grain output is not affected by natural disasters in the next several years, and it will slow growth. The analysis shows th

5、at the agricultural technology and the major natural disasters have a serious impact on the grain output in China. In order to ensure food production problems to develop agricultural technology and to do a major disaster prevention.Keyword: Time series ; Grain output ; ARIMA model.II目錄摘要II目錄III1引言11

6、.1課題背景11.2本課題研究的意義12關于我國的糧食產量問題22.1國內糧食產量的現狀22.2研究糧食產量的方法23幾種時間序列預測方法簡介23.1自回歸（AR）模型33.2移動平均（MA）模型43.3自回歸移動平均（ARMA）模型53.4差分自回歸滑動平均(ARIMA)模型54.數據的分析及模型建立74.1數據分析74.2數據平穩化84.3模型的定階94.4模型有效性檢驗124.5模型預測135結論13參 考 文 獻15致 謝17III本科生畢業設計（論文）1引言1.1課題背景我國的民生問題與糧食產量密切相關，國家經濟的可持續性發展以及政治局面的穩定等方面都受到糧食產量的影響。我國有13億

7、的人口，為世界總人口數比例的百分之二十一，但是可用耕地面積卻很少，僅占全球可用的耕地面積比例為百分之七，糧食產量問題是關系到經濟發展的安全，而且是關系到民生和國計的重大的戰略物資問題，每時每刻都不能有一絲一毫的松懈。解決十幾億人口的糧食問題是目前的頭等大事。我國的糧食產量問題受到很多方面要素的影響，并且無法從中找出任何規律。而且有的專家指出，從中長期的發展來看，因為人口的不斷增加，耕地面積在不斷的減少，城市化的步伐的加快，我國人民的生活水平有著顯著的提升。糧食是人們生活中不可或缺的重要必需品，而且需求越來越大，市場上的糧食的供求關系將會發展成是偏緊的關系。糧食始終人類生活中必不可少的特殊的商品

8、，所以我國的一項極為重要的戰略職責就是要保證糧食安全。1.2本課題研究的意義十幾年前，美國世界觀察所的所長布朗就曾說過，我國的人口在逐年大量的增加，但是耕地面積卻在的逐年減少。等到了2030年那年，我國人口將從13億增加到16億，還有由于耕地的減少，糧食產量下降，會導致糧食產量不滿足消耗需求。到時候，人們必然會到各處特別是到國際市場去搶購糧食，這樣的勢頭必然會帶動糧食價格的的上漲，形成糧食恐慌的現象，進而引起世界范圍的糧危機。然而事實在證明：我國的糧食種植面積每削減一個百分點，我國的糧食需求就會要求增添500萬噸的糧食的進口，而我國糧食減少或者增加百分之5的糧食進口量，同時，將導致國際大宗商品

9、市場糧食價格的下跌或上漲約30%，從而影響到近30個發展中國家的經濟利益。于是可以看出中國的糧食產量問題對我國甚至于全世界的糧食價格都有著極為重大的影響。據統計，耕地面積在全球土地總面積中的比重中，印度達百分之56，日本占地百分之12，美國達百分之20。據統計，當前我國的人均耕地僅為12畝（該數字比實際值小，還有待進一步的核查），我國的耕地面積不足世界人均耕地的33.33%。1980年這之后的十五年，我國的耕地減少了近65010000畝，這樣的大區域面積減少相當于縮小了近一個安徽省那么大地域的耕地。在1990年，我國的糧食的播種面積就降低到了十七億畝，在此后越來越少，到了1994年的時候我國糧

10、食的播種面積甚至降到了十六點四億畝，這是第一次比國家制訂的糧食安全警戒線要低。1994年以后，我國的糧食播種面積就一直沒有達到國家制定的安全警戒線。2000年全國糧食播種面積比上年減少了7035萬畝、201年的全國糧食播種面積比上年減少了3575萬畝，還有2002年我國的糧食播種面積比上年減少了3284萬畝。包括受災被損壞和農業布局的調整而占用了的一部分耕地之外，各種各類的建設工程也占用了大批的耕地，這些都導致了可用耕地急速減少。還有退耕還林、退耕還草等這些政府大規模組織的這類的項目，包括了明顯加大種植其他非農作物，可用的耕地面積不斷減少等都直接對糧食播種面積產生了影響，導致可用的糧食播種面積

11、的減少。近幾年來，中國糧食生產的單位面積產量的提高受到以下因素的影響，種植面積波動，農業基礎設施不足，以及自然災害的頻繁發生等。為了提高糧食生產產量和促進糧食的生產的根本就是先提供一套可以鼓勵糧食生產的政策措施,只有提高了糧食種植的效益和增加糧農的收入，才可以更加有效。但是，由于糧食生產的不確定系統的許多因素的制約，未來中國的糧食產量將如何改變，是否能實現國家糧食安全的目標已經成為一個非常有意義的話題。因此，為了政策的調整方向和保障糧食安全，更加有效的分析預測我國的糧食產量問題尤其重要。2關于我國的糧食產量問題2.1國內糧食產量的現狀就當前看來，統計動力學生長模擬法以及遙感技術預測法還有天氣預

12、測法等方法是目前國外較為受歡迎的幾種糧食預測方法。在世界預測糧食生產的問題有很多研究，尤其中國糧食產量預測更是是受到眾多研究者的關注。目前我國學者對于中國糧食產量的預測用的比較多的是時間序列模型的三大類預測模型。其中，自回歸移動平均模型、差分自回歸滑動模型大多用來是處理非平穩的時間序列常用的時間序列模型。另外，自回歸模型（AR模型）、移動平均模型（MA模型）也是時間序列模型中的兩類。雙對數模型是一種被廣發應用的回歸模型，和雙對數模型一樣回歸模型中應用的比較多的還有線性回歸模型。 2.2研究糧食產量的方法從大多的資料看來，基本上時間序列大多不是平穩的，但是假如直接對非平穩的時間序列做回歸分析處理

13、的話，就有可能會造成“偽回歸”的現象。“偽回歸”的現象是指：由于對序列進行操作是，使原本并不存在依賴關系的變量卻出現了依賴關系。本文是根據我國糧食產量的時間序列，做出時間序列的散點圖并分析數據的特點，再根據該序列的自相關函數圖以及偏自相關函數圖，然后判斷這個時間序列的平穩性。然后，平穩化非平穩時間序列的數據，把它當做一個隨機序列，然后根據該序列的特征，建立出一個與之相應的時間序列模型。最后對其進行檢驗，判斷其有效性，再利用該模型進行分析預測我國的糧食產量問題。3幾種時間序列預測方法簡介就現狀而言，數量分析方法應用的較多的方法是時間序列分析。現象是隨著時間而變化的，時間序列分析就是尋找出其中得規

14、律，然后描述出來。一切客觀現象都在不斷的發展變化中，所以對于現象在發展中的變的規律，不光要其相互關聯以及從內部結構以去認識，并且應當隨著時間演化的過程去鉆研，這就需要運用到時間序列分析方法。時間序列分析是處理動態數據數據的統計方法。數理統計方法和隨機過程理論，都是在時間序列分析的基礎上，探索一些隨機數據的統計規律，為解決實際問題的。傳統的統計分析是基于數據序列是獨立的，與時間序列分析的隨機數據序列之間的相關性。后者又可以看成隨機過程統計的一部分。一個單變量時間序列的建立模型和預測的主要方法是使用時間序列建模方法。單變量的時間序列模型(比如說ARIMA模型,季節性ARIMA延長的銷售季節ARIM

15、A和TFN的外生變量) 被廣泛應用于建模和預測糧食產量在不同時間尺度下的常態假設。包含自回歸和移動平均線的組件使得AR(1)模型非常靈活的和受歡迎的時間序列建模范圍廣泛的應用程序。此外,ARIMA模型有很多優勢,如處理序列相關性、與時間相關的變化和預測能力比其他類似的模型。時間序列是一系列的數字序列的時間序列。時間序列分析就是指利用這組序列，使用數理統計方法進行處理，以此來預測該時間序列的現象在未來一段時間內的發展形勢。在時間序列，數據或數據點在時間依賴性，這是數據值的變化，但可能不是時間的嚴格函數。定量預測包括了時間序列分析方法，時間序列分析的基本原理如下：一、對事物發展的連續性的識別。利用

16、過去的數據，可以推斷出事物的發展。二、考慮到事物發展的隨機性。任何事物的發展可能是偶然因素的影響，所以歷史數據應采用加權平均法分析處理。時趨勢變化、周期性變化、隨機變化這是時間序列預測所反映的三類實際變化規律。要使用時間序列預測出一個現象的發展方向，需要進行如下步驟：首先，通過對這一現象的歷史數據的觀測研究，我們發現隨著時間的變化而變化的現象的規律，并預測未來的現象。隨著社會的不斷發展，在影響經濟生活的許多不確定因素的影響越來越大，這必然會吸引人們的注意。時間序列分析的模型有如下：AR模型、MA模型、ARMA模型和ARIMA模型。時間序列分析預測法跟回歸分析之間的差異就是時間序列分析能夠依據一

17、個變量的值，而不添加任何輔助信息來預測自己的變化。3.1自回歸（AR）模型假設時間序列是其前期值和隨機項的線性函數，那么就能夠表示為=+ （1） 那么稱這個時間序列為自回歸序列，式（1）為自回歸模型，記做AR(p)。實參數，是自回歸系數，是該模型的待估參數。隨機項是互相獨立的白噪聲序列，而且該隨機是一個項服從均值是0、方差是的正態分布。該隨機項和滯后變量，是不相關的。不是一般性，在（1）中假設時間序列均值為0.若，則令，可將寫成（1）式的形式。是k步滯后算子，就是，那么模型（1）就可以表示為 =+ （2） 令=1-模型可簡寫為： = （3）那么要使AR（p）過程平穩的前提條件是滯后多項式的所有

18、根都必須在單位圓外，并且的根不小于等1。其中AR（p）過程中經常使用到的是AR（1）過程和AR（2）過程， 保持其平穩性的條件是特征方程 根的絕對值必須大于1，滿足。3.2移動平均（MA）模型假設時間序列是它的當前跟前期的隨機誤差的線性函數，那么可以表示為 （4） 就稱這個時間序列是移動平均序列，（2）式是q階移動平均模型，記為MA（q）模型。實參數，是移動平均系數，是模型的待估系數。引入滯后算子，并令則模型（4）可簡寫為 （5）移動平均過程，即MA過程是無條件的平穩的。然而期望AR過程和MA過程能互相表出，也就是說過程是可逆的。所以要求滯后多項式的所有根均必須在單位圓外，經過推導可得出 （6

19、）其中，其他權重可地推得到。可以說式（6）是MA模型的逆轉形式，這相當于無限階數的AR過程。3.3自回歸移動平均（ARMA）模型假設時間序列是其當前和前期的隨機誤差項和前期值的線性函數，那么能夠表示為： （7）那么這個時間序列為自回歸平均序列，（7）式是（p,q）階的自回歸移動平均模型，記為ARMA(p,q)。，為自回歸系數，為移動平均系數，都是模型的待定參數。那么就可以引入滯后算子B，模型（7）可簡記為 （8）ARMA（p,q）過程的平穩條件就是滯后多項式的根都在單位圓外。可逆條件是的根都在單根圓外。如果=0，那么方程的平穩隨機序列是p階自回歸模型，記為AR(p)模型。很明顯的可以看出，AR

20、(p)模型和MA(q)模型全都是ARMA(p,q)模型的特例。3.4差分自回歸滑動平均(ARIMA)模型差分自回歸滑動平均模型，處理非平穩的時間序列，將其轉化為平穩的時間序列。繼而將因變量對其的滯后值跟進行回歸，同時把隨機誤差項的現值也進行進行回歸，然后建立的模型。將ARIMA模型進行平穩化處理是，依據對時間序列進行回歸處理，然后根據回歸中所含的部份的差異，可以分為以下的幾類：移動平均過程（MA（q）過程）、自回歸過程（AR（p）過程）和自回歸移動平均過程（ARMA過程）以及ARIMA過程。差分自回歸滑動平均模型ARIMA(p，d，q)模型中，每個參數都有不同的意義。其中，AR指的是“自回歸”

21、，p為自回歸的項數；MA指的是“滑動平均”，q為指滑動平均的項數，d是指使它成為平穩序列而做的差分的次數（即階數）。可以這么說，對看做是ARMA(p,q)模型進行擴展處理可以得到ARIMA(p,d,q)模型。ARIMA(p,d,q)模型可以表示為： = （9）其中L是滯后算子 ，。ARIMA模型預測的基本程序數據平穩化處理：第一步，檢驗時間序列數據是平穩的還是非平穩的。并根據其散點圖或者折線圖的分析對該序列進行初步判斷。ADF單位根檢驗一般用來來精確判斷這個序列的平穩性。重復上述過程，直到其成為平穩的序列。這個時候進行微分處理了多少次，這個次數就是ARIMA(p,d,q)模型當中的階數d。從理

22、論上而言，足夠多次的差分運算可以充分地提取序列中的非平穩確定性信息。但應指出是，差分運算取出的階數不應取得太多。由于差分運算是一種信息的提取和處理，所以在實際操作中要適當的選取差分的次數（即差分的階數），以避免出現過差分的現象。一般情況下，差分次數不超過2次。（1）檢驗時間序列數據是否是平穩性的。并根據其散點圖或者折線圖的分析對該序列進行初步判斷。ADF單位根檢驗一般用來確定這個序列是否是平穩的。（2）平穩化處理非平穩序列。假如該時間序列的數據序列不是平穩的，并且還有一定的變化趨勢，那么就要對該數據做差分處理，假如該數據處理時存在異方差，那么就需要技術處理該數據，一直處理到最后得出的數據所做出

23、的自相關函數的值跟偏相關函數的值都是接無限近于零的。（3）處理一個時間序列時，依據既定的識別規則，創建與之對應的時間序列模型。對一個平穩序列而言，如果偏相關函數被截斷，但自相關函數序列顯示拖尾現象，因此可以認為該序列符合AR模型的規則；如果這個序列的偏相關函數是拖尾的，但是該序列的自相關函數卻是顯示截尾現象，那么就可斷定，此序列符合MA模型規則；同時，如果該序列的偏相關函數和自相關函數都顯示拖尾的現象，那么這個序列就符合ARMA模型的規則。一個時間序列的自相關系數趨于零，或偏相關系數趨于零的序列，顯示是拖尾現象，這一過程往往會有很多不同的形式，例如，幾何衰減序列，又比如說有正弦波形式的衰減；而

24、截尾就是指對于一個時間序列而言，從某一階數后的自相關系數為0，又或者該序列的偏相關系數為0。（4）進行參數估計。（5）診斷該時間序列的殘差序列是否為白噪聲（可以通過假設檢驗）。（6）預測分析是通過以檢驗的模型進行的。參數的判斷中AIC準則：AIC信息準則就是Akaike information criterion，它是一種衡量統計模型擬合優良性的準則，又由于它是日本統計學家赤池弘次所創立和發展的，所以又稱為赤池信息量準則。它可以權衡所估計模型的復雜度和此模型擬合數據的優良性。參數估計是指確定該時間模型的模型階數后，然后對其對與之對應的ARMA模型做參數估計。本文采是采用最小二乘法進行參數估計，

25、但要提醒的是，對移動平均模型（MA模型）進行參數估計比較困難，所以我們應當選取較低階數的自回歸移動模型。避免出現過度差分而引起的誤差過大的現象。模型試驗和模型參數估計的識別后，估計的結果應該被診斷和測試，以找出合適的模型。若不合適，應該知道下一步作何種修改。這一階段主要檢驗擬合的模型是否合理。首先對該模型的參數做出顯著性檢驗，其次對該模型是否為白噪聲進行檢驗。參數估計的意義是Q檢驗，t檢驗完成零假設。該模型的誤差項是一個白噪聲過程。4.數據的分析及模型建立很多因素都會影響到我國的糧食生產，比環境、人力和土地資源、農業科學技術以及發展戰略方向等，這些因素之間又存在著錯綜復雜的關系。所以很難分析和

26、利用結構的因果模型，預測糧食產量的問題。將往年的糧食產量整理做時間序列，并找出過去數據的變化規律，然后建立預測模型，并用此來預測我國未來的糧食產量的發展變化，這有著非常重要的意義。下面以我國20002014年我國糧食總產量的數據(見下表)為例，并用時間序列對下面數據進行數據分析處理的過程做出介紹。 4.1數據分析 20002014年我國糧食產量（單位：萬噸）如下表：年份糧食產量年份糧食產量200046,217.52200852,870.92200145,263.67200953,082.08200245,705.75201054,647.71200343,069.53201157,120.85

27、200446,946.95201258,957.97200548,402.19201360,193.84200649,804.23201460,709.90200750,160.28注：上表的數據取自中國統計局官網在ARMA模型中，一個平穩的時間序列可以由一個均值為0的平穩隨機過程發展出來，就是指在任何時間段這個序列的隨機性不變，在圖形上的表現是，在時間序列的所有樣本點都在一定水平的隨機波動。對于非平穩時間序列，需要預先對時間序列進行平穩化處理。建立時間序列模型之前需要檢驗序列的平穩性，只有平穩序列才能建立時間序列模型。時間序列的ADF檢驗是通過Eviews軟件進行數據分析，可以判斷其穩定性。

28、根據測試結果，如果得出的測試值的絕對值大于臨界值的絕對值，那么就可以判斷出，該序列就是平穩序列；否則，為非平穩數據。利用eviews軟件做圖，畫出我國糧食總產量的曲線圖，如下：從這個圖，我們可以清晰地看到，糧食產量總體呈上升趨勢，但在2002年到2004年有了下跌的出現，這是對非平穩序列的直觀表現形勢。4.2數據平穩化糧食產量的時間序列是一階差分的單位根檢驗，其結果是： 圖：糧食產量時間序列單根檢驗 由上圖的測試結果可以看出，單根為0.609816明顯大于10%水平以下的臨界值，所以該時間序列存在單位根，故而時間序列是非平穩的時間序列。4.3模型的定階 ARIMA模型的參數的確定有不同的定階方

29、法，自相關和偏相關函數定階法、FPS準則以及AIC和SC準則是三種廣泛使用的方法。下面是對我國糧食產量時間模型的定階步驟，如下： 圖：糧食產量一階差分單位根檢驗結果圖 通過圖，可以看出糧食產量異界差分后得到的序列在某一常數附近波動，可初步識別序列已平穩。并且該序列的ADF檢驗值為-5.382102，分別小于1%、5%和10%的檢驗水平的臨界值，所以它通過了ADF檢驗，可以看出該序列是平穩序列。在此，可以看出一階差分就已可以消除序列的非平穩性，因此在ARIMA(p，d，q)模型中可以確定：d=1。通過觀察平穩時間序列的自相關函數的性質，相關函數的性質決定階數。做出糧食產量一階差分序列D(x)的自

30、相關偏相關圖，如下：分別取幾個p，q的值進行檢驗，所得到的AIC值，然后根據AIC值做出判斷。p，q值AIC值p=1,q=117.25400p=1,q=216.50130p=1,q=316.49318p=2,q=117.28908p=2,q=216.61116p=2,q=316.02550p=3,q=116.23902p=3,q=216.04907p=3,q=316.39271從上表，其中p=2，q=3時，AIC值最小。以此，可以得出：p=2 q=3所以，我們可以確定模型參數，建立模型ARIMA(2,1,3)4.4模型有效性檢驗 通過利用沒使用過的時間序列的觀測值判斷模型的預測能力，這是模型的

31、有效性檢驗。利用一部分已有的歷史數據，對其進行回歸并做出預測。最后將回歸預測處理后得到的結果跟實際值做比較，這樣可以簡單對模型預測的有效性進行檢驗。 下面是20102014年的糧食產量估計值和真實值得比較及相對誤差表：20102014年估計值與真實值及相對誤差年份真實值（萬噸）估計值（萬噸）誤差（%）201054,647.7155348.001.3201157,120.8556757.450.64201258,957.9758241.201.2201360,193.8459712.600.80201460,709.9061158.11 0.74 由上表可以看出相對誤差最高為1.3%，都小于5%

32、，可以判斷出預測結果比較準確，能夠基本擬合實際值。4.5模型預測 由此模型對2015-2017年糧食產量進行預測，結果如下：表：20152017年我國糧食產量預測年份201520162017預測產量（萬噸）62607.9264066.7665524.09由上圖預測的結果我們可以看出，我國的糧食總產量在將來的幾年內將呈現增長趨勢。5結論 正常情況下，時間序列模型只適合進行短期預測，如果進行長期預測，就可能出現錯誤的預測結果。所以模型只能對短期未來幾年，中國的糧食產量預測。其次，根據采集的歷史數據可以看出，我國糧食產量數據特點是基本上是呈上升趨勢，但在2000年到2002年卻出現了下降的情況，根據

33、查詢，原因是糧食種植總面積有所減少，國家當時的退耕還林、退耕還草的項目的實施。同時農民當時收入較低，各地的農民為了獲得更高的利益增加收入而播種了利潤更高的經濟作物。導致當時糧食總產量的下降，以及突然發生的自然災害，遭受到了嚴重干旱的襲擊，這也是全年糧食總產量可能下降的原因。因此，發展現代農業，根本出路在科技。為了確保國家未來的糧食安全，進行農業科技發展，鼓勵加強創新力度。大力支持農業的發展。在中國的西部，是一個農業大省，通過不斷的努力，加強農業科技創新，培養了一批抗寒、抗倒伏獲得作物，這些作物適合本地生長。加強污染管理力度，特別是土地污染情況，進行土壤改良，提高土地利用率；大力發展農業裝備優化

34、，這樣可以大大解放了勞動力，提高了生產率；國家還應鼓勵農業科技的發展，培養了一大批農業科技人才，提高農業生產的軟實力。要使中國糧食產量呈現持續增長，除了科技進步，還要加強預防自然災害。參 考 文 獻1 劉羅曼時間序列分析中指數平滑法的應用J沈陽師范大學學報(自然科學版)，2009，27(4)：416-4182 齊雨藻，黃偉建，邱漩鴻大鵬灣夜光藻種群動態的時間序列模型J暨南大學學報(自然科學與醫學版)，1991，12(3)：961033 王燕應用時間序列分析(第三版)M北京：中國人民大學出版社20124 龍建輝，賀向陽港口投資、港口物流與經濟增長 來自寧波市19852010年的時間序列證據J經濟

35、師，2012(4)：2l22l65 李瑞瑩，康銳基于ARMA模型的故障率預測方法研究J系統工程與電子技術，2009，30(8)：15881591.6 賈朝龍，徐維祥，王福田，等基于GM (1，1)與AR模型的軌道不平順狀態預測J北京交通大學學報：自然科學版，2012，36(3)：52567 Jiang Y，Feng LTransport-related resource and environmental issues in China. Special EditionJTransportation in China，2006，7(8)：98 Hou Q，An XQ, Guo JPAn eva

36、luation of resident exposure to respirable particulate matter and health economic loss in Beijing during Beijing 2008 Olympic GamesJSci. Total Environ，2010，408(19)：239China State Environmental Protection AdministrationJChina State Environmental Protection Administration，2012，21(5)：95-9710 Anderson HR，Air pollution and mortality： a historyJ Atmos. Environ. 2009；43(1)：14215211 Yu GH，Chen HL, Wen WC (2002) A distribution-free method for forecasting non-Gaussian time seriesJStoch Environ Res Risk A 16(2)：10111112 Yu GH，Huang CC (2001) A distribution free plotting pos