В данном практикуме я создавала гомологичную модель комплекса лизоцима Meretrix lusoria (P86383) на основе структуры лизоцима радужной форели в комплексе с олигосахаридом (1LMP).

Подготовка Modeller

import modeller 
import _modeller
import modeller.automodel 

env=modeller.environ()
env.io.hetatm=True

Cоздание выравнивания только белков

alignm = modeller.alignment(env)
alignm.append(file='P86383.fasta', align_codes='all', alignment_format='FASTA')
mdl = modeller.model(env, file='1lmp.pdb', model_segment=('FIRST:'+'A', 'LAST:'+'A'))
alignm.append_model(mdl, atom_files='1lmp.pdb', align_codes='1lmp')
alignm.salign()
alignm.write(file='no_lig.ali', alignment_format='PIR')

Получилось следующее выравнивание:

>P1;sp|P86383|LYS_MERLU
sequence::: :: :::-1.00:-1.00
FAGGIVSQRCLSCICKMESGCRNVGCKMDMGSLSCGYFQIKEAYWIDCGRPGSSWKSCAASSYC--------ASL
CVQNYMKRYAKWAGCPLRCEGFAREHNGGPRGCKKGSTIGYWNRLQKISGCHGVQ*

>P1;1lmp
structureX:1lmp.pdb:1:A:+132:B:MOL_ID  1; MOLECULE  LYSOZYME; CHAIN  A; SYNONYM  MUCOPEPTIDE N-ACETYLMURAMYLHYDROLASE; EC  3.2.1.17:MOL_ID  1; ORGANISM_SCIENTIFIC  ONCORHYNCHUS MYKISS; ORGANISM_COMMON  RAINBOW TROUT; ORGANISM_TAXID  8022; ORGAN  KIDNEY: 2.00: 0.16
KVYDRCELARALKASGMDGYAGNSLPNWVCLSKWESSYNTQATNRNTDGSTDYGIFQINSRYWCDDGRTPGAKNV
CGIRCSQLLTDDLTVAIRCAKRVVLDPNGIGAWVAWRLHCQNQDLRSYVAGCGV*

В выравнивании отсутсвует лиганд, поэтому в построенной по такому выравнианию модели он также будет отсутствовать.

Учёт лиганда

Лиганд ("...", поскольку не белковой природы) был добавлен к последовательности белка вручную в качестве второй цепи. В случае необходимости это можно было бы автоматизировать, но для одной структуры легче так. Также был изменён интервал model_segment, чтобы он из структуры он также был прочитан.

alignm = modeller.alignment(env)
alignm.append_sequence('FAGGIVSQRCLSCICKMESGCRNVGCKMDMGSLSCGYFQIKEAYWIDCGRPGSSWKSCAASSYCASLCVQNYMKRYAKWAGCPLRCEGFAREHNGGPRGCKKGSTIGYWNRLQKISGCHGVQ/...')
mdl = modeller.model(env, file='1lmp.pdb', model_segment=('FIRST:'+'A', 'LAST:'+'B'))
alignm.append_model(mdl, atom_files='1lmp.pdb', align_codes='1lmp')
alignm[0].code = 'LYS_MERLU_lig'
alignm.salign()
alignm.write(file='all_in_one.ali', alignment_format='PIR')

Далее производилось моделирование но этому выравнианию.

## Выбираем объект для моделирования 
s = alignm[0]
pdb = alignm[1]

## Создаем объект automodel
a = modeller.automodel.automodel(env, alnfile='all_in_one.ali', knowns= pdb.code, sequence = s.code )

a.name='mod'+s.code
a.starting_model = 1
a.ending_model = 2
a.make()

Результат представлен на рисунке ниже. Хорошо видна разница в области делеции (петля в модели на месте продолжения спирали в исходной структуре).

Сдвиг лиганда

Подвинем лиганд, наложив ограничение на расстояние между парой атомов

class mymodel(modeller.automodel.automodel):
        
      def special_restraints(self, aln):
        rsr = self.restraints
        at = self.atoms
        
        # unpick all old restraints, connecting chain A and B
        for aA in self.chains['A'].atoms:
          for aB in self.chains['B'].atoms:
             rsr.unpick(aA, aB)

        
        # add new restraint
        for a1, a2 in [('NH2:50:A','O4:124:B')]:
            rsr.add(modeller.forms.gaussian(
                                      group=modeller.physical.xy_distance, 
                                      feature=modeller.features.distance(at[a1], at[a2]),
                                      mean=3.0, stdev=0.1))

a = mymodel(env, alnfile='all_in_one.ali', knowns= pdb.code, sequence = s.code )
a.name='move' + s.code
a.starting_model = 2
a.ending_model = 2
a.make()
        

Результат показан на рисунке ниже. Лиганд действительно сдвинулся, однако поскольку были разрушены все прочие ограничения (иначе новое не перевешивало и он оставался на месте), вместо ожидаемого небольшого сдвига вправо это дало сильный сдвиг вправо-вверх и изменением конформации Arg50.

Пептид как лиганд

Для того, чтобы заменить лиганд на аланиновый трипептид я просто вручную изменила моделруюмую последовательность, поставив "ААА" вместо "..." и переименовала три атома оригинального лиганда в CA, чтобы Modeller мог построить соответсвие.

alignm = modeller.alignment(env)
alignm.append_sequence('FAGGIVSQRCLSCICKMESGCRNVGCKMDMGSLSCGYFQIKEAYWIDCGRPGSSWKSCAASSYCASLCVQNYMKRYAKWAGCPLRCEGFAREHNGGPRGCKKGSTIGYWNRLQKISGCHGVQ/AAA')
mdl = modeller.model(env, file='1lmp_CA.pdb', model_segment=('FIRST:'+'A', 'LAST:'+'B'))
alignm.append_model(mdl, atom_files='1lmp_CA.pdb', align_codes='1lmp_CA')
alignm[0].code = 'LYS_MERLU_lig'
alignm.salign()
alignm.write(file='all_in_one_AAA.ali', alignment_format='PIR')

s = alignm[0]
pdb = alignm[1]

a = modeller.automodel.automodel(env, alnfile='all_in_one_AAA.ali', knowns= pdb.code, sequence = s.code )

a.name='mod'+s.code
a.starting_model = 3
a.ending_model = 3
a.make()

Результат представлен на рисунке ниже. Пептид действительно был расположен примерно на место оригинального лиганда. Скор (molpdf) такой модели составил 921 при скоре с оригинальным лигандом 1002 — на удивление небольшая разница. Возможно, она объясняется тем, что у структуры с оригинальным лигандом и так низкое качество в силу малого сходства последовательности белков (см. выравнивание).

Мутация

Для введения мутации я воспользовалась скриптом. Мутировала в исходном 1lmp остаток Val109 (гидрофобный на поверхности, не взаимодействует с лигандом) на Gln109 (гидрофильный, потенциально мог бы взаимодействовать с лигандом). Результат представлен на рисунке ниже. Взаимодействие с лигандом моделлер не смог предположить, поскольку лиганд непептидной природы рассматривает просто как твёрдое тело без моделирования взаимодействий помимо унаследованных из оригинального комплекса ограничений.