text: user8600: I have a problem in R. I have the following dataSet (the first three rows shown) (the 5th number is the methylation level in its region):
  chr1    15864   15866   FALSE   894 +
  chr1    534241  534243  FALSE   921 -
  chr1    710096  710098  FALSE   729 +
  
  and then I'm trying to get the sequences from every row in the dataSet (427.000 rows in total):
  library(BSgenome.Hsapiens.UCSC.hg19)
  genome <- BSgenome.Hsapiens.UCSC.hg19
  chr<-as.matrix(as.character(mydataSet[,1]))
  start<-as.matrix(as.integer(as.character(mydataSet[,2]))-50)
  end<-as.matrix(as.integer(as.character(mydataSet[,3]))+50)
  Sequences_for_50_plus_bases<-getSeq(genome,chr,start=start,end=end)
  
  I'd like to find a new fuction where I can have this dataSet, and find for its rows (which represent sequences in the human genome) if they refer to a promoter, or an enhancer or any other DNA regions' information, because I'd like to compare the differences between the regions' methylation levels.
  Is there a way to do this?
  
  user8870: You can get specific sequences efficiently using the following packages:
  library(GenomicFeatures)
  library(Biostrings)
  library(tidyverse)
  library(data.table)
  library(dtplyr)
  library(BSgenome.Hsapiens.UCSC.hg19)
  library(TxDb.Hsapiens.UCSC.hg19.knownGene)
  
  Methylation_data <- tribble(
    ~seqnames,~start,~end,~Logical,~"Methylation Level",~strand,
    "chr1",15864,15866,FALSE,894,"+",
    "chr1",534241,534243,FALSE,921,"-",
    "chr1",710096,710098,FALSE,729,"+") # Arrange your dataset
  
  Methylation_GRanges <- GRanges(  # Construct a GRanges objects
    seqnames = Methylation_data %>% pull(seqnames), # First pull in the chromosome locations
    IRanges(start = (Methylation_data %>% pull(start))-50, # then the start position minus 50 nt
            end = (Methylation_data %>% pull(end))+50), # the end position plus 50 nt
    strand = Methylation_data %>% pull(strand) # the strand your sequence is located on
  )
  > getSeq(Hsapiens,Methylation_GRanges)
  A DNAStringSet instance of length 3
  width seq
  [1]   103 GCTGCTGCTTCTCCAGCTTTCGCTCCTTCATGCTGCGCAGC...CTTGGCGGATGGACTCTAGCAGAGTGGCCAGCCACCGGAG
  [2]   103 ACAGTTCCAATGTAATCAGAGAGAACATCACACACACACCA...ACCTGAGCAGCACTCTGCAAAGCTGTCAAGGCGGTGAAGC
  [3]   103 GATAAAATAAAGCTTAGATTGGAAAAAATATTTAAGATTCT...GGAAGCTGAGTAATTGTATGTTCAAATACTTGCAAAACAT
  
  To write the sequences to a file, use writeXStringSet()
  For the second part of what you need, GenomicFeatures has a number of functions that obtain the annotated sequences of the whole genome by its type:
  transcriptsBy(x, by=c("gene", "exon", "cds"), ...)
  exonsBy(x, by=c("tx", "gene"), ...)
  cdsBy(x, by=c("tx", "gene"), ...)
  intronsByTranscript(x, ...)
  fiveUTRsByTranscript(x, ...)
  threeUTRsByTranscript(x, ...)
  
  And this is how you use one of them to find whether your sequences overlap a particular type of DNA:
  Find_overlap <- function(Ref_Data,
                           Overlap_Message){
  
    Ref_Data <- as.data.table(Ref_Data) # Convert to data.table for efficiency
  
    Methylation_data <- tribble(
      ~seqnames,~start,~end,~Logical,~"Methylation Level",~strand,
      "chr1",15864,15866,FALSE,894,"+",
      "chr1",534241,534243,FALSE,921,"-",
      "chr1",710096,710098,FALSE,729,"+") %>%
      mutate(Overlap = "-",
             start = start-50,
             end = end+50) %>%
      as.data.table()
  
    Methylation_data[Ref_Data, 
                     Overlap := Overlap_Message, 
                     on = .(start >= start, 
                            end <=  end)] # wherever the range overlaps, assign an identifier to it
  }
  
  
  Methylated_in_5UTR <- fiveUTRsByTranscript(TxDb.Hsapiens.UCSC.hg19.knownGene) %>% 
    Find_overlap(., "5' UTR") %>% 
    print
  
  seqnames  start    end Logical Methylation Level strand Overlap
  1:     chr1  15814  15916   FALSE               894      +       -
  2:     chr1 534191 534293   FALSE               921      -       -
  3:     chr1 710046 710148   FALSE               729      +  5' UTR
  ```
  
  
original_text: 8600: <p>I have a problem in R. I have the following dataSet (the first three rows shown) (the 5th number is the methylation level in its region):</p>
  
  <pre><code>chr1    15864   15866   FALSE   894 +
  chr1    534241  534243  FALSE   921 -
  chr1    710096  710098  FALSE   729 +
  </code></pre>
  
  <p>and then I'm trying to get the sequences from every row in the dataSet (427.000 rows in total):</p>
  
  <pre><code>library(BSgenome.Hsapiens.UCSC.hg19)
  genome &lt;- BSgenome.Hsapiens.UCSC.hg19
  chr&lt;-as.matrix(as.character(mydataSet[,1]))
  start&lt;-as.matrix(as.integer(as.character(mydataSet[,2]))-50)
  end&lt;-as.matrix(as.integer(as.character(mydataSet[,3]))+50)
  Sequences_for_50_plus_bases&lt;-getSeq(genome,chr,start=start,end=end)
  </code></pre>
  
  <p>I'd like to find a new fuction where I can have this dataSet, and find for its rows (which represent sequences in the human genome) if they refer to a promoter, or an enhancer or any other DNA regions' information, because I'd like to compare the differences between the regions' methylation levels.</p>
  
  <p>Is there a way to do this?</p>
   8870: <p>You can get specific sequences efficiently using the following packages:</p>
  
  <pre><code>library(GenomicFeatures)
  library(Biostrings)
  library(tidyverse)
  library(data.table)
  library(dtplyr)
  library(BSgenome.Hsapiens.UCSC.hg19)
  library(TxDb.Hsapiens.UCSC.hg19.knownGene)
  
  Methylation_data &lt;- tribble(
    ~seqnames,~start,~end,~Logical,~"Methylation Level",~strand,
    "chr1",15864,15866,FALSE,894,"+",
    "chr1",534241,534243,FALSE,921,"-",
    "chr1",710096,710098,FALSE,729,"+") # Arrange your dataset
  
  Methylation_GRanges &lt;- GRanges(  # Construct a GRanges objects
    seqnames = Methylation_data %&gt;% pull(seqnames), # First pull in the chromosome locations
    IRanges(start = (Methylation_data %&gt;% pull(start))-50, # then the start position minus 50 nt
            end = (Methylation_data %&gt;% pull(end))+50), # the end position plus 50 nt
    strand = Methylation_data %&gt;% pull(strand) # the strand your sequence is located on
  )
  &gt; getSeq(Hsapiens,Methylation_GRanges)
  A DNAStringSet instance of length 3
  width seq
  [1]   103 GCTGCTGCTTCTCCAGCTTTCGCTCCTTCATGCTGCGCAGC...CTTGGCGGATGGACTCTAGCAGAGTGGCCAGCCACCGGAG
  [2]   103 ACAGTTCCAATGTAATCAGAGAGAACATCACACACACACCA...ACCTGAGCAGCACTCTGCAAAGCTGTCAAGGCGGTGAAGC
  [3]   103 GATAAAATAAAGCTTAGATTGGAAAAAATATTTAAGATTCT...GGAAGCTGAGTAATTGTATGTTCAAATACTTGCAAAACAT
  </code></pre>
  
  <p>To write the sequences to a file, use <code>writeXStringSet()</code></p>
  
  <p>For the second part of what you need, <code>GenomicFeatures</code> has a number of functions that obtain the annotated sequences of the whole genome by its type:</p>
  
  <pre><code>transcriptsBy(x, by=c("gene", "exon", "cds"), ...)
  exonsBy(x, by=c("tx", "gene"), ...)
  cdsBy(x, by=c("tx", "gene"), ...)
  intronsByTranscript(x, ...)
  fiveUTRsByTranscript(x, ...)
  threeUTRsByTranscript(x, ...)
  </code></pre>
  
  <p>And this is how you use one of them to find whether your sequences overlap a particular type of DNA:</p>
  
  <pre><code>Find_overlap &lt;- function(Ref_Data,
                           Overlap_Message){
  
    Ref_Data &lt;- as.data.table(Ref_Data) # Convert to data.table for efficiency
  
    Methylation_data &lt;- tribble(
      ~seqnames,~start,~end,~Logical,~"Methylation Level",~strand,
      "chr1",15864,15866,FALSE,894,"+",
      "chr1",534241,534243,FALSE,921,"-",
      "chr1",710096,710098,FALSE,729,"+") %&gt;%
      mutate(Overlap = "-",
             start = start-50,
             end = end+50) %&gt;%
      as.data.table()
  
    Methylation_data[Ref_Data, 
                     Overlap := Overlap_Message, 
                     on = .(start &gt;= start, 
                            end &lt;=  end)] # wherever the range overlaps, assign an identifier to it
  }
  
  
  Methylated_in_5UTR &lt;- fiveUTRsByTranscript(TxDb.Hsapiens.UCSC.hg19.knownGene) %&gt;% 
    Find_overlap(., "5' UTR") %&gt;% 
    print
  </code></pre>
  
  <pre><code>seqnames  start    end Logical Methylation Level strand Overlap
  1:     chr1  15814  15916   FALSE               894      +       -
  2:     chr1 534191 534293   FALSE               921      -       -
  3:     chr1 710046 710148   FALSE               729      +  5' UTR
  <span class="math-container">```</span>
  </code></pre>
  