extract_bams: Detect allele-specific methylation from a bam file

In markrobinsonuzh/DAMEfinder: Finds DAMEs - Differential Allelicly MEthylated regions

Detect allele-specific methylation from a bam file

Description

The function takes a bam (from bismark) and vcf file for each sample. For each SNP contained in the vcfile it calculates the proportion of methylated reads for each CpG site at each allele. At the end it returns (saves to working directory) a GRanges list, where each GRanges contains all the CpG sites overlapping the reads containing a specific SNP.

Usage

extract_bams(
  bamFiles,
  vcfFiles,
  sampleNames,
  referenceFile,
  coverage = 4,
  cores = 1,
  verbose = TRUE
)

Arguments

bamFiles	List of bam files.
vcfFiles	List of vcf files.
sampleNames	Names of files in the list.
referenceFile	fasta file used to generate the bam files. Or DNAStringSet with DNA sequence.
coverage	Minimum number of reads covering a CpG site on each allele. Default = 2.
cores	Number of cores to use. See package parallel for description of core. Default = 1.
verbose	Default = TRUE

Value

A list of GRanges for each sample. Each list is saved in a separate .rds file.

Examples

DATA_PATH_DIR <- system.file('extdata', '.', package = 'DAMEfinder')
get_data_path <- function(file_name) file.path(DATA_PATH_DIR, file_name)
bamFiles <- get_data_path('NORM1_chr19_trim.bam')
vcfFiles <- get_data_path('NORM1.chr19.trim.vcf')
sampleNames <- 'NORM1'

#referenceFile 
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

GRanges_list <- extract_bams(bamFiles, vcfFiles, sampleNames, dna)

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