1简介

的minfi包提供了分析Illumina甲基化阵列的工具，特别是450k和EPIC(也称为850k)阵列。我们对旧的27k数组有部分支持。

在这个包中处理的任务包括预处理，QC评估，感兴趣的甲基化位点的识别和绘图功能。分析这些类型的阵列是我们和其他小组正在进行的研究。

这个包的输入数据是IDAT文件，在标准化之前表示两个不同的颜色通道。这是最完整的数据类型，因为IDAT文件包括控制探针上的测量值。可以将Genome Studio文件与此包中包含的数据结构一起使用，但只有一些功能可用，因为Genome Studio输出不包含控制探针信息。此外，通常Genome Studio的输出是使用在Genome Studio中实现的方法进行标准化的，这些方法通常被认为是较差的。

1．1引用minfi套餐

MINFI包包含由不同的非重叠作者在多个手稿中描述的方法。这使得引用包有点困难，所以这里有一个指南。

  • 如果你在出版物中使用MINFI，请注明出处(Aryee et al. 2014)．本刊物包括使用性别估计的详细资料getSex ()，质量控制使用getQC ()，分位数归一化使用preprocessQuantile ()，使用bumphunter ()块查找使用blockFinder ()．
  • 如果您正在使用MINFI来分析EPIC或27k阵列，请引用(Fortin, Triche Jr.， and Hansen 2017)．该出版物包括关于convertArray ()而且combineArrays ()，扩展NOOB工作在单样本模式以及使用estimateCellCounts ()使用来自450k数组的参考数据来估计EPIC数据的单元类型组成。
  • 如果你正在使用preprocessQuantile ()，这将是体贴的也引用(Touleimat and Tost 2012)，因为本文描述的方法本质上与preprocessQuantile ()．
  • 如果你正在使用bumphunter ()，你也应该引用原始的bump hunter出版物(Jaffe et al. 2012)．
  • 如果您正在使用中实现的SWAN规范化preprocessSWAN ()请引用(Maksimovic, Gordon, and Oshlack 2012)．
  • 如果您正在使用noob背景校正，如在preprocessNoob ()，请引用(Triche et al. 2013)．
  • 如果您正在使用中实现的函数规范化preprocessFunnorm ()，请引用(Fortin et al. 2014)．中默认的preprocessFunnorm ()就是做新手背景校正。如使用，请同时注明(Triche et al. 2013)．
  • 如果您正在评估A/B区，如在隔间()而且extractAB ()，请引用(Fortin and Hansen 2015)．

如果你使用Bibtex，你可以得到引用在这种格式通过

toBibtex(引用(“minfi”))

1．2术语

文献往往有点不具体。DNA甲基化微阵列的术语。

对于450k微阵列，每个样品都在一个阵列上测量，在两个不同的颜色通道(红色和绿色)上。每个阵列测量大约45万个CpG位置。每个CpG都与两个测量相关联:甲基化测量和“非”甲基化测量。这两个值可以用两种方法之一来测量:使用“类型I”设计或“类型II设计”。使用I型设计测量的cpg使用单一颜色测量，在同一颜色通道中使用两个不同的探针提供甲基化和非甲基化测量。使用II型设计测量的cpg使用单个探针测量，两种不同的颜色提供甲基化和非甲基化的测量值。实际上，这意味着在这个数组上有不探针和CpG位置之间的一对一对应关系。因此，我们试图精确地描述这一点，当我们提到单碱基基因组位点时，我们指的是“位点”(或“CpG”)，我们将其与“探针”区分开来。上一代27k甲基化阵列仅使用Type I设计，而EPIC阵列同时使用Type I和Type II。

样本之间DNA甲基化的差异既可以是在单个CpG上，称为差异甲基化位置(DMP)，也可以是在区域水平上，称为差异甲基化区域(DMR)。

物理上，每个样本都是在一个“阵列”上测量的。对于450k的设计，在一个物理“幻灯片”上有12个数组(以6 × 2的网格组织)。幻灯片被组织成最多包含8个幻灯片(96个数组)的“板”。EPIC阵列每个幻灯片有8个阵列，每个板有64个阵列。

1.3依赖关系

该文档具有以下依赖项

库(minfi)库(minfiData)

2minfi类

