Constructing “Proper” ROCs from Ordinal Response Data Using Weighted Power Functions

Abstract

Background. Receiver operating characteristic (ROC) analysis is the standard method for describing the accuracy of diagnostic systems where the decision task involves distinguishing between 2 mutually exclusive possibilities. The popular binormal curve-fitting model usually produces ROCs that are improper in that they do not have the ever-decreasing slope required by signal detection theory. Not infrequently, binormal ROCs have visible hooks that falsely imply worse-than-chance diagnostic differentiation where the curve lies below the no-information diagonal. In this article, we present and evaluate a 2-parameter, weighted power function (WPF) model that always results in a proper ROC curve with a positive, monotonically decreasing slope. Methods. We used a computer simulation study to compare results from binormal and WPF models. Results. The WPF model produces ROC curves that are less biased and closer to the true values than are curves obtained using the binormal model. The better performance of the WPF model follows from its design constraint as a necessarily proper ROC. Conclusions. The WPF model fits a broader variety of data sets than previously published power function models while maintaining straightforward relationships among the original decision variable, specific operating points, ROC curve contours, and model parameters. Compared with other proper ROC models, the WPF model is distinctive in its simplicity, and it avoids the flaws of the conventional binormal ROC model.

Keywords

receiver operating characteristic ROC analysis proper ROC weighted power function

Receiver operating characteristic (ROC) analysis is the standard method for evaluating and describing the accuracy of diagnostic systems with continuous test results or ordinal category data.^1–3 In contexts that involve distinguishing between 2 mutually exclusive possibilities (e.g., diseased v. nondiseased), the ROC graph describes the performance of the diagnostic system as the decision threshold is moved across all possible test values or decision points.⁴

By custom, ROC graphs depict the true positive rate (tpr = sensitivity) as a function of the false positive rate (fpr = 1 − specificity). Using line segments, one can draw a piece-wise ROC graph by connecting the (fpr, tpr) pairs that correspond to divisions between ordinal categories (e.g., 5-category diagnostic classes ranging from 1 = “strongly indicative of presence” to 5 = “strongly indicative of absence”) or from the empirical (fpr, tpr) pairs that correspond to points along a continuous response scale (e.g., reported values of a blood level).⁵ When, however, one believes that empirical data arise from some underlying, continuous response variable, one may wish to describe this belief by fitting a continuous, smooth ROC curve to the data. Beyond incorporating the just-mentioned prior belief into one’s interpretation of the data, basing inferences about accuracy on the parameters of a smoothed curve may help one to avoid inferential errors based on locations of particular (fpr, tpr) pairs that merely represent random numerical fluctuations in the data.

Of the many suggested curve-fitting schemes,^6,7 the most commonly used is the so-called binormal model, which assumes that on some monotonic transformation of the original decision variable, the data from the signal-plus-noise population (in medical applications, the “diseased” subjects) and noise-only populations (the “nondiseased” subjects) will conform to overlapping Gaussian distributions with different means and variances. The binormal model converts a linear-space ROC curve to a straight line on normal-deviate axes, so that

z (tpr) = A + Bz (fpr),

where A and B are the intercept and slope of the ROC line, and z(fpr) and z(tpr) are z-transforms of fpr and tpr.⁸ Note that $\frac{A}{B}$ is the difference between the means of the diseased (D+) and nondiseased (D−) populations on the transformed scale (measured in units of σ _D ₋, the standard deviation σ of the D− population), and B equals the ratio $\frac{σ_{D -}}{σ_{D +}}$ on the transformed scale. An equivalent statement of Equation 1a is

probit (tpr) = Φ^{- 1} (tpr) = A + B Φ^{- 1} (fpr),

where Φ(•) is the cumulative normal distribution function and Φ⁻¹(•) is its inverse. The parameters A and B are related to AUC_z, the total area under the binormal ROC curve, as follows:

A U C_{z} = Φ (\frac{A}{\sqrt{1 + B^{2}}}) .

The binormal ROC model provides a flexible, intuitively appealing, 2-parameter summary of diagnostic accuracy, and it accommodates data from a wide variety of practical contexts even when empirical data distributions deviate from normality.⁹ But the binormal model has one major drawback: Except when B = 1 (i.e., unless the ROC curve is symmetric with respect to the negative diagonal in the ROC square), the binormal model produces ROCs that at some point change from convex to concave. Binormal ROCs with B ≠ 1 also have a point located at z(fpr) = r = $\frac{A}{1 - B}$ where the curve falls below the no-information diagonal (i.e., the line defined by tpr = fpr) before tailing back toward a corner of the ROC square.

Often, the portion of the curve that falls below the no-information diagonal lies so close to a corner of the ROC square that it is imperceptible. Sometimes, however, curve-fitting maximum likelihood estimation (MLE) algorithms¹⁰ that implement the binormal model generate ROC curves that are “eccentric”¹¹ and have readily discernible hooks (“rd-hooks”) when plotted on linear coordinate axes. ROC curves with rd-hooks usually occur when −2 < r < 2. Hillis and Berbaum refer to r as the “mean-to-sigma ratio”¹² because $r = \frac{μ_{D +} - μ_{D -}}{σ_{D +} - σ_{D -}}$ ; they also treat r as a “measure of improperness” because binormal ROCs with smaller absolute values of r have more obvious hooks. ROCs with rd-hooks in the upper-right corner of ROC space occur when 0 < r < 2, and ROCs have rd-hooks in the lower-left corner when −2 < r < 0.

ROCs with Discernible Hooks

ROC curves with rd-hooks call attention to problems with the binormal model that have theoretical and practical importance. Any binormal ROC for which B ≠ 1 is “improper” in the sense that it violates a principle of signal detection theory, that is, rational use of diagnostic information produces a ROC curve with an ever-decreasing likelihood ratio as the decision variable.¹³ Using the likelihood ratio (or posterior probability) as the decision variable in a 2-choice classification task necessarily produces a proper ROC curve with an always-positive but monotonically decreasing slope.¹⁴ Binormal ROCs wrongly imply suboptimal or “improper” use of diagnostic information for the portion of the decision axis that lies beyond the curve’s inflection point, which is located where z(fpr) = $\frac{A B}{1 - A B^{2}}$ , and below-chance diagnostic differentiation where the curve lies below the no-information diagonal.¹² These regions can be especially large when the ROC has an rd-hook.

