量化跨膜序列中密码子的频率 - 应用功能？
Quantifying frequency of codons in a transmembrane sequence - apply function?
我正在尝试查看某些蛋白质跨膜结构域内的密码子使用情况。
为此，我有 TM 域的序列，我想在这些序列中搜索某些密码子出现的频率（频率）。
理想情况下，我想将新列添加到现有数据框中，其中包含每个基因每个密码子的计数。比如这个假设数据：
Gene ID
TM_domain_Seq
AAA
CAC
GGA
ENSG00000003989
TGGAGCCTCGCTC
0
0
1
ENSG00000003989
TGGAGCCTCGCTC
0
0
1
ENSG00000003989
TGGAGCCTCGCTC
0
0
1
ENSG00000003989
TGGAGCCTCGCTC
0
0
1
ENSG00000003989
TGGAGCCTCGCTC
0
0
1
我尝试了以下方法 - 创建一个函数来计算特定密码子出现的频率，并将其应用于每个 TM 序列。 我遇到的问题是如何将每个密码子的新列添加到我的数据框中，以及如何将密码子频率放入其中。
我试过 for 循环，但它们花费的时间太长
amino_search <- function(seq) {
count <- str_count(seq, pattern = codons)
return(count)
}
codon_search <- function(TMseq) {
High_cor$Newcol <- unlist(lapply(TMseq, amino_search))
}
如有任何帮助，我们将不胜感激。谢谢！
创建可能组合的向量，然后使用 str_count:
comb <- expand.grid(replicate(3, c("A", "T", "G", "C"), simplify = FALSE)) |>
apply(MARGIN = 1, FUN = paste, collapse = "")
#apply(X = _, 1, FUN = paste, collapse = "") #with the new placeholder
df[, comb] <- t(sapply(df$TM_domain_Seq, stringr::str_count, comb))
如果您只需要 in-frame 个密码子，一种方法是每三个字符添加一个 space：
gsub('(.{3})', '\1 ', df$TM_domain_Seq[1])
#[1] "TGG AGC CTC GCT C"
df[, comb] <- t(sapply(gsub('(.{3})', '\1 ', df$TM_domain_Seq), stringr::str_count, comb))
输出
# A tibble: 5 × 66
Gene_ID TM_domain_Seq AAA CAC GGA TAA GAA CAA ATA TTA GTA CTA AGA TGA CGA ACA TCA
<chr> <chr> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int>
1 ENSG00… TGGAGCCTCGCTC 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0
2 ENSG00… TGGAGCCTCGCTC 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0
3 ENSG00… TGGAGCCTCGCTC 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0
4 ENSG00… TGGAGCCTCGCTC 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0
5 ENSG00… TGGAGCCTCGCTC 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0
# … with 49 more variables: GCA <int>, CCA <int>, AAT <int>, TAT <int>, GAT <int>, CAT <int>, ATT <int>,
# TTT <int>, GTT <int>, CTT <int>, AGT <int>, TGT <int>, GGT <int>, CGT <int>, ACT <int>, TCT <int>,
# GCT <int>, CCT <int>, AAG <int>, TAG <int>, GAG <int>, CAG <int>, ATG <int>, TTG <int>, GTG <int>,
# CTG <int>, AGG <int>, TGG <int>, GGG <int>, CGG <int>, ACG <int>, TCG <int>, GCG <int>, CCG <int>,
# AAC <int>, TAC <int>, GAC <int>, ATC <int>, TTC <int>, GTC <int>, CTC <int>, AGC <int>, TGC <int>,
# GGC <int>, CGC <int>, ACC <int>, TCC <int>, GCC <int>, CCC <int>
将问题拆分为 sub-problems，分别解决它们，然后组合解决方案。
第一个子问题是：如何获得给定 (in-frame) 序列的密码子频率？答案是要么使用 pre-made 解决方案（例如 Bioconductor 的 Biostrings:: trinucleotideFrequency(…, steps = 3L)），要么像下面这样快速而肮脏的东西：
codon_frequencies = function (seq) {
# Take care of incomplete codon at end.
len = nchar(seq) - (nchar(seq) %% 3L)
start = seq(1L, len, by = 3L)
substring(seq, start, start + 2L) |> table()
}
试一试：
codon_frequencies('TGGAGCCTCGCTC')
#
# AGC CTC GCT TGG
# 1 1 1 1
… 顺便说一句，你的序列有碎片密码子是故意的吗？如果是这样，您确定它们总是以完整的密码子开始吗？
好的。下一步是为 table 中的每个基因 ID 调用此函数并收集结果。在这一点上，计数 table 可以转换为整洁的数据框这一事实对我们有所帮助：
data.frame(codon_frequencies('TGGAGCCTCGCTC'))
# Var1 Freq
# 1 AGC 1
# 2 CTC 1
# 3 GCT 1
# 4 TGG 1
为了我们的目的，这是一种方便的格式，因为它使 table 操作更容易（尤其是在使用整洁的数据格式时，我在下面使用包 'dplyr'、' tidyr' 和 'purrr'):
df |>
group_by(`Gene ID`) |>
summarize(map_dfr(TM_domain_Seq, ~ data.frame(codon_frequencies(.x))))
# # A tibble: 20 × 3
# # Groups: Gene ID [1]
# `Gene ID` Var1 Freq
# <chr> <fct> <int>
# 1 ENSG00000003989 AGC 1
# 2 ENSG00000003989 CTC 1
# 3 ENSG00000003989 GCT 1
# 4 ENSG00000003989 TGG 1
# …
至此我们可能可以收工了：这是一种方便使用的格式。但是，如果您愿意，也可以将数据转换为宽格式：
… |>
pivot_wider(
id_cols = `Gene ID`,
names_from = Var1,
values_from = Freq,
values_fill = 0L # Otherwise missing codons will be `NA`
)
# # A tibble: 5 × 11
# # Groups: Gene ID [5]
# `Gene ID` AGC CTC GCA TGG TGA GCG AGT GAT TAC GCT
# <chr> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int>
# 1 ENSG00000003981 1 1 1 1 0 0 0 0 0 0
# 2 ENSG00000003982 1 1 0 1 1 0 0 0 0 0
# 3 ENSG00000003983 1 1 0 1 0 1 0 0 0 0
# 4 ENSG00000003984 0 0 0 0 0 0 1 1 1 0
# 5 ENSG00000003989 1 1 0 1 0 0 0 0 0 1
（这是使用了一些不同的玩具数据。）
最后，如果您想要 所有 密码子的列，即使是那些不存在于您的数据中的密码子，您可以对 codon_frequencies 函数进行小的修改：
all_codons = c('A', 'C', 'G', 'T') %>% expand.grid(., ., .) |> apply(1L, paste, collapse = '')
codon_frequencies = function (seq, all = FALSE) {
# Take care of incomplete codon at end.
len = nchar(seq) - (nchar(seq) %% 3L)
start = seq(1L, len, by = 3L)
codons = substring(seq, start, start + 2L)
table(if (all) factor(codons, levels = all_codons) else codons)
}
然后在上面的代码中调用为codon_frequencies(.x, all = TRUE)。 pivot_wider 不再需要 values_fill = 0L 参数。
综合起来：
df |>
group_by(`Gene ID`) |>
summarize(
map_dfr(TM_domain_Seq, ~ data.frame(codon_frequencies(.x, all = TRUE))),
.groups = 'drop'
) |>
pivot_wider(
id_cols = `Gene ID`,
names_from = Var1,
values_from = Freq
)
我正在尝试查看某些蛋白质跨膜结构域内的密码子使用情况。
为此，我有 TM 域的序列，我想在这些序列中搜索某些密码子出现的频率（频率）。
理想情况下，我想将新列添加到现有数据框中，其中包含每个基因每个密码子的计数。比如这个假设数据：
| Gene ID | TM_domain_Seq | AAA | CAC | GGA |
|---|---|---|---|---|
| ENSG00000003989 | TGGAGCCTCGCTC | 0 | 0 | 1 |
| ENSG00000003989 | TGGAGCCTCGCTC | 0 | 0 | 1 |
| ENSG00000003989 | TGGAGCCTCGCTC | 0 | 0 | 1 |
| ENSG00000003989 | TGGAGCCTCGCTC | 0 | 0 | 1 |
| ENSG00000003989 | TGGAGCCTCGCTC | 0 | 0 | 1 |
我尝试了以下方法 - 创建一个函数来计算特定密码子出现的频率，并将其应用于每个 TM 序列。 我遇到的问题是如何将每个密码子的新列添加到我的数据框中，以及如何将密码子频率放入其中。
我试过 for 循环，但它们花费的时间太长
amino_search <- function(seq) {
count <- str_count(seq, pattern = codons)
return(count)
}
codon_search <- function(TMseq) {
High_cor$Newcol <- unlist(lapply(TMseq, amino_search))
}
如有任何帮助，我们将不胜感激。谢谢！
创建可能组合的向量，然后使用 str_count:
comb <- expand.grid(replicate(3, c("A", "T", "G", "C"), simplify = FALSE)) |>
apply(MARGIN = 1, FUN = paste, collapse = "")
#apply(X = _, 1, FUN = paste, collapse = "") #with the new placeholder
df[, comb] <- t(sapply(df$TM_domain_Seq, stringr::str_count, comb))
如果您只需要 in-frame 个密码子，一种方法是每三个字符添加一个 space：
gsub('(.{3})', '\1 ', df$TM_domain_Seq[1])
#[1] "TGG AGC CTC GCT C"
df[, comb] <- t(sapply(gsub('(.{3})', '\1 ', df$TM_domain_Seq), stringr::str_count, comb))
输出
# A tibble: 5 × 66
Gene_ID TM_domain_Seq AAA CAC GGA TAA GAA CAA ATA TTA GTA CTA AGA TGA CGA ACA TCA
<chr> <chr> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int>
1 ENSG00… TGGAGCCTCGCTC 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0
2 ENSG00… TGGAGCCTCGCTC 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0
3 ENSG00… TGGAGCCTCGCTC 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0
4 ENSG00… TGGAGCCTCGCTC 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0
5 ENSG00… TGGAGCCTCGCTC 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0
# … with 49 more variables: GCA <int>, CCA <int>, AAT <int>, TAT <int>, GAT <int>, CAT <int>, ATT <int>,
# TTT <int>, GTT <int>, CTT <int>, AGT <int>, TGT <int>, GGT <int>, CGT <int>, ACT <int>, TCT <int>,
# GCT <int>, CCT <int>, AAG <int>, TAG <int>, GAG <int>, CAG <int>, ATG <int>, TTG <int>, GTG <int>,
# CTG <int>, AGG <int>, TGG <int>, GGG <int>, CGG <int>, ACG <int>, TCG <int>, GCG <int>, CCG <int>,
# AAC <int>, TAC <int>, GAC <int>, ATC <int>, TTC <int>, GTC <int>, CTC <int>, AGC <int>, TGC <int>,
# GGC <int>, CGC <int>, ACC <int>, TCC <int>, GCC <int>, CCC <int>
将问题拆分为 sub-problems，分别解决它们，然后组合解决方案。
第一个子问题是：如何获得给定 (in-frame) 序列的密码子频率？答案是要么使用 pre-made 解决方案（例如 Bioconductor 的 Biostrings:: trinucleotideFrequency(…, steps = 3L)），要么像下面这样快速而肮脏的东西：
codon_frequencies = function (seq) {
# Take care of incomplete codon at end.
len = nchar(seq) - (nchar(seq) %% 3L)
start = seq(1L, len, by = 3L)
substring(seq, start, start + 2L) |> table()
}
试一试：
codon_frequencies('TGGAGCCTCGCTC')
#
# AGC CTC GCT TGG
# 1 1 1 1
… 顺便说一句，你的序列有碎片密码子是故意的吗？如果是这样，您确定它们总是以完整的密码子开始吗？
好的。下一步是为 table 中的每个基因 ID 调用此函数并收集结果。在这一点上，计数 table 可以转换为整洁的数据框这一事实对我们有所帮助：
data.frame(codon_frequencies('TGGAGCCTCGCTC'))
# Var1 Freq
# 1 AGC 1
# 2 CTC 1
# 3 GCT 1
# 4 TGG 1
为了我们的目的，这是一种方便的格式，因为它使 table 操作更容易（尤其是在使用整洁的数据格式时，我在下面使用包 'dplyr'、' tidyr' 和 'purrr'):
df |>
group_by(`Gene ID`) |>
summarize(map_dfr(TM_domain_Seq, ~ data.frame(codon_frequencies(.x))))
# # A tibble: 20 × 3
# # Groups: Gene ID [1]
# `Gene ID` Var1 Freq
# <chr> <fct> <int>
# 1 ENSG00000003989 AGC 1
# 2 ENSG00000003989 CTC 1
# 3 ENSG00000003989 GCT 1
# 4 ENSG00000003989 TGG 1
# …
至此我们可能可以收工了：这是一种方便使用的格式。但是，如果您愿意，也可以将数据转换为宽格式：
… |>
pivot_wider(
id_cols = `Gene ID`,
names_from = Var1,
values_from = Freq,
values_fill = 0L # Otherwise missing codons will be `NA`
)
# # A tibble: 5 × 11
# # Groups: Gene ID [5]
# `Gene ID` AGC CTC GCA TGG TGA GCG AGT GAT TAC GCT
# <chr> <int> <int> <int> <int> <int> <int> <int> <int> <int> <int>
# 1 ENSG00000003981 1 1 1 1 0 0 0 0 0 0
# 2 ENSG00000003982 1 1 0 1 1 0 0 0 0 0
# 3 ENSG00000003983 1 1 0 1 0 1 0 0 0 0
# 4 ENSG00000003984 0 0 0 0 0 0 1 1 1 0
# 5 ENSG00000003989 1 1 0 1 0 0 0 0 0 1
（这是使用了一些不同的玩具数据。）
最后，如果您想要 所有 密码子的列，即使是那些不存在于您的数据中的密码子，您可以对 codon_frequencies 函数进行小的修改：
all_codons = c('A', 'C', 'G', 'T') %>% expand.grid(., ., .) |> apply(1L, paste, collapse = '')
codon_frequencies = function (seq, all = FALSE) {
# Take care of incomplete codon at end.
len = nchar(seq) - (nchar(seq) %% 3L)
start = seq(1L, len, by = 3L)
codons = substring(seq, start, start + 2L)
table(if (all) factor(codons, levels = all_codons) else codons)
}
然后在上面的代码中调用为codon_frequencies(.x, all = TRUE)。 pivot_wider 不再需要 values_fill = 0L 参数。
综合起来：
df |>
group_by(`Gene ID`) |>
summarize(
map_dfr(TM_domain_Seq, ~ data.frame(codon_frequencies(.x, all = TRUE))),
.groups = 'drop'
) |>
pivot_wider(
id_cols = `Gene ID`,
names_from = Var1,
values_from = Freq
)