偏最小二乘回歸方法

偏最小二乘回歸方法1偏最小二乘回歸方法仲任）背景介紹在經(jīng)濟管理、教育學、農(nóng)業(yè)、社會科學、工程技術、醫(yī)學和生物學中，多元線性回歸分析是一種普遍應用的統(tǒng)計分析與預測技術。多元線性回歸中 ，一般采用最小二乘方法（OrdinaryLeastSquares:OLS）估計回歸系數(shù)，以使殘差平方和達到最小，但當自變量之間存在多重相關性時，最小二乘估計方法往往失效。而這種變量之間多重相關性問題在多元線性回歸分析中危害非常嚴重，但又普遍存在.為消除這種影響，常采用主成分分析（principalComponentsAnalysis:PCA）的方法，但采用主成分分析提取的主成分，雖然能較好地概括自變量系統(tǒng)中的信息，卻帶進了許多無用的噪聲，從而對因變量缺乏解釋能力。最小偏二乘回歸方法（PartialLeastSquaresRegression：PLS）就是應這種實際需要而產(chǎn)生和發(fā)展的一種有廣泛適用性的多元統(tǒng)計分析方法。它于1983年由S。Wold和C.Albano等人首次提出并成功地應用在化學領域。近十年來，偏最小二乘回歸方法在理論、方法和應用方面都得到了迅速的發(fā)展,己經(jīng)廣泛地應用在許多領域，如生物信息學、機器學習和文本分類等領域.偏最小二乘回歸方法主要的研究焦點是多因變量對多自變量的回歸建模，它與普通多元回歸方法在思路上的主要區(qū)別是它在回歸建模過程中采用了信息綜合與篩選技術。它不再是直接考慮因變量集合與自變量集合的回歸建模，而是在變量系統(tǒng)中提取若干對系統(tǒng)具有最佳解釋能力的新綜合變量（又稱成分），然后對它們進行回歸建模。偏最小二乘回歸可以將建模類型的預測分析方法與非模型式的數(shù)據(jù)內(nèi)涵分析方法有機地結合起來，可以同時實現(xiàn)回歸建模、數(shù)據(jù)結構簡化（主成分分析）以及兩組變量間的相關性分析（典型性關分析），即集多元線性回歸分析、典型相關分析和主成分分析的基本功能為一體。下面將簡單地敘述偏最小二乘回歸的基本原理。2偏最小二乘法的工作目標2.1偏最小二乘法的工作目標在一般的多元線性回歸模型中，如果有一組因變量丫二藥廠?。鸵唤M自變量X=｛xj，…，xp｝，當數(shù)據(jù)總體能夠滿足高斯一馬爾科夫假設條件時，根據(jù)最小二乘法，有Y=X（XtX）-1XtY『將是Y的一個很好的估計量。從這個公式容易看出，由于（XtX）必須是可逆矩陣,所以當X中的變量存在嚴重的多重相關性時，或者在X中的樣本點數(shù)與變量個數(shù)相比顯然過少時，這個最小二乘估計都會失效并將引發(fā)一系列應用方面的困難。考慮到這個問題，偏最小二乘回歸分析提出了采用成分提取的方法。在主成分分析中，對于單張數(shù)據(jù)表X,為了找到能最好地概括原數(shù)據(jù)的綜合變量，在X中提取了第一主成分F1，使得F1中所包含的原數(shù)據(jù)變異信息可達到最大，即Var(F1)fmax在典型相關分析中，為了從整體上研究兩個數(shù)據(jù)表之間的相關關系分別在X和Y中提取了典型成分F1和G1，它們滿足r(F1,G1)fmaxFiTFi=1GiTGi=1在能夠達到相關度最大的綜合變量F1和G1之間，如果存在明顯的相關關系，則可以認為，在兩個數(shù)據(jù)表之間亦存在相關關系。提取成分的做法在數(shù)據(jù)分析的方法中十分常見，除主成分、典型成分以外，常見到的還有Fisher判別法中的判別成分。實際上，如果F是X數(shù)據(jù)表的某種成分，則意味著F是X中變量的某一線性組合F=Xa，而F作為一個綜合變量，它在X中所綜合提取的信息，將滿足我們特殊的分析需要。2.2偏最小二乘回歸分析的建模方法設有q個因變量｛"，...，yq｝和p個自變量｛勺…，xp｝，為了研究因變量與自變量的統(tǒng)計關系，觀測n個樣本點，由此構成了自變量與因變量的數(shù)據(jù)表X4x1，?.,xp]n^p和YMy】".，yq】n*q。偏最小二乘法回歸分別在X與Y中提取出t1和u1(也就是說,t1是％，...,xp的線性組合，u1是y】,...，yq的線性組合)。在提取這兩個成分時，為了回歸分析的需要，有下列兩個要求：t1和u1應盡可能大地攜帶它們各自數(shù)據(jù)表中的變異信息t1和u1的相關程度能達到最大這兩個要求表明,^和u1應盡可能好地代表數(shù)據(jù)表X和Y,同時自變量的成分t1對因變量的成分u1又有最強的解釋能力。在第一個成分t1和u1被提取后，偏最小二乘法回歸分別實施X對t1的回歸以及Y對t1的回歸。如果方程達到了滿意的精度，則算法終止；否則將利用X被t1解釋后的殘余信息以及Y被t1解釋后的殘余信息進行第二輪的成分提取。如此遞推，直到能達到一個較為滿意的精度為止。若最終對X共提取了m個成分t1，.，tm，偏最小二乘法回歸將通過實施Yk對t1".,tm的回歸，然后再表達成Yk關于原變量X1,…，xp的回歸方程冰=1,…，q。3計算方法推導3.1普遍采用的計算推導過程為了數(shù)學推導方便起見，首先將數(shù)據(jù)做標準化處理。X經(jīng)標準化處理后的數(shù)據(jù)矩陣記為Eo=(E°i，…，E°p)或丫經(jīng)過標準化處理后的數(shù)據(jù)矩陣記為F=(F,…，F(xiàn))。0 010Pnzp 001 0qn*q第一步，記t1是E0的第一個成分，t1=E0w1,w1是E0的第一個軸，它是一個單位向量，即I|w1I|=1；記U1是F0的第一個成分,u1=F0c1，C1是F0的第一個軸，它是一個單位向量，即I|CJ1=1O如果要t1，%能分別很好德代表X與Y中的數(shù)據(jù)變異信息，根據(jù)主成分分析原理應該有Var(tJ—maxVar(u】)—max另一方面,由于回歸建模的需要，又要求t1對％有最大的解釋能力,由典型相關分析的思路，t1與u1的相關度應達到最大值，即r(t1，u】)一max因此綜合起來，在偏最小二乘回歸中，我們要求t1與％協(xié)方差達到最大，即Cov(t1,u1)=^Var(t1)Var(u1)r(t1,ul)^max即求解下列優(yōu)化問題max<E0w1，F(xiàn)0C1>TOC\o"1-5"\h\zw1tw1=1 (3-1)C1Tc1=1因此，將在|IwjI=1和||C1I|=1的約束條件下，去求(w/EoTFoCj)的最大值.此種情況下我們就可以用拉格朗日算法求其最優(yōu)解，記s=wtEtFC—入(wTw-1)-入(cTC—1)sVV100C1 1、VV1TVV1 2 1’C1對s分別求關于w「&、％、入2的偏導，并令之為零，有ds礦=E°tF0C1-2入1w1=0 (3—2)1dsE=F0TE0w1-2入2C1=0 (3—3)1dsd^=—(wjw1—1)=0 (3-4)1

