Resonance Assignment

Introduction

Assignment in the CCPN data model is different from that in most programs (and, we believe, better).  Instead of assigning peaks or chemical shifts directly to atoms, all connections go via the Resonance object that serves as a connector. In this way you can connect peaks, shifts and spin systems to each other before you know which atom is involved, and you can (re)assign all the peaks and shifts to different atoms in a single operation.

• • Resonances are data model objects serve as stand-ins for NMR signals. A Resonance corresponds to the signal from atom or group of atoms, but you do not need to know which group it is. NMR measurements, such as shifts, or derived quantities like pKAs all belong to one or more Resonances. Resonances also serve to assign peaks, to form spin systems groups, and to identify the origin of structural restraints.  The link from Resonances to the atoms of a molecular system is not direct, but passes through intermediate objects: ResonanceSets and AtomSets.

  • AtomSets represent group of atoms that are in fast exchange and therefore the whole set is assigned, rather than the atoms individually (e.g. CH3 groups).

  • ResonanceSets connect Resonances to AtomSets; when assigning prochiral groups, up to two Resonances may connect to up to two AtomSets.  The data model does not use pseudatoms; all Atoms are real and the assignment mechanism takes care of the ambiguities.

  • Shifts are objects that represent the measurement of chemical shift values from NMR data. Shifts are distinct from Resonances, but they link to a Resonance, so that a single Resonance can have many shift measurements, e.g. under different conditions. In the data model Shifts are just one kind of AbstractMeasurement; there are others kinds of measurement like JCoupling, RDC, T1, T2, HExchRate, etc. All such measurements are contained within measurement lists. For example Shifts are contained in ShiftLists and RDCs in RdcLists.  A ShiftList contains the set of chemical shifts that are relevant for a particular set of Experiments; each experiment can have at most one ShiftList (etc.).

  • PeakDims are objects that represent the different dimensions of a given Peak and it is these dimensions that are assigned to Resonances; e.g. an HSQC peak has two dimensions, one of which may be assigned to a 15N resonance and the other to a 1H resonance . Each peak dimension is linked to Resonances through one or more contributions (PeakDimContribs). The different contributions allow the assignment of a single dimension to multiple, distinct Resonances; in other words an ambiguous peak assignment, which would be common place on NOE spectra. Assignment contributions to a peak may be weighted and also grouped into subsets (PeakContrib object). The mechanism is quite complex since it has to accommodate multiple assignments, assignment of multiplet components (e.g. for resolved couplings) and projection spectra.

  • ResonanceGroups are used to group Resonances into spin systems and store partial assignments. ResonanceGroups correspond to Resonances that belong to the same residue. They can be annotated with a molecule type, residue type, residue linking type and residue descriptor (e.g. protonation state) and linked to a definite residue. This allows different levels of intermediate assignment information, from the vague to a specific residue. ResonanceGroups can also carry probabilistic information (e.g. the probability of being a particular kind of amino acid) and sequence linking information. Such linking allows ResonanceGroups to be connected into sequential segments before it is known which residues they correspond to.

Setup

The following examples assume that you have loaded the demonstration project, so that we can look at pre-existing data. If you have not already done so, make sure that you load such a CCPN project. This can be done as follows, remembering to change the location of the project directory to that it is appropriate to your system:

from memops.general.Io import loadProject

rootProject = loadProject('/home/user/myProjDirName')

Next we will find some objects to work with the below examples:

# There may be a pause on issuing this command

# while the NMR data is loaded from disk

nmrProject = rootProject.currentNmrProject

experiment = nmrProject.findFirstExperiment(name='HSQC')

spectrum = experiment.findFirstDataSource()

peakList = spectrum.findFirstPeakList()

molSystem = rootProject.findFirstMolSystem()

chain = molSystem.findFirstChain()

residue = chain.findFirstResidue()

Resonances

Resonances are contained directly within the NmrProject. New Resonances must have an isotopeCode; other attributes are optional or set automatically. Note that we are using a CCP project object that we either made earlier or loaded.

nmrProject = rootProject.currentNmrProject

resonance = nmrProject.newResonance(isotopeCode='13C')

Looping over all Resonances in a given NMR project is done in the usual way. The following code will print a tab-separated list of simple attributes for all resonances in serial number order:

print "Serial\tIsotopeCode\tName\tAssigned?"

