tbereau / auto_martini

Automatic MARTINI parametrization of small organic molecules
GNU General Public License v2.0
61 stars 22 forks source link

Auto_martini not working #12

Closed sizeineb closed 6 months ago

sizeineb commented 6 years ago

Hi Tristan,

First of all, I would like to thank you for your program Auto_martini.

My problem is that I'm trying to use this script to coarse grain a small organic molecule: PK11195 (C21H21ClN2O) but it's not working.

Molecule: https://pubchem.ncbi.nlm.nih.gov/compound/pk_11195#section=Top

After executing this command line: python auto_martini.py --sdf PK11195.sdf --mol PK195 --verbose

I have hundreds of lines that are displayed in my terminal and it continues for more then 15min. I'm getting no errors and at the same time no clear final results.

My molecule contains 25 heavy atoms. From one of your previous comments, I read that the code currently only supports up to 25 heavy atoms so I think that it can handle my molecule.

I don't see what could be the problem! Any help, please?

tbereau commented 6 years ago

Adding the --verbose flag will output all the intermediate steps, so it's normal you get so much output.

25 heavy atoms will take a really long time. You might need to leave it running overnight. Another option would be to cut it down into two fragments and parametrize those. You then combine them by hand.

sizeineb commented 6 years ago

Hi Tristan,

Thank you very much for your fast reply.

I didn't know that calculations can take a lot of time. Meanwhile, I tried to cut my molecule into two fragments: Part 1: the chlorophenyl & the isoquinoline (cyclic part of the molecule) Part 2: the rest of the molecule

When I tried to coarse grain Part 1, I got the following error: [15:49:59] non-ring atom 1 marked aromatic Traceback (most recent call last): File "auto_martini.py", line 1226, in ringAtomsFlat, True) File "auto_martini.py", line 848, in printAtoms molFrag = genMoleculeSMI(smiFrag) File "auto_martini.py", line 93, in genMoleculeSMI mol = Chem.AddHs(mol) Boost.Python.ArgumentError: Python argument types in rdkit.Chem.rdmolops.AddHs(NoneType) did not match C++ signature: AddHs(RDKit::ROMol mol, bool explicitOnly=False, bool addCoords=False, boost::python::api::object onlyOnAtoms=None)

Coarse graining part 2 of the molecule worked well without any errors.

NB: these are the principle packages versions that I'm using: Name Version Build beautifulsoup4 4.6.0 py27h9416283_1 boost 1.56.0 py27_4 numpy 1.11.3 py27h8a80b8c_4
Python 2.7.14 h138c1fe_31 rdkit 2016.03.4 np111py27_1 rmg readline 7.0 hc1231fa_4 requests 2.18.4 py27h9b2b37c_1

I thought that the error that I'm getting is because one of the packages are outdated. I tried to update boost & rdkit using conda but I'm still getting the same error.

Do you have any idea, what could be the problem, please?

tbereau commented 6 years ago

Could you provide the smiles code of the molecule that fails? That would help to try to reproduce your error.

sizeineb commented 6 years ago

Thank you very much for your help and consideration, Tristan.

Here are the smiles code of the molecule that fails: Clc1c(c2nc(C(=O)N(C(CC)C)C)cc3c2cccc3)cccc1

And here is a 3D picture of the molecule: https://www.dropbox.com/s/rwod1xlftuvs8tu/PK11195_part2.png?dl=0

tbereau commented 6 years ago

I meant the fragment, not the whole molecule.

sizeineb commented 6 years ago

Yeah! This is just the fragment

tbereau commented 6 years ago

This is still too long. Break it further apart then. I wouldn’t recommend building a fragment bigger than ~12 heavy atoms

sizeineb commented 6 years ago

It's working. Thank you very much, Tristan.

david-martini commented 5 years ago

Hi every body. I read almost all of the comments, but I couldn't built my .itp file. There is any body to build .itp and .gro files for me. this is my smiles code (CC(C1CCC(C(O1)OC2C(CC(C(C2O)OC3C(C(C(CO3)(C)O)NC)O)N)N)N)NC). also I attached .sdf file. thanks a lot.

Structure2D_CID_3467 (1).zip

tbereau commented 5 years ago

Your molecule is too large for the algorithm. I recommend you break it into smaller fragments (2 should be fine) and then connect the fragments together.

