原文連結：http://tecdat.cn/?p=23759

原文出處：拓端資料部落公眾號

簡介

兩階段最小二乘法（2SLS）迴歸擬合的線性模型是一種常用的工具變數估計方法。

本文的主要內容是將各種標準的迴歸診斷擴充套件到2SLS。

2SLS估計的回顧

我們需要2SLS迴歸的一些基本結果來開發診斷方法，因此我們在此簡單回顧一下該方法。2SLS迴歸是由Basmann(1957)和Theil(引自Theil 1971)在20世紀50年代獨立發明的，他們採取了略微不同但又相當的方法，都在下面描述，以得出2SLS估計器。

我們想估計線性模型y=Xβ+ε，其中y是因變數的n×1觀察向量，X是迴歸因子的n×p矩陣，通常初始列1s為迴歸常數。 β是一個p×1的迴歸係數向量，需要根據資料進行估計，ε是一個n×1的誤差向量，假定其分佈為Nn(0,σ2In)，其中Nn是多變數正態分佈，0是一個n×1的零向量，In是n階單位矩陣。假設X中的一些（也許是全部）迴歸因子是內生的，即它們被認為不獨立於ε的意義。因此，β的普通最小二乘法（OLS）估計值bOLS=(X⊤X)-1X⊤y通常是有偏的，而且不一致。

現在假設我們有另一組獨立於ϵ的q工具變數（IVs）Z，其中q≥p。如果q=p，我們可以直接應用IV來估計β，但如果q>p，我們有更多的IV，而不是我們需要的。簡單地拋棄IVs是低效的，2SLS迴歸是一個通過合理的方式將IVs的數量減少到p的程式。

2SLS的第一階段通過多元普通最小二乘法對模型矩陣X中的所有迴歸變數進行迴歸，得到q×p的迴歸係數矩陣B=（Z⊤Z）-1Z⊤X，以及擬合值Xˆ=ZB。B的列相當於X的每一列對Z的單獨最小二乘迴歸產生的係數。如果X的某些列是外生的，那麼這些列也會出現在Z中，因此，XˆX^中與外生調節器有關的列只是複製了X的相應列。

由於XˆX^的列是Z的列的線性組合，它們（漸進地）與ε不相關，使它們成為估計迴歸方程的合適IV。這個IV步驟是Theil方法中2SLS的第二個階段。

作為一種替代方法，我們可以通過對XˆX^進行OLS迴歸來獲得完全相同的β的估計值b2SLS，產生b2SLS=(Xˆ⊤Xˆ)Xˆ⊤y。這就是巴斯曼的方法，也是 "2SLS "這個名字的由來。

無論我們把第二階段看成是IV估計還是OLS迴歸，我們都可以把這兩個階段合併成一個公式。

這就是sem包中的tsls()函式(Fox, Nie, and Byrnes 2020)所做的，但是從開發迴歸診斷的角度來看，通過兩個不同的OLS迴歸來計算2SLS估計值是有利的。

對2Sls迴歸異常-資料診斷

就我們所知，用2SLS擬合的迴歸模型的診斷是一個相對被忽視的話題，但Belsley, Kuh和Welsch(1980, 266-68)簡要地討論了這個問題。刪除診斷法直接評估每個案例對擬合迴歸模型的影響，方法是刪除案例，重新擬合模型，並注意到迴歸係數或其他迴歸輸出，如殘差標準差，如何變化。

對於有影響的資料，總是可以通過粗暴的計算來獲得案例刪除診斷，即用每個案例依次刪除來重新擬合模型，但這種方法效率低下，因此在大樣本中沒有吸引力。對於某些類別的統計模型，如廣義線性模型（如Pregibon 1981），對個案刪除診斷的計算要求較低的近似值是可用的，而對於線性模型，有效的 "更新 "公式是可用的（如Belsley, Kuh, and Welsch 1980所描述的），允許精確計算個案刪除診斷的。

事實證明，正如Belsley、Kuh和Welsch所指出的，Phillips（1977年，公式15和16）給出了2SLS迴歸的精確更新公式，允許有效地計算個案選擇統計。

其中，b2SLS-i是去除第ii種情況後的2SLS迴歸係數向量，以及