If improperness were a real feature of the diagnostic data, the diagnostician would simply treat the data differently to avoid inefficient use of information. Figure 1 illustrates a concrete example of this. Suppose hematocrit were used to detect the presence of blood dyscrasias, and suppose that testing of known D+ and D− populations led to the (made-up) results shown in Figure 1a, where most D+ patients are anemic (and have below-normal hematocrit levels), some D+ patients have levels in the “normal” range, and a few D+ patients are plethoric and have high levels. Using hematocrit values as the decision axis (going from lower to higher hematocrit levels) would produce the solid, obviously hooked ROC curve shown in Figure 1b. But a better approach to these data would transform the decision axis, using distance from 45.5% (the median of the D− population) as the decision criterion. This would result in the dashed, nonhooked ROC—and in superior diagnostic performance.

When ROCs have rd-hooks, they may give a misleading picture of diagnostic detection throughout much of the ROC square.¹⁵ As we demonstrate below, when rd-hooks appear in the upper right corner, the low-fpr portion of the ROC may be biased upward (i.e., it may imply higher sensitivity than is actually the case).

Commenting on data from radiologic studies, some authors suggest that as a practical matter, ROC with discernible hooks arise mainly from small data sets or from “degenerate” categorical data with poorly allocated boundaries (e.g., multiple thresholds where fpr < 0.02), such that several rating categories contain no subjects.^16,17 However, we have encountered several data sets that led to rd-hooked ROCs but did not have these limitations. Although it is difficult to specify all the circumstances that can produce binormal ROCs with rd-hooks, one common situation is that the data have a lower bound and the variance of the D− population is much smaller than the variance of the D+ population.

For example, laboratory test results sometimes yield rd-hooked ROCs if disease is associated with elevated concentrations of an analyte (values of which have a lower bound of zero) and the D− population has low levels of the analyte. This happens, for example, if one applies binormal ROC methods to the cortisol level data published in Carroll and colleagues’ seminal article on the dexamethasone suppression test (DST) in depression.^18,19 Carroll and colleagues¹⁸ present their DST findings from 368 patients as a histogram with 33 concentration categories; the histogram shows that following dexamethasone administration, cortisol levels of depressed patients tend to be less suppressed (i.e., they are higher) than the cortisol levels of nondepressed subjects. For the present analysis, we evaluated these data by regrouping them into 13 categories with boundaries chosen so that the locations of the 12 (fpr, tpr) pairs computed from the 13-category data—which appear as large circles in Figure 2—would be well-spaced, with all categories containing at least one member. Using Metz’s ROCKIT software,²⁰ we obtained binormal parameters A = 0.7075 and B = 0.6185. In this case, then, r = 1.85, and the binormal ROC curve has an rd-hook (see the dashed line in Figure 2). Because the distribution of results from depressed (D+) subjects is much wider than the results from the nondepressed (D−) group, application of the binormal model leads (implausibly and incorrectly) to the conclusion that a subset of depressed individuals have postdexamethasone cortisol levels that fall below typical levels from nondepressed individuals.

To address problems like these, investigators have developed ROC models with mathematical constraints that require monotonically decreasing slopes and thereby avoid hooked fits to data. Although these models yield proper ROCs, they have limitations. The “proper binormal” ROC model introduced by Birdsall²¹ and developed by Metz and colleagues^22,23 uses the likelihood ratio of the conventional (improper) ROC model as the new decision variable. This binormal likelihood-ratio (“binormal-LR,” adopting the terminology of Hillis and Berbaum¹²) model involves a nonmonotonic transformation of the conventional binormal decision variable that (so to speak) divides the conventional binormal decision axis into many small line segments, then reorders the line segments according to the values of their slopes such that the overall slope of the new, proper ROC always decreases. A drawback of this approach is that it produces convoluted relationships among the original decision variable, the intuitively appealing conventional binormal parameters, and the binormal-LR parameters.

A second approach to proper ROC construction assumes that diagnostic data arise from 2 latent gamma distributions. This “bigamma model”¹⁷ produces proper ROCs that rise rapidly from the left lower corner and taper slowly, that is, ROCs that, under the conventional binormal ROCs, would have B < 1. ROC curves constructed using Lehmann’s assumptions²⁴ have similar limitations. A third approach assumes “contamination”: The D+ population arises from multiple mixed distributions, one component of which is identical to the D− distribution.²⁵ This assumption is not always plausible, however, and it requires additional distributional constraints to produce ROCs with ever-decreasing positive slopes.

In this article, we present an approach to ROC curve fitting that implements the exponential features of Lehmann-family ROCs but allows a larger variety of curve shapes by weighting 2 mirror-image power functions. In the following sections, we describe the construction of a weighted power function (WPF) ROC curve and outline why it is a proper ROC with a positive, monotonically decreasing slope. We then describe and present results from a simulation study that compares outcomes from the WPF and the conventional binormal model of ROC curve fitting.

Constructing a Weighted Power Function ROC Curve

A WPF ROC is a set of points which, measured from the positive diagonal (the no-information line) in a ROC square, comprise a weighted average of 2 power functions. We now explain how one constructs a WPF ROC, an example of which is the solid black line in Figure 3.

In Figure 3, the heavy dotted line has the formula

y = f_{1} (x) = x^{a},

and the dashed line has the formula

y = f_{2} (x) = 1 - {(1 - x)}^{1 / a},

where 0 ≤ x ≤ 1 and 0 < a ≤ 1. Notice that f₁(x) and f₂(x) are mirror-images or reflections about the negative diagonal (that is, the line y = 1 – x) of the ROC square. Line λ is the no-information diagonal (that is, the line y = x), and line θ is perpendicular to λ and therefore has a slope of −1. The WPF ROC curve is the thick solid line in Figure 3.

For notational clarity, we will designate the particular (fpr, tpr) pair at decision threshold t as (x₃, y₃). The location of threshold t = (x₃, y₃) is determined by point (x₁, y₁) on y = f₁(x) and by point (x₂, y₂) on y = f₂(x). Given a value (x₁, y₁), the equation describing line θ is y = x₁ + x₁^a – x. Therefore, (x₂, y₂) is at a point along f₂(x) such that x₂ is the root of the following equation:

