2. Using a GRanges object in shiny.gosling

Call required libraries.

require(shiny.gosling)
require(shiny)
require(GenomicRanges)

Getting a sample data for the GRanges object

We will be loading the peaks data from ChipSeq dataset with the GEO accession GSM1295076

GSM1295076_CBX6_BF_ChipSeq_mergedReps_peaks.bed.gz file will be used to create a sample GRanges object.

url <- "https://rb.gy/7y3fx"
temp_file <- file.path(tempdir(), "data.gz")
download.file(url, destfile = temp_file, method = "auto", mode = "wb")
df <- read.delim(
  temp_file,
  header = FALSE,
  comment.char = "#",
  sep = ""
)
gr <- GRanges(
  seqnames = df$V1,
  ranges = IRanges(df$V2, df$V3)
)
gr
#> GRanges object with 1331 ranges and 0 metadata columns:
#>          seqnames              ranges strand
#>             <Rle>           <IRanges>  <Rle>
#>      [1]     chr1       815092-817883      *
#>      [2]     chr1     1243287-1244338      *
#>      [3]     chr1     2979976-2981228      *
#>      [4]     chr1     3566181-3567876      *
#>      [5]     chr1     3816545-3818111      *
#>      ...      ...                 ...    ...
#>   [1327]     chrX 135244782-135245821      *
#>   [1328]     chrX 139171963-139173506      *
#>   [1329]     chrX 139583953-139586126      *
#>   [1330]     chrX 139592001-139593238      *
#>   [1331]     chrY   13845133-13845777      *
#>   -------
#>   seqinfo: 24 sequences from an unspecified genome; no seqlengths

Using the GRanges object for a plot using shiny.gosling

Method 1 - Using the track_data_gr function

You can use the track_data_gr to pass the GRanges object inside a track.

Note: Make sure to run the Shiny app using the shiny::runApp() rather than interactively running the shiny::shinyApp() object.


ui <- fluidPage(
  use_gosling(clear_files = FALSE),
  goslingOutput("gosling_plot")
)

track_1 <- add_single_track(
  width = 800,
  height = 180,
  data = track_data_gr(
    gr, chromosomeField = "seqnames",
    genomicFields = c("start", "end")
  ),
  mark = "bar",
  x = visual_channel_x(
    field = "start", type = "genomic", axis = "bottom"
  ),
  xe = visual_channel_x(field = "end", type = "genomic"),
  y = visual_channel_y(
    field = "width", type = "quantitative", axis = "right"
  ),
  size = list(value = 5)
)

composed_view <- compose_view(
  layout = "linear",
  tracks = track_1
)

arranged_view <- arrange_views(
  title = "Basic Marks: bar",
  subtitle = "Tutorial Examples",
  views = composed_view
)

server <- function(input, output, session) {
  output$gosling_plot <- renderGosling({
    gosling(
      component_id = "component_1",
      arranged_view
    )
  })
}

shiny::shinyApp(ui, server)

Method 2 - Using the track_data_csv function

You can save the GRanges object as a csv file inside the www directory which can be used in the shiny.gosling plot.

if (!dir.exists("data")) {
  dir.create("data")
}
utils::write.csv(gr, "data/ChipSeqPeaks.csv", row.names = FALSE)

track_1 <- add_single_track(
  width = 800,
  height = 180,
  data = track_data_csv(
    "data/ChipSeqPeaks.csv", chromosomeField = "seqnames",
    genomicFields = c("start", "end")
  ),
  mark = "bar",
  x = visual_channel_x(
    field = "start", type = "genomic", axis = "bottom"
  ),
  xe = visual_channel_x(field = "end", type = "genomic"),
  y = visual_channel_y(
    field = "width", type = "quantitative", axis = "right"
  ),
  size = list(value = 5)
)

composed_view <- compose_view(
  layout = "linear",
  tracks = track_1
)

arranged_view <- arrange_views(
  title = "Basic Marks: bar",
  subtitle = "Tutorial Examples",
  views = composed_view
)

shiny::shinyApp(ui = fluidPage(
  use_gosling(clear_files = FALSE),
  goslingOutput("gosling_plot")
), server = function(input, output, session) {
  output$gosling_plot <- renderGosling({
    gosling(
      component_id = "component_1",
      arranged_view
    )
  })
}, options = list(height = 1000))

Session Info


sessionInfo()
#> R Under development (unstable) (2024-10-21 r87258)
#> Platform: x86_64-pc-linux-gnu
#> Running under: Ubuntu 24.04.1 LTS
#> 
#> Matrix products: default
#> BLAS:   /home/biocbuild/bbs-3.21-bioc/R/lib/libRblas.so 
#> LAPACK: /usr/lib/x86_64-linux-gnu/lapack/liblapack.so.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] stats4    stats     graphics  grDevices utils     datasets  methods  
#> [8] base     
#> 
#> other attached packages:
#> [1] GenomicRanges_1.59.0 GenomeInfoDb_1.43.0  IRanges_2.41.0      
#> [4] S4Vectors_0.45.0     BiocGenerics_0.53.0  shiny_1.9.1         
#> [7] shiny.gosling_1.3.0 
#> 
#> loaded via a namespace (and not attached):
#>  [1] httr_1.4.7              cli_3.6.3               knitr_1.48             
#>  [4] rlang_1.1.4             xfun_0.48               stringi_1.8.4          
#>  [7] UCSC.utils_1.3.0        promises_1.3.0          jsonlite_1.8.9         
#> [10] xtable_1.8-4            glue_1.8.0              shiny.react_0.4.0      
#> [13] htmltools_0.5.8.1       httpuv_1.6.15           sass_0.4.9             
#> [16] rmarkdown_2.28          evaluate_1.0.1          jquerylib_0.1.4        
#> [19] fastmap_1.2.0           yaml_2.3.10             lifecycle_1.0.4        
#> [22] compiler_4.5.0          fs_1.6.4                Rcpp_1.0.13            
#> [25] XVector_0.47.0          later_1.3.2             digest_0.6.37          
#> [28] R6_2.5.1                GenomeInfoDbData_1.2.13 magrittr_2.0.3         
#> [31] bslib_0.8.0             tools_4.5.0             mime_0.12              
#> [34] zlibbioc_1.53.0         cachem_1.1.0