• Docs »
  • Visualizing B-factor new

Visualizing B-factor new

Problem

You want to visualize the B-factor of a ligand and its environment in order to reveal regions with high flexibility at the active site. See example in Figure 1.

hover the mouse over any ligand atom to display its B-factor

Figure 1. Example of depicting the B-factor of 2A1B complex

../_images/bfactor2img-pdb1a1b1.svg

Ligand atoms with high B-factor (i.e. high flexibility) are colored red, while blue color indicates low B-factors.

The 2D molecule surface (arcs around the ligand) represents the environment of the ligand in 3D. Each arc is associated with a ligand atom that is the center of the arc. A red arc indicates that the average B-factor of the nearby protein atoms is high i.e. the environment is more flexible. When this average B-factor is calculated, only protein atoms that are closer than 4.0 Ångströms (user-adjustable parameter) to the given ligand is considered. Blue arcs represent low average B-factor values of the nearby protein atoms. The lack of arc indicates that there is no protein heavy atom closer than 4.0 Ångströms to the corresponding ligand atom.

The linear color gradient, used for coloring B-factors, is drawn at the bottom of the image. The number (that is 18.00) marked on the color gradient is the average B-factor of the whole complex. The black box shown on the color gradient indicates the range of B-factor values close at the active size i.e. the minimum and maximum B-factor values of any protein atom that is closer than 4.0 Ångströms to any ligand atoms.

Ingredients

  • OEChem TK - cheminformatics toolkit (including OEBio TK)
  • OEDepict TK - molecule depiction toolkit
  • Grapheme TK - molecule and property visualization toolkit

Note

Requires OpenEye toolkits version 2017.Jun or later.

Difficulty level

../_images/chilly6.png ../_images/chilly6.png ../_images/chilly6.png

Solution

The ConsiderBFactor function is used to calculate the average B-factor for the whole complex. When a protein-ligand complex is read from a PDB file, an OEResidue object is attached to each atom as generic data. This OEResidue object stores a number of fields specific to processing macromolecules including the the B-factor.

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
def GetAverageBFactor(mol):
    nratoms, sumbfactor = 0, 0.0

    for atom in mol.GetAtoms():
        if not OEHasResidue(atom):
            continue
        res = OEAtomGetResidue(atom)
        sumbfactor += res.GetBFactor()
        nratoms += 1

    avgbfactor = sumbfactor / nratoms
    return avgbfactor

The GetMinAndMaxBFactor function is utilized to calculate the minimum and maximum B-factors at the region of the active site. Therefore, only protein atoms that are close to the ligand are considered (as defined by the maxdistance parameter). These nearby protein atoms are identified by using the OENearestNbrs class. The user defined atom predicate, NotHydrogenOrWater, is also used to ignore any hydrogen atoms or water molecules.

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
def GetMinAndMaxBFactor(ligand, protein, maxdistance):

    minbfactor, maxbfactor = float("inf"), float("-inf")

    # ligand atoms

    for latom in ligand.GetAtoms(OEIsHeavy()):
        if not OEHasResidue(latom):
            continue
        res = OEAtomGetResidue(latom)
        minbfactor = min(minbfactor, res.GetBFactor())
        maxbfactor = max(maxbfactor, res.GetBFactor())

    # protein atoms close to ligand atoms
    considerbfactor = NotHydrogenOrWater()

    nn = OENearestNbrs(protein, maxdistance)
    for latom in ligand.GetAtoms(OEIsHeavy()):
        for neigh in nn.GetNbrs(latom):
            ratom = neigh.GetBgn()

            if considerbfactor(ratom):
                res = OEAtomGetResidue(ratom)
                minbfactor = min(minbfactor, res.GetBFactor())
                maxbfactor = max(maxbfactor, res.GetBFactor())

    return minbfactor, maxbfactor

The NotHydrogenOrWater atom predicate returns true if the given atom is not a hydrogen and not part of a water molecule.

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
class NotHydrogenOrWater(OEUnaryAtomPred):
    def __call__(self, atom):

        if atom.GetAtomicNum() == OEElemNo_H:
            return False
        if not OEHasResidue(atom):
            return False

        res = OEAtomGetResidue(atom)
        if OEGetResidueIndex(res) == OEResidueIndex_HOH:
            return False

        return True

