marrink-lab / polyply_1.0

Generate input parameters and coordinates for atomistic and coarse-grained simulations of polymers, ssDNA, and carbohydrates
Apache License 2.0
126 stars 22 forks source link

KeyError: 0 when using gen_itp #217

Closed dim-99 closed 2 years ago

dim-99 commented 2 years ago

Hi,

I have been following issue #205 trying to pegylate a protein. After using,

polyply gen_seq -f lysozyme.itp -from_file protein:molecule_0 -from_string polymer:50:1:PEO-1.0 end:1:1:OHend-1.0 -seq protein polymer end -connects 0:1:0-0 1:2:49-0 -o sequence.json -name test -label 0:"from_itp":"molecule_0-1."

I did polyply gen_itp -lib martini3 -f lysozyme.itp PEO.martini.3b.itp OH_end.itp all_links.ff -seqf sequence.json -o lysoPEG.itp -name lysoPEG

and gives me the error

INFO - step - mapping sequence to molecule Traceback (most recent call last): File "/apps/polyply/1.3.0/bin/polyply", line 256, in main() File "/apps/polyply/1.3.0/bin/polyply", line 252, in main args.func(args) File "/apps/polyply/1.3.0/lib/python3.9/site-packages/polyply/src/gen_itp.py", line 77, in gen_params meta_molecule = MapToMolecule(force_field).run_molecule(meta_molecule) File "/apps/polyply/1.3.0/lib/python3.9/site-packages/polyply/src/map_to_molecule.py", line 292, in run_molecule new_molecule = self.add_blocks(meta_molecule) File "/apps/polyply/1.3.0/lib/python3.9/site-packages/polyply/src/map_to_molecule.py", line 196, in add_blocks fragment_nodes = list(self.fragments[self.node_to_fragment[start_node]]) KeyError: 0

Could you let me know what I'm doing wrong here?

Thank you !

fgrunewald commented 2 years ago

Hi @dim-99,

Thanks for your question. The reason the command is not working is, because your mixing force-field libraries. With the current PEGylated protein tutorial, which works with the Martini3 beta force-field, you don't need the specify the -lib martini3 flag. All input parameters needed are specified via the -f flag. If you just copy the command as presented in issue #205 it should be working. Let me know if any other issues arise with that command.

Cheers!

dim-99 commented 2 years ago

Hi @fgrunewald,

Thank you for the reply I got through the error. But similar to issue #212 I have the error,

WARNING - general - Node 634 has no resid. Setting resid to 634 + 1. WARNING - general - Node 635 has no resid. Setting resid to 635 + 1. WARNING - general - Node 636 has no resid. Setting resid to 636 + 1. WARNING - general - Node 637 has no resid. Setting resid to 637 + 1. INFO - step - mapping sequence to molecule INFO - step - applying links between residues 100%|██████████████████████████████████████████████████████████████████████████████████████████████████████████████████████| 4/4 [00:00<00:00, 73.58it/s] WARNING - general - You molecule consists of 3 disjoint parts. Perhaps links were not applied correctly.

I'm using polyply/1.3 and all my links are in the .ff file. ( polyply gen_itp -f ../molecule_0.itp all_links.ff PEO.martini.3b.itp OH_end.itp -seqf sequence.json -o lysoPEG.itp -name lysoPEG )

Is there a reason or a solution to get through this?

Thanks again !

fgrunewald commented 2 years ago

Hi,

Are you using the correct links file? In the repro it's called combined_links.ff there is no all_links.ff file. Otherwise, can you post the content of the all_links.ff

Cheers!

fgrunewald commented 2 years ago

Hi @dim-99,

I think your response got lost somehow. I'm reposting it here:

It's same as combined_links. I'm attaching MEE linker to 541 LYS. My protein has 587 residues, so MEE linker has the resid 588, PEG polymer (with 50 residues) starts at 589 the OH end is at 639.

[ link ]
[ molmeta ]
by_atom_id true
[ bonds ]
638 639 1 0.28 7000 {"group": "link"}
541 588 1 0.41 2000 {"group": "link"}
588 589 1 0.39 5000 {"group": "link"}
[ angles ]
637 638 639 2 140 25 {"group": "link"}
589 588 541 2 170 50 {"group": "link"}
588 541 542 2 150 15 {"group": "link"}

I think the confusion is that the combined_links.ff should contain atom indices not resids. So if you change the the above input file to the actual indices it should work. Can you check and let me know?

Cheers, Fabian

dim-99 commented 2 years ago

Thank you Fabian. The commands worked well for the system I tried. But my actual system is a dimer and I want to attach PEG to a single residue is my dimer. Two monomers in my dimer are identical.

I used the same commands. Protein A and B are the two monomers. I only want to attach the polymer to protein A.

polyply gen_seq -f molecule_0.itp molecule_1.itp -from_file ProteinA:proteinA ProteinB:proteinB -from_string linker:1:1:MEE-1.0 polymer:50:1:PEO-1.0 end:1:1:OHend-1.0 -seq ProteinA ProteinB linker polymer end -connects 0:2:540-0 2:3:0-0 3:4:49-0 -o sequence.json -name test -label 0:"from_itp":"proteinA-1" 1:"from_itp":"proteinB-1." (runs smoothly)

polyply gen_itp -f molecule_0.itp molecule_1.itp MEE.itp all_links.ff PEO.martini.3b.itp OH_end.itp -seqf sequence.json -o goxpegall.itp -name goxpeg

My all_links.ff is identical to the combined_links.ff

[ link ] [ molmeta ] by_atom_id true [ bonds ] 2615 2616 1 0.28 7000 {"group": "link"} 1181 2565 1 0.41 2000 {"group": "link"} 2565 2566 1 0.39 5000 {"group": "link"} [ angles ] 2614 2615 2616 2 140 25 {"group": "link"} 2566 2565 1181 2 170 50 {"group": "link"} 2565 1181 1180 2 150 15 {"group": "link"}

Gives me the error,

Traceback (most recent call last): File "/apps/polyply/1.3.0/bin/polyply", line 256, in main() File "/apps/polyply/1.3.0/bin/polyply", line 252, in main args.func(args) File "/apps/polyply/1.3.0/lib/python3.9/site-packages/polyply/src/gen_itp.py", line 79, in gen_params meta_molecule = ApplyLinks().run_molecule(meta_molecule) File "/apps/polyply/1.3.0/lib/python3.9/site-packages/polyply/src/apply_links.py", line 452, in run_molecule self.apply_link_between_residues(meta_molecule, link, link_node_to_resid) File "/apps/polyply/1.3.0/lib/python3.9/site-packages/polyply/src/apply_links.py", line 369, in apply_link_between_residues link_to_mol = match_link_and_residue_atoms(meta_molecule, File "/apps/polyply/1.3.0/lib/python3.9/site-packages/polyply/src/apply_links.py", line 228, in match_link_and_residue_atoms resid = block.nodes[list(block.nodes)[0]]["resid"] IndexError: list index out of range

I have a hard time tracking the error for the dimer system.

Thanks,

Best, Dimuthu

fgrunewald commented 2 years ago

Hi @dim-99,

I'm not 100% sure what is happening here. Protein itp files can have many issues at the CG level. If I understand you correctly though you want to attach PEG to only one of your protein dimers? Can you try generating the graph and running polyply only with one protein copy? So omit proteinB in your case and just generate the tip with proteinA. Unless the dimer has to be in a single itp-file that should work. Let me know if that solves the issues.

I'll also try to have a look what's going on with the current strategy but don't promise any fast results.

Cheers, Fabian