Machine Learning With Go 第4章：迴歸

4 迴歸

之前有轉載過一篇文章：容量推薦引擎：基於吞吐量和利用率的預測縮放，裡面用到了基本的線性迴歸來預測容器的資源利用情況。後面打算學一下相關的知識，譯自：Machine Learning With Go

我們將探究的第一組機器學習技術通常被稱為迴歸(regression)，我們可以將回歸理解為一個變數(例如銷售額)的變化是如何影響到其他變數(如使用者數)的。對於機器學習技術來說，這是一個很好的開端，它們是構成其他更加複雜技術的基礎。

機器學習中的迴歸技術通常會注重評估連續值(如股票價格、溫度或疾病進展等)。下一章討論的歸類(Classification)則會注重離散值，或離散的類集合(如欺詐/非欺詐、坐下/起立/跑動，或熱狗/非熱狗等)。正如上面提到的，迴歸技術會貫徹到機器學習中，並作為歸類演演算法的一部分，但本章中，我們將會注重其基本的應用--預測連續值。

理解迴歸模型術語

正如前面提到的，迴歸本身是一個分析一個變數和另一個變數之間關係的過程，但在機器學習中還用到了一些術語來描述這些變數以及各種型別的迴歸和與迴歸有關的過程：

響應(response)或因變數(dependent variable)：這些術語可以互用，表示基於其他一個或多個變數來預測的變數，通常使用y表示
解釋變數(Explanatory variables)、自變數(independent variables)、特徵(features)、屬性(attributes)或迴歸係數(regressors)：這些術語可以互用，表示用於預測響應的變數，通常使用x或x1 , x2表示
線性迴歸：該型別的迴歸會假設因變數會線性依賴自變數(即遵循直線方程)
非線性迴歸：該型別的迴歸會假設因變數會非線性依賴自變數(如多項式或指數)
多元迴歸(Multiple regression:)：具有超過一個自變數的迴歸
擬合(Fitting)或訓練(training)：引數化模型的過程(如迴歸模型)，可以用該模型來預測一個特定的因變數
預測：使用引數模型預測因變數的過程(如迴歸模型)

部分術語會在迴歸上下文和本書的其他上下文中使用。

線性迴歸

線性迴歸是最簡單的機器學習模型之一，但不能出於某些原因而忽略該模型。正如前面提到的，它是其他模型的基礎，且有一些非常重要的優勢。

正如在本書中討論的，完整性在機器學習應用中非常重要，模型越簡單，解釋性越強，則越容易維護其完整性。此外，如果模型簡單且具有解釋性，那麼就可以幫助理解變數之間的推斷關係，並簡化開發過程中的檢查工作。

Mike Lee Williams (來自Fast Forward Labs的)說過:

未來是演演算法，可解釋模型在人類和智慧機器之間建立了一種更安全、更高效、最終更具共同作業性的關係。

線性迴歸模型是可解釋的，因此可以為資料科學提供一種安全且高效的選項。當需要搜尋一種可以預測連續變數的模型時，如果資料和相關條件具備，則應該考慮並使用線性迴歸(或多元線性迴歸)。

線性迴歸概述

線上性迴歸中，我們會嘗試使用如下線性方程，使用一個自變數x，對因變數y進行建模：

\[y = mx + b \]

這裡，m為直線的斜率，b為截距。例如，我們想要通過每天存取網站的users對sales進行建模，為了使用線性迴歸，我們會希望通過確定m和b來讓預測公司的銷售額：

\[sales = m * (number~of~users) + b \]

這樣，我們的訓練模型就是該引數化函數。通過輸入Number of Users 來預測 sales，如下：

線性迴歸的訓練或擬合需要確定m和b的值，這樣得出的公式就有預測響應的能力。有多種方式來確定m和b，但最常見的是普通最小二乘法(ordinary least squares (OLS))。

為了使用OLS來確定m和b，首先為m和b選擇一個值來建立第一條範例線(example line)。然後測量每個已知點(如訓練集)和範例線之間的垂直距離，這些距離稱為誤差(errors)或殘差(residuals)。下圖展示了評估和驗證：

下面，我們計算這些誤差平方和:

\[\frac{error^2_1+error^2_2+...+error^2_N}{N} \]

通過調整m和b來最小化誤差的平方和。換句話說，我們訓練的線性迴歸直線是平方和最小的直線。

有很多種方式可以找出誤差平方和最小的直線，如通過OLS可以找出並分析這條直線。但最常用的減少誤差平方和的優化技術稱為梯度下降法(gradient descent)。相比於分析法，這種方法更容易實現，且便於計算(如記憶體)，也更加靈活。

可以說，線性迴歸和其他迴歸的實現都利用梯度下降來擬合或訓練線性迴歸線。實際上，梯度下降法在機器學習中無處不在，由此可以產生更加複雜的模型技術，如深度學習。

梯度下降法

梯度下降法有很多變種，且在機器學習世界中無處不在。最重要的是，它們用於確定線性或邏輯迴歸等演演算法的最佳係數，同時也在更復雜的技術中發揮著重要作用(至少部分基於線性/邏輯迴歸(如神經網路))。

