GenomicData-class          package:IRanges          R Documentation

_D_a_t_a _o_n _a _G_e_n_o_m_e

_D_e_s_c_r_i_p_t_i_o_n:

     'GenomicData' extends 'RangedData' to more conveniently manipulate
     data on genomic ranges. The spaces are now called chromosomes (but
     could still refer to some other type of sequence). The annotation
     refers to the genome and there is formal treatment of an optional
     'strand' variable.

_A_c_c_e_s_s_o_r_s:

     The 'GenomicData' class essentially adds a set of convenience
     accessors on top of 'RangedData'.   In the code snippets below,
     'x' is a 'RangedData' object.


      'chrom(x)': Gets the chromosome names (a factor) over the ranges
          in 'x'.

      'genome(x)': Gets the genome for the ranges in 'x'; a simple
          wrapper around 'annotation(x)'.

      'strand(x)': Gets the strand factor in 'x', or, if 'x' is
          missing, return an empty factor with the standard levels:
          '-', '+' and '*', referring to the negative, positive and
          both strands, respectively. Any strand factor stored in
          'GenomicData' should have those levels. 'NA''s are allowed;
          the value is all 'NA' if no strand has been specified. 


_C_o_n_s_t_r_u_c_t_o_r:


      'GenomicData(ranges, ..., strand = NULL, chrom = NULL, genome =
          NULL)': Constructs a 'GenomicData' instance with the given
          'ranges' and variables in '...' (see the 'RangedData'
          constructor). If non-'NULL', 'strand' specifies the strand of
          each range. It should be a character vector or factor of
          length equal to that of 'ranges'. All values should be either
          '-', '+', '*' or 'NA'. To get these levels, call
          'levels(strand())'. 'chrom' is analogous to 'splitter' in
          'RangedData'; if non-'NULL' it should be coercible to a
          factor indicating how the ranges, variables and strand should
          be split up  across the chromosomes. The 'genome' argument
          should be a scalar string and is treated as the 'RangedData'
          annotation. See the examples.


_C_o_e_r_c_i_o_n:


      'as(from, "GenomicData")': coerces 'from' to a 'GenomicData',
          according to its class:

          _X_R_l_e The bounds of the runs become the ranges and the values
               become a column named 'score'.


_N_o_t_e:

     This may move to another package in the future, as it is not quite
     general enough for IRanges.

_A_u_t_h_o_r(_s):

     Michael Lawrence

_S_e_e _A_l_s_o:

     'RangedData', on which this class is based.

_E_x_a_m_p_l_e_s:

       range1 <- IRanges(start=c(1,2,3), end=c(5,2,8))

       ## just ranges
       gr <- GenomicData(range1) 

       ## with a genome (annotation)
       gr <- GenomicData(range1, genome = "hg18")
       genome(gr) ## "hg18"

       ## with some data
       filter <- c(1L, 0L, 1L)
       score <- c(10L, 2L, NA)
       strand <- factor(c("+", NA, "-"), levels = levels(strand()))
       gr <- GenomicData(range1, score, genome = "hg18")
       gr[["score"]]
       strand(gr) ## all NA
       gr <- GenomicData(range1, score, filt = filter, strand = strand)
       gr[["filt"]]
       strand(gr) ## equal to 'strand'

       range2 <- IRanges(start=c(15,45,20,1), end=c(15,100,80,5))
       ranges <- c(range1, range2)
       score <- c(score, c(0L, 3L, NA, 22L)) 
       chrom <- paste("chr", rep(c(1,2), c(length(range1), length(range2))), sep="")
       
       gr <- GenomicData(ranges, score, chrom = chrom, genome = "hg18")
       chrom(gr) # equal to 'chrom'
       gr[["score"]] # unlists over the chromosomes
       gr[1][["score"]] # equal to score[1:3]

