autodock vina后处理分析

• 拆分对接结果vina_split --input result.pdbqt --ligand complex/lig
  

• 使用mv命令批量修改文件名，把01-09修改成1-9，便于批量处理

for i in `seq 1 9`; do
> mv "lig0${i}.pdbqt" "lig${i}.pdbqt"
> done

• 使用Openbabel把pdbqt转成pdb

for i in `seq 1 20`; do
> obabel -ipdbqt lig${i}.pdbqt -opdb -O lig${i}.pdb
> done

• 使用cat命令把ligand.pdb和receptor.pdb合成为complex.pdb

for i in `seq 1 20`; do
> cat ../ev71_none_2mer.pdb lig${i}.pdb > complex${i}.pdb
> done

• 使用freesasa计算接触面积
  recptor: freesasa ev71_none_2mer.pdb -n 100 --depth=residue -o sasa/rec

ligand:

for i in `seq 1 20`; do
>  freesasa complex/lig${i}.pdb -n 100 -H --depth=residue -o sasa/lig${i}
> done

complex:

for i in `seq 1 20`; do
>  freesasa complex/complex${i}.pdb -n 100 -H --depth=residue -o sasa/complex${i}
> done

接触面积的计算公式为：(rec_sasa + lig_sasa - complex_sasa) / 2

# 切换到tag工作文文件夹
os.chdir('/data5_large/home/xyli/enterovirus/ev71-tag/docking/none/vina')

# 配体ev71的sasa
with open('./sasa/rec', 'r') as rec_log:
    text = rec_log.readlines()
    total = text[15]
rec_sasa = eval(total[10:].lstrip())

inter_sasa_ls = []
for i in range(1,21):
    # 依次读取每个lig的sasa
    lig_file = 'lig%d'%i
    with open(f'./sasa/{lig_file}', 'r') as lig_log:
        text = lig_log.readlines()
        total = text[15]
    lig_sasa = eval(total[10:].lstrip())

    # 依次读取每个复合物的sasa
    complex_file = 'complex%d'%i
    with open(f'./sasa/{complex_file}', 'r') as complex_log:
        text = complex_log.readlines()
        total = text[15]
    complex_sasa = eval(total[10:].lstrip())

    # 计算binding area
    inter_sasa = (rec_sasa+lig_sasa-complex_sasa)/2
    inter_sasa_ls.append(inter_sasa)