梯度下降法的一般思想是確定某些引數的變化方向和幅度，這些引數將使預測曲線朝著正確的方向移動，以優化某些度量(如誤差)。想象站在某個地方，如果要向較低的位置移動，則需要朝向下的方向移動。這基本上就是梯度下降演演算法在優化引數時所做的事情。

讓我們看一下所謂的隨機梯度下降(SGD)，這是一種增量的梯度下降，從而對這個過程有更多直覺上的瞭解。我們在第5章"分類"的邏輯迴歸實現中使用了SGD。在該範例中，我們實現了對邏輯迴歸引數的訓練或擬合，如下所示：
// logisticRegression fits a logistic regression model
// for the given data.
func logisticRegression(features *mat64.Dense, labels []float64, numSteps int, learningRate f)
// Initialize random weights.
_, numWeights := features.Dims()
weights := make([]float64, numWeights)
s := rand.NewSource(time.Now().UnixNano())
r := rand.New(s)
for idx, _ := range weights {
  weights[idx] = r.Float64()
}
// Iteratively optimize the weights.
for i := 0; i < numSteps; i++ {
  // Initialize a variable to accumulate error for this iteration.
  var sumError float64
  // Make predictions for each label and accumulate error.
  for idx, label := range labels {
      // Get the features corresponding to this label.
      featureRow := mat64.Row(nil, idx, features)
      // Calculate the error for this iteration's weights.
      pred := logistic(featureRow[0]*weights[0] + featureRow[1]*weights[1])
      predError := label - pred
      sumError += math.Pow(predError, 2)
      // Update the feature weights.
      for j := 0; j < len(featureRow); j++ {
          weights[j] += learningRate * predError * pred * (1 - pred) * featureRow[j]
      }
  }
}
return weights
}
// Iteratively optimize the weights註釋下面的迴圈實現了通過SGD來優化邏輯迴歸引數。下面選擇這部分迴圈來看下到底發生了什麼。

首先，我們使用當前權重和預測值與理想值(即實際觀察值)之間的差值來計算模型的輸出：
// Calculate the error for this iteration's weights.
pred := logistic(featureRow[0]*weights[0] + featureRow[1]*weights[1])
predError := label - pred
sumError += math.Pow(predError, 2)
根據SGD，我們將根據如下公式來計算引數(在本例中為權weights)的更新：

\[update=leaning~rate\times~gradient~of~the~parameters \]
gradient是cost函數的數學梯度。

更多參見: http://mathworld.wolfram.com/Gradient.html

然後將該更新應用到引數，如下所示：

\[parameter=parameters-update \]
在我們的邏輯迴歸模型中，計算結果如下：
// Update the feature weights.
for j := 0; j < len(featureRow); j++ {
  weights[j] += learningRate * predError * pred * (1 - pred) * featureRow[j]
}
機器學習中廣泛使用了這種型別的SGD，但在某些場景下，這種梯度下降法可能導致過擬合或陷入區域性最小值/最大值（而不是尋找全域性最優值）。

為了解決這些問題，可以使用一個梯度下降的變種，稱為批次梯度下降(batch gradient descent)。在批次梯度下降中，可以基於所有訓練資料集中的梯度來計算每個引數更新，而不針對資料集的特定觀測值或行。這種方式有助於防止過度擬合，但它也可能很慢，並且存在記憶體問題，因為需要計算每個引數相對於整個資料集的梯度。微批次梯度下降(Mini-batch gradient descent)是另一個變種，它在試圖保持批次梯度下降的某些好處的同時，更易於計算。在微批次梯度下降法中，梯度是在訓練資料集的子集上計算的，而不是在整個訓練資料集上計算的。

在邏輯迴歸的場景中，你可能看到過使用梯度上升或下降，梯度上升與梯度下降是一回事，只是cost函數的方向不同而已。更多參見: https://stats.stackexchange.com/questions/261573/using-gradient-ascent-instead-of-gradient-descent-for-logistic-regression

gonum團隊已經實現了梯度下降法：gonum.org/v1/gonum/optimize。檔案地址：https://pkg.go.dev/gonum.org/v1/gonum/optimize

線性迴歸的假設和缺點

與所有機器學習模型一樣，線性迴歸並不能適用於所有場景，它的前提是假設你的資料之間的關係是確定的：

線性關係：線性迴歸會假設因變數線性依賴自變數(線性方程)。如果這種關係不是線性的，則線性迴歸可能會表項不佳
正態性：假設變數遵循正太分佈(看起來像鐘形)。本章後面會討論這種特性以及非正態分佈變數下的取捨。
非多重共線性：多重共線性是一個特別的術語，它意味著自變數並不是真正獨立的，它們會以某種形式相互依賴
沒有自相關性：自相關性是另一個特別的術語，意味著變數依賴自身或自身的某個版本(如存在某些可預測的時序中)。
同方差性：這可能是這一組術語中最特別的一個，但它相對比較簡單，且並不需要經常關注。線性迴歸假設迴歸線周圍的資料的方差與自變數值的方差大致相同。

從技術上講，為了使用線性迴歸，需要滿足上述所有假設。但最重要的是我們需要知道資料是如何分佈的，以及它們是如何表現的。後續在使用線性迴歸的範例中會深入討論這些假設。

