利用 Amazon HPC 加速药物研发分子对接任务

背景介绍

在药物早期发现阶段，先导化合物的分子筛选是药物发现阶段靶点确认后最重要的起始环节。虚拟筛选是一个重要且十分有效的途径受到广泛的关注并逐渐在工业应用中展开。虚拟筛选的任务是从规模庞大的分子库中寻找与特定靶点位点能结合的小分子，在该过程中需要庞大的算力来满足大规模分子库的虚拟筛选。

在虚拟筛选领域应用比较多的开源软件是Auto Dock Vina, Vina作为一款开源的分子模拟软件，最主要应用于执行配体—蛋白分子对接，分子对接（molecular docking）是依据配体与受体作用的“锁-钥原理”，模拟配体小分子与受体生物大分子相互作用的一种技术方法。配体与受体相互作用是分子识别的过程，主要包括静电作用、氢键作用、疏水作用、范德华作用等。通过计算可以预测两者间的结合模式和亲和力。通常在Vina上进行一个分子对接需要1～10秒，如果需要柔性对接则需要更长时间。

在生产应用中，用户通常需要对数百万级别的分子和蛋白结构做对接，这样就需要非常长的时间周期。在传统方式在单机模式下运行一次任务需可能要数月时间，意味着单机不可能完成大规模的虚拟机筛选任务。故使用高性能的计算集群进行大规模虚拟筛选的能力对药物研发企业加速药物发现至关重要。

本文旨在介绍如何利用亚马逊云科技云服务实现药物发现阶段的大规模分子虚拟筛选任务。

主要涉及使用服务和软件介绍

Amazon ParallelCluster 是Amazon支持的开源集群管理工具，可帮助您在Amazon Web Services云中部署和管理高性能计算 (HPC) 集群。它会根据配置自动启动所需的计算资源和共享文件系统。支持Amazon Batch和Slurm等调度程序并充分利用云端的弹性和扩展性，用户只需为实际使用的计算和存储等资源付费。

Amazon FSx for Lustre是一项完全托管的服务，可为计算工作负载提供经济高效的高性能存储。Amazon FSx for Lustre可提供亚毫秒级延迟、高达每秒数百GB的吞吐量和数百万的IOPS。它提供多种部署选项和存储类型以优化成本和性能，从而让满足您的工作负载需求。Amazon FSx for Lustre文件系统还可以链接到 Amazon S3 存储桶，使您能够同时从高性能文件系统和S3 API访问和处理数据。

AutoDock Vina 是一个开源的被广泛使用的可用于进行分子对接的应用程序。

Open Babel 是一个应用非常广泛化学工具箱，在分子对接过程中用于分子格式转换、三维构象生成。

MGLtools是一个免费的分子预处理工具，用Vina做分子对接必须用这个工具来做预处理将分子生成pdbqt格式。

测试环境部署准备

云端环境准备

通过VPC创建向导创建包含公有子网和私有子网的VPC，这里我们把子网都放到同一个可用区方便后续访问Amazon FSx for Lustre共享存储。公有子网和私有子网分别配置指向互联网网关和NAT网关路由，详细配置方法可以参考链接。

登录亚马逊云科技管理控台并切换至FSx服务，根据向导创建一个1.2TiB, 60MiB/s吞吐的SSD存储类型的文件系统。

测试软件安装

准备一台Amazon Linux2实例，参考文档安装Lustre客户端并将之前创建的Amazon FSx for Lustre文件系统挂载到/fsx目录，这里我们已经提前将需要的3个测试需要的软件下载到/fsx目录，如下图所示：

这里我们将3个软件均安装到/fsx目录，主要是方便后续计算节点启动后直接使用我们已经提前预算好在共享目录的测试软件。相关软件安装脚本如下：

