methyl_circle_plot: Draw methylation circle plot

View source: R/methyl_circle_plot.R

methyl_circle_plotR Documentation

Draw methylation circle plot

Description

Draws CpG site methylation status as points, in reads containing a specific SNP. Generates one plot per bam file.

Usage

methyl_circle_plot(
  snp,
  vcfFile,
  bamFile,
  refFile,
  build = "hg19",
  dame = NULL,
  letterSize = 2.5,
  pointSize = 3,
  sampleName = "sample1",
  cpgsite = NULL,
  sampleReads = FALSE,
  numReads = 20
)

Arguments

snp

GRanges object containing SNP location.

vcfFile

vcf file.

bamFile

bismark bam file path.

refFile

fasta reference file path. Or DNAStringSet with DNA sequence.

build

genome build used. default = "hg19"

dame

(optional) GRanges object containing a region to plot.

letterSize

Size of alleles drawn in plot. Default = 2.5.

pointSize

Size of methylation circles. Default = 3.

sampleName

FIX?: this is to save the vcf file to not generate it every time you run the function.

cpgsite

(optional) GRanges object containing a single CpG site location of interest.

sampleReads

Whether a subset of reads should be plotted. Default = FALSE.

numReads

Number of reads to plot per allele, if sampleReads is TRUE. Default = 20

Value

Plot

Examples

DATA_PATH_DIR <- system.file('extdata', '.', package = 'DAMEfinder')

get_data_path <- function(file_name) file.path(DATA_PATH_DIR, file_name)
bam_files <- get_data_path('NORM1_chr19_trim.bam')
vcf_files <- get_data_path('NORM1.chr19.trim.vcf')
sample_names <- 'NORM1'
#reference_file
suppressPackageStartupMessages({library(BSgenome.Hsapiens.UCSC.hg19)})
genome <- BSgenome.Hsapiens.UCSC.hg19
seqnames(genome) <- gsub("chr","",seqnames(genome))
dna <- DNAStringSet(genome[[19]], use.names = TRUE)
names(dna) <- 19

snp <- GenomicRanges::GRanges(19, IRanges::IRanges(292082, width = 1))
methyl_circle_plot(snp = snp,
 vcfFile = vcf_files,
 bamFile = bam_files,
 refFile = dna,
 sampleName = sample_names)


markrobinsonuzh/DAMEfinder documentation built on April 7, 2023, 6:37 a.m.