VuisterLab / cing

Automated Validation of NMR Structures
http://nmr.le.ac.uk
2 stars 4 forks source link

Need to interpret the phenyl ring protons without stereospecificity for distances and dihedral violation analysis #231

Closed GoogleCodeExporter closed 9 years ago

GoogleCodeExporter commented 9 years ago
What steps will reproduce the problem?
1. Verify buguous reporting of phenyl ring proton distance violation:
http://nmr.cmbi.ru.nl/NRG-
CING/data/kn/2knr/2knr.cing/2knr/HTML/Restraints/distance_constraint_list.html#_
o9466

771 A PHE37 HD1 A LEU96 HA  1.8 5.0 4.3+-1.0    0.2+-0.8    3.7 1   
violMax: 3.67
violSd: 0.82

What is the expected output? What do you see instead?
The restraint from:
assign ( resid   37 and name HD1  )   ( resid   96 and name HA   )   4.00  2.20 
 1.00

Measured in Pymol to NOT violate in any model.

Please use labels and text to provide additional information.

Original issue reported on code.google.com by jurge...@gmail.com on 15 Mar 2010 at 1:41

GoogleCodeExporter commented 9 years ago
By the way in way of cleanup it would help if the bits of code that you added 
that are not essential to the CCPN part (i.e. the report at the end) is a 
separate method that can be deactivated.

Original comment by wfvran...@gmail.com on 10 Mar 2011 at 3:19

GoogleCodeExporter commented 9 years ago
For the violation, the offending restraint is not the one to QE, but to QD of 
same residue:

<ccp.nmr.NmrConstraint.DistanceConstraint [1, 1, 1435]>
  AL [(<ccp.molecule.MolSystem.Atom ['2knr', 'A', 7, 'HD3']>, <ccp.molecule.MolSystem.Atom ['2knr', 'A', 7, 'HD2']>)]
  OAL [(<ccp.molecule.MolSystem.Atom ['2knr', 'A', 34, 'HB3']>,), (<ccp.molecule.MolSystem.Atom ['2knr', 'A', 34, 'HB2']>,)]
  AL [(<ccp.molecule.MolSystem.Atom ['2knr', 'A', 7, 'HD3']>, <ccp.molecule.MolSystem.Atom ['2knr', 'A', 7, 'HD2']>)]
  OAL [(<ccp.molecule.MolSystem.Atom ['2knr', 'A', 34, 'HB3']>,), (<ccp.molecule.MolSystem.Atom ['2knr', 'A', 34, 'HB2']>,)]
   original
     0 2.08909106409
     1 2.09171296735
     2 1.68766185212
    -> VIOL 0.112338147881 0.0126198594693
     3 2.21512863978
     4 2.16337194768
     5 3.1644975645
     6 2.23890958008
     7 2.02327855365
     8 2.57224557731
     9 1.54235434907
    -> VIOL 0.25764565093 0.0663812814431
     10 1.91281564654
     11 2.28623973062
     12 1.90668409911
     13 1.97041321958
     14 2.52222262374
     15 2.24463698531
     16 3.14100174541
     17 1.67649558228
    -> VIOL 0.123504417716 0.0152533411953
     18 1.76416020875
    -> VIOL 0.0358397912547 0.00128449063718
     19 2.93487485941

So it's the violation in model 9 that ends up as the highest one.

Original comment by wfvran...@gmail.com on 10 Mar 2011 at 4:04

GoogleCodeExporter commented 9 years ago
Wrt previous comment 52

>-2- The max org violation 0.258 found by FC in 2knr differs from Wattos found 
0.093 Angstrom for:
FC found the 0.258 violation in model 9 (counting from 0) for restraint 1435 
which in star reads:
       1435 2 . 1 . 1 1   7   7 LYS HD2  H . . . .   4 . HD*  rr_2knr 1 
       1435 2 . 2 . 1 1  34  34 GLU HB3  H . . . .  31 . HB1  rr_2knr 1 
       1435 3 . 1 . 1 1   7   7 LYS HD3  H . . . .   4 . HD*  rr_2knr 1 
       1435 3 . 2 . 1 1  34  34 GLU HB3  H . . . .  31 . HB1  rr_2knr 1 
In the original author data this restraint was:
assign ( resid    4 and name HD*  )   ( resid   31 and name HB1  )   4.00  2.20 
 1.00
Instead of 2 distances, FC seems to average over too many distances (4 or 8?) 
which causes a too high violation
on the lower bounds to be reported.

