pKa Prediction with Macro-pKa

Tutorial Created with Software Release: 2024-1
Topics: Catalysis & Reactivity, Pharmaceutical Formulations
Methodology: Molecular Quantum Mechanics
Products Used: Jaguar, Jaguar, MS Maestro, Maestro

Tutorial files

71.2 MB
This tutorial is written for use with a 3-button mouse with a scroll wheel.
Words found in the Glossary of Terms are shown like this: Workspacethe 3D display area in the center of the main window, where molecular structures are displayed

 

Tip: You can hover over a glossary term to display its definition. You can click on an image to expand it in the page.
Abstract:

 

This tutorial teaches how to carry out DFT-based pKa calculations with the Macro-pKa workflow and to analyze the results it produces. This workflow is recommended for treating molecular systems in which multiple functional groups can be protonated or deprotonated at the same pH. It computes micro-pKa’s associated with the protonation or deprotonation of individual functional groups of the molecule, and then combines these micro-pKa’s into a final macro-pKa value which should be directly comparable to an experimentally measured pKa. The Macro-pKa workflow can also handle molecules which exist as a mixture of tautomers and which can be protonated or deprotonated multiple times in different pH intervals. In addition to predicting micro- and macro-pKa’s of molecular systems, Macro-pKa outputs tautomer populations at any given pH. Equipped with a built-in conformational search, the workflow is capable of accounting for conformational effects.

 

Tutorial Content

  1. Introduction to pKa

  1. Creating Projects

  1. Application of Macro-pKa to methylamine

  1. Application of Macro-pKa to guanine

  1. Conclusions and References

  1. Glossary of Terms

1. Introduction to pKa

pKa, acidity, and protonation forms

The equilibrium constant that describes the protonation/deprotonation equilibrium

P ⇄ D + H+ (1)

of a simple chemical compound which can exist in either protonated (P) or a deprotonated (D) form in an aqueous solution is known as Ka:

Ka = [D][H+]/[P], (2)

where the square brackets stand for concentrations (or rather activities) of the corresponding chemical species.

The negative base 10 logarithm of Ka is denoted pKa:

pKa = -log Ka. (3)

The smaller the pKa of a given compound the easier this compound yields the proton. Neutrally charged molecules HA that dissociate into anion A- and the proton are known as acids:

HA ⇄ A- + H+ (4)

Strong acids P (such as hydrochloric acid) have pKa < 0. Weak acids (such as carboxylic acids) typically have pKa of 4-6. Very weak acids (such as amides and alcohols) have pKa > 10. Neutrally charged molecules B which can accept a proton and thus form cations BH+ are known as bases. To be consistent with equations (1)-(4) it is better to write the equilibrium for bases with the cation BH+ on the left:

BH+ ⇄ B + H+ (5)

The larger the pKa characterizing the dissociation of the cation BH+ the stronger the base B. Nitrogen-containing heterocycles and aromatic amines are weak bases with pKa typically between 3 and 6. Aliphatic amines are sufficiently strong bases with pKa’s normally ranging from 9 to 11.

If the molecule possesses multiple protonatable or deprotonatable functional groups, it can be characterized by multiple pKa values. A single molecule can act as an acid and a base. If, being in a neutral state, it has a functional group that can lose a proton, it is an acid. But if it also has a functional group that accepts a proton it is also a base. Some examples of typical acids, bases, and ambivalent molecules are shown in Figure 1.

Figure 1. Examples of neutral molecules which are acids, bases, as well as acids and bases at the same time.

Computational pKa prediction is an important task in molecular modeling projects. pKa values are calculated because simulation of chemical compounds or chemical reactions in aqueous solutions requires specification of protonation forms that would be present in the corresponding experimental conditions at a given pH. For example, it is important to know whether one should model compound XH, its protonated form XH2+, or its deprotonated form X- at a given pH. And it is precisely the knowledge of predicted pKa values that can be converted into information about pH-dependent populations of protonation forms X-, XH, and XH2+.

Micro- and macro-pKas

When multiple protonated P1, P2, … and deprotonated D1, D2, … tautomers can coexist at the same pH the equilibrium constant becomes

Ka = ([D1] + [D2] + …)[H+]/([P1] + [P2] + …). (6)

