ESN笔记

组成

输入层，储备层，输出层

和RNN的联系和区别

相同点

隐藏状态及输出表达式相同。ESN是RNN的一种，前面时刻输入的信息会通过连接权重回荡在储备池中（所以叫回声状态网络），相当于RNN中该时刻的隐藏状态会再次循环输入，生成下一时刻的隐藏状态。

不同点

1、ESN中隐藏层的权重是随机初始化的，并且在训练过程中不会更新。因此，主要的训练任务是学习输出层的权重。
优点：相对于RNN（需要更新W,U,V），整个网络只需要训练输出权重，所以它的训练过程非常快
缺点：对于高维的时序数列（多个变量（或特征）随时间变化的数据序列），比如说视频帧处理，ESN不太能胜任。
2、ESN用线性回归方式更新输出矩阵
RNN用反向传播方式更新参数

构造过程

初始化

确定储备池大小

确定储备池的大小，即神经元的个数。节点数越多，拟合能力越强。由于ESN仅仅通过调整输出权值(我们最终就是利用储备池的状态信息来确定这个W)来线性拟合输出结果，所以一般ESN需要远大于一般神经网络的节点规模。

（1）包含数目较多的神经元（与经典神经网络相比）；
（2）神经元之间的连接关系随机产生；
（3）神经元之间的连接具有稀疏性。（？）

生成输入权重矩阵

型号为（神经元个数，输入特征维度数+1）

每一行表示一个神经元，每一列表示一个输入维度并且最后一列是表示偏置值

偏置值允许拟合的函数进行平移操作，能拟合复杂的数据分布，灵活地适应不同数据模式，降低欠拟合的风险

随机生成连接矩阵

这个矩阵表示了神经网络随机的内部初始结构，即哪些神经元之间是有连接的，以及连接的方向和权值，实际上就是有向图（神经网络图）的矩阵表示。

缩放矩阵

将矩阵缩放为谱半径<1

方法

1、用一个缩放因子乘初始矩阵

构建一个与初始矩阵相同维度的对角矩阵，对角线上的元素为缩放因子

2、使用特征值缩放

标准化（将特征值缩放为均值为0、方差为1的分布）：缩放矩阵的对角线元素为个征值的标准差的倒数

目的

1、稳定性：矩阵的谱半径被用来判断系统和算法的稳定性。如果系统的所有特征值的绝对值都小于1，则系统是稳定的；如果存在特征值的绝对值大于1，则系统是不稳定的。
2、收敛速度：矩阵的谱半径也可以用来估计迭代算法的收敛速度。通常情况下，矩阵的谱半径越小，迭代算法的收敛速度越快。
3、受激活函数的影响，比如sigmoid和tanh激活函数，在0与1之间区分度比较大，但是大于1之后，激活值变化不大。故将连接权重初始化为0-1之间(或-1至1)
4、如果权重矩阵的特征值大于1，那么通过多次迭代计算，激活值可能会呈指数级增长，导致激活值爆炸。

初始化结果

这些参数会影响到网络短期记忆时间的长短。
输入权值越小而内部矩阵的谱半径越接近1，网络短期记忆时间越长。

输入权重决定了输入信号对神经元活动的影响程度，而输入权值较小意味着输入信号传递给神经元的强度较弱，即输入信号对网络状态的影响会相对较小
内部矩阵的谱半径接近1时，意味着该矩阵的特征值分布在单位圆附近，即网络处于临界状态，表明网络的动态性较高。这意味着网络内部的状态不会偏向于过度衰减（消失）或过度增长（爆炸），而是能够保持在一种平衡的状态，使得信息在网络内部能够保持一定的传播和记忆。

但是，增强记忆能力的同时，这种操作也造成了网络对“快速变化”系统的建模能力下降。

网络内部的动态性主要受到内部矩阵的影响，而输入权值较小可能会导致网络对输入的响应较为迟缓。

训练

初始化储备池

方法

空转：将输入数据送入储备池中，但不会利用储备池的输出，储备池内部的状态会逐渐达到稳定状态。

目的

储备池的内部连接是随机的，最开始的输入序列得到储备池状态的噪声会比较大，所以会先使用一些数据来初始化储备池的状态，从而降低噪声的影响。确保网络的内部状态能够适应输入数据的动态特性，并且不受随机初始状态的影响。

使用线性回归确定输出权重

计算输出权重