1 - {(1 - x_{2})}^{1 / a} {- x}_{1} - x_{1}^{a} {+ x}_{2} = 0 .

Equation 3 does not have an analytic solution, but one can use successive approximations (e.g., by implementing either the secant method or Halley’s rational method) to find a root in the region of interest. (We have found that by using x₁ as the initial “guess” for the solution to Equation 3, just 2 iterations of Halley’s method yield a value of x₂ that is correct to at least 8 decimal places.) From this value, one then can compute y₂ using Equation 2b.

Finally, the values of x₃ and y₃ are weighted averages of (x₁, y₁) and (x₂, y₂), that is

x_{3} = {bx}_{1} + (1 - b) x_{2}

y_{3} = {by}_{1} + (1 - b) y_{2},

where 0 ≤ b ≤ 1, and y₁ and y₂ are calculated from Equations 2a and 2b. Note that if a, b, and x₁ are fixed, then Equations 2a, 2b, 3, 4a, and 4b determine the values of y₁, x₂, y₂, x₃, and y₃, respectively.

When b equals or is very close to 1, WPF ROC forms will be identical or very close to ROCs generated from the Lehmann family of power ROC curves. However, WPF ROCs have more flexibility than the Lehmann curves because b can take on values between 0 and 1, with values <0.5 creating ROCs that rise relatively slowly from the lower left corner and that flatten and lie close to the upper boundary of the ROC square at higher fpr values. Note, finally, that when b = 0.5, the ROC is symmetric about the negative diagonal of the ROC square.

It is easy to show that the first derivatives of f₁(x) and f₂(x) are always positive and that their second derivatives are always negative for 0 ≤ x ≤ 1. Because Equations 4a and 4b imply that points along the ROC curve are weighted averages of points in f₁(x) and f₂(x), it seems plausible that the WPF ROC always has a positive but ever-decreasing slope. In Appendix I (online), we show that this is the case.

One can implement the construction of a WPF ROC in several ways. We used an Excel spreadsheet and chart to create Figure 3, and we have developed Excel spreadsheets, WinBUGS code, and SAS code for finding the parameters of the conventional binormal and WPF models based on ordinal category ROC data. The spreadsheets use Excel’s “Solver” function to select model parameters (i.e., estimates of operating point locations and values of A and B or a and b) that either maximize the likelihood function (discussed further below) or (as an alternative) minimize the chi-square or likelihood chi-square²⁶ (procedures that yield results similar and sometimes preferable to MLE methods²⁷). WinBUGS is a free software program that allows users to implement Markov chain Monte Carlo methods to find Bayesian solutions to statistical problems.²⁸ In Appendix II (online), we provide the SAS code to help readers better understand the relationships among WPF parameters and how the x₂ values are approximated via Halley’s method.

A useful feature of the WPF ROC derives from the fact that f₁(x) and f₂(x) in Figure 3 have closed-form integrals, and for the span 0 ≤ x ≤ 1, the areas under f₁(x) and f₂(x) equal $\frac{1}{1 + a}$ . Because the WPF ROC is a weighted average of f₁(x) and f₂(x), one would expect that its area also will equal $\frac{1}{1 + a}$ , and we show this is true in Appendix I.

Likelihood Function

Having explained how one constructs a WPF ROC, we now describe its likelihood function, adopting and modifying a conceptualization scheme used by several previous authors, but following Metz and Pan²² most closely.

We begin with ordinal, K-category data or with continuous data that have been grouped for purposes of mathematical analysis into K categories. We posit an underlying latent decision axis partitioned into K categories by K−1 thresholds t_k (k = 1, 2, …, K−1), with a progressively less demanding positivity threshold as k increases or, if the outcome space is not a simple ordinal scale, we assume that one can chop up the outcome space so that a subject who is positive according to threshold criterion (k – 1) is also positive with threshold k. Thus, the directionality of the decision axis is oriented such that fpr_k > fpr_k₋₁. We denote the probability of getting an outcome (or response) in any category k as p_k for an actually D− case and as q_k for an actually D+ case and assume that the K-category data come from N_D− independent, actually D− trials and N_D+ independent, actually D+ trials. Then the ordinal categorized response data from an experiment will consist of the outcome vector v ≡ {v₁, v₂, …, v_k | ∑_k v_k = N_D−} for the actually D− cases and outcome vector w ≡ {w₁, w₂, …, w_k | ∑_k w_k = N_D+} for the actually D+ cases. Because we assume that the trials are independent, we assume the random vectors v and w follow independent multinomial distributions. Then the probability or likelihood L of obtaining any particular set of data is

L = N_{D -}! N_{D +}! \prod_{k = 1}^{K} \frac{p_{k}^{v_{k}} q_{k}^{w_{k}}}{v_{k}! w_{k}!} .

Producing Equation 5 has not required making any assumption about the functions that give rise to p_k and q_k or the relationship (if any) between these probabilities. As we have seen, both the binormal and WPF ROC models describe the relationship between p_k and q_k in terms of 2 ROC curve parameters (A and B in the binormal model, a and b in the WPF model). This means that given a particular set of categorical outcome data, the K−1 thresholds plus the 2 curve parameters, or a total of K+1 parameters, determine the value of L. Coming up with MLE estimates of the K+1 parameters involves finding the values of those parameters that maximize L.

The natural logarithm of L, ln L,

ln L = ln (N_{D -}! N_{D +}!) + \sum_{k = 1}^{K} v_{k} ln p_{k} + \sum_{k= 1}^{K} w_{k} ln q_{k} - \sum_{k = 1}^{K} ln v_{k}! - \sum_{k = 1}^{K} ln w_{k}!

is a monotonic function of L but yields values that are easier to work with than L itself. Notice that the first, fourth, and fifth terms on the right side of Equation 6 do not depend on model parameters. We can delete these terms to create LL, the “log-likelihood function,”

L L = \sum_{k = 1}^{K} v_{k} In p_{k} + \sum_{k = 1}^{K} w_{k} In q_{k} .