Each of the protonation forms P1, P2, … has charge C, while each of the deprotonation forms D1, D2, …   has charge C-1. So, the constant Ka characterizes the transition between charges C and C-1. If multiple charge transitions are possible for this molecule, it will be characterized by several Ka values, which are typically marked by subscripts 1, 2, 3 … So that Ka1 means the first known transition between charges C and C-1, Ka2 describes the transition between charges C-1 and C-2 and so on.

The constant (6) which describes the equilibrium involving multiple tautomers is known as a macro-Ka, and its logarithmic equivalent is called macro-pKa. The word “macro” means that the single constant characterizes proton exchanges among multiple protonation states P1, P2, …, D1, D2, … One can think of the so-called micro-constant Ka,ij

Ka,ij = [Dj][H+]/[Pi] (7)

which corresponds to the proton exchange among two protonation forms Pi and Dj. The micro-constant Ka,ij is typically not measurable in experiment since one cannot isolate the forms Pi and Dj from all other forms and measure the protonation equilibrium between these two forms only. The following relation connects macro-Ka with micro-pKa’s:

(8)

In the limit that there is only a single possible protonated form P1 (Nprot = 1) and a single possible deprotonated form D1 (Ndeprot = 1) there is only one possible micro-pKa which is equivalent to the macro-pKa.

Notice that the Schrödinger workflow which predicts macro-pKa’s is called Macro-pKa (with first letter capitalized). For historical reasons the workflow which predicts micro-pKa’s is called Jaguar pKa (and not Micro-pKa).

Micro- and macro-pKa calculations in practice

Prediction of micro-pKa constants can be done with the Jaguar pKa workflow. In order to launch a Jaguar pKa calculation, the user is required to select pKa (or active) atoms, i.e. define the protonation/deprotonation process manually. Micro-pKa calculations can be useful in the following scenarios:

  1. There is only one protonatable/deprotonatable functional group in the molecule and no tautomeric transformations are likely.
  2. There are several protonatable/deprotonatable functional groups in the molecule, but the user is certain that these groups participate in protonation equilibria in different pH ranges. In other words, the micro-pKa’s of these groups are “well-separated”, with the micro-pKa difference being 2-3 and larger. The user should be certain which of the functional groups with well-separated pKa’s is responsible for the pKa signal observed in the experiment, if comparison with the experiment is planned.

For a tutorial introducing the prediction of micro-pKa constants with Jaguar pKa, see pKa Predictions with Jaguar pKa.

Examples of molecules for which the use of Jaguar pKa is recommended:

This molecule contains a single functional group that can participate in protonation/deprotonation equilibria. The choice of the pyridinic N atom as the pKa atom is straightforward.

This molecule contains two functional groups that can dissociate into protons, but from chemical intuition it is obvious that the first group to dissociate will be the carboxylic group COOH (at pH about 4-5), while the phenolic group OH will remain intact at this pH. The phenolic group should dissociate at a much higher pH (9-10), after the COOH group got completely ionized. So, the COOH and OH groups are “well- separated”. They do not interfere in each other’s dissociation process, and it is also easy to pick the pKa atom.

Jaguar pKa calculations are not recommended for the following cases:

  1. There are several protonatable/deprotonatable functional groups in the molecule, and it is not clear whether their dissociations occur in clearly separated pH regions.
  2. The molecule contains many protonatable/deprotonatable groups and potentially multiple protonation states. The risk of selecting an incorrect pKa atom, and consequently working with unphysical protonation state becomes too high.
  3. The molecule is prone to tautomerism.  

The table below shows examples of molecules for which the use of Jaguar pKa is not recommended. Use the Macro-pKa workflow for these molecules instead:

This molecule contains two protonatable nitrogen atoms which are expected to be protonated in the same pH range. Two micro-pKa’s associated with these atoms will have to be inserted into formula (8) to produce the final macro-pKa. Use Macro-pKa to avoid making mistakes when handling cases with several micro-pKa’s. Macro-pKa computes all necessary micro-pKa’s and combines them into a macro-pKa automatically.