MINFI包被设计为非常灵活的方法开发人员。bob电竞体育官网这种灵活性对用户来说是有代价的;他们需要了解一些不同的数据类:

  • RGChannelSet: IDAT文件的原始数据;这些数据是在探针(而不是CpG位点)水平上组织的。该数据有两个通道:红色和绿色。
  • MethylSet:由CpG位点水平组织的数据，但未映射到基因组。该数据有两个通道:Meth(甲基化)和Unmeth(非甲基化)。
  • RatioSet:由CpG位点水平组织的数据，但未映射到基因组。数据至少具有两个通道中的一个:Beta和/或M (Beta logratio)。它可以选择包含CN通道(拷贝号)。
  • GenomicMethylSet:像MethylSet，但映射到一个基因组。
  • GenomicRatioSet:像RatioSet，但映射到基因组。

类层次结构如下所示

minfi的等级结构

为了更清楚地说明这一点，让我们看一下示例数据RGsetEx从minfiData包中。这是我们在类层次结构中移动时的特性和类的数量(代码未运行):

RGsetEx
##类:RGChannelSet ## dim: 622399 6 ##元数据(0):## assays(2):绿红## rownames(622399): 10600313 10600322…74810490 74810492 ## rowData names(0): ## colnames(6): 5723646052_R02C02 5723646052_R04C01…5723646053_R05C02 ## 5723646053_R06C02 ## colData names(13): Sample_Name Sample_Well…注释:IlluminaHumanMethylation450k注释:ilmn12.hg19
## RGsetEx: RGChannelSet, 622,399 features MsetEx <- preprocessRaw(RGsetEx) ## MsetEx: MethylSet, 485,512 features GMsetEx <- mapToGenome(MsetEx) ## GMsetEx: GenomicMethylSet, 485,512 features

注意特性的数量是如何变化的。在RGChannelSet一个特征是一个探针(它不同于CpG轨迹，请参阅术语部分)。在MethylSet每个特征现在都是一个甲基化位点，它的特征更少，因为一些位点是用两个探针测量的。最后,GenomicMethylSet尺寸和MethylSet，但原则上它可以更小，以防你使用的注释说一些探针不能映射到基因组(在这种情况下是hg19)。

最后我们可以转换成aRatioSet通过ratioConvert ()．这两个函数ratioConvert ()而且mapToGenome ()如上面的类层次结构图所示。许多预处理函数在图中经历了几个步骤，例如，如果函数需要知道探针的基因组位置(一些预处理函数以不同的方式处理在性染色体上测量的探针)。

3.读取数据

此包支持对IDAT文件的分析，其中包含汇总的珠子信息。

根据我们的经验，大多数实验室使用“样本表”CSV文件来描述实验的布局。这是基于Illumina提供的示例表文件。我们的管道假设存在这样一个文件，但是如果它不可用的话，使用例如Excel创建这样一个文件相对容易。

我们使用了一个包含6个样本的示例数据集，分布在两张幻灯片上。首先，我们获得IDAT文件的系统路径;这需要一点，因为数据来自已安装的包

baseDir <- system。file("extdata"， package = "minfiData") list.files(baseDir)
## [1] "5723646052" "5723646053" "SampleSheet.csv"

这显示了450k数据的典型布局:每个“幻灯片”(包含12个数组，参见术语)存储在一个单独的目录中，具有一个数字名称。顶层目录包含示例表文件。在幻灯片目录中，我们可以找到IDAT文件(可能还有一些JPG图像或其他文件):

list.files(文件。路径(baseDir,“5723646052”))
## [1] "5723646052_R02C02_Grn。idat”“5723646052 _r02c02_red。idat“##[3]”5723646052_R04C01_Grn。idat”“5723646052 _r04c01_red。idat“##[5]”5723646052_R05C02_Grn。idat 5723646052 _r05c02_red.idat”

每个数组的文件都有另一个数字，由一个红色和一个Grn(绿色)IDAT文件组成。注意，对于这个示例数据，每张幻灯片只包含3个数组，而不是12个。这样做是因为文件大小的限制，因为我们只需要6个数组来说明包的功能。

首先我们读一下样品单。我们提供了一个方便的函数来读取这个文件read.metharray.sheet ()．这个函数有几个吸引人的附加功能。让我们看看输出