使用线性回归方法（建立自变量（预测变量）与因变量（目标变量）之间的线性关系模型），找到最优的参数 （Wout），使得模型预测值与真实观测值 之间的残差（以L2范数表示）最小化
通过最小二乘法，计算$min||WX-Y|{|}_2^2$（还没完事）

优化：正则化方法防止过拟合

L1L2正则化详解

正则化

正则化指all减少过拟合的方法。

范数

将权重W看作高维向量空间的一个点，其到源点的距离看作欧氏距离（高维勾股定理），则记为L2范数；若为曼哈顿距离（坐标绝对值相加），则记为L1范数。

过拟合的原因

对于相同的损失函数值，可能会得到不同的输出权重值（有大有小），如果权重值过大，在变换过程中会放大测试数据中的噪声误差，造成较大的判断误差。

正则化方法

此处使用正则化方法本质就是控制参数（输出权重$W_{out}$的范数）的大小，就是给一个可行域范围（$W_{out}$范数大小），在范围内求最值（只控制$W_{out}$不管b是因为b只是控制了图像的平移位置）

可行域范围为$||W_{out}|{|}_2^2<C$
Lagrange乘数法求最值$min||WX-Y|{|}_2^2-\lambda (||W_{out}|{|}_2^2-C)$

将Lagrange函数转换为L2正则化形式的方式：
理解1：对$W_{out}$求导可消去常数C，故转换为L2正则化形式$min||WX-Y|{|}_2^2-\lambda ||W_{out}|{|}_2^2$
理解2：$\lambda$用于调整可行域与目标函数相切处的梯度大小，使得梯度相加为0。C用于调节可行域的大小，即切点的选取。C和$\lambda$的大小绑定，有C的式子，是把C作为超参数，$\lambda$随C改变；无C的式子，是把$\lambda$作为超参数。

L1正则化和L2正则化的区别

L1正则化的切点会在坐标轴上，即$W_{out}$的某些维度被置为0，即去除了某些特征的影响，实现了特征之间的去耦合化。
L2正则化就是把$W_{out}$的绝对值缩小。

正则化前后最值的的误差（？）

求解：岭回归（带L2范数惩罚项的最小二乘回归）

岭回归和L2正则化干的是同一件事

岭回归适用场景

在使用最小二乘法计算矩阵乘法结果$\beta =(X^{T}X)Xy$时，系数矩阵$X^{T}X$在大多数情况下不满秩，即参数之间存在多重共线性，特征之间有高度相关性，使用最小二乘法估计会出现回归模型中的参数估计不准确的问题（因为回归模型假设自变量之间不存在完全相关）

形式

在损失函数后添加L2范数的惩罚项（形式同L2正则化）
使用最小二乘估计时，结果表示为$\beta =(X^{T}X+\lambda I_n)Xy$

实现原理

添加罚项用来约束模型参数的大小，从而减少模型对自变量之间高度相关性的敏感度，从而缓解多重共线性。

L1正则化对应的最小二乘方式：Lasso回归

套索回归在目标函数中引入了一个额外的L1惩罚项，它通过对参数进行稀疏化（将某些参数缩小至零）来减少模型的复杂性。套索回归的一个重要特点是可以自动进行特征选择，即将某些不相关的自变量的系数收缩为零，从而剔除了一些可能导致多重共线性问题的自变量。

代码

from numpy import *
from matplotlib.pyplot import *
import scipy.linalg
import matplotlib.pyplot as plt
import numpy as np
import math

# 加载数据
# 前1000个数据用来训练，1001-2000的数据用来测试。训练数据中，前100项用来初始化储备池，以让储备池中形成良好的回声之后再开始训练。
trainLen = 1000
testLen = 1000
initLen = 100  # 前100项用来初始化储备池

# 生成ESN储层
inSize = outSize = 1  # inSize 输入维数 K
resSize = 1000  # 储备池规模 N
a = 0.3  # 可以看作储备池更新的速度，可不加，即设为1.

random.seed(42)
# 随机初始化 Win 和 W，且使随机数的范围置为 -0.5 到 0.5 之间（为了避免权重值过大或过小，从而减少训练过程中的梯度爆炸或梯度消失问题。）
Win = (random.rand(resSize, 1 + inSize) - 0.5)  # 输入矩阵 N * 1+K
W = random.rand(resSize, resSize) - 0.5  # 储备池连接矩阵 N * N

