十四、迴歸的推斷

到目前為止，我們對變數之間關係的分析純粹是描述性的。我們知道如何找到穿過散點圖的最佳直線來繪製。在所有直線中它的估計的均方誤差最小，從這個角度來看，這條線是最好的。

但是，如果我們的資料是更大總體的樣本呢？如果我們在樣本中發現了兩個變數之間的線性關係，那麼對於總體也是如此嘛？它會是完全一樣的線性關係嗎？我們可以預測一個不在我們樣本中的新的個體的響應變數嗎？

如果我們認為，散點圖反映了被繪製的兩個變數之間的基本關係，但是並沒有完全規定這種關係，那麼就會出現這樣的推理和預測問題。例如，出生體重與孕期的散點圖，顯示了我們樣本中兩個變數之間的精確關係；但是我們可能想知道，對於抽樣總體中的所有新生兒或實際中的一般新生兒，這樣的關係是否是真實的，或者說幾乎是正確的。

一如既往，推斷思維起始於仔細檢查資料的假設。一組假設被稱為模型。大致線性的散點圖中的一組隨機性的假設稱為迴歸模型。

迴歸模型

簡而言之，這樣的模型認為，兩個變數之間的底層關係是完全線性的；這條直線是我們想要識別的訊號。但是，我們無法清楚地看到這條線。我們看到的是分散在這條線上的點。在每一點上，訊號都被隨機噪聲汙染。因此，我們的推斷目標是將訊號從噪聲中分離出來。

更詳細地說，迴歸模型規定了，散點圖中的點是隨機生成的，如下所示。

• x和y之間的關係是完全線性的。我們看不到這個“真實直線”，但它是存在的。
• 散點圖通過將線上的點垂直移動，或上或下來建立，如下所示：
• 對於每個x，找到真實直線上的相應點（即訊號），然後生成噪聲或誤差。
• 誤差從誤差總體中帶放回隨機抽取，總體是均值為 0 的正態分佈。
• 建立一個點，橫座標為x，縱座標為“x處的真實高度加上誤差”。
• 最後，從散點圖中刪除真正的線，只顯示建立的點。

基於這個散點圖，我們應該如何估計真實直線？ 我們可以使其穿過散點圖的最佳直線是迴歸線。 所以迴歸線是真實直線的自然估計。

下面的模擬顯示了迴歸直線與真實直線的距離。 第一個面板顯示如何從真實直線生成散點圖。 第二個顯示我們看到的散點圖。 第三個顯示穿過散點圖的迴歸線。 第四個顯示迴歸線和真實直線。

為了執行模擬，請使用三個引數呼叫draw_and_compare函式：真實直線的斜率，真實直線的截距以及樣本量。

執行模擬幾次，用不同的斜率和截距，以及不同的樣本量。 因為所有的點都是根據模型生成的，所以如果樣本量適中，你會看到迴歸線是真實直線的一個良好估計。

# The true line,
# the points created,
# and our estimate of the true line.
# Arguments: true slope, true intercept, number of points

draw_and_compare(4, -5, 10)

實際上，我們當然不會看到真實直線。 模擬結果表明，如果迴歸模型看起來合理，並且如果我們擁有大型樣本，那麼迴歸線就是真實直線的一個良好近似。

真實斜率的推斷

我們的模擬表明，如果迴歸模型成立，並且樣本量很大，則迴歸線很可能接近真實直線。 這使我們能夠估計真實直線的斜率。

我們將使用我們熟悉的母親和她們的新生兒的樣本，來開發估計真實直線的斜率的方法。 首先，我們來看看我們是否相信，迴歸模型是一系列適當假設，用於描述出生體重和孕期之間的關係。

correlation(baby, 'Gestational Days', 'Birth Weight')
0.40754279338885108

總的來說，散點圖看起來相當均勻地分佈在這條線上，儘管一些點分佈在主雲形的周圍。 相關係數為 0.4，迴歸線斜率為正。

這是否反映真實直線斜率為正的事實？ 為了回答這個問題，讓我們看看我們能否估計真實斜率。 我們當然有了一個估計：我們的迴歸線斜率。 這大約是 0.47 盎司每天。

slope(baby, 'Gestational Days', 'Birth Weight')
0.46655687694921522

但是如果散點圖出現的方式不同，迴歸線會有所不同，可能會有不同的斜率。 我們如何計算，斜率可能有多麼不同？

我們需要點的另一個樣本，以便我們可以繪製迴歸線穿過新的散點圖，並找出其斜率。 但另一個樣本從哪裡得到呢？

你猜對了 - 我們將自舉我們的原始樣本。 這會給我們自舉的散點圖，通過它我們可以繪製迴歸線。

自舉散點圖