Original comment by jurge...@gmail.com on 11 Mar 2011 at 10:18

GoogleCodeExporter commented 9 years ago
Wrt comment 49

-1- ambi count off for 1brv A 19 CYS HB*. Wattos got 4 out of 22 you get 2:
Wim modified and for 1brv all ambi counts match Wattos. For 2knr FC seems to 
have more ambies for a few triplets, i.e.
  1   48   GLY QA         124 no    100.0    97.0    0.4    0.4    0.0   6 2 no      0.267   0   0
seems to have 2 ambi restraints by FC but manually I see none in the 6 original 
restraints:
assign ( resid   44 and name HA   )   ( resid   45 and name HA1  )   4.00  2.20 
 1.00
assign ( resid   44 and name HN   )   ( resid   45 and name HA1  )   4.00  2.20 
 2.00
assign ( resid   44 and name HN   )   ( resid   45 and name HA2  )   4.00  2.20 
 1.00
assign ( resid   42 and name HN   )   ( resid   45 and name HA2  )   4.00  2.20 
 2.00
assign ( resid   45 and name HA2  )   ( resid   46 and name HG*  )   4.00  2.20 
 2.00
assign ( resid   45 and name HA2  )   ( resid   46 and name HD*  )   4.00  2.20 
 1.00

Original comment by jurge...@gmail.com on 11 Mar 2011 at 11:13

GoogleCodeExporter commented 9 years ago
Wrt comment 49

-3- Can you add the other DNA/RNA triplets like in entry 1a4d. Currently I only 
see Q5's.
Wim got list of additional ones. Output matches Wattos exactly except for 
trivial differences in rounding.

-4- Can you split the Phe/Tyr QR like Wattos has them?
Wim did this. Nice!

Bugs: see attached logs.

-A- For 1a24 (1c2n was completely fixed!)
Error changed to:
  File "/Users/jd/workspace35/cing/python/cing/Scripts/FC/utils.py", line 282, in <module>
    if fcProcessEntry( entry_code, ccpnTgzFile, outputCcpnTgzFile, functionToRun ):
  File "/Users/jd/workspace35/cing/python/cing/Scripts/FC/utils.py", line 195, in fcProcessEntry
    swapCheck(nmrConstraintStore, structureEnsemble)
  File "/Users/jd/workspace35/cing/python/cing/Scripts/FC/utils.py", line 125, in swapCheck
    swapCheck.checkSwapsAndClean()
  File "/Users/jd/workspace35/ccpn/python/ccpnmr/format/process/stereoAssignmentSwap.py", line 1317, in checkSwapsAndClean
    resSet.addAtomSet(otherAtomSet)
  File "/Users/jd/workspace35/ccpn/python/ccp/api/nmr/NmrConstraint.py", line 14942, in addAtomSet
    + ": %s" % (self,)
ApiError: ccp.nmr.NmrConstraint.FixedResonanceSet.addAtomSet:
       value is in list already: <ccp.nmr.NmrConstraint.FixedResonanceSet [1, 193]>

-B- For 1d3z
Nicely fixed!

Original comment by jurge...@gmail.com on 11 Mar 2011 at 11:25

GoogleCodeExporter commented 9 years ago
The problem is that this restraint is calculated twice because it's effectively 
treating each HD as an HD*, so counting the distances twice in this case (if 
you look at the original restraint it is effectively the case that it's HD* 
originally).

Anyway send me the latest version of your code and I'll add code to check for 
duplicated items, and not count them if so.

Original comment by wfvran...@gmail.com on 11 Mar 2011 at 12:52

GoogleCodeExporter commented 9 years ago
Thanks Wim! Not many changes attached.

Original comment by jurge...@gmail.com on 11 Mar 2011 at 1:17

Attachments:

GoogleCodeExporter commented 9 years ago
Wrt Comment 54, it took me a while to figure this one out, but for Arg 46 (or 
49) the HG2 and HG3 as well as the HD2 and HD3 are listed separately elsewhere, 
so in CCPN the assignments are considered stereo, and the restraints:

assign ( resid   45 and name HA2  )   ( resid   46 and name HG*  )   4.00  2.20 
 2.00
assign ( resid   45 and name HA2  )   ( resid   46 and name HD*  )   4.00  2.20 
 1.00

are in effect considered to each consist of two items, for the first one one to 
the HG2, one to HG3, for the second one to HD2, one to HD3.