This molecule contains two functional groups that can dissociate into protons, but it is not clear which of the two groups – carboxylic COOH or phenolic OH is going to dissociate first. The three electron-withdrawing F-atoms in the benzene ring make the phenol part of the molecule much more acidic in comparison with plain phenol. At the same time, the dissociation of the COOH group, separated from the electron withdrawing groups by a long methylene chain, is not likely to be greatly affected. So, one can assume that the micro-pKa’s associated with the COOH and OH groups are no longer well-separated.

This molecule can be drawn in several tautomeric forms, and it is not obvious which of them is going to be the most stable. The form shown on the left is actually unstable. Picking a wrong (unstable) tautomeric form for a micro-pKa calculation may result in pKa predictions that cannot be meaningfully compared to experimentally measured pKa. Use Macro-pKa which automatically considers all reasonable tautomers.

In contrast to Jaguar pKa, the Macro-pKa workflow can be used in all scenarios. Macro-pKa does not require the user to select pKa atoms, since it will consider all possible active atoms automatically. Therefore, it should be more convenient to set up Macro-pKa calculations than Jaguar pKa calculations. Macro-pKa results are predicted viaDFT-based calculations with a machine-learned empirical correction and are directly comparable to experimentally measured pKa’s. User mistakes are also less common when setting up Macro-pKa calculations. The only tangible downside to using Macro-pKa is its potentially higher computational cost. Macro-pKa computes micro-constants for all reasonable tautomers, as opposed to Jaguar pKa, which focuses on a tautomeric equilibrium or a few equilibria selected by the user.

2. Creating Projects

At the start of the session, change the file path to your chosen Working Directorythe location where files are saved in MS Maestro to make file navigation easier. Each session in MS Maestro begins with a default Scratch Projecta temporary project in which work is not saved. Closing a scratch project removes all current work and begins a new scratch project, which is not saved. A MS Maestro project stores all your data and has a .prj extension. A project may contain numerous entries corresponding to imported structures, as well as the output of modeling-related tasks. Once a project is saved, the project is automatically saved each time a change is made.

Structures can be built in MS Maestro or can be imported using File > Import Structures (or drag-and-dropped), and are added to the Entry Lista simplified view of the Project Table that allows you to perform basic operations such as selection and inclusion and Project Tabledisplays the contents of a project and is also an interface for performing operations on selected entries, viewing properties, and organizing structures and data. The Entry Lista simplified view of the Project Table that allows you to perform basic operations such as selection and inclusion is located to the left of the Workspacethe 3D display area in the center of the main window, where molecular structures are displayed. The Project Tabledisplays the contents of a project and is also an interface for performing operations on selected entries, viewing properties, and organizing structures and data can be accessed by Ctrl+T (Cmd+T) or Window > Project Table if you would like to see an expanded view of your project data.

OR

  1. Double-click the Maestro or Materials Science icon

Note: This Jaguar workflow can be performed in Maestro or Materials Science Maestro. Use whichever interface you are comfortable with or typically use for your projects.

Figure 2-1. Change Working Directory option.

  1. Go to File > Change Working Directory
  2. Find your directory, and click Choose
  3. Pre-generated experimental files are included to be used in the analysis stage. Download the zip file here: schrodinger.com/sites/default/files/s3/release/current/Tutorials/zip/macro_pka.zip
  4. After downloading the zip file, unzip the contents in your Working Directorythe location where files are saved for ease of access throughout the tutorial

Figure 2-2. Save Project panel.

  1. Go to File > Save Project As
  2. Change the File name to macro_pka_tutorial, click Save
    • The project is now named macro_pka_tutorial.prj

3. Application of Macro-pKa to methylamine

In this section we will apply the Macro-pKa workflow to methylamine. This simple example could also be studied with the Jaguar pKa (micro-pKa) workflow, but we will use it to demonstrate the basic operation of the Macro-pKa workflow. Note that a more complex example is available in Section 4. 

Figure 3-1. Sketching, saving and naming methylamine.

There are several ways to create a 3-dimensional structure of methylamine in the workspacethe 3D display area in the center of the main window, where molecular structures are displayed, but perhaps the easiest way is through the 2D Sketcher panel.

  1. Go to Edit > 2D Sketcher…
  2. Sketch the structure of methylamine (CH3NH2), similarly to what is shown in the Figure
  3. Click Save as New
  4. Enter methylamine for Input Entry Title
  5. Click OK
  6. Close the 2D Sketcher panel