Tristan

david-martini commented 5 years ago

Hi Tristan Thanks for your reply. As you recommended I break it into 2 smaller fragments. It woks, but I cant understand why it don't build bead with S prefix in ring parts of my molecule as suggested by martini force filed!!

Thanks

On Fri, Jun 21, 2019 at 8:57 PM Tristan Bereau notifications@github.com wrote:

Your molecule is too large for the algorithm. I recommend you break it into smaller fragments (2 should be fine) and then connect the fragments together.

Tristan

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/tbereau/auto_martini/issues/12?email_source=notifications&email_token=AMMY7WY3W4P42SF74B3BAXLP3T6OJA5CNFSM4E2F25AKYY3PNVWWK3TUL52HS4DFVREXG43VMVBW63LNMVXHJKTDN5WW2ZLOORPWSZGODYI53BI#issuecomment-504487301, or mute the thread https://github.com/notifications/unsubscribe-auth/AMMY7W3L5VJ3KYC6Z7JNIATP3T6OJANCNFSM4E2F25AA .

david-martini commented 5 years ago

Hi Dear Tristan

I attached my result. Is it valid?

Thank you so muche

On Fri, Jun 21, 2019 at 8:57 PM Tristan Bereau notifications@github.com wrote:

Your molecule is too large for the algorithm. I recommend you break it into smaller fragments (2 should be fine) and then connect the fragments together.

Tristan

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/tbereau/auto_martini/issues/12?email_source=notifications&email_token=AMMY7WY3W4P42SF74B3BAXLP3T6OJA5CNFSM4E2F25AKYY3PNVWWK3TUL52HS4DFVREXG43VMVBW63LNMVXHJKTDN5WW2ZLOORPWSZGODYI53BI#issuecomment-504487301, or mute the thread https://github.com/notifications/unsubscribe-auth/AMMY7W3L5VJ3KYC6Z7JNIATP3T6OJANCNFSM4E2F25AA .

david-martini commented 5 years ago

Hi Dear Tristan This is a reminder. Would you please guide me? Thanks

On Fri, Jun 21, 2019, 8:24 PM Ghorbi jamal ghorbi.sbmu@gmail.com wrote:

Hi Dear Tristan

I attached my result. Is it valid?

Thank you so muche

On Fri, Jun 21, 2019 at 8:57 PM Tristan Bereau notifications@github.com wrote:

Your molecule is too large for the algorithm. I recommend you break it into smaller fragments (2 should be fine) and then connect the fragments together.

Tristan

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/tbereau/auto_martini/issues/12?email_source=notifications&email_token=AMMY7WY3W4P42SF74B3BAXLP3T6OJA5CNFSM4E2F25AKYY3PNVWWK3TUL52HS4DFVREXG43VMVBW63LNMVXHJKTDN5WW2ZLOORPWSZGODYI53BI#issuecomment-504487301, or mute the thread https://github.com/notifications/unsubscribe-auth/AMMY7W3L5VJ3KYC6Z7JNIATP3T6OJANCNFSM4E2F25AA .

tbereau commented 5 years ago

Hi,

I'm not sure you did attach your results. Could you try coarse-graining a benzene ring (smiles string "c1ccccc1"), and see whether you get the expected result:

[atoms]
; id    type    resnr   residue  atom    cgnr    charge  smiles
    1     SC5     1     MOL       S01     1         0   ; c1ccccc1
    2     SC5     1     MOL       S02     2         0   ; c1ccccc1
    3     SC5     1     MOL       S03     3         0   ; c1ccccc1
david-martini commented 5 years ago

Hi Dear Tristan I coarse-grained Benzene and 5-phenylvaleric acid molecules. you can see results:

5-phenylvaleric acid

;;;; GENERATED WITH auto-martini ; INPUT SMILES: OC(=O)CCCCc1ccccc1 ; Tristan Bereau (2014)

[moleculetype] ; molname nrexcl PHE 2

[atoms] ; id type resnr residu atom cgnr charge smiles

1 P1 1 PHE P01 1 0 ; [C]C([O])=O

2 C5 1 PHE C01 2 0 ; [C][C]

3 SC4 1 PHE S01 3 0 ; [C]c1[c][c][c][c][c]1

4 SC5 1 PHE S02 4 0 ; [c]1[c][c][c][c][c]1