我們可以通過對原始樣本帶放回地隨機抽樣，來模擬新樣本，它的次數與原始樣本量相同。 這些新樣本中的每一個都會給我們一個散點圖。 我們將這個稱為自舉散點圖，簡而言之，我們將呼叫整個過程來自舉散點圖。

這裡是來自樣本的原始散點圖，以及自舉重取樣過程的四個複製品。 請注意，重取樣散點圖通常比原始圖稀疏一點。 這是因為一些原始的點沒有在樣本中被選中。

估計真實斜率

我們可以多次自舉散點圖，並繪製穿過每個自舉圖的迴歸線。 每條線都有一個斜率。 我們可以簡單收集所有的斜率並繪製經驗直方圖。 回想一下，在預設情況下，sample方法帶放回地隨機抽取，次數與表中的行數相同。 也就是說，sample預設生成一個自舉樣本。

slopes = make_array()
for i in np.arange(5000):
    bootstrap_sample = baby.sample()
    bootstrap_slope = slope(bootstrap_sample, 'Gestational Days', 'Birth Weight')
    slopes = np.append(slopes, bootstrap_slope)
Table().with_column('Bootstrap Slopes', slopes).hist(bins=20)

然後，我們可以使用percentile方法，為真實直線的斜率構建約 95% 置信區間。 置信區間從 5000 個自舉斜率的第 2.5 百分位數，延伸到第 97.5 百分位數。

left = percentile(2.5, slopes)
right = percentile(97.5, slopes)
left, right
(0.38209399211893086, 0.56014757838023777)

用於自舉斜率的函式

讓我們收集我們估計斜率的方法的所有步驟，並定義函式bootstrap_slope來執行它們。 它的引數是表的名稱，預測變數和響應變數的標籤，以及自舉複製品的所需數量。 在每個複製品中，該函式自舉原始散點圖並計算所得迴歸線的斜率。 然後繪製所有生成的斜率的直方圖，並列印由斜率的“中間 95%”組成的區間。

def bootstrap_slope(table, x, y, repetitions):

    # For each repetition:
    # Bootstrap the scatter, get the slope of the regression line,
    # augment the list of generated slopes
    slopes = make_array()
    for i in np.arange(repetitions):
        bootstrap_sample = table.sample()
        bootstrap_slope = slope(bootstrap_sample, x, y)
        slopes = np.append(slopes, bootstrap_slope)

    # Find the endpoints of the 95% confidence interval for the true slope
    left = percentile(2.5, slopes)
    right = percentile(97.5, slopes)

    # Slope of the regression line from the original sample
    observed_slope = slope(table, x, y)

    # Display results
    Table().with_column('Bootstrap Slopes', slopes).hist(bins=20)
    plots.plot(make_array(left, right), make_array(0, 0), color='yellow', lw=8);
    print('Slope of regression line:', observed_slope)
    print('Approximate 95%-confidence interval for the true slope:')
    print(left, right)

當響應變數為出生體重，預測變數為孕期時，我們呼叫bootstrap_slope來找到真實斜率的置信區間，我們得到了一個區間，非常接近我們之前獲得的東西：大約 0.38 到 0.56 盎司每天。

bootstrap_slope(baby, 'Gestational Days', 'Birth Weight', 5000)
Slope of regression line: 0.466556876949
Approximate 95%-confidence interval for the true slope:
0.378663152966 0.555005146304

現在我們有一個函式，可以自動完成估計在迴歸模型中展示斜率的過程，我們也可以在其他變數上使用它。

例如，我們來看看出生體重與母親身高的關係。 更高的女性往往有更重的嬰兒嗎？

迴歸模型似乎是合理的，基於散點圖，但相關性不高。 這隻有大約 0.2。

scatter_fit(baby, 'Maternal Height', 'Birth Weight')

correlation(baby, 'Maternal Height', 'Birth Weight')
0.20370417718968034

像之前一樣，我們使用bootstrap_slope來估計迴歸模型中真實直線的斜率。

bootstrap_slope(baby, 'Maternal Height', 'Birth Weight', 5000)
Slope of regression line: 1.47801935193
Approximate 95%-confidence interval for the true slope:
1.0403083964 1.91576886223

真實斜率的 95% 置信區間，從約 1 延伸到約 1.9 盎司每英寸。

真實斜率可能為 0 嘛？

假設我們相信我們的資料遵循迴歸模型，並且我們擬合迴歸線來估計真實直線。 如果迴歸線不完全是平的，幾乎總是如此，我們將觀察到散點圖中的一些線性關聯。

但是，如果這種觀察是假的呢？ 換句話說，如果真實直線是平的 - 也就是說，這兩個變數之間沒有線性關係 - 我們觀察到的聯絡，只是由於從樣本中產生點的隨機性。

