问题描述
我下面有一个data.table,它包含DNA字符串组(ref_seq
列)和一系列突变(alt
列中的字符)及其沿字符串的相对位置(snp_position
列)。
由于 n 个突变,因此产生的字符串数为2 ^ n ,我想做的是制作一个通用函数,该函数可以生成所有在alt
字符串中的snp_position
处替换的ref_seq
字符组合,同时(最好)保留此data.frame结构
> dt
seqnames hap_start hap_end snp alt snp_position numb_snps_per_hap
1: chr1 19600274 19600443 19600324 G 50 4
2: chr1 19600274 19600443 19600378 C 104 4
3: chr1 19600274 19600443 19600389 A 115 4
4: chr1 19600274 19600443 19600396 C 122 4
5: chr10 5482730 5482899 5482790 C 60 4
6: chr10 5482730 5482899 5482830 A 100 4
7: chr10 5482730 5482899 5482839 A 109 4
8: chr10 5482730 5482899 5482843 A 113 4
ref_seq
1: CCAGGGAGAGCGGGGCAGCGCACGGCTGGTGACAGGCTCAGCTCTGCACCGTGCAGAGGGAGGATCAGGGGCCACTGTTACCTCTGCAGTCGCTGCTGCCCATTCGGAGCCCCAGGCCGCCAAGGATGGTCTCGGACTGGCCGTCGCTGTACAGGAAGGCCGTGTCTATC
2: CCAGGGAGAGCGGGGCAGCGCACGGCTGGTGACAGGCTCAGCTCTGCACCGTGCAGAGGGAGGATCAGGGGCCACTGTTACCTCTGCAGTCGCTGCTGCCCATTCGGAGCCCCAGGCCGCCAAGGATGGTCTCGGACTGGCCGTCGCTGTACAGGAAGGCCGTGTCTATC
3: CCAGGGAGAGCGGGGCAGCGCACGGCTGGTGACAGGCTCAGCTCTGCACCGTGCAGAGGGAGGATCAGGGGCCACTGTTACCTCTGCAGTCGCTGCTGCCCATTCGGAGCCCCAGGCCGCCAAGGATGGTCTCGGACTGGCCGTCGCTGTACAGGAAGGCCGTGTCTATC
4: CCAGGGAGAGCGGGGCAGCGCACGGCTGGTGACAGGCTCAGCTCTGCACCGTGCAGAGGGAGGATCAGGGGCCACTGTTACCTCTGCAGTCGCTGCTGCCCATTCGGAGCCCCAGGCCGCCAAGGATGGTCTCGGACTGGCCGTCGCTGTACAGGAAGGCCGTGTCTATC
5: TGCCCCTCTCTTGTTCATTTGCCACCACTACTATGTAAGGCATTGCCTCACTACTCTGTGTCTCATCCCTATTCATACCTGTAGTCCTTATCACCATTTTACTAGTGAGTGAGTCTGTTTCAGAGTCCCATCCTGTGTTTCTGATCTCTAGGACCAATCTCAGTGTCCAT
6: TGCCCCTCTCTTGTTCATTTGCCACCACTACTATGTAAGGCATTGCCTCACTACTCTGTGTCTCATCCCTATTCATACCTGTAGTCCTTATCACCATTTTACTAGTGAGTGAGTCTGTTTCAGAGTCCCATCCTGTGTTTCTGATCTCTAGGACCAATCTCAGTGTCCAT
7: TGCCCCTCTCTTGTTCATTTGCCACCACTACTATGTAAGGCATTGCCTCACTACTCTGTGTCTCATCCCTATTCATACCTGTAGTCCTTATCACCATTTTACTAGTGAGTGAGTCTGTTTCAGAGTCCCATCCTGTGTTTCTGATCTCTAGGACCAATCTCAGTGTCCAT
8: TGCCCCTCTCTTGTTCATTTGCCACCACTACTATGTAAGGCATTGCCTCACTACTCTGTGTCTCATCCCTATTCATACCTGTAGTCCTTATCACCATTTTACTAGTGAGTGAGTCTGTTTCAGAGTCCCATCCTGTGTTTCTGATCTCTAGGACCAATCTCAGTGTCCAT
> dput(dt)
structure(list(seqnames = c("chr1","chr1","chr10","chr10"),hap_start = c(19600274,19600274,5482730,5482730),hap_end = c(19600443,19600443,5482899,5482899
),snp = c(19600324L,19600378L,19600389L,19600396L,5482790L,5482830L,5482839L,5482843L),alt = c("G","C","A","A"),snp_position = c(50,104,115,122,60,100,109,113),numb_snps_per_hap = c(4L,4L,4L),ref_seq = c("CCAGGGAGAGCGGGGCAGCGCACGGCTGGTGACAGGCTCAGCTCTGCACCGTGCAGAGGGAGGATCAGGGGCCACTGTTACCTCTGCAGTCGCTGCTGCCCATTCGGAGCCCCAGGCCGCCAAGGATGGTCTCGGACTGGCCGTCGCTGTACAGGAAGGCCGTGTCTATC","CCAGGGAGAGCGGGGCAGCGCACGGCTGGTGACAGGCTCAGCTCTGCACCGTGCAGAGGGAGGATCAGGGGCCACTGTTACCTCTGCAGTCGCTGCTGCCCATTCGGAGCCCCAGGCCGCCAAGGATGGTCTCGGACTGGCCGTCGCTGTACAGGAAGGCCGTGTCTATC","TGCCCCTCTCTTGTTCATTTGCCACCACTACTATGTAAGGCATTGCCTCACTACTCTGTGTCTCATCCCTATTCATACCTGTAGTCCTTATCACCATTTTACTAGTGAGTGAGTCTGTTTCAGAGTCCCATCCTGTGTTTCTGATCTCTAGGACCAATCTCAGTGTCCAT","TGCCCCTCTCTTGTTCATTTGCCACCACTACTATGTAAGGCATTGCCTCACTACTCTGTGTCTCATCCCTATTCATACCTGTAGTCCTTATCACCATTTTACTAGTGAGTGAGTCTGTTTCAGAGTCCCATCCTGTGTTTCTGATCTCTAGGACCAATCTCAGTGTCCAT"
)),row.names = c(NA,-8L),class = c("data.table","data.frame"
),sorted = "seqnames",.internal.selfref = <pointer: 0x1d50f80>)
所需的输出
> dt_final
seqnames hap_start hap_end
1: chr1 19600274 19600443
2: chr1 19600274 19600443
3: chr1 19600274 19600443
4: chr1 19600274 19600443
5: chr1 19600274 19600443
6: chr1 19600274 19600443
7: chr1 19600274 19600443
8: chr1 19600274 19600443
9: chr1 19600274 19600443
10: chr1 19600274 19600443
11: chr1 19600274 19600443
12: chr1 19600274 19600443
13: chr1 19600274 19600443
14: chr1 19600274 19600443
15: chr1 19600274 19600443
16: chr1 19600274 19600443
17: chr10 5482730 5482899
18: chr10 5482730 5482899
19: chr10 5482730 5482899
20: chr10 5482730 5482899
21: chr10 5482730 5482899
22: chr10 5482730 5482899
23: chr10 5482730 5482899
24: chr10 5482730 5482899
25: chr10 5482730 5482899
26: chr10 5482730 5482899
27: chr10 5482730 5482899
28: chr10 5482730 5482899
29: chr10 5482730 5482899
30: chr10 5482730 5482899
31: chr10 5482730 5482899
32: chr10 5482730 5482899
seqnames hap_start hap_end
sequence
1: CCAGGGAGAGCGGGGCAGCGCACGGCTGGTGACAGGCTCAGCTCTGCACCGTGCAGAGGGAGGATCAGGGGCCACTGTTACCTCTGCAGTCGCTGCTGCCCATTCGGAGCCCCAAGCCGCCCAGGATGGTCTCGGACTGGCCGTCGCTGTACAGGAAGGCCGTGTCTATC
2: CCAGGGAGAGCGGGGCAGCGCACGGCTGGTGACAGGCTCAGCTCTGCACGGTGCAGAGGGAGGATCAGGGGCCACTGTTACCTCTGCAGTCGCTGCTGCCCATCCGGAGCCCCAGGCCGCCAAGGATGGTCTCGGACTGGCCGTCGCTGTACAGGAAGGCCGTGTCTATC
3: CCAGGGAGAGCGGGGCAGCGCACGGCTGGTGACAGGCTCAGCTCTGCACGGTGCAGAGGGAGGATCAGGGGCCACTGTTACCTCTGCAGTCGCTGCTGCCCATTCGGAGCCCCAAGCCGCCAAGGATGGTCTCGGACTGGCCGTCGCTGTACAGGAAGGCCGTGTCTATC
4: CCAGGGAGAGCGGGGCAGCGCACGGCTGGTGACAGGCTCAGCTCTGCACCGTGCAGAGGGAGGATCAGGGGCCACTGTTACCTCTGCAGTCGCTGCTGCCCATCCGGAGCCCCAGGCCGCCCAGGATGGTCTCGGACTGGCCGTCGCTGTACAGGAAGGCCGTGTCTATC
5: CCAGGGAGAGCGGGGCAGCGCACGGCTGGTGACAGGCTCAGCTCTGCACGGTGCAGAGGGAGGATCAGGGGCCACTGTTACCTCTGCAGTCGCTGCTGCCCATCCGGAGCCCCAAGCCGCCCAGGATGGTCTCGGACTGGCCGTCGCTGTACAGGAAGGCCGTGTCTATC
6: CCAGGGAGAGCGGGGCAGCGCACGGCTGGTGACAGGCTCAGCTCTGCACGGTGCAGAGGGAGGATCAGGGGCCACTGTTACCTCTGCAGTCGCTGCTGCCCATCCGGAGCCCCAAGCCGCCCAGGATGGTCTCGGACTGGCCGTCGCTGTACAGGAAGGCCGTGTCTATC
7: CCAGGGAGAGCGGGGCAGCGCACGGCTGGTGACAGGCTCAGCTCTGCACGGTGCAGAGGGAGGATCAGGGGCCACTGTTACCTCTGCAGTCGCTGCTGCCCATCCGGAGCCCCAAGCCGCCCAGGATGGTCTCGGACTGGCCGTCGCTGTACAGGAAGGCCGTGTCTATC
8: CCAGGGAGAGCGGGGCAGCGCACGGCTGGTGACAGGCTCAGCTCTGCACGGTGCAGAGGGAGGATCAGGGGCCACTGTTACCTCTGCAGTCGCTGCTGCCCATCCGGAGCCCCAAGCCGCCCAGGATGGTCTCGGACTGGCCGTCGCTGTACAGGAAGGCCGTGTCTATC
9: CCAGGGAGAGCGGGGCAGCGCACGGCTGGTGACAGGCTCAGCTCTGCACGGTGCAGAGGGAGGATCAGGGGCCACTGTTACCTCTGCAGTCGCTGCTGCCCATTCGGAGCCCCAGGCCGCCAAGGATGGTCTCGGACTGGCCGTCGCTGTACAGGAAGGCCGTGTCTATC
10: CCAGGGAGAGCGGGGCAGCGCACGGCTGGTGACAGGCTCAGCTCTGCACCGTGCAGAGGGAGGATCAGGGGCCACTGTTACCTCTGCAGTCGCTGCTGCCCATCCGGAGCCCCAGGCCGCCAAGGATGGTCTCGGACTGGCCGTCGCTGTACAGGAAGGCCGTGTCTATC
11: CCAGGGAGAGCGGGGCAGCGCACGGCTGGTGACAGGCTCAGCTCTGCACCGTGCAGAGGGAGGATCAGGGGCCACTGTTACCTCTGCAGTCGCTGCTGCCCATTCGGAGCCCCAAGCCGCCAAGGATGGTCTCGGACTGGCCGTCGCTGTACAGGAAGGCCGTGTCTATC
12: CCAGGGAGAGCGGGGCAGCGCACGGCTGGTGACAGGCTCAGCTCTGCACCGTGCAGAGGGAGGATCAGGGGCCACTGTTACCTCTGCAGTCGCTGCTGCCCATTCGGAGCCCCAGGCCGCCCAGGATGGTCTCGGACTGGCCGTCGCTGTACAGGAAGGCCGTGTCTATC
13: CCAGGGAGAGCGGGGCAGCGCACGGCTGGTGACAGGCTCAGCTCTGCACCGTGCAGAGGGAGGATCAGGGGCCACTGTTACCTCTGCAGTCGCTGCTGCCCATCCGGAGCCCCAAGCCGCCCAGGATGGTCTCGGACTGGCCGTCGCTGTACAGGAAGGCCGTGTCTATC
14: CCAGGGAGAGCGGGGCAGCGCACGGCTGGTGACAGGCTCAGCTCTGCACGGTGCAGAGGGAGGATCAGGGGCCACTGTTACCTCTGCAGTCGCTGCTGCCCATTCGGAGCCCCAAGCCGCCCAGGATGGTCTCGGACTGGCCGTCGCTGTACAGGAAGGCCGTGTCTATC
15: CCAGGGAGAGCGGGGCAGCGCACGGCTGGTGACAGGCTCAGCTCTGCACGGTGCAGAGGGAGGATCAGGGGCCACTGTTACCTCTGCAGTCGCTGCTGCCCATCCGGAGCCCCAGGCCGCCCAGGATGGTCTCGGACTGGCCGTCGCTGTACAGGAAGGCCGTGTCTATC
16: CCAGGGAGAGCGGGGCAGCGCACGGCTGGTGACAGGCTCAGCTCTGCACGGTGCAGAGGGAGGATCAGGGGCCACTGTTACCTCTGCAGTCGCTGCTGCCCATCCGGAGCCCCAAGCCGCCAAGGATGGTCTCGGACTGGCCGTCGCTGTACAGGAAGGCCGTGTCTATC
17: TGCCCCTCTCTTGTTCATTTGCCACCACTACTATGTAAGGCATTGCCTCACTACTCTGTGTCTCATCCCTATTCATACCTGTAGTCCTTATCACCATTTTACTAGTGAATGAATCTGTTTCAGAGTCCCATCCTGTGTTTCTGATCTCTAGGACCAATCTCAGTGTCCAT
18: TGCCCCTCTCTTGTTCATTTGCCACCACTACTATGTAAGGCATTGCCTCACTACTCTGTCTCTCATCCCTATTCATACCTGTAGTCCTTATCACCATTTAACTAGTGAGTGAGTCTGTTTCAGAGTCCCATCCTGTGTTTCTGATCTCTAGGACCAATCTCAGTGTCCAT
19: TGCCCCTCTCTTGTTCATTTGCCACCACTACTATGTAAGGCATTGCCTCACTACTCTGTCTCTCATCCCTATTCATACCTGTAGTCCTTATCACCATTTTACTAGTGAATGAGTCTGTTTCAGAGTCCCATCCTGTGTTTCTGATCTCTAGGACCAATCTCAGTGTCCAT
20: TGCCCCTCTCTTGTTCATTTGCCACCACTACTATGTAAGGCATTGCCTCACTACTCTGTGTCTCATCCCTATTCATACCTGTAGTCCTTATCACCATTTAACTAGTGAGTGAATCTGTTTCAGAGTCCCATCCTGTGTTTCTGATCTCTAGGACCAATCTCAGTGTCCAT
21: TGCCCCTCTCTTGTTCATTTGCCACCACTACTATGTAAGGCATTGCCTCACTACTCTGTCTCTCATCCCTATTCATACCTGTAGTCCTTATCACCATTTAACTAGTGAATGAATCTGTTTCAGAGTCCCATCCTGTGTTTCTGATCTCTAGGACCAATCTCAGTGTCCAT
22: TGCCCCTCTCTTGTTCATTTGCCACCACTACTATGTAAGGCATTGCCTCACTACTCTGTCTCTCATCCCTATTCATACCTGTAGTCCTTATCACCATTTAACTAGTGAATGAATCTGTTTCAGAGTCCCATCCTGTGTTTCTGATCTCTAGGACCAATCTCAGTGTCCAT
23: TGCCCCTCTCTTGTTCATTTGCCACCACTACTATGTAAGGCATTGCCTCACTACTCTGTCTCTCATCCCTATTCATACCTGTAGTCCTTATCACCATTTAACTAGTGAATGAATCTGTTTCAGAGTCCCATCCTGTGTTTCTGATCTCTAGGACCAATCTCAGTGTCCAT
24: TGCCCCTCTCTTGTTCATTTGCCACCACTACTATGTAAGGCATTGCCTCACTACTCTGTCTCTCATCCCTATTCATACCTGTAGTCCTTATCACCATTTAACTAGTGAATGAATCTGTTTCAGAGTCCCATCCTGTGTTTCTGATCTCTAGGACCAATCTCAGTGTCCAT
25: TGCCCCTCTCTTGTTCATTTGCCACCACTACTATGTAAGGCATTGCCTCACTACTCTGTCTCTCATCCCTATTCATACCTGTAGTCCTTATCACCATTTTACTAGTGAGTGAGTCTGTTTCAGAGTCCCATCCTGTGTTTCTGATCTCTAGGACCAATCTCAGTGTCCAT
26: TGCCCCTCTCTTGTTCATTTGCCACCACTACTATGTAAGGCATTGCCTCACTACTCTGTGTCTCATCCCTATTCATACCTGTAGTCCTTATCACCATTTAACTAGTGAGTGAGTCTGTTTCAGAGTCCCATCCTGTGTTTCTGATCTCTAGGACCAATCTCAGTGTCCAT
27: TGCCCCTCTCTTGTTCATTTGCCACCACTACTATGTAAGGCATTGCCTCACTACTCTGTGTCTCATCCCTATTCATACCTGTAGTCCTTATCACCATTTTACTAGTGAATGAGTCTGTTTCAGAGTCCCATCCTGTGTTTCTGATCTCTAGGACCAATCTCAGTGTCCAT
28: TGCCCCTCTCTTGTTCATTTGCCACCACTACTATGTAAGGCATTGCCTCACTACTCTGTGTCTCATCCCTATTCATACCTGTAGTCCTTATCACCATTTTACTAGTGAGTGAATCTGTTTCAGAGTCCCATCCTGTGTTTCTGATCTCTAGGACCAATCTCAGTGTCCAT
29: TGCCCCTCTCTTGTTCATTTGCCACCACTACTATGTAAGGCATTGCCTCACTACTCTGTGTCTCATCCCTATTCATACCTGTAGTCCTTATCACCATTTAACTAGTGAATGAATCTGTTTCAGAGTCCCATCCTGTGTTTCTGATCTCTAGGACCAATCTCAGTGTCCAT
30: TGCCCCTCTCTTGTTCATTTGCCACCACTACTATGTAAGGCATTGCCTCACTACTCTGTCTCTCATCCCTATTCATACCTGTAGTCCTTATCACCATTTTACTAGTGAATGAATCTGTTTCAGAGTCCCATCCTGTGTTTCTGATCTCTAGGACCAATCTCAGTGTCCAT
31: TGCCCCTCTCTTGTTCATTTGCCACCACTACTATGTAAGGCATTGCCTCACTACTCTGTCTCTCATCCCTATTCATACCTGTAGTCCTTATCACCATTTAACTAGTGAGTGAATCTGTTTCAGAGTCCCATCCTGTGTTTCTGATCTCTAGGACCAATCTCAGTGTCCAT
32: TGCCCCTCTCTTGTTCATTTGCCACCACTACTATGTAAGGCATTGCCTCACTACTCTGTCTCTCATCCCTATTCATACCTGTAGTCCTTATCACCATTTAACTAGTGAATGAGTCTGTTTCAGAGTCCCATCCTGTGTTTCTGATCTCTAGGACCAATCTCAGTGTCCAT
sequence
您知道如何在R中执行此操作吗?
谢谢
解决方法
我只使用了前两行和第二列。我还减少了refseq的大小,并使用snp_position = 2以获得更好的可视化效果。但是,一旦名称相同,它将与您的完整数据相同。
输入:
df <- data.frame(alt=c("G","A"),snp_position=c(2,2),refseq=c("CCAGGGAGAGCGGGGCAGCGC","TCAGGGAGAGCGGGGCAGCGC"))
功能:
seq_mut <- function(refseq,snp_position,Mutated_allele){
combinations <- c("A","T","C","G")
WT_allele <- substr(refseq,snp_position)
combinations <- combinations[combinations!=WT_allele]
refseq2 <- refseq
substr(refseq2,snp_position) <- Mutated_allele
combinations <- combinations[combinations!=Mutated_allele]
seq <- refseq2
for(variant in combinations){
refseq2 <- refseq
substr(refseq2,snp_position) <- variant
seq <- paste0(seq,";",refseq2)
}
return(seq)
}
与dplyr:
library(dplyr)
df <- df %>% rowwise() %>% mutate(ref_seq_mut=seq_mut(refseq,alt))
tidyr:
library(tidyr)
df <- df %>% separate(ref_seq_mut,into=c("V1","V2","V3"))%>%
pivot_longer(cols=refseq:V3,values_to="refseq") %>% select(-name)
输出:
alt snp_position refseq
<chr> <dbl> <chr>
1 G 2 CCAGGGAGAGCGGGGCAGCGC
2 G 2 CGAGGGAGAGCGGGGCAGCGC
3 G 2 CAAGGGAGAGCGGGGCAGCGC
4 G 2 CTAGGGAGAGCGGGGCAGCGC
5 A 2 TCAGGGAGAGCGGGGCAGCGC
6 A 2 TAAGGGAGAGCGGGGCAGCGC
7 A 2 TTAGGGAGAGCGGGGCAGCGC
8 A 2 TGAGGGAGAGCGGGGCAGCGC
,
这里是使用简化示例输入数据的选项。有很多步骤需要介绍太多细节,请尝试一点一点地运行它以了解每个功能。
基本上,我是在 refseq 上拆分,然后将 refseq 拆分为字母,然后用 alt 更新位置组合,然后粘贴回去
def partial_sum(lst,i,j):
return sum(lst[i:j])