for resonance in nmrProject.sortedResonances():

 print '%s\t%s\t%s\t%s' % (resonance.serial, resonance.isotopeCode,

                           resonance.name, bool(resonance.resonanceSet))

To find a chemical shift value for a Resonance, finding in this instance an arbitrary 1H resonance, we need to look at its linked Shift objects:

resonance = nmrProject.findFirstResonance(isotopeCode='1H')

shift = resonance.findFirstShift()

print shift.value

However, the above will find an abitrary Shift, not necessarily the one we want (i.e. relating to a given experiment or set of conditions). In general, to specify the correct chemical shift measurement we find a specific ShiftList and use only Shifts from that list. For example, if you need the chemical shift for a Resonance in a given Experiment, it is easy to start from the Experiment because the experiment refers to a specific ShiftList. Note that there can be only one shift (i.e. a measurement) in any ShiftList that corresponds to a given Resonance:

shiftList = experiment.shiftList

if shiftList is not None:

  shift = shiftList.findFirstMeasurement(resonance=resonance)

 print shift.value, shift.error

If you want all shifts no matter what the shift list, you do either of the following

setOfShifts = resonance.shifts

listOfShifts = resonance.sortedShifts()

Making a Resonance name that reflects the assignment state is not a simple job, given the number of situations you need to allow for (from partial annotation to full assignment). Fortunately the job has been mostly done for you by using a high-level functions imported from the Analysis library:

from ccpnmr.analysis.core import AssignmentBasic

simpleName = AssignmentBasic.getResonanceName(resonance)

print simpleName

# If you want a full, informative annotation for GUI display etc.

fullName = AssignmentBasic.makeResonanceGuiName(resonance)

print fullName

A Resonance may be assigned to atoms via an intermediate ResonanceSet assignment object; this may seem complex but is vital to express ambiguity (an alternative to using pseudoatoms). A ResonanceSet (the assignment object ) may link to more than one AtomSet and thence to Atoms, though they will in practice be from the same residue.

If the assignment is simple (e.g. protein CA) there will be only one atom set, but if the assignment is ambiguous (e.g Serine HB2/HB3) then there will be multiple atom sets. If you are looking for the atom, you need to decide how to handle these different cases (see below). Meanwhile this code will get you to the AtomSets and hence Atoms to which a Resonance may be assigned:

resonanceSet = resonance.resonanceSet

# If we have an assigment object

if resonanceSet is not None:

 atomSets = resonanceSet.atomSets

 for atomSet in atomSets:

   for atom in atomSet.atoms:

     print atom.name

Finding the assigned residue  is a simple matter of navigating further through the assignment objects, and it doesn't matter which of the linked AtomSets or Atoms we choose:

resonanceSet = resonance.resonanceSet

if resonanceSet is not None:

 residue = resonanceSet.findFirstAtomSet().findFirstAtom().residue

 print residue.seqCode, residue.ccpCode

Finally, thinking about assignment from the Peak (rather than from Atoms), you might want all peak dimensions (PeakDims) assigned to a given Resonance. The following code will work, just bear in mind that some of the peaks may have multiple assignments.

peakDims = [contrib.peakDim for contrib in resonance.peakDimContribs]

for peakDim in peakDims:

 print peakDim.value

Note that if you work within Analysis the values for the Resonance's chemical shift measurements are automatically calculated as a (weighted) average of these PPM values; i.e. derived from the positions of the assigned peak dimensions.

Chemical Shifts

The relationship between shifts, peaks, and resonances is shown in the figure:

Creating a new Shift is simple once you have a ShiftList object. The mandatory values are the shift value (in ppm) and the Resonance that has the chemical shift:

shiftList = nmrProject.newShiftList()

shift = shiftList.newShift(value=8.37, resonance=resonance)

print shift.value

Other attributes include the error (=the estimated standard deviation of the measured value) and the figure of merit. A ShiftList can be linked to several Experiments, that can be retrieved as

# for our new shift list

experiments = shiftList.experiments

print experiments

# for a previous, well-used shift list

shiftList = experiment.shiftList

for experiment in shiftList.experiments:

 print experiment.name

A ShiftList ultimately belongs to the NmrProject, just like Experiments. However, because a ShiftList is really a sub-type of MeasurementList to find a specific shift list from an NmrProject you must find measurement lists of the "ShiftList" type. Also, for the same reason, we access the shift lists measurements attribute to get at the chemical shifts; there is no "shifts" attribute.

