R中的DNA条件频率
DNA conditional frequency in R
我试图找出 R 中 2 个不同 DNA 序列中是否存在任何条件依赖性
这是我的代码,但是我遇到了错误;
Error in `[.data.frame`(data, i) : undefined columns selected
我不确定问题出在哪里,如果我将 data[i-1]==bases[b2]
括起来,我只会得到多个 unexpected}
,这是我唯一能想到的事情。
for (b1 in 1:length(bases))
{
for (b2 in 1:length(bases))
{
count = 1
for (i in 2:length(mydata1))
{
if ((mydata1[i]==bases[b1]) & mydata1[i-1]==bases[b2])
{
count = count+1
}
}
b3 = c(bases[b1], bases[b2], count)
print(b3)
}
}
_我基本上期待某些 DNA 碱基的列表,例如,我认为 DNA 序列好像是以先前的碱基为条件的;
[1] "A" "C" "002"
[1] "A" "C" "005"
[1] "A" "C" "009"
等等,通过清楚地显示 A 的条件,可以告诉我某个碱基是否对后续碱基的身份有任何影响C.
Ok 所以本质上 mydata1
(还有 mydata2
)是 DNA 序列,也就是说 "A", "G", "C" and "T"
的列表,每个都是10,000 个碱基长。
如图所示;
V1
1 T
2 C
3 G
4 G
5 T
6 G
7 G
8 G
9 C
10 A
我的任务是尝试确定该序列是否具有相互依赖的碱基,因此如果 [1] T
影响 [2] C
的存在,等等。一个一个序列是相关的,另一个不是。
如果我没理解错的话,你想计算DNA序列中每对核苷酸i,i+1出现的次数。您可以使用 R 函数 table 实现此目的;下面提供了一个示例。
# input sequence
seq <- "ACGTACTGCACAAACTAC"
# length of input sequence
length_seq <- nchar(seq, type="chars")
# first substring: from 1 to second-last
seq1 <- substr(seq, 1, (length_seq - 1))
# second substring: from 2 to last
seq2 <- substring(seq, 2, length_seq)
# split strings
seq1_split <- strsplit(seq1, "")[[1]]
seq2_split <- strsplit(seq2, "")[[1]]
# initialize vectors
first_nt <- vector(mode="character", length = (length_seq - 1))
second_nt <- vector(mode="character", length = (length_seq -1))
# fill vectors
count = 0
for (b in seq1_split)
{
count = count + 1
first_nt[count] <- b
}
count = 0
for (b in seq2_split)
{
count = count + 1
second_nt[count] <- b
}
# create matrix with character i and i+1 in each row
mat <- matrix(c(first_nt, second_nt), nrow=(length_seq - 1))
# collapse matrix
to_table <- apply(mat, 1, paste, collapse="")
# table
my_table <- table(to_table)
print(my_table)
我试图找出 R 中 2 个不同 DNA 序列中是否存在任何条件依赖性
这是我的代码,但是我遇到了错误;
Error in `[.data.frame`(data, i) : undefined columns selected
我不确定问题出在哪里,如果我将 data[i-1]==bases[b2]
括起来,我只会得到多个 unexpected}
,这是我唯一能想到的事情。
for (b1 in 1:length(bases))
{
for (b2 in 1:length(bases))
{
count = 1
for (i in 2:length(mydata1))
{
if ((mydata1[i]==bases[b1]) & mydata1[i-1]==bases[b2])
{
count = count+1
}
}
b3 = c(bases[b1], bases[b2], count)
print(b3)
}
}
_我基本上期待某些 DNA 碱基的列表,例如,我认为 DNA 序列好像是以先前的碱基为条件的;
[1] "A" "C" "002"
[1] "A" "C" "005"
[1] "A" "C" "009"
等等,通过清楚地显示 A 的条件,可以告诉我某个碱基是否对后续碱基的身份有任何影响C.
Ok 所以本质上 mydata1
(还有 mydata2
)是 DNA 序列,也就是说 "A", "G", "C" and "T"
的列表,每个都是10,000 个碱基长。
如图所示;
V1
1 T
2 C
3 G
4 G
5 T
6 G
7 G
8 G
9 C
10 A
我的任务是尝试确定该序列是否具有相互依赖的碱基,因此如果 [1] T
影响 [2] C
的存在,等等。一个一个序列是相关的,另一个不是。
如果我没理解错的话,你想计算DNA序列中每对核苷酸i,i+1出现的次数。您可以使用 R 函数 table 实现此目的;下面提供了一个示例。
# input sequence
seq <- "ACGTACTGCACAAACTAC"
# length of input sequence
length_seq <- nchar(seq, type="chars")
# first substring: from 1 to second-last
seq1 <- substr(seq, 1, (length_seq - 1))
# second substring: from 2 to last
seq2 <- substring(seq, 2, length_seq)
# split strings
seq1_split <- strsplit(seq1, "")[[1]]
seq2_split <- strsplit(seq2, "")[[1]]
# initialize vectors
first_nt <- vector(mode="character", length = (length_seq - 1))
second_nt <- vector(mode="character", length = (length_seq -1))
# fill vectors
count = 0
for (b in seq1_split)
{
count = count + 1
first_nt[count] <- b
}
count = 0
for (b in seq2_split)
{
count = count + 1
second_nt[count] <- b
}
# create matrix with character i and i+1 in each row
mat <- matrix(c(first_nt, second_nt), nrow=(length_seq - 1))
# collapse matrix
to_table <- apply(mat, 1, paste, collapse="")
# table
my_table <- table(to_table)
print(my_table)