In omorante/shinyepico: ShinyÉPICo

knitr::include_graphics("images/logo.png")

This report was created with shinyÉPICo.
© 2020 Octavio Morante-Palacios

knitr::opts_chunk$set(echo = TRUE)

knitr::knit_engines$set(asis = function(options) {
  if (options$echo && options$eval) knitr::knit_child(text = options$code)
})

```{css, echo=FALSE} div, code {overflow-x:auto !important}

## Variable and samples selection

The group variable selected was __`r params$grouping_var`__, and the donor variable was __`r params$donor_var`__ .

Finally, the selected samples added to the RGSet file were:

```r
DT::datatable(params$rval_sheet_target, extensions = 'Buttons',
    style = "bootstrap",
    rownames = FALSE,
    selection = "single",
    options = list(dom = 'Blfrtip',
    columnDefs = list(list(
      targets = match("Basename",colnames(params$rval_sheet_target)) - 1, visible = FALSE)),
     #scrollX = TRUE,
    buttons = c('copy', 'csv', 'excel', 'pdf', 'print'),
    lengthMenu = list(c(10,25,50,-1), c(10,25,50,"All"))))

The following samples were excluded for the analysis:

cat(params$rval_sheet[[params$name_var]][!(params$rval_sheet[[params$name_var]] %in% params$rval_sheet_target[[params$name_var]])], 
    sep = ", ")

Bad probes filtering

After Minfi Normalization, 'bad' probes in RGChannelSet are calculated using a standard p.value threshold (0.01). They will be filtered in the normalized data:

  detP <- minfi::detectionP(rgset)
  bad_pos <- row.names(as.data.frame(detP))[rowMeans(detP) > detectP]

QC and Array Normalization

The selected normalization method was r params$normalization_mode. In the next plots, we can see the comparison of the processed data with the raw data.

QC Signal plot

In this plot, the average methylated signal vs the average unmethylated signal of each sample is represented. Both average signals should be higher than the signal threshold (10) to be considered OK. Lower or very variable signals can reflect problems in the sample material or in the scanning, but the user is responsible for deciding whether a sample should be excluded.

params$plot_qcraw

QC Bisulfite Conversion plot

Illumina 450k and EPIC arrays have probes to determine if bisulfite conversion has been successful or not. In this graph, the minimum ratio of converted/non-converted signal of bisulfite conversion II probes for each sample is represented.

params$plot_bisulfiterawII

Density plot

Raw

params$plot_densityplotraw

Processed (r params$normalization_mode)

params$plot_densityplot

Boxplot

Raw

   params$plot_boxplotraw

Processed (r params$normalization_mode)

   params$plot_boxplot

SNPs Heatmap

The methylation array has 65 specific SNP probes. These SNP probes are intended to be used for sample tracking and sample mixups. Each SNP probe ought to have values clustered around 3 distinct values corresponding to homo-, and hetero-zygotes. Therefore, different samples of the same donor should cluster together.

Processed (r params$normalization_mode)

params$plot_snpheatmap

Sex Prediction

Depending on the average chromosome Y signal and the average chromosome X signal is possible to predict the sex of the sample donors.

params$plot_sexprediction

knitr::kable(data.frame(name=params$rval_sheet_target[[params$name_var]],sex=params$data_sexprediction))

Principal Component Analysis

Plot of PCA is shown, with the selected principal components.

params$plot_pcaplot

DT::datatable(
  format(params$table_pcaplot, digits=2),
  selection = "single",
  rownames = TRUE,
  style = "bootstrap",
  class = "table-condensed",
  options(list(
    #scrollX = TRUE,
    lengthChange = FALSE,
    paging = FALSE,
    searching = FALSE,
    info = FALSE,
    ordering = FALSE
  ))
)

Correlations

Correlating principal components with variables we can determine if Beta values are related to our variable of interest or other variables. This can also be useful to determine possible errors in the sample hybridization randomization and to select covariables to add to the linear model.

params$plot_corrplot

Not useful variables are discarded and the variable type is autodetected. The autodetected variable types were:

DT::datatable(params$table_corrplot, extensions = 'Buttons',
    style = 'bootstrap',
    rownames = FALSE,
    selection = "single",
    options = list(dom = 'Blfrtip',
    buttons = c('copy', 'csv', 'excel', 'pdf', 'print'),
    lengthMenu = list(c(10,25,50,-1),
    c(10,25,50,"All"))))

Predict sex and drop SNPs

After Minfi normalization, we obtained a GenomicRatioSet with some transformations depending on the type of normalization chosen. To understand better the different classes of Minfi packages, and their relations depending of the normalization method, please read this vignette of Minfi creators (it is also valid for Illumina EPIC arrays).

Options selected:

• Drop SNPs: r params$dropsnps r if(params$dropsnps) paste('minimum MAF filtered:', params$maf)
• Drop CpHs: r params$dropcphs
• Drop X/Y Chromosomes: r params$dropsex

Moreover, we also removed SNP loci, CpHs and X/Y chromosome positions if these options are selected:

minfi::dropLociWithSnps(gset, maf = 0))

minfi::dropMethylationLoci(gset, dropRS = TRUE, dropCH = TRUE)

gset[rownames(minfi::getAnnotation(gset))[!(minfi::getAnnotation(gset)$chr %in% c("chrX","chrY"))],]

Finally, we added sex information:

minfi::addSex(gset)

After normalization, r params$probes DNA methylation positions remained.

Model generation

