概率統(tǒng)計(jì)方法模型上

1、第6章概率統(tǒng)計(jì)方法模型在對(duì)實(shí)際問(wèn)題進(jìn)行數(shù)學(xué)建模的過(guò)程中，人們經(jīng)常遇到隨機(jī)性的不確定問(wèn)題，用傳統(tǒng)的數(shù)學(xué)建模方法難以解決。此時(shí)，就需要基于概率論和數(shù)理統(tǒng)計(jì)知識(shí)，運(yùn)用概率統(tǒng)計(jì)的方法建立數(shù)學(xué)模型，對(duì)實(shí)際問(wèn)題進(jìn)行求解，揭示事物發(fā)展的基本規(guī)律。本章詳細(xì)介紹用概率統(tǒng)計(jì)方法建模的基本思路，結(jié)合實(shí)際的案例，指出如何用隨機(jī)變量和概率分布來(lái)描述隨機(jī)不確定事件，說(shuō)明求解概率統(tǒng)計(jì)類模型的一般過(guò)程，并指出該類數(shù)學(xué)模型在社會(huì)調(diào)查、影響因素分析、發(fā)展趨勢(shì)模擬等方面的廣泛應(yīng)用。6.1概率模型與Monte Carlo模擬概率模型（1）傳染病隨機(jī)模型在各種傳染病的流行過(guò)程中，無(wú)論健康人還是病人，任何兩個(gè)人之間接觸的機(jī)會(huì)都是隨機(jī)的，

2、而且當(dāng)健康人與病人接觸時(shí)，健康人是否被傳染也是一個(gè)隨機(jī)的事件。我們通過(guò)建立傳染病隨機(jī)模型來(lái)分析這些隨機(jī)規(guī)律。假設(shè)人群總的規(guī)模為N，在總?cè)巳褐?，病人的?shù)量為m，健康人的數(shù)量為s，即滿足N=m+s。在人們的日常生活中，任意兩人之間（包括健康人和病人）接觸的概率相同，每人平均與k個(gè)人接觸。當(dāng)健康人和病人接觸時(shí)，被傳染的概率為p。在以上假設(shè)的參數(shù)中，m和s通常是已知的，k和p可以通過(guò)專家的經(jīng)驗(yàn)和統(tǒng)計(jì)數(shù)據(jù)獲得。我們分析的目的是尋找健康人群中每天平均被感染的人數(shù)與已知參數(shù)之間的關(guān)系，以及初始參數(shù)對(duì)傳染病的擴(kuò)散速度和流行趨勢(shì)的影響。我們首先以每一名健康人為研究對(duì)象，探討其每天被感染的概率，而每一名健康人被一

3、名指定病人接觸并傳染的概率等于每名健康人與指定傳染者接觸的概率乘以接觸時(shí)感染的概率。記人群中任意兩人接觸的概率為q，則對(duì)每一名健康人來(lái)說(shuō)，其每天接觸的人數(shù)服從二項(xiàng)分布，分布函數(shù)為, （6.1.1）這個(gè)分布的期望為k，即，進(jìn)而。這樣，一名健康人被一名指定病人接觸并感染的概率為.進(jìn)一步，對(duì)人群中的每一名健康人來(lái)說(shuō)，其每天不被感染的概率為，被感染的概率為. （6.1.2）所以，對(duì)人群中的所有健康人來(lái)說(shuō)，每天被感染的人數(shù)服從二項(xiàng)分布，分布函數(shù)為, （6.1.3）每天被感染的人數(shù)期望為，標(biāo)準(zhǔn)差為為了得到簡(jiǎn)明的結(jié)果，對(duì)進(jìn)行近似計(jì)算，由于通常人群的總數(shù)，且根據(jù)Talyor展開(kāi)，得，因此，. （6.1.4）通

4、過(guò)式（6.1.4）可以看出平均每天被感染的人數(shù)與s、m、p和k之間的關(guān)系。進(jìn)而可以度量平均每天被感染人數(shù)的相對(duì)誤差即 （6.1.5）由式（6.1.4）可以看出，對(duì)于健康人群來(lái)說(shuō)，每天平均被感染的人數(shù)與人群中每人每天平均接觸的人數(shù)k，健康人與病人接觸時(shí)被感染的概率p成正比。當(dāng)n，p，k都確定的情況下，時(shí)，也就是在整個(gè)人群中，病人和健康人的數(shù)量各占一半時(shí)，每天被感染的人數(shù)達(dá)到最大。為了對(duì)傳染病的傳染過(guò)程有一個(gè)直觀的了解，假設(shè)一個(gè)人口總量n=10000的人群，在日常生活中，平均每人每天接觸的人數(shù)k=18，健康人與病人接觸時(shí)被感染的概率p=10%，對(duì)于不同的m，平均每天被感染人數(shù)與相對(duì)誤差的變化趨勢(shì)如

5、所示。可見(jiàn)被感染人數(shù)隨著病人數(shù)量的增大而增大，直到病人數(shù)量占總?cè)巳簲?shù)量的一半時(shí)達(dá)到最大，隨后呈下降趨勢(shì)。隨著病人人口的增加每天被感染人數(shù)的相對(duì)誤差一直呈減少趨勢(shì)，尤為明顯的是病人數(shù)量增長(zhǎng)的前期，相對(duì)誤差急劇減少。圖6.1.1 平均每天被感染人數(shù)的趨勢(shì)圖 平均每天被感染人數(shù)的相對(duì)誤差趨勢(shì)R編程如下：crb - function(m, n=10000, p=0.1, k=18) #函數(shù) u-(m*(n-m)*p*k)/(n-1);u m-1:10000plot(1:10000,crb(m), xlab=m, ylab=平均每天被的傳染人數(shù),type=l, col=blue) crb1 - func

