2.1的QcMetric类

质量控制度量由以下描述组成(的名字在下面的代码块)，一些QC数据(qcdata)和a状态这定义了度量是否被认为具有可接受的质量(编码为真正的)，质量差(编码为假)或尚未计算(编码为NA)．单个指标可以显示为简短的文本摘要或绘制。要实现前者，可以使用默认值显示方法。

library("qcmetrics") qc <- QcMetric(name = "A test metric") qcdata(qc， "x") <- rnorm(100) qcdata(qc) ##所有可用的qcdata

##[1]“x”

Summary (qcdata(qc， "x")) ##获取x

最小第一曲，中位数，平均第三曲，最大值。## -2.2147 -0.4942 0.1139 0.1089 0.6915 2.4016

显示(qc) ##或只显示qc

## QcMetric类的对象##名称:测试度量##状态:NA ##数据:x

status(qc) <- TRUE qc

## QcMetric类的对象##名称:测试度量##状态:TRUE ##数据:x

策划QcMetric实例需要实现与手头数据相关的绘图方法。我们可以用情节替换方法来定义自定义函数。中的代码情节使用qcdata提取相关QC数据对象然后作为参数传递给情节并使用足够的可视化来呈现QC数据。

情节(qc)

## x@plot(x，…)中的警告:没有定义特定的绘图函数

Plot (qc) <- function(对象，…)Boxplot (qcdata(对象，"x")，…)

2．2的QcMetrics类

一个QcMetrics对象本质上只是单个实例的列表。也可以设置一个元数据变量列表来描述QC度量的来源。元数据可以作为QcMetadata对象的存储方式QcMetrics实例)或直接作为命名的列表．的QcMetadata本身就是列表并可以访问和设置元数据或mdata．访问时，它将返回并显示为列表．