這是一個模擬，說明了為什麼會出現這個問題。 我們將再次呼叫draw_and_compare函式，這次要求真實斜率為 0。我們的目標是，觀察我們的迴歸線是否顯示不為 0 的斜率。

請記住函式draw_and_compare的引數是真實直線的斜率和截距，以及要生成的點的數量。

draw_and_compare(0, 10, 25)

執行模擬幾次，每次保持真實直線的斜率為 0 。你會注意到，雖然真實直線的斜率為 0，但迴歸線的斜率通常不為 0。迴歸線有時會向上傾斜，有時會向下傾斜，每次都給我們錯誤的印象，即這兩個變數是相關的。

為了確定我們所看到的斜率是否真實，我們想測試以下假設：

原假設。真實直線的斜率是 0。

備選假設。真實直線的斜率不是 0。

我們很有條件來實現它。由於我們可以為真實斜率構建一個 95% 的置信區間，我們所要做的就是看區間是否包含 0。

如果沒有，那麼我們可以拒絕原假設（P 值為 5% 的截斷值）。

如果真實斜率的置信區間確實包含 0，那麼我們沒有足夠的證據來拒絕原假設。也許我們看到的斜率是假的。

我們在一個例子中使用這個方法。假設我們試圖根據母親的年齡來估計新生兒的出生體重。根據樣本，根據母親年齡估計出生體重的迴歸線的斜率為正，約為 0.08 盎司每年。

slope(baby, 'Maternal Age', 'Birth Weight')
0.085007669415825132

雖然斜率為正，但是很小。 迴歸線非常接近平的，這就產生了一個問題，真實直線是否是平的。

scatter_fit(baby, 'Maternal Age', 'Birth Weight')

我們可以使用bootstrap_slope來估計真實直線的斜率。 計算表明，真實斜率的約 95% 的自舉置信區間左端為負，右端為正 - 換句話說，區間包含 0。

bootstrap_slope(baby, 'Maternal Age', 'Birth Weight', 5000)
Slope of regression line: 0.0850076694158
Approximate 95%-confidence interval for the true slope:
-0.104335243815 0.272791852339

因為區間包含 0，所以我們不能拒絕原假設，母親年齡與新生兒出生體重之間的真實線性關係的斜率為 0。基於此分析，使用母親年齡作為預測變數，基於迴歸模型預測出生體重是不明智的。

預測區間

迴歸的主要用途之一是對新個體進行預測，這個個體不是我們原始樣本的一部分，但是與樣本個體相似。在模型的語言中，我們想要估計新值x的y。

我們的估計是真實直線在x處的高度。當然，我們不知道真實直線。我們使用我們的樣本點的迴歸線來代替。

給定值x的擬合值，是基於x值的y的迴歸估計。換句話說，給定值x的擬合值就是迴歸線在x處的高度。

假設我們試圖根據孕期天數來預測新生兒的出生體重。我們在前面的章節中看到，這些資料非常適合迴歸模型，真實直線的斜率的 95% 置信區間不包含 0。因此，我們的預測似乎是合理的。

下圖顯示了預測位於迴歸線上的位置。紅線是x = 300。

紅線與迴歸線的相交點的高度是孕期天數 300 的擬合值。

函式fitted_value計算這個高度。像函式的相關性，斜率和截距一樣，它的引數是表的名稱和x和y的列標籤。但是它也需要第四個引數，即x的值，在這個值上進行估算。

def fitted_value(table, x, y, given_x):
    a = slope(table, x, y)
    b = intercept(table, x, y)
    return a * given_x  + b

孕期天數 300 的擬合值約為 129.2 盎司。 換句話說，對於孕期為 300 天的孕婦，我們估計的新生兒體重約為 129.2 盎司。

fit_300 = fitted_value(baby, 'Gestational Days', 'Birth Weight', 300)
fit_300
129.2129241703143

預測的可變性

我們已經開發了一種方法，使用我們樣本中的資料，根據孕期天數預測新生兒的體重。 但作為資料科學家，我們知道樣本可能有所不同。 如果樣本不同，迴歸線也會不一樣，我們的預測也是。 為了看看我們的預測有多好，我們必須瞭解預測的可變性。

為此，我們必須生成新的樣本。 我們可以像上一節那樣，通過自舉散點圖來實現。 然後，我們為每個散點圖的複製品擬合迴歸線，並根據每一行進行預測。 下圖顯示了 10 條這樣的線，以及孕期天數 300 對應的出生體重預測。

lines

slope	intercept	prediction at x=300
0.503931	-21.6998	129.479
0.53227	-29.5647	130.116
0.518771	-25.363	130.268
0.430556	-1.06812	128.099
0.470229	-11.7611	129.308
0.48713	-16.5314	129.608
0.51241	-23.2954	130.428
0.52473	-27.2053	130.214
0.409943	5.22652	128.21
0.468065	-11.6967	128.723