Figure 3-2. The Jaguar - pKa panel with the Macro-pKa radio button selected.

  1. With the methylamine entry selected(1) the atoms are chosen in the Workspace. These atoms are referred to as "the selection" or "the atom selection". Workspace operations are performed on the selected atoms. (2) The entry is chosen in the Entry List (and Project Table) and the row for the entry is highlighted. Project operations are performed on all selected entries in the entry lista simplified view of the Project Table that allows you to perform basic operations such as selection and inclusion and includedthe entry is represented in the Workspace, the circle in the In column is blue in the workspacethe 3D display area in the center of the main window, where molecular structures are displayed, go to Tasks > Materials > Quantum Mechanics > Molecular Quantum Mechanics > More Molecular QM Tasks > pKa
    • Alternatively, search the Tasks widget for pKa
    • The Jaguar - pKa panel opens

 

The Jaguar - pKa panel can be used to compute both micro-pKa’s (using the Jaguar pKa workflow) and macro-pKa (using the Macro-pKa workflow). In this tutorial we are interested in Macro-pKa calculations.

 

  1. Choose the Macro-pKa radio button

At this point it is necessary to learn how to define charge transitions in the panel.

Each macro-pKa transition (observed in experiment) is associated with a transition between two neighboring charges C and C-1 in the Macro-pKa workflow. The Set range for selector gives the user the possibility to choose between two ways which define macro-pKa transitions. Using the dropdown, one can set the pKa range or the Relative Charge range.

 

Let us first discuss defining macro-pKa transitions by setting the pKa range. The default pKa range is from 2 to 12, which should cover most pKa modeling needs, although formally Macro-pKa does not have restrictions on the range of pKa values it can calculate. For a given molecule with total charge C Macro-pKa attempts to set up tautomers of this molecule with charges C+1 and C-1, and then predict the macro-pKas for the transitions C+1 to C and C to C-1. If the macro-pKa corresponding to the C+1 to C transition is below the Min pKa limit, the generation of the tautomers in the direction of increasing charges stops. Analogously, if the macro-pKa corresponding to the C to C-1 surpasses the Max pKa limit, tautomers of charge C-2 will not be generated, and the macro-pKa prediction in this direction will stop. If either of the limits is not reached, the generation of the tautomers and the calculation of the corresponding macro-pKa transitions continues in the corresponding directions, by gradually protonating or deprotonating the molecule (i.e. increasing or decreasing the total charge), until the limit is reached. If the tautomers for a certain charge cannot be generated (for example, because the molecule cannot be protonated or deprotonated further), the search in this direction stops. All the macro-pKa’s computed in the process are reported. Note that the Min pKa and Max pKa values define the limits at which the generation of the tautomers must stop, but they are not necessarily the boundaries of the lowest and highest predicted pKa values.

 

A concrete example may help clarify the functioning of the algorithm. The user sets Min pKa and Max pKa to 2 and 12, respectively. The given molecule has total charge C = 1. Tautomeric forms with charge 2, 1, and 0 will be generated. The macro-pKa corresponding to charge transition 2 to 1 is predicted to be 1.692. This value is lower than Min pKa, and therefore no tautomers of charge 3 will be generated, and no further pKa transitions will be computed in the direction of the increasing charge. The macro-pKa corresponding to charge transition 1 to 0 is predicted to be 7.889. This is lower than the Max pKa limit, and therefore the search in the direction of the diminishing charge continues. Tautomers of charge -1 are generated and the macro-pKa transition of charge 0 to -1 is predicted to be 11.073. This is again lower than the Max pKa limit, so the search in the direction of diminishing charge must continue. When trying to generate tautomers with charge -2 it turns out that no further deprotonation is possible (as the molecule ran out of hydrogen atoms that can conceivably dissociate as protons). Therefore, the search in the direction of the decreasing charge terminates at this point. The following macro-pKa transitions will be reported:

 

Macro-pKa (2 to 1) = 1.692

Macro-pKa (1 to 0) = 7.889

Macro-pKa (0 to -1) = 11.073

 