Because LL varies monotonically with L, those values of the K+1 parameters that maximize LL will also maximize L.

As we have described the binormal model, the decision axis is measured in standard deviations of the D− population, so that the thresholds t_k are in z-transformed units of fpr. The probability p_k that an outcome from an actually D− case will fall in category k is

p_{k} - {fpr}_{k} - {fpr}_{k - 1} = Φ z ({fpr}_{k}) - Φ z ({fpr}_{k - 1}) = Φ (t_{k}) - Φ (t_{k - 1})

and the probability q_k that an outcome from an actually D+ case will fall in category k is

q_{k} = {tpr}_{k} - {tpr}_{k - 1} = Φ z ({tpr}_{k}) - Φ z ({tpr}_{k - 1}) = Φ [A + B z (t_{k})] - Φ [A + B z (t_{k - 1})],

where z(t₀) ≡ −∞ and z(t_K) ≡ +∞ (so that fpr₀ = tpr₀ = 0 and fpr_K = tpr_K = 1). By substituting Equations 8a and 8b into Equation 7, we obtain LL in the binormal model for the data represented earlier as vectors v and w. Under the binormal assumption, LL is maximized when the first derivatives of LL with respect to the K+1 parameters A, B, t₁, t₂, …, t_k all equal zero. The binormal ROC software programs developed by Dorfman and Alf²⁹ and Metz and others²⁰ trace their ancestry to Fortran computer code that implemented a Newton-Raphson iteration procedure and the “method of scoring” to find solutions to the K+1 nonlinear equations produced by partial differentiation of Equation 7.³⁰ Note, however, that several numerical optimization algorithms might be used to find a maximum value of LL.

We take advantage of this last point as we consider maximizing LL under the WPF model. Recall from our earlier discussion that (fpr, tpr) pairs under the WPF model correspond to values of x₃ and y₃ calculated in Equations 4a and 4b. We now introduce the “k” subscript so that we can speak of multiple thresholds t_k represented by (x_3,k, y_3,k) pairs. In the WPF model, these thresholds simply refer to points along the two axes of the ROC square, so that

p_{k} {= fpr}_{k} - {fpr}_{k - 1} = t_{k} - t_{k - 1} = x_{3, k} - x_{3, k - 1}

q_{k} {= tpr}_{k} - {tpr}_{k - 1} = y_{3, k} - y_{3, k - 1} .

Although we cannot write an expression analogous to Equation 9a for the values of q_k, we know from Equations 2a through 4b above how to establish values for q_k in Equation 9b from a, b, and x_1,1, x_1,2, …, x_1,k. Therefore, we can treat a, b, x_1,1, x_1,2, …, x_1,k as the K+1 parameters for which we seek MLE estimates using numerical optimization methods to maximize LL in Equation 7.

Some Examples

We next illustrate how the WPF ROC solves problems associated with data sets that generate improper binormal ROCs. In Figures 4a and 4b, we compare WPF and binormal ROCs derived from hypothetical 5-category data obtained from 20 randomly selected D+ and 20 randomly selected D− subjects. In the hypothetical data represented in Figure 4a, the 40 subjects come from a diagnostic context in which the actual (“true”) thresholds lie at fpr = 0.01, 0.02, 0.05, and 0.5, and the D+ and D− populations give rise to data that, when transformed, are normally distributed such that A = 1.5 and B = 0.8. Notice that these parameters do not, by themselves, imply an rd-hooked ROC because r = $\frac{A}{1 - B}$ = 7.5. However, the probability that a single, randomly selected D− subject would be found within the first 3 rating categories is 0.05, so in a random sample with 20 D− subjects, the probability of getting a data set in which no D− subject falls in the first 3 categories is (1 − 0.05)²⁰ ≈ 0.36. Figure 4a depicts one such result. The binormal ROC curve calculated from these “degenerate” data contains a noticeable hook and deviates substantially from the “true” ROC that would be derived from the full D+ and D− populations. By contrast, the WPF ROC constructed from the same hypothetical data set, in addition to being proper, displays a much closer fit to the true ROC.

Next, consider hypothetical 5-category data obtained from 20 D+ and 20 D− randomly selected subjects who come from a diagnostic context in which the actual operating points lie at fpr = 0.02, 0.06, 0.10, and 0.5 and the distributions of the full D+ and D− populations are described by the binormal parameters A = 1.875 and B = 1.25 (i.e., the true ROC is a mirror image about the negative diagonal of the true ROC curve shown in Figure 4a). Several varieties of 20-member random samples from the D− population would have a good chance of resulting in an rd-hooked ROC; among these are samples in which 1 D− subject fell in the first rating category, fewer than 2 subjects fell in the second, and no subject fell in the third—a situation that we expect to arise approximately 26% of the time. Figure 4b depicts one such result. The binormal ROC curve calculated from these data has an rd-hook in the lower left corner, and it gives a misleading indication of diagnostic performance throughout the ROC square except where 0.07 ≤ fpr ≤ 0.11. The WPF ROC is biased upward for fpr values between 0.1 and 0.5, but it generally gives results closer to the true ROC.

For a final example, we return to Carroll and colleagues’ DST data. Here, we do not know what the “true” values are—we have only the empirical data from the sample and our knowledge of the clinical context as the bases for judgments about the relative merit of binormal and WPF models. If one quickly compared the dashed and dotted lines in Figure 2 and noted only that the vertical distances from the initial data points to the binormal ROC were less than the vertical distances to the WPF ROC, one might conclude that the binormal ROC followed the original data more closely than the WPF ROC and gave a better depiction of the DST’s performance, especially where fpr > 0.2.

More detailed investigation shows that this is not so, however. One role of ROC analysis is to estimate locations of the actual or true operating points, which have multinomial distributions such that initial empirical data can contain errors in the horizontal as well as vertical direction. In Figure 2, we plotted the locations of the operating points as estimated by the ROCKIT software as diamonds, and we used squares to show the locations of operating points obtained using WinBUGS code that implements the WPF model. We then used straight lines to connect each model’s operating point estimates to 5 of the initial data points (the circles in Figure 2) as calculated from the crude categorical data. In all cases, the operating points estimated using the binormal software lie farther from the initial data points than do the operating points estimated from the WPF model. The average distance, δ_mean, between empirical and model operating points for the k = 12 thresholds, calculated from the individual Euclidean distances in ROC space as

