get.reads: Read genomic positions of sequencing data

Description Usage Arguments Value Author(s) See Also Examples

View source: R/get.reads.R

Description

Reads a bedfile containing positions of sequenced read aligned to a reference genome and creates a RangedData object.

Usage

1
get.reads(readsfile, filetype = c("bed", "bam"), chrcol = 1, startcol = 2, endcol = 3, idcol, zerobased = TRUE, sep = "\t", skip = 1, header = FALSE, ...)

Arguments

readsfile

name of bedfile giving the positions of aligned reads

#!!

filetype

Input file type. If "bam", the .bam file is read using scanBam, where flag options isUnmappedQuery=FALSE and isSecondaryAlignment=FALSE are used. Defaults to "bed"

# !!

chrcol

In which column in the reads bedfile there is the chromosome information (chromosome information in the file should be in string format, e.g. "chrX"). Ignored if filetype = "bam".

startcol

In which column there are the starting positions of the reads. Ignored if filetype = "bam".

endcol

In which column there are the end positions of the reads. Ignored if filetype = "bam".

idcol

In which column there are read identifiers. For single-end data it is optionally. For paired-end data it is required for some functionalities. The two reads of one pair need to have the same ID. Ignored if filetype = "bam" (the ID column is automatically included then). If read IDs include "#0/1" and "#0/2" in the end (indicating read 1 and read 2 of a pair), those characters will be removed from the IDs.

zerobased

if TRUE, start coordinates in readsfile are assumed to be 0-based and are then converted to 1-based system by adding 1. If FALSE, coordinates are not shifted. In this case they should already be 1-based in readsfile. Ignored if filetype = "bam", since scanBam converts 0-based to 1-based coordinates.

sep

Column separator character, defaults to tabs. Ignored if filetype = "bam".

skip

Number of lines of the bedfile to skip before beginning to read data; defaults to 1. Ignored if filetype = "bam".

header

A logical value indicating whether the file contains the names of the variables as its first line; defaults to FALSE. Ignored if filetype = "bam".

...

Further arguments passed to read.delim. Ignored if filetype = "bam".

Value

A RangedData table holding the read positions

Author(s)

Manuela Hummel m.hummel@dkfz.de

See Also

get.targets

Examples

1
2
3
exptPath <- system.file("extdata", package="TEQC")
readsfile <- file.path(exptPath, "ExampleSet_Reads.bed")
reads <- get.reads(readsfile, idcol=4, skip=0)

Example output

Loading required package: BiocGenerics
Loading required package: parallel

Attaching package: 'BiocGenerics'

The following objects are masked from 'package:parallel':

    clusterApply, clusterApplyLB, clusterCall, clusterEvalQ,
    clusterExport, clusterMap, parApply, parCapply, parLapply,
    parLapplyLB, parRapply, parSapply, parSapplyLB

The following objects are masked from 'package:stats':

    IQR, mad, sd, var, xtabs

The following objects are masked from 'package:base':

    Filter, Find, Map, Position, Reduce, anyDuplicated, append,
    as.data.frame, cbind, colMeans, colSums, colnames, do.call,
    duplicated, eval, evalq, get, grep, grepl, intersect, is.unsorted,
    lapply, lengths, mapply, match, mget, order, paste, pmax, pmax.int,
    pmin, pmin.int, rank, rbind, rowMeans, rowSums, rownames, sapply,
    setdiff, sort, table, tapply, union, unique, unsplit, which,
    which.max, which.min

Loading required package: IRanges
Loading required package: S4Vectors
Loading required package: stats4

Attaching package: 'S4Vectors'

The following object is masked from 'package:base':

    expand.grid

Loading required package: Rsamtools
Loading required package: GenomeInfoDb
Loading required package: GenomicRanges
Loading required package: Biostrings
Loading required package: XVector

Attaching package: 'Biostrings'

The following object is masked from 'package:base':

    strsplit

Loading required package: hwriter
[1] "read 19546 sequenced reads"

TEQC documentation built on Nov. 8, 2020, 8:07 p.m.

Related to get.reads in TEQC...