回归与时间序列分析

Posted 2022-10-28 彭祥.

这篇博客是记录在数据挖掘中的回归与时间序列分析的学习过程

基础概念

回归分析应用

回归分析是一种应用极为广泛的数量分析方法。它用于分析事物之间的统计关系，侧重考察变量之间的数量变化规律，并通过回归方程的形式描述和反映这种关系，以帮助人们准确把握变量受其他一个或多个变量影响的程度，进而预测提供科学依据。
在大数据分析中，回归分析是一种预测性的建模技术，它研究的是因变量(目标)和自变量(预测器)之间的关系。这种技术通常用于预测分析、时间序列模型，以及发现变量之间的因果关系。

相关关系与函数关系

回归分析是处理多变量间相关关系的一种数学方法。相关关系不同于函数关系，后者反映变量之间的严格依存性，而前者则表现出一定程度的波动性或随机性，对自变量的每一个取值，因变量可以有多个数值与之相对应。
当自变量为非随机变量而因变量为随机变量时，它们的关系分析成为回归分析；当两者都是随机变量时，它们的关系分析成为相关分析。

回归分析的分类

回归分析包含线性回归与非线性回归，其中线性回归又包含简单线性回归和多元线性回归。下面是简单介绍：

简单线性回归分析

线性回归分析中，如果仅有一个自变量于一个因变量，且关系大致上可用一条直线表示，则称之为简单线性回归分析。
如果发现因变量Y和自变量X之间存在高度的正相关，则可以确定一条直线方程，使得所有的数据点尽可能接近这条拟合的直线。简单线性回归分析的模型可以用以下方程表示：

Y = a + b x

其中，Y为因变量，a为截距，b为相关系数，x为自变量。

多元线性回归分析

多元线性回归分析是简单线性回归分析的推广，指的是多个因变量对多个自变量的回归分析。其中最常用的是只限于一个因变量但有多个自变量的情况，也叫多重回归分析。多重回归分析的一般形式如下：

非线性回归数据分析

对于线性回归问题，样本点落在空间中的一条直线上或该直线的附近，因此可以使用一个线性函数表示自变量和因变量间的对应关系。然而在一些应用中，变量间的关系呈现曲线形式，因此无法用线性函数表示自变量和因变量间的对应关系，而需要使用非线性函数表示。
数据挖掘中常用的一些非线性回归模型：

二次曲线模型： Y=a+b₁​ X+b ₂ X²
双曲线模型： Y=a+ X/b

由于许多非线性模型是等价的，所以模型的参数化不是唯一的，这使得非线性模型的拟合和解释相比线性模型复杂得多。在非线性回归分析中估算回归参数的最通用的方法依然是最小二乘法（因为在线性回归分析中估算回归参数使用的也是最小二乘法）。

最小二乘法

最小二乘法是我们进行回归分析中估算参数的最常用的方法：求解多元线性回归残差平方和RSS最小化的参数向量?这种通过最小化真实值和预测值之间的RSS来求解参数的方法叫做最小二乘法。
当然这其实是一种损失函数，或者是代价函数，我们在后期也可以采用其他的评估指标，如均方误差MSE(mean squared error)来衡量我们的预测值和真实值的差异:

我们首先以简单的一元线性回归为例：

一元线性回归求解过程

多元线性回归分析求解过程

多元线性回归分析求解则可以借助线性代数来解决

注意：以上是在假设X^TX存在可逆阵的情况下，因为是线性多元，所以必定满足。
在这里,逆矩阵存在的充分必要条件是特征矩阵不存在多重共线性。

非线性回归求解过程

对于非线性的回归求解有两种方式，一种还是我们利用矩阵向量推导的公式进行回归分析，另一种则是通过迭代法，我们常使用的迭代方法即为梯度下降。
线性最小二乘的解是closed-form即有可逆阵 ，而非线性最小二乘没有closed-form（即 没有可逆矩阵），通常用迭代法求解。

向量矩阵形式求参数

