Datawhale AI夏令营--第二届世界科学智能大赛物质科学赛道：催化反应产率预测（2）

一、RNN建模SMILES进行反应产率预测

本次练习围绕RNN建模来进行预测，先来了解一下吧。

RNN建模SMILES是一种使用循环神经网络（‌RNN）‌对SMILES表示法进行建模的方法。‌

SMILES（‌Simplified Molecular Input Line Entry System）‌是一种用于描述化学分子的线性表示法，‌它能够将复杂的分子结构转化为简单的字符串形式，‌便于计算机处理和分析。‌RNN特别适合处理序列数据，‌因此在化学信息学中，‌RNN被广泛应用于SMILES字符串的建模和分析，‌以预测分子的各种性质或生成新的分子结构。‌

具体来说，‌RNN在SMILES建模中的应用包括：‌

• 预测分子性质：‌通过训练，‌RNN可以学习到SMILES字符串与分子性质之间的关系，‌从而能够预测新分子的性质，‌如溶解度、‌毒性等。‌
• 生成新的分子结构：‌RNN还能够根据给定的条件或约束，‌生成符合特定要求的新的SMILES字符串，‌这有助于在药物发现等领域中探索新的化学空间。‌

此外，‌RNN在SMILES建模中的应用还涉及到处理SMILES字符串中的序列依赖性问题，‌即分子结构中的原子和化学键之间的复杂相互作用和依赖关系。‌通过训练，‌RNN能够学习到这些依赖关系，‌从而更好地理解和预测分子的性质。‌

 二、代码详细解释

 具体操作步骤和第一次的相同，我在这里就不再写一次了。具体看我的第一篇。

（https://blog.csdn.net/m0_56875011/article/details/140714771）

在这次的建模还需要一个词汇表的文件（vocab_full.txt），我放在附件了。

现在开始我们来看一下代码

1、导包

--re 模块是Python的标准库之一，用于字符串的正则表达式操作。它提供了对正则表达式的全面支持，允许你进行复杂的文本搜索、替换和字符串解析等操作。

--time 模块提供了各种与时间相关的函数。它可以用来获取当前时间、日期，执行时间测量（比如，计算某个操作的执行时间），以及进行时间格式化等。

--pandas 是一个强大的Python数据分析库，提供了高性能、易用的数据结构和数据分析工具。这里通过 as pd 将其别名设置为 pd，这是 pandas 常用的别名，便于在代码中引用。

--从 typing 模块中导入了 List 和 Tuple 类型注解。这些类型注解用于提高代码的可读性和可维护性，尤其是在使用静态类型检查工具（如mypy）时。它们不会改变代码的运行时行为，但可以帮助开发者在编写代码时捕获潜在的错误。

--torch 是PyTorch库的主要模块，PyTorch是一个广泛使用的开源机器学习库，特别是在深度学习领域。它提供了强大的张量（多维数组）操作和自动微分系统，非常适合于构建和训练神经网络。

--从 torch 库中导入了 nn 模块，并将其别名设置为 nn。nn 模块包含了构建神经网络所需的所有构件，如层（layers）、激活函数（activations）等。使用 nn 可以方便地定义和训练神经网络。

--从 torch 库中导入了 optim 模块，并将其别名设置为 optim。optim 模块包含了多种优化算法，用于神经网络的参数更新。例如，SGD（随机梯度下降）、Adam等。

--从 torch.utils.data 模块中导入了 Dataset、DataLoader 和 Subset。这些工具用于数据的加载和预处理，以便在训练神经网络时使用。

--Subset 用于从现有的数据集中提取一个子集，常用于划分训练集、验证集和测试集。

DataLoader 提供了一个可迭代的对象，它将数据集封装起来，支持批量加载、打乱数据、多进程数据加载等功能。

Dataset 是一个抽象类，用于表示数据集。你需要继承这个类并实现 __getitem__ 和 __len__ 方法来创建自己的数据集。

import re
import time
import pandas as pd
from typing import List, Tuple
import torch
import torch.nn as nn
import torch.optim as optim
from torch.utils.data import Dataset, DataLoader, Subset