目标<- read.metharray.sheet(baseDir)
# # (read.metharray。找到以下CSV文件:
## [1] "/home/biocbuild/bb -3.11- bio/ R/library/minfiData/extdata/SampleSheet.csv"
目标
# # Sample_Name Sample_Well Sample_Plate Sample_Group Pool_ID人年龄性别# # 1 GroupA_3 H5 < NA > GroupA < NA > id3 83 # # 2 GroupA_2 D5 < NA > GroupA F < NA > id2 58 # # 3 GroupB_3 C6 < NA > GroupB < NA > id3 83 # # 4 GroupB_1 F7 < NA > GroupB < NA > id1 75 F # # 5 GroupA_1七国集团(G7) < NA > GroupA < NA > id1 75 F # # 6 GroupB_2 H7 < NA > GroupB F < NA > id2 58 # #状态数组幻灯片# # 1正常R02C02 5723646052 # # 2正常R04C01 5723646052 # # 3癌症R05C02 5723646052 # # 4正常R05C02 R04C02 5723646053 # # 51 /home/biocbuild/ bbbs -3.11-bioc/R/library/minfiData/extdata/5723646052/5723646052_R02C02 ## 2/ home/biocbuild/ bbbs -3.11-bioc/R/library/minfiData/extdata/5723646052/5723646052_R04C01 ## 3/ home/biocbuild/ bbbs -3.11-bioc/R/library/minfiData/extdata/5723646052/5723646052_R05C02 ## 4 /home/biocbuild/ bbbs -3.11-bioc/R/library/minfiData/extdata/5723646052/5723646052 _r04c02 ##6 /home/biocbuild/论坛-3.11-bioc/R/library/minfiData/extdata/5723646053/5723646053_R06C02 ## 6 /home/biocbuild/论坛-3.11-bioc/R/library/minfiData/extdata/5723646053/5723646053_R06C02

首先输出:这是adata.frame．它包含一列Basename描述了对应于示例的IDAT文件的位置，以及两列数组而且幻灯片．在Illumina提供的样例表中，命名了这两列Sentrix_Position而且Sentrix_ID，但我们重新命名它们。下面将提供关于此函数使用的更多细节。的Basename列往往太大，不便于显示，这里是相对简化baseDir：

sub(baseDir， ""，目标$Basename)
## [1] "/ 5723646052/5723646052_r02c02 " "/ 5723646052/5723646052_r04c01 " ## [3] "/ 5723646052/5723646052_r05c02 " " ## [5] "/ 5723646053/5723646053_r05c02 "

(这只是为了展示)。

用这个data.frame，从数据中很容易读出

RGset <- read.metharray。Exp(目标=目标)
## readChar警告(con, nchars = n):嵌入空的截断字符串## readChar警告(con, nchars = n):嵌入空的截断字符串## readChar警告(con, nchars = n):嵌入空的截断字符串## readChar警告(con, nchars = n):嵌入空的截断字符串## readChar(con, nchars = n):嵌入空的截断字符串## readChar警告(con, nchars = n):嵌入空的截断字符串## readChar(con, nchars = n):嵌入空的截断字符串#### readChar警告(con, nchars = n):嵌入空的截断字符串## readChar警告(con, nchars = n):嵌入空的截断字符串## readChar警告(con, nchars = n):嵌入空的截断字符串## readChar警告(con, nchars = n):嵌入空的截断字符串## readChar(con, nchars = n):嵌入空的截断字符串## readChar警告(con, nchars = n):嵌入空的截断字符串## readChar(con, nchars = n):嵌入空的截断字符串##

让我们看看相关的pheno数据，这实际上只是上面的targets对象中包含的信息。

RGset
##类:RGChannelSet ## dim: 622399 6 ##元数据(0):## assays(2):绿红## rownames(622399): 10600313 10600322…74810490 74810492 ## rowData names(0): ## colnames(6): 5723646052_R02C02 5723646052_R04C01…5723646053_R05C02 ## 5723646053_R06C02 ## colData names(13): Sample_Name Sample_Well…注释:IlluminaHumanMethylation450k注释:ilmn12.hg19
pd <- pData(RGset) pd[，1:4]
##数据帧，6行4列## Sample_Name Sample_Well Sample_Plate Sample_Group ## <字符> <字符> <字符> <字符> ## 5723646052_R04C01 GroupA_2 H5 NA groupa# # 5723646052_R04C01 GroupA_2 D5 NA groupa# # 5723646052_r04c02 GroupB_3 C6 NA groupb# # 5723646053_R04C02 GroupB_1 F7 NA groupb# # 5723646053_R06C02 GroupB_2 H7 NA groupb# # 5723646053_R06C02 GroupB_2 H7 NA groupb# # 5723646053_R06C02 GroupB_2 H7 NA groupb# # 5723646053_R06C02 GroupB_2

的read.metharray.exp ()函数还可以读取整个目录或目录树(使用递归设置为真正的)通过只使用函数的实参基地而且目标=零,就像

