Installation

From Bioconductor

Development version can be installed from github

CRITICAL: After a new installation, one needs to cache the associated annotation data using the following command. This needs to be done only once per SeSAMe installation.

This function caches the needed SeSAMe annotation for all the supported platforms. One optionally caches only one specific platform (e.g., EPIC) by calling

SeSAMe annotation data is managed by the sesameData package which is based on the ExperimentHub package. You can find the path of the annotation data storage on your local computer using tools::R_user_dir("ExperimentHub", which="cache").

Beta Value Extraction

The raw Infinium BeadChip data are stored in IDAT files. Each sample has two IDAT files and they correspond to the red and green signal respectively. Green and red files for the same samples should always share the same sample name prefix. For example, 204529320035_R06C01_Red.idat and 204529320035_R06C01_Grn.idat correspond to the red and green signal of one sample. In SeSAMe, we will use the common prefix, i.e. 204529320035_R06C01, to refer to that sample. Sesame recognize both the raw idat as well as gzipped idats which are common in data stored in GEO. For example, in addition to the example above, SeSAMe also recognizes GSM2178224_Red.idat.gz and GSM2178224_Grn.idat.gz.

The function readIDATpair function reads in the signal intensity data from the IDAT pairs. The function takes the common prefix as input and outputs a SigSet. Using the two examples above, one would run the following commands.

Note that SeSAMe automatically detects and matches up the green and red signal files for the same sample. We will get back to the structure of SigSet below.

DNA methylation level (aka beta values) are defined as methylated signal/(methylated signal + unmethylated signal). It can be retrieved using the getBetas function on the SigSet. The output is a named vector with probe ID as names. For example, the following commands reads in one sample and convert it to beta values.

Please note that getBetas takes a single SigSet object instead of a list of them. To have it process many SigSets, you should combine that with looping functions lapply or mclapplys, or using the openSesame pipeline (see below).

Search IDAT Prefixes

Most often we will be working with a folder that contains many IDATs. Here is where the searchIDATprefixes function comes in handy. It lets us search all the IDATs in a folder and its subfolders recursively. Combine this with the R looping functions lets you process many IDATs without having to specify all IDAT names. searchIDATprefixes returns a named vector of prefixes with associated Red and Grn files, for readIDATpair:

A simple list of “SigSet”s are returned.

Custom-made Array

If you need to deal with a custom-made array instead of the standard array (MM285, EPIC, HM450 etc) supported natively by SeSAMe, you would need to provide a manifest that describe the probe information. You should be able to obtain that from the Illumina support website. The manifest should be made into a simple data frame with a minimum of four columns: Probe_ID, M, U and col. The easiest way to make the the manifest compatible with SeSAMe is by following internal manifests for a SeSAME-supported platform. They can be retrieved with the sesameDataGet function:

The col is either G (stand for Green) or R (stand for Red) or NA ( stand for both in the case of Infinium II design). For Infinium-II probes, the M column and col column is left to be NA. For example, one can check that both M and col columns are filled with the Infinium-I probes:

The last column mask is a logical vector that defines the default masking behavior of SeSAMe for the platform (see below for discussion of NA-masking).

With the manifest, your data can be processed using

Sum Two-channels

Beta values for Infinium-I probes can also be obtained by summing up the two in-band channel and out-of-band channel. This rescues probes with SNP hitting the extension base and hence switching color channel. More details can be found in Zhou et al 2017.

For color-channel-switching probes, extra SNP allele frequencies can be obtained by summing up methylated and umethylated alleles:

NA-Masking

If you call getBetas as is, you may notice that some of the beta values show up having NA values. This NA-masking is controlled internally using a mask attribute in SigSet. E.g.,

## cg00000029 cg00000103 cg00000109 cg00000155 cg00000158 cg00000165 
##  0.8237945         NA  0.8125637  0.9152265  0.9105163  0.8196466

To check probes to be NA-masked in a SigSet, one can use the mask function

## cg07881041 cg23229610 cg03513874 cg05451842 cg14797042 cg09838562 
##      FALSE       TRUE      FALSE      FALSE      FALSE      FALSE
## [1] 105454
## [1] 105454

Please be noted that mask in SigSet does not actually remove the probe reading but only specify how SeSAMe currently views the measurement. One can add more probes to the mask with the addMask function. Other functions such as the detection p-value calculation, also modifies mask. NA-masking influences other normalization and preprocessing functions. Therefore one should set it carefully. If one does not do any explicit masking, one simply gets the masking specified in the mask column of the manifest. This defines the default masking behavior of SeSAMe. For more details of masking, one can refer to Zhou et al 2017. One can override this masking by the resetMask function

## [1] 105454
## [1] 0

The getBetas function can also ignore NA-masking while extracting beta values:

## [1] 0

We recommend two types of probe masking:

  1. Experiment-dependent masking based on signal detection p-value (Zhou et al. 2018). Probes with p-value higher than a threshold (default: 0.05) are masked (see following for detection p-value calculation using the pOOBAH method).

  2. Experiment-independent probe masking due to design issues. This is typically designated in the mask column of the manifest (see Zhou et al. 2017): This masking supports EPIC, MM285, HM450 and HM27 and is turned on by default and can also be explicitly added by the function qualityMask.

Preprocessing

Detection p-value

As mentioned above, experiment-dependent masking based on signal detection p-values is effective in excluding artifactual methylation level reading and probes with too much influence from signal background. We recommend the pOOBAH algorithm that was based on Infinium-I probe out-of-band signal for calibrating the distribution of the signal background:

## [1] 105454
## [1] 138116
## [1] 38136

Sometimes one would want to calculation detection p-value without creating masking like in the case of having to upload the p-value to GEO. In those cases one can use the return.pval option and add pvalue-based mask later.

Background Subtraction

SeSAMe implements the background subtraction based on normal-exponential deconvolution using out-of-band probes noob (Triche et al. 2013) and optionally with extra bleed-through subtraction. Signal bleed-through happens when measurement from one channel affects the measurement in the other channel. SeSAMe’s noob further removes residual background by regressing out the green-to-red and red-to-green relationship using Type-I probes.

One can use following beta to total signal intensity to check the effect of background subtraction.

Dye Bias Correction

Dye bias refers to the difference in signal intensity between the two color channel. SeSAMe offers two flavors of dye bias correction: linear scaling (dyeBiasCorr) and nonlinear scaling (dyeBiasCorrTypeINorm). Linear scaling equalize the mean of all probes from the two color channel.

Residual dye bias can be corrected using nonlinear quantile interpolation with Type-I probes.

Under this correction, Infinium-I Red probes and Infinium-I Grn probes have the same distribution of signal.

Note that linear scaling does not shift beta values of Type-I probes while nonlinear scaling does shift beta values of Type-I probes.

Channel Inference

Sometimes Infinium-I channel spec is inaccurate in the manifest. We can infer the channel using data.

## Type-I color channel reset:
## R>R: 92070
## G>G: 48959
## R>G: 69
## G>R: 933
## Red Failed: 155
## Grn Failed: 97

The openSesame Pipeline

We have discussed noob, nonlinear dye bias correction and pOOBAH. Put together, this can be implemented as follows

Here we use two cores with mc.cores=2.

Equivalently, sesame provides the openSesame pipeline

as a quickstart default. Here idat_dir is the directory containing all the IDAT files. Like readIDATpair, openSesame also works with custom-made array with a manifest file (see above):

The SigSet Class

SeSAMe design includes alight-weight full exposure of internal signal intensities (essential information for users of Illumina methylation array data, as demonstrated in Zhou et al 2018), which permits sensitive and specific joint inference on copy number and DNA methylation.

Central to the SeSAMe platform is the SigSet data structure, an S4 class with slots containing signals for six different classes of probes:

  1. II - Type-II probes;
  2. IR - Type-I Red channel probes;
  3. IG - Type-I Grn channel probes;
  4. oobG - Out-of-band Grn channel probes (matching Type-I Red channel probes in number);
  5. oobR - Out-of-band Red channel probes (matching Type-I Grn channel probes in number);
  6. ctl - control probes.

For all save control probes, signal intensities are stored as an Nx2 numeric matrix, with N representing the number of probes in the class. The two columns of the matrix represent the methylated probe intensity and the unmethylated probe intensity. (Previously, this was implemented in an R6 Reference class, SignalSet. The current S4 implementation in SigSet complies with Bioconductor guidelines, and for backwards compatibility, the signalR6toS4 function transforms a SignalSet to a SigSet.

For example, printing the SigSet directly shows its content

## SigSet EPIC 
##  - @IG Probes: 49989 - 332 4145 70 7094 599 2958 ...
##  - @IR Probes: 92294 - 183 8040 1949 6152 833 89 ...
##  - @II Probes: 724612 - 6543 1596 3133 1011 3035 2837 ...
##  - @oobG Probes: 92294 - 138 277 107 218 232 80 ...
##  - @oobR Probes: 49989 - 1013 150 81 910 448 183 ...
##  - Control Probes: 635 ...
##  - Number of Masked Probes: 105454
##               M    U
## cg07881041 6543 1011
## cg23229610 1596 3035
## cg03513874 3133 2837
## cg05451842  376 5673
## cg14797042 3571  368
## cg09838562  131 2523

Similarly, signals for Type-I probes can be accessed from sset@IR and sset@IG and out-of-band signals from sset@oobG and sset@oobR. As one can see the probe names (row names) of IR always coincide with the probe names (row names) of oobG (and vice versa). This is because the out-of-band probe signal for red channel probes is in green channel (and vice versa).

Interaction with Minfi

SigSet can be converted back and forth from Minfi RGChannelSet in multiple ways. One can sesamize a minfi RGChannelSet which returns a GenomicRatioSet. Here we are illustrating using the FlowSorted.Blood.450k object, which is distributed in the minfi::RGChannelSet.

## Sesamizing WB_105...
## Sesamizing WB_218...
## Sesamizing WB_261...
## Sesamizing PBMC_105...
## Loading required package: IlluminaHumanMethylation450kanno.ilmn12.hg19
## Loading required package: IlluminaHumanMethylation450kmanifest
## class: GenomicRatioSet 
## dim: 485512 4 
## metadata(1): SNPs
## assays(2): Beta CN
## rownames(485512): cg13869341 cg14008030 ... cg08265308 cg14273923
## rowData names(0):
## colnames(4): WB_105 WB_218 WB_261 PBMC_105
## colData names(8): Sample_Name Slide ... CellType Sex
## Annotation
##   array: IlluminaHumanMethylation450k
##   annotation: ilmn12.hg19
## Preprocessing
##   Method: SeSAMe (type I)
##   minfi version: 1.38.0
##   Manifest version: 0.6.0