作者：落園

前言的廢話：有種被統(tǒng)計(jì)學(xué)洗腦的感覺,，跟搞統(tǒng)計(jì)的接觸的有點(diǎn)太多了吧…哎，算了,，難得有點(diǎn)感悟記錄一下吧…喵,。本文試圖以一個(gè)典型的統(tǒng)計(jì)學(xué)的思維過程來敘述，呃,，可能功力不足,，湊合著看吧，也不枉我連續(xù)幾天做夢都是直方圖 密度曲線 average treatment effect之類的,。

——————————————-廢話結(jié)束了———————————–

Variance這個(gè)詞很好玩,，如果是用在統(tǒng)計(jì)的情境中就中文譯作方差，公式就是∑(x?xˉ)2/N,，寫成列向量（N*1矩陣）的形式就是Var(X)=E(X2)?E(X)2=X′X/N?(X′1/N)2,，1這里是N×1個(gè)1的列向量,。從公式的形式就可以看出是一個(gè)二階中心距，即距離中心的距離的平方（二階）和,。然后既然是二階距,，那么自然的衡量的就是某個(gè)數(shù)據(jù)集分布的離散情況。越大越散,，越小越密,。此為一維（這里指只有一個(gè)觀測維度的情形）相信這樣的定義實(shí)在是太耳熟能詳了。方差有個(gè)小弟叫做標(biāo)準(zhǔn)差,，就是方差開平方,。這個(gè)也沒啥說的，有意思的是大家習(xí)慣用字母σ來標(biāo)注他,，于是有了著名的六-西格瑪原理…好吧,，其實(shí)最有用的就是正態(tài)分布的幾個(gè)西格瑪了（隨便偷張圖）：

然后我們看簡單的二維。二維就是散點(diǎn)圖,，沒什么特別說的,。找張著名的散點(diǎn)圖來（隨便找的，懶得自己去R里面畫了,。最后還是乖乖的去R里面的畫了,，還是自己畫的好用一些，果然偷懶不容易,，唉,，我寫博客實(shí)在是太敬業(yè)了！）,。

背景知識大家可以去自己搜搜,，反正就是黃石公園某個(gè)自然形成的間歇性噴泉，每次噴發(fā)的時(shí)間和等待時(shí)間的散點(diǎn)圖,。挺簡單的對吧,。這個(gè)數(shù)據(jù)點(diǎn)也不多，大家一眼掃過去大概百余個(gè)的樣子,?？墒沁@幅圖真的很有意思，跟Variance的聯(lián)系實(shí)在是太緊密了,。

我們先來說直覺,。如果讓你用自然語言（而非統(tǒng)計(jì)或者數(shù)學(xué)公式）來講述這個(gè)圖，比如你帶這你剛上小學(xué)的孩子去黃石公園玩,，他好奇的在等待這個(gè)噴泉噴發(fā),，這個(gè)時(shí)候你應(yīng)該怎么跟他講？嗯,，可以大概說一下,，“你看基本上呢,，你等的時(shí)間越久，下一次噴發(fā)的時(shí)間就越長哦,。讓我們一起來計(jì)時(shí)~” 然后小朋友看了一眼圖,，不服的說到，“什么嘛,，明明是等了如果超過一小時(shí)（或70分鐘）,，那么下一次基本上噴發(fā)時(shí)間才會長(4-5分鐘）?！?。那么到底哪種說法更準(zhǔn)確呢？

（吐槽君：自從上次寫了樂高機(jī)器人之后,，落園的段子里面的科普對象就從同學(xué)們降低到小朋友了,，喵~）

好啦，不跟小朋友玩了,，我們得研究一下更fancy一點(diǎn)的方法,。說白了就是,，用一個(gè)什么樣的模型可以讓經(jīng)過我們模型處理的Variance盡量的小呢,？

嗯，同學(xué)們說試試回歸唄,，明顯的正相關(guān)性啊,。你都說了啊，X的增加Y也在增加（先不要理會因果即X和Y誰先誰后,，我們只說相關(guān)）,。

所以我們雖然得到了一個(gè)顯著的正相關(guān)關(guān)系，但是回歸模型的R方只有81%（當(dāng)然已經(jīng)很好了,，只是我們希望更好嘛）,。來看看對應(yīng)的殘差分布：

殘差好像挺散的（最理想的殘差就是白噪音了，說明沒有任何信息量了）,，但是隱隱約約也能看出來有兩群這樣,。于是很自然的，有同學(xué)說,，我們?nèi)ピ囋嚲垲惏　?/p>

在去試聚類以前,，先說說大名鼎鼎的K-mean算法的一些基石。

上面不是羅嗦了一堆variance在一維情形下的定義嘛,，那么推廣到二維,，也很簡單。

定義二維的中心點(diǎn)：(xˉ,yˉ),，然后就可以定義二維方差： 每個(gè)點(diǎn)到圖中心的距離的平方和,?？磮D就知道了。