RGset2 <- read.metharray.exp(文件。路径(baseDir,“5723646052”))
## readChar警告(con, nchars = n):嵌入空的截断字符串## readChar警告(con, nchars = n):嵌入空的截断字符串## readChar警告(con, nchars = n):嵌入空的截断字符串## readChar警告(con, nchars = n):嵌入空的截断字符串## readChar(con, nchars = n):嵌入空的截断字符串## readChar警告(con, nchars = n):嵌入空的截断字符串## readChar(con, nchars = n):嵌入空的截断字符串##
RGset3 <- read.metharray。exp(baseDir, recursive = TRUE)
## readChar警告(con, nchars = n):嵌入空的截断字符串## readChar警告(con, nchars = n):嵌入空的截断字符串## readChar警告(con, nchars = n):嵌入空的截断字符串## readChar警告(con, nchars = n):嵌入空的截断字符串## readChar(con, nchars = n):嵌入空的截断字符串## readChar警告(con, nchars = n):嵌入空的截断字符串## readChar(con, nchars = n):嵌入空的截断字符串#### readChar警告(con, nchars = n):嵌入空的截断字符串## readChar警告(con, nchars = n):嵌入空的截断字符串## readChar警告(con, nchars = n):嵌入空的截断字符串## readChar警告(con, nchars = n):嵌入空的截断字符串## readChar(con, nchars = n):嵌入空的截断字符串## readChar警告(con, nchars = n):嵌入空的截断字符串## readChar(con, nchars = n):嵌入空的截断字符串##

3．1读取数据的高级说明

表中唯一重要的列data.frame用于目标的参数read.metharray.exp ()函数是列名Basename．通常，这样的对象也有命名的列数组，幻灯片，和(可选)板．

我们使用了基于Illumina提供的样例表数据文件构建的表数据文件。这是一个带有头的CSV文件。在这种情况下，我们假设表型数据从第一行开始(数据)(或者没有头文件)。

使用该函数手动读取示例表也很容易read.csv ()．在这里，我们知道要跳过文件的前7行。

Targets2 <- read.csv(file. csv)targets2 . path(baseDir， "SampleSheet.csv")， stringsAsFactors = FALSE, skip = 7
# # Sample_Name Sample_Well Sample_Plate Sample_Group Pool_ID Sentrix_ID # # 1 GroupA_3 H5 NA GroupA NA 5723646052 # # 2 GroupA_2 D5 NA GroupA NA 5723646052 # # 3 GroupB_3 C6 NA GroupB NA 5723646052 # # 4 GroupB_1 F7 NA GroupB NA 5723646052 # # 5 GroupA_1 G7 NA GroupA NA 5723646053 # # 6 GroupB_2 H7 NA GroupB NA 5723646053 # # Sentrix_Position人年龄性别地位# # 1正常R02C02 id3 83 # # 2 R04C01 id2 58 # # 3 F正常R05C02 id3 83癌症# # 4 R04C02 id1 75 F # # 5 R05C02 id1 75 F正常R06C02 id2 58 F癌症

我们现在需要填充Basename列。关于可能的方法如下

targets2$Basename <- file。path(baseDir, targets2$Sentrix_ID, paste0(targets2$Sentrix_ID, targets2$Sentrix_Position))

最后，MINFI包含一个基于文件的解析器:read.metharray ()．返回对象表示样本的红色和绿色通道测量值。我们从包中得到的一个有用的函数Biobase是结合()这将组合(“添加”)两组样本。这允许用户手动建立一个RGChannelSet．

4清单/注释

4．1每个人都需要知道的

对于甲基化阵列，我们有两种类型的注释包:包含阵列设计的“manifest”包和包含甲基化位点在基因组上的位置、它们映射到哪些基因组特征以及它们是否可能重叠任何已知snp的“annotation”包。

您可以看到哪些包正在被使用

注释(RGsetEx)
## IlluminaHumanMethylation450k“ilmn12.hg19”

4.2先进的讨论

本讨论面向希望了解MINFI内部原理的包开发人员或用户。bob电竞体育官网

一组450k的数据文件最初将被读入RGChannelSet，表示原始强度为两个矩阵:一个是绿色通道，一个是红色通道。这是一个非常类似于ExpressionSet或者一个NChannelSet．的RGChannelSet是，在一起了吗IlluminaMethylationManifest对象，预处理为MethylSet．的IlluminaMethylationManifest对象包含阵列设计，并描述探针和颜色通道如何配对在一起以测量特定CpG的甲基化水平。对象还包含关于控制探测(也称为QC探测)的信息。的MethylSet包含标准化数据，基本上由两个矩阵组成，其中包含每个CpG的甲基化和非甲基化证据。只有RGChannelSet包含关于控制探测的信息。