The second way to specify macro-pKa transitions is via the relative charge range. Here the user sets the lowest and highest relative charges to be considered. Let us imagine that the total charge of the molecule is 2. If the minimum and maximum relative charges are set to -2 and 3, respectively, then the lowest and highest total charges that will be considered are 2 - 2 = 0 and 2 + 3 = 5. As a result, the following macro-pKa transitions will be calculated: macro-pKa (5 to 4), macro-pKa (4 to 3), macro-pKa (3 to 2), macro-pKa (2 to 1), and macro-pKa (1 to 0). If total charge 4, for example, cannot be reached as no reasonable protonatable atoms remain after total charge 3, the workflow will skip computation of macro-pKa (4 to 3) and macro-pKa (5 to 4).

 

Methylamine is a very simple molecule. It can obviously have only one pKa value associated with the equilibrium:

 

CH3NH3+ ⇄ CH3NH2 + H+

 

In application to methylamine, Macro-pKa can be used with either the pKa range or the relative charge range setting. Since we expect only one charge transition for this molecule (1 to 0), it is easier to launch the Macro-pKa calculation with relative charges. If one were to set pKa range (let us say, from 2 to 12), Macro-pKa would attempt to compute macro-pKa’s for (2 to 1), (1 to 0), and (0 to -1) transitions, according to the explanation given in the previous section. The (2 to 1) and (0 to -1) calculations are unnecessary and are not likely to succeed, and therefore the best way to launch Macro-pKa for methylamine (and other conceptually similar amines) is through the relative charge setting.

Figure 3-3. Defining charge transitions and conformational search settings.

  1. From the Set range for dropdown, choose Relative Charge
  2. Set the Min and Max to 0 and 1, respectively
  3. Ensure that Perform conformational searches on input structures is checked
  4. Set the Max number of conformers to use for each species to 2

 

Note: Methylamine and its protonation forms are not conformationally flexible, and conformational search here is probably superfluous. However, we recommend doing conformational search with at least 2-5 considered conformations per generated tautomer for larger and more flexible molecules.

Figure 3-4. Adjusting Job Settings and running the job.

  1. Maintain the Job name of jag_methylamine_pka
  2. Adjust the job settings () as needed
    • This job requires a CPU host. The job can be completed in about 30 minutes on a CPU host
    • The Macro-pKa calculation that we are planning to launch is not going to be CPU-intensive, and therefore it should complete quickly even if only 1-2 CPUs are used. For larger and more flexible molecules we recommend using a larger number of CPUs for faster calculations.
  3. Click Run
  4. Close the Jaguar - pKa panel
The files are also available for importing in the provided tutorial files: go to File > Import Structures, navigate to where you downloaded the tutorial files, and Open Section_03 > jag_methylamine_pKa > jag_methylamine_pka_methylamin_macro_pka > jag_methylamine_pka_methylamine_macro_pka-outmae

Figure 3-5. The entry list after incorporation with the first new entry included in the workspace.

When the methylamine calculation finishes the results are incorporated in the Project Tabledisplays the contents of a project and is also an interface for performing operations on selected entries, viewing properties, and organizing structures and data and the Workspacethe 3D display area in the center of the main window, where molecular structures are displayed.

All the micro-pKa results will be incorporated (see the group called “micro-pKa pairs”). For methylamine there is only one micro-pKa (which is equivalent to macro-pKa, as discussed in Section 1). The computed micro-pKa is given in the title of the subgroup (in our case it is 10.590). In addition, the dominant conformation for every protonation form (tautomer) is incorporated. See the group called “tautomer populations”. One can view the tautomer populations as molecular properties in the Project Tabledisplays the contents of a project and is also an interface for performing operations on selected entries, viewing properties, and organizing structures and data. The incorporated results enable the user to examine the 3D structures of the generated tautomers and their lowest energy conformations.

 

For this simple molecule Macro-pKa does not generate many structures and micro-pKa’s. But for more complex molecules, like the one that we will consider later in this tutorial, there might be many incorporated structures.

 

Take a moment to explore the output entries.

Figure 3-6. Accessing the 2D summary via the WAM.