6、tion(s, n=10000, p=0.1, k=18) #相對(duì)誤差函數(shù) miugama-(n-1-m*p*k)/(n-m)*m*p*k)0.5;miugama m-1:6000plot(1:6000,crb1(m), xlab=m, ylab=相對(duì)誤差,type=l, col=red) （2）企鵝繁殖模型企鵝的繁殖過(guò)程是一個(gè)典型的隨機(jī)不確定模型。首先，每只母企鵝下蛋的數(shù)量是隨機(jī)的，服從泊松分布，其次，每個(gè)企鵝蛋是否可以成功孵化也是不確定的。針對(duì)這一問(wèn)題，我們?cè)诤侠砑僭O(shè)的基礎(chǔ)上，建立概率模型，求企鵝后代個(gè)數(shù)的期望值。根據(jù)人們的統(tǒng)計(jì)，企鵝生蛋的個(gè)數(shù)是服從參數(shù)為的泊松分布，即 （6.1.6）而每

7、個(gè)生蛋能發(fā)育成企鵝的概率為p，且每個(gè)生蛋能否發(fā)育成企鵝是彼此獨(dú)立的隨機(jī)事件。令代表企鵝后代的個(gè)數(shù)，且有?？梢?jiàn)取非負(fù)的整數(shù)值0, 1, 2, ，對(duì)于的概率，我們利用全概率公式 (6.1.7)注意到每個(gè)生蛋發(fā)育成小企鵝是相互獨(dú)立的，且發(fā)育成小企鵝的概率為p，因此，實(shí)際上反映了有k個(gè)生蛋，每個(gè)生蛋獨(dú)立發(fā)育，恰好發(fā)育成 k 個(gè)企鵝的概率。顯然，它是一個(gè)伯努利試驗(yàn)，因而 (6.1.8)于是有得到企鵝后代的個(gè)數(shù)服從參數(shù)為的泊松分布，從而企鵝后代個(gè)數(shù)的期望值為。也說(shuō)明了企鵝后代的個(gè)數(shù)與生蛋的個(gè)數(shù)以及發(fā)育成功的概率成正比。Monte Carlo模擬Monte Carlo（蒙特卡洛）模擬，也稱統(tǒng)計(jì)模擬方法。該方

8、法是上世紀(jì)40年代，由John von Neumann（馮諾依曼），Stanislaw Ulam和Nicholas Metropolis在洛斯阿拉莫斯國(guó)家實(shí)驗(yàn)室進(jìn)行核武器計(jì)劃的工作時(shí)發(fā)明的，后來(lái)該方法的得名是由于Ulam的叔叔常在馳名世界的賭城（摩納哥的Monte Carlo）輸錢(qián)。事實(shí)上，Monte Carlo模擬是由于科學(xué)技術(shù)的發(fā)展和電子計(jì)算機(jī)的發(fā)明，而被提出的一種以概率統(tǒng)計(jì)理論為指導(dǎo)的非常重要的數(shù)值計(jì)算方法。該方法是一種使用隨機(jī)數(shù)來(lái)解決很多計(jì)算問(wèn)題的方法。目前，蒙特卡羅方法在金融工程學(xué)，宏觀經(jīng)濟(jì)學(xué)，生物醫(yī)學(xué)，計(jì)算物理學(xué)（如粒子輸運(yùn)計(jì)算、量子熱力學(xué)計(jì)算、空氣動(dòng)力學(xué)計(jì)算）等領(lǐng)域應(yīng)用廣泛。通常

9、蒙特卡羅方法可以粗略地分成兩類：一類是所求解的問(wèn)題本身具有內(nèi)在的隨機(jī)性，借助計(jì)算機(jī)的運(yùn)算能力可以直接模擬這種隨機(jī)的過(guò)程。另一種類型是所求解問(wèn)題可以轉(zhuǎn)化為某種隨機(jī)分布的特征數(shù)，比如隨機(jī)事件出現(xiàn)的概率，或者隨機(jī)變量的期望值。通過(guò)隨機(jī)抽樣的方法，以隨機(jī)事件出現(xiàn)的頻率估計(jì)其概率，或者以抽樣的數(shù)字特征估算隨機(jī)變量的數(shù)字特征，并將其作為問(wèn)題的解。這種方法多用于求解復(fù)雜的多維積分問(wèn)題。例如，我們要計(jì)算一個(gè)不規(guī)則圖形的面積，蒙特卡羅方法基于這樣的思想：假想你有一袋豆子，把豆子均勻地朝這個(gè)圖形上撒，然后數(shù)這個(gè)圖形之中有多少顆豆子，這個(gè)豆子的數(shù)目就是圖形的面積。當(dāng)你的豆子越小，撒的越多的時(shí)候，結(jié)果就越精確。借助計(jì)

10、算機(jī)程序可以生成大量均勻分布坐標(biāo)點(diǎn)，然后統(tǒng)計(jì)出圖形內(nèi)的點(diǎn)數(shù)，通過(guò)它們占總點(diǎn)數(shù)的比例和坐標(biāo)點(diǎn)生成范圍的面積就可以求出圖形面積?？梢钥闯觯琈onte Carlo得到概率模型的解是通過(guò)試驗(yàn)得到的，而不是計(jì)算出來(lái)的。也正是由于這個(gè)原因，對(duì)于那些由于計(jì)算過(guò)于復(fù)雜而難以得到解析解或者根本沒(méi)有解析解的問(wèn)題，Monte Carlo方法是一種有效的求出數(shù)值解的方法。我們利用Monte Carlo模擬的方法實(shí)現(xiàn)對(duì)圓周率的估計(jì)?？紤]邊長(zhǎng)為1的正方形，1為半徑的四分之一圓弧，如圖6.1.1所示。 Monte Carlo對(duì)的估計(jì)在邊長(zhǎng)為1的正方形內(nèi)，等概率的產(chǎn)生n個(gè)隨機(jī)點(diǎn)，。這樣和就是（0,1）上均勻分布的隨機(jī)數(shù)。當(dāng)n

