project-gemmi / gemmi

macromolecular crystallography library and utilities
https://project-gemmi.github.io/
Mozilla Public License 2.0
217 stars 44 forks source link

Normalize in sf2map isn't normalizing? #263

Closed pschmidtke closed 1 year ago

pschmidtke commented 1 year ago

Hi,

I'm trying to generate normalized 2Fo-Fc ccp4 maps with gemmi sf2map from mtz files.

For example by running :

gemmi sf2map --normalize --mapmask=4cwr.pdb --sample=5 4cwr.mtz 4cwr_2Fo_2fofc.ccp4

That's an example using the public domain structure from the RCSB: https://files.rcsb.org/download/4CWR.pdb with the mtz provided here: https://edmaps.rcsb.org/coefficients/4cwr.mtz

The resulting file doesn't seem to be normalized (mean is not 0, rms is not 1 in the ccp4 header) + when I display the results with NGL they differ from the RCSB normalied map (https://edmaps.rcsb.org/maps/4cwr_2fofc.dsn6)

You can compare also with pymol or coot but careful they normalize upon opening unless you specify otherwise.

Could you please advise if I'm doing something wrong or if that's a bug in the current sf2map implementation?

Cheers!

CV-GPhL commented 1 year ago

Remember that normalisation works on a full asymmetric unit content (or a multiple of it). So the usual sequence is

(1) create map

(2) limit to asu (or unit cell)

(3) normalise

(4) extract map region (e.g. around PDB file)

At the end of this you will have a normalised map (step 3), but the map file around your PDB file will not necessarily have a mean of 0 and rms of 1 (if calculating those from the available map values).

So the real question is what the CCP4 map header should record? I always thought that it a bit odd to record the stats as presented by the available map values: after all, one can quickly compute those on-the-fly. It would be much better to record the values of an asymmetric unit ... but that would only be possible for crystallographic maps I guess (and what to do about cryo-EM maps).

Anyway, the important thing is if the normalisation step is done at the correct stage (I'm sure it is inside gemmi) the resulting map will be normalised no matter what region you cut out.

You can compare with a series of commands like that:

fft hklin 4cwr.mtz mapout 1.map <<e LABI F1=FWT PHI=PHWT XYZL ASU e

mapmask mapin 1.map mapout 2.map <<e SCAL SIGM e

mapmask mapin 2.map mapout 3.map xyzin 4cwr.pdb <<e BORD 5 e

which will give you a map ("mapmask mapin 3.map </dev/null") with

 Minimum density .................................    -2.83795
 Maximum density .................................     9.17930
 Mean density ....................................    -0.00167
 Rms deviation from mean density .................     0.97052

Always good to check things against the old "gold-standard" (those Fortran programs have been around/debugged for a VERY long time) ;-)

What I get with

gemmi sf2map --normalize --sample=5 4cwr.mtz 4cwr_2Fo_2fofc.ccp4

is indeed giving a map header of

       Minimum density .................................    -0.80151
       Maximum density .................................     2.33769
       Mean density ....................................     0.00000
       Rms deviation from mean density .................     0.29264

However, running then

mapmask mapin 4cwr_2Fo_2fofc.ccp4 mapout 4.map < /dev/null mapmask mapin 4.map < /dev/null

shows

       Minimum density .................................    -2.73892
       Maximum density .................................     7.98834
       Mean density ....................................    -0.00000
       Rms deviation from mean density .................     1.00000

So it looks to me as if gemmi forgets to (re-)compute the map stats before writing out the file in the "--normalise" case (the actual map seems to be normalised after all).

Cheers

Clemens

On Mon, Apr 24, 2023 at 07:28:20AM -0700, Peter Schmidtke wrote:

Hi,

I'm trying to generate normalized 2Fo-Fc ccp4 maps with gemmi sf2map from mtz files.

For example by running :

gemmi sf2map --normalize --mapmask=4cwr.pdb --sample=5 4cwr.mtz 4cwr_2Fo_2fofc.ccp4

That's an example using the public domain structure from the RCSB: https://files.rcsb.org/download/4CWR.pdb with the mtz provided here: https://edmaps.rcsb.org/coefficients/4cwr.mtz

The resulting file doesn't seem to be normalized (mean is not 0, rms is not 1 in the ccp4 header) + when I display the results with NGL they differ from the RCSB normalied map (https://edmaps.rcsb.org/maps/4cwr_2fofc.dsn6)

You can compare also with pymol or coot but careful they normalize upon opening unless you specify otherwise.

Could you please advise if I'm doing something wrong or if that's a bug in the current sf2map implementation?

Cheers!

-- Reply to this email directly or view it on GitHub: https://github.com/project-gemmi/gemmi/issues/263 You are receiving this because you are subscribed to this thread.

Message ID: @.***>

--

*--------------------------------------------------------------

pschmidtke commented 1 year ago

Thanks @CV-GPhL for looking into this.

I'll play around a bit with the head on my side to check if I understood everything correctly. From a pure user experience of sf2map & with analogy to pymol, coot & other tools, I'd have expected a normalized map & updated header as an output, thus my confusion.

Thanks a lot again!

Peter

wojdyr commented 1 year ago

Remember that normalisation works on a full asymmetric unit content (or a multiple of it). So the usual sequence is

(1) create map

(2) limit to asu (or unit cell)

(3) normalise

(4) extract map region (e.g. around PDB file)

The order in gemmi-sf2map is the same, but header statistics are calculated between (2) and (3). I recollect that I was aware that the statistics in header are not updated, but I thought that the ones from before the normalization are more interesting and useful, because after the normalization it's just mean=0 and rms=1. But that's confusing, so I'll update the statistics after normalization.

BTW, you can use gemmi map file.ccp4 to see statistics from both the header and data.