δ_{m e a n} = \frac{1}{k} \sum_{k = 1}^{12} \sqrt{{(f p r_{e s t i m a t e_{k}} - f p r_{i n i t i a l_{k}})}^{2} + {(t p r_{e s t i m a t e_{k}} - t p r_{i n i t i a l_{k}})}^{2}},

is 0.047 for the WPF ROC but 0.101 for the binormal ROC. This finding accords with the LL values calculated using Equation 7: −870.1 for the binormal ROC model and −852 for the WPF model, implying a better fit for the latter. Furthermore, the WPF ROC presents a more plausible view of the DST’s diagnostic accuracy by incorporating the assumption that at any postdexamethasone cortisol level, the proportion of depressed patients above that level will be higher than the proportion of nondepressed patients.

Evaluation via Simulation Study

We tested the WPF ROC using methods modeled on the Monte Carlo computer studies conducted by Pan and Metz¹⁵ and Dorfman and colleagues.¹⁷ In summary, the approach taken in these studies was as follows:

a. Set the parameters of the “true” D+ and D− distributions.

b. Draw sets of data samples from those distributions.

c. Estimate parameters for each sample’s data using the binormal curve fitting assumptions.

d. Estimate parameters for the same sample data with the WPF curve fitting routine.

e. Compare (c) and (d) to see which results came closer to the “true” parameters.

In the first set of Monte Carlo studies, we posited 3 conventional binormal population models with ROC parameters B = 0.8 and A = 0.6, 1.5, and 2.1; we selected these values because they would produce ROCs with shapes typical of actual diagnostic data from tests with modest, moderate, and high accuracy (AUCs = 0.68, 0.88, and 0.95, respectively). We generated 1000 sets of simulated 5-category rating data for 20 actually D− subjects and 20 actually D+ subjects, under the assumption that the actual decision thresholds t_k (k = 1, 2, 3, 4) were located at fpr = 0.01, 0.02, 0.05, and 0.50. These fpr values are not as extreme as the fpr values (viz., fpr = 0.000250, 0.000328, 0.03040, 0.281140) used by Pan and Metz and by Dorfman and colleagues, and our data may therefore have yielded more generalizable results. Still, however, our thresholds were chosen such that given the size of the D− and D+ groups, a substantial fraction of simulations were likely to contain “degenerate” data sets (i.e., data sets with no false positives in the first 3 rating categories).

In the second set of studies, we posited that data arose from 3 WPF ROC population models with b = 0.75 and a = 1/2, 3/17, and 1/19—values that imply ROC shapes similar to those from the just-mentioned binormal ROCs and for which the true AUCs are 0.667, 0.85, and 0.95. Again, we generated 1000 sets of simulated 5-category rating data assuming that 4 true operating points were located at fpr = 0.01, 0.02, 0.05, and 0.50, each data set consisting of diagnostic ratings for 20 actually D− subjects and 20 actually D+ subjects.

We generated and analyzed simulated data sets using SAS 9.1. We obtained estimates of the 6 ROC model parameters—the 2 binormal and WPF accuracy indices (either A and B, or a and b), and the 4 operating points—using the nonlinear programming procedure (“proc nlp”) to maximize the log-likelihood function (LL). Our code ran 2 optimization algorithms—the Newton-Raphson method with ridging (“nrridg”) and the double dogleg method (“dbldog”)—and then selected and retained the results that yielded the largest value of LL.

AUC often serves as a summary index for evaluating and comparing test quality, but as Figure 3 illustrates, WPF and binormal ROC curves with similar areas beneath them may cross each other and can have very different tpr values through large fpr ranges. Although we computed total AUC values, our main comparisons between WPF and binormal models used 1) the calculated tpr values at fpr = 0.01, 0.02, 0.05, 0.25, and 0.50, which allowed us to evaluate accuracy at a representative set of thresholds that are likely to be of diagnostic significance, and 2) following Pan and Metz,¹⁵ partial AUC values for the fpr ranges 0.0–0.25, 0.25–0.5, 0.5–0.75, and 0.75–1.0.

We compared results obtained from the binormal and WPF models using empirical mean values and the root mean square error (RMSE) from the tpr values and areas derived from computer-simulated data. RMSE, a frequently used metric for evaluating estimators, incorporates the variance and bias of an estimator $\hat{φ}$ and is defined as

R M S E = \sqrt{\frac{\sum_{i = 1}^{N} {({\hat{φ}}_{i} - φ_{t r u e})}^{2}}{N}}

where ${\hat{φ}}_{i}$ is the i^th estimate of the parameter $φ$ , $φ_{t r u e}$ is the true value of the parameter, and N is the number of estimates (here, N = 1000).

Results

Table 1 shows the means, medians, and RMSEs for tpr values where fpr = 0.01, 0.02, 0.05, 0.25, and 0.50 from the six 1000-instance binormal and WPF ROC simulations. The “Features” column describes the assumed “true” ROC parameters and the number of rd-hooked ROCs (i.e., those for which |r| < 2) obtained under the binormal model.

Table 1

Results of Simulation Studies: True Positive Rates When the False Positive Rate Equals 0.01, 0.02, 0.05, 0.25, or 0.50

Simulation			Binormal ROC					WPF ROC
No.	Features	fpr	0.01	0.02	0.05	0.25	0.5	0.01	0.02	0.05	0.25	0.5
1	Binormal	tpr_true	0.3590	0.4431	0.5730	0.8316	0.9332	0.3590	0.4431	0.5730	0.8316	0.9332
	A = 1.5	Median	0.4546	0.5252	0.6319	0.8407	0.9165	0.3707	0.4649	0.6067	0.8519	0.9368
	B = 0.8	Mean	0.4666	0.5203	0.6181	0.8286	0.9007	0.3817	0.4620	0.5990	0.8424	0.9289
	364 hooks	RMSE	0.2856	0.2550	0.1966	0.0979	0.0875	0.2311	0.2085	0.1656	0.0897	0.0506