11、個(gè)點(diǎn)中有k個(gè)點(diǎn)落在四分之一圓內(nèi)，既有k個(gè)點(diǎn)滿足關(guān)系式：，則當(dāng)時(shí)，有如下關(guān)系：此時(shí)，圓周率的估計(jì)值為。通過(guò)R語(yǔ)言編程如下：monte-function(n) k-0 x-runif(n) # runif( )函數(shù)的作用是產(chǎn)生均勻分布的隨機(jī)數(shù) y-runif(n) for(i in 1:n) if (xi2+yi2=1) k-k+1 pi runif(10,2,3) #在(2,3)范圍內(nèi)產(chǎn)生10個(gè)隨機(jī)數(shù)1 2.597194 2.998407 2.203209 2.897273 2.403639 2.873541 2.508925 2.8786569 2.003151 2.096483（2）產(chǎn)生n個(gè)

12、均值為，標(biāo)準(zhǔn)差為的正態(tài)分布隨機(jī)數(shù)：rnorm(n, a, b)。當(dāng)a, b默認(rèn)時(shí)，為標(biāo)準(zhǔn)正態(tài)分布N(0,1)的隨機(jī)數(shù)。例如： rnorm(30,1,1)1 0.9623346 0.1193753 0.9032579 -0.6649220 1.2346502 0.75963177 -0.2019421 0.7026093 0.2871598 -0.3824490 0.7913398 0.794809213 0.7883549 2.2074363 0.5559556 3.0239753 1.6879051 2.159033719 -0.3332347 1.4590187 1.2827723 -0.

13、9300292 0.1866562 1.163462425 -1.4930723 1.8948205 1.0855420 0.1551674 1.6671424 0.6467984（3）產(chǎn)生n個(gè)參數(shù)為的隨機(jī)數(shù)：rpois(n, )。例如： rpois(30,3)1 3 2 5 4 5 1 2 2 2 0 5 2 4 1 4 4 1 6 4 6 3 4 2 7 4 2 1 9 6 2下面通過(guò)一個(gè)例子說(shuō)明如何產(chǎn)生具有一定分布律的離散型隨機(jī)變量的隨機(jī)數(shù)。 產(chǎn)生具有分布律0的離散型隨機(jī)變量X的隨機(jī)數(shù)。解 設(shè)是(0,1)上均勻分布的隨機(jī)數(shù)，令.則()是具有隨機(jī)變量X分布律的隨機(jī)數(shù)。

