methylTFRAnnotationMm10

Overview

methylTFRAnnotationMm10 provides the genome annotations that methylTFR needs to compute GC-bias-corrected transcription factor (TF) deviation scores on mm10. It contains three kinds of resources, each with its own accessor:

Accessor Returns Help page
getTFbindsites() TF binding sites (a GRangesList) ?tf_bindsites
getGCfreq() GC frequency tables, one per motif ?motif_gcfreq
getGenomeGC() GC content in 30 nt windows ?genomewide_GC

The GC tables are what make the deviation score bias-corrected. A motif whose binding sites sit in GC-rich sequence overlaps CpG islands more often than one that does not, and CpG islands are hypomethylated. Without a correction, that difference in sequence composition would look like a difference in TF activity.

The resources are several gigabytes in total, so they are hosted on AnnotationHub instead of being shipped with the package. The first call to an accessor downloads the resource; later calls read it from the local AnnotationHub cache.

Installation

if (!requireNamespace("BiocManager", quietly = TRUE)) {
    install.packages("BiocManager")
}
BiocManager::install("methylTFRAnnotationMm10")

# methylTFR, which uses these annotations
BiocManager::install("methylTFR")

Example data used in this vignette

Downloading the real resources would take too long for a vignette, so this vignette works on small synthetic objects that have exactly the structure of the real ones. Every accessor can read .rds files from a local directory instead of AnnotationHub when the option methylTFRAnnotationMm10.datadir is set. We write the synthetic objects to a temporary directory and point the option at it. With the option unset (the default), the same calls below would download the real data.

suppressPackageStartupMessages(library(GenomicRanges))
library(methylTFRAnnotationMm10)

mock_dir <- file.path(tempdir(), "methylTFRAnnotationMm10_vignette")
dir.create(mock_dir, showWarnings = FALSE)

# Binding sites: a GRangesList with one GRanges per motif. Each range is
# a motif match extended by 200 bp on either side.
mock_tfbs <- GRangesList(
    MA0001.1 = GRanges(
        c("chr1:10001-10411", "chr1:52001-52411", "chr2:3001-3411"),
        strand = c("+", "-", "+")
    ),
    MA0002.1 = GRanges(c("chr1:20001-20415", "chr3:8001-8415"),
        strand = c("-", "+")
    )
)
saveRDS(mock_tfbs, file.path(mock_dir, "altius_tf_bindsites.rds"))

# GC frequency tables: one 5-row matrix per motif. Rows are GC
# quintiles, columns are positions of a 30 nt window along the binding
# site, and every column sums to 1. We use 6 positions to keep it short.
set.seed(1)
mock_freq <- function(n_pos) {
    m <- matrix(runif(5 * n_pos), nrow = 5)
    sweep(m, 2, colSums(m), "/")
}
mock_gcfreq <- list(MA0001.1 = mock_freq(6), MA0002.1 = mock_freq(6))
saveRDS(mock_gcfreq, file.path(mock_dir, "altius_motif_gcfreq.rds"))

# Genome-wide GC: non-overlapping 30 nt windows with their GC fraction
# and GC quintile (1 = lowest, 5 = highest).
mock_gc <- GRanges("chr1", IRanges(start = seq(1, by = 30, length.out = 10),
    width = 30))
mock_gc$GC_bias <- round(runif(10, 0.2, 0.7), 2)
breaks <- quantile(mock_gc$GC_bias, probs = seq(0, 1, 0.2))
mock_gc$GC_bin <- findInterval(mock_gc$GC_bias, breaks,
    rightmost.closed = TRUE)
metadata(mock_gc) <- list(gc_breaks = breaks, bin_scope = "genome",
    step = 30L)
saveRDS(mock_gc, file.path(mock_dir, "genomewide_GC_mm10.rds"))

old_opt <- options(methylTFRAnnotationMm10.datadir = mock_dir)
list.files(mock_dir)
## [1] "altius_motif_gcfreq.rds" "altius_tf_bindsites.rds"
## [3] "genomewide_GC_mm10.rds"

