Marylee Z. Southard

Design of novel pharmaceutical products via combinatorial optimization

Abstract This paper describes a new methodology for the application of computer-aided molecular design to the design of pharmaceutical products with prespecified physical properties. For a pharmaceutical product to be effective, it must not only cause a desirable therapeutic effect, but it must also have proper values of physical properties such as solubility and density. In order to apply computer-aided molecular design to the discovery of new pharmaceuticals, it is therefore necessary to be able to first predict the physical and biological properties of a given molecule, and then optimize over an entire set of molecules to find one which matches target values on those properties. This work employs topological descriptors to predict the physical properties of pharmaceutical molecules. The descriptors used here, called connectivity indices, are easy to compute, yet contain valuable information about the internal molecular structure of a molecule. Property prediction, via connectivity indices, can be viewed as an improvement over group contribution methods, since these indices take into account molecular connectivity and internal electronic structure in addition to the identity of each group in the molecule. Thus these indices correlate well with physical properties which are important in pharmaceutical design. Furthermore, these indices are fairly simple to compute, and a proper choice of variables to describe the molecule allows the equations for these indices to be written in a linear form. The optimization problem used here combines a set of basic groups, which are defined as a non-hydrogen atom at a given valency state bonded to a given number of hydrogens, to form a candidate molecule. Each candidate molecule can then be tested by computing estimated property values and comparing those to prespecified target values. Structural constraints are also added to the problem to ensure a connected, stable molecule is generated. The set of constraints and the correlation equations for property prediction are then combined and reformulated, resulting in a mixed integer linear program (MILP). If a certain functional group is known to be required in the molecule, this requirement can also be added to the constraint set. For problems consisting of a smaller number of basic groups or property targets, commercial solvers can be employed to solve the resulting MILP. The effectiveness of this method is presented through the solution of a small example problem.

Pharmaceutical product design using combinatorial optimization

Abstract A two-step computational method for designing new molecules in medicinal chemistry is described. In the first step, topological indices are used to develop structure-based correlations for properties of interest. Zeroth and first order connectivity indices are employed to develop linear correlations for three physical properties of interest in pharmaceutical chemistry: octanol–water partition coefficient (OWPC), melting point and water solubility. These correlations are then used within an optimization framework to design molecules having the desired properties. This step involves formulating a mixed integer linear program (MILP) which includes the property correlations, structural constraints which ensure that a stable, connected molecule is formed, and an objective function which minimizes the deviation from a set of property targets. A new data structure, known as a partitioned adjacency matrix, is employed to allow the connectivity index definitions to be written linearly, such that they can be included in an MILP and solved using a standard branch-and-bound method. The connectivity of the molecule is ensured by the inclusion of network flow constraints within the formulation. Three examples show the efficacy of this approach.

Dissolution of lonizable Drugs into Unbuffered Solution: A Comprehensive Model for Mass Transport and Reaction in the Rotating Disk Geometry

A model has been developed to describe the mass transport and reaction of ionizable compounds where mass transfer is caused by convection and diffusion from a rotating disk. Dissolution rates of benzoic acid, 2-naphthoic acid, and indomethacin in aqueous solutions of high ionic strength (I = 0.5 with potassium chloride) at 25°C were investigated. The model includes the effects of diffusion, convection, and simultaneous acid/base reaction at all points in the region adjacent to the dissolving solid. The solution of the transport equations is obtained numerically with an iterative algorithm which uses (a) closure of all material balances and (b) equilibria at the solid/liquid surface as constraints. The model solution yields both the flux of the dissolving acid and the concentration profile of each component. Reduced values of all reaction rate constants are used in the region adjacent to the dissolving surface to allow convergence of the computation. Although nonequilibrium concentration values are calculated, it is shown that the theoretical dissolution rate determined as the solution of the model is insensitive to the magnitude of the rate constants as their maximum useable values are approached. Comparisons of the model results with experimentally determined fluxes show close agreement and confirm that the transport mechanisms in the model formulation are consistent with the measured values. Further, the inclusion of convection allows accurate calculations without utilization of an arbitrary boundary layer thickness. Accurate dissolution rates can be determined using this technique under a wide range of conditions, except at low pH.

Experimental Determinations of Diffusion Coefficients in Dilute Aqueous Solution Using the Method of Hydrodynamic Stability

Diffusion coefficients were experimentally determined in dilute aqueous solution at 25 ± 0.1°C, ionic strength 0.5 M, using Taylor’s method of hydrodynamic stability. The methodology described is accurate enough to show significant differences in diffusion coefficients between the various ionic forms of the same species as a function of degree of ionization. In Taylor’s method, diffusion coefficients were measured by allowing two solutions of differing solute concentration to contact in a capillary tube, forming a stable, measurable concentration gradient. The solute diffusion coefficient is a function of the gradient, the solution viscosity, the solution density, and some capillary dimensions. Viscosity was maintained constant across experiments and values of sufficient accuracy were available in the literature. Solution densities were measured with a tuning fork densimeter. Compounds studied were o-aminobenzoic acid, benzoate anion, the four forms of phosphate and citrate, and the zwitterionic forms of glycine, diglycine, and triglycine. Based on the results for the four forms of phosphate and citrate, experimental diffusivity values vary with the ionic state of the diffusant, presumably because of the altered state of hydration as charge varies. For the glycine series, the diffusivity showed an unexpected dependency on molecular weight (size).