2、 类定义

• RNNModel(nn.Module): 继承自nn.Module，是构建所有PyTorch模型的基础。

①初始化方法 __init__

• num_embed: 嵌入层的词汇表大小，即输入词汇的数量。
• input_size: 嵌入层输出的特征维度，即每个词嵌入的向量大小。
• hidden_size: RNN隐藏层的大小。
• output_size: 模型的输出大小，通常对应于分类任务中的类别数。
• num_layers: RNN层的数量。
• dropout: 在RNN层之间应用dropout的比例，用于防止过拟合。
• device: 模型运行的设备（如CPU或GPU）。

在初始化方法中：

• self.embed: 定义了一个嵌入层，将输入的词汇索引（整数）转换为固定大小的密集向量。
• self.rnn: 定义了一个双向RNN层。batch_first=True表示输入张量的第一个维度是批量大小，dropout用于层间的dropout，bidirectional=True表示RNN是双向的，因此最终隐藏状态hn会包含两个方向的信息。
• self.fc: 定义了一个全连接层序列，首先是一个线性层将双向RNN的隐藏状态（已展平）映射到output_size大小，然后通过Sigmoid激活函数，再通过一个线性层映射到单个输出（可能是二分类问题的概率），最后再次通过Sigmoid激活函数得到最终的输出。这里的设计有点不寻常，因为通常最后的输出层不需要再次通过Sigmoid，除非有特殊的理由（比如输出需要被解释为概率）。

②前向传播方法 forward

• x: 输入的词汇索引张量，形状为[bs, seq_len]，其中bs是批量大小，seq_len是序列长度。

在前向传播方法中：

• 首先，输入x通过嵌入层转换为形状为[bs, seq_len, input_size]的张量。
• 然后，这个张量通过双向RNN层。RNN层返回输出和最终隐藏状态hn，但在这个方法中只使用了hn（因为这里可能只关心序列的最终状态或整体表示）。hn的形状是[2*num_layers, bs, hidden_size]，因为RNN是双向的，所以层数是num_layers的两倍。
• hn通过transpose(0,1)操作调整维度，使其形状变为[bs, 2*num_layers, hidden_size]，这样第一个维度就是批量大小。
• hn被展平（通过reshape操作）成一个二维张量z，形状为[bs, 2*num_layers*hidden_size]，这样它就可以作为全连接层的输入。
• 最后，z通过全连接层序列self.fc，得到最终的输出output，其形状为[bs, 1]，表示每个样本的预测结果。

# 定义RNN模型
class RNNModel(nn.Module):
    def __init__(self, num_embed, input_size, hidden_size, output_size, num_layers, dropout, device):
        super(RNNModel, self).__init__()
        self.embed = nn.Embedding(num_embed, input_size)
        self.rnn = nn.RNN(input_size, hidden_size, num_layers=num_layers, 
                          batch_first=True, dropout=dropout, bidirectional=True)
        self.fc = nn.Sequential(nn.Linear(2 * num_layers * hidden_size, output_size),
                                nn.Sigmoid(),
                                nn.Linear(output_size, 1),
                                nn.Sigmoid())

    def forward(self, x):
        # x : [bs, seq_len]
        x = self.embed(x)
        # x : [bs, seq_len, input_size]
        _, hn = self.rnn(x) # hn : [2*num_layers, bs, h_dim]
        hn = hn.transpose(0,1)
        z = hn.reshape(hn.shape[0], -1) # z shape: [bs, 2*num_layers*h_dim]
        output = self.fc(z).squeeze(-1) # output shape: [bs, 1]
        return output

 3、数据处理

• 正则表达式处理：

  • 您的_regex_match方法尝试使用正则表达式来分割SMILES字符串，但这里的实现似乎有误。您希望按字符拆分SMILES字符串，并可能通过正则表达式过滤掉一些不需要的字符，但当前的实现会对每个字符分别应用正则表达式，这不太可能是您想要的。通常，SMILES字符串不需要复杂的字符拆分，除非您想过滤掉某些特定字符。
  • 如果您的目标是简单地按字符拆分SMILES字符串，那么不需要使用正则表达式。如果需要过滤字符，可以考虑在应用拆分之前先应用一个正则表达式来清理字符串。

