Crazy pipe operator

Positions of all insertion/deletion

                Qname start_pos insert_1 insert_2 insert_3 insert_4
1      SRR2426363.615   3356120  3356201       NA       NA       NA
2      SRR2426363.698   4015323  4015492       NA       NA       NA
3     SRR2426363.1034   3912080  3912120       NA       NA       NA
4     SRR2426363.4217   4382662  4382826  4382838       NA       NA
5     SRR2426363.4732    979022   979071       NA       NA       NA
6     SRR2426363.6348   4017239  4017395       NA       NA       NA

Data manipulation skills with plyr/dplyr

#---  **ply(): a*ply, l*ply, d*ply
#---  select(): focus on a subset of variables
#---  filter(): focus on a subset of rows
#---  mutate(): add new columns
#---  arrange(): re-order the rows
#---  colwise(): functional programming example
#---  [[: subsetting
#---  ->: right assignment operator
#---  ->>: assign through parents environments
#---  %>%
#---  workflow

---

Let’s do it!

---

Option 1: read data remotely

myData = readLines('https://raw.githubusercontent.com/MingChen0919/KBS_workshop_Michigan_2016/master/SRR2426363.mapped.sam')

Option 2: download data and then read data locally

curl -O https://raw.githubusercontent.com/MingChen0919/KBS_workshop_Michigan_2016/master/SRR2426363.mapped.sam
myData = readLines('SRR2426363.mapped.sam')

Do everything in one single pipe line.

• Pipeline starts with your data
• Input -> Output: No intermediate variables generated
• Pipeline is extensible

myData %>% tail(-1) %>% ## remove the first row because it contains field names
  (function(x){
    ## wrap the str_split_fixed() function so it returns a character
    str_split_fixed_2 = function(string, pattern, n){
      x = str_split_fixed(string, pattern, n)
      return(as.character(x))
    }
    ldply(x, str_split_fixed_2, '\t', 20)[, 1:11]
  }) %>% 
  (function(x){
    colnames(x) = c('QNAME', 'FLAG', 'RNAME', 'POS', 'MAPQ',
                    'CIGAR', 'RNEXT', 'PNEXT', 'TLEN', 'SEQ', 'QUAL')
    rownames(x) = paste0('seq', 1:nrow(x))
    return(x)
  }) %>% 
  mutate(FLAG = as.numeric(FLAG),
         POS = as.numeric(POS),
         MAPQ   = as.numeric(MAPQ),
         PNEXT  = as.numeric(PNEXT),
         TLEN   = as.numeric(TLEN)) %>%
  (function(x){
    filter(x, MAPQ > 30 & MAPQ < 55) %>%
      select(QNAME, RNAME, POS, MAPQ, CIGAR) %>% nrow() 
    return(x)
  }) %>% 
  (function(x){
    arrange(x, MAPQ) %>% 
      select(QNAME, RNAME, POS, MAPQ, CIGAR) %>% head() 
    return(x)
  }) %>%
  (function(x){
    filter(x, MAPQ > 30) %>% 
      select(QNAME, RNAME, POS, MAPQ, CIGAR) %>%  
      arrange(POS, MAPQ) %>% head() 
    return(x)
  }) %>% 
  (function(x){
    select(x, CIGAR) %>% 
      laply(grep, pattern='[I]') %>% 
      slice(.data = x) %>%
      (function(x){
        x ->> df_with_Insertion
        select(x, CIGAR) ->> CIGAR_with_Insertion
        select(x, QNAME) ->> QNAME 
        select(x, POS) ->> POS
        return(select(x, CIGAR))
      })
  }) %>% 
  llply(str_split, 'I') %>%
  `[[`('CIGAR') %>%
  llply(head, -1) %>%
  llply(gsub, pattern='[A-Z]', replacement='+') %>%
  (function(x){
    eval_string = function(var1) eval(parse(text=var1))  ## define a function
    eval_string_plus = function(var2) laply(var2, eval_string) ## define another function
    llply(x, eval_string_plus)
  })  %>% 
  llply(cumsum) %>%
  (function(x){
    maxN = max(lengths(x))
    newX = vector('numeric', length = maxN)
    ldply(x, `[`, 1:maxN)
  }) %>% 
  (function(x){
    `colnames<-`(x, paste0('insert_', 1:ncol(x)))
    ## return(x)
  }) %>%
  colwise(.fun=`+`)(unlist(POS)) %>%
  mutate(Qname = unlist(QNAME), start_pos = unlist(POS)) %>%
  `[`(c('Qname', 'start_pos', 'insert_1', 'insert_2', 'insert_3', 'insert_4'))

