Installation

PscanR is being prepared for Bioconductor submission. After acceptance, install it with the following command (not executed while building this vignette):

BiocManager::install("PscanR")

For pre-submission installation from source, see the package README.

Introduction

PscanR identifies transcription factor binding motifs enriched in regulatory sequences. Promoters of co-expressed or co-regulated genes are scanned with a collection of motif matrices, and their mean best-match scores are compared with those of a reference promoter population. Enriched motifs suggest candidate regulators.

PscanR provides access to precomputed backgrounds for five model organisms: human, mouse, and fruit fly (RefSeq annotation), Arabidopsis (TAIR annotation), and yeast (SGD annotation). Custom backgrounds can also be constructed for these or other organisms, as described in Building a background for PscanR.

This vignette demonstrates motif enrichment analysis of A. thaliana promoters and the functions available for examining the results.

The data

The example uses gene expression data from a study of BASIC PENTACYSTEINE (BPC) transcription factors, which control differentiation of the medial domain of the Arabidopsis gynoecium.

Caselli F, Palermiti M, Petrella R, Morlacchi VA, Dünser K, Kleine-Vehn J, Chiara M, Gregis V. A novel role of BPCs in the control of medial domain differentiation during gynoecium development in Arabidopsis thaliana. Plant Mol Biol 116(1):4 (2025). doi:10.1007/s11103-025-01662-x

The input is a two-column table containing gene identifiers and directions of expression change.

The examples below use bundled results to avoid network access and repeated genome-wide sequence extraction.

gene_sets <- read.delim(
    fixture("arabidopsis_bpc_gene_sets.tsv"),
    stringsAsFactors = FALSE
)
str(gene_sets)
#> 'data.frame':    2253 obs. of  2 variables:
#>  $ gene_id  : chr  "AT1G25145" "AT2G43020" "AT1G29951" "AT2G01422" ...
#>  $ direction: chr  "down" "down" "down" "down" ...
table(gene_sets$direction)
#> 
#> down   up 
#> 1247 1006

The rest of this vignette follows the 1247 genes with reduced expression when BPC function is impaired.

The background

Precomputed backgrounds are indexed by organism, genome assembly, promoter window, and motif collection. ps_available_bg() lists the available resources.