作為一個資料科學家或分析師，在使用線性迴歸時需要注意到線性迴歸的不足：

使用特定範圍的自變數來訓練線性迴歸模型，在預測該範圍外的資料時應該格外小心，因為你的線性迴歸直線可能並不適用(如，某些極端數值下，因變數可能並不是線性的)。
可能為兩個並無關聯的變數建立了一個線性迴歸模型。需要確保變數之間有邏輯上的關聯性。
可能會因為擬合某些特定型別資料中的異常或極端值而偏離迴歸線，如OLS。有一些方式可以讓擬合迴歸不受異常值的影響，或針對異常值展示出不同的行為，如正交最小二乘法(orthogonal least squares)或嶺迴歸(ridge regression)。

線性迴歸的例子

為了描述線性迴歸，讓我們建立第一個機器學習模型。下面資料為廣告資料，儲存格式為.csv：

$ head Advertising.csv
TV,Radio,Newspaper,Sales
230.1,37.8,69.2,22.1
44.5,39.3,45.1,10.4
17.2,45.9,69.3,9.3
151.5,41.3,58.5,18.5
180.8,10.8,58.4,12.9
8.7,48.9,75,7.2
57.5,32.8,23.5,11.8
120.2,19.6,11.6,13.2
8.6,2.1,1,4.8

該資料集包括一系列表示廣告媒體屬性(TV，Radio和Newspaper)以及對應的銷售額(Sales)，本例中我們的目標是對銷售額(因變數)和廣告支出(因變數)進行模型。

分析資料

為了構建模型(或流程)，並確保能夠對模型的結果進行檢查，首先需要對資料進行分析(所有機器學習模型的第一個步驟)。我們需要了解變數是如何分佈的，以及變數的範圍和可變性。

為了實現該目標，我們將計算第2章矩陣、概率和統計中討論的彙總資料。這裡，我們將使用github.com/go-gota/gota/tree/master/dataframe中的內建方法，一次性計算出資料集中的所有列的彙總資訊：

// Open the CSV file.
advertFile, err := os.Open("Advertising.csv")
if err != nil {
    log.Fatal(err)
}
defer advertFile.Close()
// Create a dataframe from the CSV file.
advertDF := dataframe.ReadCSV(advertFile)
// Use the Describe method to calculate summary statistics
// for all of the columns in one shot.
advertSummary := advertDF.Describe()
// Output the summary statistics to stdout.
fmt.Println(advertSummary)

編譯並執行後得到如下結果：

$ go build
$ ./myprogram
[7x5] DataFrame
   column   TV         Radio     Newspaper  Sales
0: mean     147.042500 23.264000 30.554000  14.022500
1: stddev   85.854236  14.846809 21.778621  5.217457
2: min      0.700000   0.000000  0.300000   1.600000
3: 25%      73.400000  9.900000  12.600000  10.300000
4: 50%      149.700000 22.500000 25.600000  12.900000
5: 75%      218.500000 36.500000 45.100000  17.400000
6: max      296.400000 49.600000 114.000000 27.000000
   <string> <float>    <float>   <float>    <float>

上面以表格形式列印出所有的彙總資料，包括平均值、標準偏差、最小值、最大值、25%/75%百分位和中位數(或50%百分位)。

這些值為我們提供了良好的數值參考，後續會在訓練線性迴歸模型時將看到這些數位。但缺乏直觀上的理解，為此，我們需要為每列數值建立一個直方圖：

// Open the advertising dataset file.
f, err := os.Open("Advertising.csv")
if err != nil {
    log.Fatal(err)
}
defer f.Close()
// Create a dataframe from the CSV file.
advertDF := dataframe.ReadCSV(f)
// Create a histogram for each of the columns in the dataset.
for _, colName := range advertDF.Names() {
    // Create a plotter.Values value and fill it with the
    // values from the respective column of the dataframe.
    plotVals := make(plotter.Values, advertDF.Nrow())
    for i, floatVal := range advertDF.Col(colName).Float() {
        plotVals[i] = floatVal
    }
    // Make a plot and set its title.
    p, err := plot.New()
    if err != nil {
        log.Fatal(err)
    }
    p.Title.Text = fmt.Sprintf("Histogram of a %s", colName)
    // Create a histogram of our values drawn
    // from the standard normal.
    h, err := plotter.NewHist(plotVals, 16)
    if err != nil {
        log.Fatal(err)
    }
    // Normalize the histogram.
    h.Normalize(1)
    // Add the histogram to the plot.
    p.Add(h)
    // Save the plot to a PNG file.
    if err := p.Save(4*vg.Inch, 4*vg.Inch, colName+"_hist.png"); err != nil {
        log.Fatal(err)
    }
}

本程式會為每個直方圖建立一個.png影象：

觀察上圖以及計算出的彙總資訊，下一步考慮是否符合線性迴歸的假設條件。可以看到並不是所有的變數都是正態分佈的(鐘形的)。可以看到銷售額是鐘形的，而其他則不是正態的。

