404: Not Found
Link BandarQQ
Situs Idpro
Pkv QQ Online
Situs QQ Online
PKV GAMES QQ ONLINE
Judi DominoQQ
Pkv Games Online
Situs Judi Pkv Games
Hasil Pertandingan Sepakbola
Keluaran Togel Hari Ini
Indoxxi Subtitle Indonesia
Platform QQ Online LK21
Film Semi
BandarQQ DominoQQ Terbaik
Link Situs BandarQQ
Situs QQ Online
Idpro BandarQQ Terbaik
QQ ONLINE TERPERCAYA
Agen Judi PKV Games
Judi PKV QQ Online
Agen Games PKV Online
Judi QQ Online Terpercaya
Situs Resmi Judi BandarQQ
Judi QQ Online
Permainan BandarQ Terpercaya
Situs Poker Terpercaya
Situs BandarQQ Terpercaya
MataQQ
Pytorch实现广义回归神经网络GRNN – KAKO Academy of Sciences

Python / 机器学习 / 神经网络 · 2023年 8月 4日 2

Pytorch实现广义回归神经网络GRNN

【悲报】这篇文章的代码没有批量预测,看了基本原理后可以移步新文章:Tensor,广播机制与GRNN重构 – KAKO Academy of Sciences (kakosci.com) ,非常好更新,使我的张量升维。

General regression neural network,维基百科:General regression neural network – Wikipedia

本篇涵盖了广义回归神经网络的原理及其推导,其神经网络的构建,并根据其高度可并行化的特征,使用了Pytorch来构建算法,利用GPU来加速算法的运行。

课题项目的缘故开始接触这个算法,但似乎找不到哪个包里有能直接调用的实现?不过反正也不是很复杂,直接来动手实现好啦。先根据论文[1]整理一下算法的原理和我的推导思路:

假定当前存在一个已知的联合概率密度函数(probability density function,后文简称pdf):\(f\left (x,y \right)\),其中\(x\)代表一个随机向量,\(y\)则代表一个标量的随机变量,现在的问题是:在\(x\)为特定值\(X\)的条件下,\(y\)的期望是多少?即求条件期望\(E\left [y|X \right]\)。

先上一个条件概率公式:$$f_{y|X}(y|X)=\frac{f(X,y)}{f_{X}(X)}\cdots (1)$$
其中\(f_{y|X}\)是条件概率密度,\(f\)是最开始说的pdf,而\(f_{X}\)则代表向量\(x\)的概率密度。

现在开始展开条件期望:$$E[y|X]=\int_{-\infty}^{\infty}yf_{y|X}(y|X)\mathrm{d}y$$于是根据公式(1):$$=\int_{-\infty}^{\infty}\frac{yf(X,y)}{f_{X}(X)}\mathrm{d}y=\frac{\int_{-\infty}^{\infty}yf(X,y)\mathrm{d}y}{f_{X}(X)}$$最终得到:$$E[y|X]=\frac{\int_{-\infty}^{\infty}yf(X,y)\mathrm{d}y}{\int_{-\infty}^{\infty}f(X,y)\mathrm{d}y}\cdots (2)$$

现在我们把目光放在求\(f(X,y)\)上来,于是正式进入GRNN算法的核心内容之——高斯函数,并将样本引入到算法中,利用样本来估计概率密度。设\(x\)和\(y\)的已知样本为\(X^{i}\)和\(Y^{i}\),\(i=1,2,\dots,n\),\(n\)是样本数量,设\(p\)为向量\(x\)的维数,下面我将给出一串又臭又长)-:的公式:
$$\hat{f}(X,y)=\frac{1}{n}\sum_{i=1}^{n}\prod_{j=1}^{p}[\frac{1}{\sqrt{2\pi}\sigma}exp(-\frac{(X_j-X_j^i)^2}{2\sigma^2})]\cdot\frac{1}{\sqrt{2\pi}}exp(-\frac{(Y-Y^i)^2}{2\sigma^2})$$
$$=\frac{1}{(2\pi)^{\frac{(p+1)}{2}}\sigma^{(p+1)}}\frac{1}{n}\sum_{i=1}^{n}exp(-\frac{(X-X^i)^{T}(X-X^i)}{2\sigma^2})exp(-\frac{(Y-Y^i)^2}{2\sigma^2})\cdots (3)$$
仔细看就会明白,公式(3)实际上做的事情就是先构建\(p+1\)维的联合高斯分布,\(p\)维来自\(x\)而加上的一个维度就是\(y\)。接着根据\(n\)个样本生成了\(n\)个高斯分布,取对应样本的值作为每个高斯分布的均值\(\mu\),但方差\(\sigma\)统一,实际上这个方差将是GRNN唯一的超参数,也可以称为平滑参数。最后将这\(n\)个高斯分布加起来再除以\(n\),以保证生成的概率密度函数在整个空间上的积分为1。