方法分析Affymetrix表达式数组的范例与前一段中描述的过程非常相似affy包(一个AffyBatch被预处理成ExpressionSet使用存储在CDF环境(包)中的数组设计信息)。

5质量控制

  • 使用shinyMethyl
  • minfiQC
  • getSex检查
  • mds情节
  • 地块转换探头

6预处理

MINFI中的预处理是由一系列名为preprocessXX．不同的函数有不同的输出类。

目前，我们有

  • preprocessRaw:不处理。
  • preprocessIllumina: Illumina预处理，由Genome Studio执行(由我们进行逆向工程)。
  • preprocessSWAN: SWAN归一化，在(Maksimovic, Gordon, and Oshlack 2012)．
  • preprocessQuantile:分位数归一化(适用于DNA甲基化阵列)，在(Touleimat and Tost 2012，阿尔耶等人(2014)）
  • preprocessNoob: Noob预处理，详见(Triche et al. 2013)．
  • preprocessFunnorm:函数归一化，如(Fortin et al. 2014)．

补救办法:讨论文学

7snp和其他问题

  • 单核苷酸多态性
  • 交叉反应探针
  • 移除坏探针
  • 检测P
  • 差距狩猎

8鉴别不同甲基化区域

  • DMP发现
  • 撞狩猎
  • 块发现
  • 带注释的DMRs图。

9校正细胞类型异质性

  • estimateCellTypes
  • 参考包

10其他东西

  • 霍瓦特年龄估计
  • 其他包

11Sessioninfo

参考文献

Aryee, Martin J, Andrew E Jaffe, Hector Corrada Bravo, Christine Ladd-Acosta, Andrew P Feinberg, Kasper D Hansen和Rafael A Irizarry. 2014。“Minfi:用于Infinium DNA甲基化微阵列分析的灵活而全面的Bioconductor包。”生物信息学30(10): 1363 - 9。https://doi.org/10.1093/bioinformatics/btu049．

让-菲利普·福廷，卡斯珀·D·汉森，2015。“利用表观遗传数据中的长程相关性重建Hi-C揭示的A/B区室。”基因组生物学16:180。https://doi.org/10.1186/s13059-015-0741-y．

Fortin, Jean-Philippe, Aurélie Labbe, Mathieu Lemire, Brent W Zanke, Thomas J Hudson, Elana J Frtig, Celia MT Greenwood, Kasper D Hansen. 2014。“450k甲基化阵列数据的功能规范化提高了大型癌症研究的复制能力。”基因组生物学15(11): 503。https://doi.org/10.1186/s13059-014-0503-2．

福廷，让-菲利普，小蒂莫西·特里切，卡斯珀·D·汉森，2017。Illumina人类甲基化epic阵列与Minfi的预处理、标准化和集成生物信息学33(4): 558 - 60。https://doi.org/10.1093/bioinformatics/btw691．

贾菲，安德鲁·E，彼得·村上，李华金，杰弗里·T·莱克，M·丹尼尔·法林，安德鲁·P·范伯格，拉斐尔·A·伊瑞扎里，2012。“在表观遗传流行病学研究中寻找差异甲基化区域。”国际流行病学杂志41(1): 200 - 209。https://doi.org/10.1093/ije/dyr238．

Maksimovic, Jovana, Lavinia Gordon和Alicia Oshlack, 2012。“SWAN: Illumina Infinium HumanMethylation450串珠芯片阵列归一化子集分位数。”基因组生物学13 (6): R44机身内部。https://doi.org/10.1186/gb-2012-13-6-r44．

Touleimat, Nizar, Jörg Tost, 2012。“使用子集分位数归一化进行精确DNA甲基化估计的Infinium()人类甲基化450K珠子芯片数据处理的完整管道。”表观基因组学4(3): 325 - 41。https://doi.org/10.2217/epi.12.21．

特里谢，蒂莫西J，丹尼尔J魏森伯格，大卫范登伯格，彼得W莱尔德，金伯利D齐格蒙德。2013。Illumina Infinium DNA甲基化串珠阵列的低水平加工核酸研究41 (7): e90。https://doi.org/10.1093/nar/gkt090．