ds（ds（3—5）由(3—2)~(3-5)可以推出TOC\o"1-5"\h\z2入，=2入-=w.tE-tFcC.=〈EcW_,FC〉

1 2 1001 01 01記01=2入1=2入2=w1tE0tF0c1，所以91是優(yōu)化問題的目標函數(shù)值.把式（3—2）和式（3—3）寫成E0tF°C1=?1W1 （3—6）F0TE0W1=?1C1 （3—7）將式（3-7）代入式（3-6）,有E0TF0F0TE0W1=?12W1 （3-8）由式(3—8)可知，w1是矩陣E0tF0F0tE0特征向量，對應的特征值為?12,?1是目標函數(shù)值，要求取得其最大值，所以w1是對應于矩陣E0tF0F0tE0最大特征值?12的單位特征向量。求得軸w1和c1后，即可得到成分t1=E0W1U1=F0C1然后,分別求E0和F0對t1和u1的回歸方程E=tPt+E,F=uQT+F*,F

0 11 1 0 11 1 0其中,P其中,P]=ETt1/t12，Q1=F0TuW2,向量r1=F0Tt1/t12;E1,F1*，F(xiàn)1為回歸方程的殘差矩陣.第2成分t2的提取，以E］取代E°,、取代F°,用上面的方法求第2個軸W2和第2個成分t2,有W=.1^,t=EW2時財211同樣，馬，F(xiàn)1分別對t2做回歸，得到E=tPt+E,F=trT+F1 22 2 1 2 2 2同理可推得第h成分th,h的個數(shù)可以用交叉有效性原則進行,h小于X的秩.如此計算下去，如果X的秩為A,則會有E0=t1P1T+-+tAPATF0=t1r1T+-+tArAT+FA由于",..&均可以表示成E01，…，E0P由于",..&均可以表示成E01，…，Xj=E0J的回歸方程形式YK=bkiXi+???+bkpXp+FAK k=1,。。,q3.2一種簡潔的計算推導過程3.1中介紹的推導思路是最為常見的，在3。2中將介紹一種更為簡潔的計算方法，即直接在E。，…，Em—矩陣中提取成分t1，?，tm（m〈p）。要求th能盡可能多地攜帶X中的信息，同時，th對因變量系統(tǒng)F0有最大的解釋能力。這時無需在F0中提取成分uh，并且在迭代算法中也無需使用其殘差矩陣，而始終直接用F0進行計算。這可以使計算過程大為簡化，并且對算法結論的解釋也更為方便。下面討論成分"，...,tm（m〈=A，A=R（X））的一種新原則。在3。1中推導偏最小二乘法回歸算法時，第一步的思路是在因變量F0抽取一個成分u1=F0c1，同時在自變量E0中抽取一個成分tEw，成分的抽取原則是max〈EwFC〉。1=01 01'01在這個原則下得知w1，c1，u1，t1的計算方法如下：（1） w1是矩陣E0tF0F0tE0最大特征值的特征向量，成分t1=E0w1；（2） c1是矩陣F0tE0E0tF0最大特征值的特征向量，成分u1=F0c1；在求得成分u1，t1以后，分別實施E0在t1上的回歸，并生成殘差矩陣］以及F0在t1上的回歸，得到殘差矩陣、。再以E1，F(xiàn)1取代E0,F°進行第二輪成分的提取計算，注意到成分％，...，um是不參加回歸計算的，因此是否可以考慮不提取因變量的成分呢？為此，用下述原則提取比變量中的成分t2是與3.1中介紹的方法，結果是完全等價的，即踏*網(wǎng)）由于f0K是標準化變量，所以Cov（F0K，E0w1）=^yar（EOwl）r（F0K，E0w1）因此，該優(yōu)化原則是求成分t1=E0w1，使得t1能攜帶盡可能多的E0變異，同時&對因變量F0K（k=1，…，q）的解釋能力會綜合達到最大值。由于在目標函數(shù)上配上常量（n-1）2不影響其求解，即TOC\o"1-5"\h\z，、乙， 、工\（n」）2 Cov2（F0K，E°w1）= <%E°w1〉2k=1 k=1='w.TE?tFF?,tE"_=w.tE/（'FF?,t）E?w.=w.TE?tFFtFw.1 00K0K01 1 0 0K0K01 1 00001k=1 k=1為了求w1采用拉格朗日算法求解，記s='<Fc,,， E?w>2—X （W“Tw-1）=wTE?tFFtEw-X （W“Tw—1）0K01 1 11’ 1 00001 1 1 1 /k=1對s求關于w1和XJ的偏導，并令之為零，得ds云=2E0TF0F0TE0w1-2X1w1=0 (3-9)1ds~d^=—(w1Tw1-1)=0 (3—10)1由式(3-9)可知E0TF0F0TE0w1=X1w1可見，最優(yōu)解w1應是矩陣E0tF0F0tE0的一個特征向量，將它代入目標函數(shù)，并且由式(3-10)可得工 、、、<F?,， Ertw>2=wtE^tFFtEw=wt（Xw）=X“0K， 01 1 00001 1 1 1 1k=1因此XJ矩陣E0tF0F0tE0的最大特征根，w】則是其相應的特征向量.由此可見，在新的原則下，w1仍然是對應于E0tF0F0tE0最大特征值的特征向量，而這個新的原則完全沒有提取到F0成分u1提取。也就是說，t1=E0w1提取可以不依賴對u1的提取，而這種新的原則又從新的角度說明了J的意義。從這個新的原則出發(fā)，對&，％的計算就可以省略.不過，在偏最小二乘法回歸的一些解釋技術中，由于％可以較好地概括F0中的信息，因此，它常常也是很有用.4應用舉例下面將通過兩個具體的案例分析，以進一步理解偏最小二乘回歸的工作過程和它的特點。4。1應用舉例一應用舉例一將采用Linnerud給出的關于體能訓練的數(shù)據(jù)進行典型相關分析.在這個數(shù)據(jù)系統(tǒng)中被觀測樣本點，是某健身俱樂部的20位中年男子。被觀測變量分為兩組,第一組是身體特征指標X,包括：體重、腰圍、脈搏；第二組變量是訓練結果指標Y，包括:單杠、彎曲、跳高。原始數(shù)據(jù)表見表4-1。

