一個簡單的線性擬合問題，到底有多少種做法

相信大家都做過線性擬合問題吧，其實就是給很多點，來求線性方程的斜率和截距。早在高中數學就有這類問題，我記得很清楚，如果出現在試卷中，一般出現在解答題的第二題左右，高中中的做法就是最小二乘法，代入公式，求斜率和截距，說句好聽，就是送分題。

在科學計算中，也是采用ols(普通最小二乘法)進行回歸分析。OLS 全稱ordinary least squares,是回歸分析(regression analysis)最根本的一個形式。

深度學習

這個就是三大深度學習框架的入門demo。

（1） Keras

Keras 的核心數據結構是 model，一種組織網絡層的方式。最簡單的模型是 Sequential 順序模型，它由多個網絡層線性堆疊。

import keras import numpy as np import matplotlib.pyplot as plt %matplotlib inline from pylab import mpl #正常顯示畫圖時出現的中文 mpl.rcParams['font.sans-serif']=['SimHei'] from keras import layersx = np.linspace(0,50,50) y = 3 * x + 7 + np.random.randn(50) *5 plt.scatter(x,y)

通過numpy隨機繪制了50個散點圖，如下圖所示

下面使用keras搭建最簡單的順序Sequential模型

# 順序模型 model = keras.Sequential() # 輸入輸出是一維特征 model.add(layers.Dense(1,input_dim=1)) model.summary() # 編譯模型 # 使用adam優化器，損失函數mse均方差，在compile加入精確度 model.compile(optimizer='adam',loss='mse',metrics=['acc']) # 訓練模型 epochs 參數是把數據訓練3000遍 history = model.fit(x, y, epochs=3000) plt.scatter(x, y, c='r') plt.plot(x, model.predict(x))

運行代碼，結果如下圖所示

（2）Pytorch

我們創建一個由方程$y=2x+0.2$產生的數據集，并通過torch.rand()函數制造噪音

import torch as t from matplotlib import pyplot as plt from torch.autograd import Variable from torch import nn # 創建數據集 x = Variable(t.unsqueeze(t.linspace(-1, 1, 100), dim=1)) y = Variable(x * 2 + 0.2 + t.rand(x.size())) plt.scatter(x.data.numpy(),y.data.numpy()) plt.show()

制造出來的數據集，如下圖所示

下面我們就要開始定義我們的模型，這里定義的是一個輸入層和輸出層都只有一維的模型，并定義出損失函數MSE和優化函數sgd，這里使用均方誤差作為損失函數。

class LinearRegression(t.nn.Module):def __init__(self):#繼承父類構造函數super(LinearRegression, self).__init__() #輸入和輸出的維度都是1self.linear = t.nn.Linear(1, 1) def forward(self, x):out = self.linear(x)return outmodel = LinearRegression()#實例化對象 num_epochs = 1000#迭代次數 learning_rate = 1e-2#學習率0.01 Loss = t.nn.MSELoss()#損失函數 optimizer = t.optim.SGD(model.parameters(), lr=learning_rate)#優化函數

遍歷每次epoch，計算出loss，反向傳播計算梯度，不斷的更新梯度，使用梯度下降進行優化。

for epoch in range(num_epochs):# 預測y_pred= model(x)# 計算lossloss = Loss(y_pred, y)#清空上一步參數值optimizer.zero_grad()#反向傳播loss.backward()#更新參數optimizer.step()if epoch % 200 == 0:print("[{}/{}] loss:{:.4f}".format(epoch+1, num_epochs, loss))plt.scatter(x.data.numpy(), y.data.numpy()) plt.plot(x.data.numpy(), y_pred.data.numpy(), 'r-',lw=5) plt.text(0.5, 0,'Loss=%.4f' % loss.data.item(), fontdict={'size': 20, 'color': 'red'}) plt.show() ####結果如下#### [1/1000] loss:1.7766 [201/1000] loss:0.1699 [401/1000] loss:0.0816 [601/1000] loss:0.0759 [801/1000] loss:0.0755

運行代碼，結果如下圖所示

（3）tensorflow

下面，我們采用tensorflow進行線性擬合

首先生成我們的數據，數據標簽通過真實函數加上高斯噪聲得到。

然后為了進行梯度計算,定義了變量w和偏移量b，初值都設為0.

import tensorflow as tf import numpy as np# 一些參數 learning_rate = 0.01 # 學習率 training_steps = 1000 # 訓練次數 display_step = 50 # 訓練50次輸出一次# 訓練數據 X = np.linspace(0, 1, 50).reshape((-1,1)) Y = 2* X+ 1 + np.random.normal(0,0.1, X.shape)n_samples = 50# 隨機初始化權重和偏置 W = tf.Variable(np.random.randn(), name="weight") b = tf.Variable(np.random.randn(), name="bias")# 線性回歸函數 def linear_regression(x):return W*x + b# 損失函數 def mean_square(y_pred, y_true):return tf.reduce_sum(tf.pow(y_pred-y_true, 2)) / (2 * n_samples)# 優化器采用隨機梯度下降(SGD) optimizer = tf.optimizers.SGD(learning_rate)# 計算梯度，更新參數 def run_optimization():# tf.GradientTape()梯度帶，可以查看每一次epoch的參數值with tf.GradientTape() as g:pred = linear_regression(X)loss = mean_square(pred, Y)# 計算梯度gradients = g.gradient(loss, [W, b])# 更新W，boptimizer.apply_gradients(zip(gradients, [W, b]))# 開始訓練 for step in range(1, training_steps+1):run_optimization()if step % display_step == 0:pred = linear_regression(X)loss = mean_square(pred, Y)print("step: %i, loss: %f, W: %f, b: %f" % (step, loss, W.numpy(), b.numpy()))

