r-如何修改fasta文件的标题并将其与另一个文件的信息合并



我有一个巨大的fasta文件,它的格式是:

HQ323811.1 Abies alba tRNA Leu(trnL(基因,内含子;叶绿体gggcaagggatagggtgcagagactcaatgg

然后,我的fasta文件中有另一个文件,其中包含每个生物体的分类信息。

我想获得一份最终的fasta文件,其中只包含该物种的学名,然后是分类学信息。有办法做到这一点吗?我不知道!有人能告诉我是否有教程或我可以阅读的东西来尝试做这件事吗?

谢谢你,

设一个包含的文件in.fasta

>HQ323811.1 Abies alba tRNA-Leu (trnL) gene, intron; chloroplast
GGGCAATCCTGAGCCAAATCCGGTTCATAGAGAAAAGGGTTTCTCTCCTTCTCCTAAGGA AAGGGATAGGTGCAGAGACTCAATGG

以及每行包含一个物种的文件tax.txt,例如

Abies alba cellular organisms; Eukaryota; Viridiplantae; Streptophyta; Streptophytina; Embryophyta; Tracheophyta; Euphyllophyta; Spermatophyta; Acrogymnospermae; Pinopsida; Pinidae; Conifers I; Pinales; Pinaceae; Abies

然后你可以用R编程语言运行这个脚本

library(tidyverse)
library(tidysq)
lineages <-
# Taxon might have multiple intermediate ranks e.g. Conifers I
# Can not assume just species, genus, family, order, class, phylum, kingdom
# parsing required
read_lines("tax.txt") %>%
tibble(lineage = .) %>%
mutate(
# first two words
species = lineage %>% str_extract("^[A-z0-9]+ [A-z0-9]+"),
# last word
genus = lineage %>% str_extract("[A-z0-9]+$")
)
lineages
read_fasta("in.fasta") %>%
mutate(species = name %>% str_remove("^[A-z0-9.]+") %>% str_extract("[A-z0-9]+ [A-z0-9]+")) %>%
left_join(lineages) %>%
transmute(
sq,
name = str_glue("{row_number()}; g:{genus}, s:{species}")
) %>%
{
.x <- .
write_fasta(.x$sq, .x$name, "out.fasta")
}

导致包含的文件out.fasta

>1; g:Abies, s:Abies alba
GGGCAATCCTGAGCCAAATCCGGTTCATAGAGAAAAGGGTTTCTCTCCTTCTCCTAAGGA AAGGGATAGGTGCAGAGAC
TCAATGG

您还可以创建一个包含最终标头的文件headers.txt。第一行是in.fasta的第一个序列的新名称,依此类推:

1tax=f:Pinaceae;g:Abies;s:alba

相应的R脚本将是:

library(tidyverse)
library(tidysq)
headers <- read_lines("headers.txt")
read_fasta("in.fasta") %>%
mutate(name = headers) %>%
{
.x <- .
write_fasta(.x$sq, .x$name, "out.fasta")
}

相关内容

最新更新