Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

add genomecov #411

Merged
merged 4 commits into from
Oct 11, 2023
Merged

add genomecov #411

merged 4 commits into from
Oct 11, 2023

Conversation

kriemo
Copy link
Member

@kriemo kriemo commented Oct 11, 2023

library(valr)
library(dplyr)
#> 
#> Attaching package: 'dplyr'
#> The following objects are masked from 'package:stats':
#> 
#>     filter, lag
#> The following objects are masked from 'package:base':
#> 
#>     intersect, setdiff, setequal, union
library(microbenchmark)

genome <- read_genome("https://hgdownload.soe.ucsc.edu/goldenPath/sacCer3/bigZips/sacCer3.chrom.sizes")

ivls <- bed_random(
  genome,
  length = 100,
  n = 1e5, 
  seed = 42
) 

bp <- bed_partition(ivls, .depth = n())
bp
#> # A tibble: 154,751 × 4
#>    chrom start   end .depth
#>    <chr> <int> <int>  <int>
#>  1 chrI     59   159      1
#>  2 chrI    165   207      1
#>  3 chrI    207   265      2
#>  4 chrI    265   307      1
#>  5 chrI    456   527      1
#>  6 chrI    527   534      2
#>  7 chrI    534   556      3
#>  8 chrI    556   617      2
#>  9 chrI    617   627      3
#> 10 chrI    627   634      2
#> # ℹ 154,741 more rows

bgc <- bed_genomecov(ivls, genome)
bgc
#> # A tibble: 154,751 × 4
#>    chrom start   end .depth
#>    <chr> <int> <int>  <int>
#>  1 chrI     59   159      1
#>  2 chrI    165   207      1
#>  3 chrI    207   265      2
#>  4 chrI    265   307      1
#>  5 chrI    456   527      1
#>  6 chrI    527   534      2
#>  7 chrI    534   556      3
#>  8 chrI    556   617      2
#>  9 chrI    617   627      3
#> 10 chrI    627   634      2
#> # ℹ 154,741 more rows

identical(bp, bgc)
#> [1] TRUE

# 10 million intervals
ivls <- bed_random(
  genome,
  length = 100,
  n = 1e7, 
  seed = 42
) 
microbenchmark(bed_genomecov(ivls, genome), times = 3, unit = 's')
#> Warning in microbenchmark(bed_genomecov(ivls, genome), times = 3, unit = "s"):
#> less accurate nanosecond times to avoid potential integer overflows
#> Unit: seconds
#>                         expr      min       lq     mean   median       uq
#>  bed_genomecov(ivls, genome) 1.653511 1.677021 1.860034 1.700531 1.963296
#>      max neval
#>  2.22606     3

Created on 2023-10-11 with reprex v2.0.2

@jayhesselberth
Copy link
Member

Can you run styler::style_pkg()?

@jayhesselberth
Copy link
Member

This looks great, thanks for putting it together.

@kriemo kriemo merged commit a4093ef into main Oct 11, 2023
7 checks passed
@kriemo kriemo deleted the genomecov branch October 11, 2023 17:53
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

2 participants