• 字符到索引的映射：

  • 在_pad_token_to_idx方法中，您尝试将字符转换为索引。但是，如果SMILES字符串包含不在词汇表中的字符，您会动态地向词汇表添加这些字符并分配索引。这在某些情况下可能有用，但可能会导致模型训练时的不一致性，因为训练集和测试集之间的词汇表可能不同。
  • 更好的做法是在预处理阶段就确定一个固定的词汇表，并确保所有输入数据都只包含这个词汇表中的字符。

• tokenize方法的逻辑：

  • 您的tokenize方法试图对每个SMILES字符串（假设是单个字符串，但您的方法名和数据结构暗示了多个字符串）应用_regex_match，但_regex_match的实现似乎并不适合这个目的。
  • 如果您的输入是单个SMILES字符串，那么您应该直接在这个字符串上应用分词逻辑，而不是试图将其视为字符列表。
  • 如果您的输入是多个SMILES字符串的列表，那么您应该遍历这个列表，并对每个字符串分别应用分词逻辑。

read_data 函数 ：

这个函数用于从指定路径的CSV文件中读取数据，并根据是否用于训练集（train参数控制）来返回不同的输出。

1. 读取数据：使用pandas的read_csv函数读取CSV文件，并获取特定的列（Reactant1, Reactant2, Product, Additive, Solvent, Yield）。

2. 处理数据：

   • 将Reactant1和Reactant2拼接成一个字符串，用.分隔，然后将这个字符串与Product拼接，用>分隔。这样做是为了构造一个包含反应物和产物的输入字符串，用于后续的模型输入。
   • 如果train为True，则保留Yield（产率）数据；如果为False（通常用于测试集或验证集），则将Yield全部设置为0。

3. 构造输出：将处理后的输入字符串和对应的Yield（或全0的列表）组成一个元组列表，其中每个元组包含一个输入字符串和对应的产率。

4. 可选的序列长度统计：注释部分提供了如何统计序列长度的代码示例，这对于了解数据特性和模型设计很有帮助。

# import matplotlib.pyplot as plt
## 数据处理部分
# tokenizer，鉴于SMILES的特性，这里需要自己定义tokenizer和vocab
# 这里直接将smiles str按字符拆分，并替换为词汇表中的序号
class Smiles_tokenizer():
    def __init__(self, pad_token, regex, vocab_file, max_length):
        self.pad_token = pad_token
        self.regex = regex
        self.vocab_file = vocab_file
        self.max_length = max_length

        with open(self.vocab_file, "r") as f:
            lines = f.readlines()
        lines = [line.strip("\n") for line in lines]
        vocab_dic = {}
        for index, token in enumerate(lines):
            vocab_dic[token] = index
        self.vocab_dic = vocab_dic

    def _regex_match(self, smiles):
        regex_string = r"(" + self.regex + r"|"
        regex_string += r".)"
        prog = re.compile(regex_string)

        tokenised = []
        for smi in smiles:
            tokens = prog.findall(smi)
            if len(tokens) > self.max_length:
                tokens = tokens[:self.max_length]
            tokenised.append(tokens) # 返回一个所有的字符串列表
        return tokenised
    
    def tokenize(self, smiles):
        tokens = self._regex_match(smiles)
        # 添加上表示开始和结束的token：<cls>, <end>
        tokens = [["<CLS>"] + token + ["<SEP>"] for token in tokens]
        tokens = self._pad_seqs(tokens, self.pad_token)
        token_idx = self._pad_token_to_idx(tokens)
        return tokens, token_idx

    def _pad_seqs(self, seqs, pad_token):
        pad_length = max([len(seq) for seq in seqs])
        padded = [seq + ([pad_token] * (pad_length - len(seq))) for seq in seqs]
        return padded

    def _pad_token_to_idx(self, tokens):
        idx_list = []
        for token in tokens:
            tokens_idx = []
            for i in token:
                if i in self.vocab_dic.keys():
                    tokens_idx.append(self.vocab_dic[i])
                else:
                    self.vocab_dic[i] = max(self.vocab_dic.values()) + 1
                    tokens_idx.append(self.vocab_dic[i])
            idx_list.append(tokens_idx)
        
        return idx_list