14、例如產(chǎn)生20個(gè)隨機(jī)數(shù)，編程如下：n=100; r-runif(n); x-array(0,dim=c(n); for (i in 1:n) if (ri=0.3) xi-0 else if (ri=0.6) xi-1 else xi-2 x1 2 1 2 1 0 1 2 2 2 1 0 0 2 2 0 2 0 2 2 2 2 1 2 2 0 1 2 0 1 30 1 0 0 0 2 2 1 0 1 0 0 1 2 2 1 1 0 1 2 2 0 1 2 2 0 2 2 1 1 59 2 2 0 2 1 2 0 0 1 0 1 2 2 0 0 0 2 0 2 2 2 0 2 2 0 0 2 2

15、0 88 1 2 2 1 1 1 0 2 1 1 0 0 1 一列火車從A站開(kāi)往B站，某人每天趕往B站上火車，他已經(jīng)了解到火車從A站到B站的運(yùn)行時(shí)間是服從均值為30min，標(biāo)準(zhǔn)差為2min的正態(tài)隨機(jī)變量?；疖嚧蠹s下午13:00離開(kāi)A站，此人大約13:30達(dá)到B站?；疖囯x開(kāi)A站的時(shí)刻及概率及此人到達(dá)B站的時(shí)刻及概率如下表所示?；疖囯x開(kāi)時(shí)刻13:0013:0513:1012:55概率0.700.200.100人到站時(shí)刻13:2813:3013:3213:34概率0.1問(wèn)他能趕上火車的概率是多少。利用Monte Carlo方法進(jìn)行分析。以下是求解過(guò)程的R程序：MC-function

16、(n)r1-runif(n); r2-runif(n); t2-rnorm(n,30,2)t1-array(0,dim=c(1,n); t3-t1;for(i in 1:n)if (r1i=0.7)t1i-0else if (r1i=0.9)t1i-5elset1i-10for(i in 1:n)if (r2i=0.3)t3i-28else if (r2i=0.7)t3i-30else if (r2i=0.9)t3i-32elset3i-34kt3i) k-k+1 k/n做一萬(wàn)次模擬，得到：MC(10000)1 0.6341此人能趕上火車的概率大約是0.6341。6.2報(bào)童問(wèn)題與隨機(jī)庫(kù)存模型本

17、節(jié)應(yīng)用概率統(tǒng)計(jì)知識(shí)，首先介紹基本的Newsboy模型，然后將Newsboy模型拓展，討論隨機(jī)庫(kù)存模型。報(bào)童問(wèn)題Newsboy問(wèn)題中，報(bào)童每天清晨從報(bào)社購(gòu)進(jìn)報(bào)紙，通過(guò)一天的零售后，晚上將沒(méi)有賣掉的報(bào)紙以低于購(gòu)進(jìn)價(jià)的價(jià)格退回。設(shè)進(jìn)價(jià)為c，零售價(jià)為s，剩余退回的價(jià)格為a，問(wèn)其如何確定每天購(gòu)進(jìn)的數(shù)量，使其期望獲益最大。這里滿足。從過(guò)上述假設(shè)，報(bào)童每正常賣掉一份報(bào)紙利潤(rùn)為，退回一份賠，由于需求量事先無(wú)法確定，是隨機(jī)的。若通過(guò)以往銷售的經(jīng)驗(yàn)了解到需求量的隨機(jī)規(guī)律，銷售份的概率為，。我們根據(jù)以及報(bào)紙的進(jìn)價(jià)、零售價(jià)和剩余退回價(jià)格來(lái)建立優(yōu)化模型，求解最優(yōu)的訂購(gòu)量。假設(shè)報(bào)童早晨購(gòu)進(jìn)報(bào)紙的量為n，則 或 ，所以每天

18、的收入也是不確定的。這里考慮報(bào)童在不同銷售情況下，建立每天銷售收入的期望函數(shù)，則 （6.2.1）接下來(lái)求當(dāng)n為何值時(shí)，達(dá)到最大？由于r為離散的，這里用差分的方法來(lái)求式（6.2.1）的極值。令令，且，則 （6.2.2）也就是說(shuō)，當(dāng)，a，s和c具體確定時(shí)，n即可確定。例 6.2.1 某服裝店出售某款夏季時(shí)裝。該款衣服成本100元，售價(jià)200元。如整個(gè)夏季不能售出，則必須降價(jià)為70元。設(shè)降價(jià)后一定可以售出，已知售貨量r服從泊松分布為平均出售數(shù)，根據(jù)以往經(jīng)驗(yàn)，平均出售數(shù)為120件。問(wèn)該店的訂貨量應(yīng)該為多少單位？解：由題意知：s=200，a=70，c=100代入式（6.2.2），可得編程如下：poiss

19、on-function(r) #泊松分布 lamda-120 #期望 y-(exp(-lamda)*(lamdar)/(factorial(r); y #尋找n值f-0f1-poisson(1)for (i in 1:1000) fi+1=(10/13) breakf; i可得且，所以更接近于10/13。故，最佳訂購(gòu)量應(yīng)該為126件。隨機(jī)庫(kù)存模型由于市場(chǎng)對(duì)于商品的需求是隨機(jī)變量，事前難以知道需求的準(zhǔn)確數(shù)值。此時(shí)，無(wú)論工廠或商店無(wú)法決定存貯策略，從存貯的角度來(lái)考慮，假設(shè)在一個(gè)階段開(kāi)始的時(shí)刻原有的庫(kù)存為I，如供應(yīng)不足則須承擔(dān)缺貨費(fèi)，如供應(yīng)有余，則多余的部分仍須存貯起來(lái)。有余存在這種不確定性，就需要

20、計(jì)算隨機(jī)變量的期望值，從而定出最佳的存貯量。我們考慮一個(gè)時(shí)間段落。做下列符號(hào)假設(shè)：原有的存貯量為I；存貯貨物的單價(jià)為k；訂購(gòu)一次的訂購(gòu)費(fèi)為C1，如訂貨量為Q時(shí)，所需要的訂貨費(fèi)為；單位貨物的存貯費(fèi)為，缺貨費(fèi)為；需求量為r的概率為。當(dāng)本階段開(kāi)始時(shí)，訂貨量為Q，存儲(chǔ)量達(dá)到I+Q。則本階段所需要的各種費(fèi)用由訂貨費(fèi)、存貯費(fèi)和缺貨費(fèi)構(gòu)成。訂貨費(fèi)：；存貯費(fèi)：當(dāng)需求時(shí)，未能售出的存貯部分必須付存貯費(fèi)；時(shí)，不需要付存貯費(fèi)。因此，所需要存貯費(fèi)的期望值為：。當(dāng)時(shí)，不付存貯費(fèi)及缺貨費(fèi)。缺貨費(fèi)：當(dāng)需求時(shí)，則會(huì)發(fā)生缺貨現(xiàn)象，必須付缺貨費(fèi)缺貨費(fèi)用的期望值為：綜上，在整個(gè)階段所需的訂貨費(fèi)、缺貨費(fèi)及存貯費(fèi)的期望之和為： （6.

21、2.3）為簡(jiǎn)便起見(jiàn)，記，則式（6.2.3）即為 （6.2.4）求S值使C(S)達(dá)到最小。將需求r的隨機(jī)值按大小順序排列為：，其中，（）。S只從中取值。當(dāng)S取值為時(shí)，記為，則 。與newsboy模型中求極值的方法類似，我們求的最小值。 （6.2.5）記，則 (6.2.6) 令，由于，所以，我們有 （6.2.7）式（6.2.7）右端的數(shù)值稱為臨界值，記為。我們選使不等式成立的得最小值為S，則訂貨量為。模型中還有一個(gè)問(wèn)題需要我們解決，那就是原庫(kù)存消耗到什么水平時(shí)，需要訂貨？假設(shè)這一水平是s，當(dāng)時(shí)，可以不訂貨，當(dāng)時(shí)要訂貨，使庫(kù)存達(dá)到S，訂貨量為。要想確定s，首先需要考察不等式 （6.2.8）因s也只能

22、從中取值，使式（6.2.8）成立的（）值中最小者定為s。當(dāng)時(shí)，式（6.2.8）左端缺貨費(fèi)用的期望值雖然在增加，但訂貨費(fèi)及存貯費(fèi)期望值都在減少。在最不利的情況下，如時(shí)，不等式使成立的，因此s值一定存在。例 6.2.2 某汽車零部件生產(chǎn)企業(yè)，對(duì)某型號(hào)鋼材的需求量的概率為：需求量（噸）8090100110120P(r)0.30.1已知每噸鋼材的購(gòu)價(jià)為k=7500元，訂貨費(fèi)為元，存貯費(fèi)元，缺貨費(fèi)元。求該企業(yè)最優(yōu)的存貯策略。解：（1）臨界值另外，且，因此，S=90噸為最優(yōu)訂貨量。（2）利用式（6.2.8）求s：由于S=90，式（6.2.8）右端為當(dāng)s=80時(shí)，式（6.2.8）左端為此時(shí)

23、，式（6.2.8）成立，故。可知，該企業(yè)最優(yōu)的存貯策略為每當(dāng)鋼材的庫(kù)存低于80噸，補(bǔ)充存貯使存貯量達(dá)到90噸，當(dāng)存貯量大于80噸時(shí)，不需要補(bǔ)充。6.3線性回歸模型回歸分析是應(yīng)用數(shù)理統(tǒng)計(jì)學(xué)研究問(wèn)題的一種重要的方法模型，它的目的是研究變量之間的相互關(guān)系，建立變量之間的經(jīng)驗(yàn)公式，以便達(dá)到預(yù)測(cè)和控制的目的。一元線性回歸模型在一元線性回歸分析里，我們要考察的是：隨機(jī)變量y與普通變量x之間的聯(lián)系。對(duì)于x和y，通過(guò)觀測(cè)或?qū)嶒?yàn)，得到若干對(duì)數(shù)據(jù)，。表6.3.1 GDP與人均消費(fèi)支出年份XY年份XY1978675.1359.819901602.3797.11979716.943719911727.2861.419