The easiest way to understand the results of the Macro-pKa prediction is to display a 2D summary that is available through a WAM widget.

  1. Find the WAM widget at the top of the incorporated results, click on it, and then click on View pKa Results
    • Alternatively, one can locate an html file inside the directory jag_methylamine_pka_final_results_figure and open it directly in a browser

When you click on View pKa Results… a diagram like the one shown below appears in a browser window that is associated with your Maestro session. Here we are showing only a relevant portion of the summary:

 

The full report contains more information. The diagram includes the macro-pKa predicted for methylamine. It is 10.590, which is very close to the experimental value 10.64 [M. B. Fernandez et al., J. Heter. Chem., 17, 667-672 (1980)].

 

The diagram contains the total charges and 2D structures of the associated tautomers. The populations of the tautomers in percent are shown below the 2D structures and their indices.

 

It is important to distinguish between the pH-dependent (black) and pH-independent (green) populations. pH-dependent populations add up to 100% for all the populations across all charges. These populations are called pH-dependent because, being spread across all charges, they depend on pH at which those populations are calculated. The pH-dependent populations add up to 100% at any pH but the distribution of tautomer populations of charges 0 and 1 will be different at different pH. By default, pH is set to 7.000 (the pH value is reported in the diagram). You can set a different pH using the Update pKa Settings panel and regenerate the diagram that is applicable to that pH (see Figure 3-6). The summary diagram also contains a graphical representation of pH-dependent populations as a function of pH. This makes it convenient to estimate the pH-dependent populations at any other pH without having to regenerate them.

 

The regeneration of the diagram at a different pH should take just a few seconds. No computationally expensive calculations need to be conducted in order to recalculate pH-dependent populations at a different pH. pH-independent populations are computed per charge, and they are constant at any pH. pH-dependent populations are useful if you are interested in the most populated tautomers at a given pH, whereas pH-independent populations are useful if you are interested in a relative stability of tautomers associated with a particular total charge.

4. Application of Macro-pKa to guanine

In this section we will apply the Macro-pKa workflow to guanine. This molecule is a more complex example in which the Macro-pKa workflow is essential. 

Figure 4-1. Creating a new entry.

  1. Once you are finished viewing the results for methylamine, right-click on any empty space in the workspacethe 3D display area in the center of the main window, where molecular structures are displayed and choose Create New Entry > Create Empty Entry
    • A new empty entry is added to the entry lista simplified view of the Project Table that allows you to perform basic operations such as selection and inclusion
    • A dropdown appears prompting to name the entry
  2. For Entry Title, input guanine

Figure 4-2. Drawing and saving the guanine sketch.

Again, there are several ways to create a 3-dimensional structure of guanine in the workspacethe 3D display area in the center of the main window, where molecular structures are displayed, but perhaps the easiest way is through the 2D Sketcher panel.

  1. With the new, empty guanine entry selected(1) the atoms are chosen in the Workspace. These atoms are referred to as "the selection" or "the atom selection". Workspace operations are performed on the selected atoms. (2) The entry is chosen in the Entry List (and Project Table) and the row for the entry is highlighted. Project operations are performed on all selected entries and includedthe entry is represented in the Workspace, the circle in the In column is blue, go to Edit > 2D Sketcher…
  2. Sketch the structure of guanine, similarly to what is shown in the Figure
  3. Click Update Entry
  4. Close the 2D Sketcher panel

 

Note: If preferred, the canonical SMILES string for guanine is n1c(N)[nH]c(=O)c(c12)[nH]cn2, which can be pasted into the sketcher for immediate 2D conversion via the button

Figure 4-3. The Jaguar - pKa panel with the Macro-pKa radio button selected.

  1. With the guanine entry selected(1) the atoms are chosen in the Workspace. These atoms are referred to as "the selection" or "the atom selection". Workspace operations are performed on the selected atoms. (2) The entry is chosen in the Entry List (and Project Table) and the row for the entry is highlighted. Project operations are performed on all selected entries in the entry lista simplified view of the Project Table that allows you to perform basic operations such as selection and inclusion and includedthe entry is represented in the Workspace, the circle in the In column is blue in the workspacethe 3D display area in the center of the main window, where molecular structures are displayed, go to Tasks > Materials > Quantum Mechanics > Molecular Quantum Mechanics > More Molecular QM Tasks > pKa
    • Alternatively, search the Tasks widget for pKa
    • The Jaguar - pKa panel opens
  2. Choose the Macro-pKa radio button

