Mathematical model for melting curves and heteroduplex content

In Figure 1, we show high-resolution melting curves of the DNA amplicon from three genotypes of a SNP for which the homozygous mutation is indistinguishable from that of the wild type, due to nearest-neighbor thermodynamic symmetry of the mutation. (This says that the bases immediately surrounding the mutation are identical when the strands are interchanged, e.g., TCA/AGT ↔ TGA/ACT.) Our goal is to add the right fraction of wild-type DNA to each genotype before PCR, so that after amplification, melting, and reannealing, the mixture with the homozygous mutant sample will develop heteroduplex content but the mixture with the wild-type sample will not, making the curves distinguishable. The mixture with the heterozygous sample will have its heteroduplex content reduced from its natural 50 % value, so we must be careful that its curve remains distinct from the heteroduplex-enhanced homozygous sample.

Therefore, in this section, we develop a model for the melting curve of a mixture of genotypes in terms of the melting curves of the consituent duplex types and the mixture fractions. In the case above when the homozygous mutant and wild-type curves are indistinguishable, we show that the difference among curves of different genotypes depends solely on the heteroduplex content of the mixture. In addition, the temperature of maximum separation does not change with heteroduplex content. This is in spite of the maybe surprising fact that the two heteroduplex species have different thermodynamic melting behavior, even with respect to nearest-neighbor approximation.

Finally, we model and optimize the separation of heteroduplex contents of the three genotypes in terms of mixture fraction.

We begin by listing what kind and fraction of duplexes are present when we mix wildtype and unknown DNA in fraction x of wild-type to 1−x of either wild-type, homozygous mutant, or heterozygous mutant. The DNA mixtures are amplified by PCR, dissociated, then rapidly cooled so that complementary and nearly complementary strands form duplexes independently, i.e., solely in fraction to their fractions. When x = 0, the mixture consists of just the original genotype. When x = 1, it describes a pure wild-type sample regardless of the initial unknown genotype. Table 1a lists the resulting duplex fractions.

After denaturation and cooling, fractions of four duplexes are obtained by multiplying the corresponding strand fractions from Table 1a.

Corresponding to each duplex species d_{j} will be its normalized fluorescence curve, F_{j}(T). The function F_{j}(T) refers to the curve with background fluorescence removed and accounts for factors such as lower fluorescence per duplex behavior of heteroduplexes. Assuming fluorescence is additive, the normalized melting curves for each of the genotypes are described by

W (T) = 1F_{1}(T),
M(T, x) = x^{2}F_{1}(T) + x(1 − x)F_{2}(T) + x(1 − x)F_{3}(T) + (1 − x)^{2}F_{4}(T),

H(T, x) = \frac{(1 + x)^{2}}{4} F_{1}(T) + \frac{1 − x^{2}}{4} F_{2}(T) + \frac{1 − x^{2}}{4} F_{3}(T) + \frac{(1 − x)^{2}}{4} F_{4}(T).

Our goal is to maximize our ability to distinguish these three melting curves, as measured by the minimum separation between any two of them, defined as the maximum absolute value of their difference:

max_{x∈[0,1]} min{max_{T} |M(T, x) − W (T, x)|, max_{T} |H(T, x) − W (T, x)|, max_{T} |M(T, x) − H(T, x)|}.

For this reason we consider

W (T, x) − M(T, x) = (1 − x^{2})F_{1}(T) − x(1 − x)F_{2}(T) − x(1 − x)F_{3}(T) − (1 − x)^{2}F_{4}(T)

and

W (T, x) − H(T, x) = (1 − \frac{(1 + x)^{2}}{4} )F_{1}(T) − \frac{1 − x^{2}}{4} F_{2}(T) − \frac{1 − x^{2}}{4} F_{3}(T) − \frac{(1 − x)^{2}}{4} F_{4}(T).

We will determine H(T, x) − M(T, x) as the difference of the differences above, after simplifying.

