uncoverappLib: a R shiny package containing unCOVERApp an interactive graphical application for clinical assessment of sequence coverage at the base-pair level (original) (raw)
Contents
knitr::opts_chunk$set(
collapse = TRUE,
comment = "#>"
)
Introduction
This is a package containing unCOVERApp, a shiny graphical application for clinical assessment of sequence coverage. unCOVERApp allows:
- to display interactive plots showing sequence gene coverage down to base-pair resolution and functional/ clinical annotations of sequence positions within coverage gaps (
Coverage Analysis
page). - to calculate the maximum credible population allele frequency (AF) to be applied as AF filtering threshold tailored to the model of the disease-of-interest instead of a general AF cut-off (e.g. 1 % or 0.1 %) (
Calculate AF by allele frequency app
page). - to calculate the 95 % probability of the binomial distribution to observe at least N variant-supporting reads (N is the number of successes) based on a user-defined allele fraction that is expected for the variant (which is the probability of success). Especially useful to obtain the range of variant-supporting reads that is most likely to occur at a given a depth of coverage (DoC)(which is the number of trials) for somatic variants with low allele fraction (
Binomial distribution
page).
Installation and example
Install the latest version of uncoverappLib
using BiocManager
.
uncoverappLib
requires:
- R version >= 4.1.2
- java installed
- Annotation files for clinical assessment of low-coverage positions.
install.packages("BiocManager")
BiocManager::install("uncoverappLib")
library(uncoverappLib)
Alternatively, it can be installed from GitHub using:
#library(devtools)
#install_github("Manuelaio/uncoverappLib")
#library(uncoverappLib)
When users load uncoverappLib
for the first time, the first thing to do is a download of annotation files.getAnnotationFiles()
function allows to download the annotation files from Zenodo and parse it using uncoverappLib package. The function does not return an R object but store the annotation files in a cache (sorted.bed.gz
and sorted.bed.gz.tbi
) and show
the cache path. The local cache is managed by the BiocFileCache
Bioconductor package. It is sufficient run the function getAnnotationFiles(verbose= TRUE)
one time after installing uncoverappLib package as shown below. The preprocessing time can take few minutes, therefore during running vignette, users can provide vignette= TRUE
as a parameter to download an example annotation files, as below.
library(uncoverappLib)
#>
#>
#getAnnotationFiles(verbose= TRUE, vignette= TRUE)
The preprocessing time can take few minutes.
Input file
All unCOVERApp functionalities are based on the availability of a BED-style formatted input file containing tab-separated specifications of genomic coordinates (chromosome, start position, end position), the coverage value, and the reference:alternate allele counts for each position. In the first page Preprocessing, users can prepare the input file by specifying the genes to be examined and the BAM file(s) to be inspected. Users should be able to provide:
- a text file, with .txt extension, containing HGNC official gene name(s) one per row and to be uploaded to
Load input file
box. An example file is included in extdata of uncoverappLib packages
gene.list<- system.file("extdata", "mygene.txt", package = "uncoverappLib")
- a text file, with .list extension, containing absolute paths to BAM files (one per row) and to be uploaded to
Load bam file(s) list
box.
Type the following command to load our example:
bam_example <- system.file("extdata", "example_POLG.bam", package = "uncoverappLib")
print(bam_example)
#> [1] "/tmp/RtmpmGsUUf/Rinst288daf35e4ddf1/uncoverappLib/extdata/example_POLG.bam"
write.table(bam_example, file= "./bam.list", quote= FALSE, row.names = FALSE,
col.names = FALSE)
and launch run.uncoverapp(where="browser")
command. After running run.uncoverapp(where="browser")
the shiny app appears in your deafult browser. RStudio user can define where launching uncoverapp using where
option:
browser
option will openuncoverapp
in your default browserviewer
option will openuncoverapp
in RStudio viewerwindow
option will openuncoverapp
in RStudio RStudio
If option where
is not defined uncoverapp will launch with default option of R.
In the first page Preprocessing users can load mygene.txt
inLoad input file
and bam.list
in Load bam file(s) list
. In general, a target bed can also be used instead of genes name selecting Target Bed
option in Choose the type of your input file
. Users should also specify the reference genome in Genome
box and the chromosome notation of their BAM file(s) in Chromosome Notation
box. In the BAM file, the number option refers to 1, 2, …, X,.M chromosome notation, while the chr option refers to chr1, chr2, … chrX, chrM chromosome notation. Users can specify the minimum mapping quality (MAPQ)
value in box andminimum base quality (QUAL)
value in box. Default values for both mapping and base qualities is 1. Users can download Statistical_Summary
report to obtain a coverage metrics per genes (List of genes name
) or per amplicons (Target Bed
) according to uploaded input file. The report summarizes following information: mean, median, number of positions under 20x and percentage of position above 20x.
To run the example, choose chr chromosome notation,hg19 genome reference and leave minimum mapping and base qualities to the default settings, as shown in the following screenshot of the Preprocessing page:
Figure 1: Screenshot of Preprocessing page
unCOVERApp input file generation fails if incorrect gene names are specified. An unrecognized gene name(s) table is displayed if such a case occurs. Below is a snippet of a the unCOVERApp input file generated as a result of the preprocessing step performed for the example
chr15 89859516 89859516 68 A:68
chr15 89859517 89859517 70 T:70
chr15 89859518 89859518 73 A:2;G:71
chr15 89859519 89859519 73 A:73
chr15 89859520 89859520 74 C:74
chr15 89859521 89859521 75 C:1;T:74
The preprocessing time depends on the size of the BAM file(s) and on the number of genes to investigate. In general, if many (e.g. > 50) genes are to be analyzed, we would recommend to use buildInput
function
in R console before launching the app as shown in following example. This function also return a file with .txt estention containg statistical report of each genes/amplicon Alternatively, other tools do a similar job and can be used to generate the unCOVERApp input file ( for instance:bedtools,samtools,gatk). In this case, users can load the file directly onCoverage Analysis page in Select input file
box.
Once pre-processing is done, users can move to the Coverage Analysis page and push the load prepared input file
button.
Figure 2: Screenshot of Coverage Analysis page
To assess sequence coverage of the example, the following input parameters must be specified in the sidebar of the Coverage Analysis section
Reference Genome
: reference genome (hg19 or hg38); choose hg19Gene name
and pushApply
button: write the HGNC official gene name POLGCoverage threshold
: specify coverage threshold (e.g. 20x)Sample
: sample name to be analyzedTranscript number
: transcript number. Choose 1exon number
: to zoom in a specific exon. Choose 10
Other input sections, as Chromosome
, Transcript ID
, START genomic position
,END genomic position
and Region coordinate
, are dynamically filled.
Output
unCOVERApp generates the following outputs :
- unfiltered BED file in
bed file
and the corresponding filtered dataset inLow-coverage positions
- information about POLG gene in
UCSC gene
table
Figure 3: Screenshot of output of UCSC gene table
- information about POLG exons in
UCSC exons
table
Figure 4: Screenshot of output of Exon genomic coordinate positions from UCSC table
- sequence gene coverage plot in
Gene coverage
. The plot displays the chromosome ideogram, the genomic location and gene annotations from Ensembland the transcript(s) annotation from UCSC. Processing time is few minutes. A related table shows the number of uncovered positions in each exon given a user-defined transcript number (here transcript number is 1), and the user-defined threshold coverage (here the coverage threshold is 20x). Table and plot both show the many genomic positions that display low-DoC profile in POLG.
Figure 5: Screenshot of output of gene coverage
- plot of a specific exon, choose exon 10 in sidebar Exon number , push
Make exon
and view the plot inExon coverage
. Processing time is few minutes. A related table shows the number of low-DoC positions in ClinVar which have a high impact annotation. For this output to be generated, sorted.bed.gzand sorted.bed.gz.tbi are required to be downloaded withgetAnnotationFiles()
function. Table and plot both show that 21 low-DoC genomic positions have ClinVar annotation, suggesting several clinically relevant positions that are not adequately represented in this experiment. It is possible zooming at base pair level choosing a few interval (20-30 bp) inRegion coordinates
and moving onZoom to sequence
.
Figure 6: zoom of exon 10
- dbNSFP annotation of low coverage positions can be found in
Annotations on low-coverage positions
. Functional and clinical annotations of all potential non- synonymous single-nucleotide variants across the examined low DoC sites are made available. Potential changes that have a clinical annotation, a high impact or deleterious prediction are highlighted in yellow. In the example, a low Doc site (chr15:89868687) is predicted as pathogenic and could be potentially linked to disease.
Figure 7: Example of uncovered positions annotate with dbNSFP
By clicking on the download
button, users can save the table as spreadsheet format with certain cells colored according to pre-specified thresholds for AF, CADD, MAP-CAP, SIFT, Polyphen2, ClinVar, OMIM ID, HGVSp and HGVSc, …).
In Calculate maximum credible allele frequency page, users can set allele frequency cut-offs based on specific assumptions about the genetic architecture of the disease. If not specified, variants with allele frequency > 5 % will be instead filtered out. More details are availablehere. Moreover, users may click on the ”download” button and save the resulting table as spreadsheet format.
The Binomial distribution page returns the 95 % binomial probability distribution of the variant supporting reads on the input genomic position (Genomic position
). Users should define the expected allele fraction
(the expected fraction of variant reads, probability of success) and Variant reads
(the minimum number of variant reads required by the user to support variant calling, number of successes). The comment color change according to binomial proportion intervals. If the estimated intervals , with 95% confidence, is included or higher than user-definedVariant reads
the color of comment appears blue, otherwise if it is lower the color appears red.
Session information
#> R version 4.5.0 beta (2025-04-02 r88102)
#> Platform: x86_64-pc-linux-gnu
#> Running under: Ubuntu 24.04.2 LTS
#>
#> Matrix products: default
#> BLAS: /home/biocbuild/bbs-3.22-bioc/R/lib/libRblas.so
#> LAPACK: /usr/lib/x86_64-linux-gnu/lapack/liblapack.so.3.12.0 LAPACK version 3.12.0
#>
#> locale:
#> [1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C
#> [3] LC_TIME=en_GB LC_COLLATE=C
#> [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: America/New_York
#> tzcode source: system (glibc)
#>
#> attached base packages:
#> [1] stats graphics grDevices utils datasets methods base
#>
#> other attached packages:
#> [1] uncoverappLib_1.19.0 BiocStyle_2.37.0
#>
#> loaded via a namespace (and not attached):
#> [1] later_1.4.2
#> [2] BiocIO_1.19.0
#> [3] bitops_1.0-9
#> [4] filelock_1.0.3
#> [5] tibble_3.2.1
#> [6] graph_1.87.0
#> [7] XML_3.99-0.18
#> [8] rpart_4.1.24
#> [9] lifecycle_1.0.4
#> [10] httr2_1.1.2
#> [11] Homo.sapiens_1.3.1
#> [12] processx_3.8.6
#> [13] lattice_0.22-7
#> [14] ensembldb_2.33.0
#> [15] OrganismDbi_1.51.0
#> [16] backports_1.5.0
#> [17] magrittr_2.0.3
#> [18] openxlsx_4.2.8
#> [19] Hmisc_5.2-3
#> [20] sass_0.4.10
#> [21] rmarkdown_2.29
#> [22] jquerylib_0.1.4
#> [23] yaml_2.3.10
#> [24] rlist_0.4.6.2
#> [25] shinyBS_0.61.1
#> [26] httpuv_1.6.15
#> [27] zip_2.3.2
#> [28] Gviz_1.53.0
#> [29] DBI_1.2.3
#> [30] RColorBrewer_1.1-3
#> [31] abind_1.4-8
#> [32] GenomicRanges_1.61.0
#> [33] AnnotationFilter_1.33.0
#> [34] biovizBase_1.57.0
#> [35] BiocGenerics_0.55.0
#> [36] RCurl_1.98-1.17
#> [37] nnet_7.3-20
#> [38] VariantAnnotation_1.55.0
#> [39] rappdirs_0.3.3
#> [40] GenomeInfoDbData_1.2.14
#> [41] IRanges_2.43.0
#> [42] S4Vectors_0.47.0
#> [43] condformat_0.10.1
#> [44] codetools_0.2-20
#> [45] DelayedArray_0.35.0
#> [46] DT_0.33
#> [47] xml2_1.3.8
#> [48] tidyselect_1.2.1
#> [49] UCSC.utils_1.5.0
#> [50] shinyWidgets_0.9.0
#> [51] matrixStats_1.5.0
#> [52] stats4_4.5.0
#> [53] BiocFileCache_2.17.0
#> [54] base64enc_0.1-3
#> [55] GenomicAlignments_1.45.0
#> [56] jsonlite_2.0.0
#> [57] Formula_1.2-5
#> [58] tools_4.5.0
#> [59] progress_1.2.3
#> [60] TxDb.Hsapiens.UCSC.hg19.knownGene_3.2.2
#> [61] Rcpp_1.0.14
#> [62] glue_1.8.0
#> [63] gridExtra_2.3
#> [64] SparseArray_1.9.0
#> [65] xfun_0.52
#> [66] MatrixGenerics_1.21.0
#> [67] GenomeInfoDb_1.45.0
#> [68] dplyr_1.1.4
#> [69] BiocManager_1.30.25
#> [70] fastmap_1.2.0
#> [71] latticeExtra_0.6-30
#> [72] shinyjs_2.1.0
#> [73] digest_0.6.37
#> [74] R6_2.6.1
#> [75] mime_0.13
#> [76] colorspace_2.1-1
#> [77] GO.db_3.21.0
#> [78] jpeg_0.1-11
#> [79] dichromat_2.0-0.1
#> [80] markdown_2.0
#> [81] biomaRt_2.65.0
#> [82] RSQLite_2.3.9
#> [83] generics_0.1.3
#> [84] data.table_1.17.0
#> [85] rtracklayer_1.69.0
#> [86] prettyunits_1.2.0
#> [87] httr_1.4.7
#> [88] htmlwidgets_1.6.4
#> [89] S4Arrays_1.9.0
#> [90] pkgconfig_2.0.3
#> [91] gtable_0.3.6
#> [92] blob_1.2.4
#> [93] XVector_0.49.0
#> [94] htmltools_0.5.8.1
#> [95] bookdown_0.43
#> [96] RBGL_1.85.0
#> [97] ProtGenerics_1.41.0
#> [98] scales_1.3.0
#> [99] Biobase_2.69.0
#> [100] TxDb.Hsapiens.UCSC.hg38.knownGene_3.21.0
#> [101] png_0.1-8
#> [102] EnsDb.Hsapiens.v75_2.99.0
#> [103] knitr_1.50
#> [104] rstudioapi_0.17.1
#> [105] rjson_0.2.23
#> [106] checkmate_2.3.2
#> [107] curl_6.2.2
#> [108] org.Hs.eg.db_3.21.0
#> [109] cachem_1.1.0
#> [110] stringr_1.5.1
#> [111] parallel_4.5.0
#> [112] shinycssloaders_1.1.0
#> [113] foreign_0.8-90
#> [114] AnnotationDbi_1.71.0
#> [115] restfulr_0.0.15
#> [116] pillar_1.10.2
#> [117] grid_4.5.0
#> [118] vctrs_0.6.5
#> [119] promises_1.3.2
#> [120] dbplyr_2.5.0
#> [121] EnsDb.Hsapiens.v86_2.99.0
#> [122] xtable_1.8-4
#> [123] cluster_2.1.8.1
#> [124] htmlTable_2.4.3
#> [125] evaluate_1.0.3
#> [126] GenomicFeatures_1.61.0
#> [127] cli_3.6.4
#> [128] compiler_4.5.0
#> [129] Rsamtools_2.25.0
#> [130] rlang_1.1.6
#> [131] crayon_1.5.3
#> [132] interp_1.1-6
#> [133] ps_1.9.1
#> [134] stringi_1.8.7
#> [135] deldir_2.0-4
#> [136] BiocParallel_1.43.0
#> [137] txdbmaker_1.5.0
#> [138] munsell_0.5.1
#> [139] Biostrings_2.77.0
#> [140] lazyeval_0.2.2
#> [141] Matrix_1.7-3
#> [142] BSgenome_1.77.0
#> [143] hms_1.1.3
#> [144] bit64_4.6.0-1
#> [145] ggplot2_3.5.2
#> [146] KEGGREST_1.49.0
#> [147] shiny_1.10.0
#> [148] SummarizedExperiment_1.39.0
#> [149] memoise_2.0.1
#> [150] bslib_0.9.0
#> [151] bit_4.6.0