Original comment by wfvran...@gmail.com on 11 Mar 2011 at 1:51

GoogleCodeExporter commented 9 years ago
And here's the updated code, so the violations should now be all correct.

Original comment by wfvran...@gmail.com on 11 Mar 2011 at 1:52

Attachments:

GoogleCodeExporter commented 9 years ago
Ok that violation is gone. However FC deassigns
  1   25   PHE QB          66 no    100.0    54.8    1.2    2.3    1.0  16 9 yes     2.266   9  18
but Wattos doesn't have the high violation over 1.0 Ang and maintains SSA.
       1  25 PHE QB  67 no  100.0  54.8  1.2  2.3  1.0 16  9 no  0.236   0   0 
Can you diff the attached for more examples?

Original comment by jurge...@gmail.com on 11 Mar 2011 at 2:25

Attachments:

GoogleCodeExporter commented 9 years ago
The violation comes from this constraint:

<ccp.nmr.NmrConstraint.DistanceConstraint [1, 1, 1898]>
  AL [(<ccp.molecule.MolSystem.Atom ['2knr', 'A', 20, 'HE2']>,)]
  OAL [(<ccp.molecule.MolSystem.Atom ['2knr', 'A', 25, 'HB2']>,), (<ccp.molecule.MolSystem.Atom ['2knr', 'A', 25, 'HB3']>,)]
   original
     0 5.93237642096
    -> VIOL 0.932376420963 0.869325790368
     1 6.21573125223
    -> VIOL 1.21573125223 1.47800247766
     2 5.78791585979
    -> VIOL 0.787915859789 0.620811402107
     3 5.92835769838
    -> VIOL 0.928357698385 0.861848016151
     4 7.26623685273
    -> VIOL 2.26623685273 5.13582947268
     5 5.05497220566
    -> VIOL 0.0549722056605 0.00302194339518
     6 6.03846644439
    -> VIOL 1.03846644439 1.07841255612
     7 6.00060671932
    -> VIOL 1.00060671932 1.00121380676
     8 5.03058326638
    -> VIOL 0.0305832663817 0.000935336182577
     9 5.67695490558
    -> VIOL 0.676954905581 0.45826794419
     10 6.19919172151
    -> VIOL 1.19919172151 1.43806078493
     11 5.62496551101
    -> VIOL 0.624965511005 0.390581889946
     12 5.79086867404
    -> VIOL 0.790868674042 0.625473259581
     13 6.048156744
    -> VIOL 1.048156744 1.09863255999
     14 5.75281896117
    -> VIOL 0.75281896117 0.566736388297
     15 6.0708887323
    -> VIOL 1.0708887323 1.14680267697
     16 5.78037559333
    -> VIOL 0.780375593333 0.60898606667
     17 6.0137737736
    -> VIOL 1.0137737736 1.02773726403
     18 5.72598995808
    -> VIOL 0.725989958077 0.527061419229
     19 6.21215397427
    -> VIOL 1.21215397427 1.46931725733

Original comment by wfvran...@gmail.com on 11 Mar 2011 at 5:16

GoogleCodeExporter commented 9 years ago
Anyway probably related to the changes in the Phe/Tyr I think. Too late now, 
off to my weekend! Almost there in any case...

Original comment by wfvran...@gmail.com on 11 Mar 2011 at 5:26

GoogleCodeExporter commented 9 years ago
I doublechecked and as far as I can see this restraint 1898 is one between the 
HE2 of an aromatic ring of a Phe (not deassigned) and the HB2 of another Phe, 
the distances should be correct. I don't think they're wrongly treated... you 
do handle the HE1/HE2 separately when they're separately listed?

Original comment by wfvran...@gmail.com on 14 Mar 2011 at 9:56

GoogleCodeExporter commented 9 years ago
[deleted comment]
GoogleCodeExporter commented 9 years ago
Damn, now I remember. Since depending on the chi-1 the nomenclature swaps 
already Wattos calculates both distances (either side of ring) and reports the 
lowest violation of the two for this situation. 

In this case too, it prevents that QB from being 'wrongly' deassigned. It's 
just nomenclature messing up. 

Can you special case for this problem?
Btw, I do the same for dihedrals, i.e.
http://code.google.com/p/nmrrestrntsgrid/issues/detail?id=212

Original comment by jurge...@gmail.com on 14 Mar 2011 at 10:30

GoogleCodeExporter commented 9 years ago
So for HD1/HD2 and HE1/HE2 for aromatics you always want to report the lowest 
violation?