上面这一步其实很容易理解,想象在经典的二维坐标系中,x轴上有许多点,以每一个点为中心构建方差一致的正态分布,并加在一起再除以点的个数,就得到了一个概率密度函数。最终在x轴上的点密集的地方概率密度就高,点稀疏的地方概率密度就低,而一开始统一的方差则决定了最终图像的形状,方差高则概率密度更均匀,起伏更缓,反之则更陡,故称平滑参数。于是根据得到的概率密度函数,我们就能够估计新来的点落在x轴某个位置的概率是多少。(从这里也能看出喂给GRNN的训练样本一定要覆盖的情况够全面,该算法在缺乏训练的空间上几乎没有预测能力:-)

得到概率密度\(f(X,y)\)的估计之后,就可以回到公式(2),开始计算条件期望了,此时公式(2)改成如下形式:$$\hat{Y}(X)=\frac{\int_{-\infty}^{\infty}y\hat{f}(X,y)\mathrm{d}y}{\int_{-\infty}^{\infty}\hat{f}(X,y)\mathrm{d}y}$$
带入公式(3),消去系数:$$=\frac{\int_{-\infty}^{\infty}y\sum_{i=1}^{n}exp(-\frac{(X-X^i)^{T}(X-X^i)}{2\sigma^2})exp(-\frac{(Y-Y^i)^2}{2\sigma^2})\mathrm{d}y}{\int_{-\infty}^{\infty}\sum_{i=1}^{n}exp(-\frac{(X-X^i)^{T}(X-X^i)}{2\sigma^2})exp(-\frac{(Y-Y^i)^2}{2\sigma^2})\mathrm{d}y}$$
求和的积分等于积分的求和:$$=\frac{\sum_{i=1}^{n}\int_{-\infty}^{\infty}y\cdot exp(-\frac{(X-X^i)^{T}(X-X^i)}{2\sigma^2})exp(-\frac{(Y-Y^i)^2}{2\sigma^2})\mathrm{d}y}{\sum_{i=1}^{n}\int_{-\infty}^{\infty}exp(-\frac{(X-X^i)^{T}(X-X^i)}{2\sigma^2})exp(-\frac{(Y-Y^i)^2}{2\sigma^2})\mathrm{d}y}$$
$$=\frac{\sum_{i=1}^{n}exp(-\frac{(X-X^i)^{T}(X-X^i)}{2\sigma^2})\int_{-\infty}^{\infty}y\cdot exp(-\frac{(Y-Y^i)^2}{2\sigma^2})\mathrm{d}y}{\sum_{i=1}^{n}exp(-\frac{(X-X^i)^{T}(X-X^i)}{2\sigma^2})\int_{-\infty}^{\infty}exp(-\frac{(Y-Y^i)^2}{2\sigma^2})\mathrm{d}y}$$
现在对公式分子分母后方的积分式子进行化简,具体方法就是构造高斯函数:
$$\int_{-\infty}^{\infty}exp(-\frac{(Y-Y^i)^2}{2\sigma^2})\mathrm{d}y=\sqrt{2\pi}\sigma\int_{-\infty}^{\infty}\frac{1}{\sqrt{2\pi}\sigma}exp(-\frac{(Y-Y^i)^2}{2\sigma^2})\mathrm{d}y=\sqrt{2\pi}\sigma$$
$$\int_{-\infty}^{\infty}y\cdot exp(-\frac{(Y-Y^i)^2}{2\sigma^2})\mathrm{d}y=\sqrt{2\pi}\sigma\int_{-\infty}^{\infty}y\cdot N(Y^i,\sigma^2)\mathrm{d}y=\sqrt{2\pi}\sigma\cdot Y^i$$
计算即将完成,这里我们先设\(D_i^2=(X-X^i)^{T}(X-X^i)\),即测试样本\(X\)与每个训练样本之间欧氏距离的平方。
最终得到了GRNN的公式:$$\hat{Y}(X)=\frac{\sum_{i=1}^{n}Y^{i}exp(-\frac{D_i^2}{2\sigma^2})}{\sum_{i=1}^{n}exp(-\frac{D_i^2}{2\sigma^2})}\cdots (4)$$

所以说GRNN实际上是一个只有前向传播的统计学习方法……而因为没有需要被迭代优化的参数,它的训练过程只有接收训练样本这一步骤而已,反而在预测时测试样本需要跟每一个训练样本发生作用,耗时集中在预测这一步骤。然后你可能会注意到,似乎上面的计算中\(y\)一直都是一个标量,那如果遇到多个输出维度的情况该怎么办?当然是每个输出维度单独作为一次计算,接下来的神经网络构建部分也是综合了多维度输出的情况(具体表现在求和层)。

算法的原理和推导以上就是全部了,接下来需要将公式(4)转化为神经网络的形式,不过这里就不画神经网络的示意图了,直接用矩阵的形式表达。将构建一个具有高斯层、求和层和输出层的神经网络,其中求和层与输出层将在具体实现中合并为一个层。

高斯层(即\(exp(-\frac{D_i^2}{2\sigma^2})\)):
设向量\(x\)为测试输入样本,矩阵\(t\)为训练输入样本集合,每一个输入样本的维度均为\(p\),一共有\(n\)个训练样本,Gauss即公式\(exp(-\frac{D_i^2}{2\sigma^2})\),即令测试输入样本与每个训练输入样本经过一次计算,得到行向量\(g\),流程图:

其pytorch类代码实现如下:

class GaussLayer(nn.Module):
    def __init__(self, training_inputs, sigma):
        super(GaussLayer, self).__init__()
        self.training_inputs = training_inputs
        self.sigma = sigma  # smoothing parameter

    def forward(self, x):
        out = x - self.training_inputs
        out = (out ** 2).sum(axis=1)
        out = - out / (2 * self.sigma ** 2)
        out = torch.exp(out)
        return out

求和与输出层:
向量\(g\)自身求和得到\(s_0\),对应公式(4)的分母:

训练输出样本\(y\)同样写成矩阵的形式,其中\(k\)代表单个输出样本的维度。行向量\(g\)与矩阵\(y\)的每一行经过求和层\(Sum\)进行一次内积,得到向量\(s\),向量\(s\)再除以前面求得的\(s_0\),最终得到预测输出向量\(o\),流程图:

其pytorch类代码实现如下:

class SumAndOutputLayer(nn.Module):
    def __init__(self, training_outputs):
        super(SumAndOutputLayer, self).__init__()
        self.training_outputs = training_outputs

    def forward(self, x):
        trans = self.training_outputs.T
        s0 = x.sum()
        out = (x * trans).sum(axis=1)  # Summation Layer
        out = out / s0  # Output Layer
        return out

可以看到,矩阵化后的流程中,绝大部分步骤都符合Pytorch的广播机制,因此非常适合通过自定义层来实现其整个运算过程、再将模型和数据放在GPU上运行,可以极大地改善GRNN在预测时的耗时问题。

下面将给出完整的代码,其中GRNN类中包含了数据归一化和还原的步骤。

需要注意的是,pytorch实现中使用的变量类型与上文所使用的向量、矩阵的结构有所差异。文中不论是列向量还是行向量,在代码中均表现为一维张量;矩阵用二维张量表示,文中每列代表一个样本,但代码中每行代表一个样本;而注意超参\(\sigma\)在调用类之前也需要先转化为零维张量。再者,当训练输入样本或训练输出样本(特别是输出常有这种情况)为一维,所有样本放在一个一维张量下时,必须先使用a = a.unsqueeze(axis=1)进行升维。

import torch
from torch import nn


class GaussLayer(nn.Module):
    def __init__(self, training_inputs, sigma):
        super(GaussLayer, self).__init__()
        self.training_inputs = training_inputs
        self.sigma = sigma  # smoothing parameter

    def forward(self, x):
        out = x - self.training_inputs
        out = (out ** 2).sum(axis=1)
        out = - out / (2 * self.sigma ** 2)
        out = torch.exp(out)
        return out


class SumAndOutputLayer(nn.Module):
    def __init__(self, training_outputs):
        super(SumAndOutputLayer, self).__init__()
        self.training_outputs = training_outputs

    def forward(self, x):
        trans = self.training_outputs.T
        s0 = x.sum()
        out = (x * trans).sum(axis=1)  # Summation Layer
        out = out / s0  # Output Layer
        return out


class CudaGeneralRegNN:
    def __init__(self):
        self.t = None  # training samples (INPUT)
        self.y = None  # training samples (OUTPUT)
        self.sigma = None  # smoothing parameter
        self.t_mean = None  # mean of each feature (INPUT)
        self.t_std = None  # std of each feature (INPUT)
        self.y_mean = None  # mean of each feature (OUTPUT)
        self.y_std = None  # std of each feature (OUTPUT)
        self.net = None  # General Regression Neural Network
        self.device = None  # cpu/cuda

    def fit(self, t_samples, y_samples, sigma):
        self.device = torch.device('cuda') if torch.cuda.is_available() else torch.device('cpu')
        self.t = t_samples
        self.y = y_samples
        self.sigma = sigma
        self.t_mean = self.t.mean(axis=0)
        self.t_std = self.t.std(axis=0)
        self.y_mean = self.y.mean(axis=0)
        self.y_std = self.y.std(axis=0)

        # Normalization
        self.t = (self.t - self.t_mean) / self.t_std
        self.y = (self.y - self.y_mean) / self.y_std

        self.t = self.t.to(device=self.device)
        self.y = self.y.to(device=self.device)
        self.sigma = self.sigma.to(device=self.device)

        self.net = nn.Sequential(
            GaussLayer(self.t, self.sigma),
            SumAndOutputLayer(self.y)
        )
        self.net = self.net.to(device=self.device)

    def predict(self, x):
        x = (x - self.t_mean) / self.t_std
        x = x.to(device=self.device)

        out = self.net(x)

        out = out.to(device=torch.device('cpu'))
        out = out * self.y_std + self.y_mean

        return out

耶整理完了,其实论文[1]只给出了一些比较粗略的公式,很多部分是我自己推导的,疲れた…

参考文献:
[1]D. F. Specht, “A general regression neural network,” in IEEE Transactions on Neural Networks, vol. 2, no. 6, pp. 568-576, Nov. 1991, doi: 10.1109/72.97934.