ps_available_bg(keyword = "TAIR9", details = TRUE)
#>     schema_version    status latest organism assembly upstream downstream
#> 31               1 validated   TRUE       at    TAIR9      200         50
#> 32               1 validated   TRUE       at    TAIR9      450         50
#> 33               1 validated   TRUE       at    TAIR9      500          0
#> 34               1 validated   TRUE       at    TAIR9      950         50
#> 35               1 validated   TRUE       at    TAIR9     1000          0
#> 66               1 validated   TRUE       at    TAIR9      200         50
#> 67               1 validated   TRUE       at    TAIR9      450         50
#> 68               1 validated   TRUE       at    TAIR9      500          0
#> 69               1 validated   TRUE       at    TAIR9      950         50
#> 70               1 validated   TRUE       at    TAIR9     1000          0
#> 101              1 validated   TRUE       at    TAIR9      200         50
#> 102              1 validated   TRUE       at    TAIR9      450         50
#> 103              1 validated   TRUE       at    TAIR9      500          0
#> 104              1 validated   TRUE       at    TAIR9      950         50
#> 105              1 validated   TRUE       at    TAIR9     1000          0
#>     jaspar_release tax_group background_version
#> 31            2020    plants                  2
#> 32            2020    plants                  2
#> 33            2020    plants                  2
#> 34            2020    plants                  2
#> 35            2020    plants                  2
#> 66            2022    plants                  2
#> 67            2022    plants                  2
#> 68            2022    plants                  2
#> 69            2022    plants                  2
#> 70            2022    plants                  2
#> 101           2024    plants                  2
#> 102           2024    plants                  2
#> 103           2024    plants                  2
#> 104           2024    plants                  2
#> 105           2024    plants                  2
#>                                            artifact
#> 31  backgrounds/J2020_TAIR9_200u_50d_TAIR.psbg2.txt
#> 32  backgrounds/J2020_TAIR9_450u_50d_TAIR.psbg2.txt
#> 33   backgrounds/J2020_TAIR9_500u_0d_TAIR.psbg2.txt
#> 34  backgrounds/J2020_TAIR9_950u_50d_TAIR.psbg2.txt
#> 35  backgrounds/J2020_TAIR9_1000u_0d_TAIR.psbg2.txt
#> 66  backgrounds/J2022_TAIR9_200u_50d_TAIR.psbg2.txt
#> 67  backgrounds/J2022_TAIR9_450u_50d_TAIR.psbg2.txt
#> 68   backgrounds/J2022_TAIR9_500u_0d_TAIR.psbg2.txt
#> 69  backgrounds/J2022_TAIR9_950u_50d_TAIR.psbg2.txt
#> 70  backgrounds/J2022_TAIR9_1000u_0d_TAIR.psbg2.txt
#> 101 backgrounds/J2024_TAIR9_200u_50d_TAIR.psbg2.txt
#> 102 backgrounds/J2024_TAIR9_450u_50d_TAIR.psbg2.txt
#> 103  backgrounds/J2024_TAIR9_500u_0d_TAIR.psbg2.txt
#> 104 backgrounds/J2024_TAIR9_950u_50d_TAIR.psbg2.txt
#> 105 backgrounds/J2024_TAIR9_1000u_0d_TAIR.psbg2.txt
#>                                                      artifact_sha256
#> 31  4b72e41054e30fa15ca4a047a571691835872c46f98f3244fb621c2329aaaaf9
#> 32  063bac7565c2a3183ca0f25d5e56667cec29643a66af47e6bcf6d63b95015de3
#> 33  ec5d910fbe72b6cbc16d610ce3995330fa4746069cc1fab84de63f6431c5ceab
#> 34  bab5e89332bd11d47fefc36a2a817a0df989bb8fde9bae60aa216eb1fa271669
#> 35  c1e5ce982c9747856160ad03529ec1b8a1230a12496083e29dfb497a87470f6f
#> 66  2d985c18a41f1fdc0af8e4d6effbb81c39ab3abd349f1b1d8f8bc5415af64447
#> 67  50d39d9263c7c80bf0ec5503a9c4e36a34084657ef55198f5cc9a5dd11b5cbea
#> 68  dfb462f45e234ba55e92d284fa8693c06e10116c0b73eb20fcc4a2628abcfa47
#> 69  c19383cf66d8544a6c55a691b3edb00990de67ee2d5bca66b03285a8a118ee85
#> 70  e9a9ca242f894f902018315ca8ca74e6d0f2633d410989218977654ea8103319
#> 101 1459bf5f91f07b0ce76f569df5082e1c9e88600c71a0bf4d95fe4d0ff766f1f9
#> 102 9897d45487fd166de17119735bf8cb33c82a627e1ce16d97306bc7519048ada2
#> 103 c9c74733abc74cc86859977c28bcfc6ddb86481a9be408ba5dffe85d5a1c7430
#> 104 cdb2f37941f1a2afe53a5372d6d010a89fab96115c0b3b7af7afa3e02c547808
#> 105 e9ff5282cd0ccffd939f8b91e45a3b626705243f65ff894e43851511b250475f
#>     motif_count promoter_count
#> 31          530          36331
#> 32          530          36368
#> 33          530          36371
#> 34          530          36406
#> 35          530          36405
#> 66          656          36331
#> 67          656          36368
#> 68          656          36371
#> 69          656          36406
#> 70          656          36405
#> 101         805          36331
#> 102         805          36368
#> 103         805          36371
#> 104         805          36406
#> 105         805          36405
#>                                                      annotation_hash
#> 31  4ea50fc5006b5a92fc39cdfcb2c72305ed0aebc2d39c9be1ca89efce88a2f9ac
#> 32  4ea50fc5006b5a92fc39cdfcb2c72305ed0aebc2d39c9be1ca89efce88a2f9ac
#> 33  4ea50fc5006b5a92fc39cdfcb2c72305ed0aebc2d39c9be1ca89efce88a2f9ac
#> 34  4ea50fc5006b5a92fc39cdfcb2c72305ed0aebc2d39c9be1ca89efce88a2f9ac
#> 35  4ea50fc5006b5a92fc39cdfcb2c72305ed0aebc2d39c9be1ca89efce88a2f9ac
#> 66  4ea50fc5006b5a92fc39cdfcb2c72305ed0aebc2d39c9be1ca89efce88a2f9ac
#> 67  4ea50fc5006b5a92fc39cdfcb2c72305ed0aebc2d39c9be1ca89efce88a2f9ac
#> 68  4ea50fc5006b5a92fc39cdfcb2c72305ed0aebc2d39c9be1ca89efce88a2f9ac
#> 69  4ea50fc5006b5a92fc39cdfcb2c72305ed0aebc2d39c9be1ca89efce88a2f9ac
#> 70  4ea50fc5006b5a92fc39cdfcb2c72305ed0aebc2d39c9be1ca89efce88a2f9ac
#> 101 4ea50fc5006b5a92fc39cdfcb2c72305ed0aebc2d39c9be1ca89efce88a2f9ac
#> 102 4ea50fc5006b5a92fc39cdfcb2c72305ed0aebc2d39c9be1ca89efce88a2f9ac
#> 103 4ea50fc5006b5a92fc39cdfcb2c72305ed0aebc2d39c9be1ca89efce88a2f9ac
#> 104 4ea50fc5006b5a92fc39cdfcb2c72305ed0aebc2d39c9be1ca89efce88a2f9ac
#> 105 4ea50fc5006b5a92fc39cdfcb2c72305ed0aebc2d39c9be1ca89efce88a2f9ac
#>                                                        promoter_hash
#> 31  dbd7d029ae8be7cd91214ce1972458c7175720e74f9866e72f6faa6a909a04af
#> 32  386c50763de2bc8dde69442d05e399a4322903d64bdb3b9f30c1731c4e98fb66
#> 33  394a22087789940987c669e35cfaf2677b3c99b5880a0c12e1954401653118b6
#> 34  7ab3e5cc19986b1c6cd9fd0b7a670a57259bdcaf2155cbbb4d5681478fbb2946
#> 35  cff3236aa5cf9d592417cf95f63898e129f6567975eda5e600ef3f2c074d202b
#> 66  dbd7d029ae8be7cd91214ce1972458c7175720e74f9866e72f6faa6a909a04af
#> 67  386c50763de2bc8dde69442d05e399a4322903d64bdb3b9f30c1731c4e98fb66
#> 68  394a22087789940987c669e35cfaf2677b3c99b5880a0c12e1954401653118b6
#> 69  7ab3e5cc19986b1c6cd9fd0b7a670a57259bdcaf2155cbbb4d5681478fbb2946
#> 70  cff3236aa5cf9d592417cf95f63898e129f6567975eda5e600ef3f2c074d202b
#> 101 dbd7d029ae8be7cd91214ce1972458c7175720e74f9866e72f6faa6a909a04af
#> 102 386c50763de2bc8dde69442d05e399a4322903d64bdb3b9f30c1731c4e98fb66
#> 103 394a22087789940987c669e35cfaf2677b3c99b5880a0c12e1954401653118b6
#> 104 7ab3e5cc19986b1c6cd9fd0b7a670a57259bdcaf2155cbbb4d5681478fbb2946
#> 105 cff3236aa5cf9d592417cf95f63898e129f6567975eda5e600ef3f2c074d202b
#>                                                           motif_hash
#> 31  d7efd20729fe014f4d95b4b8bf30c61e332b6d806f3faa06ff20e092e412bc09
#> 32  d7efd20729fe014f4d95b4b8bf30c61e332b6d806f3faa06ff20e092e412bc09
#> 33  d7efd20729fe014f4d95b4b8bf30c61e332b6d806f3faa06ff20e092e412bc09
#> 34  d7efd20729fe014f4d95b4b8bf30c61e332b6d806f3faa06ff20e092e412bc09
#> 35  d7efd20729fe014f4d95b4b8bf30c61e332b6d806f3faa06ff20e092e412bc09
#> 66  dc6d29707c6c8f92262d741e7da057d9b3c4bb59fd3e64df7c16d9acbbe7455f
#> 67  dc6d29707c6c8f92262d741e7da057d9b3c4bb59fd3e64df7c16d9acbbe7455f
#> 68  dc6d29707c6c8f92262d741e7da057d9b3c4bb59fd3e64df7c16d9acbbe7455f
#> 69  dc6d29707c6c8f92262d741e7da057d9b3c4bb59fd3e64df7c16d9acbbe7455f
#> 70  dc6d29707c6c8f92262d741e7da057d9b3c4bb59fd3e64df7c16d9acbbe7455f
#> 101 ce3adfb94e3a19e3db37bc1c88a99a09933845a7a1f7eafad9f7a5a00fae58d4
#> 102 ce3adfb94e3a19e3db37bc1c88a99a09933845a7a1f7eafad9f7a5a00fae58d4
#> 103 ce3adfb94e3a19e3db37bc1c88a99a09933845a7a1f7eafad9f7a5a00fae58d4
#> 104 ce3adfb94e3a19e3db37bc1c88a99a09933845a7a1f7eafad9f7a5a00fae58d4
#> 105 ce3adfb94e3a19e3db37bc1c88a99a09933845a7a1f7eafad9f7a5a00fae58d4
#>                                          annotation_snapshot
#> 31  annotation_snapshots/TAIR9_200u_50d_dbd7d029ae8be7cd.rds
#> 32  annotation_snapshots/TAIR9_450u_50d_386c50763de2bc8d.rds
#> 33   annotation_snapshots/TAIR9_500u_0d_394a220877899409.rds
#> 34  annotation_snapshots/TAIR9_950u_50d_7ab3e5cc19986b1c.rds
#> 35  annotation_snapshots/TAIR9_1000u_0d_cff3236aa5cf9d59.rds
#> 66  annotation_snapshots/TAIR9_200u_50d_dbd7d029ae8be7cd.rds
#> 67  annotation_snapshots/TAIR9_450u_50d_386c50763de2bc8d.rds
#> 68   annotation_snapshots/TAIR9_500u_0d_394a220877899409.rds
#> 69  annotation_snapshots/TAIR9_950u_50d_7ab3e5cc19986b1c.rds
#> 70  annotation_snapshots/TAIR9_1000u_0d_cff3236aa5cf9d59.rds
#> 101 annotation_snapshots/TAIR9_200u_50d_dbd7d029ae8be7cd.rds
#> 102 annotation_snapshots/TAIR9_450u_50d_386c50763de2bc8d.rds
#> 103  annotation_snapshots/TAIR9_500u_0d_394a220877899409.rds
#> 104 annotation_snapshots/TAIR9_950u_50d_7ab3e5cc19986b1c.rds
#> 105 annotation_snapshots/TAIR9_1000u_0d_cff3236aa5cf9d59.rds
#>                 annotation_source annotation_retrieved_at      scoring_spec
#> 31  raw_data/TAIR9_GFF3_genes.gff 2026-08-04 16:51:07 UTC pscan-best-hit-v1
#> 32  raw_data/TAIR9_GFF3_genes.gff 2026-08-04 16:57:34 UTC pscan-best-hit-v1
#> 33  raw_data/TAIR9_GFF3_genes.gff 2026-08-04 17:07:25 UTC pscan-best-hit-v1
#> 34  raw_data/TAIR9_GFF3_genes.gff 2026-08-04 17:19:56 UTC pscan-best-hit-v1
#> 35  raw_data/TAIR9_GFF3_genes.gff 2026-08-04 17:37:09 UTC pscan-best-hit-v1
#> 66  raw_data/TAIR9_GFF3_genes.gff 2026-08-04 16:52:57 UTC pscan-best-hit-v1
#> 67  raw_data/TAIR9_GFF3_genes.gff 2026-08-04 17:00:59 UTC pscan-best-hit-v1
#> 68  raw_data/TAIR9_GFF3_genes.gff 2026-08-04 17:11:01 UTC pscan-best-hit-v1
#> 69  raw_data/TAIR9_GFF3_genes.gff 2026-08-04 17:26:26 UTC pscan-best-hit-v1
#> 70  raw_data/TAIR9_GFF3_genes.gff 2026-08-04 17:43:23 UTC pscan-best-hit-v1
#> 101 raw_data/TAIR9_GFF3_genes.gff 2026-08-04 16:54:45 UTC pscan-best-hit-v1
#> 102 raw_data/TAIR9_GFF3_genes.gff 2026-08-04 17:04:26 UTC pscan-best-hit-v1
#> 103 raw_data/TAIR9_GFF3_genes.gff 2026-08-04 17:14:44 UTC pscan-best-hit-v1
#> 104 raw_data/TAIR9_GFF3_genes.gff 2026-08-04 17:32:18 UTC pscan-best-hit-v1
#> 105 raw_data/TAIR9_GFF3_genes.gff 2026-08-04 17:49:51 UTC pscan-best-hit-v1
#>     pscanr_version                            pscanr_commit bioc_version
#> 31          0.99.0 47b59b6a3135a032e02fbd4eb5608397c5cd9005         3.23
#> 32          0.99.0 47b59b6a3135a032e02fbd4eb5608397c5cd9005         3.23
#> 33          0.99.0 47b59b6a3135a032e02fbd4eb5608397c5cd9005         3.23
#> 34          0.99.0 47b59b6a3135a032e02fbd4eb5608397c5cd9005         3.23
#> 35          0.99.0 47b59b6a3135a032e02fbd4eb5608397c5cd9005         3.23
#> 66          0.99.0 47b59b6a3135a032e02fbd4eb5608397c5cd9005         3.23
#> 67          0.99.0 47b59b6a3135a032e02fbd4eb5608397c5cd9005         3.23
#> 68          0.99.0 47b59b6a3135a032e02fbd4eb5608397c5cd9005         3.23
#> 69          0.99.0 47b59b6a3135a032e02fbd4eb5608397c5cd9005         3.23
#> 70          0.99.0 47b59b6a3135a032e02fbd4eb5608397c5cd9005         3.23
#> 101         0.99.0 47b59b6a3135a032e02fbd4eb5608397c5cd9005         3.23
#> 102         0.99.0 47b59b6a3135a032e02fbd4eb5608397c5cd9005         3.23
#> 103         0.99.0 47b59b6a3135a032e02fbd4eb5608397c5cd9005         3.23
#> 104         0.99.0 47b59b6a3135a032e02fbd4eb5608397c5cd9005         3.23
#> 105         0.99.0 47b59b6a3135a032e02fbd4eb5608397c5cd9005         3.23
#>                generated_at
#> 31  2026-08-04 16:51:07 UTC
#> 32  2026-08-04 16:57:34 UTC
#> 33  2026-08-04 17:07:25 UTC
#> 34  2026-08-04 17:19:56 UTC
#> 35  2026-08-04 17:37:09 UTC
#> 66  2026-08-04 16:52:57 UTC
#> 67  2026-08-04 17:00:59 UTC
#> 68  2026-08-04 17:11:01 UTC
#> 69  2026-08-04 17:26:26 UTC
#> 70  2026-08-04 17:43:23 UTC
#> 101 2026-08-04 16:54:45 UTC
#> 102 2026-08-04 17:04:26 UTC
#> 103 2026-08-04 17:14:44 UTC
#> 104 2026-08-04 17:32:18 UTC
#> 105 2026-08-04 17:49:51 UTC