# 安装 AutoDock Vina
cd /fsx
mv vina_1.2.3_linux_x86_64 vina
sudo chmod +x vina
./vina --help
# 安装Open Babel 参考https://open-babel.readthedocs.io/en/latest/Installation/install.html#compiling-open-babel
cd /fsx
tar -zxvf openbabel-openbabel-2-4-0.tar.gz
mkdir -p /fsx/build
cd build
cmake ../openbabel-openbabel-2-4-0 -DCMAKE_INSTALL_PREFIX=/fsx/openbabel -DCMAKE_BUILD_TYPE=DEBUG
make -j4
make install
# 安装MGLTools
cd /fsx
chmod +x mgltools_Linux-x86_64_1.5.7_install
./mgltools_Linux-x86_64_1.5.7_install
cd /fsx/mgltools/mgltools_x86_64Linux2_1.5.7
./install -d /fsx/mgltools

测试数据准备

（1）用户提供sdf格式的分子库，这里我们已经提前下载到/fsx目录。

（2）和受体蛋白质靶点（pdbqt格式）

（3）对接grid space（config.txt）

config.txt文件示例：

#对接中心位置
center_x = 35.371
center_y = 17.826
center_z = 40.325
#对接范围
size_x = 20.0
size_y = 20.0
size_z = 20.0

注：受体蛋白靶点文件和对接位置文件提前拷贝到/fsx/acceptor目录。

安装并配置Amazon ParallelCluster

用户根据文档安装Amazon ParallelCluster, 用户可以通过命令行方式创建配置文件或者直接手动编辑配置文件config.yaml, 参考配置如下：

Region: cn-northwest-1
Image:
  Os: alinux2
HeadNode:
  InstanceType: m5.xlarge
  LocalStorage:
    RootVolume:
      Size: 50
      VolumeType: gp3
  Networking:
    SubnetId: subnet-022255c7910df393d
  Ssh:
    KeyName: zhy
  Dcv:
    Enabled: true
  Iam:
    AdditionalIamPolicies:
      - Policy: arn:aws-cn:iam::aws:policy/AdministratorAccess
Scheduling:
  Scheduler: slurm
  SlurmSettings:
    Dns:
      DisableManagedDns: false
    ScaledownIdletime: 10
  SlurmQueues:
  - Name: m5
    CapacityType: SPOT
    ComputeResources:
    - Name: m5xlarge
      InstanceType: m5.xlarge
      MinCount: 0
      MaxCount: 5
    Networking:
      SubnetIds:
      - subnet-013c82ea802f587ea
    Iam:
      AdditionalIamPolicies:
        - Policy: arn:aws-cn:iam::aws:policy/AmazonS3FullAccess
  - Name: spe
    CapacityType: SPOT
    ComputeResources:
    - Name: c524xlarge
      InstanceType: c5.24xlarge
      MinCount: 0
      MaxCount: 50
    Networking:
      SubnetIds:
      - subnet-013c82ea802f587ea
    Iam:
      AdditionalIamPolicies:
        - Policy: arn:aws-cn:iam::aws:policy/AmazonS3FullAccess
SharedStorage:
  - Name: fsx
    MountDir: /fsx
    StorageType: FsxLustre
    FsxLustreSettings:
      FileSystemId: fs-0fcc8186fed908c3e
Monitoring:
  DetailedMonitoring: true
  Logs:
    CloudWatch:
      Enabled: true
      RetentionInDays: 7
      DeletionPolicy:  Delete
  Dashboards:
    CloudWatch:
      Enabled: true

通过命令行创建集群

pcluster create-cluster --cluster-configuration config.yaml --cluster-name vina --rollback-on-failure false

查看集群状态

pcluster list-clusters

用户登录集群Head Node节点

pcluster ssh –cluster-name vina-i zhy.pem

任务流程介绍

批量虚拟筛选任务主要涉及到以下工作流程：

分子库拆分：我们将下载的比较大的的分子库拆分成独立的分子文件，这里一个分子库可能可能包含百万级别的分子数量。

分子构象生成：通过Open Babel工具为每个分子生成三维构象（sdf格式）。

格式转换：通过工具MGLTools将每个分子的三维构象从sdf转换为pdbqt格式用于后续分子对接。

AutoDock Vina分子对接: 一个分子对接任务的输入为分子的某一个三维构象（pdbqt格式），pdbqt格式的受体蛋白文件和靶点对接坐标信息（grid space）。

分子库拆分

用户可以通过命令行模式下执行 bash 1-Prepare_sdf.sh 命令，将当前目录下的分子库拆分成单个的小分子并保存至/fsx/sdf目录。这里我们将单个分子库(sdf文件)提交一个job到m5队列中。1-Prepare_sdf.sh 内容如下：