import matplotlib.pyplot as plt plt.plot(X, Y, 'ro', label='Original data') plt.plot(X, np.array(W * X + b), label='Fitted line') plt.legend() plt.show()

上面就是深度學習的做法，你們會看到一列蒙蔽，下面我來演示機器學習做法

機器學習

import numpy as np import matplotlib.pyplot as plt %matplotlib inline x = np.linspace(0,30,50) y = x+ 2*np.random.rand(50) plt.figure(figsize=(10,8)) from sklearn.linear_model import LinearRegression #導入線性回歸 model = LinearRegression() #初始化模型 x1 = x.reshape(-1,1) # 將行變列 得到x坐標 y1 = y.reshape(-1,1) # 將行變列 得到y坐標 model.fit(x1,y1) #訓練數據 plt.scatter(x,y) x_test = np.linspace(0,40).reshape(-1,1) plt.plot(x_test,model.predict(x_test)) model.coef_ #array([[1.00116024]]) 斜率 model.intercept_ # array([0.86175551]) 截距

Numpy

import numpy as np x = np.linspace(0,30,50) y = x+ 1 + np.random.normal(0,0.1, 50) z1 = np.polyfit(x,y,1) #一次多項式擬合，相當于線性擬合 z1 # [1.00895356, 0.71872268] p1 = np.poly1d(z1) p1 # array([1.00032794, 0.9799152 ])

Scipy

from scipy.stats import linregress x = np.linspace(0,30,50) y = x + 2 +np.random.normal(0) slope, intercept, r_value, p_value, std_err = linregress(x, y) print("slope: %f intercept: %f" % (slope, intercept)) # slope: 1.000000 intercept: 1.692957 print("R-squared: %f" % r_value**2) #R-squared: 1.000000 plt.figure(figsize=(10,8)) plt.plot(x, y, 'o', label='original data') plt.plot(x, intercept + slope*x, 'r', label='fitted line') plt.legend() plt.show()

Statsmodels

# 線性模型 import statsmodels.api as sm import numpy as np x = np.linspace(0,10,100) y = 3*x + np.random.randn()+ 10 # Fit and summarize OLS model X = sm.add_constant(x) mod = sm.OLS(y,X) result = mod.fit() print(result.params) print(result.summary()) [9.65615842 3. ]OLS Regression Results ============================================================================== Dep. Variable: y R-squared: 1.000 Model: OLS Adj. R-squared: 1.000 Method: Least Squares F-statistic: 7.546e+31 Date: Thu, 25 Jul 2019 Prob (F-statistic): 0.00 Time: 21:10:18 Log-Likelihood: 3082.0 No. Observations: 100 AIC: -6160. Df Residuals: 98 BIC: -6155. Df Model: 1 Covariance Type: nonrobust ==============================================================================coef std err t P>|t| [0.025 0.975] ------------------------------------------------------------------------------ const 9.6562 2e-15 4.83e+15 0.000 9.656 9.656 x1 3.0000 3.45e-16 8.69e+15 0.000 3.000 3.000 ============================================================================== Omnibus: 4.067 Durbin-Watson: 0.161 Prob(Omnibus): 0.131 Jarque-Bera (JB): 4.001 Skew: 0.446 Prob(JB): 0.135 Kurtosis: 2.593 Cond. No. 11.7 ==============================================================================

R

R語言通過lm()函數創建預測變量與響應變量之間的關系線性回歸模型。下面是我的筆記關于體重和身高的回歸

> x <- c(151, 174, 138, 186, 128, 136, 179, 163, 152, 131) > y <- c(63, 81, 56, 91, 47, 57, 76, 72, 62, 48) > lm(y~x) Call: lm(formula = y ~ x)Coefficients: (Intercept) x -38.4551 0.6746 > summary(lm(y~x))Call: lm(formula = y ~ x)Residuals:Min 1Q Median 3Q Max -6.3002 -1.6629 0.0412 1.8944 3.9775 Coefficients:Estimate Std. Error t value Pr(>|t|) (Intercept) -38.45509 8.04901 -4.778 0.00139 ** x 0.67461 0.05191 12.997 1.16e-06 *** --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1Residual standard error: 3.253 on 8 degrees of freedom Multiple R-squared: 0.9548, Adjusted R-squared: 0.9491 F-statistic: 168.9 on 1 and 8 DF, p-value: 1.164e-06

