5.7 计量分析7：时间序列分析

例子

代码部分

# lec07
# 计量实验5-代码部分 自相关的诊断与处理

# 导入程序包
import pandas as pd
import statsmodels.api as sm
import matplotlib.pyplot as plt
from statsmodels.tsa.stattools import adfuller
from statsmodels.tsa.stattools import coint
from statsmodels.stats.diagnostic import acorr_breusch_godfrey
import ssl

# 禁用SSL证书校验
ssl._create_default_https_context = ssl._create_unverified_context

# 导入(在线)数据
df = pd.read_excel(r"https://cdn.seit2019.xyz/data/econometrics/exp3.1.xlsx")

# 展示数据
print("数据展示如下：")
print(df, "\n")

# 简要统计描述
print("简要统计描述：")
print(df.describe(), "\n")

# 名义值转换为实际值
df01 = pd.DataFrame()
df01["cons"] = df["cons"].div(df["cpi"].div(100))
df01["dpi"] = df["dpi"].div(df["cpi"].div(100))

# 描述性统计分析
print("数据的简要统计描述：")
print(df01.describe(), "\n")

# 折线图
# 设置绘图格式,支持汉字
plt.rcParams['font.sans-serif'] = ['SimHei']
plt.rcParams['axes.unicode_minus'] = False
# 定义变量
x = df["year"]
y01 = df["cons"]
y02 = df["dpi"]
y1 = df01["cons"]
y2 = df01["dpi"]
# 绘图：农村居民人均消费与支出
plt.plot(x, y01, c="blue", label="农村居民人均消费")
plt.plot(x, y02, c="green", label="农村居民人均纯收入")
plt.title("1985-2012农村居民人均消费与纯收入(按当年价格)", fontsize="x-large", fontweight="bold")
plt.legend()
plt.show()
# 绘图：实际农村居民人均消费与支出
plt.plot(x, y1, c="blue", label="农村居民人均消费")
plt.plot(x, y2, c="green", label="农村居民人均纯收入")
plt.title("1985-2012实际农村居民人均消费与纯收入(按1978价格)", fontsize="x-large", fontweight="bold")
plt.legend()
plt.show()

# XY散点图
# 定义变量
x = df01["dpi"]
y = df01["cons"]
plt.scatter(x, y, c="green")
plt.title("农村居民人均消费与纯收入散点图", fontsize="x-large", fontweight="bold")
plt.xlabel("农村居民人均纯收入(元)", fontsize="large", fontweight="bold")
plt.ylabel("农村居民人均消费(元)", fontsize="large", fontweight="bold")
plt.show()

# 相关系数
print("相关系数：")
print(df01.corr(), "\n")

# 平稳性检验—单位根检验
# 单位根检验：y
adf_stat = adfuller(y, maxlag=0)
adf_y = adf_stat[0]
adf_pvalue_y = adf_stat[1]
# 单位根检验：y的一阶差分
y_diff1 = y.diff(1)
y_diff1 = y_diff1.dropna()
adf_stat = adfuller(y_diff1, maxlag=0)
adf_y_diff1 = adf_stat[0]
adf_pvalue_y_diff1 = adf_stat[1]
# 单位根检验：y的二阶差分
y_diff2 = y_diff1.diff(1)
print(y_diff2)
y_diff2 = y_diff2.dropna()
adf_stat = adfuller(y_diff2, maxlag=0)
adf_y_diff2 = adf_stat[0]
adf_pvalue_y_diff2 = adf_stat[1]
# 单位根检验：x
adf_stat = adfuller(x, maxlag=0)
adf_x = adf_stat[0]
adf_pvalue_x = adf_stat[1]
# 单位根检验：x的一阶差分
x_diff1 = x.diff(1)
x_diff1 = x_diff1.dropna()
adf_stat = adfuller(x_diff1, maxlag=0)
adf_x_diff1 = adf_stat[0]
adf_pvalue_x_diff1 = adf_stat[1]
# 单位根检验：x的二阶差分
x_diff2 = x_diff1.diff(1)
x_diff2 = x_diff2.dropna()
adf_stat = adfuller(x_diff2, maxlag=0)
adf_x_diff2 = adf_stat[0]
adf_pvalue_x_diff2 = adf_stat[1]
# 平稳性性-单位根检验结果汇总
df_adf = pd.DataFrame(columns=["变量", "ADF值", "P值", "平稳性"])
df_adf["变量"] = ["cons", "cons_diff1", "cons_diff2", "dpi", "dpi_diff1", "dpi_diff2"]
df_adf["ADF值"] = [adf_y, adf_y_diff1, adf_y_diff2, adf_x, adf_x_diff1, adf_x_diff2]
df_adf["P值"] = [adf_pvalue_y, adf_pvalue_y_diff1, adf_pvalue_y_diff2, adf_pvalue_x, adf_pvalue_x_diff1, adf_pvalue_x_diff2]
adf_p = df_adf["P值"]
for i in range(6):
    if df_adf.iloc[i-1, 2] < 0.05:
        df_adf.iloc[i-1, 3] = "平稳"
    else:
        df_adf.iloc[i-1, 3] = "不平稳"
