多元統(tǒng)計分析第五章 主成分與經(jīng)驗正交分解.doc

第5章主成分分析與經(jīng)驗正交分解5.1主成分分析數(shù)學模型當存在若干個隨機變量時，尋求它們的少量線性組合（即主成分），用以解釋這些隨機變量，是很必要的。首先我們看一個例子。例5.1 為了調查學生的身材狀況，可以測量他們的身高()、體重()、胸圍()和坐高()?？墒怯眠@4個指標表達學生身材狀況不方便。但若用=3.6356+3.3242+2.4770+2.1650表示學生身體魁梧程度;用=-3.9739+1.3582+3.7323-1.5729表示學生胖瘦程度。則這兩個指標(,)很好概括了4個指標（-）。例5.1中，學生不同，身高()、體重()、胸圍()和坐高()不同；（, , , ）是4維隨機向量；,是他們的2個線性組合，,能很好表示, , , 的特性。類似的問題在許多地方出現(xiàn)：可觀測的隨機變量很多,需要選出所有所有隨機變量的少數(shù)線性組合，使之盡可能刻劃全部隨機變量的特性，選出的線性組合就是諸多隨機變量的主成分，又稱為主分量。尋求隨機向量主成分，并加以解釋，稱為主成分分析，又稱為主分量分析。主成分分析在許多學科中都有應用，細節(jié)可參看張堯廷（1991）、Richard(2003),主成分分析在氣象等科學中稱為PCA方法，見吳洪寶（2005）。主成分分析的數(shù)學模型是：對于隨機向量X，想選一些常數(shù)向量，用盡可能多反映隨機向量X的主要信息。也即盡量大。但是的模可以無限增大，從而使無限變大，這是我們不希望的；于是限定模的大小，而改變各分量的比例，使最大；通常取的模為1最方便。 定義5.1 設隨機向量二階矩存在，若常數(shù)向量，在條件1下使最大，則稱是X的第一主成分或第一主分量。由定義可見，盡可能多地反映原來p個隨機變量變化的信息。但是一個主成分往往不能完全反映隨機向量特色，必須建立其它主成分,它們也應當最能反映隨機向量變化，而且他們應當與第一主成分不相關（不包含的信息）。定義5.2若常數(shù)向量c=在條件l，下，使最大，則稱是 X的第二主成分；若常數(shù)向量c=在條件l，下，使最大，則稱是 X的第三主成分；。當隨機向量方差已知時，定理5.1給出主成分的計算公式。定理5.1 設隨機向量方差存在為。特征值從大到小為，對應的彼此正交單位特征向量為。則X的第j個主成分為 與X的內積，即 (6.1)證明：任取p維單位向量c,必有。于是，而在條件下，當，即時，最大，所以X的第一主成分是與X的內積。由條件，可得，于是，從而；所以在條件1、下，當時，最大，所以X的第2個主成分為與X的內積。對第三，第四主成分同樣可證。由證明過程可見：。它稱為第i個主成分的方差貢獻，表示第i個主成分變化大小，從而反映第i個主成分提供的信息的大小。例5.2 設，且則=3.87939，=0.293128，-0.84403，-0.449099=1.6527，=0.449099，-0.293128，0.84403=0.467911，=0.84403，0.449099，-0.293128所以第一主成分就是=0.293128-0.84403 -0.449099；第二主成分就是=0.449099-0.293128+0.84403；第三主成分就是=0.84403+0.449099-0.293128。它們的方差貢獻分別是；。 定義5.3 稱為主成分的方差貢獻率；稱為前k個主成分的累計方差貢獻率；與X第k個分量的相關系數(shù)稱為因子負荷量。例5.2中，方差貢獻率分別是方差貢獻與6的商，即0.6466，0.2755和0.0780。累計方差貢獻率分別是0.6466、0.6466 +0.2755=0.9220和1。當某個主成分的方差貢獻率很小時，認為它提供的信息很少，可以略去此主成分。通常取q,使前q個主成分的累計方差貢獻率達到70%-80%,然后只考慮前q個主分量，用它們解釋隨機向量X的特性，其余主成分認為是觀測誤差等隨機因素造成的。例5.2中只要前兩個主成分就夠了。在實際問題中，X的每一分量可取不同單位，單位取小時（例如長度單位取毫米，甚至微米）該分量的方差會變大，從而在主成分中變得突出；而單位選取不應影響主成分。為了避免量綱對主成分的影響。常常將隨機變量都標準化，即令，它就是無量綱量，令再求X*的主成分，即標準化后的主成分。將代入，可求隨機向量X的主成分。容易證明定理5.2 設隨機向量X的相關陣為，特征值為,對應的彼此正交單位特征向量為，則標準化后X的第j個主成分是。因此，標準化后的主成分稱為由相關陣決定的主成分。直接由隨機向量的協(xié)方差陣算出的主成分稱為由協(xié)差陣決定的主成分。同樣一組隨機變量，用它們的協(xié)差陣和相關陣求出的主成分是不一樣的。這是因為優(yōu)化的準則（目標函數(shù)）不同：前者要求=最大，而后者要求=最大，其中。例5.3 （協(xié)差陣和相關陣決定的主成分不同）設隨機變量；其協(xié)方差陣是，特征值和特征向量是， 。因而由協(xié)方差陣決定的主成分是： ,。但隨機變量標準化后得到；其中 。X*的協(xié)差陣即X的相關陣是，其特征值和特征向量是，從而由相關陣決定的主成分是：。由于主成分由方差決定，可以略去常數(shù)，因而由相關陣得到的主成分可寫為：，可見由協(xié)方差陣與相關陣決定的主成分不同。5.2 樣本主成分及其計算5.2.1樣本主成分實際問題中隨機向量的協(xié)差陣、相關陣都是未知的，只能得到樣品。這時總用樣本協(xié)差陣與樣本相關陣代替協(xié)差陣、相關陣求主成分。定義5.4 用樣本協(xié)差陣與樣本相關陣的特征向量，計算主成分。所得的主成分稱為樣本主成分。這樣求主成分是有道理的：若總體，的特征值和正交單位特征向量是和；是的極大似然估計，即。的特征值為，相應正交單位特征向量為，則可證定理5.3 若X服從正態(tài)分布，則是的極大似然估計；是的極大似然估計。因此，若X服從正態(tài)分布，應當用第j個樣本主成分作為總體主成分的估計值。從樣本協(xié)差陣或樣本相關陣出發(fā)，做主成分分析，所得樣本主成分通常簡稱為主成分。通常取為樣本協(xié)差陣（的無偏估計），由或R算出的樣本相關陣是相同的，所產生（相關差陣決定）的主成分當然相同。而R與有相同的特征向量，R的特征值是特征值的n/(n-1)倍。因而由R與所產生的（協(xié)方差陣決定的）主成分相同。若X不一定服從正態(tài)分布，這時仍可由樣本協(xié)差陣R或相關陣出發(fā)，計算主成分。同上節(jié)指出的一樣：樣本相關陣和樣本協(xié)差陣決定的主成分是不同的。5.2.2 SAS軟件計算樣本主成分 樣本主成分的計算量很大，通常用軟件計算，以下介紹用SAS軟件計算的基本方法。SAS調用PRINCOMP過程(即主成分過程)作主成分分析。調用PRLNCOMP過程時常用兩個語句：（1） PROC PRINCOMP語句。一般形式是 PROC PRINCOMP；其功能是調用PRINCOMP過程。加選項cov指示電腦用協(xié)差陣計算樣本主成分，不加選項cov則電腦用相關陣計算主成分；加選項out=文件名，指示電腦將每個觀測的主成分得分存入一個數(shù)據(jù)集，即“文件名”所表示的數(shù)據(jù)集，加選項n=k指示電腦只計算k個主成分，不加選項n=k則電腦計算全部p個主成分。例如proc princomp data=wang1 out=wang2 n=3;指示電腦對數(shù)據(jù)集wang1中數(shù)據(jù)做主成分分析，求3個主成分，并將各次觀測的主成分得分存入數(shù)據(jù)集wang2。（2） VAR語句其功能是規(guī)定要分析的變量。例如var x1-x3 u1 v2;表示將變量x1,x2,x3,u1,v作為隨機向量進行主成分分析。 計算主成分固然重要，解釋主成分的意義更重要。下面我們介紹用SAS作主成分分析的實例，并對于算出的主成分加以解釋，希望學者對練習題中的主成分也試作解釋。例5.4北京19511976年冬季的氣溫資料如表5-1，第一列為年度，第二列為該年12月的月平均溫度。第三、四列為次年1、2月的月平均溫度。試做主成分分析。 表 5-1 北京19511976年冬季月平均氣溫yearx1x2x319511.0-2.7-4.31952-5.3-5.9-3.51953-2.0-3.4-0.81954-5.7-4.7-1.11955-0.9-3.8-3.11956-5.7-5.3-5.91957-2.1-5.0-1.619580.6-4.3-0.21959-1.7-5.72.01960-3.6-3.61.31961-3.0-3.1-0.819620.1-3.9-1.11963-2.6-3.0-5.21964-1.4-4.9-1.71965-3.9-5.7-2.51966-4.7-4.8-3.31967-6.0-5.6-4.91968-1.7-6.4-5.11969-3.4-5.6-2.01970-3.1-4.2-2.91971-3.8-4.9-3.91972-2.0-4.1-2.41973-1.7-4.2-2.01974-3.6-3.3-2.01975-2.7-3.70.11976-2.4-7.6-2.2 解 因為所有變量單位相同，可用協(xié)方差陣求主成分。以變量year Dec Jan Feb分別表示年度、12月、1月、2月的溫度。采用下列程序data temperat;/*建立數(shù)據(jù)集temperat*/input year Dec Jan Feb;/*建立變量year、Dec、Jan和Feb*/cards;/*以下為數(shù)據(jù)體*/1951 1.0 -2.7 -4.31952 -5.3 -5.9 -3.51953 -2.0 -3.4 -0.8.1974 -3.6 -3.3 -2.01975 -2.7 -3.7 0.11976 -2.4 -7.6 -2.2;/*空語句，結束數(shù)據(jù)體*/proc princomp cov;/* 用協(xié)差陣做主成分分析*/var Dec Jan Feb;/* 對變量Dec Jan Feb 作主成分分析*/run; 執(zhí)行上述程序，得到得許多表，主要的是：樣本協(xié)差陣的特陣值表（表頭是Eigenvalues）由他們可得方差貢獻，方差貢獻率及累計方差貢獻率；樣本協(xié)差陣的特征向量表（即主成分的系數(shù)表，表頭為Eigenvectors）。這些表及分析如下 Eigenvalues Eigenvalue Difference Proportion Cumulative PRIN1 4.79742 2.06927 0.552919 0.55292 PRIN2 2.72815 1.57720 0.314429 0.86735 PRIN3 1.15095 . 0.132652 1.00000 上表是樣本協(xié)差陣的特征,值表（表頭為Eigenvalues），其中PRIN1、PRIN2、PRIN3表示3個主成分，上表第2列給出樣本協(xié)差陣的特征值，第4列給出方差貢獻，第5列給出方差貢獻累計百分比。由于前兩個特陣值方差貢獻累計百分比等于0.867354，它大于0.7，所以只需取兩個主成分。 Eigenvectors PRIN1 PRIN2 PRIN3 DEC 0.643587 0.709882 -.286116 JAN 0.213039 0.192899 0.957812 FEB 0.735126 -.677390 -.027085上表是特征向量表（表頭為Eigenvectors）上表給出所考察變量樣本協(xié)差陣的特征向量（0.643587，0.213039，0.735126）、(0.709882,0.192899,-0.677390)和（-0.286116，0.957812，-0.027085）。因此第一、二、三主成分分別是=0.643587Dec+0.213039Jan+0.735126Feb,=0.709882Dec+0.192899Jan-0.677390Feb,=-0.286116Dec+0.957812Jan-0.027085Feb由于第一主成分中Dec,Feb系數(shù)是較大正數(shù)，Jan系數(shù)是較小正數(shù)，說明第一主成分主要表示冬季氣溫偏高的程度，由于1月分的系數(shù)變化較小，冬季氣溫偏高主要由12月，2月溫度的偏高形成。第二主成分Dec系數(shù)與Feb系數(shù)反號較大，反映第二主成分主要表示12月與2月溫度距平的反差，即12月溫度距平減去2月溫度距平所得值的反差。 例5.5 美國各州犯罪率情況如表5-2。試以murder（謀殺）,rape（強奸）,robbery（搶劫）,assult（斗毆）,burglary（夜盜）,larceny（偷竊）,auto（汽車犯罪）為7元隨機向量，做主成分分析。 表 5-2 美國各州犯罪率（十萬人中犯罪人數(shù)）murderraperobberyassultburglarylarcenyautoAlbama14.225.296.8278.31135.51881.9280.7Alaska10.851.696.8284.01331.73369.8753.3Arirona9.534.2138.2312.32346.14467.4439.5Arkansas8.834.2138.2312.32346.14467.4439.5Califonia11.549.4287.0358.02139.43499.8663.5Colorado6.342.0170.7292.91935.23903.2477.1Conecticat4.216.8129.5131.81346.02620.7593.2Delaware6.024.9157.0194.21682.63678.4467.0Florida10.239.6187.9449.11859.93840.5351.4Geogia11.731.1140.5256.51351.12170.2297.9Hawaii7.225.5128.064.11911.53920.4489.4Idaho5.519.439.6172.51050.82599.6237.6Illinois9.921.8211.3209.01085.02828.5528.6Indiana7.426.5123.2153.51086.22498.7377.4Iowa2.310.641.289.8812.52685.1219.9Kansas6.622.0100.7180.51270.42739.3244.3Kentaky10.119.181.1123.3872.21662.1245.4Loisana15.530.9142.9335.51165.52469.9337.7Maine2.413.538.7170.01253.12350.7246.9Maryland8.034.8292.1358.91400.03177.7428.5Masschusetts3.120.8169.1231.61532.22311.31140.1Michigan9.338.9261.9274.61522.73159.0545.5Minnesota2.719.585.985.81134.72559.3343.1Mississippi14.319.665.7189.1915.61239.9144.4Missouri9.628.3189.0233.51318.32424.2378.4Montana5.416.739.2156.8804.92773.2309.3Nebraska3.918.164.7112.7760.02316.1249.1Nevada15.849.1323.1355.02453.14212.6559.2Mew Hampashare3.210.723.276.01041.72343.9293.4New Jersey 5.621.0180.4185.11435.82774.5511.5New Maxico8.839.1109.6343.41418.73008.6259.5New York10.729.4472.6319.11728.02782.0745.8North Carolina10.617.061.3318.31154.12037.8192.1North Dakoda100.99.013.343.8446.11843.0144.7Ohio7.827.3190.5181.11216.02696.8400.4Oklahoma8.629.273.8205.01288.22228.1326.8Oregan4.939.9124.1286.91636.43506.1388.9Pennsyvania5.6 19.0130.3128.0877.51624.1333.2Rhode Island3.610.586.5201.01849.52844.1791.4South Carolina 11.933.0105.9485.31613.62342.4245.1South Dakoda2.013.517.9155.7570.51704.4147.5Tennessee10.129.7145.8203.91259.71776.5314.0Texas13.333.8152.4208.21603.12988.7397.6Utah3.520.368.8147.31171.63004.6334.5Vermont1.415.930.8101.21348.22201.0265.2Virginia9.023.392.1165.7986.22521.2226.7Wasinton4.339.6106.2224.81605.63386.9360.3West Viginia6.013.242.290.9597.41341.7163.3Wiskonsin2.812.952.263.7846.92614.2220.7Wyoming5.421.939.7173.9811.62772.2282.0 解：評估美國各州犯罪率時，用7種犯罪率為7維隨機向量，以50個州的統(tǒng)計數(shù)據(jù)為50次觀測。考慮不同犯罪的犯罪率差異很大，用相關陣計算主成分。采用程序data crime; /*建立數(shù)據(jù)集crime*/input state $ 1-15 murder rape robbery assult burglary larceny auto;/*建立變量state murder rape robbery assult burglary larceny auto。state $ 1-15表示前15列存州名。murder rape robbery assult burglary larceny auto 表7種罪的犯罪率*/cards;/*以下為數(shù)據(jù)體*/Albama 14.2 25.2 96.8 278.3 1135.5 1881.9 280.7Alaska 10.8 51.6 96.8 284.0 1331.7 3369.8 753.3Arirona 9.5 34.2 138.2 312.3 2346.1 4467.4 439.5.West Viginia 6.0 13.2 42.2 90.9 597.4 1341.7 163.3Wiskonsin 2.8 12.9 52.2 63.7 846.9 2614.2 220.7Wyoming 5.4 21.9 39.7 173.9 811.6 2772.2 282.0;proc princomp out=crimprin; /*調用PRINCOMP過程，用相關陣做主成分分析*/ var murder rape robbery assult burglary larceny auto; /*對這7個變量做分析*/run;執(zhí)行以上程序，電腦按相關陣做主成分分析；輸出主要數(shù)表有：樣本相關陣的特征值（表頭為Eigenvalues ）表，由他們可得方差貢獻及方差貢獻率；樣本相關陣的特征向量（表頭為Eigenvectors）。表及解釋如下 Eigenvalues Eigenvalue Difference Proportion Cumulative 1 4.09880 2.88618 0.5855 0.5855 2 1.21263 0.45104 0.1732 0.7588 3 0.76159 0.44819 0.1088 0.8676 4 0.31340 0.05201 0.0448 0.9123 5 0.26139 0.01824 0.0373 0.9497 6 0.24315 0.13411 0.0347 0.9844 7 0.10904 0.0156 1.0000以上是特征值表（表頭為 Eigenvalues of the Correlation Matrix），由表第5列可見，前3個特征值所占比例之和為0.86757，只要取3個主成分就夠了。 Eigenvectors Prin1 Prin2 Prin3 Prin4 Prin5 Prin6 Prin7 murder 0.303311 -.634076 0.154298 -.140310 0.629812 -.040772 0.253239 rape 0.432675 -.167388 -.210696 -.007242 -.139849 0.777963 -.339919 robbery 0.391443 0.019456 0.513231 -.607941 -.427219 -.172227 -.034897 assult 0.401331 -.335621 -.090970 0.586083 -.483781 -.308674 0.212600 burglary 0.434023 0.237752 -.222826 0.045108 0.315123 -.465481 -.622444 larceny 0.361074 0.391665 -.533152 -.276436 0.050926 0.005218 0.594083auto 0.296226 0.496972 0.571152 0.434423 0.255583 0.226870 0.183132以上是特征向量表（表頭為Eigenvectors），從第2列起，每列是1個特征向量。第1個特征向量各個分量值大體相同，近似于=0.38；所以第1主成分y1=0.303311murder+0.432675rape+0.391443robbery+0.401331assult +0.4434023burglary+0.361074larceny+0.29296226auto;表示各州犯罪程度的嚴重性。第2個特征向量各分量對應murder,rape, assult,分量值為負的，對應burglary,larceny,auto分量是正的，murder,rape, assult暴力程度重， burglary,larceny,auto暴力程度輕，因此第二主成分 y2=-0.6634076murder-0.167388rape+0.019456robbery-0.335621assult +0.237752burglary+0.391665 arceny+0.496972 auto反映暴力程度的輕重,第二主成分的值越大，暴力成分越輕。第三主成分的特性不明顯，不考慮。許多統(tǒng)計資料簡化成樣本協(xié)差陣，或樣本相關陣；這時仍可用SAS的princomp過程計算，只是在data步輸入數(shù)據(jù)時要用“_type_=”語句說明。例5.6 測量雄龜甲的長、寬、厚，并求其自然對數(shù)，得到變量；所得24只龜數(shù)據(jù)的協(xié)方差陣如下表，試作主成分分析。表5-3 龜甲數(shù)據(jù)的協(xié)方差陣由于觀測資料已被處理為協(xié)方差陣，而協(xié)方差陣是對稱的，只需要輸入下三角陣即可，協(xié)差陣乘以常數(shù)不改變特征向量和累積方差貢獻率，所以0.001不用輸入。我們采用如下程序 data turtle(type=cov);/*建立數(shù)據(jù)集file*/_type_=cov; /*數(shù)據(jù)集為協(xié)方差陣類型*/input name $ x1-x3; /*建立變量name x1 x2 x3 */cards;/*以下是數(shù)據(jù)體*/x1 11.072 . .x2 8.019 6.417 .x3 8.160 6.005 6.773;/*空語句，結束數(shù)據(jù)體*/proc princomp COV;/*用協(xié)方差陣計算3個主成分*/var x1-x3;/*對變量x1 x2 x3求主成分*/run;執(zhí)行后電腦按相關陣做主成分分析；輸出主要數(shù)表有：協(xié)方差陣的特征值（表頭為Eigenvalues），特征向量表（表頭為Eigenvectors）。解釋如下 Eigenvalues Eigenvalue Difference Proportion Cumulative PRIN1 23.3035 22.7048 0.960493 0.96049 PRIN2 0.5987 0.2389 0.024676 0.98517 PRIN3 0.3598 . 0.014831 1.00000 上表是特征值表，由表第2列可見，特征值分別是23.303496、0.5986906、0.3598188；由上表第5列可見，第1特征值占總變差的96%，所以只需1個主成分，就能解釋全部變化。 Eigenvectors PRIN1 PRIN2 PRIN3 X1 0.683103 -.158344 -.712950 X2 0.510212 -.595012 0.621002 X3 0.522546 0.787964 0.325666以上是特征向量表，由表可見，第1主成分的系數(shù)0.683103、0.510212、0.522546相差不多，所以第1主成分表示龜甲的尺寸的自然對數(shù)和，即龜甲體積的自然對數(shù)。5.3 主成分得分實際問題中常需要知道主成分的值，例如例5.5中需要知道哪個州犯罪程度嚴重，哪個州犯罪程度較輕，這就需要計算每個州第一主成分的值；需要知道哪個州暴力犯罪嚴重，哪個州暴力犯罪較輕，這就需要計算每個州第二主成分的值。將各變量值代入主成分的表達式，就能計算主成分。例如例5.4中北京氣溫的第一主成分是prin1=0.643587Dec +0.213039Jan+0.735126Feb,而1951年Dec、Jan、Feb的值分別是1.0、-2.7、-4.3；所以1951年第一主分量的值就是prin1=0.643587*1.0+0.213039*（-2.7）+0.735126*（-4.3）。但是，由于經(jīng)驗正交分解的需要和計算等原因，我們往往計算主成分得分。定義5.4 當用樣本協(xié)方差陣求主成分時，求各觀測值距平（觀測值減去其平均值），再代入主成分的公式，所得稱為（協(xié)方差陣生成的）主成分得分。例如例5.4中第一主成分是0.643587*Dec+0.213039*Jan+0.735126*Feb; Dec,Jan,Feb的樣本均值分別是-2.74，-4.59，-2.27；1951年Dec,Jan,Feb的值分別是1.0,-2.7,-4.3；所以1951年（協(xié)方差陣生成的）的第一主成分得分就是0.643587*（1.0+2.74）+0.213039*（-2.7+4.59）+0.735126*（-4.3+2.27）=1.32。定義5.5 當用樣本相關陣陣求主成分時，將各觀測標準化（觀測值減去其平均值，除以樣本標準差）再代入主成分的公式，所得稱為（相關陣生成的）主成分得分。例如例5.4用相關陣計算時，第一主成分是0.6388*Dec+0.5734*Jan+0.5129*Feb。而1951年標準化的Dec,Jan,Feb的值分別是2.013，1.613，-1.034；于是1951年的（相關陣生成的）第一主成分得分就是6388*2.013+0.5734*1.613 +0.5129*（-1.034）=1.681由主成分得分的值很容易算出主成分的值，但由于主成分得分與主成分的值差一常數(shù)，因而在比較各次觀測主成分的值時，只需比較主成分得分的值即可。SAS-PRINCOMP過程作主成分分析時，能計算主成分得分，在PROC PRINCOMP語句中加選項OUT=文件名，主成分得分的值即存在該文件中。 例5.7 對于例5.4北京19511976年冬季的氣溫資料，求（協(xié)方差陣生成的）各年主成分得分。 解 采用下列程序data temperat;input year Dec Jan Feb;cards;1951 1.0 -2.7 -4.31952 -5.3 -5.9 -3.51953 -2.0 -3.4 -0.8. .1974 -3.6 -3.3 -2.01975 -2.7 -3.7 0.11976 -2.4 -7.6 -2.2;proc princomp cov out=prin;/*各次觀測的主成分值存入數(shù)據(jù)集prin。*/var Dec Jan Feb;/* 對變量Dec Jan Feb 作主成分分析*/proc print data=prin;/* 打印數(shù)據(jù)集prin所存各次觀測的的主成分得分*/run; 執(zhí)行上述程序，與例5.4相比，增加的SAS輸出是下表，其中PRIN1、PRIN2、PRIN3分別表示第1、2、3主成分得分。 表5-4 北京冬季氣溫主成分得分表 OBS YEAR DEC JAN FEB PRIN1 PRIN2 PRIN3 1 1951 1.0 -2.7 -4.3 1.32159 4.39464 0.79664 2 1952 -5.3 -5.9 -3.5 -2.82663 -1.23681 -0.48750 3 1953 -2.0 -3.4 -0.8 1.81464 -0.24090 0.88972 4 1954 -5.7 -4.7 -1.1 -1.06412 -2.91502 0.71132 5 1955 -0.9 -3.8 -3.1 0.74659 2.02081 0.25417 6 1956 -5.7 -5.3 -5.9 -4.72054 0.22071 0.26664 7 1957 -2.1 -5.0 -1.6 0.82132 -0.07862 -0.59250 8 1958 0.6 -4.3 -0.2 3.73731 1.02475 -0.73246 9 1959 -1.7 -5.7 2.0 3.57608 -2.36830 -1.47492 10 1960 -3.6 -3.6 1.3 2.28606 -2.83781 1.09907 11 1961 -3.0 -3.1 -0.8 1.23497 -0.89291 1.46318 12 1962 0.1 -3.9 -1.1 2.83912 1.35662 -0.18190 13 1963 -2.6 -3.0 -5.2 -1.72085 2.39084 1.56369 14 1964 -1.4 -4.9 -1.7 1.21963 0.50533 -0.69429 15 1965 -3.9 -5.7 -2.5 -1.14787 -0.88178 -0.72358 16 1966 -4.7 -4.8 -3.3 -2.05911 -0.73417 0.38901 17 1967 -6.0 -5.6 -4.9 -4.24241 -0.72751 0.03805 18 1968 -1.7 -6.4 -5.1 -1.79244 2.30614 -1.95308 19 1969 -3.4 -5.6 -2.0 -0.43721 -0.84625 -0.78440 20 1970 -3.1 -4.2 -2.9 -0.60750 0.24643 0.49508 21 1971 -3.8 -4.9 -3.9 -1.94226 0.29187 0.05198 22 1972 -2.0 -4.1 -2.4 0.48932 0.70789 0.26259 23 1973 -1.7 -4.2 -2.0 0.95514 0.63061 0.07014 24 1974 -3.6 -3.3 -2.0 -0.07594 -0.54455 1.47579 25 1975 -2.7 -3.7 0.1 1.96183 -1.40534 0.77828 26 1976 -2.4 -7.6 -2.2 -0.36673 -