Uppsala Software Factory Tutorial - Local Density Correlation

Uppsala Software Factory Tutorial - Local Density Correlation

This page describes how to use the local density correlation method of Randy Read to investigate structural differences between molecules (e.g., related by NCS, or in different crystal forms, mutants, complexes, etc.), using Uppsala software. The information obtained can be used either to design an appropriate weighting scheme for NCS (or for harmonic restraints or torsion angle restraints in the case of multiple crystals), or to evaluate if "observed" structural differences are supported by the experimental data.

Reference: G.J. Kleywegt (1999). Experimental assessment of differences between related protein crystal structures. Acta Cryst. D55, 1878-1884.

1 - Calculating the local density correlation map

The local density correlation map can be calculated with COMA (in the case of NCS), or with MASKIT (in the case of multiple crystals or crystal forms). In both cases, you will need:

electron density map(s) that are as little biased as possible (best: experimental map, e.g. MIR or MAD phased; if none are available: map calculated with phases obtained from NCS-averaging; note that the actual maps must not be averaged !)
Cartesian operator(s) that relate one molecule to the other(s); best if these are obtained without any reference to an atomic model (i.e., from a superpositioning of heavy-atom or selenium sites, or from molecular replacement calculations using density from one of the molecules to find the other(s))

The procedure to calculate the local density correlation map is described in the respective program manuals.

2 - Extracting information from the map

Now you are ready to extract the average local value of the local density correlation map for every atom in a model. This can be done easily with MAPMAN, using the "PEek SPhere" command:

      
 ----- EXAMPLE ----- EXAMPLE ----- EXAMPLE ----- EXAMPLE ----- EXAMPLE -----
 MAPMAN > read m1 coma4a_r2.E ccp4
 Read header
 Input map : (coma4a_r2.E) 
 ...
 Map read into memory
 Closing BINARY CCP4 map on unit : (       1) 
 Requested dynamic range : -3.5239E-01  8.5839E-01
 Value of Prod and Plus  :  2.1061E+02      74
 Actual dynamic range    : -3.5136E-01  8.5942E-01
 MAPMAN > peek sphere
 Which map ? (coma4a_r2.E) m1
 Input PDB file ? (in.pdb) p2_mol_a.pdb
 Output PDB file ? (out.pdb) corr.pdb
 Mode (Int.|Abs.int.|Rms|Mean) ? (Rms) mean
 Radius (A) ? (3.0) 1.5
 Grid spacing (A) : (   1.330    1.327    1.345) 
 Peeking ...
 New PDB file created
 ----- EXAMPLE ----- EXAMPLE ----- EXAMPLE ----- EXAMPLE ----- EXAMPLE -----

While you're at it, you may also want to MAppage the CCP4 map for plotting in O (if you use the FM_ commands in O to draw the map, it will be rescaled in terms of a meaningless "sigma" level; I prefer the good old MAP_ commands which allow you to draw the map at a level of 0.5, 0.8, etc., so you know what kind of correlation cut-off you are looking at).

3 - Inspection with O

For visual inspection only, you can now go into O, and draw the correlation map at an appropriate cut off level (typically, somewhere between 0.4 and 0.6).

      
 ----- EXAMPLE ----- EXAMPLE ----- EXAMPLE ----- EXAMPLE ----- EXAMPLE -----
map_file coma4a_r2.map
map_obj coma
map_par 100 100 100 0.5 blue 1 0.01 1
map_act    
map_draw
 ----- EXAMPLE ----- EXAMPLE ----- EXAMPLE ----- EXAMPLE ----- EXAMPLE -----

Picture

In addition, you can read in the PDB file created with MAPMAN that now has each atom's B-factor replaced by the average local density correlation value. You can then use the PAINT_RAMP or PAINT_PROPERTY commands to show which parts of the model are expected to disobey the Cartesian (NCS) operator (these will often "stick out" of the correlation map), and which parts are expected to be esentially identical in the various molecules:

      
 ----- EXAMPLE ----- EXAMPLE ----- EXAMPLE ----- EXAMPLE ----- EXAMPLE -----
s_a_i corr.pdb corr
 Sam> File type is PDB
 Sam>  Database compressed.
 Sam> Space for    273712 atoms
 Sam> Space for     10000 residues
 Sam> Molecule CORR contained 131 residues and 1039 atoms