ps_retrieve_bg() retrieves the background statistics and associates them with the matching JASPAR matrices in a PSMatrixList. Here source = "zenodo" downloads the file from Zenodo directly; by default it is retrieved through ExperimentHub.

background <- ps_retrieve_bg(
    JASPAR_matrix = "Jaspar2020",
    org = "at",
    prom_reg = c(-1000, 0),
    assembly = "TAIR9",
    version = 2,
    source = "zenodo"
)
length(background)
#> [1] 530

This returns 530 matrices from the JASPAR 2020 CORE plant collection. Each matrix includes the background promoter count and the mean and standard deviation of its best-match scores across those promoters.

The selected window extends 1000 bp upstream of the transcription start site (TSS). Foreground promoters must use the same window as the background, together with a compatible assembly and annotation.

From genes to promoters

pscan() accepts promoter sequences as a DNAStringSet. The following procedure constructs the foreground using the background’s promoter definition.

When relevant transcripts are known, their corresponding promoters can be selected directly. Gene-level differential expression results require an additional selection rule because a gene may have multiple transcription start sites. Selecting one promoter per gene avoids weighting genes by their number of annotated transcripts. A deterministic rule makes this choice reproducible, but does not establish which promoter is active under the experimental conditions.

ps_select_promoters() maps a gene list to one promoter per mapped gene using supplied promoter identifiers or sequences and a gene-to-transcript annotation table.