2	Binormal	tpr_true	0.5944	0.6762	0.7835	0.9407	0.9821	0.5944	0.6762	0.7835	0.9407	0.9821
	A = 2.1	Median	0.7217	0.7779	0.8494	0.9477	0.9736	0.6206	0.7172	0.8268	0.9483	0.9769
	B = 0.8	Mean	0.6676	0.7241	0.8174	0.9368	0.9644	0.5547	0.6620	0.7962	0.9393	0.9718
	182 hooks	RMSE	0.2763	0.2342	0.1500	0.0555	0.0439	0.2631	0.2060	0.1266	0.0495	0.0261
3	Binormal	tpr_true	0.1036	0.1485	0.2370	0.5241	0.7257	0.1036	0.1485	0.2370	0.5241	0.7257
	A = 0.6	Median	0.2104	0.2043	0.2916	0.5391	0.6904	0.1285	0.1819	0.2918	0.6133	0.8093
	B = 0.8	Mean	0.1609	0.2426	0.3066	0.5228	0.6740	0.0960	0.1493	0.2650	0.6059	0.8095
	684 hooks	RMSE	0.2073	0.2002	0.1818	0.1421	0.1480	0.1180	0.1356	0.1697	0.1829	0.1409
4	WPF	tpr_true	0.2106	0.3570	0.6339	0.7866	0.9091	0.2106	0.3570	0.6339	0.7866	0.9091
	a = 3/17	Median	0.3802	0.4532	0.5695	0.7938	0.8934	0.2902	0.3958	0.5615	0.8244	0.9247
	b = 0.75	Mean	0.4052	0.4569	0.5521	0.7902	0.8795	0.3301	0.4107	0.5510	0.8129	0.9157
	392 hooks	RMSE	0.3270	0.2655	0.2212	0.1078	0.0968	0.2442	0.2045	0.1870	0.1031	0.0579
5	WPF	tpr_true	0.5559	0.7483	0.8337	0.9464	0.9735	0.5559	0.7483	0.8337	0.9464	0.9735
	a = 1/19	Median	0.7734	0.8206	0.8726	0.9529	0.9757	0.6763	0.7607	0.8460	0.9568	0.9788
	b = 0.75	Mean	0.6879	0.7384	0.8310	0.9445	0.9663	0.5894	0.6950	0.8199	0.9473	0.9747
	223 hooks	RMSE	0.3124	0.2489	0.1552	0.0515	0.0411	0.2575	0.2019	0.1107	0.0413	0.0205
6	WPF	tpr_true	0.0591	0.0997	0.1860	0.4871	0.7168	0.0591	0.0997	0.1860	0.4871	0.7168
	a = 1/2	Median	0.0983	0.1404	0.2247	0.4873	0.6574	0.0625	0.1114	0.2206	0.5836	0.8164
	b = 0.75	Mean	0.1596	0.1876	0.2475	0.5040	0.6871	0.0951	0.1467	0.2589	0.6006	0.8082
	668 hooks	RMSE	0.1896	0.1854	0.1723	0.1508	0.1597	0.0993	0.1261	0.1795	0.2168	0.1611

Note: RMSE = root mean square error; ROC = receiver operating characteristic; WPF = weighted power function.

Looking at results from the 6 studies, one sees that for fpr values 0.01, 0.02, and 0.05, the tpr values obtained using the WPF ROC method usually fall closer to the true values, display less upward bias, and have smaller RMSEs than do values obtained under the binormal assumption. In the studies where A = 0.6, B = 0.8 and a = 1/2, b = 0.75, the tpr values obtained under the binormal assumption where fpr = 0.25 or 0.50 look better than do those obtained using the WPF model. However, careful inspection of the data suggests that these better estimates were purchased at the price of incurring rd-hooked binormal ROCs in more than two-thirds of the simulations. That is, the binormal ROCs have a pronounced upward bias at low fpr values, flatten out in the middle of the fpr range, and dip underneath the no-information diagonal at higher values of fpr.

Table 2 is structured similarly to Table 1, but it describes what our 6 simulation studies told us about the means, medians, and RMSEs for partial ROC areas where fpr = 0.0–0.25, 0.25–0.5, 0.5–0.75, and 0.75–1.0. In all 6 studies, partial areas where 0 ≤ fpr ≤ 0.5 calculated from the WPF ROC parameters fell closer to the true values, displayed less upward bias, and had smaller RMSEs than did values obtained under the binormal assumption. Given the previously described findings for tpr values where A = 0.6, B = 0.8 and a = 1/2, b = 0.75, it is not surprising that the partial areas where 0.5 ≤ fpr ≤ 1 obtained from the binormal ROC parameters look slightly better than the results obtained with the WPF model.

Table 2

Results of Simulation Studies: Comparisons of Partial Areas under the ROC Curve for 4 False-Positive Ranges (0.0–0.25, 0.25–0.5, 0.5–0.75, 0.75–1.0) and Total ROC Area

Simulation			Binormal ROC					WPF ROC
No.	Features	fpr Range	0.0–0.25	0.25–0.5	0.5–0.75	0.75–1.0	Total	0–0.25	0.25–0.5	0.5–0.75	0.75–1	Total
1	Binormal	True area	0.1691	0.2225	0.2397	0.2478	0.8793	0.1691	0.2225	0.2397	0.2478	0.8793
	A = 1.5	Median	0.1784	0.2197	0.2340	0.2443	0.8735	0.1770	0.2253	0.2395	0.2470	0.8894
	B = 0.8	Mean	0.1766	0.2176	0.2299	0.2380	0.8622	0.1743	0.2232	0.2378	0.2465	0.8819
	364 hooks	RMSE	0.0341	0.0246	0.0189	0.0173	0.0793	0.0303	0.0170	0.0091	0.0030	0.0517
2	Binormal	True area	0.2109	0.2414	0.2475	0.2496	0.9495	0.2109	0.2414	0.2475	0.2496	0.9495
	A = 2.1	Median	0.2195	0.2401	0.2458	0.2490	0.9514	0.2161	0.2413	0.2460	0.2489	0.9523
	B = 0.8	Mean	0.2162	0.2382	0.2428	0.2456	0.9430	0.2120	0.2396	0.2451	0.2486	0.9453
	182 hooks	RMSE	0.0240	0.0119	0.0102	0.0082	0.0414	0.0207	0.0088	0.0048	0.0016	0.0308