shiftList = nmrProject.findFirstMeasurementList(className='ShiftList')

for shift in shiftList.measurements:

 print shift.value

Analysis automatically recalculates shift values and errors for Shifts so that they are a weighted average of the corresponding peak positions (an different spectra can have different weights). If you want to do the same you can use the Analysis utility function:

from ccpnmr.analysis.core import AssignmentBasic

AssignmentBasic.updateAllShifts(shiftList)

Going from a Shift object to Peaks with the same assignment is just a matter of getting the Resonance and getting the peaks from that:

peaks = set()

for peakDimContrib in resonance.peakDimContribs:

 peaks.add(peakDimContrib.peakDim.peak)

print peaks

Finding peaks with PPM positions rather than a Shift object is another matter. There happens to be no utility function that matches exactly, but we can make a function easily enough. As a more complex example, also illustrating how we can work with Python functions,  we can make a subroutine that will work:

def findMatchingPeaks(peakList, matchDim, ppm, tolerance):

 """ find list of peaks from peakList that match position

     within tolerance when comparing dimension matchDim

 """

 peaks = set()

 for peak in peakList.peaks:

   peakDim = peak.findFirstPeakDim(dim=matchDim)

   delta = abs(peakDim.value-ppm)

   if delta < tolerance:

     peaks.add(peak)

 return peaks   

# Run the function

findMatchingPeaks(peakList, 1, 8.00, 0.07)

Note that if we were working from inside Analysis we could alternatively make a function based upon the existing searchPeaks(), which will do the fetching of peaks and getDataDimFullShiftRange() which tells us how wide a region to look in. Using these utility functions will be quicker for large peak lists because they enable us to avoid looping through all of the peaks within the list. We cannot use searchPeaks() outside of Analysis because the function relies upon the fast C libraries that are tied to the programs internal peak representation.

Assigning Resonances to Peaks

The figure shows the links used for assigning peaks.

Peak assignment is a matter of creating a contribution (PeakDimContrib) from a given Resonance on the relevant peak dimension (PeakDim):

contrib = peakDim.newPeakDimContrib(resonance=resonance)

Note that issuing the above command may well fail, and if it does this is likely for a good reason; there are a number of checks that should be performed before making such a link. There are a number issues that one should consider when assigning a Peak dimension  to a Resonance:

• • Does the Resonance isotopeCode match?

  • Is the dimension already assigned to that Resonance?

  • Does the resonance correspond to an aliased position of the peak?

  • Does the resonance correspond to an existing assignment?

  • Does the assignment fit with the molecular system allowed for Experiment?

It is recommended to use the existing utility functions that take care of these (ans other) issues. So the following is the preferred way of assigning peak dimensions:

from ccpnmr.analysis.core import AssignmentBasic

contrib = AssignmentBasic.assignResToDim(peakDim, resonance)

Deassigning a peak is simple enough to be done without utility functions:

for peakDim in peak.peakDims:

 for contrib in peakDim.peakDimContribs:

   contrib.delete()

A more complex use of assignment is initialising a newly picked peak in a 15N HSQC (done in more detail in ccpnmr.analysis.popups.InitRootAssignments). Note that assignResToDim() will create a new Resonance object if there is none present. After we have made the new 1H and 15N resonances we add some preliminary assignment information: firstly we set the names of the atoms that the resonances would be assigned and secondly we put both resonances in a single ResonanceGroup (a "SpinSystem") that corresbonds to an as yet uknown residue.

from ccpnmr.analysis.core import AssignmentBasic

from ccpnmr.analysis.core.PeakBasic import pickPeak

# Make and position a new peak

peak = pickPeak(peakList, (7.123, 123.45), 'ppm', doFit=False)

peakDimH, peakDimN = peak.sortedPeakDims()

peakDimContrib = AssignmentBasic.assignResToDim(peakDimH)

resonanceH = peakDimContrib.resonance

peakDimContrib = AssignmentBasic.assignResToDim(peakDimN)

resonanceN = peakDimContrib.resonance

resonanceH.addAssignName('H')

resonanceN.addAssignName('N')

spinSystem = nmrProject.newResonanceGroup()

AssignmentBasic.addSpinSystemResonance(spinSystem, resonanceH)