我們可以使用分點陣圖(quantile-quantile (q-q) p)統計工具來確定分佈與正態分佈的接近程度，甚至通過統計測試來確定變數是否服從正態分佈的概率。但大多數情況下，通過圖表就可以得出一個大致的結論。

下一步要做出決策，但至少有一部分資料在技術上並不會擬合到我們的線性迴歸模型中，可以選擇如下一種方式進行處理：

嘗試轉換變數，使其遵循正態分佈，併線上性迴歸模型中使用這些轉換的變數。這種方式的好處是可以在模型的假設中進行操作，缺點是可能會讓模型難以理解，降低可解釋性
使用不同的資料來解決問題
線上性迴歸假設中忽略該問題，並嘗試建立該模型

可能還有其他解決問題的方式，但我的建議是首先嚐試第三種選項。由於可以快速地訓練線性迴歸模型，因此該選項並不會帶來多少壞處。如果最後得出滿意的模型，那麼就可以避免引入更多的複雜性。如果得到的模型不盡如人意，那麼此時再訴諸於其他選項。

選擇自變數

現在對我們的資料有了一些直覺上的瞭解，並且已經瞭解到資料是如何擬合線性迴歸模型的假設的。那麼現在應該選擇哪個變數作為我們的自變數來預測因變數？

最簡單的方法是通過直觀地探索因變數和選擇的所有自變數之間的相關性，特別是可以通過繪製因變數與其他每個變數的散點圖(使用pkg.go.dev/gonum.org/v1/plot)來做決定：

// Open the advertising dataset file.
f, err := os.Open("Advertising.csv")
if err != nil {
    log.Fatal(err)
}
defer f.Close()
// Create a dataframe from the CSV file.
advertDF := dataframe.ReadCSV(f)
// Extract the target column.
yVals := advertDF.Col("Sales").Float()
// Create a scatter plot for each of the features in the dataset.
for _, colName := range advertDF.Names() {
    // pts will hold the values for plotting
    pts := make(plotter.XYs, advertDF.Nrow())
    // Fill pts with data.
    for i, floatVal := range advertDF.Col(colName).Float() {
        pts[i].X = floatVal
        pts[i].Y = yVals[i]
    }
    // Create the plot.
    p, err := plot.New()
    if err != nil {
        log.Fatal(err)
    }
        p.X.Label.Text = colName
    p.Y.Label.Text = "y"
    p.Add(plotter.NewGrid())
    s, err := plotter.NewScatter(pts)
    if err != nil {
        log.Fatal(err)
    }
    s.GlyphStyle.Radius = vg.Points(3)
    // Save the plot to a PNG file.
    p.Add(s)
    if err := p.Save(4*vg.Inch, 4*vg.Inch, colName+"_scatter.png"); err != nil {
        log.Fatal(err)
    }
}

如此可以建立如下散點圖：

通過這些散點圖，我們需要推斷出哪些屬性 (TV, Radio, 和/或 Newspaper)與我們的因變數(Sales)具有線性關係。是否可以在這些散點圖上畫一條線，以符合銷售趨勢和各自的屬性？這種方法並不總是行得通，且對於一個特定的問題，並不一定可以將其關聯到所有的屬性。

上述場景中，Radio和TV與Sales呈線性關係，Newspaper可能與Sales有一定的關係，但相關性並不明顯。與TV的相關性是最明顯的，因此先選擇TV作為線性迴歸模型的自變數，線性迴歸公式如下：

\[Sales = m~TV+b \]

這裡要注意的另一件事是，變數TV可能不是嚴格等方差的(線上性迴歸的假設中討論過)。這一點需要注意(可能值得在專案中歸檔的)，下面將繼續探究是否可以建立具有預測能力的線性迴歸模型。當模型表現不佳時，需要重新審視這種假設。

建立訓練和測試集

為了避免過度擬合併保證模型的推廣，我們需要將資料集劃分為訓練集和測試集即評估和驗證(Evaluation and Validation)。這裡我們不會聚焦某個測試集，因為只需要進行一次模型訓練即可，而不會在訓練和測試之間來回迭代。但如果需要多個因變數進行驗證和/或需要迭代調整模型引數時，你可能希望建立一個保留集，儲存到模型開發過程結束後進行驗證。

我們將使用github.com/go-gota/gota/blob/master/dataframe建立訓練和測試資料集，並將它們儲存到各自的.csv檔案中。該場景中，我們使用80/20的比例來劃分訓練和測試資料：