3	Binormal	True area	0.0879	0.1579	0.2006	0.2339	0.6803	0.0879	0.1579	0.2006	0.2339	0.6803
	A = 0.6	Median	0.0984	0.1509	0.1831	0.2161	0.6463	0.1016	0.1799	0.2183	0.2417	0.7447
	B = 0.8	Mean	0.0999	0.1554	0.1859	0.2114	0.6526	0.1053	0.1804	0.2179	0.2412	0.7446
	684 hooks	RMSE	0.0373	0.0350	0.0383	0.0390	0.1123	0.0422	0.0411	0.0274	0.0101	0.1145
4	WPF	True area	0.1553	0.2145	0.2348	0.2455	0.8500	0.1553	0.2145	0.2348	0.2455	0.8500
	a = 3/17	Median	0.1625	0.2127	0.2308	0.2431	0.8434	0.1665	0.2207	0.2374	0.2464	0.8732
	b = 0.75	Mean	0.1629	0.2104	0.2258	0.2358	0.8350	0.1640	0.2182	0.2357	0.2459	0.8638
	392 hooks	RMSE	0.0367	0.0247	0.0233	0.0204	0.0785	0.0329	0.0196	0.0102	0.0031	0.0591
5	WPF	True area	0.2154	0.2406	0.2454	0.2487	0.9500	0.2154	0.2406	0.2454	0.2487	0.9500
	a = 1/19	Median	0.2230	0.2410	0.2462	0.2492	0.9557	0.2197	0.2423	0.2463	0.2489	0.9582
	b = 0.75	Mean	0.2193	0.2393	0.2429	0.2452	0.9468	0.2162	0.2408	0.2456	0.2487	0.9514
	223 hooks	RMSE	0.0230	0.0112	0.0097	0.0086	0.0400	0.0174	0.0074	0.0036	0.0010	0.0253
6	WPF	True area	0.0766	0.1522	0.2015	0.2364	0.6667	0.0766	0.1522	0.2015	0.2364	0.6667
	a = 1/2	Median	0.0864	0.1438	0.1828	0.2205	0.6256	0.0947	0.1771	0.2220	0.2442	0.7361
	b = 0.75	Mean	0.0872	0.1511	0.1868	0.2136	0.6387	0.0992	0.1789	0.2181	0.2416	0.7377
	668 hooks	RMSE	0.0362	0.0374	0.0410	0.0400	0.1131	0.0474	0.0483	0.0299	0.0097	0.1298

Note: RMSE = root mean square error; ROC = receiver operating characteristic; WPF = weighted power function.

Because total AUC is often used as a summary index of diagnostic performance (despite the limitations noted earlier), we present AUC values under the 2 methods. A general conclusion from Table 2 is that the total AUCs of binormal ROCs have a downward bias and the AUCs for the WPF ROCs are biased upward. This result implies that for binormal ROCs, large upward biases in the low fpr regions are more than cancelled out by large downward biases in high fpr regions, especially for data from low-accuracy diagnostic tests. Particularly in Simulation 6, the downward bias of binormal ROCs in high-fpr regions leads to a smaller bias in total AUC compared with results obtained under the WPF ROC model.

Figure 5 portrays the tpr results from the 6 simulation studies as box-and-whisker plots of the 2.5th percentile, lower quartile, median, upper quartile, and 97.5th percentile tpr values at the 5 previously described fpr values. The figure provides readily apprehended descriptions of and insights into data that are incompletely summarized in Table 1. In all 6 plots, the 5 fpr values are the bottom labels, the WPF results appear on the left, the binormal results appear on the right, and the mean tpr results for each fpr value are marked with a solid diamond; the results are graphed as deviations (“errors,” with units marked along the vertical axis) from the true tpr values. The key point from these figures is that in most cases, the WPF method gives results that are less dispersed and that come closer to the true values of tpr.

Discussion

This article has described the weighted power function model for constructing proper ROC curves from ordered category data and has compared results using this model with those obtained from applying the (justifiably) popular “binormal” assumptions of ROC curve fitting. Previous authors have mentioned single-parameter power function ROCs (of the variety tpr = fpr^a) as candidates for ROC forms.^7,24,31 Like power function ROCs, the WPF ROC has a positive but ever-decreasing slope in the ROC square. However, the 2-parameter WPF model described here is more flexible than a single-parameter power function and can fit a wider variety of actual data sets—in particular, data sets where tpr values in the lower left portion of the ROC square lie close to the no-information diagonal (that is, data sets for which the binormal ROC parameter B would be greater than 1). In contrast to investigators who have focused on radiologic data, psychiatrists not infrequently encounter data for which B > 1 (see, for example, references 5 and 32). Therefore, the capacity of the WPF ROC to accommodate such data seems an important criterion for judging its usefulness as a curve-fitting model.

We evaluated the WPF model using methods similar to those employed by previous authors: simulation studies that were structured to produce a large number of the degenerate data sets that pose problems for the binormal model. We found that the WPF model generally provided solutions that were better—that is, less biased and closer to the true values—than the results obtained using the binormal model, even when binormal distributions were the source of the simulated data. A big contributor to the better performance of the WPF model was that WPF ROCs have monotonically decreasing slopes and cannot have hooks, whereas binormal ROCs often compensate (so to speak) for their hooks with too-high tpr values in the range of ROC space where fpr < 0.10. Thus, the WPF model performs better because of its design constraint as a necessarily proper ROC.

Our study did not compare results from the WPF and other proposed proper ROC models, and we therefore have not provided any direct statistical evidence that our approach is superior. But because the WPF model produces proper results, it seems to be a valid option for data analyses where fitting a curve to data is important but the binormal model would yield ROCs with unacceptable hooks. The choice of a proper ROC model will depend on situation-specific factors and statisticians’ preferences. One can implement the WPF model with an Excel spreadsheet, WinBUGS, or SAS (see sample code in Appendix II), and the construction of WPF ROCs is easy to conceptualize and to depict graphically. The WPF model also establishes straightforward, easily depicted relationships among the results that comprise the decision variable for the detection method, individual operating points reflected in the data, the contours of the resulting ROC curve, and model parameters a and b.

