图神经网络与分子表征:番外——基组选择

学过高斯软件的人都知道,我们在撰写输入文件 gjf 时需要准备输入【泛函】和【基组】这两个关键词。

【泛函】敲定计算方法,【基组】则类似格点积分中的密度,与计算精度密切相关。

部分研究人员借用高斯中的一系列基组去包装输入几何信息(距离、角度和二面角),这样做一方面提高了GNN的可解释性,另一方面也实实在在的提高了模型精度。从 AI 角度看,embedding则可以看作是几何信息的升维。

具体来说:

  1. 如果模型输入仅有距离信息,则采用径向基函数去embedding。常用的有 Gaussian ,也有Bessel
  2. 如果模型输入含有距离和角度信息。在直角坐标系下,可以用 Gaussian 和 sin 函数组embedding。在球坐标系下,可以考虑 spherical Bessel functions and spherical harmonics 组合。其中 spherical harmonics 采用m=0的形式。
  3. 如果模型输入含有距离,角度和二面角信息,一般采用 spherical Bessel functions and spherical harmonics 组合。可能有其他的,但目前涉及二面角的模型较少,据我了解,Spherenet和ComENet均采用的是这种组合。

下面进行简要介绍:

Gaussian 系列基组

SchNet网络架构中使用的基组,是目前用途最广的基组之一。
我们借助 DIG 框架中 schnet 的实现,对其进行可视化:

from dig.threedgraph.method.schnet.schnet import *

import numpy as np
import math
import matplotlib.pyplot as plt

import torch


dist_test = torch.arange(0.01, 5.01, 0.01)
dist_emb = emb(num_gaussians=5)
y = dist_emb(dist_test)
y = y.T

for idx, y_plot in enumerate(y):
    x = [a_dist.detach().numpy() for a_dist in dist_test]
    y = [an_emb.detach().numpy() for an_emb in y_plot]

    plt.plot(x, y, label=f"Gaussian embedding {idx}")

plt.xlabel("x")
plt.ylabel("y")
plt.legend()
plt.show()

结果如下图所示:
在这里插入图片描述
所谓,“对几何信息进行嵌入”,指,同一个距离信息对应x轴一个点。如果高斯基组有5,则,嵌入后,该距离信息就映射到了5个口袋里,获得一组长度为5的特征向量。

此处为了清晰的可视化,仅设置 num_gaussians=5 ,在实际应用中,这一数值往往设的很高。例如,原版的 schnet 将这一数值设为 300,在 DIG 版本中,这一数值是默认的 50,而在最新的 schnetpack 中,这一数值 降为了 20.

Bessel 系列基组

与高斯基组类似,Bessel 系列基组用于 embedding 距离信息,文献里用 spherical Bessel functions 表示。

其源头可以追溯到微分方程的求解,spherical Bessel functions 是作为一系列解中的径向部分存在,也常被称为 radical Bessel functions。

最早使用 Bessel functions 的(可能不严谨)GNN大概是 DimeNet。据 DimeNet 原文报道,使用 Bessel functions 会带来一定程度的精度提升。
我们借助 DIG 框架中 DimeNet 的实现,对其进行可视化:

from dig.threedgraph.method.spherenet.features import *

import numpy as np
import math
import matplotlib.pyplot as plt

import torch

dist_test = torch.arange(0.01, 5.01, 0.01)
dist_emb = dist_emb(num_radial=5)
y = dist_emb(dist_test)
y = y.T

for idx, y_plot in enumerate(y):
    x = [a_dist.detach().numpy() for a_dist in dist_test]
    y = [an_emb.detach().numpy() for an_emb in y_plot]

    plt.plot(x, y, label=f"radical_basis_{idx}")

plt.xlabel("x")
plt.ylabel("y")
plt.legend()
plt.show()

结果如下图所示:
在这里插入图片描述

spherical harmonics 基组

spherical Bessel functions 和 spherical harmonics 不是一个基组。他俩分别对应方程特解中的径向和角度部分。
(下图为 ComENet 中的概述)在这里插入图片描述
spherical harmonics 基组常常在球极坐标系下,和 spherical Bessel functions 配套使用。
如果输入的几何信息仅有角度,没有二面角,我们将 spherical harmonics 中的 m 置零。
此时得到的是一系列二维的 embedding 矩阵。
我们借助 DIG 框架中 SphereNet 的实现,对其进行可视化(源码稍微改了改,此处仅是一些思路):

from dig.threedgraph.method.spherenet.features import *

import numpy as np
import math
import matplotlib.pyplot as plt

import torch

angle_emb = angle_emb(num_spherical=4, num_radial=4, cutoff=4)
rlist = np.arange(0, 4.01, 0.005)  # Angstroms
thetalist = np.radians(np.arange(0, 361, 0.5))  # Radians
rmesh, thetamesh = np.meshgrid(rlist, thetalist)  # Generate a mesh

