分析梳理--linux實現(xiàn)虛擬篩選

作者，Evil Genius

今天我們來實現(xiàn)虛擬篩選。

第一步，處理蛋白質(zhì)PDB文件，注意這時候的PDB分子已經(jīng)進行了去水處理（pymol），其中去除水分子的注意事項。

當(dāng)PDB文件中含有水分子時，必須將其刪除，否則會嚴重影響對接結(jié)果的準(zhǔn)確性。
水分子的存在會帶來兩個核心問題：
物理位阻：水分子會占據(jù)結(jié)合口袋的空間，導(dǎo)致配體分子無法正確進入。
能量計算錯誤：對接軟件會將水分子視為蛋白的一部分，產(chǎn)生虛假的氫鍵和排斥力。
正確處理流程
正確做法是：在生成PDBQT文件前，刪除所有水分子（通常以 HOH 或 WAT 標(biāo)識），但需保留可能對催化或結(jié)合至關(guān)重要的結(jié)構(gòu)性水分子（數(shù)量極少，通常≤3個）。
推薦流程：
使用PyMOL或Chimera等可視化軟件查看結(jié)構(gòu)，定位結(jié)合口袋。
刪除所有非關(guān)鍵水分子。
保留關(guān)鍵水分子：判斷依據(jù)可查看同源蛋白結(jié)構(gòu)文獻，或通過PyMOL觀察該水分子是否與蛋白/配體形成穩(wěn)定氫鍵網(wǎng)絡(luò)。
保存為新的PDB文件，再用于生成PDBQT。
去除水分子的前后對比
去水前

去水后

第一步：生成的蛋白質(zhì)PDB文件，加氫，加電荷，生成PDBQT

prepare_receptor.py -r protein_no_water.pdb -o protein.pdbqt -A hydrogens

但是我們最好聯(lián)合reduce使用，因為

特性 prepare_receptor.py -A hydrogens reduce
加氫范圍 添加所有氫原子（調(diào)用PyBabel） 添加所有氫原子，并優(yōu)化氫鍵網(wǎng)絡(luò)
側(cè)鏈優(yōu)化 無 自動翻轉(zhuǎn)Asn/Gln/His側(cè)鏈，形成最優(yōu)氫鍵網(wǎng)絡(luò)
輸出格式 直接輸出PDBQT（含原子類型和電荷） 輸出PDB格式，還需再轉(zhuǎn)PDBQT
適用場景 常規(guī)分子對接（一鍵完成） 高精度對接、分子動力學(xué)模擬前處理

所以完整的代碼是

reduce -FLIP protein_no_water.pdb > protein_h.pdb

# 轉(zhuǎn)換為 PDBQT 格式（必須使用 pythonsh）
prepare_receptor -r protein_h.pdb -o protein.pdbqt -A hydrogens

這個過程實現(xiàn)了

對蛋白質(zhì)分子的合并非極性氫（AutoDock/Vina 力場要求）
為所有原子分配正確的 AutoDock 原子類型（如 A、NA、OA 等）
計算并寫入 Gasteiger 部分電荷

第二步：對配體分子的處理

prepare_ligand.py  -l ligand.pdb -o ligand.pdbqt -A hydrogens

這個輸入的格式包括pdb、mol2、sdf 等格式。

批量處理也很簡單

for ligand in *.pdb; do
    name=$(basename "$ligand" .pdb)
    prepare_ligand.py -l "$ligand" -o "${name}.pdbqt" -A hydrogens
done

第三步：定義柔性殘基（還是以41位的組氨酸和145的半胱氨酸為例）

用法
prepare_flexreceptor.py -r receptor_filename -s list_of_names_of_residues_to_move

例子

prepare_flexreceptor.py -r protein.pdbqt -s HIS41_CYS145

完成柔性殘基的定義后，將分別生成rigid區(qū)域和flex區(qū)域的pdbqt文件:

proH_rigid.pdbqt 剛性部分
proH_flex.pdbqt 柔性側(cè)鏈部分