# 读数据并处理
def read_data(file_path, train=True):
    df = pd.read_csv(file_path)
    reactant1 = df["Reactant1"].tolist()
    reactant2 = df["Reactant2"].tolist()
    product = df["Product"].tolist()
    additive = df["Additive"].tolist()
    solvent = df["Solvent"].tolist()
    if train:
        react_yield = df["Yield"].tolist()
    else:
        react_yield = [0 for i in range(len(reactant1))]
    
    # 将reactant拼到一起，之间用.分开。product也拼到一起，用>分开
    input_data_list = []
    for react1, react2, prod, addi, sol in zip(reactant1, reactant2, product, additive, solvent):
        input_info = ".".join([react1, react2])
        input_info = ">".join([input_info, prod])
        input_data_list.append(input_info)
    output = [(react, y) for react, y in zip(input_data_list, react_yield)]

    # 下面的代码将reactant\additive\solvent拼到一起，之间用.分开。product也拼到一起，用>分开
    '''
    input_data_list = []
    for react1, react2, prod, addi, sol in zip(reactant1, reactant2, product, additive, solvent):
        input_info = ".".join([react1, react2, addi, sol])
        input_info = ">".join([input_info, prod])
        input_data_list.append(input_info)
    output = [(react, y) for react, y in zip(input_data_list, react_yield)]
    '''

    # # 统计seq length，序列的长度是一个重要的参考，可以使用下面的代码统计查看以下序列长度的分布
    # seq_length = [len(i[0]) for i in output]
    # seq_length_400 = [len(i[0]) for i in output if len(i[0])>200]
    # print(len(seq_length_400) / len(seq_length))
    # seq_length.sort(reverse=True)
    # plt.plot(range(len(seq_length)), seq_length)
    # plt.title("templates frequence")
    # plt.show()
    return output

ReactionDataset 类

这个类继承自torch.utils.data.Dataset（尽管这里没有显式导入，但假设它是这样使用的），用于封装处理后的数据，使其能够被PyTorch等数据加载工具使用。

1. 初始化：在__init__方法中，接收一个包含元组（每个元组包含一个输入字符串和对应的产率）的列表作为数据集。

2. 长度：__len__方法返回数据集中元素的数量。

3. 获取项：__getitem__方法允许通过索引访问数据集中的单个元素（即一个元组）。

class ReactionDataset(Dataset):
    def __init__(self, data: List[Tuple[List[str], float]]):
        self.data = data
        
    def __len__(self):
        return len(self.data)

    def __getitem__(self, idx):
        return self.data[idx]
    
def collate_fn(batch):
    REGEX = r"\[[^\]]+]|Br?|Cl?|N|O|S|P|F|I|b|c|n|o|s|p|\(|\)|\.|=|#|-|\+|\\\\|\/|:|~|@|\?|>|\*|\$|\%[0-9]{2}|[0-9]"
    tokenizer = Smiles_tokenizer("<PAD>", REGEX, "../vocab_full.txt", max_length=300)
    smi_list = []
    yield_list = []
    for i in batch:
        smi_list.append(i[0])
        yield_list.append(i[1])
    tokenizer_batch = torch.tensor(tokenizer.tokenize(smi_list)[1])
    yield_list = torch.tensor(yield_list)
return tokenizer_batch, yield_list

 4、训练数据

这段代码是一个使用PyTorch框架进行训练循环的示例，旨在训练一个RNN（循环神经网络）模型来处理某种化学反应数据。下面是对代码的详细解析：

1. 参数设置