24、80763.7464.119921949.8966.61981792.4501.919932187.91048.61982851.1533.519942436.11108.71983931.4572.819952663.71213.119841059.2635.619962889.11322.819851185.271619973111.91380.919861269.6746.519983323.11460.619871393.6788.319993529.31564.419881527836.420003789.71690.819891565.9779.7圖6.3.1 中國(guó)人均GDP與人均

25、消費(fèi)支出從散點(diǎn)圖上發(fā)現(xiàn)，觀測(cè)點(diǎn)基本在一條線附近，從而可以認(rèn)為Y與X的關(guān)系式線性的，既因變量Y主要受自變量X的影響，而這些觀測(cè)點(diǎn)與直線的偏離都是由其它一些不確定因素造成的。因此，我們做如下假定： （6.3.1）其中，式（6.3.1）被稱為總體回歸函數(shù)，和是未知參數(shù)，稱為回歸系數(shù)；表示Y隨X的變化而線性變化的部分，是隨機(jī)誤差，稱為隨機(jī)干擾項(xiàng)，反應(yīng)了未列入方程式的其他一切不確定影響因素對(duì)Y影響的總和，通常假定，且隨機(jī)誤差項(xiàng)與自變量X線性無(wú)關(guān)；稱函數(shù)為一元線性回歸函數(shù)，X為自變量，Y為因變量?？紤]到，是(X,Y)的一組觀測(cè)值，則一元線性回歸模型可表示為： （6.3.2）式中，。求回歸參數(shù)的一種思路是要

26、求圖6.3.1中的點(diǎn)與直線上的點(diǎn)偏離越小越好，若和是未知參數(shù)和的估計(jì)值，則被為回歸值或擬合值。和的最小二乘估計(jì)是指使 （6.3.3）成立，經(jīng)計(jì)算可得：， （6.3.4）式中，由此可得回歸方程為。通常取 （6.3.5）為參數(shù)的估計(jì)量，進(jìn)一步可證明，為的無(wú)偏估計(jì)，即。關(guān)于和估計(jì)的標(biāo)準(zhǔn)差分別為， （6.3.6）因此，根據(jù)數(shù)理統(tǒng)計(jì)學(xué)中區(qū)間估計(jì)的原理，。我們有回歸系數(shù)的區(qū)間估計(jì)為：，從回歸參數(shù)的估計(jì)式（6.3.4）可以知道，在回歸系數(shù)的估計(jì)中，不一定要知道Y與X是否有線性關(guān)系，但如果不存在這種關(guān)系，那么回歸方程便毫無(wú)意義。因此需要對(duì)回歸方程進(jìn)行檢驗(yàn)。在統(tǒng)計(jì)意義上是E(Y)隨X線性變化的變化率，若，則E(

27、Y)實(shí)際上并不隨X作線性變化，僅當(dāng)時(shí)，一元線性回歸方程才有意義。因此假設(shè)檢驗(yàn)為：.通常用三種檢驗(yàn)方法：（1）t檢驗(yàn)。當(dāng)成立時(shí)，統(tǒng)計(jì)量對(duì)于給定的顯著性水平，檢驗(yàn)的拒絕域?yàn)椋?）F檢驗(yàn)。當(dāng)成立時(shí)，統(tǒng)計(jì)量.對(duì)于給定的顯著性水平，檢驗(yàn)的拒絕域?yàn)?（3）相關(guān)系數(shù)檢驗(yàn)。記，則稱R為樣本相關(guān)系數(shù)，對(duì)于給定的顯著性水平，查相關(guān)系數(shù)臨界值表可得，則檢驗(yàn)的拒絕域?yàn)楫?dāng)拒絕時(shí)，認(rèn)為線性回歸方程式顯著的。另外，進(jìn)一步介紹回歸方程對(duì)樣本擬合程度的評(píng)價(jià)，所謂擬合程度，是指樣本觀測(cè)值聚集在樣本回歸線周圍的緊密程度。判斷回歸模型擬合程度優(yōu)劣最常用的數(shù)量尺度為樣本決定系數(shù)，它是建立在對(duì)總離差平方和進(jìn)行分解的基礎(chǔ)之上的。樣本決定系