// Open the advertising dataset file.
f, err := os.Open("Advertising.csv")
if err != nil {
    log.Fatal(err)
}
defer f.Close()
// Create a dataframe from the CSV file.
// The types of the columns will be inferred.
advertDF := dataframe.ReadCSV(f)
// Calculate the number of elements in each set.
trainingNum := (4 * advertDF.Nrow()) / 5
testNum := advertDF.Nrow() / 5
if trainingNum+testNum < advertDF.Nrow() {
    trainingNum++
}
// Create the subset indices.
trainingIdx := make([]int, trainingNum)
testIdx := make([]int, testNum)
// Enumerate the training indices.
for i := 0; i < trainingNum; i++ {
    trainingIdx[i] = i
}
// Enumerate the test indices.
for i := 0; i < testNum; i++ {
    testIdx[i] = trainingNum + i
}
// Create the subset dataframes.
trainingDF := advertDF.Subset(trainingIdx)
testDF := advertDF.Subset(testIdx)
// Create a map that will be used in writing the data
// to files.
setMap := map[int]dataframe.DataFrame{
    0: trainingDF,
    1: testDF,
}
// Create the respective files.
for idx, setName := range []string{"training.csv", "test.csv"} {
    // Save the filtered dataset file.
    f, err := os.Create(setName)
    if err != nil {
        log.Fatal(err)
    }
    // Create a buffered writer.
    w := bufio.NewWriter(f)
    // Write the dataframe out as a CSV.
    if err := setMap[idx].WriteCSV(w); err != nil {
        log.Fatal(err)
    }
}

上述程式碼會輸出如下訓練和測試集：

$ wc -l *.csv
    201 Advertising.csv
    41  test.csv
    161 training.csv
    403 total

這裡使用的資料並沒有經過排序。但如果需要按照響應、日期或其他方式處理資料，則最好隨機劃分訓練和測試集。如果不這麼做，訓練和測試集可能會包含特定範圍的響應，這樣響應可能會受到時間/日期等人為因素的影響。

訓練模型

下面將訓練(或擬合)我們的線性迴歸模型。這也意味著需要找到誤差平方和最小的的斜率(m)和截距(b)。為了執行訓練，我們會使用來自Sajari的包：github.com/sajari/regression。Sajari是一個重度依賴Go和機器學習的網站搜尋公司，他們在生產中使用了github.com/sajari/regression。

為了使用github.com/sajari/regression來訓練迴歸模型，需要初始化一個regression.Regression值，並設定一對標籤，然後使用被標記的訓練資料來填充regression.Regression。之後就可以簡單地使用呼叫Run()來對regression.Regression的值進行訓練，以此生成線性迴歸模型。

// Open the training dataset file.
f, err := os.Open("training.csv")
    if err != nil {
    log.Fatal(err)
}
defer f.Close()
// Create a new CSV reader reading from the opened file.
reader := csv.NewReader(f)
// Read in all of the CSV records
reader.FieldsPerRecord = 4
trainingData, err := reader.ReadAll()
if err != nil {
    log.Fatal(err)
}
// In this case we are going to try and model our Sales (y)
// by the TV feature plus an intercept. As such, let's create
// the struct needed to train a model using github.com/sajari/regression.
var r regression.Regression
r.SetObserved("Sales")
r.SetVar(0, "TV")
// Loop of records in the CSV, adding the training data to the regression value.
for i, record := range trainingData {
    // Skip the header.
    if i == 0 {
        continue
    }
    // Parse the Sales regression measure, or "y".
    yVal, err := strconv.ParseFloat(record[3], 64)
    if err != nil {
        log.Fatal(err)
    }
    // Parse the TV value.
    tvVal, err := strconv.ParseFloat(record[0], 64)
    if err != nil {
        log.Fatal(err)
    }
    // Add these points to the regression value.
    r.Train(regression.DataPoint(yVal, []float64{tvVal}))
}
// Train/fit the regression model.
r.Run()
// Output the trained model parameters.
fmt.Printf("\nRegression Formula:\n%v\n\n", r.Formula)

編譯執行上述程式碼，訓練的線性迴歸公式會列印到標準輸出：

$ go build
$ ./myprogram

Regression Formula:
Predicted = 7.07 + TV*0.05

從結果可以看到，參照的包輸出的線性迴歸的截距為7.07，斜率為0.5。這裡可以進行一些簡單的檢查，因為我們在散點圖中看到了TV和Sales之間的相關性是上升和向右的(即正相關)，這也意味著公式的斜率應該是正數。

評估訓練模型

下面需要通過評估模型的表現來檢視是否可以使用自變數TV來預測Sales。為此，需要載入測試集，使用訓練過的模型對每個測試例進行預測，然後計算第3章"評估和驗證"中討論的某個評估指標。

為此，我們使用平均絕對誤差(Mean Absolute Error (MAE))作為評估指標，這樣就可以直接對比結果和Sales，而不必太擔心異常值或極端值。

為了通過訓練的regression.Regression 值來預測Sales，只需解析測試集的值，並針對regression.Regression 的值呼叫Predict()。然後，計算預測值與觀測值的差值，得到差值的絕對值，然後將所有絕對值相加，得到MAE：

// Open the test dataset file.
f, err = os.Open("test.csv")
if err != nil {
    log.Fatal(err)
}
defer f.Close()
// Create a CSV reader reading from the opened file.
reader = csv.NewReader(f)
// Read in all of the CSV records
reader.FieldsPerRecord = 4
testData, err := reader.ReadAll()
if err != nil {
    log.Fatal(err)
}
// Loop over the test data predicting y and evaluating the prediction
// with the mean absolute error.
var mAE float64
for i, record := range testData {
    // Skip the header.
    if i == 0 {
    continue
    }
    // Parse the observed Sales, or "y".
    yObserved, err := strconv.ParseFloat(record[3], 64)
    if err != nil {
        log.Fatal(err)
    }
    // Parse the TV value.
    tvVal, err := strconv.ParseFloat(record[0], 64)
    if err != nil {
        log.Fatal(err)
    }
    // Predict y with our trained model.
    yPredicted, err := r.Predict([]float64{tvVal})
    // Add the to the mean absolute error.
    mAE += math.Abs(yObserved-yPredicted) / float64(len(testData))
}
// Output the MAE to standard out.
fmt.Printf("MAE = %0.2f\n\n", mAE)

