CODEWITHSUNDEEP
r
Max Cummins
Max Cummins

Reputation: 13

Challenge in data manipulation

I am trying to process a database of BLAST output to generate a data frame containing values for a given gene and given sample. When a gene is identified within a sample I would like the scaffold on which it was identified to be reported. If a given gene is NOT identified within a given sample I would like the cell to be filled with N/A.

sample_name         scaffold    gene_title  match_(%)
P24_ST48                64      aadA12      94.56
401B_ST5223             381     blaTEM-163  99.65
P32_ST218               91      aadA24      90.41
HOS66_ST73              9       blaACT-5    72.31
HOS16_ST38              70      blaTEM-146  99.42
HOS56_ST131             48      aadA21      91.39
Ecoli_2009_1_ST131      41      sul1        99.88
PH152_ST95              37      dfrA33      83.94
Ecoli_2009_32_STNT      16      aac(3)-Ib   100.00
PH231_ST38              59      mph(D)      89.83
P44_STNT                135     blaTEM-105  99.88
Ecoli_2011_89_ST127     29      blaTEM-158  99.65
405C_ST1178             120     aadA1       99.75
P3_STNT                 15      blaTEM-68   99.19
5A_ST34                 174     blaTEM-127  99.88
P27_ST10                211     aph(3')-Ia  100.00
4D_ST767                393     blaTEM-152  98.95
P10_STNT                23      blaTEM-17   99.07
Ecoli_2014_27_ST131     49      sul2_15     99.88
Ecoli_2013_10_ST73      23      blaTEM-2    99.19

The output table would look something like:

Sample        aadA1    aadA12    aadA24    blaTEM-163 ...
P24_ST48      N/A      64        N/A       N/A
401B_ST5223   N/A      N/A       N/A       381
...

In excel I have concatenated the sample name and gene titles and reported the scaffold number on row where this string is identified using VLOOKUP - I have tried many different ways in R and am going around in circles.

Now trying to process +700 genes and +450 samples, the list of gene-sample combinations is getting somewhat laborious for excel to manage and I must find another solution with my collection of samples growing increasingly large.

Any help would be greatly appreciated.

Cheers,

Max

Upvotes: 0

Views: 57

Answers (2)

akrun
akrun

Reputation: 887088

We can use dcast from data.table

library(data.table)
dcast(setDT(df1), sample_name + match_... ~ gene_title, value.var = 'scaffold')
#      sample_name match_... aac(3)-Ib aadA1 aadA12 ...
#1:         401B_ST5223     99.65        NA    NA   
#2:         405C_ST1178     99.75        NA   120   
#3:            4D_ST767     98.95        NA    NA   
#4:             5A_ST34     99.88        NA    NA

Upvotes: 0

Pierre Lapointe
Pierre Lapointe

Reputation: 16277

Here's how to do that with spread from tidyr

library(tidyr)
df1%>%
  spread(key = gene_title,value = scaffold)

         sample_name match_... aac(3)-Ib aadA1 aadA12 ...
1        401B_ST5223     99.65        NA    NA     NA
2        405C_ST1178     99.75        NA   120     NA
3           4D_ST767     98.95        NA    NA     NA
4            5A_ST34     99.88        NA    NA     NA
5 Ecoli_2009_1_ST131     99.88        NA    NA     NA
...

Data

df1 <- read.table(text="sample_name         scaffold    gene_title  match_(%)
P24_ST48                64      aadA12      94.56
                  401B_ST5223             381     blaTEM-163  99.65
                  P32_ST218               91      aadA24      90.41
                  HOS66_ST73              9       blaACT-5    72.31
                  HOS16_ST38              70      blaTEM-146  99.42
                  HOS56_ST131             48      aadA21      91.39
                  Ecoli_2009_1_ST131      41      sul1        99.88
                  PH152_ST95              37      dfrA33      83.94
                  Ecoli_2009_32_STNT      16      aac(3)-Ib   100.00
                  PH231_ST38              59      mph(D)      89.83
                  P44_STNT                135     blaTEM-105  99.88
                  Ecoli_2011_89_ST127     29      blaTEM-158  99.65
                  405C_ST1178             120     aadA1       99.75
                  P3_STNT                 15      blaTEM-68   99.19
                  5A_ST34                 174     blaTEM-127  99.88
                  P27_ST10                211     aph(3')-Ia  100.00
                  4D_ST767                393     blaTEM-152  98.95
                  P10_STNT                23      blaTEM-17   99.07
                  Ecoli_2014_27_ST131     49      sul2_15     99.88
                  Ecoli_2013_10_ST73      23      blaTEM-2    99.19",
                        header=TRUE,stringsAsFactors=FALSE)

Upvotes: 1

Related Questions