• N: 设置为10，这里看起来是打算从数据集中选取前10个样本进行训练，但后面代码并未直接使用这个变量来限制数据加载，这可能是个遗留或未完成的逻辑。
• NUM_EMBED: 嵌入层的维度，设为294，这通常用于将分类变量（如词汇表中的单词）转换为固定大小的密集向量。
• INPUT_SIZE: 输入序列的长度，设为300。
• HIDDEN_SIZE: RNN隐藏层的维度，设为512。
• OUTPUT_SIZE: 输出层的维度，也设为512，这里可能需要根据实际任务（如回归任务的输出维度）进行调整。
• NUM_LAYERS: RNN的层数，设为10，这是一个相对较大的值，可能需要调整以避免过拟合或训练时间过长。
• DROPOUT: Dropout比率，用于防止过拟合，设为0.2。
• CLIP: 梯度裁剪的值，设为1，用于防止梯度爆炸。
• N_EPOCHS: 训练轮数，设为100。
• LR: 学习率，设为0.0001。

2. 数据加载和预处理

• 使用read_data函数读取数据（该函数未在代码中定义，假设它返回一个适合训练的数据列表）。
• 创建一个ReactionDataset实例来封装数据。
• 尝试使用Subset来选取数据的一个子集（这里实际上只创建了一个索引列表，但没有在后续中使用这个子集，可能是个错误或遗漏）。
• 使用DataLoader来加载数据，设置批量大小为128，启用打乱，并指定了一个collate_fn（未在代码中定义，通常用于自定义数据批次的处理方式）。

3. 模型定义和训练

• 定义了一个RNNModel类（未在代码中给出），它接受上述超参数并创建RNN模型。
• 将模型移至指定的设备（GPU或CPU）。
• 使用Adam优化器进行训练。
• 选择L1Loss（MAE）作为损失函数。
• 在每个epoch中，遍历训练数据，计算损失，进行反向传播，优化参数，并裁剪梯度。
• 打印每个epoch的平均损失，并在损失降低时保存模型。

4. 运行时间和性能

• 记录训练开始和结束的时间，计算并打印总运行时间。

def train():
    ## super param
    N = 10  #int / int(len(dataset) * 1)  # 或者你可以设置为数据集大小的一定比例，如 int(len(dataset) * 0.1)
    NUM_EMBED = 294 # nn.Embedding()
    INPUT_SIZE = 300 # src length
    HIDDEN_SIZE = 512
    OUTPUT_SIZE = 512
    NUM_LAYERS = 10
    DROPOUT = 0.2
    CLIP = 1 # CLIP value
    N_EPOCHS = 100
    LR = 0.0001
    
    start_time = time.time()  # 开始计时
    device = torch.device('cuda' if torch.cuda.is_available() else 'cpu')
    # device = 'cpu'
    data = read_data("../dataset/round1_train_data.csv")
    dataset = ReactionDataset(data)
    subset_indices = list(range(N))
    subset_dataset = Subset(dataset, subset_indices)
    train_loader = DataLoader(dataset, batch_size=128, shuffle=True, collate_fn=collate_fn)

    model = RNNModel(NUM_EMBED, INPUT_SIZE, HIDDEN_SIZE, OUTPUT_SIZE, NUM_LAYERS, DROPOUT, device).to(device)
    model.train()
    
    optimizer = optim.Adam(model.parameters(), lr=LR)
    # criterion = nn.MSELoss() # MSE
    criterion = nn.L1Loss() # MAE

    best_loss = 10
    for epoch in range(N_EPOCHS):
        epoch_loss = 0
        for i, (src, y) in enumerate(train_loader):
            src, y = src.to(device), y.to(device)
            optimizer.zero_grad()
            output = model(src)
            loss = criterion(output, y)
            loss.backward()
            torch.nn.utils.clip_grad_norm_(model.parameters(), CLIP)
            optimizer.step()
            epoch_loss += loss.item()
            loss_in_a_epoch = epoch_loss / len(train_loader)
        print(f'Epoch: {epoch+1:02} | Train Loss: {loss_in_a_epoch:.3f}')
        if loss_in_a_epoch < best_loss:
            # 在训练循环结束后保存模型
            torch.save(model.state_dict(), '../model/RNN.pth')
    end_time = time.time()  # 结束计时
    # 计算并打印运行时间
    elapsed_time_minute = (end_time - start_time)/60
    print(f"Total running time: {elapsed_time_minute:.2f} minutes")