Guanine, in comparison with methylamine, is a more appropriate target for the Macro-pKa workflow. The macro-pKa of methylamine could be predicted with Jaguar pKa because the micro-pKa of methylamine coincides with its macro-pKa. Not so for guanine, which (at least conceivably) can exist in multiple tautomers in most of its charge forms. For guanine, there are multiple micro-pKa’s associated with every charge transition. These micro-pKa’s values need to be computed and combined using formula (6) from Section 1. These operations, difficult to accomplish “manually”, are automatically performed by the Macro-pKa workflow.

 

Some of the tautomeric complexity of guanine is illustrated below:

For complex tautomeric equilibria, such as those shown for guanine, even experienced scientists might not always be able to guess the most stable tautomeric form on the basis of “chemical intuition” or experience. However, for our calculations it is not important to know which tautomer is the most stable, since the results of Macro-pKa do not depend (up to numerical noise) on the starting tautomer molecule.

 

In Section 3, when selecting settings for the Macro-pKa calculation on methylamine, we opted for the relative charge setting. For guanine we will choose the pKa range setting. Unless the user is an expert in pKa of the heterocycle which makes up the core of the guanine molecule, it might be difficult to guess which macro-pKa of guanine will be associated with a particular charge transition. For example, if the user attempts to calculate the pKa value closest to 7, it is not immediately clear if the user must select (1 to 0), (0 to -1), or (-1 to -2) charge transition. When dealing with relatively unfamiliar types of molecules, it might be beneficial to compute the pKa in a broad titration range (such as 2 to 12) to better understand the protonation equilibria characterizing this type of molecule, and avoid missing the pKa closest to the pH of interest.

Figure 4-4. Defining charge transitions and conformational search settings.

  1. From the Set range for dropdown, choose pKa
  2. Set the Min and Max to 2.0 and 12.0, respectively
  3. Ensure that Perform conformational searches on input structures is checked
  4. Set the Max number of conformers to use for each species to 2

Note: A larger molecule typically requires more conformations to consider, although in practice one rarely goes above 5 conformations to keep computational cost within reasonable limits.

Figure 4-5. Naming and running the job.

  1. Change the Job name to jag_guanine_pka

If you wish to run the job, adjust the job settings () as needed. This job takes many hours to complete, so output files are provided.

  1. Close the Jaguar - pKa panel

Figure 4-6. The entry list after importing.

  1. Go to File > Import Structures
  2. Navigate in the provided tutorial files to the jag_guanine_pka directory and Open the Section_04 > jag_guanine_pka > jag_guanine_pka_guanine_macro_pka > jag_guanine_pka_guanine_macro_pka-out.mae file.
    • A new entry group containing 225 entries is incorporated in the Project Tabledisplays the contents of a project and is also an interface for performing operations on selected entries, viewing properties, and organizing structures and data and Workspacethe 3D display area in the center of the main window, where molecular structures are displayed

As in case of methylamine, the incorporated results contain all the micro-pKa’s computed by the Macro-pKa workflow. Each micro-pKa is represented by the most stable conformations of the corresponding tautomers. There are also the most stable conformations of each tautomer in a separate group. The associated pH-dependent and pH-independent populations are available through the Project Tabledisplays the contents of a project and is also an interface for performing operations on selected entries, viewing properties, and organizing structures and data.

 

Even though the information associated with 3D structures should be useful, it is very detailed and might be difficult to analyze if a quick summary of the results is desirable. The easiest way to view the summary of the results in a single image is again through the WAM widget, as we did in case of methylamine. Alternatively, you can locate the .html file inside the jag_guanine_pka_guanine_macro_pka directory in the provided files and open it directly in a browser. The relevant portions of the diagram that will be opened in a browser are represented below.

 

Let us first review the predicted macro-pKa values. There are three transitions:

 

Macro-pKa (1 to 0) = 2.98

Macro-pKa (0 to -1) = 9.57

Macro-pKa (-1 to -2) = 11.97

 

 

