--- title: "seqArchRplus facilitates downstream analysis of clusters of promoter sequence architectures" author: - name: Sarvesh Nikumbh affiliation: MRC London Institute of Medical Sciences and Imperial College London date: "`r Sys.Date()`" package: seqArchRplus output: BiocStyle::html_document: toc: true BiocStyle::pdf_document: default bibliography: seqArchRplus.bib vignette: > %\VignetteIndexEntry{seqArchRplus facilitates downstream analysis of clusters of promoter sequence architectures} %\VignetteEngine{knitr::rmarkdown} \usepackage[utf8]{inputenc} --- ```{r setup, include=FALSE} knitr::opts_chunk$set( collapse = TRUE, comment = "#>" ) ``` # Introduction Analysis of the promoterome of any organism entails analyzing data from an experiment like CAGE (Cap analysis of gene expression) [@Kodzius:2006] that provides information on genome-wide transcription start sites at single nucleotide resolution. These promoters can then be further studied to identify different classes [@Carninci:2006] among them based on different attributes, for instance, their shape (broad or focused/sharp promoters), gene function (tissue specific vs. housekeeping) etc. These different promoter classes harbour a variety of promoter architectures orchestrated by different proteins together with sequence elements at near-fixed positions in the promoter that determine the position of the transcription start site. The different promoter architectures are known to be used differentially by genes either in different tissues or at different developmental timepoints. For instance, [@haberle2014shifting] have shown a dynamic interchange of promoter architecture within the same genomic region between maternal and zygotic stages of zebrafish development. Identifying and studying these different promoter architectures further is thus vital. While the R/Bioconductor package `r Biocpkg("seqArchR")` enables de novo identification of clusters of promoter sequence architectures, this package, `r Biocpkg("seqArchRplus")`, enables performing various steps in their downstream bioinformatic analyses and produce publication-ready plots (building on various other Bioconductor packages). The many steps in the downstream analyses of promoter sequence architectures enabled by `r Biocpkg("seqArchRplus")` are: - Curate the final set of clusters from `r Biocpkg("seqArchR")`. See the [`r Biocpkg("seqArchR")` vignette](http://bioconductor.org/packages/release/bioc/vignettes/seqArchR/inst/doc/seqArchR.html) , or the preprint to understand in detail why this may be required - Order the sequence architectures by the interquantile widths (IQWs) of the tag clusters (_aka_ promoter shape). See [`CAGEr` vignette](https://www.bioconductor.org/packages/release/bioc/vignettes/CAGEr/inst/doc/CAGEexp.html) for more information on tag clusters and their IQWs - Visualize distributions of IQW, TPM (Tags per million) and conservation scores (or other when available) per cluster - Visualize the percentages of annotations for genome-wide CAGE-derived transcription tag clusters for each architecture-based clusters - Ease of comparison across samples/stages: Visualize the above plots as (publication ready) combined panels - Many per cluster visualizations including: - sequence logos of cluster architectures - strand-separated sequence logos of architectures - distributions of promoters on different chromosomes and strands - GO terms enriched for each cluster/architecture - Produce BED track files of `r Biocpkg("seqArchR")` clusters for visualization in a genome browser or IGV Examples of most of these capabilities are provided in this vignette. We hope that using `r Biocpkg("seqArchRplus")` (together with `r Biocpkg("seqArchR")`) will be useful in swiftly, but comprehensively, analyzing promoters identified using CAGE for any organism, and enable insights and hypotheses generation with ease. So far, we have tested `r Biocpkg("seqArchR")` and `r Biocpkg("seqArchRplus")` on promoters in drosophila from [@schor2017promoter] and modENCODE [@chen2014comparative], zebrafish [@nepal2013dynamic], mice (Fantom [@Consortium:2014hz]), humans (ENCODE), and also in plants (barley and maize). We already have a plan for a feature in the future version of `r Biocpkg("seqArchRplus")`: generate HTML reports that help you navigate this wealth of information with relative ease. # Installation The latest stable version of seqArchjRplus can be installed from Bioconductor as shown below. ```{r seqArchRplus-install, echo=TRUE, eval=FALSE} ## When available on Bioconductor if (!require("BiocManager", quietly = TRUE)) { install.packages("BiocManager") } BiocManager::install("seqArchRplus") ``` In case of any errors or help required, please consider looking up: [https://github.com/snikumbh/seqArchRplus](https://github.com/snikumbh/seqArchRplus) and file a [new issue](https://github.com/snikumbh/seqArchRplus/issues/new). # seqArchRplus for downstream analysis of promoter sequences ```{r setup-two, echo=TRUE,include=TRUE,results="hide",message=FALSE,warning=FALSE} # Load seqArchRplus library(seqArchRplus) library(seqArchR) library(Biostrings) library(TxDb.Dmelanogaster.UCSC.dm6.ensGene) library(org.Dm.eg.db) library(ChIPseeker) # Set seed for reproducibility set.seed(1234) ``` ## Setting up To begin using `r Biocpkg("seqArchRplus")`, you require a set of promoters sequences (as a DNAStringSet object) and their clustering information (as a simple list storing the sequence IDs belonging to each cluster). With this, you can already visualize the cluster-wise sequence logos, distribution of chromosome and strand locations, GO term-enrichments, motif occurrence heatmaps. When you have additional information such as the IQWs (shape information), `r Biocpkg("seqArchRplus")` can use these to order the clusters in the visualizations. This visualization can be supplemented with the distribution of per-cluster TPM values and/or conservation scores etc. when provided. If your workflow involves CAGEr (for pre-processing raw CAGE data) and `r Biocpkg("seqArchR")` (for clusters), they can be seamlessly used to utilize the full scope of `r Biocpkg("seqArchRplus")`. In this case, it will be helpful if you have - the CAGEr object or information on the tag clusters from [https://bioconductor.org/packages/release/bioc/html/CAGEr.html](CAGEr), specifically, the width of the tag clusters, total TPM value of the cluster and that of the dominant CTSS in the cluster - the `r Biocpkg("seqArchR")` result object. ### Input data In this vignette, we use the example data provided with the package. The raw CAGE data for different developmental timepoints in Drosophila melanogaster (Schor et al. 2017) were pre-processed with CAGEr to obtain promoter sequences. These were then processed with `r Biocpkg("seqArchR")` to cluster the promoter DNA sequences. Only a subset of the complete data is provided with the package and used here to enable demonstration of package utility. In particular, the data included is for the timepoint 10-12h after egg laying. Specifically, the following files are provided: 1. promoter sequences stored as a gzipped FASTA file: two kinds, with 45 and 200 bp flanks around the TSS (filenames: `example_promoters45.fa.gz` and `example_promoters200.fa.gz`) 2. clustering information (filename: `example_clust_info.rds`) 3. GRanges object obtained from CAGEr which holds information on the tag clusters (promoters) with additional information such as the position of dominant TSS, TPM value for the dominant TSS, position of quantile boundaries (filename: `example_tc_gr.rds`) 4. BED file with same information as the GRanges object (filename: `example_info_df.rds`) 5. `r Biocpkg("seqArchR")` result object (filename: `seqArchR_result.rds`) ```{r prepare-fetch-data, eval=TRUE} ## 1. Raw DNA sequences raw_seqs <- Biostrings::readDNAStringSet( filepath = system.file("extdata", "example_promoters45.fa.gz", package = "seqArchRplus", mustWork = TRUE ) ) ## 2. Clustering information (arbitrary order/unordered) unord_clusts <- readRDS(system.file("extdata", "example_clust_info.rds", package = "seqArchRplus", mustWork = TRUE )) ## 3. GRanges object tc_gr <- readRDS(system.file("extdata", "example_tc_gr.rds", package = "seqArchRplus", mustWork = TRUE )) ## 4. BED file bed_info_fname <- system.file("extdata", "example_info_df.bed.gz", package = "seqArchRplus", mustWork = TRUE ) info_df <- read.delim(file = bed_info_fname, sep = "\t", header = TRUE) ## 5. seqArchR result seqArchR_result <- readRDS(system.file("extdata", "seqArchR_result.rds", package = "seqArchRplus", mustWork = TRUE)) ## **NOTE** Only for the example seqArchR result object provided with this ## package: ## ## Any seqArchR result object has the raw DNA sequences as a part of it. ## See details here: https://snikumbh.github.io/seqArchR/reference/seqArchR.html ## ## But, in order to reduce the size of the example seqArchR result object ## provided with this (seqArchRplus) package, the `rawSeqs` element of the ## result is missing. This is reassigned below for the purposes of this ## vignette. ## ## The seqArchR result was obtained by processing the DNA sequences in ## `example_promoters45.fa.gz`. Thus we reassign them to ## `seqArchR_result$rawSeqs` here ## if(!("rawSeqs" %in% names(seqArchR_result))) seqArchR_result$rawSeqs <- raw_seqs ``` If you are not using the result from `r Biocpkg("seqArchR")` or are already using a curated set of clusters, you can skip the following subsection and jump to section \@ref(generating-plots) for a demonstration of generating all plots. ### Curating raw clusters from `seqArchR` result {#curate-seqArchR-result} Raw clusters from `r Biocpkg("seqArchR")` result (say, the final iteration) often require curation. The `r Biocpkg("seqArchR")` result does have a final `clustSol` (the clustering solution) where the clusters from the final iteration are collated. `r Biocpkg("seqArchR")` uses hierarchical clustering for this purpose. However, hierarchical clustering, with a chosen agglomeration and distance method, does not necessarily output the best suitable collated set of clusters. Therefore, some small amount of manual curation (re-assignments) may be required to reach the ideal collated set of clusters as the final solution. This is achieved with the help of the `seqArchRplus::curate_clusters` function. The basic idea of the curate_clusters function is available in `help(seqArchRplus::curate_clusters)`. A more elaborate description is as follows. It takes as input an agglomeration method and also a distance method (see more arguments in \code{help(curate_clusters)}). On the first call to the function, a plot with the associated dendrogram resulting from the hierarchical clustering is shown together with the per cluster sequence logos. This should help in identifying if the chosen agglomeration and distance methods worked well (i.e., are clusters with similar sequence logos also together in the dendrogram?). If this is satisfactory (i.e., chosen agglomeration and distance method have done well already and only few curations will be required), count the number of overall clusters that you can visually see. Generally, if (estimated) 3-4 (even 5-6) clusters out of 30-40 require curation, I would consider it as a satisfactory outcome from the hierarchical clustering output. Now, in the second call to `curate_clusters`, specify the number of clusters based on your count. View the resulting clustering in the new plot produced as output. This time the dendrogram shows this clustering with colors and grey-scaled boxes drawn around the clusters. Now you can exactly note which clusters need curation. Imagine a hypothetical scenario where a total of 16 raw clusters are obtained from the final `r Biocpkg("seqArchR")` iteration. Collating them into 4 clusters using hierarchical clustering results in: - collated cluster 1: 1, 7, 9 - collated cluster 2: 2, 3, 4, 8, 10 - collated cluster 3: 11, 12, 13 - collated cluster 4: 5, 6, 14, 15, 16 But, due to very minor differences in the sequence logos, you may want to move raw cluster 9 from collated cluster 1 to collated cluster 3 (that of raw clusters 11-13). To do this, set the `need_change` argument of `curate_clusters` to `c(9)` and `change_to` argument to `c(11)`. Any one of the destination cluster members can be specified in `change_to`. Similarly, any other such curations can be added. Specify all such curations together in one go like this: - `need_change <- list(c(9), c(4, 7), c(16))` and - `change_to <- list(11, 14, 0)` Here, (a) raw cluster 9 is moved to the collated cluster containing raw cluster 11; (b) 4 and 7 are moved to that of 14; and (c) 16 is moved to a totally new cluster of itself. The result is the re-assigned clusters. Below, this procedure of curation is re-iterated with the help of actual function calls and resulting figures using a reduced example data available with the package. ```{r curation-call-1, echo=TRUE, eval=TRUE, warning=FALSE, fig.cap="Figure showing a combined panel of dendrogram + sequence logos at the first step of curation", fig.wide=TRUE, fig.height=3.5} sn <- "RAL28_10_to_12h" use_aggl <- "complete" use_dist <- "euclid" ## get seqArchR clusters custom curated seqArchR_curated_clusts <- curate_clusters(sname = sn, use_aggl = use_aggl, use_dist = use_dist, seqArchR_result = seqArchR_result, iter = 5, pos_lab = NULL, regularize = TRUE, topn = 50, use_cutk = 2, final = FALSE, dir_path = NULL ) seqArchR_curated_clusts$curation_plot ``` The motive of Figure \@ref(fig:curation-call-1) is to look at the clusters ordered by the hierarchical clustering step. This way, counting the tentative number of clusters (K) is easier. In Figure \@ref(fig:curation-call-1), we can see that the first and the second clusters can be collated together. Therefore, to demonstrate curation, we will set K=5 clusters and perform this minor curation as shown below. ```{r curation-set-k, echo=TRUE, eval=TRUE, warning=FALSE, fig.cap="Figure showing a combined panel of dendrogram + sequence logos at the second step of curation", fig.wide=TRUE, fig.height=3.5} ## Let us set K=5 for this example, and combine clusters 1 and 2 into one. set_cutk <- 5 ## Form the lists `need_change` and `change_to` for re-assignments need_change <- list(c(2)) change_to <- list(c(1)) seqArchR_curated_clusts <- curate_clusters(sname = sn, use_aggl = use_aggl, use_dist = use_dist, seqArchR_result = seqArchR_result, iter = 5, regularize = TRUE, topn = 50, use_cutk = set_cutk, #*** need_change = need_change, change_to = change_to, final = FALSE, #*** dir_path = NULL ) seqArchR_curated_clusts$curation_plot ``` Notice that each dendrogram branch (and leaf node) has been assigned a different color because we set K=5. We are now ready to make the final call to `curate_clusters` and obtain the curated clusters' list. ```{r curation-final-call, echo=TRUE, eval=TRUE, warning=FALSE, fig.cap="Figure showing three panels combined: dendrogram + original sequence logos + collated cluster sequence logos at the final step of curation", fig.wide=TRUE, fig.height=4, fig.width=15} ## Satisfied with the re-assignments? now set final = TRUE seqArchR_curated_clusts <- curate_clusters(sname = sn, use_aggl = use_aggl, use_dist = use_dist, seqArchR_result = seqArchR_result, iter = 5, pos_lab = NULL, regularize = FALSE, topn = 50, use_cutk = set_cutk, #*** need_change = need_change, change_to = change_to, final = TRUE, #*** dir_path = NULL ) seqArchR_curated_clusts$curation_plot ``` Notice that the re-assigned clusters (1 and 2) now have different colors than the corresponding dendrogram branches. ```{r show-clusts} str(seqArchR_curated_clusts$clusters_list) ``` This fetches us clusters in _arbitrary_ order, i.e., not related to the IQWs of the clusters. [See below for the function call `order_clusters_iqw()` that orders these clusters by their median/mean interquantile widths (section \@ref(iqw-ord)).] This curated set of clusters can now be used for further downstream analyses. ## Generating individual plots via individual functions {#generating-plots} The following subsections demonstrate how the `r Biocpkg("seqArchRplus")` functions can be used to generate the various visualizations. ### Visualize architectures for each promoter cluster as sequence logos One can visualize the architecture for each cluster of promoters as a sequence logo using the function `per_cluster_seqlogos`. When the argument `one_plot` is set to TRUE, the function returns a single plot where sequence logos for all clusters are arranged one below the other. When set to FALSE, a list of sequence logo plots is returned instead. In the later case, each plot has a title with information on the number of sequences in the cluster. ```{r seqlogos-one-plot, eval=TRUE, attr.source='.numberLines', warning=FALSE} ## seqlogos_oneplot_pl <- per_cluster_seqlogos( sname = "RAL28_10_to_12h", seqs = raw_seqs, clusts = unord_clusts, pos_lab = -45:45, bits_yax = "auto", strand_sep = FALSE, one_plot = TRUE, txt_size = 12, dir_path = NULL ) seqlogos_oneplot_pl ``` The sequence logos can also be obtained as a list instead of already combined into one plot by simply setting the argument `one_plot = FALSE`. This gives you liberty to use them (either individually or otherwise) together with any other plots as is suitable. ```{r seqlogos-list, eval=FALSE, attr.source='.numberLines', warning=FALSE} ## Obtain the sequence logos as a list for combining later by setting the ## 'one_plot' argument to FALSE seqlogos_list_pl <- per_cluster_seqlogos( sname = "RAL28_10_to_12h", seqs = raw_seqs, clusts = unord_clusts, pos_lab = -45:45, bits_yax = "auto", strand_sep = FALSE, one_plot = FALSE, dir_path = NULL, txt_size = 12 ) ``` Obtaining strand-separated sequence logos is demonstrated in section \@ref(strand-seqlogos). ### Visualize the distribution of IQW and TPM values for each cluster {#iqw-ord} The function `iqw_tpm_plots` orders the clusters by their IQW values (when quantitative information on promoter shape is available, for instance, the interquantile widths of the CAGEr tag clusters) and visualizes the distribution of the IQW and TPM values for each cluster as boxplots. The input clusters are ordered from sharp (top) to broad (bottom). ```{r iqw-tpm-plot, eval=TRUE, attr.source='.numberLines'} ## iqw_tpm_pl <- iqw_tpm_plots(sname = "RAL28_10_to_12h", dir_path = NULL, info_df = info_df, iqw = TRUE, tpm = TRUE, cons = FALSE, clusts = unord_clusts, txt_size = 15 ) iqw_tpm_pl ``` ### Fetch clusters ordered by promoter shape The cluster order visualized in the IQW_TPM plots above can be fetched using the function `order_cluster_iqw`. ```{r order-clusters, eval=TRUE, attr.source='.numberLines'} ## get clusters ordered by median IQW values seqArchR_clusts_ord <- order_clusters_iqw(sname = "RAL28_10_to_12h", clusts = unord_clusts, info_df = info_df, order_by_median = TRUE ) str(unord_clusts) str(seqArchR_clusts_ord) ``` ### Visualize sequence clusters as images colored by the nucleotides The function `seqs_acgt_image` enables visualizing a set of sequences as an image colored by the nucleotides. Here, one can choose to order sequences by the clusters (see argument `seqs_ord`). ```{r acgt-image, eval=TRUE, attr.source='.numberLines'} ## seqs_acgt_image(sname = "RAL28_10_to_12h", seqs = raw_seqs, seqs_ord = unlist(seqArchR_clusts_ord), pos_lab = -45:45, dir_path = NULL ) ``` ### Visualize motif occurrences as heatmaps Dinucleotides are often deemed important within promoter sequences. For instance, the common dinucleotides at TSSs are PyPu (a pyrimidine-purine pair). Similarly, W-boxes (WW) that occur at ~30 bp upstream of the TSS are an important TSS-determining feature in Zebrafish maternal stage promoters. Occurrences of various dinucleotides can be visualized as heatmaps using the function `plot_motif_heatmaps` (based on `r Biocpkg("heatmaps")`)) or `plot_motif_heatmaps2` (based on `r Biocpkg("seqPattern")`. We demonstrate only `plot_motif_heatmaps` here. ```{r motif-heatmaps, eval=TRUE, attr.source='.numberLines'} ## Get larger flank FASTA sequences (larger than those used for seqArchR) use_seqs <- Biostrings::readDNAStringSet(filepath = system.file("extdata", "example_promoters200.fa.gz", package = "seqArchRplus", mustWork = TRUE ) ) plot_motif_heatmaps(sname = "RAL28_10_to_12h", seqs = use_seqs, flanks = c(200), clusts = seqArchR_clusts_ord, motifs = c("WW", "SS"), dir_path = NULL, fheight = 80, fwidth = 160 ) ``` ### Per cluster annotations' plot In each cluster, the promoters (more specifically, the transcription start sites) can be annotated by their genomic feature, i.e., whether the start site is located within a promoter region of a gene, the 5' UTR, 3' UTR, exon, intron or is intergenic. `r Biocpkg("seqArchRplus")` enables visualizing these annotation proportions per cluster as stacked bar plots via the function `per_cluster_annotations`. When the argument `one_plot` is set to TRUE, this function returns all stacked bars combined within a single plot (see example below). Otherwise, a list of individual/per-cluster plots is returned. ```{r annotation-one-plot, eval=TRUE, attr.source='.numberLines'} ## annotations_oneplot_pl <- per_cluster_annotations( sname = "RAL28_10_to_12h", clusts = seqArchR_clusts_ord, tc_gr = tc_gr, cager_obj = NULL, qLow = 0.1, qUp = 0.9, txdb_obj = TxDb.Dmelanogaster.UCSC.dm6.ensGene, tss_region = c(-500, 100), orgdb_obj = NULL, dir_path = NULL, one_plot = TRUE, txt_size = 12 ) annotations_oneplot_pl ## Obtain the per cluster annotations as a list for combining later by setting ## the 'one_plot' argument to FALSE # annotations_list_pl <- per_cluster_annotations( # sname = "RAL28_10_to_12h", # clusts = unord_clusts, # tc_gr = tc_gr, # cager_obj = NULL, # qLow = 0.1, qUp = 0.9, # txdb_obj = TxDb.Dmelanogaster.UCSC.dm6.ensGene, # tss_region = c(-500, 100), # orgdb_obj = NULL, dir_path = NULL, # one_plot = FALSE, # txt_size = 12 # ) ``` ### Per cluster strand-wise sequence logos {#strand-seqlogos} Similar to the sequence logos of all promoters in a cluster (section XY), one can separate promoters in a cluster that are on the positive strand vs. the negative strand, and visualize them separately as shown below. Note that this is the same function `per_cluster_seqlogos` with the argument `strand_sep` set to TRUE. ```{r stranded-seqlogos, eval=TRUE, attr.source='.numberLines', fig.height=10, warning=FALSE} ## Obtain strand-separated sequence logos corresponding to each cluster stranded_seqlogos_pl <- per_cluster_seqlogos( sname = "RAL28_10_to_12h", seqs = raw_seqs, clusts = seqArchR_clusts_ord, pos_lab = -45:45, bits_yax = "auto", info_df = info_df, strand_sep = TRUE, #** one_plot = FALSE, #** dir_path = NULL, txt_size = 12 ) one_plot <- cowplot::plot_grid(plotlist = stranded_seqlogos_pl, ncol = 1) one_plot ``` ### Per cluster chromosome- and strand-wise distribution of promoters The `r Biocpkg("seqArchRplus")` function `per_cluster_strand_dist` can enable visualizing as barplots how promoters in each cluster are distributed on different chromosomes and strands. ```{r stranded-chrom-prop, eval=TRUE, attr.source='.numberLines', fig.height=7} pair_colrs <- RColorBrewer::brewer.pal(n = 5, name = "Set3") ## Obtain strand-separated sequence logos corresponding to each cluster stranded_dist_pl <- per_cluster_strand_dist( sname = "RAL28_10_to_12h", clusts = seqArchR_clusts_ord, info_df = info_df, dir_path = NULL, txt_size = 12 ) one_plot <- cowplot::plot_grid(plotlist = stranded_dist_pl, ncol = 1) one_plot ``` ### Per cluster GO term enrichments GO terms enriched per cluster can be obtained using the `r Biocpkg("seqArchRplus")` `per_cluster_go_term_enrichments()`. This function returns a list of plots showing the GO terms enriched. ```{r go-enrichment, eval=TRUE, attr.source='.numberLines', fig.width=15, fig.height=15} go_enrich_pl <- per_cluster_go_term_enrichments( sname = "RAL28_10_to_12h", clusts = seqArchR_clusts_ord, tc_gr = tc_gr, txdb_obj = TxDb.Dmelanogaster.UCSC.dm6.ensGene, dir_path = NULL, one_file = FALSE, #*** tss_region = c(-500,100), orgdb_obj = "org.Dm.eg.db" ) one_plot <- cowplot::plot_grid(plotlist = go_enrich_pl, ncol = 2) one_plot ``` ### Large panels combining multiple plots ```{r large-panels, eval=TRUE, attr.source='.numberLines', fig.wide=TRUE, fig.height=8, fig.width=20, fig.cap="Multiple plots combined as panels to geenrate one large figure."} ## Remove y-axis text/tick labels annot_new_pl <- annotations_oneplot_pl + ggplot2::ylab(NULL) + ggplot2::theme(axis.text.y = ggplot2::element_blank()) panel_pl <- cowplot::plot_grid(iqw_tpm_pl, seqlogos_oneplot_pl, annot_new_pl, ncol = 3, axis = "tblh", align = "h") panel_pl ``` In a multi-sample study where the samples designate either different timepoints or cell types, all such (large) figures can be generated for all sampleas to produce a list of such plots. ### Clusters as track BED files for visualization in the genome browser We recommend viewing the promoters in the genome browser or IGV to observe additional details such as the position of the dominant TSS etc. The relevant `r Biocpkg("seqArchRplus")` function that writes to disk this information as browser track BED files is `write_seqArchR_cluster_track_bed`. ```{r track-bed, eval=TRUE, attr.source='.numberLines'} ## write_seqArchR_cluster_track_bed( sname = "RAL28_10_to_12h", clusts = seqArchR_clusts_ord, info_df = info_df, use_q_bound = FALSE, one_zip_all = TRUE, org_name = "Dmelanogaster", dir_path = tempdir(), include_in_report = FALSE, strand_sep = FALSE ) ``` ### Some additional notes Depending upon the use case, one may use any ordering for the input clusters for all functions. Looking at their annotations, if some clusters happen to to occur predominantly in 5' UTR or 3' UTR or exons or introns etc., such clusters can be sequestered from further downstream analyses if suitable for the end goal. ## Generating all plots at once You can let `r Biocpkg("seqArchRplus")` generate all plots using default settings like so. ```{r generate-all-plots, eval=FALSE} generate_all_plots( sname = "RAL28_10_to_12h", bed_info_fname = bed_info_fname, seqArchR_clusts = unord_clusts, raw_seqs = raw_seqs, tc_gr = tc_gr, use_q_bound = FALSE, order_by_iqw = TRUE, use_median_iqw = TRUE, iqw = TRUE, tpm = TRUE, cons = FALSE, pos_lab = -45:45, txdb_obj = TxDb.Dmelanogaster.UCSC.dm6.ensGene, org_name = "Dmelanogaster", qLow = 0.1, qUp = 0.9, tss_region = c(-500, 100), raw_seqs_mh = use_seqs, motifs = c("WW", "SS", "TATAA", "CG"), motif_heatmaps_flanks = c(50, 100, 200), dir_path = tempdir(), txt_size = 25 ) ``` # Conclusion Promoter sequences can be identified using CAGE or variations of CAGE experiment. Clusters of promoter sequences (identified by `r Biocpkg("seqArchR")` or otherwise) can be further analyzed with `r Biocpkg("seqArchRplus")`. This vignette demonstrated many of the downstream analyses steps for studying these promoter sequence clusters/architectures. # Session Info ```{r session_info, include=TRUE, echo=TRUE, results='markup'} sessionInfo() ``` # References