if __name__ == '__main__':
    train()

5、打印生成文件，比赛官网提交文件，得到分数。

这段代码定义了一个函数 predicit_and_make_submit_file，旨在使用训练好的RNN模型对测试数据进行预测，并将预测结果保存到指定的输出文件中，以便提交。下面是对这段代码的详细解析：

1. 函数定义和参数：

   • model_file: 字符串，指向包含最佳模型状态的文件的路径。
   • output_file: 字符串，指向输出文件的路径，该文件将包含预测结果。

2. 模型参数：

   • 定义了模型架构所需的参数，如嵌入层大小（NUM_EMBED）、输入序列长度（INPUT_SIZE）、隐藏层大小（HIDDEN_SIZE）、输出层大小（OUTPUT_SIZE）、RNN层数（NUM_LAYERS）和dropout比率（DROPOUT）。

3. 设备选择：

   • 使用torch.device根据CUDA的可用性选择设备（GPU或CPU）。

4. 数据加载：

   • 使用read_data函数（未在代码中定义，假设它返回测试数据集）读取测试数据。
   • 创建ReactionDataset实例来封装测试数据。
   • 使用DataLoader以批量方式加载测试数据，设置批量大小为64，不打乱数据（因为测试数据应按顺序处理），并指定collate_fn（未在代码中定义，但通常用于处理数据批次）。

5. 模型加载和预测：

   • 创建RNNModel实例，并将其移动到选定的设备上。
   • 加载最佳模型的状态字典（从model_file中加载）。
   • 将模型设置为评估模式（model.eval()），这通常会关闭dropout等训练特定的功能。
   • 遍历测试数据加载器，对于每个批次的数据，执行前向传播以获得预测结果，并将结果从张量转换为列表。

6. 结果处理：

   • 初始化一个字符串列表ans_str_lst，用于存储最终的预测结果字符串。首先添加一个标题行。
   • 遍历预测结果列表，为每个预测结果生成一个包含test编号和预测值的字符串，并将其添加到ans_str_lst中。注意，这里假设y（即output中的元素）是一个标量（单个预测值），这可能与模型的实际输出不符，具体取决于RNNModel的定义。

7. 保存结果：

   • 将ans_str_lst中的字符串使用换行符连接成一个长字符串，并将其写入到output_file指定的文件中。

# 生成结果文件
def predicit_and_make_submit_file(model_file, output_file):
    NUM_EMBED = 294
    INPUT_SIZE = 300
    HIDDEN_SIZE = 512
    OUTPUT_SIZE = 512
    NUM_LAYERS = 10
    DROPOUT = 0.2
    device = torch.device('cuda' if torch.cuda.is_available() else 'cpu')
    test_data = read_data("../dataset/round1_test_data.csv", train=False)
    test_dataset = ReactionDataset(test_data)
    test_loader = DataLoader(test_dataset, batch_size=64, shuffle=False, collate_fn=collate_fn) 

    model = RNNModel(NUM_EMBED, INPUT_SIZE, HIDDEN_SIZE, OUTPUT_SIZE, NUM_LAYERS, DROPOUT, device).to(device)
    # 加载最佳模型
    model.load_state_dict(torch.load(model_file))
    model.eval()
    output_list = []
    for i, (src, y) in enumerate(test_loader):
        src, y = src.to(device), y.to(device)
        with torch.no_grad():
            output = model(src)
            output_list += output.detach().tolist()
    ans_str_lst = ['rxnid,Yield']
    for idx,y in enumerate(output_list):
        ans_str_lst.append(f'test{idx+1},{y:.4f}')
    with open(output_file,'w') as fw:
        fw.writelines('\n'.join(ans_str_lst))

    print("done!!!")