在X^TX存在可逆阵情况下可以使用该方法：

同时我们可以使用程序来实现拟合过程：

import numpy as np
import matplotlib.pyplot as plt
import numpy.linalg as lg
t = np.arange(1, 17, 1)
y = np.array([4, 6.4, 8, 8.8, 9.22, 9.5, 9.7, 9.86, 10, 10.20, 10.32, 10.42, 10.5, 10.55, 10.58, 10.6])
plt.figure()
plt.plot(t, y, 'k*')
# y=at^2+bt+c
A = np.c_[t ** 2, t, np.ones(t.shape)]
w = lg.inv(A.T.dot(A)).dot(A.T).dot(y)
plt.plot(t, w[0] * t ** 2 + w[1] * t + w[2])
plt.show()

迭代法求解回归分析问题

迭代法，即在每一步update未知量逐渐逼近解，可以用于各种各样的问题（包括最小二乘），比如求的不是误差的最小平方和而是最小立方和。
梯度下降是迭代法的一种，可以用于求解最小二乘问题（线性和非线性都可以）。
在python的工具库中，已经对这些方法进行封装：

# 使用带有Scikit-Learn的随机梯度下降执行线性回归
from sklearn.linear_model import SGDRegressor
# 运行1000个轮次，或者直到一个轮次期间损失下降小于1e-3，默认学习调度为0.1，不使用任何正则化
sgd_reg = SGDRegressor(max_iter = 1000, tol = 1e-3, penalty = None, eta0 = 0.1)
sgd_reg.fit(X, Y.ravel()) # ravel()对多维数据进行扁平化操作,返回的是一个数组的非副本视图
sgd_reg.intercept_, sgd_reg.coef_

回归分析解决的问题

回归分析是研究变量间函数关系的一种方法。变量之间的关系可以表示为方程的形式。
其能够进行因素分析，确定影响某一变量的若干变量（因素）中，何者为主要，何者为次要，以及它们之间的关系，比如研究房屋面积与价格，或者考虑更多因素影响价格（多元）的情况。

逻辑回归

逻辑回归也被称为广义线性回归模型，它与线性回归模型的形式基本上相同，最大的区别就在于它们的因变量不同，如果是连续的，就是多重线性回归；如果是二项分布，就是Logistic回归。
Logistic回归虽然名字里带"回归"，但它实际上是一种分类方法，主要用于二分类问题(即输出只有两种，分布代表两个类别)。逻辑回归就是这样的一个过程；面对一个回归或者分类问题，建立代价函数，然后通过优化方法迭代求解出最优的模型参数，然后测试检验这个求解的模型的好坏。它的优点有：速度快，适合二分类问题；简单、易于理解，可以直接看到各个特征的权重；能容易地更新模型吸收新的数据。它的缺点有：对数据和场景的适应能力有局限性，不如决策树算法适应性强。
逻辑回归的用途主要有以下3个方面。

寻找危险因素：寻找某一疾病的危险因素等；
预测：根据模型，预测在不同的自变量情况下，发生某种疾病或某种情况的概率有多大；
判别：实际上跟预测有些类似，也是根据模型，判断某人属于某些疾病或属于某种情况的概率有多大。
路径回归的常规步骤：寻找h函数(即预测函数)，构造J函数(损失函数)，想办法使得J函数最小并求得回归参数。

构造预测函数

二分类问题的概率与自变量之间的关系图形往往是一个S型曲线，采用sigmoid函数实现，函数形式为：

函数图像：

对于线性边界的情况，边界形式

最佳参数：

构造预测函数为：

sigmoid的函数输出是介于(0,1)之间的，中间值是0.5，公式的含义就很好理解了，因为输出是介于（0，1）之间，也就表明了数据属于某一类别的概率，例如，h θ ( x ) < 0.5 则说明当前数据属于A类；h θ ( x ) > 0.5则说明当前数据属于B类。所以可以将sigmoid函数看成样本数据的概率密度函数。
函数的值有特殊含义，它表示结果取1的概率。