---

Let’s break it down!

---

##==== Step 1: split each line by tab and return a tabular data structure  ====
myData %>% tail(-1) %>% ## remove the first row because it contains field names
  (function(x){
    ## wrap the str_split_fixed() function so it returns a character
    str_split_fixed_2 = function(string, pattern, n){
      x = str_split_fixed(string, pattern, n)
      return(as.character(x))
    }
    ldply(x, str_split_fixed_2, '\t', 20)[, 1:11]
  }) %>% 

##==== Step 2:  set column names =======
  (function(x){
    colnames(x) = c('QNAME', 'FLAG', 'RNAME', 'POS', 'MAPQ',
                    'CIGAR', 'RNEXT', 'PNEXT', 'TLEN', 'SEQ', 'QUAL')
    rownames(x) = paste0('seq', 1:nrow(x))
    return(x)
  }) %>% 

##==== Step:3 convert strings to numbers ====
  mutate(FLAG = as.numeric(FLAG),
         POS = as.numeric(POS),
         MAPQ   = as.numeric(MAPQ),
         PNEXT  = as.numeric(PNEXT),
         TLEN   = as.numeric(TLEN)) %>%

##==== Exploring  ====
  (function(x){
    filter(x, MAPQ > 30 & MAPQ < 55) %>%
      select(QNAME, RNAME, POS, MAPQ, CIGAR) %>% nrow() 
    return(x)
  }) %>% 

##==== Exploring ====
  (function(x){
    arrange(x, MAPQ) %>% 
      select(QNAME, RNAME, POS, MAPQ, CIGAR) %>% head() 
    return(x)
  }) %>%

##==== Exploring ====
  (function(x){
    filter(x, MAPQ > 30) %>% 
      select(QNAME, RNAME, POS, MAPQ, CIGAR) %>%  
      arrange(POS, MAPQ) %>% head() 
    return(x)
  }) %>% 

##==== Step 4: select rows that have insertions ====
  (function(x){
    select(x, CIGAR) %>% 
      laply(grep, pattern='[I]') %>% 
      slice(.data = x) %>%
      (function(x){
        x ->> df_with_Insertion
        select(x, CIGAR) ->> CIGAR_with_Insertion
        select(x, QNAME) ->> QNAME 
        select(x, POS) ->> POS
        return(select(x, CIGAR))
      })
  }) %>%

##==== Step 5: split CIGAR by insertion letter 'I'  ====
  llply(str_split, 'I') %>%

##==== Step 6: Extract data from a nested list ====
  `[[`('CIGAR') %>%

##==== Step 7: delete the last element ====
  llply(head, -1) %>%

##==== Step 8: build calculation expression by replace remaining CIGAR letters with '+' ====
  llply(gsub, pattern='[A-Z]', replacement='+') %>%

##==== Step 9: evalate the expression ====
 (function(x){
    eval_string = function(var1) eval(parse(text=var1))  ## define a function
    eval_string_plus = function(var2) laply(var2, eval_string) ## define another function
    llply(x, eval_string_plus)
  })  %>% 

##==== Step 10: get positions ====
 llply(cumsum) %>%
  (function(x){
    maxN = max(lengths(x))
    newX = vector('numeric', length = maxN)
    ldply(x, `[`, 1:maxN)
  }) %>% 

##==== Step 10: column names and row names ====
  (function(x){
    `colnames<-`(x, paste0('insert_', 1:ncol(x)))
    ## return(x)
  }) %>%

##==== Step 11: add starting position =======
colwise(.fun=`+`)(unlist(POS)) %>%

##==== Step 11: restructure the data frame to get final results ====
  mutate(Qname = unlist(QNAME), start_pos = unlist(POS)) %>%
  `[`(c('Qname', 'start_pos', 'insert_1', 'insert_2', 'insert_3', 'insert_4'))