qcm <- QcMetrics(qcdata = list(qc)

类“QcMetrics”的对象，包含1个QC指标。##和没有元数据变量。

元数据(qcm) <- list(author = " professor - Who"， lab = "Big lab") qcm

类“QcMetrics”的对象，包含1个QC指标。##和两个元数据变量。

mdata(药物)

# # $ # #作者[1]“教授”# # # # # # $实验室[1]“大实验室”

元数据可以在同一接口下更新。如果传递了新的命名项，则通过添加新元素来更新元数据。如果已存在命名项，则更新其值。

metadata(qcm) <- list(author = " professor - Who"， lab = "Cabin lab"， University = " University -ity") mdata(qcm)

# # $ # #作者[1]“教授”# # # # # # $实验室[1]“机舱实验室”# # # # # #大学美元[1]“Universe-ity”

的QcMetrics然后可以传递到qcReport方法生成报告，如下面的详细信息所述。

3．1微阵列退化

Microarray降级部分已被删除，因为它所依赖的软件包已弃用。

3.2蛋白质组学原始数据

为了说明蛋白质组学数据的简单QC分析，我们将下载数据集PXD00001mzXML格式的ProteomeXchange存储库(Pedrioli等2004)．从质谱运行的MS2光谱然后读入Rand存储为MSnExp使用readMSData函数从MSnbase包(Gatto and Lilley 2012)．

library("RforProteomics") msfile <- getPXD000001mzXML() library("MSnbase") exp <- readMSData(msfile, verbose = FALSE)

为了节省时间，这个代码块已经被预先计算，并且是子集(1 / 3)经验值实例随包一起分发。加载数据

负载(执行(“extdata /经验。Rda "， package = "qcmetrics"))

的QcMetrics将包括3个项目，即用MS2光谱前体强度构建的色谱图、说明MS空间中前体电荷的图和说明MS2光谱用于识别的适用性的m/z δ图(见plotMzDelta ?或(Foster et al. 2011))．

qc1 <- QcMetric(name ="色谱图")x <- rtime(exp) y <- precursorIntensity(exp) o <- order(x) qcdata(qc1， "x") <- x[o] qcdata(qc1， "y") <- y[o] plot(qcdata(object， "x")， qcdata(object， "y")， col ="深灰色"，type ="l"， xlab ="保留时间"，ylab ="前体强度")

qc2 <- QcMetric(name = "MS space") qcdata(qc2， "p2d") <- plot2d(exp, z = "charge"， plot = FALSE) plot(qc2) <- function(object) {require("ggplot2") print(qcdata(object， "p2d"))}

qc3 <- QcMetric(name = "m/z delta plot") qcdata(qc3， "pmz") <- plotMzDelta(exp, plot = FALSE, verbose = FALSE) plot(qc3) <- function(对象)suppressWarnings(print(qcdata(对象，"pmz")))

请注意，我们不会将原始数据存储在上述任何实例中，而是始终预先计算必要的数据或图，然后将其存储为qcdata．如果原始数据需要在多个QcMetric实例，我们可以重用相同的实例qcdata环境避免不必要的拷贝使用Qcdata (qc2) <- qcenv(qc1)并通过自定义实现不同的视图情节方法。

现在让我们把这三项合并成aQcMetrics对象中的MIAPE信息，用自定义元数据装饰它MSnExp对象并生成报告。

protqcm <- QcMetrics(qcdata = list(qc1, qc2, qc3)) metadata(protqcm) <- list(data = "PXD000001"， instrument = experimentData(exp)@instrumentModel, source = experimentData(exp)@ionSource, analyser = experimentData(exp)@analyser, detector = experimentData(exp)@detectorType, manufacurer = experimentData(exp)@instrumentManufacturer)

汇总表的状态列为空，因为我们还没有设置QC项的状态。

qcReport(protqcm, reportname = "protqc")

完整的pdf报告可在以下网址下载:

browseURL (example_reports(“protqc”))

3．3处理N15标签数据

在本节中，我们将描述一组N15代谢标记QC指标(Krijgsveld et al. 2003)．数据是磷富集N15标记拟南芥样品按(格林等，2013)．数据是用内部工具处理的，可作为MSnSet实例。简单地说，用Mascot引擎搜索MS2光谱，用Mascot Percolator调整识别分数。然后在调查扫描中搜索重和轻对，并根据肽序列和重对成员的同位素包膜来估计N15的结合公司特征变量)。最后整合重肽和轻肽同位素包络区，得到未标记和N15定量数据。的psmobject为psm(肽谱匹配)提供后验误差概率< 0.05的数据，可以唯一地与蛋白质匹配。

我们首先载入MSnbase软件包(支持MSnSet的数据结构)和示例数据qcmetrics包中。我们将利用ggplot2策划方案。

library("ggplot2") library("MSnbase") data(n15psm) psm

## MSnSet (storageMode: lockedEnvironment) ## assayData: 1772个特征，2个样本##元素名称:exprs ##协议数据:无##表型数据:无## featureData ## featureNames: 3 5…4499 (1772 total) ## fvarLabels: Protein_Accession Protein_Description…inc (21 total) ## fvarMetadata: labelDescription ## experimentData: use 'experimentData(object)' ## pubMedIds: 23681576 ##注释:## - - - -处理信息- - - - ##子集[22540,2][1999,2]Tue Sep 17 01:34:09 2013 ##删除超过0个NAs的特性:Tue Sep 17 01:34:09 2013 ##删除featureData的级别Tue Sep 17 01:34:09 2013 ## MSnbase版本:1.9.7

第一个QC项目检查N15掺入率，可在公司特征变量。我们还定义了一个中值合并率阈值tr等于97.5，用于设置QC状态。

##纳入率QC度量qcinc <- QcMetric(name = "15N纳入率")qcdata(qcinc， "inc") <- fData(psm)$inc qcdata(qcinc， "tr") <- 97.5 status(qcinc) <- median(qcdata(qcinc， "inc")) > qcdata(qcinc， "tr")

接下来，我们实现一个自定义显示方法，该方法打印变量分布的5个汇总值。

show(qcinc) <- function(object) {qcshow(object, qcdata = FALSE) cat(" QC阈值:"，qcdata(object， "tr")， "\n") cat("并入率\n") print(summary(qcdata(object， "inc")) invisible(NULL)}

然后我们定义度规情节函数，表示PSM的合并率的分布为箱线图，显示所有的个别率为抖动点，并表示tr阈值用红色虚线表示。

plot(qcinc) <- function(object) {inc <- qcdata(object， "inc") tr <- qcdata(object， "tr") lab <- "纳入率" dd <- data.frame(inc = qcdata(qcinc， "inc")) p <- ggplot(dd, aes(factor("" ")， inc))) + geom_jitter(color = "#4582B370"， size = 3) + geom_boxplot(fill = "#FFFFFFD0"， color = "#000000"，异常值。size = 0) + geom_hline(yintercept = tr, color = "red"， linetype = "虚线"，size = 1) +实验室(x = ""， y = "掺入率")p}

高质量的N15实验的特点是掺入率高，这允许反褶积重和轻肽同位素包络和准确定量。

第二个指标检查psm的log2折叠变化，独特的修饰肽，独特的肽序列(不考虑修饰)和蛋白质。计算这些各自的数据集combineFeatures函数(见combineFeatures ?详情)。

fData(psm)$modseq <- ## pep seq + PTM paste(fData(psm)$Peptide_Sequence, fData(psm)$ variable_modification, sep = "+") pep <- combineFeatures(psm, as.character(fData(psm)$Peptide_Sequence)， "median"， verbose = FALSE) modpep <- combineFeatures(psm, fData(psm)$modseq， "median"， verbose = FALSE) prot <- combineFeatures(psm, as.character(fData(psm)$Protein_Accession)， "median"， verbose = FALSE)

然后计算所有特征的log2折叠变化，并存储为我们下一个QC项目的QC数据。我们还存储了一对值explfc它定义了一个区间，在这个区间中，我们期望PSM log2的折叠变化中值。

##计算log fold-change qclfc <- QcMetric(name = "Log2 fold-changes") qcdata(qclfc， "lfc.psm") <- Log2 (exprs(psm)[，"unlabelled"] / exprs(psm)[，" N15"]) qcdata(qclfc， "lfc.pep") <- Log2 (exprs(modpep)[，"unlabelled"] / exprs(modpep)[，" N15"]) qcdata(qclfc， "lfc.modpep") <- Log2 (exprs(modpep)[，"unlabelled"] / exprs(modpep)[，" N15"]) qcdata(qclfc， "lfc.prot") <- Log2 (exprs(prot)[，"unlabelled"] / exprs(modpep)[，" N15"]) qcdata(qclfc， "lfc.prot") <- Log2 (exprs(prot)[，"unlabelled"] / exprs(prot)[，" N15"]) qcdata(qclfc， "explfc .prot") <- c(-0.5, 0.5) status(qclfc) <- median(qcdata(qclfc， "explfc .prot")) qcdata(qclfc， "explfc .prot") <- c(-0.5, 0.5) status(qclfc) <- median(qcdata(qclfc， ""lfc.psm")) > qcdata(qclfc， "explfc")[1] & median(qcdata(qclfc， "lfc.psm")) < qcdata(qclfc， "explfc")[2]

与前面一样，我们提供了一个自定义显示方法，该方法显示四个折叠更改的汇总值。的情节函数说明了各自的log2折叠变化密度和预期中位数PSM折叠变化范围(红色矩形)。预期的0 log2折叠变化显示为黑色虚线，观测到的PSM中值显示为蓝色虚线。

show(qclfc) <- function(object) {qcshow(object, qcdata = FALSE) ## default cat(" QC阈值:"，qcdata(object， "explfc")， "\n") cat(" * PSM log2 fold-changes\n") print(summary(qcdata(object， "lfc.psm")) cat(" * Modified peptide log2 fold-changes\n") print(summary(qcdata(object， "lfc.modpep")) cat(" * peptide log2 fold-changes\n") print(summary(qcdata(object， "lfc.pep")) cat(" * Protein log2 fold-changes\n") print(summary(qcdata(object， "lfc.pep")) cat(" * Protein log2 fold-changes\n") print(summary(qcdata(object， "lfc.pep")) cat(" * Protein log2 fold-changes\n") print(summary(qcdata(object， "lfc.pep"))“lfc.prot”)))无形的情节(NULL)} (qclfc) < -函数(对象){x < - qcdata(对象,“explfc”)情节(密度(qcdata(对象,“lfc.psm”),主要= ",子= " ",坳=“红色”,ylab = " ", lwd = 2, xlab =表达式(日志[2]~叠化))行(密度(qcdata(对象,“lfc.modpep”)),坳=“钢蓝色的”,lwd = 2)线(密度(qcdata(对象,“lfc.pep”)),坳=“蓝色”,lwd = 2)线(密度(qcdata(对象,“lfc.prot”)),坳=“橙色”)abline (h = 0,坳=“灰色”)abline (v = 0, lty =“点缀”)矩形(x [1], 1,x[2]， 1, col = "#EE000030"， border = NA) abline(v = median(qcdata(object， "lfc.psm"))， lty = "破折号"，col = "蓝色")legend("topright"， c("PSM"， "多肽"，"修饰肽"，"蛋白质")，col = c("红色"，"钢蓝"，"蓝色"，"橙色")，lwd = 2, bty = "n")}

一个高质量的实验应该有一个以0为中心的紧密分布。较大的偏差将表明未完全纳入，使用的轻材料和重材料各自数量的误差，广泛的分布将反映数据的巨大可变性。

我们的最后一个QC项目检查了在实验中已经确定的特征的数量。我们还研究了在PSM水平上观察到多少肽(有或没有考虑修饰)以及每个蛋白质的独特肽的数量。在这里，我们不指定任何期望值，因为观察到的特征的数量是实验特定的;QC状态保持为NA．

##功能数量qcnb <- QcMetric(name = "功能数量")qcdata(qcnb， "count") <- c(PSM = nrow(PSM)， ModPep = nrow(ModPep)， Pep = nrow(Pep)， Prot = nrow(Prot)) qcdata(qcnb， "peptab") <- table(fData(PSM)$Peptide_Sequence) qcdata(qcnb， "modpeptab") <- table(fData(PSM)$modseq) qcdata(qcnb， " upppe .per. Prot ") <- fData(PSM)$Number_Of_Unique_Peptides

计数按new显示显示并被绘制成柱状图情节方法。

显示(qcnb) < -函数(对象){qcshow(对象,qcdata = FALSE)打印(qcdata(对象,“数”))}情节(qcnb) < -函数(对象){par (mar = c(5、4、2、1)布局(矩阵(c(1、2、1,3,4),ncol = 3)) barplot (qcdata(对象,“数”),水平的= TRUE,拉斯维加斯= 2)barplot(表(qcdata(对象,“modpeptab”)),xlab =“改良肽”)barplot(表(qcdata(对象,“peptab”)),xlab =“肽”)barplot(表(qcdata(对象,“upep.per.prot”)),xlab =“独特的肽/蛋白”)}

在下面的代码块中，我们将3个QC项合并为一个QcMetrics类中提取的元数据生成报告psmMSnSet实例。

(qcdata = list(qcinc, qclfc, qcnb))

qcReport(n15qcm, reportname = "n15qcreport"， title = expinfo(experimentData(psm))["title"]， author = expinfo(experimentData(psm))["contact"]， clean = FALSE)

一旦确定了一组适当的质量度量，就可以生成QcMetrics可以对实例进行包装以实现自动化。

我们为这个示例提供了这样一个包装器函数n15qc功能完全自动化上述管道。特征变量列的名称和前两个QC项的阈值作为参数提供。如果没有给出报告名称，则使用带有日期和时间的自定义标题，以避免覆盖现有的报告。

完整的pdf报告可通过

browseURL (example_reports(“n15qc”))

4．1自定义报告

4.2新的报表类型

报告功能是一种功能

• 转换适当的QC项目部分(例如Qc2Tex2上面描述的函数)。

• 可选地将QC项目部分包括到附加页眉和页脚中，可以直接写入这些部分，也可以将这些部分插入到适当的模板中。中提供的报告功能qcmetrics可在qcReport ?：reporting_tex对于tex类型，报告\ _pdf型pdf这些函数应该使用相同的参数qcReport尽可能的。

• 一旦写入到报表源文件，最终的报表类型就生成了。针织是用来将Rnw源转换为tex编译成pdf使用工具::texi2pdf．Rmd内容直接写入一个文件，该文件编织并转换为html使用knit2html(这叫markdownTOHTML)．

新reporting_abc函数可以直接调用或传递给qcReport使用记者论点。

4.３质量控制方案

而在章节中给出的例子3.采用灵活快速的方法设计QC管道原型，生产管道需要更稳健的机制。R打包机制非常适合这种情况，因为它提供了版本控制、文档、单元测试以及易于分发和安装的工具。

虽然包开发的详细描述超出了本文档的范围，但是提供一个QC包开发的概述是有意义的。使用包装器函数，可以使用它创建包结构

包中。skeleton("N15QC"， list = "N15QC")

的描述文件将需要更新。的包qcmetrics,MSnbas将需要指定为依赖项进口:中导入的名称空间文件。文档文件N15QC /人/ N15QC。理查德·道金斯(可选的)将需要更新。问