simpleRNASeq,BamFileList,RnaSeqParam-method {easyRNASeq}R Documentation

simpleRNASeq method

Description

This function is a wrapper around the more low level functionalities of the package. It is the simplest way to get a RangedSummarizedExperiment object from a set of bam files. RangedSummarizedExperiment are containers meant to hold any Next-Generation Sequencing experiment results and metadata. The simpleRNASeq method replaces the easyRNASeq function to simplify the usability. It does the following:

Usage

## S4 method for signature 'BamFileList,RnaSeqParam'
simpleRNASeq(
  bamFiles = BamFileList(),
  param = RnaSeqParam(),
  nnodes = 1,
  verbose = TRUE,
  override = FALSE
)

Arguments

bamFiles

a BamFileList object

param

RnaSeqParam a RnaSeqParam object that describes the RNA-Seq experimental setup.

nnodes

The number of CPU cores to use in parallel

verbose

a logical to be report progress or not.

override

Should the provided parameters override the detected ones

Value

returns a RangedSummarizedExperiment object.

Author(s)

Nicolas Delhomme

See Also

Examples


# the data
tdir <- tutorialData()
annot <- fetchData("Drosophila_melanogaster.BDGP5.77.with-chr.gtf.gz")

 # create the BamFileList, get the BAM and BAI index files from the Bioc cache
 filenames <- dir(tdir,pattern="[A,T].*\\.bam$",full.names=TRUE)
 indexnames <- sapply(paste0(sub(".*_","",basename(filenames)),".bai"),fetchData)
 bamFiles <- getBamFileList(filenames,indexnames)

  # create the AnnotParam
  annotParam <- AnnotParam(annot,type="gtf")

  # create the RnaSeqParam
  rnaSeqParam <- RnaSeqParam(annotParam=annotParam)

  # get a RangedSummarizedExperiment containing the counts table
  sexp <- simpleRNASeq(
    bamFiles=bamFiles,
    param=rnaSeqParam,
    verbose=TRUE
  )

  # get the counts
  assays(sexp)$exons


[Package easyRNASeq version 2.28.0 Index]