When we compare these predicted pKa values with the reported experimental pKa values [see References],

 

pKa1 = 3.0

pKa2 = 9.3

pKa3 = 12.6

we find good agreement for pKa1, pKa2 , and pKa3 (error within 1.0 pH units). Next, from the pH-dependent numerical values in the above Figure, and from a plot in the below Figure, we see that at pH 7 the workflow predicts the most populated tautomers to be 4, 11, and 18, with 12 only slightly less populated than 11.  For this analysis we were only looking at tautomers with pH values between 2 and 12 as that was what was originally asked for, even though the panel gives a wider range of values by design. The pH-independent values (in green in the previous image) show that the tautomer from which we started our calculation (form 11) is indeed the most stable tautomer.

 

 

The situation around the relative stability of tautomers with total charge +2 (tautomers 1 and 2) illustrates why it is important to report pH-independent populations in addition to pH-dependent ones. If pH-independent populations were not reported, it would be actually impossible to know the relative stability of tautomers 1 and 2 from pH-dependent populations alone, as they are both reported as 0.00%. In contrast, the reported pH-independent populations of tautomers 1 and 2 which are 60.33%, 39.59%, respectively, give a clear picture of the relative stabilities of these tautomers.

5. Conclusion and References

This tutorial provides an introduction to practical DFT-based prediction of pKa of organic molecules. It goes over the differences between micro-pKa and macro-pKa values and explains which of the two Schrödinger workflows, Jaguar pKa or Macro-pKa, is applicable to typical cases of pKa prediction. Then, the tutorial focuses on the Macro-pKa workflow in application to two molecules. The first molecule, methylamine, is conceptually very simple with regards to its protonation equilibria, and both Jaguar pKa and Macro-pKa workflows can be used to predict its pKa. The pKa prediction of the second molecule, guanine, would be practically very difficult and error-prone if one used the Jaguar pKa workflow. This is because guanine possesses numerous tautomers and its protonation equilibria are characterized by multiple macro-pKa values. The tutorial shows that an application of the Macro-pKa workflow to guanine automates the process of generating the molecule’s tautomers, as well as the process of computing its micro- and macro-pKa’s. The macro-pKa’s predicted by the Macro-pKa workflow for guanine are in excellent or in good agreement with experimentally measured values.

For further learning:

For introductory content, focused on navigating the Schrödinger Materials Science interface, an Introduction to Materials Science Maestro tutorial is available. Please visit the materials science training website for access to 70+ tutorials. For scientific inquiries or technical troubleshooting, submit a ticket to our Technical Support Scientists at help@schrodinger.com. 

For self-paced, asynchronous, online courses in Materials Science modeling, including access to Schrödinger software, please visit the Schrödinger Online Learning portal on our website.

For some related practice, proceed to explore other relevant tutorials:

For further reading:
  • Visit the panel help documentation here.
  • Quantum chemical pKa prediction for complex organic molecules. DOI:10.1002/qua.25561
  • Multiconformation, Density Functional Theory-Based pKa Prediction in Application to Large, Flexible Organic Molecules with Diverse Functional Groups. DOI:10.1021/acs.jctc.6b00805
  • The Synthesis and Determination of Acidic Ionization Constants of Certain 5-Substituted 2-Aminopyrrolo[2,3-d]pyrimidin-4-ones and Methylated Analogs. https://doi.org/10.1002/jhet.5570330341

6. Glossary of Terms

Entry List - a simplified view of the Project Table that allows you to perform basic operations such as selection and inclusion

Included - the entry is represented in the Workspace, the circle in the In column is blue

Project Table - displays the contents of a project and is also an interface for performing operations on selected entries, viewing properties, and organizing structures and data

Scratch Project - a temporary project in which work is not saved. Closing a scratch project removes all current work and begins a new scratch project

Selected - (1) the atoms are chosen in the Workspace. These atoms are referred to as "the selection" or "the atom selection". Workspace operations are performed on the selected atoms. (2) The entry is chosen in the Entry List (and Project Table) and the row for the entry is highlighted. Project operations are performed on all selected entries

Working Directory - the location where files are saved

Workspace - the 3D display area in the center of the main window, where molecular structures are displayed