藍(lán)色的就是中心點(diǎn),。這里我們就不羅嗦什么均值容易受極值影響之類的了,，那些也是看菜下料的，我們的數(shù)據(jù)看起來還好,。（突然間為什么有種牛郎織女鵲橋相會的即視感…原來古人觀星也是有異曲同工之妙呀,，天空就是一個(gè)大大的散點(diǎn)圖——勿噴，我保證下面不跑題了）

對于一個(gè)線性回歸模型來講,，我們看的就是殘差項(xiàng)的方差——?dú)埐铐?xiàng)方差越大,，表示他們分布的越散，那模型捕捉到的信息就少,。

對于聚類呢,，也可以看相應(yīng)的方差：每個(gè)類里面的點(diǎn)到類中心的距離平方和 -> K-means。K-means雖然是通過迭代來實(shí)現(xiàn)的,，但他的原理大致也是讓二維的二階中心距最?。ㄓ擅恳淮蔚膿p失函數(shù)來決定）。一眼掃過去大概可以分成牛郎織女兩堆星星,，那么我們就聚兩類好了,。顯然不用跑程序我們也知道，聚成兩類之后的組內(nèi)方差和肯定比直接跟中心點(diǎn)算一個(gè)方差要小,。

聚成兩類之后呢,，我們類似定義的殘差就是兩類中每個(gè)點(diǎn)距離其中心點(diǎn)的Y軸距離（因?yàn)闀苯影阎行狞c(diǎn)作為每類的預(yù)測值）。還是畫個(gè)殘差圖看看,。

紅色是K-MEANS給出的殘差,，藍(lán)色是回歸給出的殘差。貌似這兩個(gè)長得還是挺像的,，也是左右兩群,，雖然每群中兩者長得都不太一樣…這個(gè)時(shí)候我們就可以回到和小朋友的對話了：你們說的都有道理，都有8成的準(zhǔn)確率,，誰也沒比誰更好很多,。

于是我們嘗試在每組內(nèi)再做回歸？會有效果么,？

見效寥寥…引入聚類后回歸模型的R方從81%升到了84%,，才3個(gè)百分點(diǎn)。這主要是在每一類里面,，我們很難找到殘差的規(guī)律了,，所以這樣只是通過組別信息的增加減少了組間方差，而其實(shí)從上圖我們也可以看出每個(gè)組內(nèi)的殘差方差還是很大的,、每條回歸擬合線的斜率從最初的10降低到6和4,，每個(gè)子回歸的R方只有10%不到,，能給予的信息已經(jīng)很少了，所以整體模型只是增加了一點(diǎn)點(diǎn)準(zhǔn)確性,。但是無論如何也比直接玩回歸的效果要好（之所以用k-means而不是簡單粗暴的用類似x>3.5這樣來分成兩類,，是因?yàn)閗-means是按照其損失函數(shù)優(yōu)化過的、給出的是最優(yōu)的兩類聚類結(jié)果）,。

問題來了,，為什么我只聚成兩類而不是三類五類甚至更多呢？主要是怕過擬合,，數(shù)據(jù)才200多個(gè)點(diǎn),，聚太多的話很容易過擬合了。大數(shù)據(jù)的情況下可以嘗試其他辦法,。

好了,，廢話了這么多，其實(shí)統(tǒng)計(jì)學(xué)家們已經(jīng)不僅僅玩這些了,。其實(shí)我們上面的殘差圖傳達(dá)了一個(gè)很重要的信息就是殘差的分布不是一個(gè)白噪聲（或者說不是均值為0,、方差為常數(shù)的正態(tài)分布），稱之為異方差(heterogeneity),。異方差有很多很多情形,，最簡單的就是隨著X的增加而增加。還是網(wǎng)上隨便找了個(gè)圖：

異方差的存在使得我們模型的估計(jì)（或者基于訓(xùn)練數(shù)據(jù)預(yù)測）精度有限（先不考慮validate- test那些）,，所以統(tǒng)計(jì)建模常年在跟殘差項(xiàng)的分布作斗爭——反正看到你有什么規(guī)律,，我就可以提取出來更多的信息放到我的模型中。計(jì)量經(jīng)濟(jì)學(xué)家喜歡說殘差項(xiàng)是垃圾桶,，但其實(shí)也是金礦——沒辦法呀，資源有限,，不能太浪費(fèi),。

然后是一些補(bǔ)充的廢話時(shí)間。

殘差和模型：只要模型可以做預(yù)測,，那么殘差就是存在的,，跟模型是線性還是樹形之類無關(guān)。殘差的定義僅僅是真實(shí)的y減去預(yù)測的y（基于train數(shù)據(jù)）,。

