| Title: | Analysis of OpenArray PCR Data |
|---|---|
| Description: | Provides a suite of R functions to analyze gene expression experiments on the OpenArray real-time PCR platform. OAtools fits logistic regressions to fluorescence curves to distinguish between real amplification and false positives. OAtools supports data import, analysis, and visualization through plots and a dynamic HTML report. |
| Authors: | Aidan Shea [aut, cre] (ORCID: <https://orcid.org/0009-0001-9256-7439>) |
| Maintainer: | Aidan Shea <[email protected]> |
| License: | GPL (>= 3) |
| Version: | 1.1.0 |
| Built: | 2026-05-30 09:50:24 UTC |
| Source: | https://github.com/bioc/OAtools |
Helper method to buildApp() responsible for construction of the Shiny web server
.buildServer(input, output, session).buildServer(input, output, session)
input |
reactive values from UI widgets |
output |
render functions for results |
session |
environment for user session |
A Shiny web server
Helper method to buildApp() responsible for construction of the graphical user interface
.buildUI().buildUI()
A Shiny UI
A thin R wrapper for the python3 function fit_curve(), which attempts to optimize 5-parameter logistic regressions to PCR fluorescence curves.
.runFitCurve(pcr_data, linear_threshold).runFitCurve(pcr_data, linear_threshold)
pcr_data |
A
|
linear_threshold |
numeric value specifying the minimum overall change-in-fluorescence over the PCR reaction required for the optimizer to attempt fitting a logistic model |
The model as a nested list
Data Import Server
.serverDataImport(id).serverDataImport(id)
id |
identifier unique to |
The Data Import server component
Curve Fitting Server
.serverFitCurves(id, se).serverFitCurves(id, se)
id |
identifier unique to |
se |
SummarizedExperiment with OpenArray assay data |
The Fit Curves server component
Graphics Server
.serverGraphics(id, se).serverGraphics(id, se)
id |
identifier unique to |
se |
SummarizedExperiment with OpenArray assay data |
The Graphics server component
Reporting Server
.serverReporting(id, se_fit).serverReporting(id, se_fit)
id |
identifier unique to |
se_fit |
SummarizedExperiment with OpenArray assay data |
The Reporting server component
Data Import UI
.uiDataImport(id).uiDataImport(id)
id |
identifier unique to |
The Data Import page
Curve Fitting UI
.uiFitCurves(id).uiFitCurves(id)
id |
identifier unique to |
The Fit Curves page
Graphics UI
.uiGraphics(id).uiGraphics(id)
id |
identifier unique to |
The Graphics page
Reporting UI
.uiReporting(id).uiReporting(id)
id |
identifier unique to |
The Reporting page
Builds a Shiny application for interactively running OAtools
buildApp()buildApp()
Temporarily changes the maximum file size for upload to 30MB, restoring the original setting upon application exit.
A Shiny web application
app <- buildApp()app <- buildApp()
Computes 5-parameter logistic models optimized to PCR data contained in the
specified assay matrix, either fluo_normalized or fluo_reporter. The
former refers to the base-line adjusted, normalized fluorescence from the
Amplification Data tab of the original Excel output. The latter refers to
the multicomponent fluorescence from the Multicomponent Data tab.
The expected input for this function is a SummarizedExperiment object
containing OpenArray PCR data, which can be generated by calling the
excelToSE() function of the package on the raw Excel output from
QuantStudio software.
computeModels(se, assay_name, linear_threshold = 400)computeModels(se, assay_name, linear_threshold = 400)
se |
OpenArray PCR data as a SummarizedExperiment object |
assay_name |
character value specifying the assay matrix from which to
load PCR data for curve-fitting, either |
linear_threshold |
numeric value specifying the minimum overall change-in-fluorescence over the PCR reaction required for the optimizer to attempt fitting a logistic model |
Under the hood, this function invokes .runFitCurve(), a thin wrapper which
calls python3 code to fit models with scipy.optimize, on each PCR reaction
of the OpenArray plate separately. The computed models are added to the
SummarizedExperiment container in metadata, named as assay_name +
_models. For example, fluo_normalized_models would be created in the
experiment metadata if computeModels(assay_name = "fluo_normalized") is
called.
OpenArray PCR data as a SummarizedExperiment object with information from the computed model stored as metadata.
path <- system.file( "extdata", "oa_gene_expression_1.xlsx", package = "OAtools" ) se <- excelToSE(excel_path = path) |> computeModels(assay_name = "fluo_normalized") |> computeModels(assay_name = "fluo_reporter")path <- system.file( "extdata", "oa_gene_expression_1.xlsx", package = "OAtools" ) se <- excelToSE(excel_path = path) |> computeModels(assay_name = "fluo_normalized") |> computeModels(assay_name = "fluo_reporter")
Assigns positive or negative PCR results to each reaction depending on the equation of the model optimized to observed multicomponent fluorescence values. The reporter dye fluorescence without normalization is used to distinguish between real amplification and false positives.
determinePCRResults(se, key_path)determinePCRResults(se, key_path)
se |
a SummarizedExperiment object containing OpenArray qPCR Data |
key_path |
file path to the target-threshold key |
The key is an Excel file storing values used by determinePCRResults() to
categorize PCR curves as positive or negative. In the key, each gene is
associated with threshold values for Crt, change-in-fluorescence, and slope
at the reaction midpoint.
For more specifics, input ?target_threshold_key into the R console.
a SummarizedExperiment
data(example_se) key_path = system.file( "extdata", "target_threshold_key.xlsx", package = "OAtools" ) se <- example_se |> determinePCRResults(key_path = key_path)data(example_se) key_path = system.file( "extdata", "target_threshold_key.xlsx", package = "OAtools" ) se <- example_se |> determinePCRResults(key_path = key_path)
A sample of OpenArray gene expression data from respiratory tract microbiota profiling experiments conducted at the UW Virology research lab on human nasal swabs. This file (.rda) stores the PCR data in a SummarizedExperiment container and is intended for use with package examples and unit testing.
a SummarizedExperiment object with the following:
information associated with each PCR well
cycle numbers
matrices of fluorescence values by PCR well and cycle
Load this object into the environment with: data(example_se)
This object contains PCR data imported from the initial Excel QuantStudio output. Logistic and linear models were fit to the PCR data and stored as metadata.
The object can be reproduced from the package example data by running the following commands in the R console:
path <- system.file(
"extdata",
"oa_gene_expression_1.xlsx",
package = "OAtools"
)
se <- excelToSE(excel_path = path) |>
computeModels(assay_name = "fluo_normalized") |>
computeModels(assay_name = "fluo_reporter")
Transforms raw gene expression run data exported from the OpenArray QuantStudio 12K Flex Software from .xlsx format into an instance of the SummarizedExperiment class from Bioconductor.
excelToSE(excel_path, header_rows = 17, skip = 19)excelToSE(excel_path, header_rows = 17, skip = 19)
excel_path |
file path to the Excel document containing the PCR data |
header_rows |
number of rows of run metadata to read in from the header |
skip |
number of rows to skip when reading fluorescence data or results |
OpenArray PCR data as a SummarizedExperiment object
path = system.file( "extdata", "oa_gene_expression_1.xlsx", package = "OAtools" ) se <- excelToSE(excel_path = path)path = system.file( "extdata", "oa_gene_expression_1.xlsx", package = "OAtools" ) se <- excelToSE(excel_path = path)
Knits an HTML report summarizing the OpenArray experiment and saves to the specified directory.
generateReport(se, path = tempdir(), model_results = FALSE)generateReport(se, path = tempdir(), model_results = FALSE)
se |
a SummarizedExperiment containing OpenArray qPCR data |
path |
intended outfile path, defaults to a temporary directory |
model_results |
boolean value indicating whether to include the
|
An HTML Report summarizing the OpenArray experiment
data(example_se) generateReport(se = example_se)data(example_se) generateReport(se = example_se)
The first of two Excel files storing run data from separate OpenArray gene expression experiments. The context behind the experiment was respiratory tract microbiota profiling on human nasopharyngeal swabs from patients experiencing respiratory syndromes.
An Excel file (.xlsx) with three tabs:
normalized fluorescence by cycle number
spectral contribution of the reporter dye by cycle number
PCR results, PCR well metadata, and QC metrics
Access this file with:
system.file("extdata", "oa_gene_expression_1.xlsx", package = "OAtools")
This file was exported from QuantStudio 12K Flex Software after a gene expression run, then filtered down to include 12 samples and 4 genes for file size concerns.
The second of two Excel files storing run data from separate OpenArray gene expression experiments. The context behind the experiment was respiratory tract microbiota profiling on human nasopharyngeal swabs on patients experiencing respiratory syndromes.
An Excel file (.xlsx) with three tabs:
normalized fluorescence by cycle number
spectral contribution of the reporter dye by cycle number
PCR results, sample metadata, and QC metrics
Access this file with:
system.file("extdata", "oa_gene_expression_2.xlsx", package = "OAtools")
This file was exported from QuantStudio 12K Flex Software after a gene expression run, then filtered down to include 12 samples and 4 genes for file size concerns.
Generates a box and whisker plot visualizing the distribution of Crt values measured on an OpenArray plate by Gene.
plotCrt(se)plotCrt(se)
se |
a SummarizedExperiment object containing OpenArray qPCR data |
a ggplot2 figure
data(example_se) plotCrt(example_se)data(example_se) plotCrt(example_se)
Juxtaposes the fluorescence values predicted by the model optimized to the measured fluorescence vs. cycle data for a particular well.
plotModel( se, well_id, assay_name, include_mdpt_tangent = FALSE, include_coldata_annotation = FALSE )plotModel( se, well_id, assay_name, include_mdpt_tangent = FALSE, include_coldata_annotation = FALSE )
se |
a SummarizedExperiment object containing OpenArray qPCR data |
well_id |
a character representing the name of the well to plot as listed in the assay matrix |
assay_name |
name for the assay matrix from which to pull observed
fluorescence values, either |
include_mdpt_tangent |
boolean determines whether to annotate the midpoint of the reaction and draw a tangent line to the model curve at that point |
include_coldata_annotation |
boolean determines whether to annotate the coldata onto the top left of the plot. |
a ggplot2 figure
data(example_se) plotModel( example_se, well_id = "well_2665", assay_name = "fluo_reporter", include_mdpt_tangent = TRUE, include_coldata_annotation = TRUE )data(example_se) plotModel( example_se, well_id = "well_2665", assay_name = "fluo_reporter", include_mdpt_tangent = TRUE, include_coldata_annotation = TRUE )
Generates an overview graphic summarizing qPCR results for each combination of sample and gene on an OpenArray experiment.
plotOverview(se)plotOverview(se)
se |
a SummarizedExperiment object containing OpenArray qPCR data |
a ggplot2 figure
data(example_se) plotOverview(example_se)data(example_se) plotOverview(example_se)
Generates a 3-dimensional quality control plot comparing the amplification status to the crt, Cq conf, and amplification score metrics output by QuantStudio 12K Flex Software.
plotQC(se)plotQC(se)
se |
a SummarizedExperiment object containing OpenArray qPCR data |
a plotly figure
data(example_se) plotQC(example_se)data(example_se) plotQC(example_se)
Transforms OpenArray run data contained within the SummarizedExperiment container into a qPCRBatch object. This conversion allows for convenient gene expression analyses with the NormqPCR package.
seToQPCRBatch(se)seToQPCRBatch(se)
se |
A SummarizedExperiment object with OpenArray run data |
a qPCRBatch object
path = system.file( "extdata", "oa_gene_expression_1.xlsx", package = "OAtools" ) se <- excelToSE(excel_path = path) qpcr <- seToQPCRBatch(se)path = system.file( "extdata", "oa_gene_expression_1.xlsx", package = "OAtools" ) se <- excelToSE(excel_path = path) qpcr <- seToQPCRBatch(se)
This Excel file (.xlsx) contains a table associating each assay target contained in the example gene expression run data with thresholds values. This key is optionally used to aid in concurrently interpreting data including numerous targets with dissimilar behaviors.
An Excel sheet (.xlsx) with four columns:
targetthe target assay ID
slope_thresholdminimum acceptable slope in the exponential phase for a positive result
delta_thresholdminimum acceptable overall change in fluorescence for a positive result
crt_thresholdmaximum acceptable crt for a positive result
The example key is stored in the inst/extdata/ directory of the package.
Access this file with:
system.file("extdata", "target_threshold_key.xlsx", package = "OAtools").
Written manually for use with resulting functions.