n = 1
l = 1
fig = plt.figure()
info = angle_emb(torch.tensor(rlist), torch.tensor(thetalist))
info_0 = info[n, l]
info_0 = info_0.detach().numpy()

info_0 = info_0.reshape(len(rlist), len(thetalist))
info_0 = info_0.T
fig, ax = plt.subplots(subplot_kw=dict(projection='polar'))
ax.contourf(thetamesh, rmesh, info_0, 100, cmap='RdBu')
ax.set_rticks([])
ax.set_xticks([])
plt.savefig(f'./basis/n_{n}_l_{l}.png', dpi=400)

结果如下图所示:
请添加图片描述
我们可以得到一系列能够embedding角度和距离信息的函数。
下图是DimeNet原文中的图:
在这里插入图片描述
需要注意的是,DimeNet源码中对 l=0 的径向函数进行了修改,所以无法复现 Figure 2 第一行。

我们还可以借助 scipy 进行实现,例如,下面我们对角度部分 ( spherical harmonics )进行可视化(不涉及径向部分,径向部分在 scipy.special._spherical_bessel 里):

  1. 借用plotly实现可交互的可视化
import plotly.graph_objects as go
import numpy as np
from mpl_toolkits.mplot3d import Axes3D
from scipy.special import sph_harm

# from scipy.special._spherical_bessel import


# l, m = 3, 0

for l in range(0, 4):
    for m in range(-l, l+1):
        theta = np.linspace(0, np.pi, 100)
        phi = np.linspace(0, 2 * np.pi, 100)
        theta, phi = np.meshgrid(theta, phi)
        xyz = np.array([np.sin(theta) * np.sin(phi),
                        np.sin(theta) * np.cos(phi),
                        np.cos(theta)])

        Y = sph_harm(abs(m), l, phi, theta)

        if m < 0:
            Y = np.sqrt(2) * (-1) ** m * Y.imag
        elif m > 0:
            Y = np.sqrt(2) * (-1) ** m * Y.real
        Yx, Yy, Yz = np.abs(Y) * xyz

        fig = go.Figure(data=[go.Surface(x=Yx, y=Yy, z=Yz, surfacecolor=Y.real), ])

        fig.update_layout(title=f'Y_l_{l}_m_{m}', )

        fig.write_html(rf'./pics_html/Y_l_{l}_m_{m}.html')

  1. 借用matplotlib实现静态的可视化:
import numpy as np
import matplotlib.pyplot as plt
import matplotlib.gridspec as gridspec
# The following import configures Matplotlib for 3D plotting.
from mpl_toolkits.mplot3d import Axes3D
from scipy.special import sph_harm


# plt.rc('text', usetex=True)

# Grids of polar and azimuthal angles
theta = np.linspace(0, np.pi, 100)
phi = np.linspace(0, 2*np.pi, 100)
# Create a 2-D meshgrid of (theta, phi) angles.
theta, phi = np.meshgrid(theta, phi)
# Calculate the Cartesian coordinates of each point in the mesh.
xyz = np.array([np.sin(theta) * np.sin(phi),
                np.sin(theta) * np.cos(phi),
                np.cos(theta)])

def plot_Y(ax, el, m):
    """Plot the spherical harmonic of degree el and order m on Axes ax."""

    # NB In SciPy's sph_harm function the azimuthal coordinate, theta,
    # comes before the polar coordinate, phi.
    Y = sph_harm(abs(m), el, phi, theta)

    # Linear combination of Y_l,m and Y_l,-m to create the real form.
    if m < 0:
        Y = np.sqrt(2) * (-1)**m * Y.imag
    elif m > 0:
        Y = np.sqrt(2) * (-1)**m * Y.real
    Yx, Yy, Yz = np.abs(Y) * xyz

    # Colour the plotted surface according to the sign of Y.
    cmap = plt.cm.ScalarMappable(cmap='RdBu')
    cmap.set_clim(-0.5, 0.5)

    ax.plot_surface(Yx, Yy, Yz,
                    facecolors=cmap.to_rgba(Y.real),
                    rstride=2, cstride=2)

    # Draw a set of x, y, z axes for reference.
    ax_lim = 0.5
    ax.plot([-ax_lim, ax_lim], [0,0], [0,0], c='0.5', lw=1, zorder=10)
    ax.plot([0,0], [-ax_lim, ax_lim], [0,0], c='0.5', lw=1, zorder=10)
    ax.plot([0,0], [0,0], [-ax_lim, ax_lim], c='0.5', lw=1, zorder=10)
    # Set the Axes limits and title, turn off the Axes frame.
    # ax.set_title(r'$Y_{{{},{}}}$'.format(el, m))
    ax.set_title('Y_l_{}_m_{}'.format(el, m))
    ax_lim = 0.5
    ax.set_xlim(-ax_lim, ax_lim)
    ax.set_ylim(-ax_lim, ax_lim)
    ax.set_zlim(-ax_lim, ax_lim)
    ax.axis('off')