A numerical convective-diffusion model for dissolution of neutral compounds under laminar flow conditions

The convective-diffusion model reported earlier by Shah and Nelson (J. Pharm. Sci., 64 (1975) 1518–1520) described the dissolution of neutral compounds under laminar flow conditions. A linear velocity profile over the tablet surface was assumed to allow calculation of an analytical solution. In the present study, dissolution under laminar flow conditions was modified to include the actual parabolic velocity profile. The modified model was solved numerically using finite difference approximations. The numerical solution to the model predicted dissolution rates to be within 10% of experimental values for benzocaine, a neutral compound, and for benzoic acid and naproxen at suppressed ionization conditions. The modified model predictions were slightly better at high flow rates (100 mlmin) when compared to the earlier model. Overall, however, the Shah-Nelson simplified assumption was found to predict the experimental results as closely as the numerical solution. Using the same apparatus as described earlier by Shah and Nelson (1975), the flow conditions in the device were characterized over the entire range of flow rates considered (1.10–110 mlmin). Flow was found to be in the laminar region from Reynolds number calculations and to be fully developed, devoid of entrance effects, before reaching the tablet surface. As expected, the concentration boundary layer thickness grew with distance from the leading edge of the tablet, but stayed within the dimensions of the flow chamber even at the lowest flow rates employed (1.10 mlmin).

A Convective-Diffusion Model for Dissolution of Two Non-interacting Drug Mixtures from Co-compressed Slabs Under Laminar Hydrodynamic Conditions

A numerical convective-diffusion dissolution model has been extended to describe dissolution of two neutral non-interacting drugs co-compressed in a slab geometry. The model predicted the experimental dissolution rates of naproxen/phenytoin mixtures and hydrocortisone/nitrofurantoin mixtures quite accurately, except for phenytoin in the naproxen/phenytoin mixture at low weight proportions. A non-linear dependence of dissolution rate on weight proportion with a positive deviation from linearity was observed. An increase in flow rate increased the dissolution rate and the cube-root dependency of dissolution rate on the flow rate for a given weight proportion of the component in the slab, as proposed earlier by Shah and Nelson for pure compounds, was also observed here, suggesting that the changes in dissolution profile were caused by changes in surface area only. As expected from the model an increase in particle size of the powders used to make the slab decreased the dissolution rate. This was explained by an increase in the average length of the component resulting in a bigger ‘carryover’ of material from one section of the component in the slab to the next section of the same component, due to convection, and hence lower flux.

Dependence of Dissolution Rate on Surface Area: Is a Simple Linear Relationship Valid for Co-Compressed Drug Mixtures?

A quantitative analysis of the dependence of dissolution rate on the relative surface area occupied by two non-interacting drug mixtures from co-compressed slabs is described. The results from the experimental dissolution rates of each component from naproxen/ phenytoin co-compressed slabs under laminar flow conditions, when corrected for the area occupied by that component in the slab, contradict the stagnant layer model predictions, where dissolution rates are assumed to be directly proportional to the occupied surface area. Simulations from non-mixed co-compressates of naproxen and phenytoin indicated that dissolution rates are proportional to bL2/3, as reported for pure compounds in the laminar dissolution apparatus by Shah and Nelson. However, for a well mixed co-compressate, which differs with the non-mixed case only in the distribution of particles, this proportionality did not hold. The deviation was explained by ‘carryover’ of material from one section of the component to the next due to fluid flow, resulting in an increase in apparent effective length of the component in the slab (Leff).

Dissolution of weak acids under laminar flow and rotating disk hydrodynamic conditions: Application of a comprehensive convection–diffusion–migration–reaction transport model

A steady-state mass transfer model that incorporates convection, diffusion, ionic migration, and ionization reaction processes was extended to describe the dissolution of weak acids under laminar flow and a rotating disk hydrodynamics. The model accurately predicted the experimental dissolution rates of benzoic acid, 2-naphthoic acid, and naproxen in unbuffered and monoprotic buffers within the physiological pH range for both hydrodynamic systems. Simulations at various flow rates indicated a cube root dependency of dissolution rate on the flow rate for a given bulk pH value for the laminar hydrodynamic system, as proposed earlier by Shah and Nelson (1975. J Pharm Sci 64(9):1518-1520) for neutral compounds. The model has limitations in its ability to accurately predict the dissolution of weak acids under certain conditions that imposed steep concentration gradients, such as high pH values, and for polyprotic buffer systems that caused the numerical solution to be unstable, suggesting that alternative numerical techniques may be required to obtain a stable numerical solution at all conditions. The model presents many advantages, most notably the ability to successfully predict the complex process under physiological conditions without simplifying assumptions, and therefore accurately representing the system in a comprehensive manner.