2
$\begingroup$

I have a file with sample data like this:

1   90042011    90042031    AGTGCCACCAGGTGGGGCCG    90042017    2   5008    T:5006  C:2     90184334    90184354    GAGGCCACAAGAGGGCACAA    90184337    2   5008    C:5007  T:1     ß@  5
1   90042011    90042031    AGTGCCACCAGGTGGGGCCG    90042017    2   5008    T:5006  C:2     90184334    90184354    GAGGCCACAAGAGGGCACAA    90184336    2   5008    A:5007  G:1     ßß  2
1   90042011    90042031    AGTGCCACCAGGTGGGGCCG    90042030    2   5008    C:5007  A:1     90184334    90184354    GAGGCCACAAGAGGGCACAA    90184337    2   5008    C:5007  T:1     @@  1
1   90042011    90042031    AGTGCCACCAGGTGGGGCCG    90042030    2   5008    C:5007  A:1     90184334    90184354    GAGGCCACAAGAGGGCACAA    90184336    2   5008    A:5007  G:1     @ß  3
1   90042011    90042031    AGTGCCACCAGGTGGGGCCG    90042029    2   5008    G:5006  T:2     90184334    90184354    GAGGCCACAAGAGGGCACAA    90184337    2   5008    C:5007  T:1     @@  1
1   90042011    90042031    AGTGCCACCAGGTGGGGCCG    90042029    2   5008    G:5006  T:2     90184334    90184354    GAGGCCACAAGAGGGCACAA    90184336    2   5008    A:5007  G:1     @ß  3
1   94853396    94853416    GCCGCCACTAGATGGTGCTA    94853398    2   5008    T:5007  C:1     94969254    94969274    CAGTCCTCTAGAGGGAGCCC    94969266    2   5008    G:5006  A:2     ß@  5
2   103584283   103584303   TGGGCCACTAGGAGGCACTG    103584300   2   5008    C:5006  T:2     103841436   103841456   CATGCCACAAGAGGGCATCA    103841456   2   5008    A:5006  G:2     @ß  3
2   103584283   103584303   TGGGCCACTAGGAGGCACTG    103584300   2   5008    C:5006  T:2     103841436   103841456   CATGCCACAAGAGGGCATCA    103841445   2   5008    C:5006  T:2     @@  1
2   103584283   103584303   TGGGCCACTAGGAGGCACTG    103584284   2   5008    A:5007  G:1     103841436   103841456   CATGCCACAAGAGGGCATCA    103841456   2   5008    A:5006  G:2     ßß  2
2   103584283   103584303   TGGGCCACTAGGAGGCACTG    103584284   2   5008    A:5007  G:1     103841436   103841456   CATGCCACAAGAGGGCATCA    103841445   2   5008    C:5006  T:2     ß@  5

As seen, some lines have repeating values for $2 and $3, line 1-6 and line 8-11 in this case.

$20 of my file contains codes of two characters representing each line, in the format "@@","ßß", "@ß","@N","ß@","ßN","N@","Nß" or "NN". Each line can have one of these 9 possible codes.

I want to write a script, which splits the characters of $20 (@,ß or N) of the lines which are repeated (1-6 & 8-11 here), and reduces them into one line in the following manner:

Count the number of @, ß or N in the first and second positions respectively of the repeated lines. For the first and second characters of each repeated line, if:

  1. @ occurs most, print @
  2. ß occurs most, print ß
  3. N occurs most, print N
  4. @ occurs equal to ß, print N
  5. @ occurs equal to N, print @
  6. ß occurs equal to N, print ß

In this sample dataset with repeated lines, the output file generated should be something like this, only columns $1, $2, $3, $4, $11,$12,$13, and a new column with the new "@","ß","N" codes after said counting is reported:

1   90042011    90042031    AGTGCCACCAGGTGGGGCCG    90184334    90184354    GAGGCCACAAGAGGGCACAA    @N  
1   94853396    94853416    GCCGCCACTAGATGGTGCTA    94969254    94969274    CAGTCCTCTAGAGGGAGCCC    ß@  
2   103584283   103584303   TGGGCCACTAGGAGGCACTG    103841436   103841456   CATGCCACAAGAGGGCATCA    NN  

In the first repeated case(1-6), for the first character @ occurs four times and ß occurs twice and for the second character @ and ß occur equal number of times, so the final code is @N. In the second repeated case(8-11), @ and ß for both character positions occur equally, so the final code is NN.

$\endgroup$
2
$\begingroup$

OK. This is a bit convoluted, but it works. There might be simpler solutions...

First we'll write a function that takes a data frame and a character position (ie; first or second position in the string), groups the data frame by the second and third columns, counts occurrences of each character in the specified position of the 18th column, and returns a vector according to your ranking scheme (this doesn't include counts of N since it doesn't occur in your example data, but it's a simple extension of the case_when line to add it):

library(dplyr)
library(tidyr)
library(stringr)

summarise_characters <- function(df, pos) {
  df %>% mutate(character=str_sub(X18, pos, pos)) %>% 
  group_by(X2, X3, character) %>% 
  summarise(n=n()) %>% 
  ungroup() %>% 
  spread(character, n, fill=0) %>% 
  mutate(pos1=case_when(`@` > ß ~ '@', 
                        `@` == ß ~ 'N', 
                        `@` < ß ~ 'ß')) %>%
  `$`(pos1)
}

Next we remove duplicate lines from the data frame, then create a new column, pasting together the output of our function called on the first and second character:

df %>% group_by(X2, X3) %>%  
  slice(1) %>%
  ungroup() %>%
  mutate(counts=paste0(summarise_characters(df, 1), summarise_characters(df, 2))) %>% 
  select(X1, X2, X3, X4, X11, X12, X13, counts)

# A tibble: 3 x 8
     X1        X2        X3                   X4       X11                  X12       X13 counts
  <int>     <int>     <int>                <chr>     <int>                <chr>     <int>  <chr>
1     1  90042011  90042031 AGTGCCACCAGGTGGGGCCG  90184354 GAGGCCACAAGAGGGCACAA  90184337     @N
2     1  94853396  94853416 GCCGCCACTAGATGGTGCTA  94969274 CAGTCCTCTAGAGGGAGCCC  94969266     ß@
3     2 103584283 103584303 TGGGCCACTAGGAGGCACTG 103841456 CATGCCACAAGAGGGCATCA 103841456     NN
| improve this answer | |
$\endgroup$
  • $\begingroup$ This looks good. One question, what is the position input in the first function? Because I want to do this to all such repeating lines in my dataset $\endgroup$ – rishi Jan 18 '18 at 13:01
  • $\begingroup$ position in the character string in column 18. As I understand it you want to use the same logic on both the first and second position in this field $\endgroup$ – heathobrien Jan 18 '18 at 13:07
  • $\begingroup$ It gives me this error: Error in mutate_impl(.data, dots) : C stack usage 7970884 is too close to the limit $\endgroup$ – rishi Jan 18 '18 at 13:19
  • $\begingroup$ Interesting. What are the names of the columns in your data frame? $\endgroup$ – heathobrien Jan 18 '18 at 13:25
  • $\begingroup$ It's V1,V2,V3. But I replaced the Xs in your code with Vs already. And changed 18 to 20 since it was actually $20 as mentioned in my question. $\endgroup$ – rishi Jan 18 '18 at 13:27

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.