python 计量 python计量经济学建模 案例

目录

• 多元线性回归
• 因变量预测
• 化为线性的非线性实例
• 虚拟变量问题
• 受约束回归
• 多重共线性
• 异方差性
• 内生解释变量问题（待更）

多元线性回归

模型假设：
假设中国2013年各地区人均现金消费支出与工资性收入、其他收入之间的关系为:
Y$=\beta_0$+$\beta_1X_1$+$\beta_2X_2$+$\mu$
通过$python$的$statsmodels$库对数据进行回归计算：

import statsmodels.api as sm
import seaborn as sns
import numpy as np
import matplotlib.pyplot as plt
import pandas as pd
from sklearn import model_selection

data = pd.read_excel(r'./计量经济学数据.xlsx', sheet_name='Sheet1')
fit = sm.formula.ols(formula='现金消费支出Y ~ 工资性收入X1 + 其他收入X2', data=data).fit()
print(fit.summary())

sns.lmplot(x='工资性收入X1', y='现金消费支出Y', data=data, ci=None)
plt.rcParams['font.sans-serif'] = ['SimHei']
plt.rcParams['axes.unicode_minus'] = False
plt.show()
sns.pairplot(data.loc[:, ['现金消费支出Y', '工资性收入X1', '其他收入X2']])
# 显示图形
plt.show()

模型检验：
$H_0:\beta_j=0$
$H_1:\beta_j$不全部为零

拟合优度检验:
从回归估计来看,模型拟合较好，可决系数$R^{2}=0.922.$
F检验:
F值为166.6，查表得$F_{\alpha}(k,n-k-1)=3.34$，其中k=2，n=31，显然有$F>F_{\alpha}(k,n-k-1)$，表明模型的线性关系在5%的显著水平下显著成立.所以拒绝原假设。
t检验:
$\left|t_1\right|=8.448,\left|t_2\right|=5.772,t_{\alpha/2}(n-k-1)=2.048$
由于$\left|t\right|>t_{\alpha/2}(n-k-1)$，所以拒绝零假设.
综上可得中国2013年各地区人均现金消费支出与工资性收入、其他收入之间的关系为：
$Y=2599.1455+0.4865X_1+0.6017X_2$
$\beta_1<\beta_2$即其他收入对人均现金消费支出的贡献率要大于工资性收入。(此处有疑问，难道大家都靠外快过活吗？)

因变量预测

有时建立完模型并对其进行检验后，还需观察实际值和预测值具体情况，以确定模型的可用性。

data4 = pd.read_excel(r'./计量经济学数据.xlsx', sheet_name='Sheet1')
train, test = model_selection.train_test_split(data4, test_size=0.2, random_state=1234)
fit4 = sm.formula.ols(formula='现金消费支出Y ~ 工资性收入X1 + 其他收入X2', data=train).fit()
test_X = test.drop(labels='现金消费支出Y', axis=1)
pred = fit4.predict(exog=test_X)
print('对比预测值和实际值：\n', pd.DataFrame({'prediction': pred, 'real': test.现金消费支出Y}))

对比预测值和实际值：
    prediction    real
7   13874.648201  14161.7
10  25068.272118  23257.2
4   16645.508042  19249.1
1   21539.239415  21711.9
29  15077.077324  15321.1
8   28477.482744  28155.0
3   15073.999588  13166.2

由预测值和实际值对比可以看出，有的预测值和实际值相差比较大，但总体上来说预测值与实际值比较接近，也就一定程度上说明了这个模型的可用性。

化为线性的非线性实例

模型假设：

由Cobb-Dauglas生产函数$Y=AK^{\beta1}L^{\beta2}$，A代表既定的工程技术水平，$\beta_1$、$\beta_2$分别为资本与劳动投入的产出弹性，当$\beta_1+\beta_2=1时$，当大于1或小于1时，表明规模收益递增或递减。为了便于比较，下面将会对此模型进行线性变换，即假设2010年中国制造业各行业的总产出及要素投入的关系为：
$Y=\beta_0+\beta_1\log K+\beta_2\log L+\mu$

data2 = pd.read_excel(r'./计量经济学数据.xlsx', sheet_name='Sheet2')
fit2 = sm.formula.ols(formula='np.log(工业总产值) ~ np.log(资本投入) + np.log(年均从业人员)', data=data2).fit()
fit2.summary()

sns.pairplot(data2.loc[:, ['工业总产值', '资本投入', '年均从业人员']])
plt.show()

模型检验：