Identifier suffixes have different meanings across annotation systems. The suffix in the RefSeq accession NM_000546.6 denotes a sequence version and is removed for matching. In the TAIR identifier AT1G01110.2, the suffix denotes a transcript variant and is preserved. PscanR detects and reports the scheme used for matching.

The following minimal example illustrates TAIR transcript selection:

promoters <- c("AT1G01110.2", "AT1G01110.1", "AT1G01160.1")
annotation <- data.frame(
    GENE = sub("\\..*$", "", promoters),
    TX = promoters
)

selection <- ps_select_promoters(
    genes = c("AT1G01110", "AT1G01160"),
    promoter_ids = promoters,
    annotation = annotation,
    gene_col = "GENE", transcript_col = "TX",
    mode = "representative"
)
#> Detected transcript scheme: "tair" (3/3 identifiers matched)
#>   AT1G01110.2 -> transcript AT1G01110.2
#>   the trailing .N is a SPLICE VARIANT and is preserved, not stripped

The returned mapping records the selected promoter and the selection rule:

Selected promoters and the rules used for their selection.
gene	promoter_id	scheme	suffix_role	transcript_class	selection_source	n_candidates	decided_by
AT1G01110	AT1G01110.1	tair	variant	primary	TAIR primary model	2	rule
AT1G01160	AT1G01160.1	tair	variant	primary	TAIR primary model	1	rule

For Arabidopsis the rule is to prefer the .1 primary gene model and fall back to the lowest available variant.

This dataset requires promoter sequences matching the background selected above: TAIR9 annotation, 1000 bp upstream of the TSS.

Annotations derived from newer gene models should not be assumed to reproduce TAIR9 promoters, even when transcript identifiers are shared. Changes to annotated transcript ends can shift the TSS and the extracted sequence. This example therefore uses the TAIR9 annotation specified by the background.

The following code illustrates promoter extraction and selection. It is not evaluated during vignette compilation because it requires an external annotation file and the txdbmaker, GenomicFeatures, GenomicRanges, GenomeInfoDb, Biostrings, and BSgenome.Athaliana.TAIR.TAIR9 packages.

TAIR9_GFF3_genes.gff can be downloaded from the TAIR9 GFF3 release directory. Set tair9_gff to the local path of the downloaded file.


down_genes <- gene_sets$gene_id[gene_sets$direction == "down"]

# Use the TAIR9 annotation used to construct the background.

tair9_gff <- "TAIR9_GFF3_genes.gff" # path to the downloaded annotation

