偏最小二乘回歸方法(PLS).docx

偏最小二乘回歸方法1 偏最小二乘回歸方法(PLS)背景介紹在經(jīng)濟(jì)管理、教育學(xué)、農(nóng)業(yè)、社會(huì)科學(xué)、工程技術(shù)、醫(yī)學(xué)和生物學(xué)中，多元線性回歸分析是一種普遍應(yīng)用的統(tǒng)計(jì)分析與預(yù)測(cè)技術(shù)。多元線性回歸中，一般采用最小二乘方法(Ordinary Least Squares :OLS)估計(jì)回歸系數(shù)，以使殘差平方和達(dá)到最小，但當(dāng)自變量之間存在多重相關(guān)性時(shí)，最小二乘估計(jì)方法往往失效。而這種變量之間多重相關(guān)性問(wèn)題在多元線性回歸分析中危害非常嚴(yán)重，但又普遍存在。為消除這種影響，常采用主成分分析(principal Components Analysis :PCA)的方法，但采用主成分分析提取的主成分，雖然能較好地概括自變量系統(tǒng)中的信息，卻帶進(jìn)了許多無(wú)用的噪聲，從而對(duì)因變量缺乏解釋能力。最小偏二乘回歸方法(Partial Least Squares Regression：PLS)就是應(yīng)這種實(shí)際需要而產(chǎn)生和發(fā)展的一種有廣泛適用性的多元統(tǒng)計(jì)分析方法。它于1983年由S.Wold和C.Albano等人首次提出并成功地應(yīng)用在化學(xué)領(lǐng)域。近十年來(lái)，偏最小二乘回歸方法在理論、方法和應(yīng)用方面都得到了迅速的發(fā)展，己經(jīng)廣泛地應(yīng)用在許多領(lǐng)域，如生物信息學(xué)、機(jī)器學(xué)習(xí)和文本分類等領(lǐng)域。偏最小二乘回歸方法主要的研究焦點(diǎn)是多因變量對(duì)多自變量的回歸建模，它與普通多元回歸方法在思路上的主要區(qū)別是它在回歸建模過(guò)程中采用了信息綜合與篩選技術(shù)。它不再是直接考慮因變量集合與自變量集合的回歸建模，而是在變量系統(tǒng)中提取若干對(duì)系統(tǒng)具有最佳解釋能力的新綜合變量(又稱成分)，然后對(duì)它們進(jìn)行回歸建模。偏最小二乘回歸可以將建模類型的預(yù)測(cè)分析方法與非模型式的數(shù)據(jù)內(nèi)涵分析方法有機(jī)地結(jié)合起來(lái)，可以同時(shí)實(shí)現(xiàn)回歸建模、數(shù)據(jù)結(jié)構(gòu)簡(jiǎn)化(主成分分析)以及兩組變量間的相關(guān)性分析(典型性關(guān)分析)，即集多元線性回歸分析、典型相關(guān)分析和主成分分析的基本功能為一體。下面將簡(jiǎn)單地?cái)⑹銎钚《嘶貧w的基本原理。2 偏最小二乘法的工作目標(biāo)2.1 偏最小二乘法的工作目標(biāo)在一般的多元線性回歸模型中，如果有一組因變量Y=y1,yq和一組自變量X=x1,xp，當(dāng)數(shù)據(jù)總體能夠滿足高斯馬爾科夫假設(shè)條件時(shí)，根據(jù)最小二乘法，有=X（XTX）-1XTY將是Y的一個(gè)很好的估計(jì)量。從這個(gè)公式容易看出，由于（XTX）必須是可逆矩陣，所以當(dāng)X中的變量存在嚴(yán)重的多重相關(guān)性時(shí)，或者在X中的樣本點(diǎn)數(shù)與變量個(gè)數(shù)相比顯然過(guò)少時(shí)，這個(gè)最小二乘估計(jì)都會(huì)失效并將引發(fā)一系列應(yīng)用方面的困難?？紤]到這個(gè)問(wèn)題，偏最小二乘回歸分析提出了采用成分提取的方法。在主成分分析中，對(duì)于單張數(shù)據(jù)表X，為了找到能最好地概括原數(shù)據(jù)的綜合變量，在X中提取了第一主成分F1,使得F1中所包含的原數(shù)據(jù)變異信息可達(dá)到最大，即Var(F1)max在典型相關(guān)分析中，為了從整體上研究?jī)蓚€(gè)數(shù)據(jù)表之間的相關(guān)關(guān)系，分別在X和Y中提取了典型成分F1和G1，它們滿足r(F1,G1) maxF1T F1=1G1T G1=1在能夠達(dá)到相關(guān)度最大的綜合變量F1和G1之間，如果存在明顯的相關(guān)關(guān)系，則可以認(rèn)為，在兩個(gè)數(shù)據(jù)表之間亦存在相關(guān)關(guān)系。提取成分的做法在數(shù)據(jù)分析的方法中十分常見(jiàn)，除主成分、典型成分以外，常見(jiàn)到的還有Fisher判別法中的判別成分。實(shí)際上，如果F是X數(shù)據(jù)表的某種成分，則意味著F是X中變量的某一線性組合F=Xa,而F作為一個(gè)綜合變量，它在X中所綜合提取的信息，將滿足我們特殊的分析需要。2.2 偏最小二乘回歸分析的建模方法設(shè)有q個(gè)因變量y1,yq和p個(gè)自變量x1,xp,為了研究因變量與自變量的統(tǒng)計(jì)關(guān)系，觀測(cè)n個(gè)樣本點(diǎn)，由此構(gòu)成了自變量與因變量的數(shù)據(jù)表X=【x1,xp】n*p和Y=【y1,yq】n*q。偏最小二乘法回歸分別在X與Y中提取出t1和u1(也就是說(shuō)，t1是x1,xp的線性組合，u1是y1,yq的線性組合)。在提取這兩個(gè)成分時(shí)，為了回歸分析的需要，有下列兩個(gè)要求：（1） t1和u1應(yīng)盡可能大地?cái)y帶它們各自數(shù)據(jù)表中的變異信息（2） t1和u1的相關(guān)程度能達(dá)到最大這兩個(gè)要求表明，t1和u1應(yīng)盡可能好地代表數(shù)據(jù)表X和Y，同時(shí)自變量的成分t1對(duì)因變量的成分u1又有最強(qiáng)的解釋能力。在第一個(gè)成分t1和u1被提取后，偏最小二乘法回歸分別實(shí)施X對(duì)t1的回歸以及Y對(duì)t1的回歸。如果方程達(dá)到了滿意的精度，則算法終止；否則，將利用X被t1解釋后的殘余信息以及Y被t1解釋后的殘余信息進(jìn)行第二輪的成分提取。如此遞推，直到能達(dá)到一個(gè)較為滿意的精度為止。若最終對(duì)X共提取了m個(gè)成分t1,tm，偏最小二乘法回歸將通過(guò)實(shí)施YK對(duì)t1,tm的回歸，然后再表達(dá)成YK關(guān)于原變量x1,xp的回歸方程，k=1,q。3 計(jì)算方法推導(dǎo)3.1 普遍采用的計(jì)算推導(dǎo)過(guò)程為了數(shù)學(xué)推導(dǎo)方便起見(jiàn)，首先將數(shù)據(jù)做標(biāo)準(zhǔn)化處理。X經(jīng)標(biāo)準(zhǔn)化處理后的數(shù)據(jù)矩陣記為E0=(E01,E0P)n*p,Y經(jīng)過(guò)標(biāo)準(zhǔn)化處理后的數(shù)據(jù)矩陣記為F0=(F01,F0q)n*q。第一步，記t1是E0的第一個(gè)成分，t1=E0w1, w1是E0的第一個(gè)軸，它是一個(gè)單位向量，即|w1|=1；記u1是F0的第一個(gè)成分，u1=F0c1, c1是F0的第一個(gè)軸，它是一個(gè)單位向量，即|c1|=1。如果要t1,u1能分別很好德代表X與Y中的數(shù)據(jù)變異信息，根據(jù)主成分分析原理，應(yīng)該有Var(t1)maxVar(u1)max另一方面，由于回歸建模的需要，又要求t1對(duì)u1有最大的解釋能力，由典型相關(guān)分析的思路，t1與u1的相關(guān)度應(yīng)達(dá)到最大值，即r(t1,u1)max因此綜合起來(lái)，在偏最小二乘回歸中，我們要求t1與u1協(xié)方差達(dá)到最大，即Cov(t1,u1)=Var(t1)Var(u1)r(t1,u1)max即求解下列優(yōu)化問(wèn)題maxw1T w1=1 （3-1）c1T c1=1因此，將在|w1|=1和|c1|=1的約束條件下，去求（w1TE0TF0c1）的最大值。此種情況下我們就可以用拉格朗日算法求其最優(yōu)解，記s= w1TE0TF0c1-1(w1T w1-1)- 2（c1T c1-1）對(duì)s分別求關(guān)于w1、c1、1、2的偏導(dǎo)，并令之為零，有 E0TF0c1-21 w1=0 （3-2） F0TE0w1-22 c1=0 （3-3） -( w1T w1-1)=0 （3-4） -( c1T c1-1)=0 （3-5）由（3-2）（3-5）可以推出21=22= w1TE0TF0c1=記1=21=22= w1TE0TF0c1,所以1是優(yōu)化問(wèn)題的目標(biāo)函數(shù)值。把式（3-2）和式（3-3）寫成 E0TF0c1= 1 w1 （3-6） F0TE0w1= 1 c1 （3-7） 將式（3-7）代入式（3-6），有 E0TF0F0TE0w1= 12 w1 (3-8)由式（3-8）可知，w1是矩陣E0TF0F0TE0特征向量，對(duì)應(yīng)的特征值為12，1是目標(biāo)函數(shù)值，要求取得其最大值，所以w1是對(duì)應(yīng)于矩陣E0TF0F0TE0最大特征值12的單位特征向量。求得軸w1和c1后，即可得到成分t1=E0w1u1=F0c1然后，分別求E0和F0對(duì)t1和u1的回歸方程其中，向量；E1，F(xiàn)1*，F(xiàn)1為回歸方程的殘差矩陣。 第2成分t2的提取，以E1取代E0 , F1取代F0 , 用上面的方法求第2個(gè)軸W2和第2個(gè)成分t2 ,有 同樣，E1 , F1分別對(duì)t2做回歸, 得到同理可推得第h 成分th , h 的個(gè)數(shù)可以用交叉有效性原則進(jìn)行, h 小于X 的秩。如此計(jì)算下去，如果X的秩為A，則會(huì)有E0=t1P1T+tAPATF0= t1r1T+tArAT+FA由于t1,tA均可以表示成E01,E0P的線性組合，因此，上式可以還原成YK= F0K關(guān)于XJ=E0J的回歸方程形式Y(jié)K=bk1X1+ bkPXP+FAK k=1,.,q3.2一種簡(jiǎn)潔的計(jì)算推導(dǎo)過(guò)程3.1中介紹的推導(dǎo)思路是最為常見(jiàn)的，在3.2中將介紹一種更為簡(jiǎn)潔的計(jì)算方法，即直接在E0,Em-1矩陣中提取成分t1,tm(mp)。要求th能盡可能多地?cái)y帶X中的信息，同時(shí)，th對(duì)因變量系統(tǒng)F0有最大的解釋能力。這時(shí)無(wú)需在F0中提取成分uh，并且在迭代算法中也無(wú)需使用其殘差矩陣，而始終直接用F0進(jìn)行計(jì)算。這可以使計(jì)算過(guò)程大為簡(jiǎn)化，并且對(duì)算法結(jié)論的解釋也更為方便。下面討論成分t1,tm(m=A,A=R(X)的一種新原則。在3.1中推導(dǎo)偏最小二乘法回歸算法時(shí)，第一步的思路是在因變量F0抽取一個(gè)成分u1=F0c1，同時(shí)在自變量E0中抽取一個(gè)成分t1=E0w1,成分的抽取原則是max。在這個(gè)原則下得知w1，c1，u1，t1的計(jì)算方法如下：（1）w1是矩陣E0TF0F0TE0最大特征值的特征向量，成分t1=E0w1；（2）c1是矩陣F0TE0E0TF0最大特征值的特征向量，成分u1=F0c1；在求得成分u1，t1以后，分別實(shí)施E0在t1上的回歸，并生成殘差矩陣E1，以及F0在t1上的回歸，得到殘差矩陣F1。再以E1，F(xiàn)1取代E0，F(xiàn)0進(jìn)行第二輪成分的提取計(jì)算，注意到成分u1,um是不參加回歸計(jì)算的，因此是否可以考慮不提取因變量的成分呢？為此，用下述原則提取比變量中的成分t2是與3.1中介紹的方法，結(jié)果是完全等價(jià)的，即由于F0K是標(biāo)準(zhǔn)化變量，所以Cov(F0K, E0w1)=Var(EOw1) r(F0K, E0w1)因此，該優(yōu)化原則是求成分t1=E0w1，使得t1能攜帶盡可能多的E0變異，同時(shí)，t1對(duì)因變量F0K(k=1,q)的解釋能力會(huì)綜合達(dá)到最大值。由于在目標(biāo)函數(shù)上配上常量（n-1）2不影響其求解，即（n-1）2Cov2(F0K, E0w1)= 2= w1T E0TF0KF0KTE0w1= w1T E0T(F0KF0KT)E0w1= w1T E0TF0F0TE0w1為了求w1采用拉格朗日算法求解，記s= 2-1(w1T w1-1)= w1T E0TF0F0TE0w1-1(w1T w1-1)對(duì)s求關(guān)于w1和1的偏導(dǎo)，并令之為零，得 2 E0TF0F0TE0w1-21 w1=0 （3-9） -( w1T w1-1)=0 （3-10）由式（3-9）可知E0TF0F0TE0w1=1 w1可見(jiàn)，最優(yōu)解w1應(yīng)是矩陣E0TF0F0TE0的一個(gè)特征向量，將它代入目標(biāo)函數(shù)，并且由式（3-10）可得2= w1T E0TF0F0TE0w1= w1T （1 w1）=1因此1矩陣E0TF0F0TE0的最大特征根，w1則是其相應(yīng)的特征向量。由此可見(jiàn)，在新的原則下，w1仍然是對(duì)應(yīng)于E0TF0F0TE0最大特征值的特征向量，而這個(gè)新的原則完全沒(méi)有提取到F0成分u1提取。也就是說(shuō)，t1=E0w1提取可以不依賴對(duì)u1的提取，而這種新的原則又從新的角度說(shuō)明了t1的意義。從這個(gè)新的原則出發(fā)，對(duì)c1,u1的計(jì)算就可以省略。不過(guò)，在偏最小二乘法回歸的一些解釋技術(shù)中，由于u1可以較好地概括F0中的信息，因此，它常常也是很有用。4 應(yīng)用舉例下面將通過(guò)兩個(gè)具體的案例分析, 以進(jìn)一步理解偏最小二乘回歸的工作過(guò)程和它的特點(diǎn)。4.1 應(yīng)用舉例一應(yīng)用舉例一將采用Linnerud給出的關(guān)于體能訓(xùn)練的數(shù)據(jù)進(jìn)行典型相關(guān)分析。在這個(gè)數(shù)據(jù)系統(tǒng)中被觀測(cè)樣本點(diǎn)，是某健身俱樂(lè)部的20位中年男子。被觀測(cè)變量分為兩組，第一組是身體特征指標(biāo)X，包括：體重、腰圍、脈搏；第二組變量是訓(xùn)練結(jié)果指標(biāo)Y，包括：?jiǎn)胃堋澢?、跳高。原始?shù)據(jù)表見(jiàn)表4-1。表4-1 原始數(shù)據(jù)表 在簡(jiǎn)化算法中，對(duì)于h=1,2,3時(shí)，有h=2=（n-1）2Cov2(F0K, th)計(jì)算可得： 1/(n-1)2=1.2724262/(n-1)2=0.0387633/(n-1)2=0.026655而成分th的方差，uk的方差以及th與 uk相關(guān)系數(shù)的平方r2(th, uk)在表4-2中列出。表4-2 Var(th), Var(uk)和r2(th, uk)hVar(th)Var(uk)r2(th, uk)12.02522.05030.306620.43811.90710.046430.23551.15050.0983記第h個(gè)軸是wk，第h個(gè)成分tk為tk=Eh-1wh （h=1,2,3）其中tk亦可以表示成原自變量E0的線性組合，即tk=E0wh*則wh*=（1-wjpjT）wh。表4-3給出wh*與wh的取值。表4-3 wh*與wh的取值自變量w1w2w3w1*w2*w3*X1-0.589890.46879-0.65747-0.589890.36793-0.93459X2-0.77134-0.568010.28706-0.77134-0.699890.80231X30.23888-0.67647-0.696660.23888-0.63562-0.22282在利用Eh-1對(duì)th進(jìn)行回歸時(shí)，有回歸系數(shù)向量pk，h=1,2,3,見(jiàn)表4-4。表4-4 回歸系數(shù)pkp1P2P3-0.6659-0.0197-0.6574-0.6760-0.35460.28700.3589-1.1942-0.6966成分tk=Eh-1wh的取值見(jiàn)表4-5。表4-5 tk取值表NOt1t2t31-0.6430.591-0.1312-0.7700.1670.1343-0.907-0.5210.04840.688-0.6800.3465-0.4871.133-0.1826-0.229-0.0720.0257-1.404-0.077-0.57280.744-0.211-0.03291.715-0.655-1.557101.1630.1670.333110.3650.7010.201120.7430.6980.002131.187-0.7570.36614-4.390-0.7600.25515-0.8230.974-0.08316-0.749-0.521-0.66717-0.393-0.2030.564181.1990.7830.092191.0480.3730.319201.942-1.1290.568通過(guò)交叉驗(yàn)證的方法可得，之取一個(gè)成分t1時(shí)，擬合方程的預(yù)測(cè)性為最佳，不過(guò)為了后面作圖和解釋的方便起見(jiàn)，我們?nèi)蓚€(gè)成分t1，t2擬合預(yù)測(cè)模型。yk=r1kt1+ r2kt2 k=1,2,3由于成分th可以寫成自變量xj的函數(shù)，即有th=wh1*x1+ wh2*x2+ wh3*x3由此可得兩個(gè)成分t1，t2所建立的偏最小二乘回歸模型為yk=r1k（w11*x1+ w12*x2+ w13*x3）+ r2k（w21*x1+ w22*x2+ w23*x3） =（r1k w11*+ r2k w21*）x1+（r1k w12*+ r2k w22*）x2+（r1k w13*+ r2k w23*）x3回歸系數(shù)的計(jì)算結(jié)果見(jiàn)表4-6。表4-6 回歸系數(shù)rkKr1r2r310.34160.33630.477220.41600.29070.455430.14290.0651-0.2125所以，有F01=-0.077E01-0.499 E02-0.132 E03F02=-0.138E01-0.524 E02-0.085 E03F01=-0.060E01-0.156 E02-0.007 E03將標(biāo)準(zhǔn)化變量Fok(k=1,2,3)和Eoj(j=1,2,3)分別還原成原始變量，yk(k=1,2,3)以及xj(j=1,2,3),則回歸方程為：Y1*=47.02-0.0166x1-0.824x2-0.097x3Y2*=612.57-0.351x1-10.52x2-0.741x3Y3*=183.98-0.125x1-2.497x2-0.052x3為了快速直觀地觀察出各個(gè)自變量在解釋Yk時(shí)的作用，可以繪制回歸系數(shù)圖，見(jiàn)圖4-1圖4-1 回歸系數(shù)的直方圖從回歸系數(shù)圖中可以立刻觀察到，腰圍變量在解釋三個(gè)回歸方程時(shí)起到了極為重要的作用，然而，與單杠及彎曲相比，跳高成績(jī)的回歸方程顯然不夠理想，三個(gè)自變量對(duì)它的解釋能力均很低。因此有必要考慮對(duì)自變量做適當(dāng)?shù)恼{(diào)整。為了考察這三個(gè)回歸方程的模型精度，我們以（yik*,yik）為坐標(biāo)值，對(duì)所有的樣本點(diǎn)繪制預(yù)測(cè)圖。yik*是第k個(gè)變量，第i個(gè)樣本點(diǎn)（yik）的預(yù)測(cè)值。在預(yù)測(cè)直方圖上，如果所有樣本點(diǎn)都能在圖的對(duì)角線附近