AssignmentBasic.addSpinSystemResonance(spinSystem, resonanceN)

print AssignmentBasic.makePeakDimAnnotation(peakDimH)

print AssignmentBasic.makePeakDimAnnotation(peakDimN)

As a final example consider finding all 15N HSQC-NOESY peaks that match a 15N HSQC amide peak and transferring the assignments. Finding the best utility functions again requires some thought, but the following will work. Again the advantage of using utility functions is that issues like correction for aliasing, and correct matching of dimensions are already taken care of. Note that the below example assumes the "peakList" variable is still the previous 15N HSQC  peak list.

from ccpnmr.analysis.core import PeakBasic

peak = peakList.findFirstPeak()

# Find a destination peak list

noeExperiment = nmrProject.findFirstExperiment(name='3dNOESY')

noePeakList = noeExperiment.findFirstDataSource().findFirstPeakList()

peakDimH, peakDimN = peak.sortedPeakDims()

peaksFitH = set(PeakBasic.findMatchingPeaks(noePeakList, peakDimH))

peaksFitN = set(PeakBasic.findMatchingPeaks(noePeakList, peakDimN))

usePeaks = peaksFitH.intersection(peaksFitN)

AssignmentBasic.propagatePeakAssignments(usePeaks, refPeak=peak)

Assigning Resonances to Atoms

From the top, the figure shows:

• • A stereospecific atom, e.g. Ser HB2 or Ser HA:

  • one Resonance, one ResonanceSet, one AtomSet, one Atom.

  • A stereospecific group of atoms, e.g. The Ala CH3, HB*:

  • in one Resonance, one ResonanceSet, one AtomSet, many Atoms.

  • Non specific prochiral atoms, e.g. Ser HBa, HBb:

  • two Resonances, one ResonanceSet, two AtomSets, one Atom each.

  • If only one resonance is observed the other Resonance object is not created.

  • Non specific prochiral atom groups, e.g. Val CH3a, CH3b:

  • one or two Resonance(s), one ResonanceSet, two AtomSets, three Atoms each.

  • If only one resonance is observed (as in the diagram) the other Resonance object is not created.

AtomSets are contained directly inside the NmrProject. Except for nuclei in fast exchange, each Atom has its own private AtomSet.  CcpNmr Analysis creates all necessary AtomSets from the reference information in the chemical compound templates (ChemComps). However, they are simple enough to create once you know which atoms to combine:

# This will fail with an existing project

# where the atom sets are already made

atoms = residue.findAllAtoms(name='CA')

atomSet = nmrProject.newAtomSet(atoms=atoms)

To look up existing AtomSets it is a simple matter to find them from the appropriate Atoms you wish to assign to, noting that when we are fetching a methyl AtomSet it doesn't matter which of the three atoms we use:

# Get residues

ser = chain.findFirstResidue(ccpCode='Ser')

ala = chain.findFirstResidue(ccpCode='Ala')

val = chain.findFirstResidue(ccpCode='Val')

# Get atomSets

set_SerHb2 = ser.findFirstAtom(name='HB2').atomSet

set_SerHb3 = ser.findFirstAtom(name='HB3').atomSet

set_AlaHb = ala.findFirstAtom(name='HB1').atomSet

set_ValHg1 = val.findFirstAtom(name='HG11').atomSet

set_ValHg2 = val.findFirstAtom(name='HG21').atomSet

Once you have the right AtomSets, atom assignment is a matter of creating one or more ResonanceSets. The examples from the figure would work like this:

from ccpnmr.analysis.core import AssignmentBasic

getName = AssignmentBasic.makeResonanceGuiName

# Get some new demo resonances to work with

resonance1 = nmrProject.newResonance(isotopeCode='1H')

resonance2 = nmrProject.newResonance(isotopeCode='1H')

resonance3 = nmrProject.newResonance(isotopeCode='1H')

resonance4 = nmrProject.newResonance(isotopeCode='1H')

resonance5 = nmrProject.newResonance(isotopeCode='1H')

# Do the atom assignment

nmrProject.newResonanceSet(resonances=[resonance1],

                          atomSets=[set_SerHb2])

nmrProject.newResonanceSet(resonances=[resonance2],

                          atomSets=[set_AlaHb])

nmrProject.newResonanceSet(resonances=[resonance3, resonance4],

                          atomSets=[set_SerHb2, set_SerHb3])