txdb <- txdbmaker::makeTxDbFromGFF(
    tair9_gff, format = "gff3",
    dataSource = "TAIR9", organism = "Arabidopsis thaliana"
)
genome <- BSgenome.Athaliana.TAIR.TAIR9::BSgenome.Athaliana.TAIR.TAIR9

# The background covers the five nuclear chromosomes.
tx <- GenomicFeatures::transcripts(txdb, columns = "tx_name")
tx <- tx[as.character(GenomicRanges::seqnames(tx)) %in% paste0("Chr", 1:5)]
names(tx) <- GenomicRanges::mcols(tx)$tx_name

# promoters() takes the region upstream of each 5' end, so it is strand-aware
# and anchored on the transcription start site.
ranges <- GenomicRanges::promoters(tx, upstream = 1000, downstream = 0)

# Exclude out-of-bounds ranges before sequence extraction; getSeq() otherwise
# raises an error. One transcript is affected here.
inside <- GenomicRanges::start(ranges) >= 1L &
    GenomicRanges::end(ranges) <= GenomeInfoDb::seqlengths(genome)[
        as.character(GenomicRanges::seqnames(ranges))
    ]
ranges <- ranges[inside]

promoters_tair9 <- Biostrings::getSeq(genome, ranges)
names(promoters_tair9) <- names(ranges)

# Exclude one additional promoter with more than 50% ambiguous N bases.
promoters_tair9 <- promoters_tair9[
    Biostrings::letterFrequency(promoters_tair9, "N", as.prob = TRUE) <= 0.5
]

annotation <- data.frame(
    GENE = sub("\\..*$", "", names(promoters_tair9)),
    TX = names(promoters_tair9)
)

down_mapping <- ps_select_promoters(
    genes = down_genes,
    promoter_sequences = promoters_tair9,
    annotation = annotation,
    gene_col = "GENE", transcript_col = "TX",
    mode = "representative",
    return = "mapping"
)

down_promoters <- promoters_tair9[down_mapping$promoter_id]

The mapping and scan results are bundled in a precomputed object, allowing the remaining examples to run without repeating promoter extraction or scanning.

precomputed <- readRDS(fixture("arabidopsis_bpc_results.rds"))
down_mapping <- precomputed$down_mapping

ps_selection_summary() reports the numbers of input and mapped genes, including those with multiple candidate promoters.

ps_selection_summary(down_mapping)$overall
#>   scheme genes_requested genes_mapped genes_unmapped genes_with_choice
#> 1   tair            1247         1242              5               210
#>   genes_decided_by_tie_break candidates_excluded_protein
#> 1                          0                           0

Unmapped genes have no eligible promoter in the supplied reference. In this dataset, annotation coverage accounts for several missing genes, as discussed in the appendix. selection_source records the rule used to select each promoter:

ps_selection_summary(down_mapping)$by_source
#>         selection_source promoters
#> 1     TAIR primary model      1239
#> 2 TAIR alternative model         3

A small executable scan

The following example uses twelve genes from the same downregulated list and all their candidate TAIR9 promoters, extracted over -1000/0 from the archived annotation snapshot. Selection preserves splice-variant suffixes. Two motifs keep execution short; this demonstration does not replace the full analysis, and its FDR adjustment covers only those two motifs.

small_example <- readRDS(fixture("promoter_selection_example.rds"))
small_sequences <- ps_select_promoters(
    genes = small_example$genes,
    promoter_sequences = small_example$promoter_sequences,
    annotation = small_example$annotation, gene_col = "GENE",
    transcript_col = "TX", scheme = "tair", mode = "representative",
    return = "sequences"
)
small_background <- precomputed$results_subset[1:2]
small_results <- pscan(small_sequences, small_background)
ps_results_table(small_results)
#>          NAME    BG_AVG   BG_STDEV    FG_AVG   ZSCORE    P.VALUE        FDR
#> MA1404.1 BPC1 0.7246550 0.05827205 0.7596819 2.082244 0.01866008 0.02789806
#> MA1402.1 BPC6 0.7385571 0.06169966 0.7726231 1.912625 0.02789806 0.02789806

Running the scan

pscan() takes foreground promoter sequences and a background PSMatrixList, and returns a PSMatrixList containing the per-promoter results for each motif.

The BPPARAM argument configures execution through BiocParallel:

SerialParam() uses the current R process and simplifies package checks and debugging.
SnowParam() uses independent R workers connected through sockets. It supports Windows and Unix-like systems, with additional startup and data-transfer overhead.
MulticoreParam() uses forked workers on supported Unix-like systems. It can reduce data-transfer overhead, but suitability depends on the execution environment. Forked execution is unavailable on Windows.

The following scan is not evaluated during vignette compilation to limit computation time. It can be run after the preceding data-preparation steps; runtime depends on the hardware and execution strategy.


BPPARAM <- BiocParallel::SerialParam() # Execute in the current R process.

# Portable parallel execution with four independent R workers:
# BPPARAM <- BiocParallel::SnowParam(workers = 4, type = "SOCK")

# Forked execution on supported Unix-like systems:
# BPPARAM <- BiocParallel::MulticoreParam(workers = 4)

results <- pscan(
    down_promoters, background,
    BPPARAM = BPPARAM
)

ps_results_table() converts the result object to a motif enrichment table.

enrichment <- ps_results_table(results)

The following objects contain the saved scan results. Subsequent examples operate on these bundled results.

# Saved enrichment table for downregulated promoters.
enrichment <- precomputed$enrichment_down
dim(enrichment)
#> [1] 530   7

# A motif subset of the PSMatrixList retains the per-promoter results.
results_subset <- precomputed$results_subset
results_subset
#> PSMatrixList of length 14
#> names(14): MA1402.1 MA1404.1 MA1403.1 MA0931.1 ... MA0956.1 MA1081.1 MA0563.1

