Description Usage Arguments Details Value Author(s) See Also Examples
View source: R/callPeaksMultivariate.R
Fit a HMM to multiple ChIP-seq samples to determine the combinatorial state of genomic regions. Input is a list of uniHMM
s generated by callPeaksUnivariate
.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 |
hmms |
A list of |
use.states |
A data.frame with combinatorial states which are used in the multivariate HMM, generated by function |
max.states |
Maximum number of combinatorial states to use in the multivariate HMM. The states are ordered by occurrence as determined from the combination of univariate state calls. |
per.chrom |
If |
chromosomes |
A vector specifying the chromosomes to use from the models in |
eps |
Convergence threshold for the Baum-Welch algorithm. |
keep.posteriors |
If set to |
num.threads |
Number of threads to use. Setting this to >1 may give increased performance. |
max.time |
The maximum running time in seconds for the Baum-Welch algorithm. If this time is reached, the Baum-Welch will terminate after the current iteration finishes. The default |
max.iter |
The maximum number of iterations for the Baum-Welch algorithm. The default |
keep.densities |
If set to |
verbosity |
Verbosity level for the fitting procedure. 0 - No output, 1 - Iterations are printed. |
temp.savedir |
A directory where to store intermediate results if |
Emission distributions from the univariate HMMs are used with a Gaussian copula to generate a multivariate emission distribution for each combinatorial state. This multivariate distribution is then kept fixed and the transition probabilities are fitted with a Baum-Welch. Please refer to our manuscript at http://dx.doi.org/10.1101/038612 for a detailed description of the method.
A multiHMM
object.
Aaron Taudt, Maria Colome Tatche
multiHMM
, callPeaksUnivariate
, callPeaksReplicates
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 | # Get example BAM files for 2 different marks in hypertensive rat
file.path <- system.file("extdata","euratrans", package='chromstaRData')
files <- list.files(file.path, full.names=TRUE, pattern='SHR.*bam$')[c(1:2,6)]
# Construct experiment structure
exp <- data.frame(file=files, mark=c("H3K27me3","H3K27me3","H3K4me3"),
condition=rep("SHR",3), replicate=c(1:2,1), pairedEndReads=FALSE,
controlFiles=NA)
states <- stateBrewer(exp, mode='combinatorial')
# Bin the data
data(rn4_chrominfo)
binned.data <- list()
for (file in files) {
binned.data[[basename(file)]] <- binReads(file, binsizes=1000, stepsizes=500,
experiment.table=exp,
assembly=rn4_chrominfo, chromosomes='chr12')
}
# Obtain the univariate fits
models <- list()
for (i1 in 1:length(binned.data)) {
models[[i1]] <- callPeaksUnivariate(binned.data[[i1]], max.time=60, eps=1)
}
# Call multivariate peaks
multimodel <- callPeaksMultivariate(models, use.states=states, eps=1, max.time=60)
# Check some plots
heatmapTransitionProbs(multimodel)
heatmapCountCorrelation(multimodel)
|
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.