lm(y~x)返回的是 截距和斜率，即身高和體重的線性回歸方程：y = 0.6746x -38.4551

SPSS

SPSS也可以做回歸分析。

在菜單欄中選擇【分析】，在下拉菜單中選擇【回歸】，右側彈出子菜單中選擇【散線性】，彈出線性回歸窗口，將“cholesterol“變量移至右側的【因變量】框中，”time”變量移至右側的【自變量】框選擇因變量和自變量，如下圖所示

tep1：在菜單欄中選擇【分析】，在下拉菜單中選擇【一般線性模型】中選擇【多變量】，彈出“多變量”窗口，將“身高”、“體重”選入右側【因變量】框中，將“班級”選入【固定因子】框中，如下圖所示

Step2：點擊【事后比較】在【多變量：實測平均值的事后比較】窗口中，將“class”選入右側【事后檢驗】中，并選中【圖基】，點擊【繼續】，如下圖所示。

Stata

SPSS出來了，我會放過Stata?

. set obs 10 number of observations (_N) was 0, now 10 . gen x=_n . gen y= x+runiform() . list+---------------+| x y ||---------------|1. | 1 1.348872 |2. | 2 2.266886 |3. | 3 3.136646 |4. | 4 4.028557 |5. | 5 5.868933 ||---------------|6. | 6 6.350855 |7. | 7 7.071105 |8. | 8 8.323368 |9. | 9 9.555103 |10. | 10 10.87599 |+---------------+. reg y xSource | SS df MS Number of obs = 10 -------------+---------------------------------- F(1, 8) = 1107.40Model | 89.9664664 1 89.9664664 Prob > F = 0.0000Residual | .649928434 8 .081241054 R-squared = 0.9928 -------------+---------------------------------- Adj R-squared = 0.9919Total | 90.6163949 9 10.0684883 Root MSE = .28503------------------------------------------------------------------------------y | Coef. Std. Err. t P>|t| [95% Conf. Interval] -------------+----------------------------------------------------------------x | 1.044271 .0313806 33.28 0.000 .9719076 1.116635_cons | .1391392 .1947113 0.71 0.495 -.3098658 .5881443 ------------------------------------------------------------------------------

在結果界面中， Coef為1.044271表示回歸斜率，_cons為.1391392表示回歸截距，R-squared和Adj R-squared分別為0.9928和0.9919，說明回歸方程擬合效果很好。

如果我們需要繪制對應圖形，只需要使用命令avplots就可以繪制折線圖，如下圖所示。

Excel

其實，最簡單的完成這個回歸分析的當然是Excel

現在需要求得當價格取得多少時，獲得最大的利潤。只需要CTRL+A全選數據，選擇插入X-Y散點圖，繪制帶平滑線和數據標記的散點圖。點擊坐標軸，右鍵，選擇坐標軸選項，設置邊界最小值和最大值。點擊散點圖，添加趨勢線。在趨勢線選項中，存在指數，線性，對數，多項式，乘冪和移動平均的回歸選擇，如下圖所示。

同時點擊顯示公式和R平方值，不斷地換回歸關系的選項，使R平方值靠近1，即擬合程度越高。這里選擇二階多項式。具體方程：y = -4.6716x2 + 286.75x - 2988.4，R2 = 0.9201，如下圖所示。

還有powerbi也是可以做回歸分析。

總結

一種簡單的線性擬合竟然有十幾種做法，對此，你會選擇哪一種下，是我肯定選擇Excel，20秒搞定。

假如有個需求就是要線性擬合，有的數據分析專家馬上跑DL，花了一個下午搞定了，卻不知道我用Excel，SPSS，20秒，點一點搞定，而且還你的好。

工具，要選擇利器。我記得有個廣告說：如何不加班通過Python完成上百份Excel報表的合并，你還不學Python。

我只想說一句：你學過Excel嗎？

我們點擊數據標簽下邊的【新建查詢】→【從文件】→【從文件夾】，如下圖：

然后進入PowerQuery界面瞬間搞定，干嘛一定要寫Python代碼？

再比如常見的因子和主成分分析，excel是搞不定，然后趕緊用Python做因子和主成分分析，需要嗎？

SPSS，Stata就是因子分析的好東西

再比如，之前做時間序列，結果我去用keras回滾3個月，就開始寫代碼，結果卻不知道經常使用的Stata 就可以做時間序列，突然覺得自己很蠢，白白浪費我時間

在stata中，幾乎所有的機器模型都有。

當你，打開jupyter notebook的時候，你可以想一想excel能不能處理，比如，換個數，那就直接在excel換唄，還裝什么導入pandas。

在做數據分析的時候，請不要敲代碼，因為根本不需要，excel，powerbi，SPSS，Stata絕對能搞定，但是敲代碼都是高手，其實我看就是簡單想搞復雜的“高手”。

高手請不要噴我，我知道你很厲害！

總結

以上是生活随笔為你收集整理的一个简单的线性拟合问题，到底有多少种做法的全部內容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網站內容還不錯，歡迎將生活随笔推薦給好友。