Skip to content

Receptor Preparation error #393

@NitPo

Description

@NitPo

Hi i'm having a problem preparing a receptor with meeko i get this error while using polymer object to convert a cleaned pdb (using prody function below) gotten from PDB.

I suspect the error is due to the failed residue matching but i can't find anything about fixing it in github/meeko docs

No template matched for residue_key='A:836'

tried 6 templates for residue_key='A:836'excess_H_ok=False

LYS heavy_miss=4 heavy_excess=0 H_excess=[] bond_miss=set() bond_excess=set()

NLYS heavy_miss=4 heavy_excess=0 H_excess=[] bond_miss=set() bond_excess={4}

CLYS heavy_miss=5 heavy_excess=0 H_excess=[] bond_miss=set() bond_excess={1}

LYN heavy_miss=4 heavy_excess=0 H_excess=[] bond_miss=set() bond_excess=set()

NLYN heavy_miss=4 heavy_excess=0 H_excess=[] bond_miss=set() bond_excess={4}

CLYN heavy_miss=5 heavy_excess=0 H_excess=[] bond_miss=set() bond_excess={1}

matched with excess inter-residue bond(s): A:836

No template matched for residue_key='A:847'

tried 6 templates for residue_key='A:847'excess_H_ok=False

LYS heavy_miss=4 heavy_excess=0 H_excess=[] bond_miss=set() bond_excess=set()

NLYS heavy_miss=4 heavy_excess=0 H_excess=[] bond_miss=set() bond_excess={4}

CLYS heavy_miss=5 heavy_excess=0 H_excess=[] bond_miss=set() bond_excess={1}

LYN heavy_miss=4 heavy_excess=0 H_excess=[] bond_miss=set() bond_excess=set()

NLYN heavy_miss=4 heavy_excess=0 H_excess=[] bond_miss=set() bond_excess={4}

CLYN heavy_miss=5 heavy_excess=0 H_excess=[] bond_miss=set() bond_excess={1}

matched with excess inter-residue bond(s): A:847

No template matched for residue_key='A:848'

tried 3 templates for residue_key='A:848'excess_H_ok=False

ASN heavy_miss=3 heavy_excess=0 H_excess=[] bond_miss=set() bond_excess=set()

NASN heavy_miss=3 heavy_excess=0 H_excess=[] bond_miss=set() bond_excess={4}

CASN heavy_miss=4 heavy_excess=0 H_excess=[] bond_miss=set() bond_excess={1}

matched with excess inter-residue bond(s): A:848

No template matched for residue_key='A:893'

tried 6 templates for residue_key='A:893'excess_H_ok=False

GLU heavy_miss=4 heavy_excess=0 H_excess=[] bond_miss=set() bond_excess=set()

NGLU heavy_miss=4 heavy_excess=0 H_excess=[] bond_miss=set() bond_excess={4}

CGLU heavy_miss=5 heavy_excess=0 H_excess=[] bond_miss=set() bond_excess={1}

GLH heavy_miss=4 heavy_excess=0 H_excess=[] bond_miss=set() bond_excess=set()

NGLH heavy_miss=4 heavy_excess=0 H_excess=[] bond_miss=set() bond_excess={4}

CGLH heavy_miss=5 heavy_excess=0 H_excess=[] bond_miss=set() bond_excess={1}

matched with excess inter-residue bond(s): A:893

  • Template matching failed for: ['A:836', 'A:847', 'A:848', 'A:893'] Ignored due to allow_bad_res.

multiprocessing.pool.RemoteTraceback:

"""

Traceback (most recent call last):

File "/opt/conda/envs/Screener/lib/python3.9/multiprocessing/pool.py", line 125, in worker

result = (True, func(*args, **kwds))

File "/opt/conda/envs/Screener/lib/python3.9/multiprocessing/pool.py", line 51, in starmapstar

return list(itertools.starmap(args[0], args[1]))

File "/home/screener/./scripts/screener.py", line 456, in prepare_target

receptor=PDBQTReceptor.from_pdbqt_filename(file_path)

File "/opt/conda/envs/Screener/lib/python3.9/site-packages/meeko/receptor_pdbqt.py", line 127, in from_pdbqt_filename

receptor = cls(pdbqt_string, skip_typing)

File "/opt/conda/envs/Screener/lib/python3.9/site-packages/meeko/receptor_pdbqt.py", line 117, in init

self._atoms, self._atom_annotations = _read_receptor_pdbqt_string(pdbqt_string, skip_typing)

File "/opt/conda/envs/Screener/lib/python3.9/site-packages/meeko/receptor_pdbqt.py", line 76, in _read_receptor_pdbqt_string

atom_annotations[atom_property_definitions[atom_type]].append(idx)

KeyError: 'O'

#this function use prody (prody_fix function) to clean pdb then call the preparation function prepare_target

def process_targets(self,
        source_receptor_files:Iterable=None,
        from_PDB:Iterable=None,
        n_proc:int=None,
        hydrate:bool=False,
        preparation:Callable=prepare_target,
        flexres:dict=None,
        map_mode:str="o4e",
        charge_type="geisteger",
        config:dict=None,
        ph:float=7.0,
        backend:str="prody",
        *args,
        **kwargs
        ):

        if (source_receptor_files is None) and (from_PDB is None):
            raise Exception("Missing input receptor(s)")
        pass

        if n_proc is None:
            n_proc=mp.cpu_count()-2 #use max cores #remove -2 for HPC
        
        
        receptor_ids=set()
        receptor_files=set()

        if source_receptor_files is not None:
            receptor_ids=set(source_receptor_files) #get unique files
            receptor_files=set([str(path.join(self.source_dir,Id)) for Id in receptor_ids]) #find absolute path
        
        if from_PDB is not None:
            pdb_entries = Screener.from_protein_data_bank(self,from_PDB) #Download pdb files prom protetein databank and returns paths
            from_PDB=set(from_PDB)
            receptor_ids=receptor_ids | from_PDB #add downloaded protein ids and merge with source files using union operator other methods fail
            src_receptor_files=receptor_files | pdb_entries #add downloaded files paths as per ids using union operator other methods fail
           
        for receptor_id in receptor_ids:
            makedirs(path.join(self.targets_workpath,receptor_id),exist_ok=True)
        #clean up and add hydrogens to recepetor before preparation




        pdb_cleaner={"pdbfix":Screener.pdb_fix,"prody":Screener.prody_fix}
        PDB=pdb_cleaner.get(backend)
        if PDB is None:
            raise ValueError(f"Invalid backend\nSupported pdbfix,prody")
        receptor_files=[PDB(file,self.targets_workpath,Id,kwargs.get(Id)) for file,Id in zip(src_receptor_files,receptor_ids)]
        #add H using reduce doi.org/10.1006/jmbi.1998.2401
        H_receptor_files = [files.replace("_clean.pdb", "_H.pdb") for files in receptor_files]
        print("Adding H to receptor")
        phil_params = [
            ["--overwrite",
            "--quiet", 
            clean_pdb,
            f"approach=add",
            f"add_flip_movers=True",
            f"output.filename={H_file}"]
            for clean_pdb, H_file in zip(receptor_files, H_receptor_files)
        ]
        for params in phil_params:
            run_program(program_class=reduce2.Program,args=params) #da errore che significa che non si può usare in multiprocessing
        
        tasks = [(self,files,ids,hydrate) for files,ids in zip(receptor_files,receptor_ids)]  #calculate task paramenters
        
        
        start = time.perf_counter()
        with Pool(n_proc) as pool:
            print("Start target preparation")
            results=pool.starmap(preparation, tasks)

        #sostituire questo con dataframe in shared memory
        #print(results)
        data={}
        for column in results[0].keys():
            data[column] = tuple(result[column] for result in results)

        self.targets_data=pd.DataFrame(data) #
        end = time.perf_counter()
        
        print(f"Terminated in {end-start:.3f}s,avg:{(end-start)/len(results):.3f} s/receptor")

```preparation function:

    def prepare_target(self,
        file_path:str,
        receptor_id:str,
        hydrate:bool=True,
        charge_type="geisteger",
        flexres:list=None,
        config:dict=None,
        backend:str="prody",
        ph:float=7.0):


       
        
        if config is not None:
            preparation_config=config.get("preparation_config")
            blunt_ends=config.get("blunt_ends")
            wanted_altloc=config.get("wanted_altloc")
            delete_residues=config.get("delete_residues")
            allow_bad_res=config.get("allow_bad_res")
            set_template=config.get("set_template")
            charge_type=config.get("charge_type")
        file_format=file_path.split(".")[-1]


        #fare in modo che reduce non dia più errore trasferendo questo il pezzo della scelta del converter  in process

        
        
        receptor_pdb_path = path.join(self.targets_workpath,receptor_id,f"{receptor_id}_H.pdb")
        outpath=path.join(self.targets_workpath,receptor_id,f"{receptor_id}.pdbqt")
        
        
        
        if not path.exists(outpath): #checks if pdbqt is already prepared, if yes skips it
            

            #set target preparation parameters
            templates = ResidueChemTemplates.create_from_defaults() #create from defaults for now
            if config is not None:
                mk_prep = MoleculePreparation.from_config(config)
            else:
                mk_prep = MoleculePreparation(hydrate=hydrate)
        
            #load pdb with H
            if backend.lower() == "prody":
                target=self.prody_load(receptor_pdb_path)
                polymer=Polymer.from_prody(
                target,
                templates,  # residue_templates, padders, ambiguous,
                mk_prep,
                #set_template,
                #delete_residues,
                allow_bad_res=True,
                #blunt_ends=True,
                #wanted_altloc=wanted_altloc,
                default_altloc="A"
            )
        
            elif backend.lower() == "file":
                with open(receptor_pdb_path,"r") as pdb_file:
                    pdb_string=pdb_file.read()
                
                polymer=Polymer.from_pdb_string(
                pdb_string,
                templates,  # residue_templates, padders, ambiguous,
                mk_prep,
                #set_template,
                #delete_residues,
                #allow_bad_res,
                #blunt_ends=blunt_ends,
                #wanted_altloc=wanted_altloc,
                #default_altloc=args.default_altloc,
            )
            else:
                raise ValueError(f"Invalid back end:{backend}")


            #pdb_string=mmCIF_to_pdb()
           
            
            
            pdbqt_tuple = PDBQTWriterLegacy.write_from_polymer(polymer)
            Screener.write_pdbqt(dir=self.targets_workpath,filename=f"{receptor_id}",pdbqt=pdbqt_tuple)
            receptor=PDBQTReceptor.from_pdbqt_filename(file_path)
        else:
            receptor=PDBQTReceptor.from_pdbqt_filename(file_path)
        
        receptor.assign_type_chrages(charge_type)
        

        if flexres is not None:
            flex_residues = set()
            for string in flexres:
                for res_id in parse_residue_string(string):
                    flex_residues.add(res_id)
            for res_id in flex_residues:
                receptor.flexibilize_sidechain(res_id, mk_prep)
        
        

        pdbqt_tuple = PDBQTWriterLegacy.write_from_polymer(receptor)
        rigid_pdbqt, flex_pdbqt_dict = pdbqt_tuple
        rigid=Screener.write_pdbqt(dir=self.target_workpath,filename=f"{receptor_id}_rigid",pdbqt=rigid_pdbqt)
        if flexres:
            flex=Screener.write_pdbqt(dir=self.target_workpath,filename=f"{receptor_id}_flex",pdbqt=flex_pdbqt)
        pass
pdb_fix function:

def prody_fix(filepath:str,
Dir:str,Id:str,hydrate:bool=False,prody_string:str=None,*args):
        if prody_string is None:
            prody_string="chain A not hetero"
        if hydrate:
            prody_string += "and water"
        base_pdb=parsePDB(filepath)
        clean_pdb=f"{Dir}/{Id}/{Id}_clean.pdb" 
        print(checkNonstandardResidues(base_pdb))
        if checkNonstandardResidues(base_pdb):
            print(checkNonstandardResidues(base_pdb))
            #remove std residue
            std_only=base_pdb.select(prody_string)
            writePDB(clean_pdb,std_only)
        else:
            protein=base_pdb.select(prody_string)
            writePDB(clean_pdb,protein)
        
        
        
        return clean_pdb
Incriminated PDBQT: http://pastebin.com/ntVuQrCP

source pdb: https://pastebin.com/ipS44fhV

Cleaned pdb (with H) : https://pastebin.com/4121Pnd1

Sorry for my bad grammar,Eng is not my first language.

Metadata

Metadata

Assignees

No one assigned

    Labels

    No labels
    No labels

    Type

    No type

    Projects

    No projects

    Milestone

    No milestone

    Relationships

    None yet

    Development

    No branches or pull requests

    Issue actions