2016-09-12 17 views
1

私はRの新作で、私は最初のバイオインフォマティクスの割り当てを行っています。だから馬鹿げた間違いがあれば、私を許してください。私はデータフレームをループするのが難しいです。私がしたいことはファイルからdnaの配列を読み、それをアミノ酸コドンに変換することです。この問題は、私が3核塩基をコドンに変換しようとしているときに起こります。条件が満たされず常に偽である場合

私は、bases.csvファイルが3核酸塩基のコドンで次のようになっています。内容は以下のとおりです。

a, b, c, amino 
A, G, G, R 
A, G, A, R 
A, G, C, S 
A, G, T, S 
A, A, G, K 
A, A, A, K 
A, A, C, N 
A, A, T, N 
A, C, G, T 
A, C, A, T 
A, C, C, T 
A, C, T, T 
A, T, G, M 
A, T, A, I 
A, T, C, I 
A, T, T, I 
C, G, G, R 
C, G, A, R 
C, G, C, R 
C, G, T, R 
C, A, G, Q 
C, A, A, Q 
C, A, C, H 
C, A, T, H 
C, C, G, P 
C, C, A, P 
C, C, C, P 
C, C, T, P 
C, T, G, L 
C, T, A, L 
C, T, C, L 
C, T, T, L 
T, G, G, W 
T, G, C, C 
T, G, T, C 
T, A, C, Y 
T, A, T, Y 
T, C, G, S 
T, C, A, S 
T, C, C, S 
T, C, T, S 
T, T, G, L 
T, T, A, L 
T, T, C, F 
T, T, T, F 
G, G, G, G 
G, G, A, G 
G, G, C, G 
G, G, T, G 
G, A, G, E 
G, A, A, E 
G, A, C, D 
G, A, T, D 
G, C, G, A 
G, C, A, A 
G, C, C, A 
G, C, T, A 
G, T, G, V 
G, T, A, V 
G, T, C, V 
G, T, T, V 

EDIT: 完全なソースコード:

lines <- c(readLines("demo.txt")) 
bases <- read.csv(file="bases.csv", header=TRUE, sep=",") 


gene_start <- FALSE 
gene <- "" 
amino <- "" 

convert_to_amino <- function(a, b, c) { 
    index <- 1 
    while (index <= nrow(bases)) { 
     if ((bases[index, 'a'] == a) && (bases[index, 'b'] == b) && (bases[index, 'c'] == c)) { 
      return (bases[index, 'amino']) 
     } 

     index <- index + 1 

    } 


} 



lines <- strsplit(lines, "")[[1]] 
i <- 0 
while ((i + 2 < length(lines))) { 

    if (gene_start == FALSE) { 

     if (lines[i] == 'A' && lines[i + 1] == 'T' && lines[i + 2] == 'G') { 
      gene_start <- TRUE 
      print ("gene found") 
      amino <- convert_to_amino(lines[i], lines[i+1], lines[i+2]) 
      gene <- paste(gene, amino, sep="") 
      i <- i + 3 
     } 

     i <- i + 1 

    } 

    else if (gene_start == TRUE) { 

     if ((lines[i] == 'T' && lines[i + 1] == 'A' && lines[i + 2] == 'A') || (lines[i] == 'T' && lines[i + 1] == 'A' && lines[i + 2] == 'G') || (lines[i] == 'T' && lines[i + 1] == 'G' && lines[i + 2] == 'A')) { 

      gene_start <- FALSE 
      print (gene) 
      gene <- "" 
     } 

     else { 
      amino <- convert_to_amino(lines[i], lines[i+1], lines[i+2]) 
      gene <- paste(gene, amino, sep="") 
     } 

     i <- i + 3 
    } 

    else 
     i <- i + 1 
} 

私はここに達成したいことはabc拠点の組み合わせは、データフレームに存在しているかどうかを確認しています。それがアミノ酸コードであればそれを割り当てます。 しかし、出力によれば、ifの条件は決して満たされません。ここで何が間違っているのかの指針は参考になります。

+0

どのようにデータフレームにCSVファイルを変換していますか? –

+0

'bases < - read.csv(file =" bases.csv "、header = TRUE、sep ="、 ")' – MrPyCharm

+0

データフレームと関数をローカルにロードしましたが、コードは関数呼び出し ' convert_to_amino( 'A'、 'G'、 'G') '...あなたはどのように関数を使用していますか? –

答えて

2

私はうまくいく方法を見つけましたが、それは不器用です。
最初にbasesを4列のデータフレームにする代わりに、2列のデータフレームで作業しました。

ref=bases 
ref$amino=as.character(ref$amino) 
ref$codon=paste0(ref$a,ref$b,ref$c) 
ref$codon=gsub(" ","",ref$codon) 
ref=ref[,c("amino","codon")] 

は、それは次のようになります。

amino codon 
1  R AGG 
2  R AGA 
3  S AGC 
4  S AGT 
5  K AAG 

demoと:CATGTTTCCACTTACAGATCCTTCAAAAAGAGTGTTTCAAAACTGCTCTATGA(あなたのサンプルのみ)

demo="CATGTTTCCACTTACAGATCCTTCAAAAAGAGTGTTTCAAAACTGCTCTATGA" 
demo<- strsplit(demo, "(?<=.{3})", perl = TRUE)[[1]] 

これは(3つの文字、コドンのビットに変身最後のビットはランダムにサンプルの長さを選んだので2だけです)

> demo 
[1] "CAT" "GTT" "TCC" "ACT" "TAC" "AGA" "TCC" "TTC" "AAA" "AAG" "AGT" "GTT" "TCA" "AAA" "CTG" "CTC" "TAT" "GA" 

Iは、次いで、各々がrefから参照してコドン関連付ける:

sapply(demo,function(x)ref$amino[x==ref$codon]) 

(にのみ試料)が得られる:

$CAT 
[1] " H" 

$GTT 
[1] " V" 

$TCC 
[1] " S" 

$ACT 
[1] " T" 

これにはフォーマットはさらに再加工を必要とするかもしれないリストです。結果は参考文献と一致します。

1

私たちは、その後、塩基とマージ3列のデータフレームにシーケンスを変換することができます:

# dummy input 
x <- "CATGTTTCCACTTACAGATCCTTCAAAAAGAGTGTTTCAAAACTGCTCTATGAAAAGGAATGTTCAACTCTGTGAGTTAAATAAAAGCAT" 
nchar(x) 
# [1] 90 

# convert input to a dataframe with 3 columns matching bases columns 
xdf <- data.frame(t(matrix(unlist(strsplit(x, "")), nrow = 3))) 
colnames(xdf) <- colnames(bases)[1:3] 
xdf$ix <- seq(nrow(xdf)) 
head(xdf) 
# a b c ix 
# 1 C A T 1 
# 2 G T T 2 
# 3 T C C 3 
# 4 A C T 4 
# 5 T A C 5 
# 6 A G A 6 

# merge to get amino column 
xdfAmino <- merge(xdf, bases, all.x = TRUE) 

# mark non-matches with "_" 
xdfAmino$amino[is.na(xdfAmino$amino)] <- "_" 

# result 
paste(xdfAmino$amino[order(xdfAmino$ix)], collapse = "") 
#[1] "HVSTYRSFKKSVSKLLYEKECSTL_VK_KH" 
関連する問題