這裡，yi是第i個案例的因變數值，x⊤ixi⊤是模型矩陣X的第i行，z⊤izi⊤是工具變數模型矩陣Z的第i行。

Belsley, Kuh和Welsch特別研究了（用我們的符號）dfbetai=b2SLS-b2SLS-i的值。他們還討論了殘差標準差s-i的刪除值。

然後，Belsley、Kuh和Welsch計算它們對擬合值(和迴歸係數)影響的綜合度量dffits為

其中（如前）x⊤ixi⊤是模型矩陣X的第i行，XˆX^是第二階段迴歸變數的模型矩陣。

讓

代表將y轉換為擬合值的n×n矩陣，yˆ=H∗y。在OLS迴歸中，類似的量是hat矩陣H=X（X⊤X）-1X⊤。Belsley, Kuh和Welsch指出，H∗與H不同，它不是一個正交投影矩陣，將y正交地投影到X的列所跨越的子空間上。特別是，儘管H∗和H一樣，是等值的（H∗=H∗H∗），並且trace(H∗)=ptrace(H∗)=p，但H∗和H不同，是不對稱的，因此它的對角線元素不能被當作槓桿的總結性措施，也就是說，不能被當作hat值。

Belsley, Kuh和Welsch建議簡單地使用第二階段迴歸的hat值。這些是H2=Xˆ(Xˆ⊤Xˆ)-1Xˆ⊤的對角線條目hi=hii。我們在下面討論一些替代方案。

除了hatvalues、dfbeta、s-i和dfits之外，還計算cook距離Di，這基本上是dfits的一個稍有不同的比例版本，它使用總體殘差標準差s來代替刪除的標準差s-i。

因為它們具有相等的方差，並且在正態線性模型下近似於t分佈，所以 studentized殘差對於檢測異常值和解決正態分佈誤差的假設非常有用。studentized殘差與OLS迴歸相類似，定義為

其中ei=yi-x⊤ib 2SLS是第i種情況的因變數殘差。

如前所述，Belsley, Kuh, and Welsch (1980)建議使用第二階段迴歸的 hatvalues。這是一個合理的選擇，但是它有可能遺漏那些在第一階段有高槓杆率但在第二階段迴歸中沒有的案例。讓h(1)i代表第一階段的hatvalues，h(2)i代表第二階段的hatvalues。如果模型包括一個截距，兩組hatvalues都以1/n和1為界，但第一階段的平均hatvalues是q/n，而第二階段的平均hatvalues是p/n。為了使兩個階段的hatvalues具有可比性，我們將每個hatvalues除以其平均值，h(1∗)i=h(1)iq/n；h(2∗)i=h(2)ip/n。然後我們可以把兩階段的hatvalue定義為每種情況下兩者中較大的一個，hi=(p/n)×max(h(1∗)i,h(2∗)i)，或者定義為它們的幾何平均。

異常資料診斷

標準的R迴歸模型通用方法，包括anova()（用於模型比較），predicted()用於計算預測值，model.matrix()（用於模型或第一或第二階段的迴歸），print()，residuals()（有幾種），summary()，update()，和vcov()。

例子

資料在Kmenta（1986年）中用來說明（通過2SLS和其他方法）對線性聯立方程計量經濟學模型的估計。這些資料代表了經濟從1922年到1941年的年度時間序列，有以下變數。

• Q,人均食品消費
• P，食品價格與一般消費價格的比率
• D, 可支配收入
• F, 前一年農民收到的價格與一般消費價格的比率
• A, 年為單位的時間

該資料集很小，我們可以對其進行檢查。

估計以下兩個方程式模型，第一個方程式代表需求，第二個代表供應。

變數D、F和A被視為外生變數，當然常數迴歸因子（一列1）也是如此，而兩個結構方程中的P是內生解釋變數。由於有四個工具變數可用，第一個結構方程有三個係數，是過度識別的，而第二個結構方程有四個係數，是剛剛識別的。

外生變數的數值是真實的，而內生變數的數值是由Kmenta根據模型生成（即模擬）的，引數的假設值如下。

解決內生變數P和Q的結構方程，可以得到模型的簡化形式