Then, the linear model was generated using r params$limma_voi as variable of interest, r if(length(params$limma_covar) == 0) "None" else paste(params$limma_covar, collapse=", ") as covariable(s), and r if(length(params$limma_inter) == 0) "None" else paste(params$limma_inter, collapse=", ") as interaction(s).

The resulting design matrix was:

DT::datatable(
  params$rval_design,
  selection = "single",
  rownames = TRUE,
  style = "bootstrap",
  options(list(
    #scrollX = TRUE,
    lengthChange = FALSE,
    searching = FALSE
  ))
)

Optionally, ArrayWeights can be activated to weigh the value of each sample depending on their adjust to the linear model:

      if (weights_option){
        try({weights = limma::arrayWeights(MValues, design = design})
      }
        else { weights = NULL}

      fit = limma::lmFit(MValues, design, weights = weights)

Where weights is the option specified in the form. In this case, weights is r params$limma_arrayweights .

The relation variance/average of the fit model was:

  params$plot_plotSA

With this sample groups, the possible contrasts were:

cat(params$rval_contrasts, sep = "\n")

Using the calculated model, the t-statistics and p.values are obtained automatically for each contrast:

  mcontrast = makeContrasts(contrasts=contrast, levels = design)
  fitting = contrasts.fit(fit, mcontrast)
  fitting = eBayes(fitting, trend=trend, robust=robust)
  toptable = topTable(fitting, coef= 1,  adjust.method = "fdr", number =  Inf)

Where trend and robust are the options specified in the form. In this case, trend is r params$limma_ebayes_trend and robust is r params$limma_ebayes_robust .

Filtered results (DMPs)

For each contrast ( r paste(params$dmrs_contrasts, collapse=", ") ) , the statistics for each CpG were generated and, after, filtered with these criteria:

• Max. FDR: r params$max_fdr
• Min. Δβ: r params$min_deltabeta
• Max. pvalue: r params$max_pvalue

Table with all contrasts (DMPs)

knitr::kable(params$table_dmps)

Custom heatmap (DMPs)

• Column dendogram: r params$Colv
• Scale: r params$scale
• Cluster Algorithm: r params$clusteralg
• Distance: r params$distance

The data used for the representation are r if(!params$removebatch) "the normalized beta values" else "the normalized beta values after removing covariable effects on the data (with the limma removeBatchEffect function), preserving the variable of interest".

In the custom heatmap the selected contrasts to plot were: r paste(params$contrasts2plot, collapse=", ") and the groups were: r paste(params$groups2plot, collapse=", ") .

    if(!is.null(params$filteredlist2heatmap)) create_heatmap(params$filteredlist2heatmap, factorgroups =  factor(params$rval_voi[params$rval_voi %in% params$groups2plot], levels = params$groups2plot), groups2plot = params$rval_voi %in% params$groups2plot, Colv = as.logical(params$Colv), clusteralg = params$clusteralg, distance = params$distance, scale = params$scale, static=TRUE, ColSideColors = TRUE, RowSideColors = params$rval_dendrogram) else cat("Differences are not in the plotting range (<12000, >1)")

DMPs in heatmap: r if(is.null(nrow(params$filteredlist2heatmap))) "None" else nrow(params$filteredlist2heatmap)

```{asis filt_res, echo=TRUE, eval=!is.null(params$dmrs_contrasts)}

Filtered results (DMRs)

For each contrast ( r paste(params$rval_contrasts, collapse=", ") ) , the statistics for each CpG were generated and, after, filtered with these criteria:

• Max. FDR: r params$dmrs_max_fdr
• Min. Δβ: r params$dmrs_min_deltabeta
• Max. pvalue: r params$dmrs_max_pvalue

### Table with all contrasts (DMRs)

```r
if(!is.null(params$dmrs_contrasts)) knitr::kable(params$table_dmrs) else cat("DMRs were not calculated in this analysis.")

```{asis custom_heatmap, echo=TRUE, eval=!is.null(params$dmrs_contrasts)}

Custom heatmap (DMRs)

• Column dendogram: r params$dmrs_Colv
• Scale: r params$dmrs_scale
• Cluster Algorithm: r params$dmrs_clusteralg
• Distance: r params$dmrs_distance

The data used for the representation are r if(!params$dmrs_removebatch) "the normalized beta values" else "the normalized beta values after removing covariable effects on the data (with the limma removeBatchEffect function), preserving the variable of interest".

In the custom heatmap the selected contrasts to plot were: r paste(params$dmrs_contrasts2plot, collapse=", ") and the groups were: r paste(params$dmrs_groups2plot, collapse=", ") .

```r

if(!is.null(params$filteredmcsea2heatmap) & !is.null(params$dmrs_contrasts)) create_heatmap(params$filteredmcsea2heatmap, factorgroups =  factor(params$rval_voi[params$rval_voi %in% params$dmrs_groups2plot], levels = params$dmrs_groups2plot), groups2plot = params$rval_voi %in% params$dmrs_groups2plot, Colv = as.logical(params$dmrs_Colv), clusteralg = params$dmrs_clusteralg, distance = params$dmrs_distance, scale = params$dmrs_scale, static=TRUE, ColSideColors = TRUE, RowSideColors = params$dmrs_rval_dendrogram) else cat("Differences are not in the plotting range (<12000, >1)")

DMRs in heatmap: r if(is.null(nrow(params$filteredmcsea2heatmap))) "None" else nrow(params$filteredmcsea2heatmap)

System and analysis information

This report was generated at r Sys.time().

The session information was the following:

print(sessionInfo())

omorante/shinyepico documentation built on May 11, 2023, 7:21 p.m.