# 对W进行缩放，以满足稀疏的要求。
# 方案 1 - 直接缩放 (快且有脏数据, 特定储层): W *= 0.135
# 方案 2 - 归一化并设置谱半径 (正确, 慢):
print('计算谱半径...')
rhoW = max(abs(linalg.eig(W)[0]))  # linalg.eig(W)[0]:特征值   linalg.eig(W)[1]:特征向量
W *= 1.25 / rhoW    # 乘一个缩放因子1.25 / rhoW，权重矩阵的谱范数（即最大特征值的绝对值）归一化为 1.25（？）

# 为设计（收集状态）矩阵分配内存
X = zeros((1+inSize+resSize,trainLen-initLen))  # 储备池的状态矩阵x(t)：每一列是每个时刻的储备池状态。后面会转置
#状态矩阵每一行的意义 1：偏置值；insize:每一时刻输入值；ressize：每一时刻各神经元的状态
# 直接设置相应的目标矩阵
data_array = data.split('\n')  # 将读入数据（字符串类型）按回车分割为数组
Yt = data_array[initLen + 1:trainLen + 1]# 输入矩阵:每一行是一个时刻的输入
print(Yt)
# 输入所有的训练数据，然后得到每一时刻的输入值和储备池状态。
x = zeros((resSize, 1))
for t in range(trainLen):
    u = Yt[t]
    print(u)
    # 计算每个t的状态值
    # a为更新率，控制了新旧信息的相对重要性，较大的 a 值会更加强调新信息，而较小的 a 值则更加保留先前状态的信息，即控制状态更新速度
    x = (1 - a) * x + a * tanh(dot(Win, vstack((1, u))) + dot(W, x))  # vstack((1, u)):将偏置量1加入输入序列
    if t >= initLen:  # 空转100次后，开始记录储备池状态
        X[:, t - initLen] = vstack((1, u, x))[:, 0]  # 对初始化之后的时刻赋值该时刻的状态；vstack将每个t的偏置值、输入值、状态值堆叠成列赋值给X的每一列

# 使用Wout根据输入值和储备池状态去拟合目标值，这是一个简单的线性回归问题，这里使用的是岭回归(Ridge Regression)。
reg = 1e-8  # 正则化系数lambda
# 对X转置
X_T = X.T
# Wout:  1 * 1+K+N
Wout = dot(dot(Yt, X_T), linalg.inv(dot(X, X_T) + \
                                    reg * eye(
    1 + inSize + resSize)))  # linalg.inv矩阵求逆；numpy.eye()生成对角矩阵，规模:1+inSize+resSize，默认对角线全1，其余全0
# Wout = dot( Yt, linalg.pinv(X) )

# 使用训练数据进行前向处理得到结果
# run the trained ESN in a generative mode. no need to initialize here,
# because x is initialized with training data and we continue from there.
Y = zeros((outSize, testLen))
u = data[trainLen]
for t in range(testLen):
    x = (1 - a) * x + a * tanh(dot(Win, vstack((1, u))) + dot(W, x))
    y = dot(Wout, vstack((1, u, x)))  # 输出矩阵(1 * 1+K+N)*此刻状态矩阵(1+K+N * 1)=此刻预测值
    Y[:, t] = y  # t时刻的预测值   Y: 1 * testLen
    # 生成模型
    u = y
    # 预测模型
    # u = data[trainLen+t+1]

# 计算第一个errorLen时间步长的MSE
errorLen = 500
mse = sum(square(data[trainLen + 1:trainLen + errorLen + 1] - Y[0, 0: errorLen])) / errorLen
print('MSE = {0}'.format(str(mse)))

# 绘制测试集的真实数据和预测数据
figure(1).clear()
plot(data[trainLen + 1:trainLen + testLen + 1], 'g')
plot(Y.T, 'b')
title('Target and generated signals $y(n)$ starting at $n=0$')
legend(['Target signal', 'Free-running predicted signal'])

# 绘制储备池中前200个时刻状态(x(t))的前20个储层结点值
figure(2).clear()
plot(X[0:20, 0:200].T)
title('Some reservoir activations $\mathbf{x}(n)$')

# 绘制在各神经元上的输出权重值
figure(3).clear()
# bar(np.arange(1 + inSize + resSize), Wout.T, 8)
plot(np.arange(1 + inSize + resSize), Wout.T)

title('Output weights $\mathbf{W}^{out}$')
show()