1315434&4142151051424746501505015193364567031162033762122101201717637544602518157325211Z30ao156335415225732?13S33211D均帽178*SO3S.4OS&IO9*小I4&.5S7U.3024.693^3507-215-296Z.575=2S表4—1原始數(shù)據(jù)表Ha*■體重■圖*3uwJfi單狂>2尊曲KA119】3G50516260218937522LIOSO3I9i3A12101Id4162356212105375】時3S461315558618236郁410142721L3881D1388167345U612540917?3174】520040101&43356172S1250113450J712C3812 |166335213210115在簡化算法中，對于h=1,2,3時，有N \E,入h= <F0k，Eh—wi〉2=(n-1)2 Cov2(F0k，th)k=1k=1計算可得：入1/(n—1)2=1.272426入2/(n-1)2=0.038763入3/(n-1)2=0.026655而成分th的方差，uk的方差以及th與uk相關系數(shù)的平方r2(th，uk)在表4-2中列出.表4-2Var(th),Var(uk)和r2(th，uk)hVar(th)Var(uk)r2(th，uk)

12。02522。05030。306620。43811.90710.046430.23551。15050.0983記第h個軸是wk，第h個成分tk為tk=Eh—iwh(h=1，2,3)其中tk亦可以表示成原自變量E0的線性組合，即L=E0Wh*則wh*=h (1-w.p.T)wh。表4—3給出wh*與wh的取值。j=1表4-3wh*與wh的取值自變量w1w2w3w1*w?*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對th進行回歸時，有回歸系數(shù)向量pk，h=1,2,3,見表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的取值見表4-5。表4—5tk取值表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.025