Available motif sets

availableMotifSets() lists the motif sets that can be passed to getTFbindsites() and getGCfreq():

availableMotifSets()
## [1] "altius"     "cisbpv2"    "jaspar2020"

Binding sites: getTFbindsites()

getTFbindsites() returns a GRangesList named by motif ID, with one GRanges of binding sites per motif. Each range covers the motif match plus 200 bp on either side: the footprint window over which methylTFR reads methylation. All ranges of a motif therefore have the same width.

tf_bindsites <- getTFbindsites("altius")
tf_bindsites
## GRangesList object of length 2:
## $MA0001.1
## GRanges object with 3 ranges and 0 metadata columns:
##       seqnames      ranges strand
##          <Rle>   <IRanges>  <Rle>
##   [1]     chr1 10001-10411      +
##   [2]     chr1 52001-52411      -
##   [3]     chr2   3001-3411      +
##   -------
##   seqinfo: 3 sequences from an unspecified genome; no seqlengths
## 
## $MA0002.1
## GRanges object with 2 ranges and 0 metadata columns:
##       seqnames      ranges strand
##          <Rle>   <IRanges>  <Rle>
##   [1]     chr1 20001-20415      -
##   [2]     chr3   8001-8415      +
##   -------
##   seqinfo: 3 sequences from an unspecified genome; no seqlengths

The number of binding sites per motif, and the width of each motif’s window:

lengths(tf_bindsites)
## MA0001.1 MA0002.1 
##        3        2
vapply(width(tf_bindsites), unique, integer(1))
## MA0001.1 MA0002.1 
##      411      415

GC frequency tables: getGCfreq()

getGCfreq() returns a named list with one matrix per motif, using the same motif IDs as getTFbindsites(). Each matrix has five rows, one per genome-wide GC quintile (row 1 lowest GC, row 5 highest). Each column is one position of a 30 nt window sliding along the binding site. An entry is the fraction of the motif’s binding sites whose window at that position falls into that GC quintile, so every column sums to one. In the real data a matrix has several hundred columns.

gcfreqs <- getGCfreq("altius")
names(gcfreqs)
## [1] "MA0001.1" "MA0002.1"
dim(gcfreqs$MA0001.1)
## [1] 5 6
round(gcfreqs$MA0001.1, 2)
##      [,1] [,2] [,3] [,4] [,5] [,6]
## [1,] 0.11 0.28 0.09 0.15 0.43 0.19
## [2,] 0.16 0.30 0.08 0.21 0.10 0.01
## [3,] 0.25 0.21 0.31 0.29 0.30 0.19
## [4,] 0.39 0.20 0.17 0.11 0.06 0.44
## [5,] 0.09 0.02 0.35 0.23 0.12 0.17
colSums(gcfreqs$MA0001.1)
## [1] 1 1 1 1 1 1

methylTFR combines these fractions with the average methylation it observes in each GC quintile to get the methylation a motif would show from GC content alone, and reports the deviation from that expectation.

Genome-wide GC distribution: getGenomeGC()

getGenomeGC() returns a GRanges of non-overlapping 30 nt windows. GC_bias is the GC fraction of each window and GC_bin its genome-wide GC quintile. methylTFR overlaps each methylation call with these windows to assign it to a GC bin. The package serves a single assembly, so no argument is needed; getGenomeGC("mm10") is accepted too.

gc_dist <- getGenomeGC()
gc_dist
## GRanges object with 10 ranges and 2 metadata columns:
##        seqnames    ranges strand |   GC_bias    GC_bin
##           <Rle> <IRanges>  <Rle> | <numeric> <integer>
##    [1]     chr1      1-30      * |      0.66         5
##    [2]     chr1     31-60      * |      0.35         2
##    [3]     chr1     61-90      * |      0.43         3
##    [4]     chr1    91-120      * |      0.37         2
##    [5]     chr1   121-150      * |      0.53         4
##    [6]     chr1   151-180      * |      0.33         1
##    [7]     chr1   181-210      * |      0.44         3
##    [8]     chr1   211-240      * |      0.58         4
##    [9]     chr1   241-270      * |      0.24         1
##   [10]     chr1   271-300      * |      0.64         5
##   -------
##   seqinfo: 1 sequence from an unspecified genome; no seqlengths
table(gc_dist$GC_bin)
## 
## 1 2 3 4 5 
## 2 2 2 2 2