Kmenta獨立地從N(0,1)中抽出20個δ1和δ2的值，然後設定ν1=2δ1和

結構方程估計如下（比較Kmenta 1986, 686）。

預設情況下，summary()會輸出2SLS迴歸的三個 "診斷 "測試的結果。這些測試不是本文的重點，所以我們只對它們進行簡單的評論。

• 一個好的工具變數與一個或多個解釋變數高度相關，同時與誤差保持不相關。如果一個內生的迴歸者與工具變數只有微弱的關係，那麼它的係數將被不精確地估計。在弱工具的診斷測試中，我們希望有一個大的測試統計量和小的p值，Kmenta模型中的兩個迴歸方程就是如此。

• 應用於2SLS迴歸中，Wu-Hausman檢驗是對內生性的一種檢驗。如果所有的迴歸者都是外生的，那麼OLS和2SLS的估計都是一致的，並且OLS的估計更有效，但是如果一個或多個迴歸者是內生的，那麼OLS的估計就不一致了。大的檢驗統計量和小的p值，就像在這個例子中一樣，表明OLS估計器是不一致的，因此，2SLS估計器是首選。

• Sargan檢驗是對過度識別的檢驗。也就是說，在一個過度識別的迴歸方程中，如Kmenta的需求方程中，工具變數比要估計的係數多，工具變數有可能提供關於係數值的衝突資訊。因此，大的檢驗統計量和小的Sargan檢驗的pp值表明，該模型被錯誤地指定了。在這個例子中，儘管我們知道（通過資料的構建方式）需求方程是正確的，但我們還是偶然得到了一個適度小的pp值0.084。Sargan檢驗不適用於剛剛確定的迴歸方程，其工具變數和係數的數量相等，如Kmenta的供給方程。

lm "類物件的幾個方法與產生的物件正常工作。例如，物件的plot()方法呼叫了相應的 "lm "方法併產生了可解釋的圖，這裡是Kmenta模型中需求方程的2SLS擬合。

1. par(mfrow=c(2, 2))
2. plot(deq)

然而，在這種情況下，我們更喜歡描述的這些診斷圖的版本。

資料表現良好。例如，在第一個結構方程中，學生化殘差的QQ圖和hatvalues、學生化殘差和庫克cook距離的 "影響圖 "都是不明顯的，除了幾個高槓杆但在一起的案例。

qqPlot

influence

影響圖中的圓圈面積與Cook's D成正比，水平線畫在學生化殘差標度的0和±2處（rstudent=2處的水平線不在圖中），垂直線在2×h¯和3×h¯處。

為了產生一個更有趣的例子，我們將把高槓杆的第20種情況（即1941年）的QQ值從Q20=106.232改為Q20=95，這個值完全在資料中QQ的範圍內，但與其他資料不一致。

然後重複對第一個結構方程的2SLS擬合，並將結果與未被破壞的資料進行比較，發現迴歸係數有很大變化。

compareCoefs(deq, deq1)

有問題的第20個案例（1941年）通過異常資料迴歸診斷法清楚地顯示出來。

qqPlot(deq1)

outlier

influence Plot

avPlots(deq1)

去掉第20種情況，產生的估計係數接近於未被破壞的資料的係數。

compareCoefs(deq, deq1, deq1.20)

估計係數的標準誤差比原來大，因為我們現在有19個而不是20個案例，也因為解釋變數的變化減少了。

發現hatvaues的三種定義是否對這個例子產生了實際的影響，是有一定的意義的。三種hatvalues的散點圖矩陣表明，它們都產生了類似的結果。

最後，讓我們驗證一下刪除診斷的計算結果是否正確。

非線性診斷法

Cook(1993)和Cook and Croos-Dabrera(1998)系統地探討了成分、殘差圖作為非線性診斷的理論屬性。按照這些作者的說法，並把重點放在解釋變數x1上，讓我們假設響應y與x1的部分關係可能是非線性的，由部分迴歸函式f(x1)表示，而y與其他xxs的部分關係是線性的，因此，資料的準確模型是。

我們不知道f()，所以改用工作模型來擬合

