trichelab / biscuiteer Goto Github PK
View Code? Open in Web Editor NEWtools for analysing Biscuit (https://huishenlab.github.io/biscuit/) output
tools for analysing Biscuit (https://huishenlab.github.io/biscuit/) output
The flogit and fexpit functions are not invertable as written similar to compartmap usage of these functions. Just noting here as well.
Old flogit:
#' Calculate squeezed logit function
#'
#' Helper function for calculating squeezed logit
#'
#' @param x A vector of values between 0 and 1 inclusive
#' @param sqz The amount by which to 'squeeze' (DEFAULT: 0.000001)
#'
#' @return A vector of values between -Inf and +Inf
#'
#' @importFrom gtools logit
#'
flogit <- function(x,
sqz = 0.000001) {
x[ which(x < sqz) ] <- sqz
x[ which(x > (1 - sqz)) ] <- (1 - sqz)
gtools::logit(x)
}
New flogit as written by @ttriche in compartmap
#' Helper function: squeezed logit
#'
#' @param p a vector of values between 0 and 1 inclusive
#' @param sqz the amount by which to 'squeeze', default is .000001
#'
#' @return a vector of values between -Inf and +Inf
#'
#' @examples
#'
#' set.seed(1234)
#' p <- runif(n=1000)
#' summary(p)
#'
#' sqz <- 1 / (10**6)
#' x <- flogit(p, sqz=sqz)
#' summary(x)
#'
#' all( abs(p - fexpit(x, sqz=sqz)) < sqz )
#' all( abs(p - fexpit(flogit(p, sqz=sqz), sqz=sqz)) < sqz )
#'
#' @export
flogit <- function(p, sqz=0.000001) {
midpt <- 0.5
deflate <- 1 - (sqz * midpt)
if (any(p > 1 | p < 0)) stop("Values of p outside (0,1) detected.")
squoze <- ((p - midpt) * deflate) + midpt
return( log( squoze / (1 - squoze)) )
}
Will file PR to change.
When combining single sample bsseq objects using either biscuiteer::unionize() or bsseq::combineList() with a list of individual bsseq objects, the trans slot becomes very large (e.g. 120 gb). Further, this is propagated to the SummarizedExperiment object when saving as an HDF5-backed SE.
Looking at hansenlab/bsseq#53 (comment), it appears that this is just a function to do a smooth. Within the bsseq::combineList(), the trans argument is given as follows at construction of the new, combined bsseq object:
ans_trans <- function() NULL
The worry here is that for some reason, we are in fact copying data around, leading loss of all gains made by making it an HDF5-backed object.
One workaround is to reset the trans slot to be the same as whatever you're combining it with. There is an explicit check in bsseq::combineList() but not in biscuiteer::unionize() for this. Alternatively, we can set both to the above and combine.
Importing single samples that are proper bed files without a VCF:
#grab file names
my_samp_names <- list.files("path/to/bed/files") %>% gsub("*.bed", "", .)
#import using lapply
my_bsseq_list <- lapply(my_samp_names, function(x) read.biscuit(BEDFile = "path/to/my/bed", sampleNames = x, genome = "hg19")
returns objects M and Cov matrices instead of DelayedMatrix. Ultimately, this also returns straight beta values from biscuit-derived bed files (vcf2bed) instead of methylated read counts.
Will check to see which version of biscuiteer this was introduced.
Using chr11p15 (for reducing computation time), the hypoMeth and hyperMeth arrays have different sizes. This becomes an issue when trying to calculate the ratio between them, as R can't conform the arrays to one another.
Should the ratio be included? Or, can it be removed and left up to the user to calculate the ratio on their own if they decide they want to do that?
Hi. I'm running reverse dependency checks on parallelly, I noticed that your package lacks Suggests: BiocStyle
, cf. https://github.com/HenrikBengtsson/parallelly/blob/develop/revdep/problems.md#biscuiteer. Without, it you get:
*** checking re-building of vignette outputs ... WARNING
Error(s) in re-building vignettes:
...
--- re-building ‘biscuiteer.Rmd’ using rmarkdown
Error: processing vignette 'biscuiteer.Rmd' failed with diagnostics:
there is no package called ‘BiocStyle’
--- failed re-building ‘biscuiteer.Rmd’
SUMMARY: processing the following file failed:
‘biscuiteer.Rmd’
Error: Vignette re-building failed.
Execution halted
Working on my first analysis w/ BISCUIT/biscuiteer, but I’ve encountered some issues handling the data. I have 20 gzip/tabix’d VCFs (15-20Gb each) with accompanying bed.gz files. Biscuiteer seems to be working just fine with small/toy datasets. However, I’ve been having issues merging all these samples into a single bsseq
object. I think part of the issue is simply due to the large sample number and the amount of data for each sample. I have attempted to solve this issue with two approaches that have failed thus far:
biscuiteer::readBiscuit()
for each sample individually and then use biscuiteer::unionize()
to get a single object.biscuiteer::readBiscuit()
Do you have any advice for a better/ideal approach in this situation?
thanks in advance!
When calling atRegions(bsseq, regions, mappings, nm, ...), you can pass it either dropNA or impute, which will then be passed on to summarizeBsSeqOver. When trying to pass it dropNA, atRegions will fail because of an incorrect number of dimensions (see error message below). The issue is arising because summarizeBsSeqOver will return a numeric, rather than a matrix, when only one sample is given. Should summarizeBsSeqOver always return a matrix, and then let atRegions reduce this to a numeric if there is only one sample?
Code:
tcga_bed <- system.file("extdata", "TCGA_BLCA_A13J_chr11p15_merged.bed.gz",
package = "biscuiteer")
tcga_vcf <- system.file("extdata", "TCGA_BLCA_A13J_header_only.vcf.gz",
package = "biscuiteer")
bisc <- read.biscuit(BEDfile = tcga_bed, VCFfile = tcga_vcf,
merged = TRUE, genome = "hg38", verbose = TRUE)
reg <- Ranges(seqnames = rep("chr11",5),
strand = rep("*",5),
ranges = IRanges(start = c(0,2.8e6,1.17e7,1.38e7,1.69e7),
end= c(2.8e6,1.17e7,1.38e7,1.69e7,2.2e7))
)
regions <- atRegions(bsseq = bisc, regions = reg, dropNA = TRUE)
Error message:
Error in regional[as.character(regions), ] :
incorrect number of dimensions
dim(myBS)
25555777 7
saveRDS(myBS, file = "biscuiteer.RDS") # 7.6 Gb on disk
test <- myBS[1,]
dim(test)
1 7
saveRDS(test, file = "test.RDS") # 7.4 Gb on disk
Hi,
Trying to install biscuiteer with
BiocManager::install("trichelab/biscuiteer",update = F)
(update is False becuase I'm running on cluster, and I don't have root permission to update)
And I'm getting the error:
"this R is version 3.5.1, package 'biscuiteer' requires R >= 3.6.0"
Is there something I can do, giving that I don't have root permission to update R?
A declarative, efficient, and flexible JavaScript library for building user interfaces.
🖖 Vue.js is a progressive, incrementally-adoptable JavaScript framework for building UI on the web.
TypeScript is a superset of JavaScript that compiles to clean JavaScript output.
An Open Source Machine Learning Framework for Everyone
The Web framework for perfectionists with deadlines.
A PHP framework for web artisans
Bring data to life with SVG, Canvas and HTML. 📊📈🎉
JavaScript (JS) is a lightweight interpreted programming language with first-class functions.
Some thing interesting about web. New door for the world.
A server is a program made to process requests and deliver data to clients.
Machine learning is a way of modeling and interpreting data that allows a piece of software to respond intelligently.
Some thing interesting about visualization, use data art
Some thing interesting about game, make everyone happy.
We are working to build community through open source technology. NB: members must have two-factor auth.
Open source projects and samples from Microsoft.
Google ❤️ Open Source for everyone.
Alibaba Open Source for everyone
Data-Driven Documents codes.
China tencent open source team.