OEGrapheme TK’s 2D molecule surface depiction style is used to visualize the B-factor of the protein around the ligand. The SetAverageBFactorOfNearbyProteinAtoms function below loops over the ligand atoms and calculate the average B-factor of nearby protein atoms. This average B-factor will later be used to determine the color of the arc drawn by the BFactorArcFxn class.

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
def SetAverageBFactorOfNearbyProteinAtoms(ligand, protein, itag, maxdistance):

    considerbfactor = NotHydrogenOrWater()

    nn = OENearestNbrs(protein, maxdistance)
    for latom in ligand.GetAtoms(OEIsHeavy()):
        sumbfactor = 0.0
        neighs = []
        for neigh in nn.GetNbrs(latom):
            ratom = neigh.GetBgn()
            if considerbfactor(ratom):
                res = OEAtomGetResidue(ratom)
                sumbfactor += res.GetBFactor()
                neighs.append(ratom)

        avgbfactor = 0.0
        if len(neighs) > 0:
            avgbfactor = sumbfactor / len(neighs)
        latom.SetDoubleData(itag, avgbfactor)

The linear color gradient used in the depiction is initialized by the GetBFactorColorGradient function. The color gradient does not depend on the depicted complex thereby allowing to easily compare sets of B-factor images.

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
def GetBFactorColorGradient():

    colorg = OELinearColorGradient()
    colorg.AddStop(OEColorStop(0.0, OEDarkBlue))
    colorg.AddStop(OEColorStop(10.0, OELightBlue))
    colorg.AddStop(OEColorStop(25.0, OEYellowTint))
    colorg.AddStop(OEColorStop(50.0, OERed))
    colorg.AddStop(OEColorStop(100.0, OEDarkRose))

    return colorg

The ColorLigandAtomByBFactor class is used to draw the colored circles below the ligand atoms indicating their B-factor. The color gradient used in this class are generated by the GetBFactorColorGradient function.

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
class ColorLigandAtomByBFactor(OEAtomGlyphBase):
    def __init__(self, colorg):
        OEAtomGlyphBase.__init__(self)
        self.colorg = colorg

    def RenderGlyph(self, disp, atom):
        adisp = disp.GetAtomDisplay(atom)
        if adisp is None or not adisp.IsVisible():
            return False

        if not OEHasResidue(atom):
            return False

        res = OEAtomGetResidue(atom)
        bfactor = res.GetBFactor()
        color = self.colorg.GetColorAt(bfactor)

        pen = OEPen(color, color, OEFill_On, 1.0)
        radius = disp.GetScale() / 3.0

        layer = disp.GetLayer(OELayerPosition_Below)
        OEDrawCircle(layer, OECircleStyle_Default, adisp.GetCoords(), radius, pen)
        return True

    def CreateCopy(self):
        return ColorLigandAtomByBFactor(self.colorg).__disown__()

The BFactorArcFxn class draws the arc segments around the depicted ligand. The color of each arc represents the average B-factor of the protein atoms around the corresponding ligand atom. These values are calculated by the SetAverageBFactorOfNearbyProteinAtoms function and attached to the ligand atoms as generic data. The color gradient used in this class are generated by the GetBFactorColorGradient function.

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
class BFactorArcFxn(OESurfaceArcFxnBase):
    def __init__(self, colorg, itag):
        OESurfaceArcFxnBase.__init__(self)
        self.colorg = colorg
        self.itag = itag

    def __call__(self, image, arc):
        adisp = arc.GetAtomDisplay()
        if adisp is None or not adisp.IsVisible():
            return False

        atom = adisp.GetAtom()
        if atom is None:
            return False

        avgresiduebfactor = atom.GetDoubleData(self.itag)
        if avgresiduebfactor == 0.0:
            return True
        color = self.colorg.GetColorAt(avgresiduebfactor)

        pen = OEPen(color, color, OEFill_Off, 5.0)

        center = arc.GetCenter()
        bAngle, eAngle = arc.GetBgnAngle(), arc.GetEndAngle()
        radius = arc.GetRadius()

        OEDrawDefaultSurfaceArc(image, center, bAngle, eAngle, radius, pen)

        return True

    def CreateCopy(self):
        return BFactorArcFxn(self.colorg, self.itag).__disown__()