在我們的例子中，通過2SLS迴歸，得到估計的迴歸係數a′,b′1,b′2,...,b′k。Cook和Croos-Dabrera的工作表明，只要迴歸估計是一致的，XXS是線性相關的。部分殘差b′1x1+e可以被繪製出來，並對x1x1進行平滑處理，以顯示f()的估計。其中e=y-(a′+b′1x1+b′2x2+⋯b′kxk)是因變數殘差。在實踐中，如果x1和其他xxs之間有很強的非線性關係，或者y與另一個與x1相關的x有非線性關係，那麼分量加殘差圖就會被分解為f()的準確表示。

Fox和Weisberg（2018）將成分加殘差圖擴充套件到更復雜的迴歸模型，例如可以包括互動作用，將偏殘差新增到預測變數效應圖中。這些圖也可以應用於由2SLS迴歸擬合的線性模型。

診斷非線性：一個例子

我們再一次轉向Kmenta的資料和模型的需求方程來說明成分殘差圖，資料再一次表現良好。為一個加法迴歸方程中的所有數字解釋變數構建了分量殘差圖。比如說。

crPlots(deq, smooth=list(span=1))

我們在圖中為區域性加權迴歸loess平滑（Cleveland, Grosse, and Shyu 1992）設定了一個較大的跨度，因為在資料集中只有n=20個案例。跨度的預設值是2/3。在每個面板中，紅線給出的loess平滑度與藍線給出的最小二乘線緊密匹配，藍線代表的是解釋變數方向的擬合迴歸面，左邊是P，右邊是D。因此，兩種偏關係似乎都是線性的。

CERES圖(Cook 1993)，是成分加殘差圖的一個版本，它使用平滑器而不是線性迴歸，因此對預測因子之間的非線性關係更加穩定。

ceresPlots(deq, smooth=list(span=1))

在當前的例子中，這是一個加性模型，我們得到的圖形與之前的基本相同，只是y軸的縮放比例不同。

plot(predictorEffects)

預測效應圖中的藍色陰影區域代表擬合的部分迴歸線周圍95%的置信度包絡。

然而，假設我們對資料擬合了錯誤的模型。

1. deq2 <- update(deq, . ~ I((P - 85)^4/10^5) + D)
2. crPlots(deq2, smooth=list(span=1))

因為max(P)/min(P)=113.49/86.50=1.3的比率不比1大多少，所以我們在把變數提高到4次方之前，從P中減去一個比min(P)略小的數字，以引起擬合部分迴歸曲線中的非線性。變換後的P的成分加殘差圖清楚地反映了由此產生的缺乏擬合，而D的圖仍然是合理的線性。

帶有部分殘差的預測器效應圖顯示了對同一情況的不同看法，它將P而不是轉換後的P放在橫軸上，並揭示了擬合的非線性部分迴歸函式未能捕獲資料的線性模式。

1. plot(predictorEffects(deq2, residuals=TRUE),
2. partial.residuals=list(span=1))

藍線代表擬合模型，紅線代表平滑後的偏殘差；兩條線之間的差異表明缺乏擬合。

非恆定誤差方差

標準的最小二乘法非恆定方差（"異方差"）診斷法可以直接延伸到2SLS迴歸中。例如，我們可以繪製殘差與擬合值的對比圖，以發現前者的變異性隨著後者的水平而變化（通常是增加）的趨勢。對於的模型中的需求方程。

1. plot(fitted(deq), rstudent(deq))
2. abline(h=0)

似乎沒有問題。

Fox(2016)提出了這個圖形的一個變體，改編了Tukey的擴散水平圖(Tukey 1977)，繪製了絕對螺距殘差的對數與擬合值的對數，假設後者是正數。如果該圖的擬合線的斜率為b，則方差穩定化的功率轉換由yλ=y1-b給出。因此，如果b>0，建議的變換是順著Tukey的冪和根的階梯，例如，λ=1-b=1/2代表平方根變換，λ=1-b=0代表對數變換，以此類推。對於模型，我們有

這表明分佈有隨水平提高的輕微趨勢。λ=-2.45的轉換似乎很強，直到我們注意到QQ的值離0很遠，而且最大和最小的值Qmax/Qmin=106.23/92.42=1.15的比率接近1，所以Q-2.45幾乎是QQ的線性轉換，也就是說，實際上根本沒有轉換。