print("平稳性检验的汇总结果：")
print(df_adf.to_string(index=False), "\n")

# 回归分析
# 定义变量
cons = y_diff2
dpi = x_diff2
# 添加常数项
X = sm.add_constant(dpi)
# 回归
OLS_model = sm.OLS(cons, X).fit()
print("回归结果：")
print(OLS_model.summary(), "\n")

# 残差自相关检验
# Breusch-Godfrey test
stats_bg = acorr_breusch_godfrey(OLS_model, nlags=1)
print("残差的Breusch-Godfrey检验结果：")
if stats_bg[1] < 0.05:
    print("残差序列存在一阶自相关。")
else:
    print("残差序列不存在一阶自相关。")

# 广义差分
# 定义变量
y = OLS_model.resid
x = y.shift(1)
# 辅助回归
OLS_aux = sm.OLS(y, x, missing="drop").fit()
print(OLS_aux.summary(), "\n")
rho = OLS_aux.params[0].round(4)
print("广义差分的参数rho值为：", rho, "\n")
# 广义差分
cons_gdiff = cons - rho*(cons.shift(1))
dpi_gdiff = dpi - rho*(dpi.shift(1))
# 回归分析：广义差分数据
X = sm.add_constant(dpi_gdiff)
OLS_Gdiff = sm.OLS(cons_gdiff, X, missing="drop").fit()
print("广义差分后的回归结果：")
print(OLS_Gdiff.summary(), "\n")

# 残差自相关检验
# Breusch-Godfrey test
print("广义差分模型残差检验结果：")
stats_bg = acorr_breusch_godfrey(OLS_Gdiff, nlags=1)
if stats_bg[1] < 0.05:
    print("残差序列存在一阶自相关。", "\n")
else:
    print("残差序列不存在一阶自相关。", "\n")

# 协整检验
# 检验一阶差分后的数据是否具有协整关系
stats_coint = coint(y_diff1, x_diff1)
print("协整检验结果：")
if stats_coint[1] < 0.05:
    print("经过一阶差分后的cons和dpi之间存在协整关系。", "\n")
else:
    print("经过一阶差分后的cons和dpi之间不存在协整关系。", "\n")

# 基于协整关系的回归
# 定义变量
cons = y_diff1
dpi = x_diff1
X = sm.add_constant(dpi)
OLS_coint = sm.OLS(cons, X).fit()
print("基于协整关系的回归：")
print(OLS_coint.summary(), "\n")
# 残差检验
# 残差平稳性检验
print("协整关系模型残差检验结果：")
e_coint = OLS_coint.resid
adf_e_coint = adfuller(e_coint, maxlag=0)
if adf_e_coint[1] < 0.05:
    print("残差序列平稳。")
else:
    print("残差序列不平稳。")
# 残差自相关检验
e_coint_bg = acorr_breusch_godfrey(OLS_coint, nlags=1)
if e_coint_bg[1] < 0.05:
    print("残差序列存在一阶自相关。")
else:
    print("残差序列不存在一阶自相关。")