# fig = plt.figure(figsize=plt.figaspect(1.))

for l in range(0, 4):
    for m in range(-l, l+1):
        fig = plt.figure()
        ax = fig.add_subplot(projection='3d')
        plot_Y(ax, l, m)
        plt.savefig('./pics_png/Y_l_{}_m_{}.png'.format(l, m))

静态效果如下:
请添加图片描述
OK,至此,GNN中常用的基组(至少我所了解到的)介绍完了。
一般来说,仅涉及距离信息的架构常常采用 gaussian 基组。

如果要用 spherical harmonics 这种涉及角度的基组,一般需要将几何坐标转到球极坐标下,而这将导致网络适应等变架构时遇到困难。

当然,还有使用 tensor field 做基组的,这块我还了解的少,但看起来好像也是套的 spherical harmonics 。

本文来自互联网用户投稿,该文观点仅代表作者本人,不代表本站立场。本站仅提供信息存储空间服务,不拥有所有权,不承担相关法律责任。如若转载,请注明出处:/a/91008.html

如若内容造成侵权/违法违规/事实不符,请联系我们进行投诉反馈qq邮箱809451989@qq.com,一经查实,立即删除!

相关文章

Linux线程 --- 生产者消费者模型(C语言)

在学习完线程相关的概念之后&#xff0c;本节来认识一下Linux多线程相关的一个重要模型----“ 生产者消费者模型” 本文参考&#xff1a; Linux多线程生产者与消费者_红娃子的博客-CSDN博客 Linux多线程——生产者消费者模型_linux多线程生产者与消费者_两片空白的博客-CSDN博客…

Discuz!论坛发帖标题字数限制80字符可以修改吗?修改发帖标题字数的方法

Discuz!论坛发帖标题字数限制80字符修改方法 1.数据库修改2.修改JS验证字符数文件3.修改模板中写死的字符限制数4.修改函数验证文件5.修改语言包文件6.更新缓存 Discuz X3.4论坛网站帖子标题字数限制80字符&#xff0c;当我们想使用长标题的时候就得一删再删&#xff0c;实在是…

JAVA switch case 穿透问题

1&#xff0c;前提 其实开发中很少会用到switch &#xff0c;一般更倾向于if-else&#xff0c; 但是最近接手的项目&#xff0c;前人写的代码都用switch &#xff0c; 但是我一直以来对switch 的理解就跟if一样&#xff0c; 然后项目运用的时候才发现这玩意居然还有穿透问题 …

nginx(七十八)nginx配置http2

一 ngx_http_v2模块 1、本文不讲解HTTP2的知识2、只讲解nginx中如何配置HTTP2 ① 前置条件 1、openssl的版本必须在1.0.2e及以上2、开启https加密,目前http2.0只支持开启了https的网站编译选项&#xff1a;--with-http_ssl_module --with-http_v2_module 特点&#xff1a…

[管理与领导-51]:IT基层管理者 - 8项核心技能 - 6 - 流程

前言&#xff1a; 管理者存在的价值就是制定目标&#xff0c;即目标管理、通过团队&#xff08;他人&#xff09;拿到结果。 要想通过他人拿到结果&#xff1a; &#xff08;1&#xff09;目标&#xff1a;制定符合SMART原则的符合业务需求的目标&#xff0c;团队跳一跳就可以…

Vue3 [Day11]

Vue3的优势 create-vue搭建Vue3项目 node -v npm init vuelatest npm installVue3项目目录和关键文件 Vetur插件是Vue2的 Volarr插件是Vue3的 main.js import ./assets/main.css// new Vue() 创建一个应用实例 > createApp() // createRouter() createStore() // 将创建实…

基于JAYA算法优化的BP神经网络(预测应用) - 附代码

基于JAYA算法优化的BP神经网络&#xff08;预测应用&#xff09; - 附代码 文章目录 基于JAYA算法优化的BP神经网络&#xff08;预测应用&#xff09; - 附代码1.数据介绍2.JAYA优化BP神经网络2.1 BP神经网络参数设置2.2 JAYA算法应用 4.测试结果&#xff1a;5.Matlab代码 摘要…

容器技术,1. Docker,2. Kubernetes(K8s):

目录 容器技术 1. Docker&#xff1a; 2. Kubernetes&#xff08;K8s&#xff09;&#xff1a; Docker和Kubernetes 容器的主要应用场景有哪些&#xff1f; 容器技术 有效的将单个操作系统的资源划分到孤立的组中&#xff0c;以便更好的在孤立的组之间平衡有冲突的资源使…