$H_0:\beta_j=0$
$H_1:\beta_j$不全部为零
拟合优度检验:
从回归估计来看,模型拟合较好，可决系数$R^{2}=0.941.$
F检验:
F值为286.3，查表得$F_{\alpha}(k,n-k-1)=3.26$，其中$k=2$，$n=39$，显然有$F>F_{\alpha}(k,n-k-1)$，表明模型的线性关系在5%的显著水平下显著成立.所以拒绝原假设。
t检验:
$\left|t_1\right|=8.344,\left|t_2\right|=3.395,t_{\alpha/2}(n-k-1)=2.036$
由于$\left|t\right|>t_{\alpha/2}(n-k-1)$，所以拒绝零假设.
综上可得2010年中国制造业各行业的总产出及要素投入的关系为：$Y=1.8003+0.6778\log K+0.2911\log L，0.6778+0.2911=0.9689$，以上结果表明，在2010年，中国工业总产出关于资本投入的产出弹性为0.6778，表明当其他因素不变时，工业的资本每增加1%，总产出将增加0.6778%，同样地，当其他因素不变时，劳动力投入每增长1%，总产出将增加0.2911%，可见，资本投入的增加对工业总产出的增长起到了更大的作用。

虚拟变量问题

在一些数据中，通常会有一些变量无法通过量化来进行处理，但是这些变量往往对模型结果产生较大的影响，所以，这类因素是无法被丢弃的，因此引入了“虚拟变量”，又叫做哑变量，来进行“量化处理”。下面我们将会以城镇居民为基准线对2013年中国农村与城镇居民家庭人均工资收入、其他收入和生活消费支出进行模型建立。
假设模型为：
$Y=\alpha_0+\alpha_1X_1+\alpha_2X_2+C$

data3 = pd.read_excel(r'./计量经济学数据.xlsx', sheet_name='Sheet3')
fit3 = sm.formula.ols(formula='生活消费 ~ 工资收入 + 其他收入 + C(农村or城镇)', data=data3).fit()
fit3.summary()

模型检验：

$H_0:\beta_j=0$
$H_1:\beta_j$不全部为零
拟合优度检验:
从回归估计来看,模型拟合较好，可决系数$R^{2}=0.975.$
F检验:
F值为758.1，查表得$F_{\alpha}(k,n-k-1)=4.16$，其中k=3，n=62，显然有$F>F_{\alpha}(k,n-k-1)$，表明模型的线性关系在5%的显著水平下显著成立.所以拒绝零假设。
t检验:
$\left|t_1\right|=0.291,\left|t_2\right|=13.978,,\left|t_3\right|=7.666，t_{\alpha/2}(n-k-1)=2.010$
$由于\left|t\right|>t_{\alpha/2}(n-k-1)$，所以拒绝零假设.
综上可得2013年中国农村与城镇居民家庭人均工资收入、其他收入和生活消费支出的关系为：
$Y=1783.7377+0.5477X_1+0.5589X_2+140.8608$城镇居民，
以上结果表明，当其他因素不变时，中国城镇居民平均消费支出比农村居民平均消费水平多140.8608元。

受约束回归

在建立回归模型时，有时根据经济理论需要对自变量之间的关系进行约束，比如两个回归系数$\beta_1$和$\beta_2$之间的约束条件使得$\beta_1+\beta_2=1$或者使得$\beta_1=\beta_2$，此时称为此回归模型为受约束回归。

首先建立无约束回归模型
即：$ln(Q)=\beta_0+\beta_1ln(X/P_0)+\beta_2(P_1/P)+\beta_3(P_2/P)+\beta_4P_{01}+\beta_5P_{02}+\beta_6P_{03}$

import statsmodels.api as sm
data4 = pd.read_excel(r'./计量经济学数据.xlsx', sheet_name='Sheet4')
Q = data4["蛋类消费量Q（千克）"]
X = data4["人均消费支出X（元）"]
P0 = data4["居民消费价格指数P0"]
P = data4["蛋类P（价格指数）"]
P1 = data4["肉禽类P1（价格指数）"]
P2 = data4["水产类P2（价格指数）"]
P01 = data4["粮食P3（价格指数）"]
P02 = data4["油脂P4（价格指数）"]
P03 = data4["蔬菜P5（价格指数）"]
df = pd.DataFrame({"log(X/P0)":np.log(X/P0),
                  "P1/P":P1/P,
                 "P2/P":P2/P,
                  "P01":P01,
                  "P02":P02,
                  "P03":P03},)
df = sm.add_constant(df)
fit = sm.OLS(np.log(Q),df).fit()
fit.summary()