第四步：設(shè)置對接盒子

示例如下,contig.txt
receptor = receptor.pdbqt
center_x = 15.0
center_y = 20.0
center_z = 10.0
size_x = 20.0
size_y = 20.0
size_z = 20.0
exhaustiveness = 32
num_modes = 9
energy_range = 3.0

其中參數(shù)的意義

exhaustiveness = 32 (搜索 exhaustive-ness，即搜索徹底性)
直譯：窮舉性、徹底性。
作用：這是最重要的參數(shù)，控制著對接搜索的強度或努力程度。
具體含義：值越大，程序在配體-受體結(jié)合口袋中進行的全局搜索次數(shù)就越多。這能增加找到真正全局最優(yōu)結(jié)合構(gòu)象的概率，但同時會顯著增加計算時間（幾乎與設(shè)定值成正比）。
取值范圍與建議：
默認值：8。
快速篩選/測試：1 - 4（速度極快，用于排除明顯不結(jié)合的分子，但可能漏掉最優(yōu)解）。
常規(guī)對接：8 - 32（在精度和速度之間取得良好平衡。8是官方默認值，32已經(jīng)是相當(dāng)高的精度）。
高精度/最終確認：32 - 64或更高（用于已篩選出的幾個關(guān)鍵分子，追求最可靠的結(jié)果）。
核心原則：exhaustiveness 增加一倍，運行時間也大約增加一倍。建議先用小值（如4或8）快速測試流程，正式運行時再提高。

energy_range = 3.0 (Energy Range，即能量范圍)

直譯：能量范圍。
作用：與 num_modes 配合使用，用于控制輸出構(gòu)象之間的能量差異。
具體含義：程序只會輸出那些結(jié)合能（Affinity）與最優(yōu)構(gòu)象（排名第一）的差值 ≤ energy_range 的構(gòu)象。超過這個能量差的構(gòu)象，即使排名再靠前，也不會被輸出。
取值范圍與建議：
默認值：3.0 (單位：kcal/mol)。
原理：這是一個過濾機制。例如，如果最優(yōu)結(jié)合能為 -9.0 kcal/mol，設(shè)置 energy_range = 3.0，那么 Vina 只會輸出結(jié)合能優(yōu)于 -6.0 kcal/mol 的構(gòu)象。即使你設(shè)置了 num_modes = 9，但排名第9的結(jié)合能是 -5.8 kcal/mol（差值3.2 > 3.0），它也不會被輸出，實際輸出可能少于9個。
建議：
默認 3.0 即可，足以覆蓋大多數(shù)有意義的構(gòu)象變化。
如果只想看能量最低的幾個極近似構(gòu)象，可以減小到 1.0 或 1.5。
如果想看到更多樣化（但能量可能更高）的結(jié)合模式，可以增大到 4.0 或 5.0。

第五步：執(zhí)行分子對接

vina --config config.txt --receptor proH_rigid.pdbqt --flex proH_flex.pdbqt --ligand ligand.pdbqt --out out.pdbqt

批量對接呢？

#!/bin/bash

# 剛性受體和柔性殘基文件
RECEPTOR_RIGID="proH_rigid.pdbqt"
RECEPTOR_FLEX="proH_flex.pdbqt"
CONFIG_FILE="config.txt"

# 配體文件夾和輸出目錄
LIGAND_DIR="./ligands"
OUTPUT_BASE_DIR="./results"

# 創(chuàng)建輸出目錄
mkdir -p $OUTPUT_BASE_DIR