时间序列分析

回归分析中如果自变量是时间，则将按时间顺序发生的离散型观测数据序列 X_t（t=1,2,3,4）称为时间序列，根据时间序列，揭示内在系统的统计特性和发展规律的统计方法，称为时间序列分析。

时间序列分析和回归分析的差别

时间序列分析方法明确强调变量值顺序的重要性，而其它统计分析方法则不必如此;
时间序列各观察值之间存在一定的依存关系， 而其它统计分析一般要求每一变量各 自独立;
时间序列分析根据序列自身的变化规律来预测未来，而其它统计分析则根据某一变量与其它变量间的因果关系来预测该变量的未来;
时间序列是一-组随机变 量的一次样本实现，而其它统计分析的样本值一般是 对同一随机变量进行N次独立重复实验的结果;
二者建模思路不同。

为何有如此差别？

回归分析对数据的假设：独立性

在回归分析中，我们假设数据是相互独立的。这种独立性体现在两个方面：一方面，自变量(X)是固定的，已被观测到的值，另一方面，每个因变量(y)的误差项是独立同分布，对于线性回归模型来说，误差项是独立同分布的正态分布，并且满足均值为0，方差恒定。

这种数据的独立性的具体表现就是：在回归分析中，数据顺序可以任意交换。在建模的时候，你可以随机选取数据循序进行模型训练，也可以随机选取一部分数据进行训练集和验证集的拆分。也正因为如此，在验证集中，每个预测值的误差都是相对恒定的：不会存在误差的积累，导致预测准确度越来越低。

时间序列对数据的假设：相关性

但对于时间序列分析而言，我们必须假设而且利用数据的相关性。核心的原因是我们没有其他任何的外部数据，只能利用现有的数据走向来预测未来。因此，我们需要假设每个数据点之间有相关性，并且通过建模找到对应的相关性，利用它去预测未来的数据走向。这也是为什么经典的时间序列分析(ARIMA)会用ACF（自相关系数）和PACF（偏自相关系数）来观察数据之间的相关性。

时间序列对相关性的假设直接违背了回归分析的独立性假设。在多段时间序列预测中，一方面，对于未来预测的自变量可能无法真实的观察到，另一方面，随着预测越来越远，误差会逐渐积累：你对于长远未来的预测应该会比近期预测更不确定。因此，时间序列分析需要采用一种完全不同的视角，用不同的模型去进行分析研究。

常用模型

这些都是统计学上使用的模型，而我们作为计算机系如果要实现时间序列预测的话，可以选用RNN,LSTM等神经网络模型来实现。

AR模型（自回归模型）和线性回归模型的“相似”&区别

时间序列分析中一个基础模型就是AR(Auto-Regressive)模型。它利用过去的数据点来预测未来。举例而言，AR(1)模型利用当前时刻的数据点预测未来的值，它们的数学关系可以被表示为：

它的表达形式的确和线性回归模型非常类似，甚至连一般的AR(n)模型都和线性回归有很高的相似性。唯一的差别就是等式右边的自变量(X)变成了过去的因变量(y)。而正是因为这一点微小的差异，导致两者的解完全不同。在AR模型中，由于模型自变量成为了过去的因变量，使得自变量与过去的误差之间有相关性。而这种相关性使得利用线性模型得到的AR模型（自回归模型）的解会是有偏估计(biased)。
对于上述结论的实际证明需要引入过多的概念。在此我们只对AR(1)模型作为一个特例来分析。不失一般性，我们可以通过平移数据将AR(1)模型表示成如下的形式：

对于这类模型，线性回归会给出以下的估计值：

对于一般的线性回归模型而言，由于所有的自变量都会被视为已经观测到的真实值。所以当我们取均值的时候，我们可以把分母当作已知，通过过去观测值和未来误差无关的性质得到无偏的结论。