由Breusch和Pagan（1979）提出的最小二乘迴歸中的非恆定誤差方差的普通分數測試，是基於模型的

其中函式g()是未指定的，變數z1,...,zs是誤差方差的預測因子。在最常見的應用中，由Cook和Weisberg（1983）獨立提出，有一個zz，即迴歸的擬合值yˆ，儘管使用初級迴歸中的迴歸者x作為zs也很常見。測試是通過將標準化殘差的平方e2i/σˆ2迴歸到zs上實現的，其中σˆ2=∑e2i/n。然後，在誤差方差不變的無效假設下，該輔助迴歸的迴歸平方和除以2的漸近分佈為χ2s。

Breusch-Pagan/Cook-Weisberg檢驗很容易適應2SLS迴歸。對於需求方程。

在這裡，第一個檢驗是針對擬合值的，第二個更一般的檢驗是針對需求方程中的解釋變數的；這兩個檢驗的p值都很大，表明沒有什麼證據反對恆定方差的假說。

2SLS迴歸中對非恆定方差的補救方法與最小二乘迴歸中的補救方法相似。

• 我們已經提出，如果誤差方差隨著響應水平的提高（或降低），並且因變數是正的，那麼我們就可以通過對因變數進行冪變換來穩定誤差方差。

• 另外，如果我們知道誤差的方差與一個常數成正比，那麼我們可以對2SLS估計使用反方差加權。

• 最後，我們可以在2SLS中使用係數協方差矩陣的估計（或自舉法：例如，見Davison和Hinkley 1997）來修正非恆定誤差方差的標準誤差，就像Huber（1967）和White（1980；也見Long和Ervin 2000）提出的最小二乘法迴歸那樣。

對於例子，非恆定誤差方差的證據是輕微的，標準誤差與傳統的2SLS標準誤差相似，甚至比它略小。

我們將修改資料以反映非恆定誤差方差，像最初那樣從還原形式方程中重新生成資料，將內生變數P和Q表示為外生變數D、F和A的函式，以及還原形式誤差ν1和ν2。

3. w <- 3*(EQ - min(EQ) + 0.1)/(max(EQ) - min(EQ))
4. v1 <- v1*w # 非恆定方差
5. Q <- EQ + v1
6. P <- EP + v2

將取樣的簡化形式誤差v1與Q的期望值作圖，顯示出明顯的異質性模式。

 plot(EQ, v1))

然後將需求方程與新的資料集重新匹配，我們得到

而非恆定的誤差方差在診斷中得到了明顯的反映；例如：

ncvTest(deq2)

從建議的QQ的冪變換的極值λ=-23，反映了（正如我們之前指出的）max(Q)/min(Q)並沒有比1大多少。

在我們的例子中，標準誤差與傳統的標準誤差沒有太大區別。

如前所述，bootstrapping提供了一種替代標準誤差的方法，作為對非恆定誤差方差的修正，實現了個案再抽樣的bootstrapping，並返回一個適合與boot包中的函式一起使用的 "boot "類物件（Davison和Hinkley 1997；Canty和Ripley 2020）。比如說。

b.deq2 <- Boot(deq2)

在這個例子中，bootstrap的標準誤差比傳統的標準誤差大。

Boostrap置信區間也可以從Boot()返回的物件中計算出來，預設情況下報告BCa（偏差校正，加速）區間（。

confint(b.deq2)

加權的2SLS迴歸

假設我們修改了迴歸模型y=Xβ+ε，那麼現在Nn(0,σ2W-1)其中W=diag{wi}是一個已知反方差權重的n×n對角矩陣；即V(εi)=σ2/wi。如前所述，X的某些列可能與誤差ε相關，但我們有足夠的工具變數Z，它們與誤差無關。

那麼加權的2SLS估計是

或者，我們可以將2SLS的兩個階段視為加權最小二乘法（WLS）問題，在每個階段都要最小化加權殘差平方和。

Phillips的2SLS迴歸的更新公式也可以針對加權情況進行修改，但是更簡單的方法（是將加權2SLS問題轉換成非加權問題，通過使用W1/2=diag{wi--√}，W的Cholesky平方根將資料轉換為恆定方差。W的平方根特別簡單，因為W是對角線的。然後在Phillips的更新公式中，我們用y∗=W1/2y代替y，用X∗=W1/2X代替X，用Z∗=W1/2Z代替Z。

