Analyze counts for sequential windows in a single region

Given a SummarizedExperiment with modified and total base counts in assays "Nmod" and "Nvalid", perform a pairwise statistical test for differential modification.

Usage

getDifferentiallyModifiedWindows(
  se,
  assayNameMod = "Nmod",
  assayNameValid = "Nvalid",
  groupCol = "group",
  verbose = FALSE
)

Arguments

se

SummarizedExperiment, for example returned by quantifyWindowsInRegion. It is expected to at least contain assays for modified and total counts (given by assayNameMod and assayNameValid) and a colData column that defines the groups (given by groupCol).

assayNameMod, assayNameValid

Character scalars that give the assay names in se containing the modified and total counts, respectively.

groupCol

Character scalar giving the column in colData(se) that defines the groups of samples to be compared.

verbose

Logical scalar. If TRUE, report on progress.

Value

The GRanges object constructed from the topTags output obtained for the statistical analysis, with an additional column named "dirNegLog10PValue", calculated as the sign of the logFC multiplied with the -log10(PValue).

Author

Panagiotis Papasaikas, Sebastien Smallwood, Charlotte Soneson, Michael Stadler

Examples

modbamfiles <- system.file("extdata",
                           c("6mA_1_10reads.bam", "6mA_1_10reads.bam",
                             "6mA_2_10reads.bam", "6mA_2_10reads.bam"),
                           package = "footprintR")
se <- quantifyWindowsInRegion(bamfiles = modbamfiles,
                              region = "chr1:6940000-6955000", modbase = "a",
                              BPPARAM = BiocParallel::SerialParam())
se$group <- c("group1", "group1", "group2", "group2")
gr <- getDifferentiallyModifiedWindows(se, groupCol = "group")
class(gr)
#> [1] "GRanges"
#> attr(,"package")
#> [1] "GenomicRanges"
head(gr)
#> GRanges object with 6 ranges and 6 metadata columns:
#>       seqnames          ranges strand |     logFC    logCPM        LR
#>          <Rle>       <IRanges>  <Rle> | <numeric> <numeric> <numeric>
#>   [1]     chr1 6940001-6940024      * |  6.922263   13.0368 25.338333
#>   [2]     chr1 6940013-6940036      * |  7.025231   13.4291 29.825436
#>   [3]     chr1 6940025-6940048      * |  5.122833   13.7647  9.165793
#>   [4]     chr1 6940037-6940060      * |  0.586798   13.7916  0.166427
#>   [5]     chr1 6940049-6940072      * |  1.214905   13.6515  1.284355
#>   [6]     chr1 6940061-6940084      * | -1.164122   13.5286  1.256727
#>            PValue         FDR dirNegLog10PValue
#>         <numeric>   <numeric>         <numeric>
#>   [1] 4.81053e-07 3.22306e-05          6.317807
#>   [2] 4.72749e-08 6.33483e-06          7.325370
#>   [3] 2.46581e-03 6.60836e-02          2.608041
#>   [4] 6.83307e-01 1.00000e+00          0.165384
#>   [5] 2.57091e-01 1.00000e+00          0.589913
#>   [6] 2.62272e-01 1.00000e+00         -0.581249
#>   -------
#>   seqinfo: 1 sequence from an unspecified genome; no seqlengths