文章目录
1 前序
常见的两类方法:
- 基于相关原则的回归分析法
- 基于惯性原则的时间序列法
回归分析的基本步骤:
- (1)重点考察一个特定的变量(因变量),而把其他变量(自变量)看作是影响这一变量的因素,并通过适当的数学模型将变量间的关系表达出来;
- (2)利用样本数据建立模型的估计方程;
- (3)对模型进行显著性检验;
- (4)通过一个或几个自变量的取值来估计或预测因变量的取值。
预测方法的选择(仅供参考):
数据模式 | 预测方法 | 对数据的要求 | 预测期 |
---|---|---|---|
平稳序列 | 移动平均 | 数据个数与移动平均步长相等 | 非常短 |
平稳序列 | 简单指数平滑 | 5个以上 | 短期 |
线性趋势 | Holt指数平滑 | 5个以上 | 短期至中期 |
线性趋势 | 一元线性回归 | 10个以上 | 短期至中期 |
非线性趋势 | 指数模型 | 10个以上 | 短期至中期 |
非线性趋势 | 多项式函数 | 10个以上 | 短期至中期 |
趋势和季节成分 | Winter指数平滑 | 至少有四个周期的季节或月份数据 | 短期至中期 |
趋势和季节成分 | 季节性多元回归 | 至少有四个周期的季节或月份数据 | 短期、中期、长期 |
趋势、季节成分和循环成分 | 分解预测 | 至少有四个周期的季节或月份数据 | 短期、中期、长期 |
预测方法的评估:一种预测方法的好坏取决于预测误差的大小。预测误差是预测值与实际值的差距,度量方法有:
- 平均误差(Mean Error)
- 平均绝对误差(Mean Absolute Deviation)
- 均方误差(Mean Square Error,MSE)(常用)
- 平均百分比误差(Mean Percentage Error)
- 平均绝对百分比误差(Mean Absolute Percentage Error)
2 预测方法及案例
2.1 回归分析
2.1.1 含有哑变量的线性回归分析案例
为研究员工月工资收入与工作年限和性别之间的关系,从某公司职员中随机抽取男女各4名,他们的月工资收入与工作年限和性别之间的关系表如下:
月工资收入(元) | 工作年限 | 性别 |
---|---|---|
2900 | 2 | 男 |
3000 | 6 | 女 |
4800 | 8 | 男 |
1800 | 3 | 女 |
2900 | 2 | 男 |
4900 | 7 | 男 |
4200 | 9 | 女 |
4800 | 8 | 女 |
令 y y y表示月工资收入, x 1 x_1 x1表示工作年限, x 2 x_2 x2表示性别,性别作为哑变量引入时,回归方程如下: y = β 0 + β 1 x 1 + β 2 x 2 y=\beta_0+\beta_1 x_1 + \beta_2 x_2 y=β0+β1x1+β2x2,于是我们可以得到:
- 女( x 2 = 0 x_2=0 x2=0): y 女 性 = β 0 + β 1 x 1 y_{女性}=\beta_0+\beta_1 x_1 y女性=β0+β1x1
- 男( x 2 = 1 x_2=1 x2=1): y 男 性 = ( β 0 + β 2 ) + β 1 x 1 y_{男性}=(\beta_0+\beta_2)+\beta_1 x_1 y男性=(β0+β2)+β1x1
其中各参数的含义如下:
- β 0 \beta_0 β0的含义是女性职工的基本月工资收入
- ( β 0 + β 2 ) (\beta_0+\beta_2) (β0+β2)的含义是男性职工的基本月工资收入
- β 1 \beta_1 β1的含义是工作年限每增加1年,男性或女性工资的平均增加值
- β 2 \beta_2 β2的含义是男性职工的月工资收入与女性职工的月工资收入之间的差值,即 y 男 性 − y 女 性 = ( β 0 + β 2 ) + β 1 x 1 − β 0 + β 1 x 1 = β 2 y_{男性}-y_{女性}=(\beta_0+\beta_2)+\beta_1 x_1-\beta_0+\beta_1 x_1=\beta_2 y男性−y女性=(β0+β2)+β1x1−β0+β1x1=β2
python实现代码如下:
import pandas as pd
import numpy as np
import statsmodels.api as sm
data = pd.DataFrame({
'月工资收入':[2900,3000,4800,1800,2900,4900,4200,4800],
'工作年限':[2,6,8,3,2,7,9,8],
'性别':['男','女','男','女','男','男','女','女']
})
# 哑变量处理
dummy_variables = pd.get_dummies(data=data['性别'].values)
X = np.column_stack(tup=(data['工作年限'].values,dummy_variables.values))
X = sm.add_constant(data=X) # 加一列常数项
y = data['月工资收入'].values
# 用最小二乘法拟合回归方程
linear_model = sm.OLS(endog=y,exog=X)
ols_result = linear_model.fit()
# 输出拟合结果
print(ols_result.summary())
'''
OLS Regression Results
==============================================================================
Dep. Variable: y R-squared: 0.901
Model: OLS Adj. R-squared: 0.862
Method: Least Squares F-statistic: 22.78
Date: Sun, 22 May 2022 Prob (F-statistic): 0.00307
Time: 17:02:55 Log-Likelihood: -58.036
No. Observations: 8 AIC: 122.1
Df Residuals: 5 BIC: 122.3
Df Model: 2
Covariance Type: nonrobust
==============================================================================
coef std err t P>|t| [0.025 0.975]
------------------------------------------------------------------------------
const 950.7246 247.685 3.838 0.012 314.029 1587.420
x1 397.5845 60.183 6.606 0.001 242.879 552.290
x2 -85.0242 231.138 -0.368 0.728 -679.184 509.135
x3 1035.7488 172.207 6.015 0.002 593.076 1478.421
==============================================================================
Omnibus: 4.593 Durbin-Watson: 1.536
Prob(Omnibus): 0.101 Jarque-Bera (JB): 1.483
Skew: 1.049 Prob(JB): 0.477
Kurtosis: 3.219 Cond. No. 7.55e+16
==============================================================================
Notes:
[1] Standard Errors assume that the covariance matrix of the errors is correctly specified.
[2] The smallest eigenvalue is 5.64e-32. This might indicate that there are
strong multicollinearity problems or that the design matrix is singular.
从上面结果可以看到,拟合优度R方为0.901,调整R方为0.862,模型拟合显著性检验的P值为0.00307 < 0.05,说明模型拟合效果还是可以的,
再看到参数拟合结果,x2的参数p值为0.728 > 0.05,因此下面我们剔除x2后再做一次拟合
'''
X1 = X[:,[0,1,3]]
# 用最小二乘法拟合回归方程
linear_model1 = sm.OLS(endog=y,exog=X1)
ols_result1 = linear_model1.fit()
# 输出拟合结果
print(ols_result1.summary())
'''
OLS Regression Results
==============================================================================
Dep. Variable: y R-squared: 0.901
Model: OLS Adj. R-squared: 0.862
Method: Least Squares F-statistic: 22.78
Date: Sun, 22 May 2022 Prob (F-statistic): 0.00307
Time: 17:18:21 Log-Likelihood: -58.036
No. Observations: 8 AIC: 122.1
Df Residuals: 5 BIC: 122.3
Df Model: 2
Covariance Type: nonrobust
==============================================================================
coef std err t P>|t| [0.025 0.975]
------------------------------------------------------------------------------
const 865.7005 447.091 1.936 0.111 -283.583 2014.984
x1 397.5845 60.183 6.606 0.001 242.879 552.290
x2 1120.7729 323.747 3.462 0.018 288.554 1952.992
==============================================================================
Omnibus: 4.593 Durbin-Watson: 1.536
Prob(Omnibus): 0.101 Jarque-Bera (JB): 1.483
Skew: 1.049 Prob(JB): 0.477
Kurtosis: 3.219 Cond. No. 20.8
==============================================================================
Notes:
[1] Standard Errors assume that the covariance matrix of the errors is correctly specified.
可以看到,模型整体的拟合是没有变的,在显著性水平为0.05的条件下,x1和x2的参数显著性通过检验,但常数项的显著性不明显,因此下面我们剔除常数项再做一次拟合
'''
X2 = X[:,[1,3]]
# 用最小二乘法拟合回归方程
linear_model2 = sm.OLS(endog=y,exog=X2)
ols_result2 = linear_model2.fit()
# 输出拟合结果
print(ols_result2.summary())
'''
OLS Regression Results
=======================================================================================
Dep. Variable: y R-squared (uncentered): 0.986
Model: OLS Adj. R-squared (uncentered): 0.981
Method: Least Squares F-statistic: 210.6
Date: Sun, 22 May 2022 Prob (F-statistic): 2.77e-06
Time: 17:22:47 Log-Likelihood: -60.274
No. Observations: 8 AIC: 124.5
Df Residuals: 6 BIC: 124.7
Df Model: 2
Covariance Type: nonrobust
==============================================================================
coef std err t P>|t| [0.025 0.975]
------------------------------------------------------------------------------
x1 499.5470 35.188 14.197 0.000 413.446 585.648
x2 1502.1518 310.270 4.841 0.003 742.948 2261.356
==============================================================================
Omnibus: 0.202 Durbin-Watson: 1.923
Prob(Omnibus): 0.904 Jarque-Bera (JB): 0.253
Skew: -0.258 Prob(JB): 0.881
Kurtosis: 2.297 Cond. No. 10.5
==============================================================================
Notes:
[1] R² is computed without centering (uncentered) since the model does not contain a constant.
[2] Standard Errors assume that the covariance matrix of the errors is correctly specified.
可以看大,模型显著性和参数显著性都得到了提升,最终我们模型估计结果为:y = 499.5470*x1 + 1502.1518*x2
'''
2.1.2 自变量之间有交互作用的回归分析案例
某牙膏制造企业要求销售部门根据市场调查,找出公司生产的牙膏销售量与销售价格、广告投入等之间的关系,从而预测出在不同价格和广告费用下的销售量。为此,销售部收集了过去30个销售周期公司生产的牙膏的销售量、销售价格、广告费用,以及同期其他厂家生产的同类牙膏的平均销售价格,数据如下:
销售周期 | 公司销售价格(元) | 其他厂家平均价格(元) | 广告费用(百万元) | 价格差(元) | 销售量(百万支) |
---|---|---|---|---|---|
1 | 3.85 | 3.8 | 5.5 | -0.05 | 7.38 |
2 | 3.75 | 4 | 6.75 | 0.25 | 8.51 |
3 | 3.7 | 4.3 | 7.25 | 0.6 | 9.52 |
…… | …… | …… | …… | …… | …… |
试根据这些数据建立一个数学模型,分析牙膏销售量与其他因素的关系,为制定价格策略和广告投入策略提供数据依据
在购买同类产品的牙膏时,顾客会在意不同品牌之间的价格差异,而不是价格本身。设牙膏销售量为y,价格差为x1,广告费为x2,其他厂家平均价格为x3,公司销售价格为x4,x1=x3-x4。
回归分析过程如下:
- (1)计算相关系数矩阵,观察自变量y与各因变量的相关程度,以及各因变量之间是否存在相关性;
- (2)绘制自变量与因变量的散点图,通过散点图初步判断模型是线性还是非线性;
- (3)建立模型,并对模型及其参数进行评估;
- (4)模型改进。在初始假定中,价格差和广告费对牙膏销售量的影响是相互独立的,而据直觉和经验可以猜想,价格差与广告费的交互作用也会对牙膏销售量有影响。
- (5)模型最终确定。
python代码如下:
import pandas as pd
import numpy as np
import statsmodels.api as sm
from matplotlib import pyplot as plt
plt.rcParams['font.sans-serif'] = ['SimHei'] # 解决中文显示问题
plt.rcParams['axes.unicode_minus'] = False # 解决坐标轴负数的负号显示问题
data = pd.read_excel(r'G:\牙膏销售数据表.xlsx')
print(data.info())
# 1、计算相关系数矩阵
print(data[['公司销售价格(元)','其他厂家平均价格(元)','广告费用(百万元)','价格差(元)','销售量(百万支)']].corr())
'''
公司销售价格(元) 其他厂家平均价格(元) 广告费用(百万元) 价格差(元) 销售量(百万支)
公司销售价格(元) 1.000000 0.078367 -0.468793 -0.322067 -0.469220
其他厂家平均价格(元) 0.078367 1.000000 0.604540 0.918566 0.740948
广告费用(百万元) -0.468793 0.604540 1.000000 0.759964 0.875954
价格差(元) -0.322067 0.918566 0.759964 1.000000 0.889672
销售量(百万支) -0.469220 0.740948 0.875954 0.889672 1.000000
可以看到,销售量与价格差、广告费用的相关系数分别是89%、88%,相关度很高,另外价格差和广告费用的相关度也不低(76%左右)
'''
# 2、绘制自变量与因变量的散点图
fig1 = plt.figure()
plt.scatter(x=data['广告费用(百万元)'].values,y=data['销售量(百万支)'].values,marker='*')
plt.xlabel('广告费用(百万元)')
plt.ylabel('销售量(百万支)')
fig2= plt.figure()
plt.scatter(x=data['价格差(元)'].values,y=data['销售量(百万支)'].values,marker='o')
plt.xlabel('价格差(元)')
plt.ylabel('销售量(百万支)')
plt.show()
'''
可以看到,价格差对销售量的走势是近似线性的,而广告费对销售量的走势有一点弯曲,
因此在假设价格差和广告费对牙膏销售量的影响是相互独立的前提下,准对价格差和销售量建立线性回归方程,准对广告费和销售量建立非线性回归方程(二次项)
'''
# 3、建立模型,并对模型及其参数进行估计
# 3.1、价格差和销售量的线性回归方程估计
X1 = sm.add_constant(data=data['价格差(元)'].values) # X1由1列常数项和价格差组成
y = data['销售量(百万支)']
# 用最小二乘法估计方程
model1 = sm.OLS(endog=y,exog=X1)
result1 = model1.fit()
print(result1.summary())
'''
OLS Regression Results
==============================================================================
Dep. Variable: 销售量(百万支) R-squared: 0.792
Model: OLS Adj. R-squared: 0.784
Method: Least Squares F-statistic: 106.3
Date: Fri, 27 May 2022 Prob (F-statistic): 4.88e-11
Time: 22:26:27 Log-Likelihood: -7.0261
No. Observations: 30 AIC: 18.05
Df Residuals: 28 BIC: 20.85
Df Model: 1
Covariance Type: nonrobust
==============================================================================
coef std err t P>|t| [0.025 0.975]
------------------------------------------------------------------------------
const 7.8141 0.080 97.818 0.000 7.650 7.978
x1 2.6652 0.258 10.310 0.000 2.136 3.195
==============================================================================
Omnibus: 5.481 Durbin-Watson: 2.414
Prob(Omnibus): 0.065 Jarque-Bera (JB): 4.092
Skew: 0.883 Prob(JB): 0.129
Kurtosis: 3.391 Cond. No. 4.69
==============================================================================
Notes:
[1] Standard Errors assume that the covariance matrix of the errors is correctly specified.
可以看到,模型估计的拟合优度为0.792,在假设显著性水平α=0.05的条件下,模型显著性检验和参数显著性检验均通过,得到的回归方程为:y = 7.8141 + 2.6652*x1
'''
# 3.2、广告费和销售量的非线性回归方程估计
# 多项式处理
from sklearn.preprocessing import PolynomialFeatures
X2 = PolynomialFeatures(degree=2).fit_transform(X=data[['广告费用(百万元)']])
# 用最小二乘法估计方程
model2 = sm.OLS(endog=y,exog=X2)
result2 = model2.fit()
print(result2.summary())
'''
OLS Regression Results
==============================================================================
Dep. Variable: 销售量(百万支) R-squared: 0.838
Model: OLS Adj. R-squared: 0.826
Method: Least Squares F-statistic: 69.81
Date: Fri, 27 May 2022 Prob (F-statistic): 2.14e-11
Time: 22:45:32 Log-Likelihood: -3.2455
No. Observations: 30 AIC: 12.49
Df Residuals: 27 BIC: 16.69
Df Model: 2
Covariance Type: nonrobust
==============================================================================
coef std err t P>|t| [0.025 0.975]
------------------------------------------------------------------------------
const 25.1091 6.863 3.659 0.001 11.028 39.190
x1 -6.5589 2.217 -2.958 0.006 -11.109 -2.009
x2 0.6101 0.178 3.432 0.002 0.245 0.975
==============================================================================
Omnibus: 0.063 Durbin-Watson: 1.523
Prob(Omnibus): 0.969 Jarque-Bera (JB): 0.224
Skew: -0.090 Prob(JB): 0.894
Kurtosis: 2.617 Cond. No. 5.98e+03
==============================================================================
Notes:
[1] Standard Errors assume that the covariance matrix of the errors is correctly specified.
[2] The condition number is large, 5.98e+03. This might indicate that there are
strong multicollinearity or other numerical problems.
可以看到,模型估计的拟合优度为0.838,在假设显著性水平α=0.05的条件下,模型显著性检验和参数显著性检验均通过,得到的回归方程为:y = 25.1091 - 6.5589*x2 + 0.6101*(x2)**2
'''
# 4、模型改进:加入交互项和广告费的二次项
data['价格差*广告费用'] = data['价格差(元)'] * data['广告费用(百万元)']
X3 = np.column_stack(tup=(X1,X2[:,[1,2]],data['价格差*广告费用'].values))
# 用最小二乘法估计方程
model3 = sm.OLS(endog=y,exog=X3)
result3 = model3.fit()
print(result3.summary())
'''
OLS Regression Results
==============================================================================
Dep. Variable: 销售量(百万支) R-squared: 0.921
Model: OLS Adj. R-squared: 0.908
Method: Least Squares F-statistic: 72.78
Date: Fri, 27 May 2022 Prob (F-statistic): 2.11e-13
Time: 23:11:51 Log-Likelihood: 7.5137
No. Observations: 30 AIC: -5.027
Df Residuals: 25 BIC: 1.979
Df Model: 4
Covariance Type: nonrobust
==============================================================================
coef std err t P>|t| [0.025 0.975]
------------------------------------------------------------------------------
const 29.1133 7.483 3.890 0.001 13.701 44.525
x1 11.1342 4.446 2.504 0.019 1.978 20.291
x2 -7.6080 2.469 -3.081 0.005 -12.693 -2.523
x3 0.6712 0.203 3.312 0.003 0.254 1.089
x4 -1.4777 0.667 -2.215 0.036 -2.852 -0.104
==============================================================================
Omnibus: 0.242 Durbin-Watson: 1.512
Prob(Omnibus): 0.886 Jarque-Bera (JB): 0.148
Skew: -0.153 Prob(JB): 0.929
Kurtosis: 2.843 Cond. No. 9.81e+03
==============================================================================
Notes:
[1] Standard Errors assume that the covariance matrix of the errors is correctly specified.
[2] The condition number is large, 9.81e+03. This might indicate that there are
strong multicollinearity or other numerical problems.
可以看到,模型估计的拟合优度提升至0.921,在假设显著性水平α=0.05的条件下,模型显著性检验和参数显著性检验均通过,最终确定回归方程为:
y = 29.1133 + 11.1342*x1 - 7.6080*x2 + 0.6712*(x2)**2 - 1.4777*x1*x2
'''
2.1.3 非线性回归分析——预测第三产业国内生产总值案例
某省为了研究第三产业在本省宏观经济发展中的运行情况,对影响第三产业的资本要素、劳动力要素及科技进步要素这三项主要因素进行了统计分析,并运用道格拉斯生产函数建立了基本数学模型。原始数据如下:
年份 | 第三产业国内生产总值 | 资本投入 | 从业人员 |
---|---|---|---|
1992 | 448.96 | 524085 | 470.08 |
1993 | 611.23 | 1068889 | 480.77 |
1994 | 834.93 | 1632884 | 529.08 |
…… | …… | …… | …… |
2002 | 3120 | 10029357 | 903.14 |
现需要预测当资本投入为11738245、劳动力投入为987.37时,第三产业国内生产总值是多少。
道格拉斯生产函数对应数学模型为 Y = A K α L β Y=AK^{\alpha}L^{\beta} Y=AKα