The quintile boundaries are stored with the object, so that the genome table and the motif tables always use the same scale:

metadata(gc_dist)$gc_breaks
##    0%   20%   40%   60%   80%  100% 
## 0.240 0.346 0.406 0.476 0.592 0.660

Using the annotations with methylTFR

The three objects are passed to methylTFR::run_methyltfr(). With real data this looks like:

library(methylTFR)
deviations <- run_methyltfr(
    "samples.tsv",               # sample annotation
    "samples_dir",               # directory with methylation calls
    tf_bindsites = getTFbindsites("jaspar2020"),
    gcfreqs      = getGCfreq("jaspar2020"),
    gc_dist      = getGenomeGC(),
    filetype     = "EPP"
)

See the methylTFR vignettes for a complete analysis.

Reading the real data from AnnotationHub

Once the local-directory option is removed, the accessors go to AnnotationHub. The records can also be browsed directly:

library(AnnotationHub)
ah <- AnnotationHub()
query(ah, "methylTFRAnnotationMm10")

If you have downloaded the .rds files by hand, set options(methylTFRAnnotationMm10.datadir = "/path/to/files") (or the environment variable METHYL_TFRANNOTATION_Mm10_DIR) to read them from disk.

Cleaning up

Restore the option and remove the synthetic files:

options(old_opt)
unlink(mock_dir, recursive = TRUE)
getOption("methylTFRAnnotationMm10.datadir")
## NULL

Session info

sessionInfo()
## R version 4.6.1 (2026-06-24)
## Platform: x86_64-pc-linux-gnu
## Running under: Ubuntu 26.04.1 LTS
## 
## Matrix products: default
## BLAS:   /usr/lib/x86_64-linux-gnu/openblas-pthread/libblas.so.3 
## LAPACK: /usr/lib/x86_64-linux-gnu/openblas-pthread/libopenblasp-r0.3.32.so;  LAPACK version 3.12.0
## 
## locale:
##  [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C              
##  [3] LC_TIME=en_US.UTF-8        LC_COLLATE=en_US.UTF-8    
##  [5] LC_MONETARY=en_US.UTF-8    LC_MESSAGES=en_US.UTF-8   
##  [7] LC_PAPER=en_US.UTF-8       LC_NAME=C                 
##  [9] LC_ADDRESS=C               LC_TELEPHONE=C            
## [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C       
## 
## time zone: Etc/UTC
## tzcode source: system (glibc)
## 
## attached base packages:
## [1] stats4    stats     graphics  grDevices utils     datasets  methods  
## [8] base     
## 
## other attached packages:
## [1] methylTFRAnnotationMm10_0.99.12 GenomicRanges_1.65.4           
## [3] Seqinfo_1.3.2                   IRanges_2.47.5                 
## [5] S4Vectors_0.51.10               BiocGenerics_0.59.12           
## [7] generics_0.1.4                  BiocStyle_2.41.0               
## 
## loaded via a namespace (and not attached):
##  [1] cli_3.6.6           knitr_1.52          rlang_1.3.0        
##  [4] xfun_0.61           otel_0.2.0          jsonlite_2.0.0     
##  [7] buildtools_1.0.0    htmltools_0.5.9     maketools_1.3.2    
## [10] sys_3.4.3           sass_0.4.10         rmarkdown_2.32     
## [13] evaluate_1.0.5      jquerylib_0.1.4     fastmap_1.2.0      
## [16] yaml_2.3.12         lifecycle_1.0.5     BiocManager_1.30.27
## [19] compiler_4.6.1      digest_0.6.39       R6_2.6.1           
## [22] bslib_0.12.0        tools_4.6.1         cachem_1.1.0