nmrProject.newResonanceSet(resonances=[resonance5],

                          atomSets=[set_ValHg1, set_ValHg2])

print getName(resonance1)

print getName(resonance2)

print getName(resonance3)

print getName(resonance4)

print getName(resonance5)

Again, you are safer using a utility function that will take care of the bookkeeping relevant to modifying pre-existing assignments. Once you have the list of atomSets, cases 1, 2, and 4, from the figure are all handled by:

AssignmentBasic.assignAtomsToRes([set_AlaHb], resonance2)

Only case 3 is different, because we wish to link the two resonances by the same ResonanceSet (assignment object):

atomSets =[set_SerHb2, set_SerHb3]

newResonanceSet = AssignmentBasic.assignAtomsToRes(atomSets, resonance3)

AssignmentBasic.assignAtomsToRes(atomSets, resonance4, newResonanceSet)

Getting the Atom(s) that a Resonance is assigned to is in itself simple. The bigger question is what you want to do in the cases where the Resonance corresponds to several AtomSets, maybe with several Atoms each. The code you need simply finds any ResonanceSet and then loops through the linked AtomSets to get at the atoms:

resonanceSet = resonance5.resonanceSet

atoms = []

if resonanceSet:

 for atomSet in resonanceSet.sortedAtomSets():

   atoms.extend(atomSet.atoms)

print [atom.name for atom in atoms]

Going the other way - from an Atom to the  Resonances it is assigned to - is a very similar exercise, except that the resonanceSets have no physical meaning and so are not worth recording:

atom = val.findFirstAtom(name='HG21')

atomSet = atom.atomSet

resonances = []

if atomSet is not None:

 for resonanceSet in atomSet.resonanceSets:

   resonances.extend(resonanceSet.resonances)

# Print out some PPM values

for resonance in resonances:

 for shift in resonance.shifts:

   print shift.value

Spin Systems: Resonance Groups

ResonanceGroups are groups of resonances that should correspond to atoms from a single residue. Typically this will contain backbone resonances, and then side chain resonances will be added as an assignment grows. Sometimes side chain resonances may be placed into separate spin systems, until such time as they can be merged into an existing backbone spin system. ResonanceGroups are contained directly inside the NmrProject and can be created like this, as they have no mandatory attributes except for the automatic serial number

resonanceGroup = nmrProject.newResonanceGroup()

You can set the resonance by e.g.

resonanceGroup.addResonance(resonanceH)

# or

listOfResonances = [resonanceH, resonanceN]

resonanceGroup.resonances = listOfResonance

and get it back in the normal ways:

setOfResonances = resonanceGroup.resonances

# or

listOfResonances = resonanceGroup.sortedResonances()

print setOfResonances

print listOfResonances

You can assign ResonanceGroups to molecule and residue types:

ccpCode = residue.ccpCode

print ccpCode

resonanceGroup.molType = 'protein'

resonanceGroup.ccpCode = ccpCode

Or you can assign a specific residue directly:

resonanceGroup.residue = residue

It is also possible to set residue probabilities, residue type probabilities, sequential neighbour probabilities and the probability of a resonance belonging within a ResonanceGroup. One example should suffice. Say that you have four ResonanceGroups, resGroupA, resGroupsB, resGroupsC and resGroupD, and your experiments show you that  the sequential order must be either

resGroupB-resGroupA-resGroupC or

resGroupB-resGroupA-resGroupD.

You would then do:

# Setup some demo spin systems for this example

resGroupA = nmrProject.newResonanceGroup()

resGroupB = nmrProject.newResonanceGroup()

resGroupC = nmrProject.newResonanceGroup()

resGroupD = nmrProject.newResonanceGroup()

resGroupA.newResonanceGroupProb(linkType='sequential', sequenceOffset=-1,

                               isSelected=True, possibility=resGroupB)

resGroupA.newResonanceGroupProb(linkType='sequential', sequenceOffset=1,

                               weight=3, possibility=resGroupC)

resGroupA.newResonanceGroupProb(linkType='sequential', sequenceOffset=1,

                               weight=1, possibility=resGroupD)

Spin system resGroupB is the only possibility for i-1 link, and that possibility is marked as selected (i.e. correct). For the i+1 link there are two possibilities, neither is selected as correct, and the link to resGroupC is three times more likely than the link to resGroupD.