The WPF model overcomes the susceptibility to rd-hooks that can arise from the conventional binormal ROC model, is more flexible than the Lehmann curve and bigamma ROC models, and seems simpler than the binormal-LR model. The WPF ROC curve-fitting code we have developed offers investigators and researchers convenient tools for analyzing data from a wide variety of medical diagnostic technologies and methods.

Footnotes

This publication was supported by an Institutional Clinical and Translational Science Award, NIH/NCATS Grant No. 8UL1TR000077-04. Its contents are solely the responsibility of the authors and do not necessarily represent the official views of the NIH.

Supplementary material for this article is available on the Medical Decision Making Web site at .

References

Pepe

. The Statistical Evaluation of Medical Tests for Classification and Prediction. New York: Oxford University Press; 2003.

Zhou

X-H

Obuchowski

McClish

. Statistical Methods in Diagnostic Medicine. New York: Wiley-Interscience; 2002.

Wagner

Metz

Campbell

. Assessment of medical imaging systems and computer aids: a tutorial review. Acad Radiol. 2007;14:723–48.

Pesce

Metz

Berbaum

. On the convexity of ROC curves estimated from radiological test results. Acad Radiol. 2010;17:960–8.

Mossman

. Resampling techniques in the analysis of non-binormal ROCs. Med Decis Making. 1995;15:358–66.

Swets

. Indices of discrimination or diagnostic accuracy: their ROCs and implied models. Psychol Bull. 1986;99:100–17.

Swets

. Form of empirical ROCs in discrimination and diagnostic tasks: implications for theory and measurement of performance. Psychol Bull. 1986;99:181–98.

Swets

Pickett

Whitehead

. Assessment of diagnostic technologies. Science. 1979;205:753–9.

Hanley

. The robustness of the “binormal” assumptions used in fitting ROC curves. Med Decis Making. 1988;8:197–203.

10.

Dorfman

Alf

. Maximum-likelihood estimation of parameters of signal-detection theory and determination of confidence intervals—rating method data. J Math Psychol. 1969;6:487–96.

11.

Somoza

. Eccentric diagnostic tests: redefining sensitivity and specificity. Med Decis Making. 1996;16:15–23.

12.

Hillis

Berbaum

. Using the mean-to-sigma ratio as a measure of the improperness of binormal ROC curves. Acad Radiol. 2011;18:143–54.

13.

Green

Swets

. Signal Detection Theory and Psychophysics. New York: Wiley; 1966.

14.

Egan

. Signal Detection Theory and ROC Analysis. New York: Academic Press; 1975.

15.

Pan

Metz

: The “proper” binormal model: parametric ROC curve estimation with degenerate data. Acad Radiol. 1997;4:380–9.

16.

Metz

. Some practical issues of experimental design and data analysis in radiological ROC studies. Invest Radiol. 1989;24:234–45.

17.

Dorfman

Berbaum

Metz

. Proper receiver operating characteristic analysis: the bigamma model. Acad Radiol. 1997;4:138–49.

18.

Carroll

Feinberg

Greden

. A specific laboratory test for the diagnosis of melancholia: standardization, validation, and clinical utility. Arch Gen Psychiatry. 1981;38:15–22.

19.

Mossman

Somoza

: Maximizing diagnostic information from the dexamethasone suppression test: an approach to criterion selection using receiver operating characteristic analysis. Arch Gen Psychiatry. 1989;46:653–60.

20.

Metz

Herman

Roe

. Statistical comparison of two ROC-curve estimates obtained from partially-paired datasets. Med Decis Making. 1998;18:110-21. ROCKIT is available from: URL: http://www-radiology.uchicago.edu/krl/KRL_ROC/software_index6.htm.

21.

Birdsall

. The theory of signal detectability: ROC curves and their character. Doctoral thesis. Ann Arbor, MI: University of Michigan; 1966.

22.

Metz

Pan

. “Proper” binormal ROC curves: theory and maximum likelihood estimation. J Math Psychol. 1999;43:1–33.

23.

Pesce

Metz

. Reliable and computationally efficient maximum likelihood estimation of “proper” binormal ROC curves. Acad Radiol. 2007;14:814–29.

24.

Gönen

Heller

. Lehmann family of ROC curves. Med Decis Making. 2010;30:509–17.

25.

Dorfman

Berbaum

. A contaminated binormal model for ROC data, part II: a formal model. Acad Radiol. 2000;7:427–37.

26.

Berkson

. Minimum chi-square, not maximum likelihood! Ann Stat. 1980;8:457–87.

27.

Keating

Mason

Sen

. Pitman’s Measure of Closeness: A Comparison of Statistical Estimators. Philadelphia, PA: Society for Industrial and Applied Mechanics; 1993.

28.

Lunn

Spiegelhalter

Thomas

Best

. The BUGS project: evolution, critique and future directions. Stat Med. 2009;28:3049–67.

29.

Dorfman

Alf

. Maximum likelihood estimation of parameters of signal detection theory and determination of confidence intervals—rating method data. J Math Psychol. 1969;6:487–96.

30.

Schönemann

Tucker

. A maximum likelihood solution for the method of successive intervals allowing for unequal stimulus dispersions. Psychometrika. 1967;32:403–17.

31.

Egan

Greenberg

Schulman

. Operating characteristics, signal detectability, and the method of free response. J Acoust Soc Am. 1961;33:993–1007.

32.

Barzman

Brackenbury

Sonnier

. Brief Rating of Aggression by Children and Adolescents (BRACHA): development of a tool to assess risk of inpatients’ aggressive behavior. J Am Acad Psychiatry Law. 2011;39:170–9.

Supplementary Material

Please find the following supplemental material available below.

For Open Access articles published under a Creative Commons License, all supplemental material carries the same license as the article it is associated with.

For non-Open Access articles published, all supplemental material carries a non-exclusive license, and permission requests for re-use of supplemental material or any part of supplemental material shall be sent directly to the copyright owner as specified in the copyright notice associated with the article.

0.00 MB

0.07 MB

0.01 MB