28、數(shù)的公式為：的取值范圍為。由公式可以看出當(dāng)所有的樣本點(diǎn)都位于回歸直線上時(shí)，說(shuō)明總離差可以完全由所估計(jì)的樣本回歸直線來(lái)解釋。根據(jù)圖6.3.1中，人均消費(fèi)支出與人均GDP的數(shù)據(jù)，建立回歸方程模型，進(jìn)行參數(shù)估計(jì)并做相應(yīng)的檢驗(yàn)，編程如下：x-c(675.1,716.9,763.7,792.4,851.1,931.4,1059.2,1185.2,1269.6,1393.6,1527,1565.9,1602.3,1727.2,1949.8,2187.9,2436.1,2663.7,2889.1,3111.9,3323.1,3529.3,3789.7);y-c(359.8,437,464.1,501.9,5

29、33.5,572.8,635.6,716,746.5,788.3,836.4,779.7,797.1,861.4,966.6,1048.6,1108.7,1213.1,1322.8,1380.9,1460.6,1564.4,1690.8);regression-lm(y1+x) #作線性模型summary(regression)運(yùn)行結(jié)果見(jiàn)圖6.3.2。圖6.3.2 回歸結(jié)果第一部分（call）列出了相應(yīng)的回歸模型的公式，其中y1+x表示，第二部分（Residuals）列出的是殘差的最小值點(diǎn)，1/4分位點(diǎn)，中位數(shù)點(diǎn)，3/4分位點(diǎn)和最大值點(diǎn)。在計(jì)算結(jié)果的第三部分（Coefficients）中，Es

30、timate表示回歸方程參數(shù)的估計(jì)，即和；Std. Error表示回歸參數(shù)的標(biāo)準(zhǔn)差，即和；t value為t值，即，表示p值，即，并且有顯著性程度的標(biāo)記。在計(jì)算結(jié)果的第四部分，Residual standard error表示殘差的標(biāo)準(zhǔn)差，即式（6.3.5）中的，自由度為n-2。Multiple R-Squared為相關(guān)系數(shù)的平方，即，F(xiàn)-statistic表示F統(tǒng)計(jì)量，即，其自由度為(1, n-2)，p-value為p值，即概率值。從計(jì)算結(jié)果可以看出，回歸方程通過(guò)了回歸參數(shù)的顯著性檢驗(yàn)和回歸方程的檢驗(yàn)，因此得到的回歸方程為：Y=196+0.3881X.進(jìn)一步，可以通過(guò)R軟件求預(yù)測(cè)值和預(yù)測(cè)區(qū)間

31、。例如求回歸方程中X=2000時(shí)的預(yù)測(cè)區(qū)間。這里即為求X=2000時(shí)的預(yù)測(cè)值，和置信程度為的置信區(qū)間。編程如下：new-data.frame(x=2000)lm.pred-predict( regression, new, interval=prediction, level=0.95)lm.pred fit lwr upr1 972.2582 893.0552 1051.461因此當(dāng)X=2000時(shí)，得到相應(yīng)的預(yù)測(cè)值為972.2582，預(yù)測(cè)區(qū)間為893.0552, 1051.461。多元線性回歸模型一元線性回歸是一個(gè)主要影響因素作為自變量來(lái)解釋因變量的變化。然而，在現(xiàn)實(shí)問(wèn)題研究中，一種現(xiàn)象常常

32、是與多個(gè)因素相聯(lián)系的，此時(shí)就需要用兩個(gè)或兩個(gè)以上的影響因素作為自變量來(lái)解釋因變量的變化，由多個(gè)自變量的最優(yōu)組合共同來(lái)預(yù)測(cè)或估計(jì)因變量，比只用一個(gè)自變量進(jìn)行預(yù)測(cè)或估計(jì)更有效，更符合實(shí)際。在回歸分析中，如果有兩個(gè)或兩個(gè)以上的自變量，就稱為多元回歸或多重回歸。在實(shí)際的應(yīng)用中，多元線性回歸比一元線性回歸用途更廣且實(shí)用意義更大。在建立多元線性回歸模型時(shí)，隨機(jī)變量與一般變量的多元線性回歸模型為：其中是個(gè)未知參數(shù)，稱為回歸常數(shù)項(xiàng)，稱為回歸系數(shù)；稱為被解釋變量（因變量），是個(gè)可以精確測(cè)量并可控制的一般變量，稱為解釋變量（自變量）。為隨機(jī)擾動(dòng)項(xiàng)，代表主觀或客觀原因造成的不可觀測(cè)的隨機(jī)誤差，它是一個(gè)隨機(jī)變量通常假

33、定滿足。則多元線性總體回歸方程為：系數(shù)表示在其它自變量不變的情況下，自變量變動(dòng)一個(gè)單位時(shí)引起的因變量的平均變動(dòng)單位，其它回歸系數(shù)的含義類似。（1）樣本回歸模型建立設(shè)，是隨機(jī)變量與一般變量的n次獨(dú)立觀測(cè)值，則此時(shí)多元線性模型可表示為： （）其中，獨(dú)立同分布。多元線性回歸樣本方程為：,式中為的估計(jì)值。為了方便起見(jiàn)，令，則式（6.3.7）可改寫(xiě)為： （6.3.8）且滿足，?；貧w方程可改寫(xiě)為： （6.3.9）多元線性回歸方程中同歸系數(shù)的估計(jì)采用最小二乘法。若對(duì)于變量記殘差平方和為，根據(jù)微積分中求極小值原理，可知?dú)埐钇椒胶痛嬖谧钚≈?，欲使達(dá)到最小，用對(duì)的偏導(dǎo)數(shù)使其值等于零。加以整理后可以得到個(gè)方程式：通

34、過(guò)求解這一方程組便可求出的估計(jì)值。得 （）取為殘差向量，取 （6.3.11）為的估計(jì)，也稱為的最小二乘估計(jì)?？梢宰C明：.進(jìn)一步可以證明的方差估計(jì)為：.相應(yīng)的的標(biāo)準(zhǔn)差為其中是對(duì)角線上第i個(gè)元素。（2）顯著性檢驗(yàn)與一元線性回歸分析不同，在多元線性回歸分析中，很難用圖形來(lái)判斷E(y)是否隨作線性變化，因而顯著性檢驗(yàn)尤為重要。對(duì)多元線性回歸方程的擬合程度進(jìn)行測(cè)定、檢驗(yàn)回歸方程和回歸系數(shù)的顯著性。 擬合優(yōu)度檢驗(yàn)測(cè)定多元線性回歸的擬合程度，使用多重判定系數(shù)，其定義為：式中SSR為回歸平方和，SSE為殘差平方和，SST為總離差平方和。當(dāng)?shù)闹捣秶鸀椋浇咏?，回歸平面擬合程度越高；反之越接近0，回歸平面擬合程