mol corr
pai_ramp
 Paint> Colour-ramp a property in molecule CORR
 Paint> Property [residue_irc] : atom_b
 Paint> Minimum and maximum value of property [-0.0440 0.8447] : 
 Paint> First colour [red] : 
 Paint>                  red
 Paint> Second colour [blue] : 
 Paint>                 blue
obj ramp zo ; end
 Mol> No connectivity Db for CORR
 Mol>  Nothing marked for deletion, so no compression.
  131O      not found
 mol connectivity is        1521
 ----- EXAMPLE ----- EXAMPLE ----- EXAMPLE ----- EXAMPLE ----- EXAMPLE -----

Picture

If you like, you can use PAINT_PROP instead of PAINT_RAMP, as follows:

      
 ----- EXAMPLE ----- EXAMPLE ----- EXAMPLE ----- EXAMPLE ----- EXAMPLE -----
pai_prop atom_b > -999 red
 Paint>                  red
pai_prop atom_b > 0 orange
 Paint>               orange
pai_prop atom_b > 0.25 yellow
 Paint>               yellow
pai_prop atom_b > 0.5 green
 Paint>                green
pai_prop atom_b > 0.75 blue
 Paint>                 blue
obj prop zo ; end
 ----- EXAMPLE ----- EXAMPLE ----- EXAMPLE ----- EXAMPLE ----- EXAMPLE -----

Picture

4 - Plotting

If you want to generate line plots of, for instance, CA-CA distance and local density correlation values for the main chain, then you can use LSQMAN (to generate plot files of CA-CA distance as a function residue number; see the manual), MOLEMAN2 (to generate plot files of average main-chain density correlation as a function of residue number), and O2D (to convert the plot files to PostScript files). If you are handy with an editor, you can combine the two separate plot files into one.

In LSQMAN, you could do as shown below. Make sure that you use exactly the same RT operator as the one that you used in the calculation of the local density correlation map !!!

      
 ----- EXAMPLE ----- EXAMPLE ----- EXAMPLE ----- EXAMPLE ----- EXAMPLE -----
 LSQMAN > re m1 1ytt.pdb
 Cell : (  65.508   72.216   45.035   90.000   90.000   90.000) 
 ...
 LSQMAN > old m1 m1 rt_a_to_b.o
 Operator bringing : (M1) 
 on top of         : (M1) 
 Opening O datablock : (rt_a_to_b.o) 
 ...
 Rotation angle                             176.185
 Operator looks okay !
 LSQMAN > atom ca
 Nr of atom types : (       1) 
 Type : ( CA) 
 LSQMAN > dist m1 b1-999 m1 a1 cacadist.plt
 WARNING - mol1 == mol2 !
 Central-atom distance plot
 Central atom type : ( CA) 
 Plot of M1 B1-999
 And     M1 A1
 Nr of residues matched : (        112) 
 Average distance : (   0.658) 
 Minimum distance : (   0.028) 
 Maximum distance : (   7.729) 
 Plot file written
 ----- EXAMPLE ----- EXAMPLE ----- EXAMPLE ----- EXAMPLE ----- EXAMPLE -----

In MOLEMAN2, you could do as shown below to get a plot of average local correlation per residue for all atoms, and one for the main-chain atoms.

      
 ----- EXAMPLE ----- EXAMPLE ----- EXAMPLE ----- EXAMPLE ----- EXAMPLE -----
 MOLEMAN2 > read mola_corr.pdb
 ...
 MOLEMAN2 > sel all
 Select ALL atoms
 Selection history : (ALL |) 
 Nr of selected atoms : (        892) 
 MOLEMAN2 > bf plot all.plt aver
 Residues with selected atoms : (        115) 
 MOLEMAN2 > sel and class main
 AND atom selection
 With atoms for which : (CLASS) 
 Equals : (MAIN) 
 Selection history : (ALL | AND CLass = Main |) 
 Nr of selected atoms : (        460) 
 MOLEMAN2 > bf plot main.plt aver
 Residues with selected atoms : (        115) 
 ----- EXAMPLE ----- EXAMPLE ----- EXAMPLE ----- EXAMPLE ----- EXAMPLE -----

Latest update at 18 November, 1999.