Original comment by wfvran...@gmail.com on 14 Mar 2011 at 10:44

GoogleCodeExporter commented 9 years ago
Yes, just for phe/tyr because of their symmetry.

Original comment by jurge...@gmail.com on 14 Mar 2011 at 10:50

GoogleCodeExporter commented 9 years ago
And do you do that violation calculation over all models or for each model 
individually? So can the same constraint item 'swap' within a set models? 
Anyway hope this is the last one.

Original comment by wfvran...@gmail.com on 14 Mar 2011 at 11:20

GoogleCodeExporter commented 9 years ago
For each model individually. Thanks Wim!

Original comment by jurge...@gmail.com on 14 Mar 2011 at 12:00

GoogleCodeExporter commented 9 years ago
OK so in that case what do you do if you have two items in one distance 
constraint, e.g.:

Ala 13 HB - Phe 22 HE1
Ala 13 HB - Phe 22 HE2

you make sure that if you're checking the first item and you use for model 1 
HE1 but for model 2 HE2, for the second item you use HE2 for model 1 and HE1 
for model 2?

Original comment by wfvran...@gmail.com on 14 Mar 2011 at 12:03

GoogleCodeExporter commented 9 years ago
If I understand correctly yes. Always aiming for lowest violation per model 
even.

Original comment by jurge...@gmail.com on 14 Mar 2011 at 12:13

GoogleCodeExporter commented 9 years ago
OK should be working as requested now, new version attached.

Original comment by wfvran...@gmail.com on 14 Mar 2011 at 2:19

Attachments:

GoogleCodeExporter commented 9 years ago
I need to do this in CING proper as well. Now still wrong at e.g.:
http://nmr.cmbi.ru.nl/NRG-CING/data/kn/2knr/2knr.cing/2knr/HTML/Molecule/A/PHE22
/index.html#_top

I'm looking into your new version Wim instead of Wattos:
       1  25 PHE QB  67 no  100.0  54.8  1.2  2.3  1.0 16  9 no  0.236   0   0 
you now get:
  1   25   PHE QB          66 no    100.0    99.9    0.5    0.5    0.0  16   9 no      0.055   0   0
so the max violation dropped to 0.055 whereas Wattos was finding 0.236.

According to above link to the NRG-CING residue page this can only be for 
restraint:
1897    A PHE17 HE2 A PHE22 HB2 1.8 5.0 5.9+-0.5    0.9+-0.5    2.3 18  violMax: 2.27 
violSd: 0.45
Now, I'll measure the related distances for them...

I get to reproduce your numbers for model 6:
6.0000  5.1177  5.0550  5.055   0.0549999999999997
I can debug Wattos.

Sorry it's taking so long to finalize man.

Original comment by jurge...@gmail.com on 14 Mar 2011 at 4:06

GoogleCodeExporter commented 9 years ago
Sorry but I changed my mind and suggest we go back to your code from 
comment 59 without the special treatment for the phe/tyr. The reason is 
that the violations might be masked but not corrected and then they will 
show up in CING and other software tools. It's not a big issue because we 
rarely see SSA phe/tyr ring protons. Wattos was never doing this either, I
was mistaken there.

For now I suggest we stick with this. I would like to commit our changes to 
stereoAssignmentSwap.py
to sf.net. Mostly your comment 59 version. Do you want me to tailor the 
reporting
further? Unfortunately, I can't commit:
cvs ci -m "Updating to comment 59 code from Wim documented in 
http://code.google.com/p/cing/issues/detail?id=231#c55" -l 
"/ccpn/python/ccpnmr/format/process/stereoAssignmentSwap.py"
    **** Access denied: Insufficient Karma (jurgenfd|ccpn/python/ccpnmr/format/process|)
Can you give me some Karma? I attach the code so you can do the commit already.

Can you still fix the comment 49 bug for 1a24?

Original comment by jurge...@gmail.com on 15 Mar 2011 at 12:32

Attachments:

GoogleCodeExporter commented 9 years ago
Fix done (plus another one I found in the test set), new code checked into SF, 
so just use that version. Hope this is fine now, that was painful at the end!

Original comment by wfvran...@gmail.com on 15 Mar 2011 at 2:52

GoogleCodeExporter commented 9 years ago
Looking good from this end as well. I'm closing this issue. 
Thanks for sticking with it!

Original comment by jurge...@gmail.com on 15 Mar 2011 at 3:24