The subset includes BPC motifs relevant to the original study, the highest-ranked G-box-associated motifs, and two unrelated motifs for comparison.

TFBSTools::name(results_subset)
#> MA1402.1 MA1404.1 MA1403.1 MA0931.1 MA0941.1 MA0570.2 MA0551.1 MA1338.1 
#>   "BPC6"   "BPC1"   "BPC5"   "ABI5"   "ABF2"   "ABF1"    "HY5"  "AREB3" 
#> MA1061.1 MA0561.1 MA0552.1 MA0956.1 MA1081.1 MA0563.1 
#>    "SPT"   "PIF4"   "PIF1"   "BEE2" "WRKY25"   "SEP3"

The enrichment table contains one row per motif, indexed by JASPAR matrix identifier:

head(enrichment[order(-enrichment$ZSCORE), ], 6)
#>           NAME    BG_AVG   BG_STDEV    FG_AVG    ZSCORE      P.VALUE
#> MA0931.1  ABI5 0.8447239 0.06394459 0.8635191 10.354506 1.996220e-25
#> MA0570.2  ABF1 0.8083704 0.05914895 0.8256749 10.306202 3.303015e-25
#> MA0941.1  ABF2 0.8467722 0.05148996 0.8617719 10.262302 5.209659e-25
#> MA0551.1   HY5 0.7509377 0.05883005 0.7675842  9.968032 1.051777e-23
#> MA1061.1   SPT 0.8637590 0.05122061 0.8781784  9.917211 1.752059e-23
#> MA1338.1 AREB3 0.7903007 0.05242647 0.8049847  9.866820 2.898665e-23
#>                   FDR
#> MA0931.1 8.752990e-23
#> MA0570.2 8.752990e-23
#> MA0941.1 9.203730e-23
#> MA0551.1 1.393604e-21
#> MA1061.1 1.857182e-21
#> MA1338.1 2.560487e-21

G-box-associated motifs from the bHLH and bZIP families occupy the highest ranks. ABI5, ABF1, and ABF2 are bZIP factors whose motifs contain an ACGT core.

The background parameters remain associated with the result and can be retrieved from it at any point with ps_get_bg_table():

ps_get_bg_table(results_subset)
#>          BG_SIZE   BG_MEAN   BG_STDEV
#> MA1402.1   36405 0.7385571 0.06169966
#> MA1404.1   36405 0.7246550 0.05827205
#> MA1403.1   36405 0.6158596 0.06436089
#> MA0931.1   36405 0.8447239 0.06394459
#> MA0941.1   36405 0.8467722 0.05148996
#> MA0570.2   36405 0.8083704 0.05914895
#> MA0551.1   36405 0.7509377 0.05883005
#> MA1338.1   36405 0.7903007 0.05242647
#> MA1061.1   36405 0.8637590 0.05122061
#> MA0561.1   36405 0.8575087 0.06116661
#> MA0552.1   36405 0.8159746 0.05672559
#> MA0956.1   36405 0.8939195 0.04176545
#> MA1081.1   36405 0.8917160 0.05135195
#> MA0563.1   36405 0.9071071 0.03229260

BG_SIZE is the number of Arabidopsis promoters the background was computed over; BG_MEAN and BG_STDEV describe how that motif scores across them.

For each motif, the foreground mean score is compared with the background mean using an enrichment z-score. The table reports the corresponding upper-tail p-value and its Benjamini–Hochberg adjustment across the motif collection in the FDR column.

ps_motif_barplot() displays the highest-ranked motifs. The n argument limits their number to at most 20 by default. With group = "class", bar colours indicate transcription factor structural classes extracted from the JASPAR metadata by ps_motif_class().


table(ps_motif_class(results_subset))
#> 
#> Basic helix-loop-helix factors (bHLH)   Basic leucine zipper factors (bZIP) 
#>                                     4                                     5 
#>                      MADS box factors                                 Other 
#>                                     1                                     3 
#>                                  WRKY 
#>                                     1

ps_motif_barplot(results_subset, group = "class")

Within the subset the G-box binders rank above the three BPC motifs, as they do in the full table above.

The function returns a ggplot object that can be customised with a title, theme, or additional layers:

ps_motif_barplot(results_subset, group = "class") +
    ggplot2::labs(title = "Downregulated promoters")

Inspecting the result object

ps_results_table() summarises results across promoters. The PSMatrixList also retains per-promoter results, which can be retrieved using accessor functions.

A PSMatrixList can be indexed by JASPAR matrix identifier; MA0931.1 selects the ABI5 entry:

abi5 <- results_subset[["MA0931.1"]]
TFBSTools::name(abi5)
#> [1] "ABI5"

ps_zscore(abi5)
#>        z 
#> 10.35451
ps_pvalue(abi5)
#> [1] 1.99622e-25
c(background = ps_bg_avg(abi5), foreground = ps_fg_avg(abi5))
#> background foreground 
#>  0.8447239  0.8635191
ps_hits_size(abi5)
#> [1] 1241

ps_hits_table() returns each promoter’s best motif match, including its score, position, strand, and matched oligonucleotide.

head(ps_hits_table(abi5), 5)
#>             SCORE POS STRAND      OLIGO
#> AT2G29490.1     1  17      + TGACACGTGG
#> AT3G49830.1     1 104      + TGACACGTGG
#> AT1G12420.1     1 313      + TGACACGTGG
#> AT3G27870.1     1 387      - CCACGTGTCA
#> AT3G07350.1     1 392      + TGACACGTGG

The ABI5 matches illustrate sequence variation around the G-box core CACGTG. Oligonucleotides are reported in the forward orientation of the supplied promoter sequence, including matches assigned to the reverse strand.

Motif score patterns across promoters

ps_z_table() returns a matrix of standardised best-hit scores, with promoters as rows and motifs as columns. These per-promoter values differ from the motif-level enrichment z-scores, which compare foreground and background means.