編譯並執行評估程式，得到如下結果：

$ go build
$ ./myprogram

Regression Formula:
Predicted = 7.07 + TV*0.05

MAE = 3.01

那麼如果MAE = 3.01，我們怎麼知道該值是好的還是壞的？這也是為什麼擁有一個良好的資料心智模型很重要的原因。我們已經計算了銷售額的平均值、範圍和標準差。平均銷售額為14.02，標準差為5.21。這樣我們的MAE小於銷售額數值的標準差，約為平均值的20%，說明我們的模型具有一定的預測能力。

恭喜，我們構建了第一個具有預測能力的機器學習模型。

為了更直觀地瞭解模型的執行狀況，可以藉助圖形來幫助視覺化線性迴歸線(利用gonum.org/v1/plot)。首先，建立一個可以執行預測的函數(而需要使用 github.com/sajari/regression包，相當於樁函數)，通過這種方式可以提供輕量級記憶體訓練模型：

// predict uses our trained regression model to made a prediction.
func predict(tv float64) float64 {
    return 7.07 + tv*0.05
}

然後建立視覺化的迴歸線：

// Open the advertising dataset file.
f, err := os.Open("Advertising.csv")
if err != nil {
    log.Fatal(err)
}
defer f.Close()
// Create a dataframe from the CSV file.
advertDF := dataframe.ReadCSV(f)
// Extract the target column.
yVals := advertDF.Col("Sales").Float()
// pts will hold the values for plotting.
pts := make(plotter.XYs, advertDF.Nrow())
// ptsPred will hold the predicted values for plotting.
ptsPred := make(plotter.XYs, advertDF.Nrow())
// Fill pts with data.
for i, floatVal := range advertDF.Col("TV").Float() {
    pts[i].X = floatVal
    pts[i].Y = yVals[i]
    ptsPred[i].X = floatVal
    ptsPred[i].Y = predict(floatVal)
}
// Create the plot.
p, err := plot.New()
if err != nil {
    log.Fatal(err)
}
p.X.Label.Text = "TV"
p.Y.Label.Text = "Sales"
p.Add(plotter.NewGrid())
// Add the scatter plot points for the observations.
s, err := plotter.NewScatter(pts)
if err != nil {
    log.Fatal(err)
}
s.GlyphStyle.Radius = vg.Points(3)
// Add the line plot points for the predictions.
l, err := plotter.NewLine(ptsPred)
if err != nil {
    log.Fatal(err)
}
l.LineStyle.Width = vg.Points(1)
l.LineStyle.Dashes = []vg.Length{vg.Points(5), vg.Points(5)}
// Save the plot to a PNG file.
p.Add(s, l)
if err := p.Save(4*vg.Inch, 4*vg.Inch, "regression_line.png"); err != nil {
    log.Fatal(err)
}

編譯和執行後得到如下圖：

可以看到，我們訓練的迴歸線與實際的資料趨勢相匹配。

多元線性迴歸

線性迴歸並不侷限於依賴單個自變數的簡單線性公式。多元線性迴歸與前面討論的類似，但具有多個自變數(x₁,x₂等)。這種場景下的直線方程如下：

\[y=m_1x_1+m_1x_2+...+m_Nx_N+b \]

這裡x作為自變數，m作為與自變數相關的斜率，此外還有一個截距b。

多元線性迴歸相對比較難以視覺化和思考，因為它不再是一條可以在二維中視覺化的直線。而是一條在二維、三維或多維的線性曲面。但它使用了很多在一元線性迴歸中用過的技術。

多元線性迴歸具有與一元線性迴歸相同的假設，但需要注意的是與之相關的陷阱：

過擬合：通過為模型新增越來越多的自變數，會增加模型的複雜度，並存在過擬合的風險。可以使用之前推薦的技術：正則化(regularization)來解決這種問題。正則化在模型中建立一個懲罰項，該懲罰項是一個與模型複雜性有關的函數，有助於控制這種影響。
相對比例(Relative Scale)：在某些場景下，如果其中某個自變數的比例與另一個自變數的比例相差幾個數量級，那麼較大的變數可能會抵消較小變數帶來的影響，因此可能需要考慮規範化變數。

記住以上兩點，下面嘗試將Sales模型從一元線性迴歸模型擴充套件到多元線性迴歸模型。回顧一下前面章節中的散點圖，可以看到Radio似乎也與Sales呈線性關係，因此可以嘗試建立一個多元線性迴歸模型，如下：

\[Sales=m_1TV+m_2Radio+b \]

使用github.com/sajari/regression時，需要在regression.Regression中標記其他變數，並確保這些值在訓練資料點中成對出現：