DepictBFactor is the main function that does the depiction of the ligand and draws molecule surface around it to represent its environment.

  1. First, the depiction coordinates of the ligand are generated by calling the OEPrepareDepictionFrom3D function (lines 5). Then the OEPrepareDepictionOptions class and OEPrepareDepiction function of the OEDepict TK are utilized to orient the generated 2D coordinates horizontally (lines 7-10).
  2. The BFactorArcFxn arc drawing functor is added to each ligand (see also step 7)
  3. Scaling of the depicted ligand is reduced by using the OEGetMoleculeSurfaceScale function to make sure that arcs representing the molecule surface of the ligand will not be drawn outside the image (lines 12-15).
  4. The molecule display for the ligand is constructed (lines 19).
  5. If the output image file format is svg, then the OEDrawSVGHoverText function can be used to associate each display atom with a text of its B-factor value. The generated svg image will be interactive i.e. the B-factor texts will be displayed when the mouse is hovered over the ligand atoms in the image (lines 21-31).
  6. The OEAddGlyph function is used to draw a glyph (in this case a circle colored by the B-factor) underneath each ligand atom as defined in the ColorLigandAtomByBFactor class (lines 32-33).
  7. The OEDraw2DSurface function is called that will invoke the arc drawing functor set on each ligand atom in step 2. resulting in the rendering of the the molecule into the display (lines 35).
  8. Finally, the OERenderMolecule function renders the molecule display into the image (lines 37).
 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
def DepictBFactor(image, ligand, opts, colorg, itag, issvg):

    # prepare ligand for depiction

    OEPrepareDepictionFrom3D(ligand)

    clearcoords, suppressH = False, False
    popts = OEPrepareDepictionOptions(clearcoords, suppressH)
    popts.SetDepictOrientation(OEDepictOrientation_Horizontal)
    OEPrepareDepiction(ligand, popts)

    arcfxn = BFactorArcFxn(colorg, itag)
    for atom in ligand.GetAtoms():
        OESetSurfaceArcFxn(ligand, atom, arcfxn)
    opts.SetScale(OEGetMoleculeSurfaceScale(ligand, opts))

    # render ligand and visualize BFactor

    disp = OE2DMolDisplay(ligand, opts)

    if issvg:
        font = OEFont(OEFontFamily_Default, OEFontStyle_Default,
                      14, OEAlignment_Center, OEBlack)
        for adisp in disp.GetAtomDisplays():
            atom = adisp.GetAtom()
            if not OEHasResidue(atom):
                continue
            res = OEAtomGetResidue(atom)
            hovertext = "bfactor=%.1f" % res.GetBFactor()
            OEDrawSVGHoverText(disp, adisp, hovertext, font)

    colorbfactor = ColorLigandAtomByBFactor(colorg)
    OEAddGlyph(disp, colorbfactor, OEIsTrueAtom())

    OEDraw2DSurface(disp)

    OERenderMolecule(image, disp)

Download code

bfactor2img.py

Usage:

prompt > python3 bfactor2img.py complex.pdb bfactor.png

Discussion

The B-factor (temperature factor) is a measure of how much an atom vibrates around the position specified in the PDB file. Atoms at side-chain termini are expected to exhibit more freedom of movement than main-chain atoms. Similarly, chain ligand atoms exposed to solvent can have more freedom of movement (see 1A4K example below).

B-factors can indicate not only the mobility of the atoms but they can also reveal where there are errors in model building. Visualization helps to make instant judgment on the quality of model.

PDB: 1A4K PDB: 1ATL

../_images/bfactor2img-pdb1a4k.svg ../_images/bfactor2img-pdb1atl.svg

See also in OEChem TK manual

Theory

API

See also in OEDepict TK manual

Theory

API

See also in Grapheme TK manual

Theory

API