5 SC5 1 PHE S03 5 0 ; [c]1[c][c][c][c][c]1

[bonds] ; i j funct length force.c. 1 2 1 0.25 1250 2 3 1 0.31 1250

[constraints] ; i j funct length 3 4 1 0.24 3 5 1 0.24 4 5 1 0.24

[angles] ; i j k funct angle force.c. 1 2 3 2 119.9 25.0 2 3 4 2 124.2 25.0 2 3 5 2 115.5 25.0

[dihedrals] ; i j k l funct angle force.c. 2 3 4 5 2 -101.9 10.0 2 3 5 4 2 116.2 10.0

Benzene ; INPUT SMILES: c1ccccc1 [moleculetype] ; molname nrexcl BENZ 2

[atoms] ; id type resnr residu atom cgnr charge smiles

1 SC5 1 BENZ S01 1 0 ; [c]1[c][c][c][c][c]1

2 SC5 1 BENZ S02 2 0 ; [c]1[c][c][c][c][c]1

3 SC5 1 BENZ S03 3 0 ; [c]1[c][c][c][c][c]1

[constraints] ; i j funct length 1 2 1 0.24 1 3 1 0.24 2 3 1 0.24

Thanks

On Tue, Jun 25, 2019 at 10:06 AM Tristan Bereau notifications@github.com wrote:

Hi,

I'm not sure you did attach your results. Could you try coarse-graining a benzene ring (smiles string "c1ccccc1"), and see whether you get the expected result:

[atoms] ; id type resnr residue atom cgnr charge smiles 1 SC5 1 MOL S01 1 0 ; c1ccccc1 2 SC5 1 MOL S02 2 0 ; c1ccccc1 3 SC5 1 MOL S03 3 0 ; c1ccccc1

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/tbereau/auto_martini/issues/12?email_source=notifications&email_token=AMMY7W5MTWVGEWZRFH7H7LDP4GVFJA5CNFSM4E2F25AKYY3PNVWWK3TUL52HS4DFVREXG43VMVBW63LNMVXHJKTDN5WW2ZLOORPWSZGODYPB6RA#issuecomment-505290564, or mute the thread https://github.com/notifications/unsubscribe-auth/AMMY7W7QKJRFYE3YXBYSKVLP4GVFJANCNFSM4E2F25AA .

tbereau commented 5 years ago

That looks good! They both yield the expected parametrizations. Can you point to a more specific problem, possibly with a minimally-working example?

david-martini commented 5 years ago

Hi Tristan Thanks for your reply. As I said previously, I break my molecule into 2 smaller fragments. It works, but I can't understand why it don't build bead with S prefix in ring parts of my molecule as suggested by martini force filed!! Also, the beads in the ring parts seem strange. I attached my result. Is it valid? Thank you so much

*Fragment 1 GENERATED WITH auto-martini ; INPUT SMILES: O(C1C(O)CC(N)CC1N)C1OC(CCC1N)C(NC)C ; Tristan Bereau (2014)

[moleculetype] ; molname nrexcl GENT 2 [atoms] ; id type resnr residu atom cgnr charge smiles 1 P3 1 GENT P01 1 0 ; [O][C][C][O] 2 P1 1 GENT P02 2 0 ; [C][C][N] 3 P1 1 GENT P03 3 0 ; [C][C][N] 4 P1 1 GENT P04 4 0 ; [C]O[C] 5 C5 1 GENT C01 5 0 ; [C][C] 6 P2 1 GENT P05 6 0 ; [C][N] 7 Nd 1 GENT N01 7 0 ; [C][C][N][C]

[bonds] ; i j funct length force.c. 1 2 1 0.31 1250 1 3 1 0.24 1250 1 4 1 0.23 1250 2 3 1 0.26 1250 4 5 1 0.24 1250 4 6 1 0.24 1250 4 7 1 0.24 1250 5 6 1 0.25 1250

[angles] ; i j k funct angle force.c. 1 2 3 2 50.0 25.0 1 4 5 2 94.3 25.0 1 4 6 2 61.6 25.0 1 4 7 2 121.5 25.0 2 1 4 2 134.6 25.0 3 1 4 2 119.0 25.0 4 5 6 2 58.9 25.0 5 4 7 2 64.3 25.0 6 4 7 2 127.0 25.0

