7. Averaging of carbamoylsarcosine hydrolase The averaging of the electron density of the carbamoylsarcosin hydrolase is divided in three parts: averaging with 1 mask, 2 half masks and finally with 4 masks - one for each molecule. Here, only those aspects done with MAIN are described. The first part describes the averaging inside 1 mask. In one mask only proper averaging can be applied. By splitting density into 2 halves, averaging could be performed applying proper electron density averaging inside each half of the asymmetric unit followed by the improper averaging between the two halves. When the model was complete enough, the density could be divided into 4 parts, one for each molecule. These 4 parts could be then averaged by applying only improper symmetry averaging. 7.1. One mask There are only two uranium positions per asymmetric unit and they are not related by a local symmetry operation. The density in the areas around the uranium atom positions in the solvent flattened electron density map derived from the heavy atom phases (U and Hg) has large density fluctuations in range between -60.0 and 20.0 sigma units. These areas disturbed the self-correlation calculations of the masked region. Therefore the density 4.5A around the uranium positions should be erased (set to 0.0). This changed interval of the electron density map to -10.0 to 10.0 sigma units. This procedure is stored in the ERASE.COM file. ERASE.COM: $ main Read the solvent flattened density map CSURHG_SOLV.RMAP (the whole cell) with 130, 120 and 70 grid points per unit cell cell lengths. > read file csurhg_solv.rmap map number Transform all the grid points 140, 130, and 80 grid points around the current center (0.0, 0.0, 0.0) in map 1, whose density is in the range between -100.0 to -20.0, to the density points. The command word CUT disables the extrapolations of the density outside the region defined by the map boundaries. > make point initialize map 1 140 130 80 range -100. -20. cut Convert the generated points to the atoms and mask the area around them. The masked grid points obtain values 9999.0 (since it is a real*4 map) and > make atom points > make map 1 atom 4.5 select all end are set to 0.0. > make map 1 set 9000. 10000. 0.0 The resulting electron density map is then written to the file CSURHG_SOLV_ERASE.RMAP. > write file csurhg_solv_erase.rmap map 1 number > quit 7.1.2. Conversion of X-CONTOUR to MAIN mask The conversion of an X-CONTOUR mask to a MAIN mask is done via a protein density file. The empty region in X-CONTOUR has values 0.0 and the masked region 1.0. When a PROTEIN density file is read by MAIN its density values are scaled so that 1.0 in MAIN corresponds to 1.0 sigma of the map. In order to convert the X-CONTOUR mask to MAIN mask, the empty region has to be set to -9999.0 and the masked region to 9999.0 in the case of real*4 maps, while for a character*1 map conversion with PROMAP is sufficient. MAKE_MAIN_MASK.COM: $ main Read the PROTEIN mask file. > read file X-contour_mask.diy map protein Modify the empty and masked regions to meet MAIN requirements. > make map 1 set -10 0 -9999.0 > make map 1 set 0 100 9999.0 Write the MAIN mask file. > write file cs_mask.rmap map 1 > quit 7.1.2. Patterson of the mask In order to do the test for consistency of the mask, the density inside the masked region has to be transferred to a P1 cell and its Patterson density should be calculated in order to do a self-rotation search that shows if the asymmetric unit included in the mask is complete. The masked region from the file CS_MASK.RMAP has to be filled with density from CSURHG_SF_ERASE.RMAP and transformed to a larger P1 cell where no contacts with the crystal symmetry related masked regions are possible. $ main > read file cs_mask.rmap map number > read file csurhg_sf.rmap map number > make map 1 from 2 copy > write file mask_dens.rmap map 1 > quit The simplest way how to modify the size of a map unit cell is to use a text editor. The map files in MAIN format are ASCII with record length of 80 characters. So the header of the file MASK_DENS.RMAP that includes data about the cell constants, number of grid points per cell lengths, origin of the map and its last point and symmetry operations should be modified by a text editor. The cell length and number of grid divisions per cell lengths were multiplied approximately by factor 1.5, so that the grid sizes remain unchanged. The density file header 0.13622E+03 0.12229E+03 0.70870E+02 0.90000E+02 0.91820E+02 0.90000E+02 130 120 70 -10 40 0 86 81 86 4 1 1.000 0.000 0.000 0.000 1.000 0.000 0.000 0.000 1.000 0.000 0.000 0.000 2 -1.000 0.000 0.000 0.000 1.000 0.000 0.000 0.000 -1.000 0.000 0.000 0.000 3 1.000 0.000 0.000 0.000 1.000 0.000 0.000 0.000 1.000 0.500 0.500 0.000 4 -1.000 0.000 0.000 0.000 1.000 0.000 0.000 0.000 -1.000 0.500 0.500 0.000 should be changed to 0.20433E+03 0.18343E+03 0.10529E+03 0.90000E+02 0.91834E+02 0.90000E+02 195 180 104 0 0 0 86 81 86 1 1 1.000 0.000 0.000 0.000 1.000 0.000 0.000 0.000 1.000 0.000 0.000 0.000 The multiplication of 70 by 1.5 yields 105 and not 104 as above. 105 grid divisions could not be applied since the number of grid points in a y density layer should be even, as required by the P1SF. So either 195 or 105 had to be changed to an even number. The cell lengths were respectively increased with factors (1.5., 1.5, 104./70.); the cell angles remained unchanged. The beginning of the new (P1) density was moved to origin (0,0,0) and the number of symmetry operations was reduced to one (X,Y,Z). 7.1.3. New X-CONTOUR mask After the approximate center of the tetramer was found, it was evident that the first mask is not placed correctly since the center of rotations was placed in the 'lower' third of the mask. Therefore the solvent flattened density has to be rotated around all three 2-folds, in the hope that the areas that are related by the local symmetry operations will obtain higher density than the background. The following procedure differs from the usual procedure (ROTA_MAPS and AVER_MAPS) only in the definition of the mask. In this case the whole area of the map (P1_100.RMAP) is masked (maps 2, 3 and 4). The maps 2, 3 and 4 are created by rotation of density in map 1. After that their densities are averaged and the stored in the file P1_AVER.RMAP. MAIN_MAP has to be used, since the usual MAIN doesn't have enough memory for such large real*4 maps. AVER_P1.COM: $ main_map Read the large area map from which is suppose to be larger than the asymmetric unit. > read file p1_100.rmap map numbers Make 3 maps of the same size as the map 1 and assign all grid points as masked area. > make map 2 from 1 initialize 9999 > make map 3 from 1 initialize 9999 > make map 4 from 1 initialize 9999 Rotate density of map 1 to map 2, 3 and 4 by all 3 twofolds. > set matrix 2 polar 91. 82. 180. > make map 2 from 1 rotate matrix 2 cent coordinate 30.41 93.3 49.85 - > translate vector 8.60891 -3.61775 1.45716 > set matrix 2 polar 45. 171. 180. > make map 3 from 1 rotate matrix 2 cent coordinate 30.41 93.3 49.85 - > translate vector 2.30575 2.49675 1.96312 > set matrix 2 polar 45. 352. 180. > make map 4 from 1 rotate matrix 2 cent coordinate 30.41 93.3 49.85 - > translate vector 5.96678 -5.85973 2.16550 Set the points that are transformed to empty points by the previous step to 0.0. Some points are transformed to empty because they fell out of map 1 when the local symmetry operations have been applied. > make map 2 set -10000 -9000 0.0 > make map 3 set -10000 -9000 0.0 > make map 4 set -10000 -9000 0.0 Add all maps to map 1 > make map 1 from 2 add > make map 1 from 3 add > make map 1 from 4 add Divide map 1 by 4 and > scale 0.25 map 1 save the results in the MAIN (NUMBER) and Lyn Ten Eyck (FFT) format files. > write file p1_100_aver.rmap map 1 number > write file p1_100_aver.lte map 1 fft > quit The proper 222 averaging: The following procedure in the DO_ALL.COM differs from the ones described by cathepsin B: - the initial density in this case is the solvent flattened density with modified region around the uranium positions, - the mask was created by X-CONTOUR and not derived from an atomic model. The solvent flattened density is first copied to the file CS.RMAP. The first averaged density is stored in file AVER_1.RMAP in order to keep the density before it is iteratively modified by fourier transforms of averaged electron density and structure factor calculations. Therefore, the first @ROTA_MAPS and @AVER_MAPS command procedures are excluded from the LOOP. $ setup protein $ cd pic$exa:[aver.cs.222] $ tmp :== aaa $ cop [--.dat]uruohg_sf_i.rmap cs.rmap $ ii = 0 $ @rota_maps ! transforms density by local symmetry $ @aver_maps ! averages the density in place of molecules $ copy aver.rmap aver_1.rmap $ goto here $ loop: $ @rota_maps ! transforms density by local symmetry $ @aver_maps ! averages the density in place of the mask $ here: $ @cell ! fills the cell from the averaged densities applying $ ! crystal symmetry operations $ @fftfc $ @fourier $ ii = ii + 1 $ if ii .lt. 8 then goto loop The rotational and translational parameters for averaging were found by a search procedure applying program PROTEIN. They consist of a matrix defined by polar angles, rotation around a center and additional translation vector. At the beginning the local symmetry operations were applied as they came out of rotational and translational search done with PROTEIN. The density didn't seem to be good, so we start to realize that the search didn't result in the ideal tetrahedral local symmetry. Since perfect tetrahedral symmetry for the tetramer was expected, we thought that the density wasn't good enough to find it. So the rotational matrices were adjusted to it by slightly modifying the psi and phi angles and translational search was repeated again. (The rotational axes used had to be perpendicular to each other.) Since the density still wasn't good enough we tried to check also the deviation from ideality of the translational component. The check of deviations from the ideal 222 symmetry was done as described below: A network of points should be created. The starting point coordinates are set to (31.0, 91.0, 48.0), since that is approximately the center that was used for the density autocorrelation calculations done by PROTEIN. The 3-dimensional network is built by increasing each coordinate by 1.0 five times in all directions. The created points are converted to atoms, because the connectivities can be defined only on an atom array. Atoms created from points are dummies and between dummy atoms no covalent bond can be calculated, so they should be renamed to real atoms as carbons (CX) for example. MAIN> make points planar 31. 91. 48. 1. 1. 1. 5 5 5 MAIN> make atom point MAIN> rename select all end atom CX The atoms are then written to the file T.COOR in order to be read again four times (for each local symmetry operation once). MAIN> write file t.coor coordinate symbol MAIN> read file t.coor coordinate MAIN> read file t.coor coordinate MAIN> read file t.coor coordinate MAIN> read file t.coor coordinate With each read file statement a new segment is created, so that they all can be differentiated by their segment number. MAIN> calculate bond select segment number 1 end distance 1.1 MAIN> calculate bond select segment number 2 end distance 1.1 MAIN> calculate bond select segment number 3 end distance 1.1 MAIN> calculate bond select segment number 4 end distance 1.1 The atoms in segments obtain colors in the SET_COLOR module, the segment number 1 starting at 10 and then increased by 60 for each next segment. MAIN> set col SET_COL> select all end segment 10 60 SET_COL> exit Local symmetry operations are applied on the segments by running the commands in the MAKE_LOC.COM. MAIN> @make_loc MAKE_LOC.COM: The coordinates of the segment 1 are copied to the segments 2, 3 and 4. > copy coordinate select segment number 1 end select segment number 2 end > copy coordinate select segment number 1 end select segment number 3 end > copy coordinate select segment number 1 end select segment number 4 end Segments 2, 3 and 4 are one after another transformed to their local symmetry positions. > set matrix 2 polar 90. 82. 180. > rotate atom select segment number 2 end matrix 2 - > center coordinates 34.30 91.84 51.00 > translate atom select segment number 2 end - > vector 0.75766 -0.73374 -0.10119 > set matrix 2 polar 45. 172. 180. > rotate atom select segment number 3 end matrix 2 - > center coordinates 34.30 91.84 51.00 > translate atom select segment number 3 end - > vector -0.2768 -0.0200 -0.09995 > set matrix 2 polar 45. 352. 180. > rotate atom select segment number 4 end matrix 2 - > center coordinates 34.30 91.84 51.00 > translate atom select segment number 4 end - > vector 0.43366 -0.42802 0.20238 > return Attach to the picture system and set the image center to the arithmetic middle of the selected atoms in the segments numbered from 1 to 4. The resulting center is stored to the file CENT.DAT by the show command. MAIN> image initialize center calculate select segment number 1 4 end MAIN> show file cent.dat center Create an image of the pattern of covalent bonds. The colors are taken from the set, that was specified above in the SET_COLOR module. MAIN> image set bond The coordinates of the segment number 2 (related by the local symmetry operation to the origin segment number 1) are then copied to the segment number 1. The local symmetry operations are this time applied on already once transformed coordinates. The procedure is then repeated so that each by local symmetry related atomic set obtained also the back transformation. In the ideal case only 4 cubes should be seen. The rotational matrices were already set so they just have to be applied. MAIN> copy coordinate select segment number 2 end select segment number 1 end MAIN> @make_loc MAIN> image set bond MAIN> cop coordinate select segment number 5 end select segment number 1 end MAIN> @make_loc MAIN> cop coordinate select segment number 3 end select segment number 1 end MAIN> @make_loc MAIN> image set bond MAIN> cop coordinate select segment number 5 end select segment number 1 end MAIN> @make_loc MAIN> cop coordinate select segment number 4 end select segment number 1 end MAIN> @make_loc MAIN> image set bond The four displayed tetramers showed that only two of the local symmetry operations obtained are approximately consistent. We decided to introduce ideal tetrahedral symmetry by choosing the center of rotation for all 3 local 2-fold axes by taking the arithmetic middle of the first tetramer as already described above and stored the center coordinates (34.6325, 91.3915, 51.1034) to the file CENT.DAT. MAIN> show file cent.dat center The following ROTA_MAPS.COM applies the ideal tetrahedral symmetry operations. It differs from the ROTA_MAPS.COM procedures described at cathepsin B by the definition of matrix, that is defined directly from the polar angles and not from a command file (SET_TRAN_nn_TO_nn.COM) and the center of rotation is here explicitly defined. ROTA_MAPS.COM $ main > read file cs.diz map protein > read file [-.mask]cs_add_mask.rmap map number > make map 2 from 1 copy > show map > write file mol_01.rmap map 2 > set matrix 2 polar 90. 82. 180. > read file [-.mask]cs_add_mask.rmap map number over 2 > make map 2 from 1 rotate matrix 2 cent coordinate 34.6325 91.3915 51.1034 - > translate vector 0.0 0.0 0.0 > sh map > write file mol_02.rmap map 2 > set matrix 2 polar 45. 172. 180. > read file [-.mask]cs_add_mask.rmap map number over 2 > make map 2 from 1 rotate matrix 2 cent coordinate 34.6325 91.3915 51.1034 - > translate vector 0.0 0.0 0.0 > write file mol_03.rmap map 2 > set matrix 2 polar 45. 352. 180. > read file [-.mask]cs_add_mask.rmap map number over 2 > make map 2 from 1 rotate matrix 2 cent coordinate 34.6325 91.3915 51.1034 - > translate vector 0.0 0.0 0.0 > write file mol_04.rmap map 2 > quit 7.2. Two halves When it became evident that the averaging with ideal tetrahedral symmetry doesn't produce density with sufficient recognizable secondary structure for preparation of four separate masks, we tried to recognize in the solvent flattened density borders between the molecules. However, with the usual 3-dimensional electron density presentations based on contouring at a specific density level, it is possible only to display relatively small parts of the cell and it is not possible to get rid of the noise in the density. (I define noise in electron density as small, scattered clouds of electron density above the contouring level.) The following procedure prepares a line presentation of the solvent flattened density. The procedure is based on search of peaks that are later joined with the lines. The number of peaks is further reduced so that the close lying ones are joined. The parts with only a few connected peaks can be removed as they are assumed to be noise and the larger clouds of continuous electron density are kept. The density has to be divided into boxes in which the peaks are searched. The required operations can be performed only on atomic data array so the search in the density map should be done in small enough boxes, so that number of peaks didn't exceed the boundaries of atomic data array. CELL.COM: $ set def pic$exa:[aver.cs.2half.cell] $ main Read the mask filled with density. > read file mask_dens.rmap map number Create a grid of points started at coordinates (-10., 60., 0. ) with corresponding increments (20., 20., 20.) applied (5 x 5 x 5) times. the points cover the whole masked density region in grid coordinates. > point plane -10. 60. 0. 20. 20. 20. 5 5 5 Since the required transformations can not be applied to points they were converted to atoms. > make atom points The created atoms are marked in a key keep. This selection is defined not to delete them later. > key keep select all end The atomic coordinates were divided by the number of grid per cell length and orthogonalized ( converted from grid coordinates to angstroems in orthogonal space): > scale 1. / 130. 1. / 120. 1. / 70. coordinate select all end > set coordinate select all end orthogonal Initialize the variables iatom and nat_old. 'iatom' is the index of the current atom used as the center for search of peaks nat_old stores the current atom number. 'natom' is the variable defined by the program and keeps the current atomic number: > set variable iatom = 0 > set variable nat_old = natom SHOW the list of presented variables: > show variables OPEN a new file MERGE.XPL and connects it to the logical unit number 10. The open unit statement is used since the data written by the write unit nn are appended to the end of the file with logical unit number nn, while the usual write file statement writes data to a new file and closes it afterwards: > open unit 10 file merge.xpl write_enable Transfer the command input to the file CELL_LOOP.COM: > @cell_loop The CELL_LOOP.COM searches for density peaks around each atom, transforms them to atoms, reduces their number according to distance and connectivity criteria, and writes atoms to the file MERGE.XPL. CELL_LOOP.COM: Increase the variable iatom by 1: > set variable iatom = iatom + 1 If all of atomic centers have been passed, return command input to the previous command input level: > if iatom .gt. natom return SET CENTER to the atom number iatom: > set cent atom iatom SHOW the CENTER coordinates: > show center Search for peaks in density map 1 in a box of 20x20x20 grid points about the current center by the weighting algorithm. Only grid points between 0.8 and 30.0 sigma are considered. No extrapolation of the density outside the defined region is allowed (CUT). > make points initialize map 1 10 10 10 weigh 0.8 30. cut Show the number of found points (peaks): > show points If no peaks are found the search jumps to the next center by rewinding CELL_LOOP.COM. The next command sentence is again the the first one in the file: > if npoin .le. 0 rewind file Create atoms from points: > make atom points Rename the newly created atoms to OX and connects them by covalent bonds if they are closer than 1.9A. The atoms created from points are dummies (X). The program doesn't calculate any covalent bond between the dummy atoms. > rename select .not. keep end atom OX > calculate bond distance 1.9 select .not. keep end Rearrange the atoms so that a network connected via calculated covalent bonds forms one segment: > make segments select .not. keep end Delete all the segments that have fewer than 8 connected atoms: > delete atoms select segment range 1 8 .and. .not. keep end Reduce the cross links: > make atom merge select .not. keep end distance 1.1 > make atom merge select .not. keep end distance 1.3 > make atom merge select .not. keep end distance 1.5 Make sure that no bond is longer than 2.0. If it is, a new atom is introduced in the middle of it and connected to both neighbors. > make atom bond 2.0 select .not. keep end Again rearranges atoms so that a network forms a segment, since the new created atoms with MAKE BOND ... statement are inserted as new segments. > make segments select .not. keep end If more than 10 atoms remain in the atom array, then they are written to the file assigned to unit 10 in the X-PLOR format. > if nat_old + 10 .lt. natom write unit 10 select .not. keep end coordinate xplor Delete all but the keep atoms (The keep atoms are the centers for the search. > delete atom select .not. keep end Rewind the current command input file. > rewind file The peaks stored in the file MERGE.XPL as atoms were later displayed and divided into to two halves (UP-DOWN). The procedure was done interactively at the display. It started by running the command file SHOW_2HALF.COM $ main MAIN> @show_2half SHOW_2HALF.COM: First a 3D-cross is created, oriented in direction of the local axes, translated to the center of the ideal tetrahedral symmetry and assigned to the key local (local axes). The cross should help to recognize the borders between the molecules. > read coordinates > * > C0 0. 0. 0. > CX 30. 0. 0. > CX -30. 0. 0. > CZ 0. 0. 30. > CZ 0. 0. -30. > CY 0. 30. 0. > CY 0. -30. 0. > key local select all end > calculate bond select local end distance 31. > set matrix 2 initialize z-ax -53. > rotate atom select all end matrix 2 center coordinates 0.0 0.0 0.0 > translate atom select all end vector 34.223 91.836 50.80 Here the MERGE.XPL file is read to obtain the coordinates of the peaks converted to atom in the procedure CELL.COM. The density atoms are assigned to KEY merge and connected via covalent bond when the distance between two of them is shorter than 2.0A. > read file merge.xpl coordinate xplor > key merge select .not. local end > calculate bonds distance 2.0 residues select merge end The picture system is attached and center of the image is defined. > image initialize center coordinates 34.223 91.836 50.80 Both atom selections merge and local are then displayed with two colors, one for the cross (110-red) and the other (300-light blue) for the density pattern. > image select local end color 110 bond > image select merge end color 300 bond > return From the image it is evident that two of the local axes define a plane that intersects the displayed density in two well defined separate regions. From three of the 'local' cross atoms the plane equation is derived. MAIN> set plan atom 2 1 6 Parts of the density are assigned to keys 'up' and 'down' according to their distance to the previously defined plane: 'down' for negative and 'up' for positive distances. (The distance can be negative since it calculated as a scalar product between the atom coordinates and the norm of the plane plus constant.) MAIN> key down select plan 1 -200. 0. .and. merge end MAIN> key up select .not. down .and. merge end The image of local axes and of key merge are stored on the display list as the background objects in order to save the time when using the RE_DRAW command. MAIN> image from erase The colors of the 'up' and 'down' atoms are set and the densities are displayed according to the atomic color set. MAIN> set col SET_COL> select down end col 140 ! orange SET_COL> select up end col 320 ! blue SET_COL> exit MAIN> image select merge end set bond The displayed density has shown that there are several connected groups of atoms (residues) that have atoms in both parts (UP and DOWN). Since contacts between two proteins are certainly not in a plane and since we assumed that continuous parts of density belong either to the 'up' or 'down' pair of molecules, we tried to include the whole residue only in one part. To the key 'border' all atoms with color not equal to color 140 which are covalently linked to atoms WITH color 140 ('down') are assigned. With the next statement the whole residues that have a color border are assigned to the KEY 'mix'. MAIN> key border select by bond col 140 .and. .not. col 140 end MAIN> key mix select by residue residue end The images of the 'up' and 'down' are deleted and only the 'mix' selection is displayed. This distinguishes 'mix' residues on the display from the background objects because of the overlapping colors: MAIN> image era MAIN> image select mix end set bond MAIN> image dial The mixed residues are then individually assigned to the UP or DOWN part of the density by subjective criteria. That is done so that an atom of a residue is picked and than the TRANSLATE function is activated on that residue. The parts to be translated are automatically assigned to KEYS 'ob_1', 'ob_2' etc... and one by one removed from 'up' or 'down' parts as shown below: MAIN> key up select up .or. ob_1 end MAIN> key down select down .and. .not. up end MAIN> key border select by bond col 140 .and. .not. col 140 end MAIN> key mix select by residue border end MAIN> image erase MAIN> image select mix end set bond This last group of command sentences is repeated so long that there are no mixed residues left (the KEY 'mix' becomes empty): We need still to store the 'up' and 'down' atoms to the files (UP.XPL and DOWN.XPL) and to derive separate two masks: MAIN> write file up.xpl select up end coordinate xplor MAIN> write file down.xpl select down end coordinate xplor MAIN> quit The two masks can be defined from the atoms by the usual procedure. In this case all grid points 5.0A around the UP and DOWN atoms were assigned to the masks by the MAKE_MASKS.COM. The masks had an overlapping region, but at this stage we assumed that the division is justified and accurate enough. After two separate masks have been derived, the local symmetry relation between them had to be defined. Rotational and translational parameters for averaging of the split masks were derived from correlation of electron density. Inside each mask ('up' and 'down') the parameters for an ideal 2-fold rotation were optimized by the following procedure. $ main MAIN> Read the solvent flattened density and upper mask file and copy the solvent flattened density into the mask: MAIN> read file main$exa:[aver.cs.dat]uruohgko_sf_zero.rmap map number MAIN> read file [-.mask]up_mask.rmap map number MAIN> make map 2 from 1 copy Set the MAIN center. MAIN> set center coordinates 23.26 77.26 52.21 At least 1 atom is required in order to assign density peaks to it. (Otherwise no energy calculation or minimization is possible.) MAIN> read coordinates MAIN> * MAIN> C MAIN> Initialize the point array (sets the number of points to zero), derive density peaks from the density map 1 in a box of 80x80x80 grid points and less than 20.0A away from the MAIN center specified above, which values are within the density range from 3.5 to 50. sigma units. MAIN> make point initialize map 2 40 40 40 radi 20. rang 3.5 50. Turn ALL ENERGY terms OFF and turn ON the DENSITY term. The density energy should be calculated from the points (by default atoms are used) and set the density energy SCALE factor to 1.0. MAIN> energy all off density on density map 2 density - MAIN> scale 1. density points Calculate the auto correlation energy and use it further as a norm (1000.0) and set the density of all grid points in the map that are empty or have negative density to 0.0, so that only grid points with positive density are contributing to the autocorrelation calculations. MAIN> energy select all end MAIN> energy density scale 1000.0 / 20328.15 MAIN> make map 2 set -10000 0 0.0 Here the search procedure can begin by first defining the rotational matrix (2) for energy calculations further applied in the interactive translational search. MAIN> set matrix 2 polar 44. 353. 180. MAIN> energy select all end mat 2 translate 0. 0. 0. MAIN> energy select all end mat 2 translate 1. 0. 0. This has to be repeated until a correlation maximum is found. Afterwards the number of density points has to be increased (from 190 to 3000) by lowering the border of density interval from 3.5 to 2.0 sigma. MAIN> make points initialize map 2 40 40 40 radi 20. rang 2.0 50. The interactive search is followed by an automatic optimization of 2 rotational and 3 translational components using for the starting point the best interactive solution. (Rotation is always appled before translation without regarding the input order.) MAIN> minimize by body select all end translate -1.3 1 2.3 1 -5.0 1 - MAIN> polar 44. 1 353. 1 180. 0 1 or 0 after a real number tells to the program that the preceding real parameter is a variable to be optimized (1) or a constant (0). The resulting rotational and translational parameters are the written to the file UP.DAT. MAIN> show file up.dat rms The correlation energy obtained was 167. Then the same procedure has to be done for the DOWN part of the density, choosing the arithmetic center of the DOWN atoms generated by [-.MASK]CELL.COM as center the of rotation. This time the first commands are read from a file DOWN.COM. MAIN> @down DOWN.COM: > read file down_dens.rmap map number over 2 > set cent coordinate 46.2884 104.7735 51.3439 > make map 2 set -10000 0 0.0 > make points initialize map 2 40 40 40 radi 20. rang 3.5 50. > energy dens scale 1. > set matrix 2 polar 44. 353. 180. > energy select number at 1 end matrix 1 translate 0.0 0.0 0.0 > return The resulting correlation energy was 176. After this the density inside UP and DOWN was averaged according to their 2-fold symmetry ([-.AVER]ROTA_MAPS.COM, [-.AVER]AVER_MAPS.COM). The averaged densities were then used to determine the next two local symmetry operations by a similar procedure as described above. The difference was that now the density peaks from the UP part were correlated to the density in the DOWN part and vice versa. The correlation energy was this time in both cases over 300. 7.2.2. Averaging The procedure for averaging of the split mask (UP and DOWN areas), differs from the usual one by map rotations and averaging. The procedure is a combination of proper and improper averaging, where proper averaging of the UP and DOWN areas (ROTA_MAPS.COM, AVER_MAP.COM) is followed by improper averaging of the properly averaged UP and DOWN maps (ROTA_MAPS_2.COM, AVER_MAP_2.COM). The input density for the cyclic averaging procedure is the solvent flattening density with erased density in the surrounding of the positions of uranium atoms (as described before). DO_ALL.COM: $ setup protein $ cd pic$exa:[aver.cs.2half.aver] $ jj = 0 $ ii = 0 $ cycle: $ @create_p1sf $ @create_merge $ ii = 0 $ ! MAPOUT_R converts the real*4 MAIN density file into a protein density file. $ mapout_r [--.dat]uruohgko_sf_erase.rmap cs.diz cs solvent flat today $ loop: $ @rota_maps ! transforms UP and DOWN density areas by the local 2-fold $ @aver_map ! averages the UP and DOWN density areas $ @rota_maps_2 ! transforms averaged UP density area to DOWN area $ ! and vice versa by improper '2-folds' $ @aver_map_2 ! averages already averaged UP and DOWN areas with the $ ! with the transformed ones $ @cell ! fills the cell from the averaged UP and DOWN areas applying $ ! crystal symmetry operations $ @fftfc $ @fourier $ delete hb:[scr]hkl.'tmp'.* $ ii = ii + 1 $ if ii .lt. 13 then goto loop $ exit $ @create_loop $ jj = jj + 1 $ if jj .lt. 1 then goto cycle The ROTA_MAPS.COM does the density transformations of the UP and DOWN areas into themselves. The density inside each area is rotated by the 2-fold axes, defined by polar angles (the matrix 2), around the specified center including additional translational shift, that in principle only changes the rotational center. ROTA_MAPS.COM: $ main > read file cs.diz map protein Perform UP area identity transformation. > read file [-.mask]up_mask.rmap map number > make map 2 from 1 copy > write file mol_up_01.rmap map 2 Perform UP area 2-fold rotation. > read file [-.mask]up_mask.rmap map number over 2 > set matrix 2 polar 44.644 353.353 180.0 > make map 2 from 1 rotate matrix 2 cent coordinate 23.26 77.76 52.21 - > translate vector -1.369 2.298 -5.035 > write file mol_up_02.rmap map 2 Perform DOWN area identity transformation. > read file [-.mask]down_mask.rmap map number over 2 > make map 2 from 1 copy > write file mol_down_01.rmap map 2 Perform DOWN area 2-fold rotation. > read file [-.mask]down_mask.rmap map number over 2 > set matrix 2 polar 45.105 353.519 180.0 > make map 2 from 1 rotate matrix 2 cent coordinate 46.2884 104.7735 51.3439 - > translate vector 1.028 -1.232 2.894 > write file mol_down_02.rmap map 2 > quit Both pairs of UP and DOWN density are then averaged in the AVER_MAP.COM. The resulting averaged densities MOL_UP.RMAP and MOL_DOWN.RMAP are transformed by the ROTA_MAPS_2.COM to DOWN and UP areas respectively applying transformation parameters obtained as described above. ROTA_MAPS_2.COM: $ set proc/name=rota_maps_2 $ main > read file mol_down.rmap map number Transform averaged density of DOWN to UP (local symmetry operation 1). > read file [-.mask]up_mask.rmap map number > set matrix 2 polar 90.876 82.011 180.093 > make map 2 from 1 rotate matrix 2 cent coordinate 34.63 91.39 51.10 - > translate vector 0.005 -0.336 0.106 > write file mol_up_2.rmap map 2 Transform averaged density of DOWN to UP (local symmetry operation 2). > read file [-.mask]up_mask.rmap map number over 2 > set matrix 2 polar 45.163 170.874 179.392 > make map 2 from 1 rotate matrix 2 cent coordinate 34.63 91.39 51.10 - > translate vector -0.445 0.345 0.370 > write file mol_up_3.rmap map 2 Read the averaged UP density and overrides the averaged DOWN density. > read file mol_up.rmap map number over 1 Transform averaged density of UP to DOWN (local symmetry operation 1). > read file [-.mask]down_mask.rmap map number > set matrix 2 polar 90.147 81.995 180.099 > make map 2 from 1 rotate matrix 2 cent coordinate 34.63 91.39 51.10 - > translate vector -0.2 -0.3 0.25 > write file mol_down_2.rmap map 2 Transform averaged density of UP to DOWN (local symmetry operation 2). > read file [-.mask]down_mask.rmap map number over 2 > set matrix 2 polar 45.294 171.096 179.695 > make map 2 from 1 rotate matrix 2 cent coordinate 34.63 91.39 51.10 - > translate vector 0.247 -0.406 0.294 > write file mol_down_3.rmap map 2 > quit Then density maps from files MOL_UP.RMAP, MOL_UP_1.RMAP and MOL_UP_2.RMAP are averaged by the procedure AVER_MAP_2.COM into AVER_UP.RMAP. The same was done for the DOWN density maps. The next step is to build the crystal unit cell with the density maps AVER_UP.RMAP and AVER_DOWN.RMAP with CELL.COM and to transform the density to reciprocal space and back and repeat the whole cycle until the density map R-factor converges. 7.3. Four masks The electron density averaging was then applied in three ways: - averaging the solvent flattened density - averaging the density obtained from combined phases of the model and solvent flattened density - averaging the density obtained from combined phases of the model and heavy atom derivatives 7.3.1. Averaging the solvent flattened density When enough secondary structure elements were recognized and built in the UP and DOWN averaged electron density, it became possible to split the asymmetric unit in 4 parts, for one molecule each. First model was built into the UP area of the density without paying attention to which molecule a secondary structure element might belong. The aim was to interpret as much density as possible. From the resulting model the tetramer was built by applying the 2-fold for the local UP area to create the dimer and from the dimer, by applying a local symmetry operation between UP and DOWN area, the tetramer: Read 4 identical segments (for one molecule each). MAIN> read file model.rdi coordinate diamond MAIN> read file model.rdi coordinate diamond MAIN> read file model.rdi coordinate diamond MAIN> read file model.rdi coordinate diamond Transform coordinates of the segment number 2 by the UP 2-fold. MAIN> set matrix 2 polar 44.644 353.353 180.0 MAIN> rotate atom select segment number 2 end matrix 2 - MAIN> cent coordinate 23.26 77.76 52.21 MAIN> translate atom select segment number 2 end vector -1.369 2.298 -5.035 Copy the atomic coordinates of molecule 2 to the molecule 4 in order to form another dimer consisting of molecules 3 and 4. MAIN> copy coordinates select segment number 2 end - MAIN> select segment number 4 end Transform the dimer coordinates of molecules 3 and 4 by the operation UP to DOWN. MAIN> set matrix 2 polar 90.876 82.011 180.093 MAIN> rotate atom select segment number 3 4 end matrix 2 - MAIN> center coordinate 34.63 91.39 51.10 MAIN> translate atom select segment number 3 4 end - MAIN> vector 0.005 -0.336 0.106 Display the tetramer and delete the interfering parts of the model (twice interpreted areas). The parts for deletion should be chosen so that the rest of the monomer belongs to one molecule only and not to any of local symmetry related molecules. Delete atoms by the sequence criteria so that the interfering parts are removed from all 4 molecules at once as shown in the example. The command sentence, MAIN> delete atom select sequence 52 63 end deletes all the residues between residues with sequence names 52 and 63 including the residues 52 and 63. Write the remaining tetramer to the file TETRA.PDB in the PDB format: MAIN> write file tetra.pdb coordinates pdb The transformation parameter search: After starting MAIN read the coordinates of the tetramer and the solvent flattened density. MAIN> read file tetra.pdb coordinates pdb MAIN> read file pic$exa:[aver.cs.dat]uruohgko_sf_zero.rmap map number Generate density points at the atomic coordinates. By this operation each density point obtains the index of the ancestor atomic number. MAIN> make point map 1 200 200 200 atom select all end The generated density points have a density that is interpolated from the neighboring 8 grid points. Since the program doesn't include a command to change the density value of a point, the simplest solution is to write the density points to a file (POINT.DENS) MAIN> write file point.dens points density select all end and then with an editor set the column with density values to 1.00. MAIN> edt point.dens where the whole column of density values is set to 1.00. from: 1 0 0 8 36.292 58.836 52.273 0.47 36.00 57.00 51.00 1 to: 1 0 0 8 36.292 58.836 52.273 1.00 36.00 57.00 51.00 1 The file POINT.DENS should be read again by initializing the points array. MAIN> read file point.dens points density initialize Turn OFF the atomic energy terms and set the DENSITY ENERGY calculation to the density POINTS and the MAP to number 1. MAIN> energy all off density on density map 1 density point Optimize the positions of molecules according to their correlation with the density as described in UP and DOWN for each molecule separately. MAIN> energy select segment number 1 end translate 1. 0. 0. mat 1 .... MAIN> mini by body select segment number 1 end Transform the corresponding atoms by the optimized parameters just obtained. The optimized rotational matrix is stored in matrix number 8 and the translational parameters had to be included interactively (the translational components are given only as an example). MAIN> rotate atom select segment number 3 end matrix 8 MAIN> translate atom select segment number 3 end vector 0.1 -0.4 0.8 Write the modified tetramer coordinates to a file. Prepare masks from the model: When the starting model is rather small and the molecules are quite close together, the masks shouldn't be derived applying equal atomic radii for all the atoms. The atoms closer to the neighboring molecules should have smaller radii and the ones further apart, larger. So the overlap of masks from different molecules is prevented and masks can expand far away in the regions of empty space. The following MAKE_MASKS.COM procedure is applied already after several cycles of crystallographic refinement using X-PLOR, but is in essence still the same as the first one. The difference is only in the input file name. MAKE_MASKS.COM: The strategy is the following. Select the atoms of the contact region in several keys. Assign the atoms closer than 2A to any neighboring molecule atom to the key null, the ones closer then 4A to the key one and so on. The masks are then defined so that the atoms from the null key don't contribute to the masks, and the others can cover only half of the distance to their closest atom of the neighboring molecule. The procedure stops at the 6A radius which is also the maximal atomic radius applied for mask generations. $ main Read all coordinates. > read file cs_ncs_15.xpl coordinate xplor > key mol1 select segment number 1 end > key mol2 select segment number 2 end > key mol3 select segment number 3 end > key mol4 select segment number 4 end Assign atom to keys: > key null select around key mol1 distance 2.0 end > key null select null .or. around key mol2 distance 2.0 end > key null select null .or. around key mol3 distance 2.0 end > key null select null .or. around key mol4 distance 2.0 end > key two select around key mol1 distance 4.0 end > key two select two .or. around key mol2 distance 4.0 end > key two select two .or. around key mol3 distance 4.0 end > key two select two .or. around key mol4 distance 4.0 end > key three select around key mol1 distance 6.0 end > key three select three .or. around key mol2 distance 6.0 end > key three select three .or. around key mol3 distance 6.0 end > key three select three .or. around key mol4 distance 6.0 end > key four select around key mol1 distance 8.0 end > key four select four .or. around key mol2 distance 8.0 end > key four select four .or. around key mol3 distance 8.0 end > key four select four .or. around key mol4 distance 8.0 end > key five select around key mol1 distance 10.0 end > key five select five .or. around key mol2 distance 10.0 end > key five select five .or. around key mol3 distance 10.0 end > key five select five .or. around key mol4 distance 10.0 end Read the MAP header. > read file head.rmap map number > make map 2 from 1 initialize -9999 around 10 select mol1 end Mask the region around molecule 1. > make map 2 atom 1.0 select mol1 .and. two .and. .not. null end > make map 2 atom 2.0 select mol1 .and. three .and. .not. two end > make map 2 atom 4.0 select mol1 .and. four .and. .not. three end > make map 2 atom 5.0 select mol1 .and. five .and. .not. four end > make map 2 atom 6.0 select mol1 .and. .not. five end > write file mask_1.rmap map 2 number > make map 2 from 1 initialize -9999 around 10 select mol2 end Mask the region around molecule 2. > make map 2 atom 1.0 select mol2 .and. two .and. .not. null end > make map 2 atom 2.0 select mol2 .and. three .and. .not. two end > make map 2 atom 4.0 select mol2 .and. four .and. .not. three end > make map 2 atom 5.0 select mol2 .and. five .and. .not. four end > make map 2 atom 6.0 select mol2 .and. .not. five end > write file mask_2.rmap map 2 number > make map 2 from 1 initialize -9999 around 10 select mol3 end Mask the region around molecule 3. > make map 2 atom 1.0 select mol3 .and. two .and. .not. null end > make map 2 atom 2.0 select mol3 .and. three .and. .not. two end > make map 2 atom 4.0 select mol3 .and. four .and. .not. three end > make map 2 atom 5.0 select mol3 .and. five .and. .not. four end > make map 2 atom 6.0 select mol3 .and. .not. five end > write file mask_3.rmap map 2 number > make map 2 from 1 initialize -9999 around 10 select mol4 end Mask the region around molecule 4. > make map 2 atom 1.0 select mol4 .and. two .and. .not. null end > make map 2 atom 2.0 select mol4 .and. three .and. .not. two end > make map 2 atom 4.0 select mol4 .and. four .and. .not. three end > make map 2 atom 5.0 select mol4 .and. five .and. .not. four end > make map 2 atom 6.0 select mol4 .and. .not. five end > write file mask_4.rmap map 2 number > quit Averaging is done as in the case of cathepsin B (using the real*4 maps and the LONG procedure). A slight modification is introduced by the creation of 4 masks and the rotated maps are averaged already in the ROTA_MAPS.COM so that the AVER_MAPS.COM is not needed. The following file DO_ALL.COM is taken from the phase combination averaging procedure. The difference to the procedure that averages the solvent flattened density is at the beginning. In the case of the phase combinations the initial density is calculated by combining the model and heavy atom (or solvent flattened density) phases while in the case of averaging the solvent flattened density that one is simply taken as the input to the cyclic averaging procedure starting at label cycle. $ setup protein $ set def main$exa:[aver.cs.4mol.phco] $ tmp :== aaa $ jj = 0 $ ii = 0 $ !goto here $ set def [-.mask] $ @make_masks $ set def [-.phco] $ @create_p1sf $ @create_merge $ @fftfc_atom $! here: $ @phasco $ @fourier_atom $ copy cs.diz cs_atom.diz ! saves the combined phase density $ cycle: $ @create_p1sf $ @create_merge $ ii = 0 $ loop: $ @rota_maps ! transforms density from local symm. equivalent molecules $ ! to the mask of molecule 1 symmetry 4 $ @cell ! fills the cell from the averaged densities applying $ ! crystal symmetry operations $ @fftfc $ @fourier $ ii = ii + 1 $ pu scr: $ if ii .lt. 13 then goto loop $ exit $ @create_loop $ jj = jj + 1 $ if jj .lt. 1 then goto cycle The ROTA_MAPS.COM including averaging differs from the ROTA_MAPS.COM procedure applied by cathepsin B in that it operates with 3 maps and not with 2. Map number 1 is the unit cell density, map number 2 is the place where all transformed densities of each molecule are added together, and map number 3 is the density where each transformed density is stored. In the case of cathepsin B maps number 2 are written to the files and averaged separately, while here only already averaged density maps are stored. ROTA_MAPS.COM: $ set proc/name=cs_rota_maps $ main_map > read file cs.diz map protein > set variable xcenter = 0.0 > set variable ycenter = 0.0 > set variable zcenter = 0.0 molecule 1 > read file [-.mask]mask_1.rmap map number > make map 2 from 1 copy > @[-.set]set_tran_1_to_2 > read file [-.mask]mask_1.rmap map number > make map 3 from 1 rotate matrix 2 cent coordinate xcenter ycenter zcenter - > translate vector xt yt zt The transformed density (map 3) is added to the map 2. > make map 2 from 3 add > @[-.set]set_tran_1_to_3 > read file [-.mask]mask_1.rmap map number over 3 > make map 3 from 1 rotate matrix 2 cent coordinate xcenter ycenter zcenter - > translate vector xt yt zt > make map 2 from 3 add > @[-.set]set_tran_1_to_4 > read file [-.mask]mask_1.rmap map number over 3 > make map 3 from 1 rotate matrix 2 cent coordinate xcenter ycenter zcenter - > translate vector xt yt zt > make map 2 from 3 add > scale 0.25 map 2 > write file scr:aver_1.rmap map 2 number And so on for molecules 2, 3 and 4. > read file [-.mask]mask_2.rmap map number over 2 > make map 2 from 1 copy > @[-.set]set_tran_2_to_1 > read file [-.mask]mask_2.rmap map number > make map 3 from 1 rotate matrix 2 cent coordinate xcenter ycenter zcenter - > translate vector xt yt zt > make map 2 from 3 add 7.3.2. Phase combination The phase combination combines the phases from the heavy atom derivatives or from a solvent flattened density and the atomic model. Phase combination is controlled by the PHASCO.COM. PHASCO.COM: $ mir_hlc :== cs.abcd ! The heavy atom phases file $ pro_a_in :== scr:cs_atom.mrg ! Atomic model structure factors $ pro_b_out :== scr:cs.mrg ! The combined phases for FOURIER map When phases are in the form of the "G5" file they have to be converted with the COPABCD.EXE to a Hendrickson and Lattman coefficient file that can be read by PHASCO. $ ASSIGN/USER_MODE 'MIR_HLC' FOR001 $ ASSIGN/USER_MODE hb:[scr]COP.'TMP' FOR002 $ RUN DISK$HUBER:[JD.DEISENHOF]COPABCD When phases are available in the form of the structure factors, as calculated with the P1SF from an electron density file and later merged into a PROTEIN file, the required procedure is a bit longer: Command procedure to extract Hendrickson-Lattmann coefficients from a merged protein reflection file. The averaged map was backtransformed and phases merged with COPY A B to the protein_in-file. This command procedure calculates the HL-coefficients in the ITYP = 2 step and extracts them in the ITYP=8 step. $! $!**********************PHC2.COM******************** $! $ PROTEIN_IN := scr:cs.mrg $ PROTEIN_OUT := scr:cs.tmp $ HL_COEFF := cs.abcd $ ASSIGN /USER_MODE 'PROTEIN_IN' FOR011 $ ASSIGN /USER_MODE 'PROTEIN_OUT' FOR012 $ RUN disk$huber:[JD.EXE]PHASCO 2 11 12 0 20. 3.0 10 NATI 0.1 0 50. 58.3,1.60 500. 10 1 2 1 1. 0 1. 0. 1.0,0 $! $! ****************HLEXTR.COM************************** $! $ ASSIGN /USER 'PROTEIN_OUT' FOR004 $ ASSIGN /USER 'HL_COEFF' FOR007 $ RUN disk$huber:[SCHNEIDER.PHASCO]PHASCO 8 4 0 7 The phases in the form of CS.ABCD are later combined to a PROTEIN reflection file, that serves as input for the electron density calculation. $! *** RUN PHASCO: PHME $ @create_ph3 The command file CREATE_PH3.COM extracts, similar as already described by the cathepsin B, from the files information about the resolution range (LOOP.DAT, RESOL.DAT) and Fcalc scaling factors (SIN_K.DAT) and combines it with the pattern file (PH3.ORIG) into the PHASCO input file PH3.INP. The PH3.INP prepares the phase combination later done PH2.INP. CREATE_PH3.COM: $ ass/user loop.dat for001 $ ass/user sin_k.dat for002 $ ass/user resol.dat for003 $ ass/user ph3.orig for004 $ ass/user ph3.inp for007 $ r create_ph3 PH3.ORIG: 3 3 4 2 20.0 3.0 10 NATI 0.0 0 50. 18.40000 35.10000 300. 10 (3I4,4A4) $ ASSIGN/USER_MODE cs.abcd FOR002 $ ASSIGN/USER_MODE 'PRO_A_IN' FOR003 $ ASSIGN/USER_MODE scr:PHME_OUT.TMP FOR004 $ ass/user ph3.inp for005 $ RUN DISK$HUBER:[JD.DEISENHOF]PHASCO After the phases of heavy atoms are combined with atomic model phases they have also to be written in the PROTEIN format reflection file (CS.MRG). The PH2.INP is prepared similar as the PH3.INP. $! $! *** RUN PHASCO: PHC2 $! $ @create_ph2 CREATE_PH2.COM: $ ass/user loop.dat for001 $ ass/user sin_k.dat for002 $ ass/user resol.dat for003 $ ass/user ph2.orig for004 $ ass/user ph2.inp for007 $ r create_ph2 2 11 12 0 20. 3.0 10 NATI 0.0 0 50. 18.40000 35.10000 500. 10 1 2 1 1. 0 1. 0. 1.,0 $ ass/user ph2.inp for005 $ assign/user_mode SCR:phme_out.tmp for011 $ assign/user_mode 'pro_b_out' for012 $ run disk$huber:[jd.deisenhof]phasco $ delete scr:*.tmp;* 7.4. Preparing files for X-PLOR: Initially only the solvent flattened density was averaged in [AVER] using improper symmetry parameters obtained from the first model. When it seemed that the model could not be improved any more the first crystallographic refinement was done using X-PLOR. Since the model consists of several nonlinked peptide chains each segment must be stored in a separate file. In order to find the segments and write them to separate files the following procedure replaced the use of a text editor. Read the coordinates of molecule 1 and rename the segment to POLY MAIN> read file mol1.pdb coordinates pdb MAIN> rename select all end segment POLY Data from topology and parameter files for energy calculations are used MAIN> @main$u:get_top_par to create the covalent bonds between atoms according to the list of bonds stored in topology files. The whole chain becomes a network of covalent bonds. However the peptide bonds of the residues not linked during a modeling session are sometimes extremely long. According to these lengths of peptide bonds, the segments along the chain from N to C terminal could easily be identified. The distance criterion works properly when all connected residues follow each other. The covalent bonds in the length range between 2.0 and 200.0A are than written to the file BOND.LIST. MAIN> write file bond.list select bond range 2.0 200.0 end bond The list of the breaks is then introduced to the file WR_SEGM.COM with editor. The chains are selected by the sequence name criteria. WR_SEGM.COM: write file sg_01.xpl select sequence 30 52 end coordinate xplor write file sg_02.xpl select sequence 58 98 end coordinate xplor .... write file sg_14.xpl select sequence 801 810 end coordinate xplor return By running it and with the command MAIN>@wr_segm the whole list of segments, each in a separate file, is generated. From these files, with X-PLOR structure files for MOL1, MOL2, MOL3 and MOL4, segments are generated and hydrogens atoms are built to molecule 1. From the molecule 1 the remaining tetramer molecules are generated by applying the local symmetry operations with the command file READ_LOC.COM. READ_LOC.COM: $ main Molecule 1 has to be read four times. Each time another key (mol1, mol2, mol3 and mol4) is assigned to it and each key obtained another segment name (MOL1, MOL2, MOL3 and MOL4). > read file mol1.pdb coordinate xplor > key mol1 select all end > key old select all end > read file mol1.pdb coordinate xplor > key mol2 select .not. old end > key old select all end > read file mol1.pdb coordinate xplor > key mol3 select .not. old end > key old select all end > read file mol1.pdb coordinate xplor > key mol4 select .not. old end > key old select all end > rename select mol1 end segment MOL1 > rename select mol2 end segment MOL2 > rename select mol3 end segment MOL3 > rename select mol4 end segment MOL4 > @make_loc The MAKE_LOC.COM derives from coordinates of molecule 1 molecules 2, 3 and 4 by applying local symmetry operations (as created before with the RMS.COM). MAKE_LOC.COM: > copy coordinate select segment number 1 end select segment number 2 end > copy coordinate select segment number 1 end select segment number 3 end > copy coordinate select segment number 1 end select segment number 4 end Set the rotation matrix for (4) local symmetry operation and rotates atoms around the center (upper) 2 molecules. > @[-.set]set_tran_1_to_2 > rotate atom matrix 2 select segment number 2 end > translate atom vector xt yt zt select segment number 2 end > @[-.set]set_tran_1_to_3 > rotate atom matrix 2 select segment number 3 end > translate atom vector xt yt zt select segment number 3 end > @[-.set]set_tran_1_to_4 > rotate atom matrix 2 select segment number 4 end > translate atom vector xt yt zt select segment number 4 end > return After returning command back to READ_LOC.COM the four segments are written again with crystallographic weights of hydrogens set to 0.0. > set weigh select atom name H* end 0.0 > write file mol1.xpl select segment number 1 end coordinate xplor > write file mol2.xpl select segment number 2 end coordinate xplor > write file mol3.xpl select segment number 3 end coordinate xplor > write file mol4.xpl select segment number 4 end coordinate xplor > quit After crystallographical refinement the local symmetry relations have to be recalculated. The procedure described at the cathepsin B example, RMS.COM, had to be modified, since the RMS fit of two molecules didn't always converge to a satisfactorily low RMS difference value. The procedure was improved by introducing an initial guess that was close enough to the final solution. The initial guess consists of 3 rotational angles (X, Y', X'') followed by the translational component. Since the default center of rotation for the RMS fit procedure is set to the coordinate system origin (0.0, 0.0, 0.0) and the molecules are related to each other by an almost perfect 2-fold axes the translational components are relatively large. RMS.COM: $ main > read file cs_ncs_12.xpl coordinate pdb > key real select weigh 0.9 1.1 .and. name at CA end > key s1 select segment number 1 .and. real end > key s2 select segment number 2 .and. real end > key s3 select segment number 3 .and. real end > key s4 select segment number 4 .and. real end > open unit 10 file rms.dat write_enable > rms coordinate all select s1 end select s2 end - > initial_guess 98.255 -93.326 -81.834 -61.073 46.249 82.7 > show unit 10 matrices rms > rms coordinate all select s1 end select s3 end - > initial_guess 0.776 195.971 179.185 82.323 181.361 9.106 > show unit 10 matrices rms > rms coordinate all select s1 end select s4 end - > initial_guess -81.233 90.497 -98.484 117.289 134.304 109.991 > show unit 10 matrices rms > rms coordinate all select s2 end select s1 end - > initial_guess 81.835 93.326 98.972 -60.972 46.224 82.787 > show unit 10 matrices rms > rms coordinate all select s2 end select s3 end - > initial_guess 279.015 88.624 -99.285 116.362 136.299 109.739 > show unit 10 matrices rms > rms coordinate all select s2 end select s4 end - > initial_guess 180.402 163.44 -0.360 82.45 182.08 10.79 > show unit 10 matrices rms > rms coordinate all select s3 end select s1 end - > initial_guess 180.819 -195.970 -0.772 82.362 181.350 8.962 > show unit 10 matrices rms > rms coordinate all select s3 end select s2 end - > initial_guess 99.233 -89.424 -278.258 116.183 135.757 110.660 > show unit 10 matrices rms > rms coordinate all select s3 end select s4 end - > initial_guess 82.994 91.714 -262.666 -61.284 48.988 85.817 > show unit 10 matrices rms > rms coordinate all select s4 end select s1 end - > initial_guess 99.42 -88.992 82.078 115.575 136.14 111.172 > show unit 10 matrices rms > rms coordinate all select s4 end select s2 end - > initial_guess 0.078 -163.705 178.840 82.324 181.699 9.719 > show unit 10 matrices rms > rms coordinate all select s4 end select s3 end - > initial_guess 262.731 -91.822 -82.983 -61.203 49.212 85.607 > show unit 10 matrices rms > quit