z <- ps_z_table(results_subset)
dim(z)
#> [1] 1241   14
round(z[1:4, 1:5], 2)
#>             MA1402.1 MA1404.1 MA1403.1 MA0931.1 MA0941.1
#> AT2G43020.1    -1.15    -1.48    -0.77     2.01     2.15
#> AT1G29951.1     1.32     1.24     0.66    -0.86    -1.00
#> AT2G37450.1     1.65     1.11     0.42     0.73     0.11
#> AT5G39280.1    -1.08    -1.25    -1.27     0.20     0.50

ps_zscore_heatmap() displays these values to compare motif score patterns across promoters. By default, it scales each row and clusters both axes using correlation distances.

ps_zscore_heatmap(
    results_subset, FDR = 0.01,
    show_rownames = FALSE, treeheight_row = 0
)

The column dendrogram groups the three BPC motifs, consistent with their similar GA-rich sequence preferences. A separate group contains G-box-associated motifs, including HY5, ABF/AREB bZIPs, and PIF and SPT bHLHs. The row-scaled heatmap suggests that promoters with relatively high BPC scores differ from those with relatively high G-box motif scores. This heterogeneity is not apparent from the motif-level enrichment table alone.

Positional distribution of motif matches

ps_hitpos_map() displays the position of each promoter’s best motif match.

ps_hitpos_map(
    results_subset, FDR = 0.01, shift = -1000,
    show_rownames = FALSE, treeheight_row = 0
)

shift = -1000 relabels the axis so that coordinates are relative to the transcription start site, matching how the promoters were extracted.

ps_density_plot() estimates the positional distribution of the best hits for a single motif.

ps_density_plot(
    abi5,
    shift = -1000, st = "strict",
    window = c(-1000, 0), bins = 20
)

Higher-scoring ABI5 matches are concentrated near the TSS, with the estimated density peaking within approximately 100 bp upstream.

The window argument specifies the promoter interval and enables reflection boundary correction for the kernel density estimate. The boundaries correspond to valid motif start positions, accounting for motif width.

The bins argument adds a histogram of the same positions on the density scale, allowing comparison with the smoothed curve.

The BPC6 motif, with its GA-rich sequence preference, provides a comparison with the ACGT-containing ABI5 motif.

bpc6_hits <- results_subset[["MA1402.1"]]

ps_density_plot(
    bpc6_hits,
    shift = -1000, st = "strict",
    window = c(-1000, 0), bins = 20
)

Both profiles show concentration near the TSS. The ABI5 density peaks approximately 75–100 bp upstream and then decreases, whereas the BPC6 profile increases towards the final bin. BPC6 has width 21 bp, so the last valid motif start lies 20 bp before the end of the window.

The st argument selects the displayed hits: "all" includes every promoter, "loose" selects scores at or above the background mean, and "strict" selects scores at least one background standard deviation above that mean. This can highlight positional patterns among stronger matches. These display filters do not alter the motif enrichment statistics.

ps_hit_score_plot() plots best-hit score against position for each promoter and shows the st score thresholds as reference lines.

ps_hit_score_plot(
    results_subset[c("MA1402.1", "MA0931.1")],
    shift = -1000
)

Higher-scoring BPC6 hits are concentrated within the final few hundred bases, whereas lower-scoring hits are distributed across the window. Higher-scoring ABI5 hits show a broader distribution and generally occur farther upstream. The joint display distinguishes positional concentration from match strength.

Analysing a motif-selected promoter subset

The preceding analysis considers all downregulated promoters. A complementary analysis restricts the foreground to promoters containing a high-scoring match to a selected motif.

pscan_filtered() scans the promoters with one chosen matrix, keeps those scoring at least n standard deviations above that motif’s background mean, and rescans the retained promoters against the supplied collection.

The following code can be run after the preceding data-preparation steps. It is not evaluated during vignette compilation to limit computation time.

bpc6 <- background[[which(TFBSTools::ID(background) == "MA1402.1")]]

filtered <- pscan_filtered(
    down_promoters, bpc6, n = 1,
    background = background,
    BPPARAM = BiocParallel::SerialParam()
)
ps_results_table(filtered)

Filtering by BPC6 score retains 234 of the 1242 downregulated promoters.

Enrichment within the 234 promoters with a high-scoring BPC6 match.
motif_id	NAME	ZSCORE	P.VALUE
MA1402.1	BPC6	33.87	8.73e-252
MA1404.1	BPC1	30.90	6.87e-210
MA1403.1	BPC5	29.32	2.62e-189
MA1416.1	RAMOSA1	28.99	4.22e-185
MA1381.1	AT3G46070	8.75	1.07e-18
MA1372.1	STZ	4.85	6.05e-07
MA1274.1	OBP3	4.67	1.54e-06
MA1278.1	OBP1	4.12	1.91e-05
MA1268.1	AT1G69570	4.01	2.98e-05
MA0984.1	DOF5.7	3.76	8.47e-05

The GA-rich BPC and RAMOSA1 motifs have the highest ranks after selection by BPC6 score. G-box-associated bZIP and bHLH motifs are absent from the leading results shown here, which instead include C2H2 and Dof motifs. This change is consistent with different sequence properties in the BPC6-selected subset.

Filtering provides a way to examine subsets of a heterogeneous foreground. Because selection is based on a motif score, enrichment of that motif and similar motifs is expected and should not be treated as independent evidence of regulation.

Summary

The core workflow comprises four steps: select a background, map genes to promoters with a matching definition, scan the sequences, and inspect the motif enrichment table.

background <- ps_retrieve_bg(
    "Jaspar2020", "at", c(-1000, 0), "TAIR9", version = 2
)
down_promoters <- ps_select_promoters(
    genes = down_genes, promoter_sequences = promoters_tair9,
    annotation = annotation, gene_col = "GENE", transcript_col = "TX",
    mode = "representative", return = "sequences"
)
results <- pscan(down_promoters, background)
ps_results_table(results)