We are most interested in the situation where F_{1}(T) and F_{4}(T) are difficult to distinguish, or we would not be adding DNA of one type to do so. Therefore we assume F_{1}(T) = F_{4}(T), as predicted by the nearest-neighbor approximation for SNPs with the symmetry described in the previous section. This symmetry does not make F_{2} and F_{3} agree, but this is not required for our derivation. Using the subscript _{=} to indicate melting curves referring to this case, we replace F_{4} by F_{1} and combine terms to get

W_{=}(T, x) − M_{=}(T, x) = (2x − 2x^{2})F_{1}(T) − x(1 − x)F_{2}(T) − x(1 − x)F_{3}(T)

and

W_{=}(T, x) − H_{=}(T, x) = (\frac{1 − x^{2}}{2} F_{1}(T) − \frac{1 − x^{2}}{4} F_{2}(T) − \frac{1 − x^{2}}{4} F_{3}(T).

Finally, noticing that the coefficient of F_{1} in both expressions is minus twice the common coefficients of F_{2} and F_{3}, we divide it equally between F_{1} and F_{4} again to arrive at

M_{=}(T, x) − W_{=}(T, x) = m(x)F (T)

H_{=}(T, x) − W_{=}(T, x) = h(x)F (T)

and therefore

H_{=}(T, x) − M_{=}(T, x) = (h(x) − m(x))F (T),

where

F (T) = \frac{F_{1}(T) + F_{4}(T)}{2} − \frac{F_{2}(T) + F_{3}(T)}{2}

is the difference between the average homoduplex melting curve and the average heteroduplex melting curve,

m(x) = 2x(1 − x)

h(x) = \frac{1 − x^{2}}{2}

h(x) − m(x) = \frac{3x^{2} − 4x + 1}{2}.

Since there are two types of heteroduplex, m(x) and h(x) represent the total heteroduplex content of the mutant homozygous and heterozygous plus wild-type mixtures, respectively, and along with h(x) − m(x), the difference in heteroduplex content between any two mixtures. The graphs of the absolute values of the three heteroduplex content difference functions are given in Figure 2a, annotated with key features derived below. The graphs of the theoretical duplex melting curves F_{j}(T) and the corresponding F (T) are given in Figure 2b.

These expressions have two important consequences. First, they show that in the case of interest, the pointwise separation of the melting curves of different genotypes is proportional solely to the difference in heteroduplex content of the mixtures. Second, they uncouple the x (fraction of added DNA) and T (temperature) dependence of the differences among fluorescence curves of different genotypes. This allows us to optimize the fraction added independently, regardless of the specific nature of individual duplex curves contributing to the superpositions. This is perhaps surprising, since in general, different combinations of different functions result in varying shapes and locations of their extrema. Since all expressions are linear, we can obtain the difference of the negative derivative curves of the different genotypes simply by replacing F_{j} in the expressions above with their negative derivatives, and all of the results we will obtain regarding maxima and minima will still hold.

Symbolically, our original genotype separation problem reduces to

max_{x∈[0,1]} min{m(x), h(x), |h(x) − m(x)|}F,

where

F = max_{T} F (T)

represents the twice the separation of the original (no added wild-type DNA) heterozygous curve from the common original wild-type and homozygous SNP curves, since h(0) = \frac{1}{2}.

What remains is to optimize the fraction of added DNA x to maximize

min{m(x), h(x), |h(x) − m(x)|}.

In Figure 2a, this corresponds to finding the highest point on the lowest of the graphs of m, h, and |m − h|. Calculus confirms our visual intuition that the maximum of of the lowest graph can only occur where the slope of its tangent is zero, or it does not have a well-defined tangent. To find these points, we first decide which of the functions has the minimum value for each value of x. To discover which is the smaller of h(x) and m(x) as well as to simplify |h(x) − m(x)|, we locate the points where m(x) = h(x) or h(x) − m(x) = 0. Using h(x) − m(x) = \frac{3x^{2}−4x+1}{2} = \frac{1}{2}(3x − 1)(x − 1), we find that for x < \frac{1}{3}, m(x) < h(x) and |h(x) − m(x)| = h(x) − m(x). Further comparing m(x) and h(x) − m(x), we find that they are equal in this interval when h(x) − 2m(x) = \frac{7}{2}x^{2} − 4x + \frac{1}{2} = \frac{1}{2}(7x − 1)(x − 1) = 0, or x = \frac{1}{7}. So for 0 < x < \frac{1}{7}, m(x) < h(x)−m(x) < h(x), and for \frac{1}{7} < x < \frac{1}{3} h(x)−m(x) < m(x) < h(x). For \frac{1}{3} < x < 1, the order becomes h(x) − m(x) < h(x) < m(x). The only place the minimum of these functions has a horizontal tangent is halfway between its roots \frac{1}{3} and 1, i.e., at x = \frac{2}{3}. This corresponds to adding twice as much wild-type DNA as there was unknown DNA and gives a separation of \frac{1}{6}F, between the heterogygous and homozygous SNP curves, or \frac{1}{3} of the original separation between the heterozygous melting curve and the other two. For this value, the separation between the wild-type melting curve and the other two melting curves will be even greater.

The only places the graph of the minimum of the three functions does not have a tangent are where it changes form, i.e., at x = \frac{1}{7} and x = \frac{1}{3}. Comparing the values at these points, we find the optimal fraction of added DNA occurs at x = x_{∗} = \frac{1}{7}, as indicated in Figure 2a. Here, the melting curve of the mixture with homozygous DNA is halfway between the other curves at the temperature of maximum separation, where its separation from each is \frac{12}{49}F. This is only barely less than half of the original separation of \frac{1}{2}F = \frac{24}{48}F between the heterozygous melting curve and the other two genotype curves.

A simple heuristic explanation for this value, corresponding to adding one part wildtype DNA to six parts unknown sample is based upon the observation that the melting curves will be optimally separated when the homozygous mutant curve becomes equidistant from both the wild-type and heterozygous melting curves, so the heteroduplex content of the wild-type-heterozygous mixture must be exactly twice that of the wild-type-homozygous mixture: h(x) = 2m(x) The ratio of 1 part wild-type to 6 parts unknown is optimal because 6 is the unique number which can be divided in equal parts (3 + 3, the heterozygous sample strands), and when one of the parts (3, SNP strands) is multiplied by the other plus 1 (4 = 3 + 1, wild-type sample plus added wild-type strands) one obtains exactly twice the original number (6, the homozygous SNP strands) multiplied by 1 (added wild-type strands.) At the simplest level, it is because (\frac{6}{2})(\frac{6}{2} + 1) = 2(6)(1) that x_{∗} = \frac{1}{7} of the total unknown plus added DNA is the optimal fraction of wild-type DNA to add. This is visualized in the animation

http://www.math.utah.edu/∼palais/pcr/michael/spike/spike.html

When the wild-type and homozygous SNP curves are not the same, the max-min problem does not separate into individual fraction added and temperature problems. In this case, the optimal added fraction still may be characterized by the two-dimensional generalization of the above criteria, but it must be found numerically depending upon the specific curves F_{j}(T). That is, the fraction-temperature rectangle [0, 1] × [T_{1}, T_{2}] must be divided into regions on which the corresponding h(x, T) − m(x, t) is either positive or negative, and within these regions, local extrema are characterized by when the gradient of the smaller of h(x, T) and m(x, T), denoted s(x, T) or of the larger minus the smaller, denoted l(x, T) − s(x, T) is zero, or when s(x, T) = l(x, T) − s(x, T). It seems reasonable to expect that the best separation could be obtained by including the possibility of adding whichever type of homozygous SNP or wild-type curve is already closer to the heterozygous curve.