#创建所有文件目录
#sdf：分子库拆分成的小分子文件，sdf格式
#split_txt：小分子分组目录，即任务单元
#3D_sdf：分子的三维构象，sdf格式
#3D_mol2：分子的三维构象，mol2格式
#pdbqt: 分子的三维构象，pdbqt格式
#log: 日志
#!/bin/bash
cd /fsx
mkdir sdf 3D_mol2 3D_sdf split_txt pdbqt log

#将当前目录下的所有分子库（sdf文件）拆成单个分子，并保存至sdf目录
for f in *.sdf
 do
 echo $f
 b=`basename $f .sdf`
 sbatch -c 1 -p m5 --output /fsx/log/$f.out  --job-name Splitsdf --wrap "srun bash Prepare_sdf.sh $f $b"
done

Prepare_sdf.sh脚本内容如下：

#!/bin/bash
/fsx/openbabel/bin/obabel -isdf $1 -osdf -O /fsx/sdf/$2_.sdf -p 7.4 -m

分子构象生成

分子库拆分以后得到的每一个分子均保存在/fsx/sdf目录下，一个分子对应为一个sdf文件。这些分子从结构上看都是平面的，没有三维构象。由于对接需要分子的三维构象，所以需要先为每一个分子生成三维构象。这里我们通过软件Open Babel来实现。

用户在Shell模式通过运行脚本执行bash 2-Generate_3D_comformation.sh提交任务实现上述功能，脚本内容如下：

#!/bin/bash
#每$1个分子组成一个list,结果写入split_txt文件夹
find sdf | split -l $1 - /fsx/split_txt/Mol_ -d -a 6
#每一个list提交一个job.sh使用
for f in `ls /fsx/split_txt`
 do
 echo $f
 sbatch -c 16 -p spe --job-name GenerateComforation --output /fsx/log/$f.out  --wrap "srun bash Generate_3D_Conformation-Mol.sh /fsx/split_txt/$f"
done

其中Generate_3D_Conformation-Mol.sh将每一个list中的小分子生成三维构象，并保存至3D_mol2目录下。Generate_3D_Conformation-Mol.sh内容如下：

#!/bin/bash
k=`basename $1 /`
mkdir /fsx/3D_sdf/$k /fsx/3D_mol2/$k

#并发设置
tempfifo=$$.fifo
trap "exec 100>&-;exec 1000<&-;exit 0" 2
mkfifo $tempfifo
exec 1000<>$tempfifo
rm -rf $tempfifo

for ((i=0;i<100;i++))
do
   echo >&1000
done

#分子三维构象生成
for f in `cat $1`
 do
         read -u1000
         {
                 b=`basename $f .sdf`
                 echo $f
                 echo $b
                 echo Processing ligand $f $b
                 #生成三维构象，sdf格式，保存在3D_sdf目录下
                 /fsx/openbabel/bin/obabel /fsx/sdf/$b.sdf -O /fsx/3D_sdf/$k/$b.sdf --confab --rcutoff 1
                 #将每一个分子生成的多个三维构象，逐一转化为mol2格式，保存在3D_mol2目录下
                 /fsx/openbabel/bin/obabel -isdf /fsx/3D_sdf/$k/$b.sdf -omol2 -O /fsx/3D_mol2/$k/$b-.mol2 -p 7.4 --append ID -m
                 echo >&1000
 }&
 done
wait
echo $k finish!

执行 bash 2-Generate_3D_comformation.sh 200 命令后，将会每200个分子生成一个list，保存在split_txt 目录下，文件名为（例如）：Mol_000004。然后在/fsx/3D_sdf和/fsx/3D_mol2 目录下将会看到生成的各个分子的三维结构。

分子格式转换

同样我们通过执行脚本bash 3-GeneratePDBQT.sh 将3D_mol2 目录下的每一个Mol_00xxxx 文件夹提交一个任务将分子转化从 mol2 转化为pdbqt 格式，生成的 pdbqt 结果保存在/fsx/pdbqt 目录下对应的Mol_00xxxx 文件夹内。