Additional analyses examine per-promoter scores (ps_z_table()), match positions and sequences (ps_hits_table(), ps_hitpos_map(), ps_density_plot(), ps_hit_score_plot()), and promoter subsets selected by motif score (pscan_filtered()). These complementary summaries support the interpretation of candidate regulatory patterns.

Appendix: Comparison with the published analysis

The source study reports Pscan results that permit motif-level comparison. Using the comparison criterion described below, this analysis recovers 164 of the 168 motifs listed for the published downregulated set. The agreement is substantial but incomplete. The following methodological details affect the comparison.

The promoter window. The published analysis used −1000/+100 relative to the TSS. The bundled catalogue does not include that window, so this analysis uses −1000/0. A separate analysis of the same genes with the −950/+50 background gave a Spearman correlation of ρ = 0.99 across all 530 motif z-scores, although some significance classifications changed. This comparison supports similarity between the two tested windows; it does not directly measure the effect of adding the full 100 bp downstream interval used in the publication.

The annotation release. The published promoters were defined on TAIR10 gene models; this background uses TAIR9. Shared genome coordinates do not guarantee identical transcript models or TSS annotations. In the recorded reference comparison, 2244 of the study’s 2253 genes have promoters in the TAIR9 universe (99.6%). The nine missing genes include loci absent from the older annotation and a mitochondrial gene outside the nuclear reference set; they occur in both expression directions. High gene coverage does not establish promoter-level identity between annotation releases.

Multiple testing correction. PscanR reports nominal p-values and Benjamini–Hochberg adjusted p-values in the FDR column. The published methods specify Benjamini–Hochberg correction with a 0.01 threshold. However, inspection of the deposited Supplementary Table 3 in the local reproduction found adjusted values approximately equal to nominal p-values multiplied by 530, with listed values below 0.1. The overlap of 164 out of 168 therefore refers to the reproduction’s Bonferroni criterion, p * 530 < 0.1, chosen to match those deposited values; it is not the count obtained from the default PscanR FDR column.

Data provenance

The installed inst/scripts/fixture-provenance.md describes the sources, selection rules, serialization, and reproduction recipes for the bundled example data. inst/DATA_SOURCES.md records attribution and reuse terms. The large published analyses use saved results; small workflow examples are executed during the build. Reproducing a full analysis requires the matching annotation and promoter window, not just the same gene identifiers.

Session information

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] stats     graphics  grDevices utils     datasets  methods   base     
#> 
#> other attached packages:
#> [1] PscanR_0.99.0    BiocStyle_2.41.0
#> 
#> loaded via a namespace (and not attached):
#>  [1] tidyselect_1.2.1            dplyr_1.2.1                
#>  [3] farver_2.1.2                blob_1.3.0                 
#>  [5] filelock_1.0.3              Biostrings_2.81.9          
#>  [7] S7_0.2.2                    bitops_1.1-0               
#>  [9] fastmap_1.2.0               RCurl_1.98-1.20            
#> [11] BiocFileCache_3.3.0         GenomicAlignments_1.49.2   
#> [13] XML_3.99-0.25               digest_0.6.39              
#> [15] DirichletMultinomial_1.55.1 lifecycle_1.0.5            
#> [17] pwalign_1.9.1               RSQLite_3.53.3             
#> [19] magrittr_2.0.5              compiler_4.6.1             
#> [21] rlang_1.3.0                 sass_0.4.10                
#> [23] tools_4.6.1                 yaml_2.3.12                
#> [25] rtracklayer_1.73.0          knitr_1.52                 
#> [27] labeling_0.4.3              bit_4.6.0                  
#> [29] curl_8.0.0                  RColorBrewer_1.1-3         
#> [31] BiocParallel_1.47.0         withr_3.0.3                
#> [33] purrr_1.2.2                 BiocGenerics_0.59.12       
#> [35] sys_3.4.3                   grid_4.6.1                 
#> [37] stats4_4.6.1                caTools_1.18.4             
#> [39] ggplot2_4.0.3               scales_1.4.0               
#> [41] gtools_3.9.5                cli_3.6.6                  
#> [43] rmarkdown_2.32              crayon_1.5.3               
#> [45] generics_0.1.4              otel_0.2.0                 
#> [47] TFMPvalue_1.0.0             httr_1.4.9                 
#> [49] rjson_0.2.23                BiocBaseUtils_1.15.1       
#> [51] DBI_1.3.0                   cachem_1.1.0               
#> [53] parallel_4.6.1              BiocManager_1.30.27        
#> [55] XVector_0.53.0              restfulr_0.0.17            
#> [57] matrixStats_1.5.0           vctrs_0.7.3                
#> [59] Matrix_1.7-6                jsonlite_2.0.0             
#> [61] IRanges_2.47.5              S4Vectors_0.51.10          
#> [63] bit64_4.8.6                 seqLogo_1.79.0             
#> [65] maketools_1.3.2             jquerylib_0.1.4            
#> [67] glue_1.8.1                  TFBSTools_1.51.0           
#> [69] codetools_0.2-20            gtable_0.3.6               
#> [71] BiocIO_1.23.3               GenomicRanges_1.65.4       
#> [73] tibble_3.3.1                pillar_1.11.1              
#> [75] htmltools_0.5.9             Seqinfo_1.3.2              
#> [77] JASPAR2020_0.99.10          BSgenome_1.81.1            
#> [79] httr2_1.3.0                 dbplyr_2.6.0               
#> [81] R6_2.6.1                    evaluate_1.0.5             
#> [83] Biobase_2.73.2              lattice_0.23-1             
#> [85] Rsamtools_2.29.0            cigarillo_1.3.1            
#> [87] pheatmap_1.0.13             memoise_2.0.1              
#> [89] bslib_0.12.0                Rcpp_1.1.2                 
#> [91] xfun_0.61                   buildtools_1.0.0           
#> [93] pkgconfig_2.0.3

Motif enrichment analysis with PscanR: Arabidopsis BPC target promoters