　但是在时间序列下就无法得到无偏的兴致了，因为分子分母会互相干扰。因为自变量无法被视为已知，而且未来的观察值会与过去的误差项相互联系。因此，相关性使得利用线性模型得到的AR模型的解会是有偏估计(biased)。
　更直观的数据模拟也可以说明这个问题。如下图所示，左边是当参数真实值为0.9时通过数据模拟作出的平均值，可以看到真实值（黑线）和模拟值（红线）有一定的差距，但随着数据量的增大，差距在逐渐的缩小。右边是真实参数不同的时候，偏差的大小。可以看到，它们的误差一直存在，但是随着数据量的增加，误差逐渐变小。

事实上，我们会用线性回归模型去近似求解AR模型。因为虽然结果会是有偏的，但是却是一致估计。也就是说，当数据量足够大的时候，求解的值会收敛于真实值。

时间序列数据特征

长期趋势 (Trend)：时间序列在长期呈现的向上or向下的波动
季节变动 (Seasonal Effect)：随着期间变化（月，季，年）而变化
循环波动 (Cyclical Fluctuation)：沿着trend波动，不是单一方向的，而是涨落相同的交替波动。
不规则波动 (Residual/error)：除去趋势，季节波动，周期波动之外的随机波动。只含有随机波动的序列是平稳序列

时间序列性质

平稳性 (Stationarity)（转换为稳定数据）

1.概念：一个序列的统计特性不受观测时间的影响，即均值，方差，自相关系数等不因时间变化而变化。如果时间序列有季节性和趋势性，这个序列则为不平稳的。
差分：将数据转换为平稳性数据转换为稳定数据
大部分时间序列模型需要先转换为平稳模型，再进行建模。
差分是最常见的一种方法之一。其目的是：为了消除序列的水平的变化，即为消除季节性或趋势性，帮助稳定时间序列的均值。
我们在实现时间序列预测时，可以选择使用统计学方法，也可以使用深度学习的方法，这里我选用LSTM模型来进行时间序列预测。

# 所谓的差分，就是相邻的两个数相减，求得两个数之间的差，，一两个数之间的差作为一个数组，这样的数组体现了相邻两个数的变化情况
# 所谓的去趋势也就是一个说法而已，可以忽略
# 做差分，去趋势，获得差分序列
def difference(dataset, interval=1):
    diff = list()
    for i in range(interval, len(dataset)):
        value = dataset[i] - dataset[i - interval]  # 当前时间步t的值减去时间步t-interval的值
        diff.append(value)
    return Series(diff)#Series这个方法是把一个数组建立起一个一一对应的索引，参考：

转换为监督数据

def timeseries_to_supervised(data, lag=1):  # lag表示的是当前的值只与历史lag个时间步长的值有关，也就是用lag个数据预测下一个
    df = DataFrame(data)
    colums = [df.shift(i) for i in range(1, lag + 1)]  # 原始数据时间窗向后移动lag步长
    colums.append(df)  # 拼接数据
    df = concat(colums, axis=1)  # 横向拼接重塑数据，格式:input putput
    df.fillna(0, inplace=True)  # 由于数据整体向后滑动lag后，前面的lag个数据是Na形式，用0来填充
    return df

建立预测模型并训练

使用的是LSTM模型

fit LSTM来训练数据

def fit_lstm(train, batch_size, nb_epoch, neurons):
    X, y = train[:, 0:-1], train[:, -1]
    X = X.reshape(X.shape[0], 1, X.shape[1])
    model = Sequential()
    # 添加LSTM层
    model.add(LSTM(neurons, batch_input_shape=(batch_size, X.shape[1], X.shape[2]), stateful=True))
    model.add(Dense(1))  # 输出层1个node
    # 编译，损失函数mse+优化算法adam
    model.compile(loss='mean_squared_error', optimizer='adam')
    for i in range(nb_epoch):
        # 按照batch_size，一次读取batch_size个数据
        model.fit(X, y, epochs=1, batch_size=batch_size, verbose=0, shuffle=False)
        model.reset_states()
        print("当前计算次数："+str(i))
    return model