## ----style, echo = FALSE, results = 'asis', message = FALSE, warning=FALSE---- BiocStyle::markdown() library(knitr) ## ----loadingPackage, warning=FALSE, message=FALSE-------------------------- library(similaRpeak) ## ----bam_extract_coverage, collapse=TRUE----------------------------------- suppressMessages(library(GenomicAlignments)) suppressMessages(library(rtracklayer)) suppressMessages(library(Rsamtools)) bamFile01 <- system.file("extdata/align1.bam", package = "similaRpeak") region <- GRanges(Rle(c("chr1")), IRanges(start = 172081530, end = 172083529), strand= Rle(strand(c("*")))) param <- ScanBamParam(which = region) alignments01 <- readGAlignments(bamFile01, param = param) coveragesRegion01 <- coverage(alignments01)[region] coveragesRegion01 ## ----iranges_to_vector, collapse=TRUE-------------------------------------- coveragesRegion01 <- as.numeric(coveragesRegion01$chr1) length(coveragesRegion01) summary(coveragesRegion01) ## ----bam_count, collapse=TRUE---------------------------------------------- param <- ScanBamParam(flag = scanBamFlag(isUnmappedQuery=FALSE)) count01 <- countBam(bamFile01, param = param)$records coveragesRPMRegion01 <- coverage(alignments01, weight=1000000/count01)[region] coveragesRPMRegion01 <- as.numeric(coveragesRPMRegion01$chr1) length(coveragesRPMRegion01) summary(coveragesRPMRegion01) ## ----demo_profiles_loading, collapse=T------------------------------------- data(demoProfiles) str(demoProfiles) ## ----ratio_area_graph, echo=FALSE, fig.width=11, fig.height=8-------------- par(mar=c(6,4,2,2)) par(mfrow=c(2,2)) plot(demoProfiles$chr2.70360770.70361098$H3K27ac, type="l", col="blue", xlab="", ylab="", ylim=c(0, 25), main="chr2:70360770-70361098") par(new=TRUE) plot(demoProfiles$chr2.70360770.70361098$H3K4me1, type="l", col="darkgreen", xlab="Position", ylab="Coverage in reads per million (RPM)", ylim=c(0, 25)) legend("topright", c("H3K27ac","H3K4me1"), cex=1.2, col=c("blue","darkgreen"), lty=1) text(250, 17, "RATIO_AREA=1.03", cex=1.5, col="red") plot(demoProfiles$chr8.43092918.43093442$H3K27ac, type="l", col="blue", xlab="", ylab="", ylim=c(0, 1500), main="chr8:43092918-43093442") par(new=TRUE) plot(demoProfiles$chr8.43092918.43093442$H3K4me1, type="l", col="darkgreen", xlab="Position", ylab="Coverage in reads per million (RPM)", ylim=c(0, 1500)) legend("topright", c("H3K27ac","H3K4me1"), cex=1.2, col=c("blue","darkgreen"), lty=1) text(375, 1000, "RATIO_AREA=0.06", cex=1.5, col="red") plot(demoProfiles$chr3.73159773.73160145$H3K27ac, type="l", col="blue", xlab="", ylab="", ylim=c(0, 125), main="chr3:73159773-73160145") par(new=TRUE) plot(demoProfiles$chr3.73159773.73160145$H3K4me1, type="l", col="darkgreen", xlab="Position", ylab="Coverage in reads per million (RPM)", ylim=c(0, 125)) legend("topright", c("H3K27ac","H3K4me1"), cex=1.2, col=c("blue","darkgreen"), lty=1) text(285, 80, "RATIO_AREA=2.23", cex=1.5, col="red") plot(demoProfiles$chr19.27739373.27739767$H3K27ac, type="l", col="blue", xlab="", ylab="", ylim=c(0, 350), main="chr19:27739373-27739767") par(new=TRUE) plot(demoProfiles$chr19.27739373.27739767$H3K4me1, type="l", col="darkgreen", xlab="Position", ylab="Coverage in reads per million (RPM)", ylim=c(0, 350)) legend("topright", c("H3K27ac","H3K4me1"), cex=1.2, col=c("blue","darkgreen"), lty=1) text(275, 225, "RATIO_AREA=0.12", cex=1.5, col="red") ## ----ratio_area_calculation, collapse=T------------------------------------ metrics <- similarity(demoProfiles$chr2.70360770.70361098$H3K27ac, demoProfiles$chr2.70360770.70361098$H3K4me1) metrics$metric$RATIO_AREA metrics <- similarity(demoProfiles$chr8.43092918.43093442$H3K27ac, demoProfiles$chr8.43092918.43093442$H3K4me1) metrics$metric$RATIO_AREA metrics <- similarity(demoProfiles$chr3.73159773.73160145$H3K27ac, demoProfiles$chr3.73159773.73160145$H3K4me1) metrics$metric$RATIO_AREA metrics <- similarity(demoProfiles$chr19.27739373.27739767$H3K27ac, demoProfiles$chr19.27739373.27739767$H3K4me1) metrics$metric$RATIO_AREA ## ----diff_pos_max, echo=FALSE, fig.width=11, fig.height=8------------------ par(mar=c(6,4,2,2)) par(mfrow=c(2,2)) plot(demoProfiles$chr2.70360770.70361098$H3K27ac, type="l", col="blue", xlab="", ylab="", ylim=c(0, 25), main="chr2:70360770-70361098") par(new=TRUE) plot(demoProfiles$chr2.70360770.70361098$H3K4me1, type="l", col="darkgreen", xlab="Position", ylab="Coverage in reads per million (RPM)", ylim=c(0, 25)) legend("topright", c("H3K27ac","H3K4me1"), cex=1.2, col=c("blue","darkgreen"), lty=1) text(250, 17, "DIFF_POS_MAX=-20", cex=1.5, col="red") plot(demoProfiles$chr8.43092918.43093442$H3K27ac, type="l", col="blue", xlab="", ylab="", ylim=c(0, 1500), main="chr8:43092918-43093442") par(new=TRUE) plot(demoProfiles$chr8.43092918.43093442$H3K4me1, type="l", col="darkgreen", xlab="Position", ylab="Coverage in reads per million (RPM)", ylim=c(0, 1500)) legend("topright", c("H3K27ac","H3K4me1"), cex=1.2, col=c("blue","darkgreen"), lty=1) text(375, 1000, "DIFF_POS_MAX=-0.5", cex=1.5, col="red") plot(demoProfiles$chr3.73159773.73160145$H3K27ac, type="l", col="blue", xlab="", ylab="", ylim=c(0, 125), main="chr3:73159773-73160145") par(new=TRUE) plot(demoProfiles$chr3.73159773.73160145$H3K4me1, type="l", col="darkgreen", xlab="Position", ylab="Coverage in reads per million (RPM)", ylim=c(0, 125)) legend("topright", c("H3K27ac","H3K4me1"), cex=1.2, col=c("blue","darkgreen"), lty=1) text(285, 80, "DIFF_POS_MAX=2.5", cex=1.5, col="red") plot(demoProfiles$chr19.27739373.27739767$H3K27ac, type="l", col="blue", xlab="", ylab="", ylim=c(0, 350), main="chr19:27739373-27739767") par(new=TRUE) plot(demoProfiles$chr19.27739373.27739767$H3K4me1, type="l", col="darkgreen", xlab="Position", ylab="Coverage in reads per million (RPM)", ylim=c(0, 350)) legend("topright", c("H3K27ac","H3K4me1"), cex=1.2, col=c("blue","darkgreen"), lty=1) text(275, 225, "DIFF_POS_MAX=2.5", cex=1.5, col="red") ## ----diff_pos_max_calculation, collapse=T---------------------------------- metrics <- similarity(demoProfiles$chr2.70360770.70361098$H3K27ac, demoProfiles$chr2.70360770.70361098$H3K4me1) metrics$metric$DIFF_POS_MAX metrics <- similarity(demoProfiles$chr8.43092918.43093442$H3K27ac, demoProfiles$chr8.43092918.43093442$H3K4me1) metrics$metric$DIFF_POS_MAX metrics <- similarity(demoProfiles$chr3.73159773.73160145$H3K27ac, demoProfiles$chr3.73159773.73160145$H3K4me1) metrics$metric$DIFF_POS_MAX metrics <- similarity(demoProfiles$chr19.27739373.27739767$H3K27ac, demoProfiles$chr19.27739373.27739767$H3K4me1) metrics$metric$DIFF_POS_MAX ## ----ratio_max_max, echo=FALSE, fig.width=11, fig.height=8----------------- par(mar=c(6,4,2,2)) par(mfrow=c(2,2)) plot(demoProfiles$chr2.70360770.70361098$H3K27ac, type="l", col="blue", xlab="", ylab="", ylim=c(0, 25), main="chr2:70360770-70361098") par(new=TRUE) plot(demoProfiles$chr2.70360770.70361098$H3K4me1, type="l", col="darkgreen", xlab="Position", ylab="Coverage in reads per million (RPM)", ylim=c(0, 25)) legend("topright", c("H3K27ac","H3K4me1"), cex=1.2, col=c("blue","darkgreen"), lty=1) text(250, 17, "RATIO_MAX_MAX=0.95", cex=1.5, col="red") plot(demoProfiles$chr8.43092918.43093442$H3K27ac, type="l", col="blue", xlab="", ylab="", ylim=c(0, 1500), main="chr8:43092918-43093442") par(new=TRUE) plot(demoProfiles$chr8.43092918.43093442$H3K4me1, type="l", col="darkgreen", xlab="Position", ylab="Coverage in reads per million (RPM)", ylim=c(0, 1500)) legend("topright", c("H3K27ac","H3K4me1"), cex=1.2, col=c("blue","darkgreen"), lty=1) text(380, 1000, "RATIO_MAX_MAX=0.06", cex=1.5, col="red") plot(demoProfiles$chr3.73159773.73160145$H3K27ac, type="l", col="blue", xlab="", ylab="", ylim=c(0, 125), main="chr3:73159773-73160145") par(new=TRUE) plot(demoProfiles$chr3.73159773.73160145$H3K4me1, type="l", col="darkgreen", xlab="Position", ylab="Coverage in reads per million (RPM)", ylim=c(0, 125)) legend("topright", c("H3K27ac","H3K4me1"), cex=1.2, col=c("blue","darkgreen"), lty=1) text(290, 80, "RATIO_MAX_MAX=2.5", cex=1.5, col="red") plot(demoProfiles$chr19.27739373.27739767$H3K27ac, type="l", col="blue", xlab="", ylab="", ylim=c(0, 350), main="chr19:27739373-27739767") par(new=TRUE) plot(demoProfiles$chr19.27739373.27739767$H3K4me1, type="l", col="darkgreen", xlab="Position", ylab="Coverage in reads per million (RPM)", ylim=c(0, 350)) legend("topright", c("H3K27ac","H3K4me1"), cex=1.2, col=c("blue","darkgreen"), lty=1) text(280, 225, "RATIO_MAX_MAX=0.12", cex=1.5, col="red") ## ----ratio_max_max_calculation, collapse=T--------------------------------- metrics <- similarity(demoProfiles$chr2.70360770.70361098$H3K27ac, demoProfiles$chr2.70360770.70361098$H3K4me1) metrics$metric$RATIO_MAX_MAX metrics <- similarity(demoProfiles$chr8.43092918.43093442$H3K27ac, demoProfiles$chr8.43092918.43093442$H3K4me1) metrics$metric$RATIO_MAX_MAX metrics <- similarity(demoProfiles$chr3.73159773.73160145$H3K27ac, demoProfiles$chr3.73159773.73160145$H3K4me1) metrics$metric$RATIO_MAX_MAX metrics <- similarity(demoProfiles$chr19.27739373.27739767$H3K27ac, demoProfiles$chr19.27739373.27739767$H3K4me1) metrics$metric$RATIO_MAX_MAX ## ----ratio_intersect, echo=FALSE, fig.width=11, fig.height=8--------------- par(mar=c(6,4,2,2)) par(mfrow=c(2,2)) plot(demoProfiles$chr2.70360770.70361098$H3K27ac, type="l", col="blue", xlab="", ylab="", ylim=c(0, 25), main="chr2:70360770-70361098") par(new=TRUE) plot(demoProfiles$chr2.70360770.70361098$H3K4me1, type="l", col="darkgreen", xlab="Position", ylab="Coverage in reads per million (RPM)", ylim=c(0, 25)) legend("topright", c("H3K27ac","H3K4me1"), cex=1.2, col=c("blue","darkgreen"), lty=1) text(242, 17, "RATIO_INTERSECT=0.63", cex=1.5, col="red") plot(demoProfiles$chr8.43092918.43093442$H3K27ac, type="l", col="blue", xlab="", ylab="", ylim=c(0, 1500), main="chr8:43092918-43093442") par(new=TRUE) plot(demoProfiles$chr8.43092918.43093442$H3K4me1, type="l", col="darkgreen", xlab="Position", ylab="Coverage in reads per million (RPM)", ylim=c(0, 1500)) legend("topright", c("H3K27ac","H3K4me1"), cex=1.2, col=c("blue","darkgreen"), lty=1) text(375, 1000, "RATIO_INTERSECT=0.06", cex=1.5, col="red") plot(demoProfiles$chr3.73159773.73160145$H3K27ac, type="l", col="blue", xlab="", ylab="", ylim=c(0, 125), main="chr3:73159773-73160145") par(new=TRUE) plot(demoProfiles$chr3.73159773.73160145$H3K4me1, type="l", col="darkgreen", xlab="Position", ylab="Coverage in reads per million (RPM)", ylim=c(0, 125)) legend("topright", c("H3K27ac","H3K4me1"), cex=1.2, col=c("blue","darkgreen"), lty=1) text(288, 87, "RATIO_INTERSECT=", cex=1.5, col="red") text(288, 75, "0.43", cex=1.5, col="red") plot(demoProfiles$chr19.27739373.27739767$H3K27ac, type="l", col="blue", xlab="", ylab="", ylim=c(0, 350), main="chr19:27739373-27739767") par(new=TRUE) plot(demoProfiles$chr19.27739373.27739767$H3K4me1, type="l", col="darkgreen", xlab="Position", ylab="Coverage in reads per million (RPM)", ylim=c(0, 350)) legend("topright", c("H3K27ac","H3K4me1"), cex=1.2, col=c("blue","darkgreen"), lty=1) text(275, 225, "RATIO_INTERSECT=0.12", cex=1.5, col="red") ## ----ratio_intersect_calculation, collapse=T------------------------------- metrics <- similarity(demoProfiles$chr2.70360770.70361098$H3K27ac, demoProfiles$chr2.70360770.70361098$H3K4me1) metrics$metric$RATIO_INTERSECT metrics <- similarity(demoProfiles$chr8.43092918.43093442$H3K27ac, demoProfiles$chr8.43092918.43093442$H3K4me1) metrics$metric$RATIO_INTERSECT metrics <- similarity(demoProfiles$chr3.73159773.73160145$H3K27ac, demoProfiles$chr3.73159773.73160145$H3K4me1) metrics$metric$RATIO_INTERSECT metrics <- similarity(demoProfiles$chr19.27739373.27739767$H3K27ac, demoProfiles$chr19.27739373.27739767$H3K4me1) metrics$metric$RATIO_INTERSECT ## ----ratio_normalized_intersect, echo=FALSE, fig.width=11, fig.height=8---- par(mar = c(6,4,2,2)) par(mfrow = c(2,2)) plot(demoProfiles$chr2.70360770.70361098$H3K27ac* length(demoProfiles$chr2.70360770.70361098$H3K27ac)/ sum(demoProfiles$chr2.70360770.70361098$H3K27ac, na.rm=TRUE), type="l", col="blue", xlab="", ylab="", ylim=c(0, 3), main="chr2:70360770-70361098") par(new=TRUE) plot(demoProfiles$chr2.70360770.70361098$H3K4me1* length(demoProfiles$chr2.70360770.70361098$H3K4me1)/ sum(demoProfiles$chr2.70360770.70361098$H3K4me1, na.rm=TRUE), type="l", col="darkgreen", xlab="Position", ylab="Normalized Coverage (Coverage/Mean Coverage)", ylim=c(0, 3)) legend("topright", c("H3K27ac","H3K4me1"), cex=1.2, col=c("blue","darkgreen"), lty=1) text(245, 2, "RATIO_NORMALIZED_", cex=1.5, col="red") text(245, 1.75, "INTERSECT=0.63", cex=1.5, col="red") plot(demoProfiles$chr8.43092918.43093442$H3K27ac* length(demoProfiles$chr8.43092918.43093442$H3K27ac)/ sum(demoProfiles$chr8.43092918.43093442$H3K27ac, na.rm=TRUE), type="l", col="blue", xlab="", ylab="", ylim=c(0, 12), main="chr8:43092918-43093442") par(new=TRUE) plot(demoProfiles$chr8.43092918.43093442$H3K4me1* length(demoProfiles$chr8.43092918.43093442$H3K4me1)/ sum(demoProfiles$chr8.43092918.43093442$H3K4me1, na.rm=TRUE), type="l", col="darkgreen", xlab="Position", ylab="Normalized Coverage (Coverage/Mean Coverage)", ylim=c(0, 12)) legend("topright", c("H3K27ac","H3K4me1"), cex=1.2, col=c("blue","darkgreen"), lty=1) text(370, 8, "RATIO_NORMALIZED_", cex=1.5, col="red") text(370, 7, "INTERSECT=0.89", cex=1.5, col="red") plot(demoProfiles$chr3.73159773.73160145$H3K27ac* length(demoProfiles$chr3.73159773.73160145$H3K27ac)/ sum(demoProfiles$chr3.73159773.73160145$H3K27ac, na.rm=TRUE), type="l", col="blue", xlab="", ylab="", ylim=c(0, 8), main="chr3:73159773-73160145") par(new=TRUE) plot(demoProfiles$chr3.73159773.73160145$H3K4me1* length(demoProfiles$chr3.73159773.73160145$H3K4me)/ sum(demoProfiles$chr3.73159773.73160145$H3K4me, na.rm=TRUE), type="l", col="darkgreen", xlab="Position", ylab="Normalized Coverage (Coverage/Mean Coverage)", ylim=c(0, 8)) legend("topright", c("H3K27ac","H3K4me1"), cex=1.2, col=c("blue","darkgreen"), lty=1) text(288, 5.5, "RATIO_NORMALIZED_", cex=1.5, col="red") text(288, 4, "INTERSECT=0.78", cex=1.5, col="red") plot(demoProfiles$chr19.27739373.27739767$H3K27ac* length(demoProfiles$chr19.27739373.27739767$H3K27ac)/ sum(demoProfiles$chr19.27739373.27739767$H3K27ac, na.rm=TRUE), type="l", col="blue", xlab="", ylab="", ylim=c(0, 12), main="chr19:27739373-27739767") par(new=TRUE) plot(demoProfiles$chr19.27739373.27739767$H3K4me1* length(demoProfiles$chr19.27739373.27739767$H3K4me1)/ sum(demoProfiles$chr19.27739373.27739767$H3K4me1, na.rm=TRUE), type="l", col="darkgreen", xlab="Position", ylab="Normalized Coverage (Coverage/Mean Coverage)", ylim=c(0, 12)) legend("topright", c("H3K27ac","H3K4me1"), cex=1.2, col=c("blue","darkgreen"), lty=1) text(275, 7, "RATIO_NORMALIZED_", cex=1.5, col="red") text(275, 6, "INTERSECT=0.84", cex=1.5, col="red") ## ----ratio_normalized_intersect_calculation, collapse=TRUE----------------- metrics <- similarity(demoProfiles$chr2.70360770.70361098$H3K27ac, demoProfiles$chr2.70360770.70361098$H3K4me1) metrics$metric$RATIO_NORMALIZED_INTERSECT metrics <- similarity(demoProfiles$chr8.43092918.43093442$H3K27ac, demoProfiles$chr8.43092918.43093442$H3K4me1) metrics$metric$RATIO_NORMALIZED_INTERSECT metrics <- similarity(demoProfiles$chr3.73159773.73160145$H3K27ac, demoProfiles$chr3.73159773.73160145$H3K4me1) metrics$metric$RATIO_NORMALIZED_INTERSECT metrics <- similarity(demoProfiles$chr19.27739373.27739767$H3K27ac, demoProfiles$chr19.27739373.27739767$H3K4me1) metrics$metric$RATIO_NORMALIZED_INTERSECT ## ----spearman_corr_graphs, echo=FALSE, fig.width=11, fig.height=8---------- par(mar=c(6,4,2,2)) par(mfrow=c(2,2)) plot(demoProfiles$chr2.70360770.70361098$H3K27ac, type="l", col="blue", xlab="", ylab="", ylim=c(0, 25), main="chr2:70360770-70361098") par(new=TRUE) plot(demoProfiles$chr2.70360770.70361098$H3K4me1, type="l", col="darkgreen", xlab="Position", ylab="Coverage", ylim=c(0, 25)) legend("topright", c("H3K27ac","H3K4me1"), cex=1.2, col=c("blue","darkgreen"), lty=1) text(225, 17, "SPEARMAN_CORRELATION", cex=1.5, col="red") text(225, 14, "=0.06", cex=1.5, col="red") plot(demoProfiles$chr8.43092918.43093442$H3K27ac, type="l", col="blue", xlab="", ylab="", ylim=c(0, 1500), main="chr8:43092918-43093442") par(new=TRUE) plot(demoProfiles$chr8.43092918.43093442$H3K4me1, type="l", col="darkgreen", xlab="Position", ylab="Coverage", ylim=c(0, 1500)) legend("topright", c("H3K27ac","H3K4me1"), cex=1.2, col=c("blue","darkgreen"), lty=1) text(325, 1000, "SPEARMAN_CORRELATION", cex=1.5, col="red") text(325, 850, "=0.82", cex=1.5, col="red") plot(demoProfiles$chr3.73159773.73160145$H3K27ac, type="l", col="blue", xlab="", ylab="", ylim=c(0, 125), main="chr3:73159773-73160145") par(new=TRUE) plot(demoProfiles$chr3.73159773.73160145$H3K4me1, type="l", col="darkgreen", xlab="Position", ylab="Coverage", ylim=c(0, 125)) legend("topright", c("H3K27ac","H3K4me1"), cex=1.2, col=c("blue","darkgreen"), lty=1) text(288, 92, "SPEARMAN_", cex=1.5, col="red") text(288, 78, "CORRELATION=0.95", cex=1.5, col="red") plot(demoProfiles$chr19.27739373.27739767$H3K27ac, type="l", col="blue", xlab="", ylab="", ylim=c(0, 350), main="chr19:27739373-27739767") par(new=TRUE) plot(demoProfiles$chr19.27739373.27739767$H3K4me1, type="l", col="darkgreen", xlab="Position", ylab="Coverage", ylim=c(0, 350)) legend("topright", c("H3K27ac","H3K4me1"), cex=1.2, col=c("blue","darkgreen"), lty=1) text(270, 225, "SPEARMAN_CORRELATION", cex=1.5, col="red") text(270, 192, "=0.62", cex=1.5, col="red") ## ----spearman_correlation_calculation, collapse=TRUE----------------------- metrics <- similarity(demoProfiles$chr2.70360770.70361098$H3K27ac, demoProfiles$chr2.70360770.70361098$H3K4me1) metrics$metric$SPEARMAN_CORRELATION metrics <- similarity(demoProfiles$chr8.43092918.43093442$H3K27ac, demoProfiles$chr8.43092918.43093442$H3K4me1) metrics$metric$SPEARMAN_CORRELATION metrics <- similarity(demoProfiles$chr3.73159773.73160145$H3K27ac, demoProfiles$chr3.73159773.73160145$H3K4me1) metrics$metric$SPEARMAN_CORRELATION metrics <- similarity(demoProfiles$chr19.27739373.27739767$H3K27ac, demoProfiles$chr19.27739373.27739767$H3K4me1) metrics$metric$SPEARMAN_CORRELATION ## ----libraryLoad----------------------------------------------------------- suppressMessages(library(similaRpeak)) ## ----profiles, collapse=TRUE----------------------------------------------- data(chr7Profiles) str(chr7Profiles) ## ----graphProfiles, echo=FALSE, fig.align='center', fig.height=6----------- plot(chr7Profiles$chr7.61968807.61969730$H3K27ac, type="l", col="blue", xlab="", ylab="", ylim=c(0, 700), main="chr7:61968807-61969730") par(new=TRUE) plot(chr7Profiles$chr7.61968807.61969730$H3K4me1, type="l", col="darkgreen", xlab="Position", ylab="Coverage in reads per million (RPM)", ylim=c(0, 700)) legend("topleft", c("H3K27ac","H3K4me1"), cex=1.2, col=c("blue","darkgreen"), lty=1) ## ----metricCalculation----------------------------------------------------- metrics <- similarity(chr7Profiles$chr7.61968807.61969730$H3K27ac, chr7Profiles$chr7.61968807.61969730$H3K4me1, ratioAreaThreshold=5, ratioMaxMaxThreshold=2, ratioIntersectThreshold=5, ratioNormalizedIntersectThreshold=2, diffPosMaxThresholdMinValue=10, diffPosMaxThresholdMaxDiff=100, diffPosMaxTolerance=0.01) ## ----metricReturn, collapse=TRUE------------------------------------------- metrics ## ----getInfo, collapse=TRUE------------------------------------------------ metrics$areaProfile1 metrics$areaProfile2 metrics$metrics$RATIO_INTERSECT ## ----graphProfilesNorm, echo=FALSE, fig.align='center', fig.height=6------- plot(chr7Profiles$chr7.61968807.61969730$H3K27ac* length(chr7Profiles$chr7.61968807.61969730$H3K27ac)/ sum(chr7Profiles$chr7.61968807.61969730$H3K27ac, na.rm=TRUE), type="l", col="blue", xlab="", ylab="", ylim=c(0, 3.5)) par(new=TRUE) plot(chr7Profiles$chr7.61968807.61969730$H3K4me1* length(chr7Profiles$chr7.61968807.61969730$H3K4me1)/ sum(chr7Profiles$chr7.61968807.61969730$H3K4me1, na.rm=TRUE), type="l", col="darkgreen", xlab="Position", ylab="Normalized Coverage (Coverage/Mean Coverage)", ylim=c(0, 3.5)) legend("topleft", c("H3K27ac","H3K4me1"), cex=1.2, col=c("blue","darkgreen"), lty=1) ## ----factory--------------------------------------------------------------- factory = MetricFactory$new(diffPosMaxTolerance=0.04) ## ----factoryDemo, collapse=TRUE-------------------------------------------- ratio_max_max <- factory$createMetric(metricType="RATIO_MAX_MAX", profile1=demoProfiles$chr2.70360770.70361098$H3K27ac, profile2=demoProfiles$chr2.70360770.70361098$H3K4me1) ratio_max_max ratio_normalized_intersect <- factory$createMetric( metricType="RATIO_NORMALIZED_INTERSECT", profile1=demoProfiles$chr2.70360770.70361098$H3K27ac, profile2=demoProfiles$chr2.70360770.70361098$H3K4me1) ratio_normalized_intersect ## ----sessionInfo, echo=FALSE----------------------------------------------- sessionInfo()