對於修改後的資料，我們知道需求方程的誤差方差與變數w成反比。這反映了資料的構建方式。因此，加權的2SLS估計被計算為

將求和殘差與擬合值作圖，並測試非恆定誤差方差，並不表明有異方差問題，但有一個相對較大的求和殘差，約為-3，與其他數值相比有些突出。

plot(fitted(deqw), rstudent(deqw))

Bonferroni離群值測試表明，學生化殘差並不是異常大，模型是正確的。

outlierTest

共線性關係診斷

除了異常資料診斷外，Belsley, Kuh和Welsch（1980）還簡要地將他們的共線性關係診斷方法擴充套件到2SLS迴歸中。我們認為，這種將共線性關係同化為數值不穩定性的方法是有缺陷的，因為它考慮到了 "與截距的共線性關係"。也就是說，數值遠離0的迴歸者與常數迴歸者的乘積之和很大，產生了截距的大標準誤差，只是反映了截距將擬合的迴歸面遠遠超出了資料範圍的事實。

Fox和Monette(1992)描述了一種基於廣義方差膨脹因子的最小二乘法擬合的線性模型中串聯性診斷的替代方法。廣義方差膨脹因子採用了係數的估計協方差矩陣，一般適用於有線性預測因子的模型，包括由2SLS估計的線性模型。

例如，對於模型中的需求方程。

sqrt(vif(deq))

取VIF的平方根將它們放在係數標準誤差的刻度上。也就是說，P和D的係數的標準誤差比估計的係數不相關時要大23%。像這裡一樣，模型中的每個項只有一個係數時，廣義和普通方差膨脹因子是一致的。P和D的VIFs相等是兩個迴歸變數（超越迴歸常數）的情況下所特有的。

邊際/條件圖是由car軟體包中的mcPlots()函式生成的，它將新增變數圖疊加到相應的迴歸因子邊際散點圖上。這些圖使我們能夠直觀地看到由於共線性關係造成的每個係數估計精度的降低，共線性關係降低了迴歸變數相對於其邊際變化的條件變化。 例如，對於需求方程。

Plots(deq)

每個面板中的藍色點代表邊際散點圖，紅色點代表（部分）新增變數圖，箭頭顯示兩組點之間的關係。

結語

仔細的迴歸分析需要有效地觀察資料的方法。許多潛在的問題可以通過在擬合迴歸模型之前檢查資料來解決，減少（如果不是消除）擬合後診斷的必要性。毫無疑問，採用2SLS的謹慎的資料分析者一直都是這樣做的。然而，擁有允許人們對2SLS擬合的迴歸模型進行批評的方法，至少在某些情況下會建議對模型進行改進，或者對資料進行修正。

參考文獻

Basmann, R. L. 1957. “A Generalized Classical Method of Linear Estimation of Coefficients in a Structural Equation.”Econometrica25: 77–83.https://doi.org/10.2307/1907743.

Belsley, D. A., E. Kuh, and R. E. Welsch. 1980.Regression Diagnostics: Identifying Influential Data and Sources of Collinearity. New York: John Wiley & Sons.Regression Diagnostics | Wiley Series in Probability and Statistics.

Breusch, T. S., and A. R. Pagan. 1979. “A Simple Test for Heteroscedasticity and Random Coefficient Variation.”Econometrica47: 1287–94.https://doi.org/10.2307/1911963.

最受歡迎的見解

1.R語言多元Logistic邏輯迴歸 應用案例

2.面板平滑轉移回歸(PSTR)分析案例實現

3.matlab中的偏最小二乘迴歸（PLSR）和主成分迴歸（PCR）

4.R語言泊松Poisson迴歸模型分析案例

5.R語言迴歸中的Hosmer-Lemeshow擬合優度檢驗

6.r語言中對LASSO迴歸，Ridge嶺迴歸和Elastic Net模型實現

7.在R語言中實現Logistic邏輯迴歸

8.python用線性迴歸預測股票價格

9.R語言如何在生存分析與Cox迴歸中計算IDI，NRI指標