# 批量對接
for ligand_file in $LIGAND_DIR/*.pdbqt; do
    ligand_name=$(basename "$ligand_file" .pdbqt)
    output_dir="$OUTPUT_BASE_DIR/$ligand_name"
    mkdir -p $output_dir
    
    echo "正在對接: $ligand_name"
    
    vina --config $CONFIG_FILE \
         --receptor $RECEPTOR_RIGID \
         --flex $RECEPTOR_FLEX \
         --ligand $ligand_file \
         --out "$output_dir/out.pdbqt" \
         --log "$output_dir/log.txt"
done

echo "批量對接完成！"

其中使用autodock和vina力場的區(qū)別。

對比維度 AutoDock (特指AutoDock 4) AutoDock Vina
模型類型 基于物理的力場
(Physics-based Force Field)		 經(jīng)驗性 + 知識基礎(chǔ)
(Empirical + Knowledge-based)
能量項構(gòu)成 - 范德華力 (Lennard-Jones 12-6)
- 氫鍵 (方向性 12-10 勢)
- 靜電作用 (庫侖勢)
- 去溶劑化能 (基于電荷的成對項)
- 構(gòu)象熵罰 		- 范德華類勢 (引力高斯函數(shù) + 排斥項)
- 非定向氫鍵項
- 疏水項
- 構(gòu)象熵罰
特點 物理模型更完整，考慮了靜電和溶劑效應(yīng)，因此計算量更大、速度較慢 。 針對速度和準(zhǔn)確性進行了高度優(yōu)化，刪去了靜電和顯式去溶劑化項，用更簡單的經(jīng)驗項替代，使得速度比AutoDock快數(shù)十倍到上百倍

深入解讀：兩種力場的差異

雖然它們都能用來預(yù)測分子間結(jié)合模式，但背后的“方法論”完全不同：

  • AutoDock 4 系列：追求“還原”的物理模型
    AutoDock 4 的力場更像是一個精簡版的分子力學(xué)力場，試圖直接計算物理相互作用。它包含了靜電能和去溶劑化能，使得對結(jié)合自由能的預(yù)測在理論上更“完整”。這種方法的優(yōu)點是準(zhǔn)確性高，但代價是計算復(fù)雜，速度較慢。為了應(yīng)對特殊場景，它還發(fā)展出了專門針對含鋅金屬蛋白的 AutoDock4Zn 力場，以及針對多環(huán)化合物的處理方案。

  • AutoDock Vina：追求“效率”的經(jīng)驗?zāi)Ｐ?/strong>
    Vina的設(shè)計哲學(xué)從一開始就偏向虛擬篩選——即在短時間內(nèi)處理成千上萬個分子。為此，它巧妙地“舍棄”了計算最耗時的靜電和溶劑化項，轉(zhuǎn)而用幾個經(jīng)驗性的、擬合出來的高斯函數(shù)來近似描述范德華力和氫鍵作用。它獨特的疏水項也能很好地捕捉疏水效應(yīng)。這套組合拳讓Vina在保持高精度的同時，速度極快，成為目前應(yīng)用最廣泛的對接軟件之一

調(diào)用ad4力場

prepare_gpf.py -l ligand.pdbqt -r proH_rigid.pdbqt -y

autogrid4 -p proH_rigid.gpf -l proH_rigid.glg

vina --flex proH_flex.pdbqt \
     --ligand ligand.pdbqt \
     --maps proH_rigid \
     --scoring ad4 \
     --exhaustiveness 32 \
     --out out_ad4.pdbqt

生活很好，有你更好。

?著作權(quán)歸作者所有,轉(zhuǎn)載或內(nèi)容合作請聯(lián)系作者
【社區(qū)內(nèi)容提示】社區(qū)部分內(nèi)容疑似由AI輔助生成，瀏覽時請結(jié)合常識與多方信息審慎甄別。
平臺聲明：文章內(nèi)容（如有圖片或視頻亦包括在內(nèi)）由作者上傳并發(fā)布，文章內(nèi)容僅代表作者本人觀點，簡書系信息發(fā)布平臺，僅提供信息存儲服務(wù)。
禁止轉(zhuǎn)載，如需轉(zhuǎn)載請通過簡信或評論聯(lián)系作者。

相關(guān)閱讀更多精彩內(nèi)容

友情鏈接更多精彩內(nèi)容