# MCCE Mechanism

### <span class="C9DxTc ">An MCCE run is a 4-step procedure</span>

- <span class="C9DxTc ">Step 1: **M**odify PDB (file formatting)</span>
- <span class="C9DxTc ">Step 2: **C**onformer/Rotamer making </span>
- <span class="C9DxTc ">Step 3: **C**alculate energy look-up table </span>
- <span class="C9DxTc ">Step 4: **E**xtract microstates; Monte Carlo sampling of conformers at each pH or Eh </span>

<span class="C9DxTc ">MCCE program can run any steps providing the prerequisite files exist in the working directory. Files required and written out by the program are illustrated in this chart. This file flow chart shows the summary of file dependencies:</span>

<div class="CobnVe yMxPgf  aP9Z7e" id="bkmrk--2" style="text-align: justify;"></div>
[![MCCE-flowchart.png](https://mccewiki.levich.net/uploads/images/gallery/2025-04/scaled-1680-/mcce-flowchart.png)](https://mccewiki.levich.net/uploads/images/gallery/2025-04/mcce-flowchart.png)

### <span class="C9DxTc ">Step 1: Modify PDB</span>


#### <span class="C9DxTc ">Input files:</span>

- - <span class="C9DxTc ">PDB file: </span><span class="C9DxTc ">input structure file in [PDB format](https://pdb101.rcsb.org/learn/guide-to-understanding-pdb-data/pdb-overview)</span>
    - <span class="C9DxTc ">name.txt (optional) : </span><span class="C9DxTc ">residue and atom renaming rule</span>
    - <span class="C9DxTc ">list\_rot.gold (optional): </span><span class="C9DxTc ">hot residue spot definition</span>

#### <span class="C9DxTc ">Output files:</span>

- - <span class="C9DxTc ">acc.res: </span><span class="C9DxTc ">solvent accessibility of residues</span>
    - <span class="C9DxTc ">acc.atm: </span><span class="C9DxTc ">solvent accessibility of atoms</span>
    - <span class="C9DxTc ">new.tpl (not always created): </span><span class="C9DxTc ">parameter file template of unrecognized cofactors</span>
    - <span class="C9DxTc ">head1.lst (optionally used by step 2): </span><span class="C9DxTc ">summary of rotamer making policy of residues</span>
    - <span class="C9DxTc ">step1\_out.pdb (used by step 2): </span><span class="C9DxTc ">step 1 output file is in mcce extended pdb format</span>

<span class="C9DxTc ">Step 1 prepares an extended PDB file, suitable to be read into step 2. The input PDB file is in standard PDB format. It can have alternative side chain positions, but MCCE can't process alternative backbone positions. Alternative side chains are treated as side chain conformers. When side chain atoms are missing, MCCE will complete the side chain atoms at the torsion minimum. In this step several things will happen:</span>

- 1. <span class="C9DxTc ">With the instructions in the renaming rule file "name.txt", residues and atoms will be renamed so that a cofactor can be split into several independent ionizable groups (for example, heme can be divided into heme and two propionates) and several groups can be combined as one (for example, heme can be grouped with the axial ligands).</span>
    2. <span class="C9DxTc ">Identify unrecognized cofactors and interpret them as non-charged atom assemblies.</span>
    3. <span class="C9DxTc ">Complete the missing atoms in each known residue (note: If your file has individual residues with missing atoms the .tpl file for that residue will be used to add missing atoms).</span>
    4. <span class="C9DxTc ">Calculate the solvent accessible surface (SAS) area and strips off exposed water and salt (HOH, NO3 and SO4). The SAS threshold of stripping off water and salt is defined by $(H2O\_SASCUTOFF in run.prm or as user defined parameter in step1.py command).</span>
    5. <span class="C9DxTc ">Give warnings on geometry clashes between atoms not supposed to be bonded.</span>
    6. <span class="C9DxTc ">Identifies heme ligands and cys-cys disulfide bridge.</span>
    7. <span class="C9DxTc ">Extracts N terminus and C terminus of a chain.</span>

<span class="C9DxTc ">The renaming rule file "name.txt" instructs MCCE program to rename atom name, residue name, sequence number, and chain ID. Here are several sample lines in this file:</span>

```bash
# Symbol "*" in the first string is a wildcard that matchs any character.
# It means "do not replace" in the second string.
#
# The replace is accumulative in the order of appearing in this file.#
*****HEA******  *****HEM******
*****HEC******  *****HEM******
*CAA*HEM******  *****PAA******      extract PAA from heme
*CBA*HEM******  *****PAA******
*CGA*HEM******  *****PAA******
```

<span class="C9DxTc ">The line started with "#" and the line shorter than 30 characters are comment lines. For other lines, the first 30 characters should be two 14-character strings separated by exactly two spaces, and the rest of the line is comment field. A valid line instructs MCCE to replace string 1 with string 2. MCCE will match this string with position 13 to 26 of a coordinate line in the input PDB file. The symbol "\*" is the wild card that matches any character in strings. The replace action is accumulative and order sensitive. For example, The line</span>

`<span class="C9DxTc ">HETATM 1683  CAA HEC     1       1.317  -3.987  -1.685  1.00  0.00           C</span>`

<span class="C9DxTc ">will be renamed to</span>

`<span class="C9DxTc ">HETATM 1683  CAA HEM     1       1.317  -3.987  -1.685  1.00  0.00           C</span>`

<span class="C9DxTc ">then</span>

`<span class="C9DxTc ">HETATM 1683  CAA PAA     1       1.317  -3.987  -1.685  1.00  0.00           C</span>`

<span class="C9DxTc ">The output files of step 1 "acc.res" and "acc.atm" contain the solvent accessible surfaces of residues and atoms. In "acc.res", both absolute value and percentage of the solvent accessibility are listed.</span>

<span class="C9DxTc ">When an unrecognized cofactor is encountered in step 1, a parameter file will be created with name "new.tpl" and a warning message will be issued. The atom connectivity is guessed and all atoms are assumed to have charge 0. This file can be the starting point of making a parameter for a new cofactor. If the program is resumed from step 2 or 3, this file will be read in to define the molecule parameter.</span>

<span class="C9DxTc ">The file "head1.lst" lists the rotamer making policy of the residues. When $(ROT\_SPECIF) in "run.prm" is set to "t", this file will be used as instruction of rotamer making of step 2, otherwise, file "head1.lst" will be ignored. This file can be modified and will be effective if the program resumes from step 2. The rotamer making is also dependent on step 2 rotamer making level. If step 2 is running at level 1 (quick run), the rotation rotamers will not be made even though head1.lst rules say say. The behavior will change in a planned new version of MCCE.</span>

<span class="C9DxTc ">The file "step1\_out.pdb" is a formatted PDB file, which will be read in by step 2. This MCCE extended PDB format contains three more fields than a standard PDB file: charge, size and rotamer making history.</span>

<div class="GV3q8e yMxPgf  aP9Z7e" id="bkmrk--7" style="text-align: justify;"></div>### <span class="C9DxTc ">Step 2: Conformer/Rotamer making</span>


#### <span class="C9DxTc ">Input files:</span>

- - <span class="C9DxTc ">step1\_out.pdb: </span><span class="C9DxTc ">input structure file of step 2 in MCCE extended PDB format</span>
    - <span class="C9DxTc ">head1.lst (optional) : </span><span class="C9DxTc ">rotamer making policy of residues</span>

#### <span class="C9DxTc ">Output files:</span>

- - <span class="C9DxTc ">progress.log: </span><span class="C9DxTc ">progress report file, dynamically updated</span>
    - <span class="C9DxTc ">rot\_stat: </span><span class="C9DxTc ">rotamer making statistics, dynamically updated</span>
    - <span class="C9DxTc ">hvrot.pdb: </span><span class="C9DxTc ">heavy atom (without hydrogen atoms) rotamer pdb file</span>
    - <span class="C9DxTc ">head2.lst (optionally used by step 3): </span><span class="C9DxTc ">summary of rotamers made in step 2</span>
    - <span class="C9DxTc ">step2\_out.pdb (used by step 3): </span><span class="C9DxTc ">step 2 output file with multiple rotamers in extended pdb format</span>

<span class="C9DxTc ">Step 2 makes and optimizes rotamers from the structure in "step1\_out.pdb". In this step, the rotatable bonds (defined in parameter files) of each residue is rotated by the steps defined in "run.prm". Then the self Van Der Waals (VDW) potential (interaction among atoms of the same side chain excluding 1-2 and 1-3 interactions, and interaction between the side chain and backbone atoms) is calculated. Side chain rotamers with high self VDW potentials are deleted. Then the side chains are optimized with possible hydrogen bond partners. A number of repackings starting from randomized initial structures (one conformer from one residue) are performed to reduce side chain rotamers to those with low energy local packings. Ionization states are then created and protons are placed on side chains. At the end of side chain rotamer optimization, MD simulations are carried out locally to relax the structure.</span>

<span class="C9DxTc ">The input file step1\_out.pdb is the only essential file step 2 will use. You can modify this file to add, delete or edit residues without causing problems as long as the file observes MCCE extended PDB format. Sometimes a little editing of this file is necessary, for example, the terminal residues are not always correctly identified by MCCE due to "broken chains" caused by the missing residues in the PDB file.</span>

<span class="C9DxTc ">The file "head1.lst" provides residue specific rotamer making rule. It will be used only when $(ROT\_SPECIF) is set to be "t". The line of this file such as:</span>

<span class="C9DxTc ">`NTR A0003_ R t 06 S f  0.0 H t 06 M 000`</span>

<span class="C9DxTc ">is interpreted as "Rotate is true and 6 steps per bond, then swing is false (angle is 0 if any), then Hydroxyl relaxation is true and the maximum number of starting conformers per hydroxyl is 6, and the maximum total number of the conformers is not limited". If you want to investigate a specific site in details, change the step 6 to 12. But making 12 steps for many sites (&gt;30) are not recommended because it is expensive in terms of memory and CPU time.</span>

<span class="C9DxTc ">The file "progress.log" is a dynamically updated progress report file. It lists mainly the repacking progress.</span>

<span class="C9DxTc ">The file "rot\_stat" is an important file to review the rotamer making history. This file lists the number of rotamers of each residue. It is easy to tell what residues get rotamers and if the total number of rotamers is manageable (a 2000 conformer structure may need one day to run the next two steps, step 3 and 4).</span>

<span class="C9DxTc ">The file hvrot.pdb is a file at the same format as step1\_out.pdb and step2\_out.pdb but lacks hydrogen atoms. The main use of this file is to rename it to "step1\_out.pdb" and run step 2 again with "swing" rotamers instead of "rotate" rotamers. This provides a way to relax the structure by "swinging" the rotatable bond a little and reevaluate the rotamers. This feature is most for advanced users to calibrate hydrogen bond directed rotamer making algorithm and MD relaxation subroutine.</span>

<span class="C9DxTc ">The file "head2.lst" is a summary of the rotamers made in step 2, and is not used by step 3.</span>

<span class="C9DxTc ">The file "step2\_out.pdb" is in the MCCE extended PDB format, and it connects step 2 and step 3.</span>

<div class="GV3q8e yMxPgf  aP9Z7e" id="bkmrk--9" style="text-align: justify;"></div><div class="PPhIP rviiZ" id="bkmrk--4" jsname="haAclf"><div aria-describedby="h.vr6f2x2gqzdj_l" aria-disabled="false" aria-hidden="true" aria-label="Copy heading link" class="U26fgb mUbCce fKz7Od LRAOtb Znu9nd M9Bg4d" data-tooltip="Copy heading link" data-tooltip-horizontal-offset="0" data-tooltip-position="top" data-tooltip-vertical-offset="12" jsaction="click:cOuCgd; mousedown:UX7yZ; mouseup:lbsD7e; mouseenter:tfO1Yc; mouseleave:JywGue; focus:AHmuwe; blur:O22p3e; contextmenu:mg9Pef;" jscontroller="mxS5xe" jsshadow="" role="presentation"><svg class="OUGEr QdAdhf" fill="currentColor" focusable="false" height="22px" viewbox="0 0 24 24" width="22px"></svg>  
</div></div>### <span class="C9DxTc ">Step 3: Calculate energy lookup table</span>


#### <span class="C9DxTc ">Input files:</span>

- - <span class="C9DxTc ">step2\_out.pdb: </span><span class="C9DxTc ">input structure file of step 3 in MCCE extended PDB format</span>

#### <span class="C9DxTc ">Output files:</span>

- - <span class="C9DxTc ">progress.log: </span><span class="C9DxTc ">progress report file, dynamically updated</span>
    - <span class="C9DxTc ">energies (directory of energy lookup table used by step 4): </span><span class="C9DxTc ">opp files, each opp file is pairwise interaction to a conformer</span>
    - <span class="C9DxTc ">head3.lst (used by step 4): </span><span class="C9DxTc ">list of conformers and self energy of conformers. It is used by step 4.</span>
    - <span class="C9DxTc ">step3\_out.pdb: </span><span class="C9DxTc ">step 3 output file with multiple rotamers in extended pdb format</span>

<span class="C9DxTc ">Step 3 calls Poisson Boltzmann equation solver, DelPhi, to calculate reaction field energy and electrostatic pairwise interaction. The result is stored as together with Van dDer Waals interactions as one file per conformer. These files have extension "opp" and are located under directory energies. The self-energy terms (not dependent on side chains of other residues) of conformers are listed in file "head3.lst" The progress is dynamically updated is file "progress.log".</span>

<span class="C9DxTc ">The file "progress.log" reports the progress of DelPhi calculation, which can be used to estimate the total time of this step. There will be two parts of DelPhi calculations: the first is pairwise calculation and the second part is the reaction field energy calculation. It takes less time on reaction field energy calculation than on pairwise calculation.</span>

<span class="C9DxTc ">The directory "energies" holds the pairwise conformer interaction lookup table as files with extension '.opp' and starting with the conformer\_id, e.g. GLU02B0012\_002.opp.</span>

<span class="C9DxTc ">These files are header-less, but following is each column description:</span>

```
          column #:         1           2          3            4                   5                 6
          description:    conf#        name     corr_el      vdw_pwise          delphi_el        post_bdry_corr_el    (kcal/mol)
```

The file "head3.lst" contains self energy of each conformer and control flags of step 4. The flag is: "t" for fixed occupancy 0 or 1, or "f" for free to sample. The energy unit is Kcal/mol.

<span class="C9DxTc ">The file "step3\_out.pdb" is an extended pdb file with multiple conformers. The conformer number is sorted to be continuous and consistent with the conformer numbers in file "head3.lst" and step 4 output file fort.38. This file is identical to "step2\_out.pdb" if "step2\_out.pdb" is an unmodified file created by step 2.</span>

<div class="GV3q8e yMxPgf  aP9Z7e" id="bkmrk--11" style="text-align: justify;"></div><div class="PPhIP rviiZ" id="bkmrk-the-vdw_pwise-%28van-d" jsname="haAclf"><div aria-describedby="h.n88xj96o6b9i_l" aria-disabled="false" aria-hidden="true" aria-label="Copy heading link" class="U26fgb mUbCce fKz7Od LRAOtb Znu9nd M9Bg4d" data-tooltip="Copy heading link" data-tooltip-horizontal-offset="0" data-tooltip-position="top" data-tooltip-vertical-offset="12" jsaction="click:cOuCgd; mousedown:UX7yZ; mouseup:lbsD7e; mouseenter:tfO1Yc; mouseleave:JywGue; focus:AHmuwe; blur:O22p3e; contextmenu:mg9Pef;" jscontroller="mxS5xe" jsshadow="" role="presentation"><div class="VTBa7b MbhUzd" jsname="ksKsZd">The vdw_pwise (Van der Waals pairwise potential) term in opp files, vdw0 (conformer internal vdw potential) and vdw1 (conformer to backbone vdw interaction potential) can be recalculated by command vdw_pw.py.</div><svg class="OUGEr QdAdhf" fill="currentColor" focusable="false" height="22px" viewbox="0 0 24 24" width="22px"></svg></div></div>### <span class="C9DxTc ">Step 4: Extract microstates; simulate pH or Eh titration w/Monte Carlo sampling</span>


#### <span class="C9DxTc ">Input files:</span>

- - <span class="C9DxTc ">energies (directory):</span><span class="C9DxTc "> energy lookup table for pairwise interaction between conformers</span>
    - <span class="C9DxTc ">head3.lst: </span><span class="C9DxTc ">self-energy of conformers and Monte Carlo sampling flags of conformers</span>

#### <span class="C9DxTc ">Output files:</span>

- - <span class="C9DxTc ">mc\_out: </span><span class="C9DxTc ">progress of Monte Carlo sampling and energy tracing</span>
    - <span class="C9DxTc ">fort.38: </span><span class="C9DxTc ">conformer occupancies</span>

<span class="C9DxTc ">Step 4 is a titration simulation by Monte Carlo sampling. The Monte Carlo sampling is performed at specified set of pH/Eh. At each titration point, there will be several (predefined in "run.prm", the default is 6) independent samplings. Each sampling goes through annealing, reducing, and equilibration stages. Statistics of conformer occupancy is only done at equilibration statge. Yifan's Monte Carlo subroutine will check early convergence and quit sampling early to save time. The result is reported as conformer occupancy in file "fort.38".</span>

<span class="C9DxTc ">The file "mc\_out" is the progress report of Monte Carlo sampling. It contains running energy tracing which can be used to calculate the average E or enthopy of the system, or verify if the system is trapped at local energy minima. By "grep Sg mc\_out", you can find the standard deviation of independent samplings.</span>