35、度越低。 回歸方程的顯著性檢驗(yàn)（F檢驗(yàn)）所謂回歸方程的顯著性檢驗(yàn)就是檢驗(yàn)假設(shè)：所有回歸系數(shù)都等于零，即檢驗(yàn)： ；不全為0。多元線性回歸方程的顯著性檢驗(yàn)一般采用F檢驗(yàn)。F統(tǒng)計(jì)量的定義為，平均的回歸平方和與平均的殘差平方和（均方誤差）之比，對(duì)于多元線性回歸方程，在成立的條件下：式中，SSR為回歸平方和，SSE為殘差平方和，為樣本，為自變量個(gè)數(shù)。F統(tǒng)計(jì)量服從的是第一自由度為，第二自由度為的F分布。從F統(tǒng)計(jì)量的定義式可看出，如果F值較大，則說(shuō)明自變量造成的因變量的變動(dòng)遠(yuǎn)遠(yuǎn)大于隨機(jī)因素對(duì)因變量造成的影響。另外，從另一個(gè)角度來(lái)看，F(xiàn)統(tǒng)計(jì)量也可以反映回歸方程的擬合優(yōu)度。將F統(tǒng)計(jì)量的公式與的公式作結(jié)合轉(zhuǎn)換，可

36、得：可見(jiàn)，如果回歸方程的擬合優(yōu)度高，F(xiàn)統(tǒng)計(jì)量就越顯著；F統(tǒng)計(jì)量越顯著，回歸方程的擬合優(yōu)度就越高。利用F統(tǒng)計(jì)量進(jìn)行回歸方程顯著性檢驗(yàn)的步驟總結(jié)如下：Step 1 提出假設(shè)：,不全為0 Step 2 在成立條件下，計(jì)算F統(tǒng)計(jì)量由樣本觀測(cè)值計(jì)算F值。Step 3根據(jù)給定的顯著性水平確定臨界值，或者計(jì)算F值所對(duì)應(yīng)的相伴概率值p。如果（或者），就拒絕原假設(shè)，接受備擇假設(shè)，認(rèn)為所有回歸系數(shù)同時(shí)與零有顯著性差異，自變量與應(yīng)變量之間存在顯著性的線性關(guān)系，自變量的變化確實(shí)能夠反映因變量的線性變化，回歸方程顯著。如果（或者），則接受原假設(shè)，認(rèn)為所有回歸系數(shù)同時(shí)與零無(wú)顯著性差異，自變量與應(yīng)變量之間不存在顯著性的線性

37、關(guān)系，自變量的變化無(wú)法反映因變量的線性變化，回歸方程不顯著。 回歸系數(shù)顯著性檢驗(yàn)（t檢驗(yàn)）回歸方程的顯著性檢驗(yàn)是對(duì)線性回歸方程的一個(gè)整體性檢驗(yàn)。如果我們檢驗(yàn)的結(jié)果是拒絕原假設(shè)，則意味著因變量Y線性地依賴于自變量，這個(gè)回歸自變量的整體。但是，這并不排除Y并不依賴于其中某些自變量。因此，我們還要對(duì)每個(gè)自變量逐一做顯著性檢驗(yàn)，即回歸系數(shù)的顯著性檢驗(yàn)?；貧w系數(shù)的顯著性檢驗(yàn)是檢驗(yàn)各自變量對(duì)因變量的影響是否顯著，從而找出哪些自變量對(duì)的影響是重要的，哪些是不重要的。對(duì)于多元回歸方程，回歸系數(shù)的顯著性檢驗(yàn)，即檢驗(yàn)假設(shè)，在假設(shè)成立的條件下，T統(tǒng)計(jì)量式中為的對(duì)角線上第j個(gè)元素。t檢驗(yàn)步驟如下：Step 1 提出假

38、設(shè)；式中，表示零假設(shè)，表示備擇假設(shè)。如果零假設(shè)成立，則說(shuō)明對(duì)沒(méi)有顯著性的影響，反之，則說(shuō)明對(duì)有顯著性的影響；Step 2 在成立的前提下，計(jì)算回歸系數(shù)的t統(tǒng)計(jì)量；Step 3 給定的顯著性水平，確定臨界值,或者計(jì)算t值所對(duì)應(yīng)的相伴率值p的大小。應(yīng)注意的是，t檢驗(yàn)的臨界值是由顯著性水平和自由度決定的，這里進(jìn)行的檢驗(yàn)是雙側(cè)檢驗(yàn)，所以臨界值為。如果 （或者）,就拒絕原假設(shè)，接受備擇假設(shè)，認(rèn)為回歸系數(shù)與零有顯著性差異，該自變量和應(yīng)變量之間存在顯著的線性關(guān)系，它的變動(dòng)較好地解釋說(shuō)明應(yīng)變量的變動(dòng)，應(yīng)保留在回歸方程中；反之，如果 （或者）,就接受原假設(shè)，認(rèn)為回歸系數(shù)與零無(wú)顯著性差異，該自變量和應(yīng)變量之間不存

39、在顯著的線性關(guān)系，它的變動(dòng)無(wú)法較好地解釋說(shuō)明應(yīng)變量的變動(dòng)，應(yīng)剔除回歸方程。例 近年來(lái)，高校招生規(guī)模急劇擴(kuò)大，在教育產(chǎn)業(yè)化的背景下，選取1985-2003年相關(guān)數(shù)據(jù)，對(duì)影響我國(guó)高校招生人數(shù)的各因素及其影響程度的大小進(jìn)行定量分析。表 6.3.2 高校招生人數(shù)影響因素的多元線性回歸分析原始數(shù)據(jù)年份高校招生數(shù)y國(guó)家財(cái)政教育經(jīng)費(fèi) x1農(nóng)村家庭平均收入x2年份高校招生數(shù)y國(guó)家財(cái)政教育經(jīng)費(fèi)x1農(nóng)村家庭平均收入x2198546871227.9397.61995510531028.41577.7198641310270.4423.81996593981211.911926.1198739017285.9462.