#!/bin/bash
for f in `ls /fsx/split_txt`
 do
 echo $f
 sbatch -c 16 -p spe --job-name GeneratePDBQT --output /fsx/log/$f-pdbqt.out  --wrap "srun bash Generate_pdbqt.sh /fsx/split_txt/$f"
done

脚本Generate_pdbqt.sh内容如下：

#!/bin/bash
j=`basename $1 /`
mkdir /fsx/pdbqt/$j
cd /fsx/3D_mol2/$j
#并发
tempfifo=$$.fifo
trap "exec 100>&-;exec 1000<&-;exit 0" 2
mkfifo $tempfifo
exec 1000<>$tempfifo
rm -rf $tempfifo
for ((i=0;i<300;i++))
do
   echo >&1000
done

for f in ./*
 do
         read -u1000
         {
                 k=`basename $f .mol2`
                 echo $k
#使用mgltool工具将分子的mol2格式转变成pdbqt文件。 
/fsx/mgltools/mgltools_x86_64Linux2_1.5.7/bin/pythonsh /fsx/mgltools/mgltools_x86_64Linux2_1.5.7/MGLToolsPckgs/AutoDockTools/Utilities24/prepare_ligand4.py -l $f -F -A hydrogens
#pdbqt文件迁移到/fsx/pdbqt目录下 
mv $k.pdbqt ../../pdbqt/$j/
 echo >&1000
 }&
 done
wait
echo $1 finish!

批量虚拟筛选

在分子对接阶段我们通过如下指令将对接任务提交到Slurm的spe队列，参数为对接位置和靶点蛋白的pdbqt格式文件。

bash /fsx/4-docking.sh /fsx/acceptor/config.txt /fsx/acceptor/1a3w-acceptor.pdbqt

4-docking.sh内容如下：

#!/bin/bash
#docking日志
mkdir /fsx/docking_log/
#docking结果
mkdir /fsx/docking_result/
#docking分数
mkdir /fsx/docking_scores/

for f in `ls /fsx/pdbqt`
 do
#pdbqt中每一个文件夹一个job
 sbatch -c 16 -p spe --job-name docking --output /fsx/log/$f-docking.out  --wrap "srun bash docking.sh $1 $2 /fsx/pdbqt/$f"
done

docking.sh将使用Vina执行文件目录下所有分子构象（pdbqt）的对接任务，内容如下:

#!/bin/bash
j=`basename $3 /`
echo $j
#创建结果文件夹
mkdir /fsx/docking_result/$j
mkdir /fsx/docking_scores/$j
tempfifo=$$.fifo
trap "exec 100>&-;exec 1000<&-;exit 0" 2
mkfifo $tempfifo
exec 1000<>$tempfifo
rm -rf $tempfifo
for ((i=0;i<5;i++))
do
   echo >&1000
done

for f in `ls /fsx/pdbqt/$j`
 do
         read -u1000
         {
         echo $f
         m=`basename $f /`
         echo $m
         #执行对接
         /fsx/vina --config $1 --ligand /fsx/pdbqt/$j/$m --out /fsx/docking_result/$j/$m-docked.pdbqt --receptor $2  --exhaustiveness 24
         echo >&1000
 }&
 done

wait

#整合结果保存至/fsx/docking_scores/$j-sum.txt
grep "VINA RESULT" /fsx/docking_result/$j/*.pdbqt| sort -k4 -r | sed -e 's/:REMARK VINA RESULT://g' > /fsx/docking_scores/$j-sum.txt

在分子对接中我们登录其中使用的一台c5.24xlarge的SPOP计算实例，通过htop指令我们看到所有的vCPU都被充分利用，CloudWatch也看到CPU利用率接近100%，最大化的利用计算资源。

总结

通过分析结果我们可到结果输出文件，每个分子对应一个对接分数，我们重点关注对接分数较低的分子。完成虚拟筛选之后后续的分子可用于实验验证活性，若有活性则获得该蛋白靶点的先导化合物（Hit）用于后续的药物开发。通过Amazon ParallelCluster 创建的高性能集群大规模调用海量计算资源，轻松在几个小时内完成百万级别的分子对接任务，在加速药物研发周期早期发挥极为重要的作用。

亚马逊AWS官方博客