深度讲解PyG实现自定义layer的GCN
完整代码自定义layer传播方式从节点角度解读GCN原理逐行讲解代码原理initforwardmessage目前的代码讲解基本都是直接使用PyG内置的包实现固定结构的网络层。虽然我们可以通过每层使用不同的传递方式来建立不同的网络，但是却不能自定义网络层的传递方式，对于做创新性的研究工作而言是一个不足。
本篇在GCN学习：Pytorch-Geometric教程（二）的基础上，自定义了GCN的层传递方式（仍然是按照论文中的传递方式建立，但是我们以后也可以建立其他传递方式），其他代码与系列（二）的代码相同。
完整代码
import torchfrom torch_geometric.nn import MessagePassingfrom torch_geometric.utils import add_self_loops,degreefrom torch_geometric.datasets import Planetoidimport sslimport torch.nn.functional as Fclass GCNConv(MessagePassing):def __init__(self,in_channels,out_channels):super(GCNConv,self).__init__(aggr='add')self.lin=torch.nn.Linear(in_channels,out_channels)def forward(self,x,edge_index):edge_index, _ = add_self_loops(edge_index,num_nodes=x.size(0))x=self.lin(x)row,col=edge_index#计算度矩阵deg=degree(col,x.size(0),dtype=x.dtype)#度矩阵的-1/2次方deg_inv_sqrt=deg.pow(-0.5)norm=deg_inv_sqrt[row]*deg_inv_sqrt[col]return self.propagate(edge_index,x=x,norm=norm)def message(self,x_j,norm):return norm.view(-1,1)*x_jssl._create_default_https_context = ssl._create_unverified_contextdataset = Planetoid(root='Cora', name='Cora')print(dataset)print(dataset.num_node_features)print(dataset.num_classes)class Net(torch.nn.Module):def __init__(self):super(Net, self).__init__()self.conv1 = GCNConv(dataset.num_node_features, 16)self.conv2 = GCNConv(16, dataset.num_classes)def forward(self, data):x, edge_index = data.x, data.edge_indexx = self.conv1(x, edge_index)x = F.relu(x)x = F.dropout(x, training=self.training)x = self.conv2(x, edge_index)return F.log_softmax(x, dim=1)device = torch.device('cuda' if torch.cuda.is_available() else 'cpu')model = Net().to(device)data = dataset[0].to(device)optimizer = torch.optim.Adam(model.parameters(), lr=0.01, weight_decay=5e-4)model.train()for epoch in range(200):optimizer.zero_grad()out = model(data)loss = F.nll_loss(out[data.train_mask], data.y[data.train_mask])loss.backward()optimizer.step()model.eval()_, pred = model(data).max(dim=1)correct = int(pred[data.test_mask].eq(data.y[data.test_mask]).sum().item())acc = correct/int(data.test_mask.sum())print('Accuracy:{:.4f}'.format(acc))>>>Cora()14337Accuracy:0.8010Process finished with exit code 0
自定义layer传播方式
class GCNConv(MessagePassing):def __init__(self,in_channels,out_channels):#将邻居节点的特征加和super(GCNConv,self).__init__(aggr='add')self.lin=torch.nn.Linear(in_channels,out_channels)def forward(self,x,edge_index):#x的规模为[N,in_channels]#edge_index为[2,E]#A波浪实际就是在原图上将每一个节点都加入了自环edge_index, _ = add_self_loops(edge_index,num_nodes=x.size(0))#节点特征矩阵的线性变换x=self.lin(x)row,col=edge_index#计算节点的度deg=degree(col,x.size(0),dtype=x.dtype)#度的-1/2次方deg_inv_sqrt=deg.pow(-0.5)norm=deg_inv_sqrt[row]*deg_inv_sqrt[col]#进行层传播return self.propagate(edge_index,x=x,norm=norm)def message(self,x_j,norm):#x_j规模为[E,out_channels]#归一化节点特征return norm.view(-1,1)*x_j
从节点角度解读GCN原理
论文定义的网络传播公式如上。
H和W和特征矩阵和参数矩阵。我们需要构造的为A波浪和D波浪的-1/2次方。
但要理解这部分代码，我们需要牢记的是PyG中直接操作的对象是图，所以我们并不能直接以矩阵形式编写代码，而需要对图节点进行操作进行代码编写。
所以这里我们采用下面的公式，以节点为目标进行传递，图卷积层将每个节点表示为其相邻节点的聚合。下面的公式中的j属于N（i）U{i},N(i)为邻居节点，{i}即相当于添加了自环。
这一个公式转化其实困惑了我一整天，直到第二天看到这篇文章的解释才越来越理解：/p/89503068
其实从矩阵的角度来理解还不够直观，可以直接从节点的角度来理解，毕竟我们的直接操作对象就是节点：
那么我们的步骤可以拆解为以下5步：
1.增加自连接的邻接矩阵，即邻接矩阵的对角线元素为1，得到L=D−AL=D-AL=D−A;
2.对节点的特征矩阵进行线性变换，将特征变换到维度D，即X*W;
3对节点特征进行规范化， 即乘以归一化的拉普拉斯算子，即前面的归一化;
4.对邻居节点特征进行聚合操作，这里是求和，即矩阵形式公式中的A×X；
5.返回新的节点embedding；
这样我们再回头看代码，就能结合矩阵形式公式理解清楚代码的原理了！
逐行讲解代码原理
init
def __init__(self,in_channels,out_channels):super(GCNConv,self).__init__(aggr='add')self.lin=torch.nn.Linear(in_channels,out_channels)
init代码中的self.lin很好理解，就是定义了特征的线性变换函数。但上面的super函数，笔者理解它很费周折。主要费解的地方在于定义了是求和，即求和操作实在何时运行的。
forward
def forward(self,x,edge_index):print('edge_index.size',edge_index.size())print('x.size',x.size())print('x',x)edge_index, _ = add_self_loops(edge_index,num_nodes=x.size(0))print('edge_index.size:',edge_index.size())x=self.lin(x)print('x.size after lin:',x.size())row,col=edge_index#计算度矩阵deg=degree(col,x.size(0),dtype=x.dtype)print('deg.size:',deg.size())#度矩阵的-1/2次方deg_inv_sqrt=deg.pow(-0.5)print('deg_inv_sqrt:',deg_inv_sqrt)print('deg_inv_sqrt.size:',deg_inv_sqrt.size())print('deg_inv_sqrt[row]',deg_inv_sqrt[row])print('deg_inv_sqrt[row].size',deg_inv_sqrt[row].size())norm=deg_inv_sqrt[row]*deg_inv_sqrt[col]print('norm.size:',norm.size())return self.propagate(edge_index,x=x,norm=norm)>>>edge_index.size torch.Size([2, 10556])x.size torch.Size([2708, 1433])x tensor([[0., 0., 0., ..., 0., 0., 0.],[0., 0., 0., ..., 0., 0., 0.],[0., 0., 0., ..., 0., 0., 0.],...,[0., 0., 0., ..., 0., 0., 0.],[0., 0., 0., ..., 0., 0., 0.],[0., 0., 0., ..., 0., 0., 0.]])edge_index.size: torch.Size([2, 13264])x.size after lin: torch.Size([2708, 16])deg.size: torch.Size([2708])deg_inv_sqrt: tensor([0.5000, 0.5000, 0.4082, ..., 0.7071, 0.4472, 0.4472])deg_inv_sqrt.size: torch.Size([2708])deg_inv_sqrt[row] tensor([0.5000, 0.5000, 0.5000, ..., 0.7071, 0.4472, 0.4472])deg_inv_sqrt[row].size torch.Size([13264])norm.size: torch.Size([13264])
这里加的print代码都是笔者为了更好理解代码运行原理所添加的。相比于直接查看edge_index和x的值，查看其size其实会更直观展现函数效果。forward函数的每一步代码其实比较容易理解。只要记住我们是按照对节点操作的公式来实现的，即下面这个公式：
edge_index, _ = add_self_loop(edge_index,num_nodes=x,size(0))执行完成后，获得了添加了自环后的edge_index。这里代码需要按照上面的格式来写，即在edge_index后加上’, _’,否则会报错。可以看到执行完下面edge_index的规模发生了变化，增加了2708条边，正是顶点的数量。
x=self.lin(x)则很简单，就是把特征进行了线性映射。
deg=degree(col,x.size(0),dtype=x.dtype)就是计算添加了自环后的度矩阵。度矩阵的规模显然是与顶点数量相同的。deg_inv_sqrt=deg.pow(-0.5)计算度矩阵的-1/2次方。但下面的代码比较费解：norm=deg_inv_sqrt[row]*deg_inv_sqrt[col]其中norm的规模是13264，即所有边的数量。从形式来看，这一步对应于deg(i)的-1/2×deg(j)的-1/2。之前的运算已经得到了的矩阵的-1/2。那么deg_inv_sqrt[row]操作究竟是什么意思呢。从规模来看，deg_inv_sqrt为2708，但deg_inv_sqrt[row]已经变成了13264，注意这里的row并不是索引，而是edge_index，edge_index本身是一个两行矩阵：new_edge_index tensor([[ 0, 0, 0, ..., 2705, 2706, 2707], [ 633, 1862, 2582, ..., 2705, 2706, 2707]])第一行代表边起点节点，二行代表边终点节点。所以这一步的作用实际上相当于将deg_inv_sqrt按照edge_index中的两行进行扩充，按照edge_index中边的起点顺序进行排列，deg_inv_sqrt[row]和deg_inv_sqrt[col]中每项分别代表起点节点和终点节点对应值在deg_inv_sqrt中的索引。所以norm的计算即对应degi的-1/2*degj的-1/2。self.propagate(edge_index,x=x,norm=norm)对上面计算得到的edge_index,x,norm进行传播。
message
def message(self,x_j,norm):print('x_j:',x_j)print('x_j.size',x_j.size())print('norm',norm)print('norm.size',norm.size())print('norm.view.size',norm.view(-1,1).size())ww=norm.view(-1,1)*x_jprint('message_result',ww)print('message_result.size',ww.size())return ww>>>x_j: tensor([[ 0.0361, 0.0182, 0.0163, ..., 0.0196, -0.0226, -0.0110],[ 0.0361, 0.0182, 0.0163, ..., 0.0196, -0.0226, -0.0110],[ 0.0361, 0.0182, 0.0163, ..., 0.0196, -0.0226, -0.0110],...,[ 0.1105, 0.0317, -0.0821, ..., -0.0598, 0.0412, -0.0177],[ 0.0834, -0.0274, -0.0346, ..., -0.0309, -0.0100, 0.0086],[-0.0859, -0.0149, -0.0071, ..., -0.0422, -0.0136, 0.0341]],grad_fn=<IndexSelectBackward>)x_j.size torch.Size([13264, 16])norm tensor([0.2500, 0.2236, 0.2500, ..., 0.5000, 0.2000, 0.2000])norm.size torch.Size([13264])norm.view.size torch.Size([13264, 1])message_result tensor([[ 0.0090, 0.0045, 0.0041, ..., 0.0049, -0.0057, -0.0028],[ 0.0081, 0.0041, 0.0037, ..., 0.0044, -0.0051, -0.0025],[ 0.0090, 0.0045, 0.0041, ..., 0.0049, -0.0057, -0.0028],...,[ 0.0553, 0.0158, -0.0410, ..., -0.0299, 0.0206, -0.0088],[ 0.0167, -0.0055, -0.0069, ..., -0.0062, -0.0020, 0.0017],[-0.0172, -0.0030, -0.0014, ..., -0.0084, -0.0027, 0.0068]],grad_fn=<MulBackward0>)message_result.size torch.Size([13264, 16])
x_j实际就是计算公式中所有要计算的顶点，1个i可能对应多个邻居j节点，那么对于这个节点i，我们要计算1+j次。不过在message步中，我们只是将每次正则化计算即下图这部分的结果计算出来，所以共有13264个结果。
这时候我们再回头看__init__函数部分，我们就可以知道，sum是对message中同个节点i的所有邻居j及自环i的计算结果进行sum操作，得到了最终一个[2708,16]规模的x！