40、61997637491357.732090.1198835645340.7544.91998725081565.5922162198928569397.7601.51999922251815.762210.3199029649433.9686.320001284842085.6792253.4199129679482.2708.620011651972582.382366.4199233439564.978420022026113114.2382475.6199342145644.4921.620032689253453.862622.21994508648841221數(shù)據(jù)EXCEL文件如圖6

41、.3.3，首先將數(shù)據(jù)文件存為csv格式。然后編程如下：X-read.csv(d:program FilesRchengxudata2.csv,header=TRUE)y-X,2;x1-X,3;x2-X,4;lm.sol-lm(yx1+x2)summary(lm.sol)運(yùn)行結(jié)果見(jiàn)圖6.3.4。從運(yùn)行結(jié)果可以看出，回歸系數(shù)和回歸方程的檢驗(yàn)都是顯著的，因此，回歸方程為：.圖 數(shù)據(jù)文件圖 6.3.4 運(yùn)行結(jié)果（3）回歸系數(shù)的區(qū)間估計(jì)由參數(shù)的統(tǒng)計(jì)性質(zhì)可知： （）因此，的區(qū)間估計(jì)為：（4）預(yù)測(cè)當(dāng)多元線性回歸方程經(jīng)過(guò)檢驗(yàn)是顯著的，且其中每個(gè)系數(shù)均顯著不為0，則說(shuō)明回歸的結(jié)果是合理的，在此基礎(chǔ)上可用回歸方程

42、作預(yù)測(cè)。當(dāng)時(shí)，代入回歸方程可得.當(dāng)觀測(cè)值為時(shí)，的置信度為的預(yù)測(cè)區(qū)間為：.例如求回歸方程中，當(dāng) x=(900, 1300)T 時(shí)的預(yù)測(cè)區(qū)間。這里即為求x=(900, 1300)T時(shí)的預(yù)測(cè)值，和置信程度為的置信區(qū)間。編程如下：new-data.frame(x1-900,x2-1300)lm.pred-predict( lm.sol, new, interval=prediction, level=0.95)lm.pred fit lwr upr1 51668.3 28980.95 74355.66因此，當(dāng)x=(900, 1300)T時(shí)，得到相應(yīng)的預(yù)測(cè)值為51668.3，預(yù)測(cè)區(qū)間為28980.95,

43、 74355.66。逐步回歸模型逐步回歸的實(shí)質(zhì)是在建立多元回歸方程的過(guò)程中，首先按偏相關(guān)系數(shù)大小次序，將自變量逐個(gè)引入方程，并且對(duì)引入方程中的每個(gè)自變量偏相關(guān)系數(shù)進(jìn)行統(tǒng)計(jì)檢驗(yàn)，效應(yīng)顯著的自變量留在回歸方程內(nèi)。繼續(xù)遴選下一個(gè)自變量，如果效應(yīng)不顯著，停止引入新自變量。由于新自變量的引入，原已引入方程中的自變量由于變量之間的相互作用其效應(yīng)有可能變的不顯著，經(jīng)統(tǒng)計(jì)檢驗(yàn)后要隨時(shí)從方程中剔除，只保留效應(yīng)顯著的自變量，直到不再引入和剔除自變量為止，從而得到最優(yōu)的回歸方程。一般來(lái)說(shuō)，如果在一個(gè)回歸方程中忽略了對(duì)Y有顯著影響的自變量，那么所建立的方程必與實(shí)際有較大的偏離，但變量選的過(guò)多，使用就不方便，特別當(dāng)方程

44、中含有對(duì)Y影響不大的變量時(shí)，可能由于殘差平方和自由度的減小而使的估計(jì)增大，從而影響使用回歸方程作預(yù)測(cè)的精度。因此，適當(dāng)?shù)剡x擇變量以建立一個(gè)“最優(yōu)”的回歸方程是十分重要的。在多元線性逐步回歸中，“最優(yōu)”的含義是指從可供選擇的所有變量中選出對(duì)Y有顯著影響的變量建立方程，且在方程中不含對(duì)Y無(wú)顯著影響的變量。R軟件提供了較為方便的“逐步回歸”計(jì)算函數(shù)step( )，它是以信息統(tǒng)計(jì)量為準(zhǔn)則，通過(guò)選擇最小的AIC信息統(tǒng)計(jì)量，來(lái)達(dá)到刪除或增加變量的目的。我們通過(guò)一個(gè)例子來(lái)說(shuō)明如何通過(guò)R軟件來(lái)實(shí)現(xiàn)逐步回歸的過(guò)程。 某水泥在凝固時(shí)放出的熱量y（單位：卡/克）與水泥中下列四種化學(xué)成分有關(guān)：3CaO.Al2O3的成

45、分(%);：3CaO.Si2O3的成分(%);：4CaO.Al2O3.Fe2O3的成分(%);：2CaO.SiO2的成分(%)。所測(cè)定的數(shù)據(jù)如表6.3.3所示。用回歸分析建立y與四種化學(xué)成分的線性回歸模型。表 水泥中所含化學(xué)成分表試驗(yàn)序號(hào)172666078.52129155274.331156820104.34113184787.6575263395.961155922109.27371176102.78131224472.59254182293.1102147426115.911140233483.8121166912113.3131068812109.4數(shù)據(jù)來(lái)源：首先做多元線性回歸如下：x-read.csv(d:program FilesRchengxuzhubuhuigui.csv,header=F)y-x,5;x1-x,1;x2-x,2;x3-x,3;x4-x,4;lms|t|) (Intercept) 62