extract_bams: Detect allele-specific methylation from a bam file

View source: R/split_bams.R

extract_bamsR Documentation

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.