countFragmentOverlaps {FourCSeq}R Documentation

Count fragment overlaps

Description

countFragmentOverlaps counts the number of reads mapping to each fragment end in rowRanges of the FourC object.

Usage

countFragmentOverlaps(object, trim=0, minMapq=0, shift=0)

Arguments

object

A FourC object.

trim

Number of bases that should be trimmed at the start of a read. This is necessary if the read still contains the restriction enzyme sequence. Default is 0 bases.

minMapq

Minimum mapping quality required for counting the read. Default is 0. If set to negative values the filtering step is skipped.

shift

Maximum difference in starts or ends between read and fragment positions. Default is 0.

Value

Updated FourC object that contains two new assays countsLeftFragmentEnd and countsRightFragmentEnd with the count values at the respective fragment end.

Author(s)

Felix A. Klein, felix.klein@embl.de

See Also

FourC, findViewpointFragments, countFragmentOverlapsSecondCutter

Examples




metadata <- list(projectPath=tempdir(),
                 fragmentDir="re_fragments",
                 referenceGenomeFile=system.file("extdata/dm3_chr2L_1-6900.fa", 
                                                 package="FourCSeq"),
                 reSequence1="GATC",
                 reSequence2="CATG",
                 primerFile=system.file("extdata/primer.fa", 
                                        package="FourCSeq"),
                 bamFilePath=system.file("extdata/bam", package="FourCSeq"))

colData <- DataFrame(viewpoint = "testdata", 
                     condition = factor(rep(c("WE_68h", "MESO_68h", "WE_34h"),                    
                                            each=2),
                                        levels = c("WE_68h", "MESO_68h", "WE_34h")),
                     replicate = rep(c(1, 2), 
                                     3),
                     bamFile = c("CRM_ap_ApME680_WE_6-8h_1_testdata.bam", 
                                 "CRM_ap_ApME680_WE_6-8h_2_testdata.bam",       
                                 "CRM_ap_ApME680_MESO_6-8h_1_testdata.bam", 
                                 "CRM_ap_ApME680_MESO_6-8h_2_testdata.bam", 
                                 "CRM_ap_ApME680_WE_3-4h_1_testdata.bam",
                                 "CRM_ap_ApME680_WE_3-4h_2_testdata.bam"),
                     sequencingPrimer="first")

fc <- FourC(colData, metadata)
fc

fc <- addFragments(fc)

findViewpointFragments(fc) 

fc <- addViewpointFrags(fc)
fc

fc <- countFragmentOverlaps(fc, trim=4, minMapq=30)
fc


[Package FourCSeq version 1.24.0 Index]