10  All in one

In DESeq2, la funzione DESeq() offre un modo conveniente per eseguire l’intera analisi dell’espressione differenziale in un unico passaggio.

Questa funzione integra i tre passaggi chiave appena visti:

1. estimateSizeFactors(): Normalizza i dati di conteggio per tenere conto delle differenze nella dimensione della libreria.
2. estimateDispersions(): Stima la dispersione dei geni, ovvero la variabilità dell’espressione tra i replicati.
3. nbinomWaldTest(): Esegue il test di Wald per identificare i geni differenzialmente espressi.

Utilizzo di dati non trasformati e non normalizzati

È fondamentale sottolineare che la funzione DESeq() deve essere applicata ai dati di conteggio grezzi, non trasformati e non normalizzati. Questo perché DESeq2 esegue internamente la normalizzazione e la trasformazione dei dati in modo appropriato per il modello statistico.

Vantaggi:

• Semplicità: Esegue l’intera analisi in un unico passaggio.
• Accuratezza: Applica i metodi di normalizzazione e trasformazione appropriati.
• Flessibilità: Offre diverse opzioni per personalizzare l’analisi.

ddsf <- DESeq(dds)
#> estimating size factors
#> estimating dispersions
#> gene-wise dispersion estimates
#> mean-dispersion relationship
#> final dispersion estimates
#> fitting model and testing
#> -- replacing outliers and refitting for 21 genes
#> -- DESeq argument 'minReplicatesForReplace' = 7 
#> -- original counts are preserved in counts(dds)
#> estimating dispersions
#> fitting model and testing
res <- results(ddsf, contrast = c(params$mycondition, params$mynum, params$mydenom), alpha = params$mypval, pAdjustMethod = params$mypadj)

La funzione results() in DESeq2 è fondamentale per estrarre ed esplorare i risultati dell’analisi dell’espressione differenziale.

results() estrae i risultati del test di Wald per l’espressione differenziale da un oggetto DESeqDataSet dopo aver eseguito DESeq(). Restituisce un oggetto DataFrame contenente informazioni sui geni, come log2 fold change, p-value e statistiche del test.

out_res <- tibble(
  gene = res@rownames,
  baseMean = res@listData$baseMean,
  log2FoldChange = res@listData$log2FoldChange,
  lfcSE = res@listData$lfcSE,
  # stat = res@listData$stat,
  pvalue = res@listData$pvalue,
  padj = res@listData$padj
) |>
  mutate(across(where(is.double), ~ round(., digits = 2)))

datatable(out_res,
  extensions = "Buttons",
  filter = "top",
  options = list(
    pageLength = 20,
    autoWidth = TRUE,
    dom = "Bfrtip",
    buttons = c("csv", "excel", "pdf")
  ), rownames = TRUE
)
#> Warning in instance$preRenderHook(instance): It seems your data is too big
#> for client-side DataTables. You may consider server-side processing:
#> https://rstudio.github.io/DT/server.html

summary(res)
#> 
#> out of 27430 with nonzero total read count
#> adjusted p-value < 0.01
#> LFC > 0 (up)       : 595, 2.2%
#> LFC < 0 (down)     : 592, 2.2%
#> outliers [1]       : 0, 0%
#> low counts [2]     : 10104, 37%
#> (mean count < 11)
#> [1] see 'cooksCutoff' argument of ?results
#> [2] see 'independentFiltering' argument of ?results