<jupyter_start><jupyter_text>最小二乘法与线性回归 $预测值 = \sum (回归值 * 输入项)$即 $Y_i = X^Tw$使用最小二乘法计算w为找到$\min {\sum_{i=1}^m (y_i - x_i^T w)^2}$我们对其求导，其矩阵表达式为：$$\begin{aligned}s(w) & =(Y - Xw)^T(Y - Xw) \\& =(Y^T - w^TX^T)(Y - Xw) \\& =(Y_TY - Y_TXw - b^TX^TY - w^TX^TXw) \\& =\underbrace{Y^T Y}_{①} - \underbrace{2w^TX^TY}_{②} + \underbrace{b^TX^TXb}_{③}\end{aligned}$$对w进行矩阵求导,并令其等于零：$$\begin{aligned}\frac{\partial {s(w)}} {\partial {w}} = 0\end{aligned}$$①:$\frac{\partial {Y^TY}} {\partial {w}} = 0$②:$\frac{\partial {2b^TX^TY}} {\partial {w}} = 2\frac{\partial {b^T(X^TY)}} {\partial {w}} = 2X^TY$③:$\frac{\partial {b^TX^TXb}} {\partial {w}} = \frac{\partial {b^T(X^TX)b}} {\partial {w}}= [X^TX + (X^TX)^T]w = 2X^TXw$$\therefore\frac{\partial {s(w)}} {\partial {w}} = 0 - 2X^TY + 2X^TXw = 0 \\X^TXw = X^TY \\w = (X^TX)^{-1}X^TY$即$\hat{w} = (X^TX)^{-1}X^TY$  样例数据导入<jupyter_code>import numpy as np

def load_data(filename):
    data = []
    label = []
    fr = open(filename)
    for line in fr.readlines():
        line_arr = []
        curline = line.strip().split('\t')
        for i in range(len(curline) - 1):
            line_arr.append(float(curline[i]))
        data.append(line_arr)
        label.append(float(curline[-1]))
    return data, label
dataset, label = load_data(r'ex0.txt')
print(dataset[0:2])
def stand_regress(x_arr, y_arr):
    x = np.mat(x_arr)
    y = np.mat(y_arr).T
    xTx = x.T * x
    if np.linalg.det(xTx) == 0.0:
        print("this matrix cannot inverse")
        return
    ws = xTx.I * (x.T*y)
    return ws
ws = stand_regress(dataset, label)
ws
import matplotlib.pyplot as plt
x = np.mat(dataset)
y = np.mat(label)
y_hat = x*ws
fig = plt.figure()
ax = fig.add_subplot(111)
ax.scatter(x[:, 1].flatten().A[0], y.T[:, 0].flatten().A[0], s=2)
ax.plot(x[:, 1], y_hat[:, 0], c='red')<jupyter_output><empty_output><jupyter_text>结果评估利用相关系数来评估预测值和真实值之间的相关性<jupyter_code>np.corrcoef(y_hat.T, y)<jupyter_output><empty_output><jupyter_text>局部加权线性回归线性回归——欠拟合解决方法：引入偏差$\hat{w} = (X^TWX)^-1X^TWy$其中，w是一个矩阵，w的值由高斯核计算。高斯核的权重为：$$\begin{aligned}w(i,i) = exp(\frac{|x_i - x|}{-2k^2})\end{aligned}$$<jupyter_code>def lwlr(testpoint, x_arr, y_arr, k):
    x_mat = np.mat(x_arr)
    y_mat = np.mat(y_arr).T
    m = np.shape(x_mat)[0]
    weight = np.mat(np.eye((m)))
    for i in range(m):
        diff = testpoint - x_mat[i, :]
        weight[i,i] = np.exp(diff * diff.T / (-2.0 * k ** 2))
    xTx = x_mat.T * (weight*x_mat)
    if np.linalg.det(xTx) == 0.0:
        print('this matrix cannot inverse')
        return
    ws = xTx.I * (x_mat.T *(weight* y_mat))
    return testpoint * ws
def lwlrtest(test_arr, x_arr, y_arr, k=1.0):
    m = np.shape(np.mat(test_arr))[0]
    yhat = np.zeros(m)
    for i in range(m):
        yhat[i] = lwlr(test_arr[i], x_arr, y_arr, k)
    return yhat