7-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通過交叉驗證的方法可得，之取一個成分t1時，擬合方程的預測性為最佳，不過為了后面作圖和解釋的方便起見，我們?nèi)蓚€成分t1，t2擬合預測模型。y=L.t+n.Lk=1,2，3

7k1k12k2 '由于成分th可以寫成自變量*的函數(shù)，即有t=w大x+w大x+w*x

hwh1x1wh2x2wh3x3由此可得兩個成分t],t2所建立的偏最小二乘回歸模型為vr（w*x+w*x+w*x）+r（w*x+w*x+w*x）TOC\o"1-5"\h\zyk='1k、w11Jvv12x2vv13 3 ’2ktw21x1w22x2w23 3=（r,w一大+ r. w一大） x+（r,w一大+ r.w*）k.+（L. w_*+ r.w一大） x?表4-6回歸系數(shù)表4-6回歸系數(shù)r.Kr1r2r310。34160.33630.477220。41600。29070。455430.14290.0651—0.2125回歸系數(shù)的計算結果見表4-6。k所以，有FO1=—0。077E01-0.499E^0。132E03F02=-0.138E01-0.524e°20.085E03F0i=—°。060E01-0。156Es0。007E03將標準化變量Fok(k=1,2,3)和Eoj(j=1,2,3)分別還原成原始變量，yk(k=1,2,3)以及*(j=1,2,3)，則回歸方程為：Y「=47.02—0。0166x1—0O824x2—0.097x3Y2*=612。57-0.351x1-10.52x2—0.741x3Y3*=183。98-0.125x]-2。497x2-0。052x3為了快速直觀地觀察出各個自變量在解釋Yk時的作用，可以繪制回歸系數(shù)圖，見圖4-1圖4-1回歸系數(shù)的直方圖單杠 ■曲 眺高從回歸系數(shù)圖中可以立刻觀察到，腰圍變量在解釋三個回歸方程時起到了極為重要的作用，然而，與單杠及彎曲相比，跳高成績的回歸方程顯然不夠理想，三個自變量對它的解釋能力均很低.因此有必要考慮對自變量做適當?shù)恼{(diào)整。為了考察這三個回歸方程的模型精度，我們以(yik*,yik)為坐標值，對所有的樣本點繪制預測圖。yik*是第k個變量，第i個樣本點(yik)的預測值.在預測直方圖上，如果所有樣本點都能在圖的對角線附近均勻分布，則方程的擬合值與原值差異很小，這個方程的擬合效果就是滿意的。體能訓練的預測圖如4-2所示。