Kotlin协程flow发送时间间隔debounce

Kotlin协程flow发送时间间隔debounce debounce的作用是让连续发射的数据之间间隔起来。典型的应用场景是搜索引擎里面的关键词输入&#xff0c;当用户输入字符时候&#xff0c;有时候&#xff0c;并不希望用户每输入任何一个单字就触发一次后台真正的查询&#xff0c;而是希望…

《Dive into Deep Learning》

《Dive into Deep Learning》&#xff1a;https://d2l.ai/ Interactive deep learning book with code, math, and discussionsImplemented with PyTorch, NumPy/MXNet, JAX, and TensorFlowAdopted at 500 universities from 70 countries 《动手学深度学习》中文版&#xff1…

深度学习10:Attention 机制

目录 Attention 的本质是什么 Attention 的3大优点 Attention 的原理 Attention 的 N 种类型 Attention 的本质是什么 Attention&#xff08;注意力&#xff09;机制如果浅层的理解&#xff0c;跟他的名字非常匹配。他的核心逻辑就是「从关注全部到关注重点」。 Attention…

基于FPGA的Lorenz混沌系统verilog开发,含testbench和matlab辅助测试程序

目录 1.算法运行效果图预览 2.算法运行软件版本 3.部分核心程序 4.算法理论概述 5.算法完整程序工程 1.算法运行效果图预览 将vivado的仿真结果导入到matlab显示三维混沌效果&#xff1a; 2.算法运行软件版本 vivado2019.2 matlab2022a 3.部分核心程序 testbench如下所…

Docker安装ES+kibana8.9.1

参考&#xff1a;基于Docker安装Elasticsearch【保姆级教程、内含图解】_docker elasticsearch_Acloasia的博客-CSDN博客 创建网络 docker network create es-net 基于Docker安装Elasticsearch 拉取镜像 docker pull elasticsearch:8.9.1 挂载文件 mkdir -p /usr/local/e…

stm32之USART(总结)

串行通信 UART串口内部结构示意图 普中科技的详细介绍 中断知识补充 代码 #ifndef __USART_H #define __USART_H #include "stdio.h" #include "stm32f10x_usart.h" #define USART1_REC_LEN 200 //定义最大接收字节数 200extern u8 USART1_RX_BUF[US…

wx.request配置服务器域名,只能包含英文大小写字母、数字,解决办法

前言.小程序服务器域名配置常见错误及解决方法 1.配置入口&#xff1a; 小程序后台->-开发->开发设置->服务器域名 2.常见错误及原因分析&#xff1a; 3.实战中出现的错误 4.解决办法&#xff1a;应把域名后边的路径去掉&#xff0c;只写域名即可

【安卓】自定义View实现画板涂鸦等功能

一、实现效果 二、代码 1、MainActivity.class package com.lsl.mydrawingboarddemo;import androidx.appcompat.app.AppCompatActivity; import androidx.core.content.ContextCompat;import android.os.Bundle; import android.os.Handler; import android.view.View; impo…

4.17 如何基于 UDP 协议实现可靠传输?

目录 QUIC 是如何实现可靠传输的&#xff1f; Packet Header QUIC Frame Header QUIC 是如何解决 TCP 队头阻塞问题的&#xff1f; 什么是TCP对头阻塞问题&#xff1a; HTTP/2 的队头阻塞: 没有队头阻塞的 QUIC QUIC 是如何做流量控制的&#xff1f; QUIC 实现流量控制…

探索未来世界,解密区块链奥秘!

你是否曾好奇&#xff0c;区块链是如何影响着我们的生活与未来&#xff1f;想要轻松了解这个引领着技术革命的概念吗&#xff1f;那么这本令人着迷的新书《区块链导论》绝对值得你拥有&#xff01; 内容丰富多彩&#xff0c;让你轻松掌握&#xff1a; **1章&#xff1a;区块链…

MyBatis与Spring整合以及AOP和PageHelper分页插件整合

目录 前言 一、MyBatis与Spring整合的好处以及两者之间的关系 1.好处 2.关系 二、MyBatis和Spring集成 1.导入pom.xml 2.编写配置文件 3.利用mybatis逆向工程生成模型层代码 三、常用注解 四、AOP整合pageHelper分页插件 创建一个切面 测试 前言 MyBatis是一个开源的…

Unity 图片资源的适配

前言 最近小编做Unity项目时&#xff0c;发现在资源处理这方面和Android有所不同&#xff1b;例如&#xff1a;Android的资源文件夹res下会有着mipmap-mdpi&#xff0c;mipmap-hdpi&#xff0c;mipmap-xhdpi&#xff0c;mipmap-xxhdpi&#xff0c;mipmap-xxxhdpi这五个文件夹&a…