模型检验：
$H_0:\beta_j=0$
$H_1:\beta_j$不全部为零
拟合优度检验:
从回归估计来看，调整的可决系数$R^{2}=0.409$，但是此回归结果并不作为预测模型来进行预测，所以可以不必过分关注可决系数.
F检验:
F值为 4.641，$P=0.00290$，表明模型的线性关系在5%的显著水平下显著成立，所以拒绝零假设。
t检验:
以上结果除变量log(X/P0) 和P2/P在5%的显著水平下拒绝原假设，其他变量均无法通过t检验，在其他条件不变的情况下，农村人均消费支出会明显增加蛋类消费量，同时，当水产类价格上升速度大于蛋类产品时，会刺激农村消费者倾向于消费更多的蛋类产品，即在农村消费者的消费倾向中，水产品类与蛋类产品有一定的替代作用。

建立受约束回归模型
约束条件为$H_0:\beta_2=\beta_4=\beta_5=\beta_6=0$，即回归模型为$ln(Q)=\beta_0+\beta_1ln(X/P_0)+\beta_3(P2/P)$

从以上结果可以看出，在约束条件下，线性关系检验（F检验）和回归系数检验（t检验）在5%的显著水平下更加显著，拒绝原假设的理由更加充分，即更加印证了无约束回归所说明的结论。

多重共线性

检测方法以及解决办法详见下面链接。
链接: 多重共线性.

异方差性

异方差性：对于模型$Y_t=\beta_0+\beta_1X_{i1}+\beta_2X_{i2}+...+\beta_nX_{in}+\mu_i，i=1,2...,n$，其中残差$\mu_i$的方差随着$X_i$的变化而变化，而不是一个常数，这就说明模型存在异方差性，当然这只是模型假设的一个理想条件。
异方差性的后果：（1）参数估计量非有效（2）变量的显著性检验失去意义（3）模型预测失效
检验方法：（1）图示法（2）布罗施-帕甘检验（3）怀特检验

下面我们通过参考书上的案例来进行具体分析。
模型假设为$ln(Y)=\beta_0+\beta_1ln(X_1)+\beta_2ln(X_2)+\mu$

data7 = pd.read_excel(r'./计量经济学数据.xlsx', sheet_name='Sheet7')
X = np.log(data7.loc[:,[ '从事农业经营的纯收入X1', '其他来源纯收入X2']])
Y = np.log(data7["人均消费支出Y"])
X = sm.add_constant(X)
fit = sm.OLS(Y,X).fit()
fit.summary()

从拟合结果来看，很明显常数项、变量$X_1$和$X_2$都无法通过$t$检验，从拟合参数来看，变量$X_2$对因变量的解释性更大，所以如果存在异方差性，则更有可能变量$X_2$引起的，下面我们来进行异方差性检验。

1. 图示法

计算残差并绘制异方差性检验图，即$e^2$与$X_2$的散点图。

plt.figure(figsize=(12,8))
plt.scatter(X["其他来源纯收入X2"],np.power(fit.predict(X)-Y,2))
fit.predict(X)
plt.title("异方差性检验图",fontsize=16)
plt.show()

从图中可以看出，拟合模型中可能存在异方差性，而且存在着递增型异方差性，但是从图形并不能很直接的说明问题，存在一定的主观因素在里面，所以下面另外两种方法更为直接地检验一下。

2. 布罗施-帕甘检验

将原模型普通最小二乘法估计的残差项的平方后关于变量$X_2$做回归。
模型假设为$e_2=\beta_0+\beta_2X_2$

e2=np.power(fit.predict(X)-Y,2)
X2 = sm.add_constant(X["其他来源纯收入X2"])
fit2 = sm.OLS(e2,X2).fit()
fit2.summary()

从上图可以得出残差平方与变量$X_2$的拟合方程$e^2=-141+0.0236X_2$，同时在5%的显著水平都通过了$t$检验，同时在5%显著水平下也通过了$F$检验。综上即可说明在5%的显著水平下拒绝原模型随机干扰项方差相同的原假设，即原模型存在异方差性。

3. 怀特检验

记$e_i^2$为对原模型进行普通最小二乘回归得到的残差平方，将其与$ln(X_1)、ln(X_2)$及其平方项与交叉项作辅助回归。
模型假设为$e_2=\beta_0+\beta_1ln(X_1)+\beta_2ln(X_2)+\beta_3(ln(X_1))^2+\beta_4(ln(X_1))^2+\beta_5ln(X_1)ln(X_2)$

X["从事农业经营的纯收入(X1)^2"]=np.power(X["从事农业经营的纯收入X1"],2)
X["其他来源纯收入(X2)^2"]=np.power(X["其他来源纯收入X2"],2)
X["X1*X2"]=X["从事农业经营的纯收入X1"]*X["其他来源纯收入X2"]
X3=X
fit3 = sm.OLS(e2,X3).fit()
fit3.summary()

从上图可以清楚的看到，在5%的显著水平下，拟合方程在线性关系检验上是显著的，也就是拒绝了原拟合方程同方差的原假设，即原拟合方程存在异方差性。

内生解释变量问题（待更）