*Fragment 2 GENERATED WITH auto-martini ; INPUT SMILES: OC1OCC(O)(C(NC)C1O)C ; Tristan Bereau (2014)

[moleculetype] ; molname nrexcl GENT 2

[atoms] ; id type resnr residu atom cgnr charge smiles 1 Nda 1 GENT N01 1 0 ; [C]C[O] 2 P1 1 GENT P01 2 0 ; [C][N][C] 3 P4 1 GENT P02 3 0 ; [O][C]C[O]

[bonds] ; i j funct length force.c. 1 2 1 0.25 1250 1 3 1 0.26 1250 2 3 1 0.25 1250

[angles] ; i j k funct angle force.c. 1 2 3 2 61.7 25.0

tbereau commented 5 years ago

For fragment 2 I'm getting:

[atoms]
; id    type    resnr   residue  atom    cgnr    charge  smiles
    1     SP1     1     MOL       S01     1         0   ; OC1CCCCO1
    2     SNda    1     MOL       S02     2         0   ; CC1(O)CCCOC1
    3     P1      1     MOL       P01     3         0   ; CNC
    4     SP1     1     MOL       S03     4         0   ; OC1CCCOC1

Are you using the refactor branch?

david-martini commented 5 years ago

Hi Tristan Its really different!! No, I didn't use refactor branch. How to use it? Thanks

tbereau commented 5 years ago

You need to checkout that branch:

git checkout refactor

Best, Tristan

On Thu, Jun 27, 2019 at 9:44 AM david-martini notifications@github.com wrote:

Hi Tristan Its really different!! No, I didn't use refactor branch. How to use it? Thanks

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/tbereau/auto_martini/issues/12?email_source=notifications&email_token=AAGE2MXFZ2IPJYHITSC2HGLP4RVXXA5CNFSM4E2F25AKYY3PNVWWK3TUL52HS4DFVREXG43VMVBW63LNMVXHJKTDN5WW2ZLOORPWSZGODYWIA4I#issuecomment-506232945, or mute the thread https://github.com/notifications/unsubscribe-auth/AAGE2MRGWYHYP6MNXO7DRBDP4RVXXANCNFSM4E2F25AA .

david-martini commented 5 years ago

Hi I run this comment (git checkout refactor) in terminal, but I get error: fatal: Not a git repository (or any of the parent directories): .git

tbereau commented 5 years ago

You need to use the git repository, not the zip version. First clone the repo from the github page.

david-martini commented 5 years ago

Hi Tristan I clone repository from the github page. It works. Thanks for your help and thank auto-martini.

GENERATED WITH auto_martini.py ; None ; Tristan Bereau (2014)

[moleculetype] ; molname nrexcl GEN 2

[atoms] ; id type resnr residue atom cgnr charge smiles 1 SP1 1 GEN S01 1 0 ; [O][C]1[C][C][C][C]O1 2 SNda 1 GEN S02 2 0 ; [C]C1([O])[C][C][C]O[C]1 3 P1 1 GEN P01 3 0 ; [C][N][C] 4 SP1 1 GEN S03 4 0 ; [O][C]1[C][C][C]O[C]1

david-martini commented 5 years ago

Hi Dear Tristan

I clone repository from the github page and use refactor branch. Why auto_martini.py can't build the .map file!? I need map file. How do I build it?

Thanks

tbereau commented 5 years ago

Hi David,

it's simply not a feature that's implemented as of now, but you're welcome to contribute and push it to the repo!

Best, Tristan

david-martini commented 5 years ago

Hi Dear Tristan

Thanks for your reply. It will be my pleasure to assist you.

Best wishes, David

georopon commented 4 years ago

Hi i've this problem con NH elements azole.

WARNING:main:ALOGPS can't predict fragment: CNCNH+C

tbereau commented 4 years ago

ALOGPS is unfortunately not able to make a prediction for all chemical groups. In case this happens, your best bet is to maybe look up what the expected water/octanol partitioning is for that group, and manually assign the appropriate bead type.

georopon commented 4 years ago

fPREP option is best for this case?

tbereau commented 4 years ago

the --fpred option will use an atom-based prediction algorithm. The quality may be a lot lower. You can try it, but I would recommend you check the accuracy of the bead typing manually afterwards.