Automated post-run analysis of arrayed quantitative PCR amplification curves using machine learning

Introduction

Quantitative PCR has become an important approach for detection of pathogens across etiologic studies of common infectious diseases syndromes, including diarrhea, pneumonia, and acute febrile illness, sepsis, and meningitis.^1–5 Because of the wide range of pathogens that can be implicated for each of these syndromes, platforms that provide simultaneous detection of a number of pathogens are preferred, and the use of quantitative PCR has helped differentiate etiologic detections from background carriage^6,7 as well as to interpret drug resistance and pathogen type from the stool metagenome.^8,9 The TaqMan Array Card (TAC) is one such platform and has the advantage of being a closed system with arrayed singleplex or duplex target detection of both DNA and RNA targets, as long as the assays can be adapted to universal cycling conditions.⁶

A limitation of the use of these array cards is the need for manual post-run adjustment of target baselines. While the software that accompanies the instrument can perform an automated analysis of baselines and returns a cycle threshold (Ct) for each amplification curve, these are subject to error and thus require manual inspection and correction. Because each sample lane can include duplex qPCR assays in each of 48 sub-wells, with eight samples per card, there can theoretically be up to 768 individual PCR amplifications per card. Large studies may include thousands of samples and thus in excess of a million PCR amplifications to analyze manually. Further, while a protocol is used to guide the post-run analysis, it is intrinsically subject to interpretation and can result in variability in Ct calls between individual analysts.

A fully automated alternative to manual post-run analysis would save time and effort and improve reproducibility. We thus sought to use machine learning approaches to derive models for prediction of both PCR amplification and, when amplification was observed, Ct values using TAC data from two previously-published enteric infectious diseases studies.^10,11 We then performed out-of-sample validation and compared the prediction model to post-run analyses performed by laboratories as part of an external quality assessment (EQA).¹⁰

Methods

Results

In 5-fold cross-validation, we first identified the optimal hyperparameter values for the two models ( Table 1). In internal validation on the 20% testing set, we achieve an accuracy of 0.996, a sensitivity of 0.997, a specificity of 0.993, a PPV of 0.998, and NPV of 0.991 ( Table 2). Of 33,043 reaction in the testing set, 125 (0.38%) were misclassified compared to the manual gold standard ( Figure 1). Our model misclassified 53 curves as not amplified, most of which had exponential growth at high cycle counts (mean=35.7, SD=4.86, 76% with a ground truth Ct ≥ 35). For the 72 curves our model misclassified as amplified, the mean predicted Ct was 33.2 (SD=5.33). We achieved an MAE of 0.590 when predicting the Ct value on the internal testing set. Excluding the reactions that the model incorrectly classified as true amplifications, the MAE reduced to 0.544 ( Figure 2). Model performance was similar between dyes and for singleplex vs. duplex reactions. Of the 125 misclassified reactions, 81 (64.8%) were flagged based on a continuous model prediction value between 0.1 and 0.9, in comparison to 124 (3.8%) of 32918 correctly classified reactions. The model also outperformed the automated relative threshold included in the machine software, which had an accuracy of 99.2% and a MAE of 1.00.

Figure 1. Plots of curves that were misclassified (off-diagonal) and 20 randomly-selected correctly classified curves (main-diagonal), with the cycle on the x-axis and normalized fluorescence on the y-axis in internal validation.

Figure 2. Manual (gold standard) Ct vs. automated (model-predicted) Ct in internal validation.

In external validation, we achieved an accuracy of 0.997, a sensitivity of 0.998, a specificity of 0.991, a PPV of 0.997, an NPV of 0.994, and a MAE of 0.611 when classifying amplification and Ct values ( Table 2 and Figure 3). The MAE reduced to 0.579 when we excluded misclassified targets. Five of 1472 curves (0.34%) were incorrectly classified by the automated analysis compared to the manual gold standard ( Figure 4). For the three curves identified as true amplification by the manual gold standard but which the automated analysis missed, the Ct values approached 40 (mean=38.7, SD=0.59). The automated analysis had an accuracy for reaction classification better than 14 of 17 laboratory scientists ( Table 3). With a MAE for Ct prediction of 0.611, the automated analysis outperformed four of the 17 labs.

Figure 3. qPCR Ct predicted by either the automated (model-based) analysis or manual analyses performed by laboratory scientists as part of an EQA (Y axis), all shown in comparison to the manual gold standard (X axis).

Manual analysis plots are ordered by site accuracy from lowest to highest.

Figure 4. Amplification curves for incorrectly classified reactions in the external validation set.

The vertical red line indicates the Ct predicted by the automated analysis and the vertical black line indicates the Ct assigned by the manual gold standard. In either case, a value of 40 was assigned if no amplification was identified.

Discussion

Using a large sample of qPCRs that were previously assessed using the background threshold method, we were able to accurately predict both presence of amplification and the qPCR Ct directly from normalized fluorescence values. In external validation, the accuracy was superb and as good as the best-performing laboratory scientists in a multisite EQA analysis. It also outperformed the automated algorithms included in the commercial qPCR software. This approach will allow for a substantial reduction in hands-on laboratory time, laboratory training and EQA, and provide an automated, fully reproducible method for the post-run analysis of arrayed qPCR.

Previous approaches to automating qPCR amplification and cycle threshold assessment have focused on first characterizing the curves by extracting the structure and key features, for example presence of a hook effect.^13,14 However, this requires additional computation and abstraction of the fluorescence data. Our approach minimizes data manipulation and is computationally efficient. We chose the tree-based algorithm, XGBoost, due to its ability to recognize complex interactions between predictors and nonlinear relationships between the predictors and outcomes, in addition to its efficiency.¹⁵ All of the data used in this analysis came from qPCR on the stool metagenome. This is a convenient source of data for training and validating these models because of the high rate of target amplification (about one quarter of all qPCR reactions) and because of the wide range in pathogen quantity present. However, we would expect that these models would work equally well on qPCR reactions derived from other sample types. While the accuracy of amplification prediction in external validation was excellent, the cycle threshold error was higher than most of the labs, however the absolute difference in error was very small. We believe that this is a more than acceptable trade-off given the other clear advantages to this approach.

This automated post-run pipeline will be of particular value for multisite studies, removing an important source of between-site variation, as well as for studies where quick turnaround of qPCR data is valuable such as for clinical decision-making. The prediction models can be operationalized e.g. with a web portal, allowing for full-automated and reproducible amplification prediction, cycle threshold calling, data cleaning, and data aggregation. As with most evaluations of diagnostics, the absence of a true gold standard is a limitation, however in this case the goal was to consistently reproduce data generated using the best-available manual approach. In summary, we were able to accurately predict qPCR amplification and prediction by training machine learning models directly on normalized fluorescence data, which will allow for a substantial improvement in the efficiency and transparency of data generation.