x_arr, y_arr = load_data('ex0.txt')
print(y_arr[0])
lwlr(x_arr[0], x_arr, y_arr, 1.0)
lwlr(x_arr[0], x_arr, y_arr, 0.001)
yhat1 = lwlrtest(x_arr, x_arr, y_arr, 1.0)
yhat2 = lwlrtest(x_arr, x_arr, y_arr, 0.01)
yhat3 = lwlrtest(x_arr, x_arr, y_arr, 0.003)
x_mat = np.mat(x_arr)
str_index = x_mat[:,1].argsort(0)
# print(x_mat)
x_sort = x_mat[str_index][:,0,:]
print(x_sort)<jupyter_output>[[1.       0.014855]
 [1.       0.015371]
 [1.       0.033859]
 [1.       0.038326]
 [1.       0.040486]
 [1.       0.045353]
 [1.       0.047744]
 [1.       0.058014]
 [1.       0.066379]
 [1.       0.067119]
 [1.       0.067154]
 [1.       0.067732]
 [1.       0.068854]
 [1.       0.070095]
 [1.       0.070237]
 [1.       0.073224]
 [1.       0.085031]
 [1.       0.106352]
 [1.       0.110521]
 [1.       0.115487]
 [1.       0.116089]
 [1.       0.116163]
 [1.       0.116598]
 [1.       0.124931]
 [1.       0.126762]
 [1.       0.132791]
 [1.       0.133295]
 [1.       0.138306]
 [1.       0.147846]
 [1.       0.148364]
 [1.       0.15335 ]
 [1.       0.156383]
 [1.       0.156433]
 [1.       0.162176]
 [1.       0.166819]
 [1.       0.171846]
 [1.       0.175724]
 [1.       0.178449]
 [1.       0.181951]
 [1.       0.182026]
 [1.       0.189912]
 [1.       0.19683 ]
 [1.       0.199526]
 [1.       0.201473]
 [1.       0.204337]
 [1.       0.212575]
 [1.       0.218023]
 [1.       0.[...]<jupyter_text>k=1.0的回归结果<jupyter_code>fig = plt.figure()
ax = fig.add_subplot(111)
ax.plot(x_sort[:, 1], yhat1[str_index])
ax.scatter(x_mat[:, 1].flatten().A[0], np.mat(y_arr).T.flatten().A[0], s=2 ,c='red')<jupyter_output><empty_output><jupyter_text>k=0.01的回归结果<jupyter_code>fig = plt.figure()
ax = fig.add_subplot(111)
ax.plot(x_sort[:, 1], yhat2[str_index])
ax.scatter(x_mat[:, 1].flatten().A[0], np.mat(y_arr).T.flatten().A[0], s=2 ,c='red')<jupyter_output><empty_output><jupyter_text>k=0.003的回归结果<jupyter_code>fig = plt.figure()
ax = fig.add_subplot(111)
ax.plot(x_sort[:, 1], yhat3[str_index])
ax.scatter(x_mat[:, 1].flatten().A[0], np.mat(y_arr).T.flatten().A[0], s=2 ,c='red')<jupyter_output><empty_output><jupyter_text>预测鲍鱼年龄<jupyter_code>def ress_error(y_arr, yhat_arr):
    return ((y_arr - yhat_arr)**2).sum()
abx,aby = load_data('abalone.txt')
yhat_ab1 = lwlrtest(abx[0:99], abx[0:99], aby[0:99], 0.1)
yhat_ab2 = lwlrtest(abx[0:99], abx[0:99], aby[0:99], 1.0)
yhat_ab3 = lwlrtest(abx[0:99], abx[0:99], aby[0:99], 10.0)
ress_error(aby[0:99], yhat_ab1.T)
ress_error(aby[0:99], yhat_ab2.T)
ress_error(aby[0:99], yhat_ab3.T)
yhat_test1 = lwlrtest(abx[100:199], abx[0:99], aby[0:99], 0.1)
ress_error(aby[100:199], yhat_test1.T)
yhat_test2 = lwlrtest(abx[100:199], abx[0:99], aby[0:99], 1)
ress_error(aby[100:199], yhat_test2.T)
yhat_test3 = lwlrtest(abx[100:199], abx[0:99], aby[0:99], 10)
ress_error(aby[100:199], yhat_test3.T)<jupyter_output><empty_output><jupyter_text>k=10 的测试误差最小 但训练误差最大<jupyter_code>ws = stand_regress(abx[0:99], aby[0:99])
yhat_regress = np.mat(abx[100:199])*ws
ress_error(aby[100:199], yhat_regress.T.A)<jupyter_output><empty_output><jupyter_text>岭回归 如果矩阵的特征值比样本点多（$n>m$），那么$X$非满秩，$X^TX$为奇异矩阵无法求逆。岭回归就是在$X^TX$的基础上加上$\lambda I$使得矩阵非奇异，进而求逆。此时，回归系数的计算公式为：$$\begin{aligned}\hat{w} = (X^TX + \lambda I)^{-1}X^Ty\end{aligned}$$<jupyter_code>def ridge_regress(x_mat, y_mat, lam=0.2):
    xTx = x_mat.T * x_mat
    demon = xTx + np.eye(np.shape(x_mat)[1]) * lam
    if np.linalg.det(demon) == 0.0:
        print("this matrix cannot inverse")
        return 
    ws = demon.I * x_mat.T * y_mat
    return ws
def ridge_test(x_arr, y_arr):
    x_mat = np.mat(x_arr)
    y_mat = np.mat(y_arr).T
    y_mean = np.mean(y_mat, 0)
    y_mat = y_mat - y_mean
    x_mean = np.mean(x_mat, 0)
    x_var = np.var(x_mat, 0)
    x_mat = (x_mat - x_mean) / x_var
    num_test = 30
    w_mat = np.zeros((num_test, np.shape(x_mat)[1]))
    for i in range(num_test):
        ws = ridge_regress(x_mat, y_mat, np.exp(i - 10))
        w_mat[i, :] = ws.T
    return w_mat
ridge_weight = ridge_test(abx, aby)
fig = plt.figure()
ax = fig.add_subplot(111)
ax.plot(ridge_weight)<jupyter_output><empty_output><jupyter_text>前向逐步回归基于贪心，所有权重设为1，之后逐渐增加一个很小的值，不断迭代来选择误差最小的w<jupyter_code>def regularize(in_mat):
    out_mat = in_mat.copy()
    out_mean = np.mean(out_mat, 0)
    out_var = np.var(out_mat, 0)
    out_mat = (out_mat - out_mean) / out_var
    return out_mat
def stage_wise(x_arr, y_arr, eps=0.01,num=100):
    x_mat = np.mat(x_arr)
    y_mat = np.mat(y_arr).T
    y_mean = np.mean(y_mat, 0)
    y_mat = y_mat - y_mean
    x_mat = regularize(x_mat)
    m, n = np.shape(x_mat)
    return_mat = np.zeros((num, n))
    ws = np.zeros((n,1))
    ws_test = ws.copy()
    ws_max = ws.copy()
    # count = 0
    for i in range(num):
        print(ws.T)
        lowest_err = np.inf
        for j in range(n):
            for sign in [-1,1]:
                # print(sign)
                ws_test = ws.copy()
                ws_test[j] += eps*sign
                y_test = x_mat*ws_test
                ress_err = ress_error(y_mat.A, y_test.A)
                if ress_err < lowest_err:
                    lowest_err = ress_err
                    ws_max = ws_test
                # count += 1
        ws = ws_max.copy()
        return_mat[i, :] = ws.T
    return return_mat
stage_wise(abx, aby, 0.01, 200)
return_mat = stage_wise(abx, aby, 0.001, 5000)
fig = plt.figure()
ax = fig.add_subplot(111)
ax.plot(return_mat)
plt.show()
x_mat = np.mat(abx)
y_mat = np.mat(aby).T
x_mat = regularize(x_mat)
y_mean = np.mean(y_mat, 0)
y_mat = y_mat - y_mean
stand_regress(x_mat, y_mat.T).T<jupyter_output><empty_output><jupyter_text>示例：预测乐高玩具的价格  Step1：收集数据——使用Google购物的API<jupyter_code>from time import sleep
import json
import urllib.request
def searchForSet(retX, retY, setNum, yr, numPce, origPrc):
    sleep(10)
    myAPIstr = '<KEY>'
    searchURL = 'https://www.googleapis.com/shopping/search/v1/public/products?key=%s&country=US&q=lego+%d&alt=json' % (myAPIstr, setNum)
    pg = urllib.request.urlopen(searchURL)
    retDict = json.loads(pg.read())
    for i in range(len(retDict['items'])):
        try:
            currItem = retDict['items'][i]
            if currItem['product']['condition'] == 'new':
                newFlag = 1
            else: newFlag = 0
            listOfInv = currItem['product']['inventories']
            for item in listOfInv:
                sellingPrice = item['price']
                if  sellingPrice > origPrc * 0.5:
                    print("%d\t%d\t%d\t%f\t%f" % (yr,numPce,newFlag,origPrc, sellingPrice))
                    retX.append([yr, numPce, newFlag, origPrc])
                    retY.append(sellingPrice)
        except: print('problem with item %d' % i)
    
def setDataCollect(retX, retY):
    searchForSet(retX, retY, 8288, 2006, 800, 49.99)
    searchForSet(retX, retY, 10030, 2002, 3096, 269.99)
    searchForSet(retX, retY, 10179, 2007, 5195, 499.99)
    searchForSet(retX, retY, 10181, 2007, 3428, 199.99)
    searchForSet(retX, retY, 10189, 2008, 5922, 299.99)
    searchForSet(retX, retY, 10196, 2009, 3263, 249.99)
    
lgx = []
lgy = []
setDataCollect(lgx, lgy)<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.18775181305398872, 'fraction_numerical': 0.07846494762288476, 'mean_word_length': 3.0329000812347684, 'pattern_counts': {'":': 0, '<': 46, '<?xml version=': 0, '>': 47, 'https://': 1, 'lorem ipsum': 0, 'www.': 1, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '11', 'dir': 'jupyter-structured-clean-dedup', 'id': '275608', 'license': 'MIT', 'n_tokens_mistral': 5150, 'n_tokens_neox': 4779, 'n_words': 739, 'path': 'linearRegression/linearRegress.ipynb', 'repo_name': 'Dream-Fu/ML_In_Action', 'size': '9962'}