// Open the training dataset file.
f, err := os.Open("training.csv")
if err != nil {
    log.Fatal(err)
}
defer f.Close()
// Create a new CSV reader reading from the opened file.
reader := csv.NewReader(f)
// Read in all of the CSV records
reader.FieldsPerRecord = 4
trainingData, err := reader.ReadAll()
if err != nil {
    log.Fatal(err)
}
// In this case we are going to try and model our Sales
// by the TV and Radio features plus an intercept.
var r regression.Regression
r.SetObserved("Sales")
r.SetVar(0, "TV")
r.SetVar(1, "Radio")
// Loop over the CSV records adding the training data.
for i, record := range trainingData {
    // Skip the header.
    if i == 0 {
        continue
    }
    // Parse the Sales.
    yVal, err := strconv.ParseFloat(record[3], 64)
    if err != nil {
        log.Fatal(err)
    }
    // Parse the TV value.
    tvVal, err := strconv.ParseFloat(record[0], 64)
    if err != nil {
        log.Fatal(err)
    }
    // Parse the Radio value.
    radioVal, err := strconv.ParseFloat(record[1], 64)
    if err != nil {
        log.Fatal(err)
    }
    // Add these points to the regression value.
    r.Train(regression.DataPoint(yVal, []float64{tvVal, radioVal}))
}
// Train/fit the regression model.
r.Run()
// Output the trained model parameters.
fmt.Printf("\nRegression Formula:\n%v\n\n", r.Formula)

編譯並執行，得到如下回歸公式：

$ go build
$ ./myprogram

Regression Formula:
Predicted = 2.93 + TV*0.05 + Radio*0.18

可以看到，迴歸公式增加了一個額外的自變數項。截距值也發生了變化。

可以用與一元線性迴歸模型類似的方式，使用Predict()方法來測試該模型：

// Open the test dataset file.
f, err = os.Open("test.csv")
if err != nil {
    log.Fatal(err)
}
defer f.Close()
// Create a CSV reader reading from the opened file.
reader = csv.NewReader(f)
// Read in all of the CSV records
reader.FieldsPerRecord = 4
testData, err := reader.ReadAll()
if err != nil {
    log.Fatal(err)
}
// Loop over the test data predicting y and evaluating the prediction
// with the mean absolute error.
var mAE float64
for i, record := range testData {
    // Skip the header.
    if i == 0 {
        continue
    }
    // Parse the Sales.
    yObserved, err := strconv.ParseFloat(record[3], 64)
    if err != nil {
        log.Fatal(err)
    }
    // Parse the TV value.
    tvVal, err := strconv.ParseFloat(record[0], 64)
    if err != nil {
        log.Fatal(err)
    }
    // Parse the Radio value.
    radioVal, err := strconv.ParseFloat(record[1], 64)
    if err != nil {
        log.Fatal(err)
    }
    // Predict y with our trained model.
    yPredicted, err := r.Predict([]float64{tvVal, radioVal})
    // Add the to the mean absolute error.
    mAE += math.Abs(yObserved-yPredicted) / float64(len(testData))
}
// Output the MAE to standard out.
fmt.Printf("MAE = %0.2f\n\n", mAE)

執行該程式將為我們的新多元迴歸模型得出如下MAE：

$ go build
$ ./myprogram

Regression Formula:
Predicted = 2.93 + TV*0.05 + Radio*0.18
MAE = 1.26

新的多元迴歸模型提高了MAE值。現在，我們可以根據廣告支出來預測Sales了。你還可以嘗試將Newspaper新增到模型。

注意，模型複雜性增加的同時，也會犧牲掉簡易性，並增加過擬合的風險，因此只考慮當新增的複雜效能夠提升模型的表現、並帶來更大的價值時。

非線性以及其他型別的迴歸

雖然本章節主要關注線性迴歸，但不會僅限於使用線性方程來執行迴歸。你可以使用一個或多個非線性(如冪、指數或其他變換)自變數來為因變數建模。例如，我們可以通過一系列TV項來為Sales建模：

\[Sales=m_1TV+m_2TV^2+m_3TV^3+...+b \]

注意，增加複雜性的同時也增加了過擬合的風險。

為了實現非線性迴歸，不能使用github.com/sajari/regression(僅限於線性迴歸)，但可以使用github.com/go-hep/hep/tree/main/fit 來擬合或訓練特定的非線性模型。在Go社群中有很多人已經或正在開發非線性模型工具。

除了OLS外還有其他線性迴歸技術，可以幫助克服最小二乘線性迴歸中的一些假設和弱點。包括嶺迴歸和套索迴歸(lasso regression)。這兩種技術使用懲罰迴歸係數來減輕自變數的多重共線性和非正態性帶來的影響。

github.com/berkmancenter/ridge中實現了Go語言的嶺迴歸。與 github.com/sajari/regression不同，我們的自變數和因變數資料是通過gonum矩陣輸入github.com/berkmancenter/ridge的。為了說明該方法，我們首先構造一個包含廣告支出 (TV, Radio, 和Newspaper)的矩陣，以及包含Sales資料的矩陣。注意在github.com/berkmancenter/ridge中，如果想在模型中有一個截距，則需要為截距的輸入自變數矩陣顯式地新增一列，該列中的每個值僅為1.0。

