收入时间序列——之模型探索篇

前文《收入时间序列——之数学理解篇》已经梳理了时序分析所具备的基本数学原理，现在开始着手探索收入数据的内在规律，主要提出以下几个问题并给予解答。

1. 收入时间序列是平稳的吗？(偏)自相关情况如何？
2. 收入时间序列有趋势或季节性效应吗？
3. 用什么模型拟合效果比较好？
4. 有无异方差特征？用GARCH模型效果如何？
5. 收入增长率的平稳性情况如何？

下面开始具体探讨过程。

一. 线性模型时序分析步骤

当我们拿到一个时间序列（多半是非平稳的），一般作如下分析流程：

这里稍微留意一下，平稳性检验未通过后，差分运算只是其中一种方法，做完差分再回去通过平稳性检验后，一般会对原序列拟合  模型（图中是对差分后的序列拟合  模型），这里的  就是差分阶数，注意不是差分步数，差分步数是消除季节效应的。

二. 准备程序

为方便调用，将经常调用的绘图和假设检验程序事先定义好：

时序图

def draw_series(timeSeries, title):
    fig = plt.figure(figsize=(12,8))
    fig.add_subplot()
    plt.plot(timeSeries)
    plt.title(title)

ADF检验

# Dickey-Fuller检测:
def ADF_test(timeseries):
    dftest = adfuller(timeseries, autolag='AIC')
    print(dftest)
    # 尝试三种类型
    dftest_ct = adfuller(timeseries, regression='ct', autolag='AIC')
    dftest_c = adfuller(timeseries, regression='c', autolag='AIC')
    dftest_nc = adfuller(timeseries, regression='nc', autolag='AIC')
    print(dftest_ct,dftest_c,dftest_nc)

LB检验

# 白噪声检验
def LB_test(timeseries):
    [[lb], [p]] = acorr_ljungbox(timeseries, lags=1)
    if p < 0.05:
        print(u"原始序列为非白噪声序列，对应的p值为：%s" % p)
    else:
        print(u"原始序列为白噪声序列，对应的p值为：%s" % p)

# 带延迟参数
def LB_test_detail(timeseries, lag):
    acf, q, p = sm.tsa.acf(timeseries, nlags=lag, qstat=True)  ## 计算自相关系数 及p-value
    out = np.c_[range(1, lag+1), acf[1:], q, p]
    output = pd.DataFrame(out, columns=['lag', "AC", "Q", "P-value"])
    output = output.set_index('lag')
    print(output)

ACF/PACF图

# 自相关和偏相关图——全部序列
def draw_acf_pacf(timeseries):
    fig = plt.figure(figsize=(12,8))
    ax1 = fig.add_subplot(211)
    plot_acf(timeseries, ax=ax1)
    ax2 = fig.add_subplot(212)
    plot_pacf(timeseries, ax=ax2)

# 自相关和偏相关图——延迟参数
def draw_acf_pacf_lag(timeseries, lags):
    fig = plt.figure(figsize=(12,8))
    ax1 = fig.add_subplot(211)
    plot_acf(timeseries, lags=lags, ax=ax1)
    ax2 = fig.add_subplot(212)
    plot_pacf(timeseries, lags=lags, ax=ax2)

三. 数据处理

# 导入收入表结构
f = open('某店每日收入指标跟踪表.csv', encoding='UTF-8')
# 去掉EXCEL表头
data = pd.read_csv(f, header=[0,1,2])

# 特征选取前三列并重新命名
data = data.iloc[:,0:3]
data.columns = ['日期属性','账期','本日收入']

# 按账期正时序排序
data['账期'] = pd.to_datetime(data['账期'],infer_datetime_format=True)
data = data.sort_values(by=['账期'])
data = data.reset_index(drop=True)  #设定一个新的索引列

# 异常数据处理
# 用对应日期属性的收入中位数替换收入小于等于0的样本
for index in data.ix[data['本日收入'] <= 0].index:
    for name, group in data['本日收入'].groupby(data['日期属性']):
        if data.ix[index, 0] == name:
            data.ix[index, 2] = group.median()

# 检验是否还有异常数据
print("非正值情况：\n", (data['本日收入'] <= 0).sum(axis=0).sort_values(ascending=False))

# 重新将账期设为索引
data.set_index('账期', inplace=True)
income_data  = data.ix[:,'本日收入']

读入数据，并只保留时间和收入，日期属性是为计算中位数替代异常值，输出的时间序列数据定义为 income_data，是一个 Series 类型。

四. 收入平稳判断

# 收入时序图
draw_series(income_data, '本日收入')

显示上下幅度变化较大，貌似还有一点上升趋势，直观上不易观察出是否平稳：

选取前100天时序，看上去很平稳：

draw_series(income_data[:100], '本日收入')

取全部时序的ACF和PACF图，显示ACF虽然收敛至0但收敛很缓慢；PACF基本在零值附近：

# 自相关图/偏自相关图
draw_acf_pacf(income_data)

再取延迟40阶，很明显ACF除1阶延迟有正相关外，还在6/7/8、13/14/15、20/21/22等延迟阶处呈周期性显著相关；PACF在延迟8阶后基本收敛至2倍标准差以内，但也出现在15阶和22阶显著大于2倍标准差的情况。总体来说，自相关具有非常明显的周期特点：

draw_acf_pacf_lag(income_data, 40)

上图具体分析就是：本日收入和昨天、上周/上上周等的环比当天及前后各一天有显著正相关性。举个例子：今天是周三，和昨天周二、上周二/三/四、上上周二/三/四、以及上上上周二/三/四等收入正相关，可以粗略理解为它涨我也涨，它跌我也跌，有一种线性关系，正相关系数在0.5-0.7左右。随着延迟阶数增大，周环比前后那两天的正相关性大约在一个月以后不显著，唯独周环比的正相关性大约在延迟250天后才不显著。注意这家门店是个老店，说明周环比的正相关性会在8-9个月时间范围内一直显著存在，有长期记忆。我们分别取100阶和365阶延迟的ACF图再看下：

draw_acf_pacf_lag(income_data, 100)