一：當訓練集為1維時

#進行資料分析所需庫，可以看做是對numpy工具的補充
import pandas as pd
import numpy as np

#應該把Seaborn視為matplotlib的補充,作圖所用工具，在大多數情況下使用seaborn就能做出很具有吸引力的圖，而使用matplotlib就能製作具有更多特色的圖
import seaborn as sns
import matplotlib.pyplot as plt
#設定繪畫的圖示格式和顏色
sns.set(context="notebook", style="whitegrid", palette="dark
 
")
#讀取資料並賦予列名，此時df共有兩列，列名分別為population和profit
df = pd.read_csv('ex1data1.txt', names=['population', 'profit'])


def get_x(df):
    #ones:格式化為一個m列的全為1的向量
    ones = pd.DataFrame({'ones':np.ones(len(df))})
    #將date格式化為ones向量右邊加上df矩陣的矩陣
    data = pd.concat([ones, df], axis=1)
    #返回data矩陣的前兩列，返回m * 2的矩陣
 

    return data.iloc[:, :-1].as_matrix() # 這個操作返回 ndarray,不是矩陣


def get_y(df):
    #返回df矩陣的最後一列，返回m × 1的向量
    return np.array(df.iloc[:, -1])


def normalize_feature(df):
    # 特徵縮放,對df中的兩列資料分別進行特徵縮放 ：（x - x平均值）/ x方差
    return df.apply(lambda column: (column - column.mean())/column.std())


X = get_x(df)
 
print(X.shape, type(X)) #X為m × 2的矩陣
y = get_y(df)
print(y.shape, type(y)) #y為m × 1的向量

theta = np.zeros(X.shape[1])#X.shape[1]=2,代表特徵數n X.shape(0) = 97;X.shape(1) = 2;將theta初始化為2行的零向量（列向量）


def lr_cost(theta, X, y):
#     """
#     計算theta固定時，此時的代價函式值
#     X: R(m*n), m 樣本數, n 特徵數
#     y: R(m)
#     theta : R(n), 線性迴歸的引數
#     """
    m = X.shape[0]#m為樣本數

    #X矩陣和theta矩陣的點積所得矩陣（m × 1） - 矩陣y（m × 1）
    inner = X @ theta - y  # R(m*1)，X @ theta等價於X.dot(theta)

    #square_sum為inner的每個元素平方之和，即二範數
    square_sum = inner.T @ inner

    #得到此時代價函式的值
    cost = square_sum / (2 * m)

    return cost


#將代價函式對theta求導，即求梯度
def gradient(theta, X, y):
    #m為樣本數
    m = X.shape[0]

    #這裡的inner為代價函式對theta求導的結果
    inner = X.T @ (X @ theta - y)  # (m,n).T @ (m, 1) -> (n, 1)，X @ theta等價於X.dot(theta)

    return inner / m

#梯度下降函式
def batch_gradient_decent(theta, X, y, epoch, alpha=0.01):
#   擬合線性迴歸，返回引數和代價
#     epoch: 批處理的輪數
#     alpha: theta移動的步長
#     """

    #得到初始theta時的代價
    cost_data = [lr_cost(theta, X, y)]

    # 拷貝一份，不和原來的theta混淆
    _theta = theta.copy()

    #開始迭代epoch次
    for _ in range(epoch):
        #根據當前梯度更新theta
        _theta = _theta - alpha * gradient(_theta, X, y)
        #記錄此時的代價
        cost_data.append(lr_cost(_theta, X, y))

    return _theta, cost_data


#迭代500次求最小代價和所對應的theta
epoch = 500
final_theta, cost_data = batch_gradient_decent(theta, X, y, epoch)

# 計算最終的代價
lr_cost(final_theta, X, y)

#畫出代價函式值變化圖
#可以看到從第二輪代價資料變換很大，接下來平穩了
ax = sns.tsplot(cost_data, time=np.arange(epoch+1))
ax.set_xlabel('epoch')
ax.set_ylabel('cost')
plt.show()


b = final_theta[0] # intercept，Y軸上的截距
m = final_theta[1] # slope，斜率

#畫出原資料點和線性迴歸的最終結果
plt.scatter(df.population, df.profit, label="Training data")
plt.plot(df.population, df.population*m + b, label="Prediction")
plt.legend(loc=2)
plt.show()

二：當訓練集為多維時

當訓練集為多維時只需將上述程式碼稍作改動，完整程式碼如下：

#進行資料分析所需庫，可以看做是對numpy工具的補充
import pandas as pd
import numpy as np

#應該把Seaborn視為matplotlib的補充,作圖所用工具，在大多數情況下使用seaborn就能做出很具有吸引力的圖，而使用matplotlib就能製作具有更多特色的圖
import seaborn as sns
import matplotlib.pyplot as plt

#設定繪畫的圖示格式和顏色
sns.set(context="notebook", style="whitegrid", palette="dark")
#讀取資料並賦予列名，此時df共有兩列，列名分別為population和profit
df = pd.read_csv('ex1data2.txt', names=['square', 'bedrooms', 'price'])
df.head()



def get_x(df):
    #ones:格式化為一個m列的全為1的向量
    ones = pd.DataFrame({'ones':np.ones(len(df))})
    #將date格式化為ones向量右邊加上df矩陣的矩陣
    data = pd.concat([ones, df], axis=1)
    #返回data矩陣的前兩列，返回m * 2的矩陣
    return data.iloc[:, :-1].as_matrix() # 這個操作返回 ndarray,不是矩陣


def get_y(df):
    #返回df矩陣的最後一列，返回m × 1的向量
    return np.array(df.iloc[:, -1])


def normalize_feature(df):
    # 特徵縮放,對df中的兩列資料分別進行特徵縮放 ：（x - x平均值）/ x方差
    return df.apply(lambda column: (column - column.mean())/column.std())

data = normalize_feature(df) #特徵縮放

X = get_x(data)
print(X.shape, type(X)) #X為m × 2的矩陣
y = get_y(data)
print(y.shape, type(y)) #y為m × 1的向量

theta = np.zeros(X.shape[1])#X.shape[1]=2,代表特徵數n X.shape(0) = 97;X.shape(1) = 2;將theta初始化為2行的零向量（列向量）


def lr_cost(theta, X, y):
#     
#     計算theta固定時，此時的代價函式值
#     X: R(m*n), m 樣本數, n 特徵數
#     y: R(m)
#     theta : R(n), 線性迴歸的引數
#     
    m = X.shape[0]#m為樣本數

    #X矩陣和theta矩陣的點積所得矩陣（m × 1） - 矩陣y（m × 1）
    inner = X @ theta - y  # R(m*1)，X @ theta等價於X.dot(theta)

    #square_sum為inner的每個元素平方之和，即二範數
    square_sum = inner.T @ inner

    #得到此時代價函式的值
    cost = square_sum / (2 * m)

    return cost


#將代價函式對theta求導，即求梯度
def gradient(theta, X, y):
    #m為樣本數
    m = X.shape[0]

    #這裡的inner為代價函式對theta求導的結果
    inner = X.T @ (X @ theta - y)  # (m,n).T @ (m, 1) -> (n, 1)，X @ theta等價於X.dot(theta)

    return inner / m


#梯度下降函式
def batch_gradient_decent(theta, X, y, epoch, alpha=0.01):
#   擬合線性迴歸，返回引數和代價
#     epoch: 批處理的輪數
#     alpha: theta移動的步長
#     

    #得到初始theta時的代價
    cost_data = [lr_cost(theta, X, y)]

    # 拷貝一份，不和原來的theta混淆
    _theta = theta.copy()

    #開始迭代epoch次
    for _ in range(epoch):
        #根據當前梯度更新theta
        _theta = _theta - alpha * gradient(_theta, X, y)
        #記錄此時的代價
        cost_data.append(lr_cost(_theta, X, y))

    return _theta, cost_data


base = np.logspace(-1, -5, num=4)
candidate = np.sort(np.concatenate((base, base*3)))

epoch = 50

fig, ax = plt.subplots(figsize=(16, 9))

for alpha in candidate:
    _, cost_data = batch_gradient_decent(theta, X, y, epoch, alpha=alpha)
    ax.plot(np.arange(epoch+1), cost_data, label=alpha)

ax.set_xlabel('epoch', fontsize=18)
ax.set_ylabel('cost', fontsize=18)
ax.legend(bbox_to_anchor=(1.05, 1), loc=2, borderaxespad=0.)
ax.set_title('learning rate', fontsize=18)
plt.show()
'''
#畫出代價函式值變化圖
#可以看到從第二輪代價資料變換很大，接下來平穩了
ax = sns.tsplot(cost_data, time=np.arange(epoch+1))
ax.set_xlabel('epoch')
ax.set_ylabel('cost')
plt.show()
'''