化合物から水素原子を引き抜いたラジカルの組み合わせ

ある化合物がCYPにより水酸化されて代謝されると分かっている場合に、全ての水素原子の乖離エネルギーを求めてラジカル生成エネルギーの低い順に代謝部位のあたりをつけたりしますが、毎回手で水素を消去した組み合わせを作るのは面倒なのでopenbabelでやってみます。

molオブジェクトのコピーをとるメソッドが見つからなかったのでcloneという関数を定義してます。

import openbabel as ob
import sys

def clone(mol):
   new_mol = ob.OBMol()
   for atom in ob.OBMolAtomIter(mol):
       new_mol.AddAtom(atom)

   for bond in ob.OBMolBondIter(mol):
       new_mol.AddBond(bond)
   return new_mol

def abs_h(mol):
   radicals = []
   mol.AddHydrogens()
   h_indexs = [atom.GetIdx() for atom in ob.OBMolAtomIter(mol) if
atom.GetType() == 'H']

   for h_index in h_indexs:
       new_mol = clone(mol)
       new_mol.DeleteAtom(new_mol.GetAtom(h_index))
       new_mol.SetTotalSpinMultiplicity(2)
       new_mol.SetTotalCharge(0)
       radicals.append(new_mol)
   return radicals


if __name__ == '__main__':
   input = sys.argv[1]

   obc = ob.OBConversion()
   obc.SetInAndOutFormats("smi", "mol")

   mol = ob.OBMol()
   next = obc.ReadFile(mol,input)

   while next:
       radicals = abs_h(mol)
       for r in radicals:
           print obc.WriteString(r)

       mol = ob.OBMol()
       next = obc.Read(mol)

あとはGamessでエネルギー計算すればいいですね。