// Open the training dataset file.
f, err := os.Open("training.csv")
if err != nil {
    log.Fatal(err)
}
defer f.Close()
// Create a new CSV reader reading from the opened file.
reader := csv.NewReader(f)
reader.FieldsPerRecord = 4
// Read in all of the CSV records
rawCSVData, err := reader.ReadAll()
if err != nil {
    log.Fatal(err)
}
// featureData will hold all the float values that will eventually be
// used to form our matrix of features.
featureData := make([]float64, 4*len(rawCSVData))
yData := make([]float64, len(rawCSVData))
// featureIndex and yIndex will track the current index of the matrix values.
var featureIndex int
var yIndex int
// Sequentially move the rows into a slice of floats.
for idx, record := range rawCSVData {
    // Skip the header row.
    if idx == 0 {
        continue
    }
    // Loop over the float columns.
    for i, val := range record {
        // Convert the value to a float.
        valParsed, err := strconv.ParseFloat(val, 64)
        if err != nil {
            log.Fatal(err)
        }
        if i < 3 {
        // Add an intercept to the model.
        if i == 0 {
            featureData[featureIndex] = 1
            featureIndex++
        }
        // Add the float value to the slice of feature floats.
        featureData[featureIndex] = valParsed
        featureIndex++
        }
        if i == 3 {
            // Add the float value to the slice of y floats.
            yData[yIndex] = valParsed
            yIndex++
        }
    }
}
// Form the matrices that will be input to our regression.
features := mat64.NewDense(len(rawCSVData), 4, featureData)
y := mat64.NewVector(len(rawCSVData), yData)

下面使用自變數和因變數矩陣建立一個新的ridge.RidgeRegression值，然後呼叫Regress() 方法來訓練模型，最後列印訓練的迴歸公式：

// Create a new RidgeRegression value, where 1.0 is the
// penalty value.
r := ridge.New(features, y, 1.0)
// Train our regression model.
r.Regress()
// Print our regression formula.
c1 := r.Coefficients.At(0, 0)
c2 := r.Coefficients.At(1, 0)
c3 := r.Coefficients.At(2, 0)
c4 := r.Coefficients.At(3, 0)
fmt.Printf("\nRegression formula:\n")
fmt.Printf("y = %0.3f + %0.3f TV + %0.3f Radio + %0.3f Newspaper\n\n", c1, c2, c3, c4)

編譯並執行程式，得出如下回歸公式：

$ go build
$ ./myprogram
Regression formula:
y = 3.038 + 0.047 TV + 0.177 Radio + 0.001 Newspaper

可以看到TV和Radio的係數與最小二乘迴歸得到的結果類似，但略微不同。另外可以看到新增了一個Newspaper項。

可以通過建立自己的預測函數來測試嶺迴歸公式：

// predict uses our trained regression model to made a prediction based on a
// TV, Radio, and Newspaper value.
func predict(tv, radio, newspaper float64) float64 {
    return 3.038 + tv*0.047 + 0.177*radio + 0.001*newspaper
}

然後使用該predict函數和測試資料來測試嶺迴歸公式：

// Open the test dataset file.
f, err := os.Open("test.csv")
if err != nil {
    log.Fatal(err)
}
defer f.Close()
// Create a new CSV reader reading from the opened file.
reader := csv.NewReader(f)
// Read in all of the CSV records
reader.FieldsPerRecord = 4
testData, err := reader.ReadAll()
if err != nil {
    log.Fatal(err)
}
// Loop over the holdout data predicting y and evaluating the prediction
// with the mean absolute error.
var mAE float64
for i, record := range testData {
    // Skip the header.
    if i == 0 {
        continue
    }
    // Parse the Sales.
    yObserved, err := strconv.ParseFloat(record[3], 64)
    if err != nil {
        log.Fatal(err)
    }
    // Parse the TV value.
    tvVal, err := strconv.ParseFloat(record[0], 64)
    if err != nil {
        log.Fatal(err)
    }
    // Parse the Radio value.
    radioVal, err := strconv.ParseFloat(record[1], 64)
    if err != nil {
        log.Fatal(err)
    }
    // Parse the Newspaper value.
    newspaperVal, err := strconv.ParseFloat(record[2], 64)
    if err != nil {
        log.Fatal(err)
    }
    // Predict y with our trained model.
    yPredicted := predict(tvVal, radioVal, newspaperVal)
    // Add the to the mean absolute error.
    mAE += math.Abs(yObserved-yPredicted) / float64(len(testData))
}
// Output the MAE to standard out.
fmt.Printf("\nMAE = %0.2f\n\n", mAE)

編譯並執行，等到新的MAE：

$ go build
$ ./myprogram

MAE = 1.26

注意在模型中新增Newspaper並不會提高MAE，因此在這種場景下並不適合新增Newspaper項，因為此時提高了複雜度，但並沒有顯著影響到模型。

新增到模型中的任何複雜性或複雜度都應該有其可衡量的理由。使用一個複雜的模型，這看起來很有趣，但同時也會讓人頭疼。