每一行的預測都不相同。 下表顯示了 10 條線的斜率、截距以及預測。

自舉預測區間

如果我們增加重取樣過程的重複次數，我們可以生成預測的經驗直方圖。這將允許我們建立預測區間，使用為斜率建立自舉置信區間時的相同的百分比方法。

讓我們定義一個名為bootstrap_prediction的函式來實現。該函式有五個引數：

• 表的名稱
• 預測變數和響應變數的列標籤
• 用於預測的x的值
• 所需的自舉重複次數

在每次重複中，函式將自舉原始散點圖，並基於x的指定值查詢y的預測值。具體來說，它呼叫我們在本節前面定義的函式fitted_value，來尋找指定x處的擬合值。

最後，繪製所有預測值的經驗直方圖，並列印由預測值的“中間 95%”組成的區間。它還列印基於穿過原始散點圖的迴歸線的預測值。

# Bootstrap prediction of variable y at new_x
# Data contained in table; prediction by regression of y based on x
# repetitions = number of bootstrap replications of the original scatter plot

def bootstrap_prediction(table, x, y, new_x, repetitions):

    # For each repetition:
    # Bootstrap the scatter; 
    # get the regression prediction at new_x; 
    # augment the predictions list
    predictions = make_array()
    for i in np.arange(repetitions):
        bootstrap_sample = table.sample()
        bootstrap_prediction = fitted_value(bootstrap_sample, x, y, new_x)
        predictions = np.append(predictions, bootstrap_prediction)

    # Find the ends of the approximate 95% prediction interval
    left = percentile(2.5, predictions)
    right = percentile(97.5, predictions)

    # Prediction based on original sample
    original = fitted_value(table, x, y, new_x)

    # Display results
    Table().with_column('Prediction', predictions).hist(bins=20)
    plots.xlabel('predictions at x='+str(new_x))
    plots.plot(make_array(left, right), make_array(0, 0), color='yellow', lw=8);
    print('Height of regression line at x='+str(new_x)+':', original)
    print('Approximate 95%-confidence interval:')
    print(left, right)
bootstrap_prediction(baby, 'Gestational Days', 'Birth Weight', 300, 5000)
Height of regression line at x=300: 129.21292417
Approximate 95%-confidence interval:
127.300774171 131.361729528

上圖顯示了基於 5000 次重複的自舉過程，孕期天數 300 的預測出生體重的自舉經驗直方圖。經驗分佈大致是正泰的。

我們已經通過預測的“中間 95%”，即預測的第 2.5 百分位數到第 97.5 百分位數的區間，構建了分數的約 95% 的預測區間。 區間範圍從大約 127 到大約 131。基於原始樣本的預測是大約 129，接近區間的中心。

改變預測變數的值的效果

下圖顯示了孕期天數 285 的 5,000 次自舉預測的直方圖。 基於原始樣本的預測是約 122 盎司，區間範圍從約 121 盎司到約 123 盎司。

bootstrap_prediction(baby, 'Gestational Days', 'Birth Weight', 285, 5000)
Height of regression line at x=285: 122.214571016
Approximate 95%-confidence interval:
121.177089926 123.291373304

請注意，這個區間比孕婦天數 300 的預測區間更窄。 讓我們來調查其原因。

孕婦天數均值約為 279 天：

np.mean(baby.column('Gestational Days'))
279.10136286201021

所以 285 比 300 更接近分佈的中心。 通常，基於自舉樣本的迴歸線，在預測變數的分佈中心附近彼此更接近。 因此，所有的預測值也更接近。 這解釋了預測區間的寬度更窄。

你可以在下面的圖中看到這一點，它顯示了 10 個自舉複製品中每一個的x = 285和x = 300的預測值。 通常情況下，直線在x = 300處比x = 285處相距更遠，因此x = 300的預測更加可變。

注意事項

我們在本章中進行的所有預測和測試，都假設迴歸模型是成立的。 具體來說，這些方法假設，散點圖中的點由直線上的點產生，然後通過新增隨機正態噪聲將它們推離直線。

如果散點圖看起來不像那樣，那麼模型可能不適用於資料。 如果模型不成立，那麼假設模型為真的計算是無效的。

因此，在開始基於模型進行預測，或者對模型引數進行假設檢驗之前，我們首先要確定迴歸模型是否適用於我們的資料。 一個簡單的方法就是，按照我們在本節所做的操作，即繪製兩個變數的散點圖，看看它看起來是否大致線性，並均勻分佈在一條線上。 我們還應該使用殘差圖，執行我們在前一節中開發的診斷。