3)

(。單杠成績回歸方程頊街圈「(frXffiJSBETS方程預H圖j(QH離成皺回歸方覆BJSJ圖*+4。2應用舉例二這是Cornell在1990年采用的一個化工方面的例子.此后，偏最小二乘的提出者S。Wold等人多次引用，成為單因變量偏最小二乘回歸的一個經(jīng)典案例。該例中，有個自變量x1~x7,因變量記為y，如表4—7所示：表4-7自變量和應變量對照表x1-一直接蒸餾成分X2—-重整汽池；xX3—-原油熱裂化油X4--原油催化裂化油；X5——聚合物X6--烷基化物x7-天然香精y-一原辛烷值表4-8給出了12種混合物關于這8個變量的觀測數(shù)據(jù)。要求建立y對x1~x7，的回歸方程，以確定7種構成元素x1~x7對y的影響。表4-812種混合物關于8個變量的觀測數(shù)據(jù)表表化工案例的原始數(shù)據(jù)NoJiX2明了4Jt5XTy10*000.230.000,000.000,740.039&720.0G0.100.000.00。?120.740.0497.830.000.000.000.100.120.740.0496.640.000.490.000.000.120.370.0292,05。，000.000.000.620?120.180.0886.6|-60.000.620*000.000?000.370.0191.270.170.270.100.380,000.000.0881.980.170.190.100.380.020.06。,0883.190.170.210.100.380.000.060.0882.4100.170.150.100.380.020.100.0883.2110.210.360.120.250?000.000.0681.4120.000.000,000.550.000.370.088&1

這8個變量的相關系數(shù)矩陣見表4-9。從相關系數(shù)矩陣中可以看出，在自變量之間存在嚴重的多重相關性，例如r（x1，x3）=0.999，r（x平x7）=0.92，r（x1,x6）=—0.80.實際上，這7個自變量之間有如下關系:x1+x2+...+x7=1表4—98個變量的相關系數(shù)矩陣H?…)Xix3心玦X?yJi0.100.9990.37一(h55—0.800.60—0.84JT30_w—o.54—0,29—0,19—CK5-9—a.07。?37—0-55一6800.61—(k84工4-0-21-0.64-ot?1心0.46—0,27—0.66G.18x7—0.74第一步C—O.8373鹽]——O,070 —O.838Ox3——CL7067jt4-f-O.4938jrjj十O,9851 —O*741Im門/L916=—O.437Ox.一O.037O*》一O.4373j3一O.368 +2577xn-4-O.514Ix^—O+386Rx7做，在C上的計笄:結果為去i士O.96104F=133回歸方程是y =O-482Ot.=一O.2106Xj—O.017 —O,2108jc5—O.177Sx4十O,1242jt5+O.2^78jrt一O.186Ax7第二步，還可以抽取成分r”從而對上面的回歸模型做進一步改善r記Xu『…是"…土在>上回歸后的殘差向量，有[Q.03] If — 0.131 5x[3+0.032 6xu — 0.024 5X|0.070 lxl5 + 0*099 9xtfi+CL049 3x!72/0. 193 =CO,031 7X —(—0.23672)—0.131 5X (