殘差和損失函數(shù)：對于最小二乘法來說,，迭代解數(shù)值解的話，損失函數(shù)就是殘差平方和了（截距項(xiàng)的存在保證殘差的均值為0,，所以等同于殘差的方差）,。同理，大部分機(jī)器學(xué)習(xí)算法的損失函數(shù)都跟殘差項(xiàng)相關(guān),，用到二階的話基本就是在最小化殘差的方差,。從這個(gè)角度而言,，統(tǒng)計(jì)建模都是在跟殘差作斗爭。這里和專注于因果推斷的計(jì)量經(jīng)濟(jì)學(xué)家就很不同了——前幾年大家炒得火熱的trade off between consistency and efficiency就大約來自各自的堅(jiān)持,。

殘差,、頻率派與貝葉斯派：上面說的其實(shí)更多是頻率派的東西，我好像沒提過什么 posterior 之類的對吧,。不過從我的理解來看上述的分析在貝葉斯的框架下很多也是等價(jià)的,，所以這個(gè)不牽扯。下面說的有頻率也有貝葉斯派（的觀點(diǎn)）,，大家湊活著理解吧,。

殘差、異方差和穩(wěn)健標(biāo)準(zhǔn)誤：上面說了異方差是我們又愛又恨的——異方差說明我們的模型不夠好,，但是如果能看到異方差規(guī)律的話也給我們進(jìn)一步改進(jìn)的機(jī)會,。在上個(gè)世紀(jì)為了大約的減少異方差的影響，就出來了一個(gè)穩(wěn)健標(biāo)準(zhǔn)誤(robust standard error)的概念——Heteroscedasticity-consistent standard errors,，最著名的就是Eicker–Huber–White了,。Eicker–Huber–White主要是改變了三明治矩陣那兒——常用R里面sandwich()之類的話不會陌生，有肉有面包,。Eicker–Huber–White簡單的說就是用殘差的估計(jì)值放進(jìn)去改變meat（普通最小二乘法假設(shè)meat那里有個(gè)常量）,，以期系數(shù)估計(jì)值的方差會少受影響。但是吧,，Eicker–Huber–White的問題就是沒有試圖去進(jìn)一步分析殘差里面可能的patterns…所謂的廣義最小二乘,、可行廣義最小二乘（FGLS）之類的，同理,，試圖去找到一個(gè)殘差簡單的線性或者非線性的規(guī)律,，然后改進(jìn)模型中從而減小殘差的方差）。另一類有意思的Eicker–Huber–White就是引入聚類的概念——在每一類中殘差為白噪音,，類別間允許異方差,，稱之為cluster robust standard error，只是大家記得啊cluster的數(shù)量趨近于無窮的時(shí)候這貨才是consistent的,，不要像我一樣忘了這個(gè)假設(shè)…

殘差項(xiàng)和降維：一般說來,，引入的跟Y相關(guān)的維度越多，殘差含有的噪音就會越小,。但是出于過擬合的考慮,，我們需要降維（加入regularization即可）。此外很多時(shí)候模型的目的不僅僅是預(yù)測,，而是提取顯著特征,，那么對于特征的選擇也是需要從模型的角度最大化信息量（比如AIC或者BIC準(zhǔn)則，或者直接regularization看懲罰增加的時(shí)候各個(gè)X的衰減速度，然后選擇一個(gè)懲罰程度對應(yīng)的特征估計(jì)值,，但這里并不一定是validate給出的最優(yōu)懲罰程度,，目的不一樣嘛）。

殘差和實(shí)驗(yàn)數(shù)據(jù)分析：異方差的存在不會影響我們對于average treatment effect估計(jì)的一致性——一致性是從實(shí)驗(yàn)設(shè)計(jì)的隨機(jī)劃分來的,，跟估計(jì)精度無關(guān),。可是大家對于簡單的實(shí)驗(yàn)組對照組的對比越來越不滿,，尤其是在實(shí)驗(yàn)對象很多而且彼此之間異質(zhì)性很強(qiáng)的情況下,，單單一個(gè)ATE給出的平均意義上的效果已經(jīng)滿足不了大家的胃口了。于是大家有了heterogenous treatment effect這個(gè)說法——實(shí)驗(yàn)效果本身也是異質(zhì)的,。簡單的來說,，如果除了隨機(jī)分組之外我們還能觀測到另外一個(gè)維度（稱之為控制變量或者covariates），那么分位數(shù)回歸（quantile regression）就可以給出隨著那個(gè)維度的變化實(shí)驗(yàn)效果的變化,?？墒蔷S度很多的時(shí)候呢？先降維,？

說了這么久,，看到heterogenous第一反應(yīng)顯然應(yīng)該是variance啊，要不拿什么來衡量heterogenous呢,？不過其實(shí)吧,，我是想說，有一類很好玩的paper就是試圖看我們還能不能從treatment effect里面找到比ATE更多的信息量,。比如下面這篇：

Ding, P., Feller, A., and Miratrix, L. W. (2015 ). Randomization Inference for Treatment Effect Variation. http:///RzTsAnl

好吧,，到此為止吧，已經(jīng)越來越學(xué)術(shù)了…本來就是隨便記點(diǎn)東西的,。

哦,，還有code…懶得保存R了，直接貼過來吧,。