str="D10_vs_D4"
str=$(echo "$str" | sed 's/_vs_//')
echo "$str"
2. seqkit stats C1_1.fq.gz -a 更多信息 -T tab分割
命令：fasta文件按照500个序列进行**拆分**
seqkit split -i -s 500 zea.pep.fasta -O zea 按照500个序列分割文件 `-i`: **这是一个标志** 输入fasta文件使用
3. 提取序列展示开头两行的序列 :: seqkit head -n 2 fa
4. 提取ID seqkit seq -n
tr|C0PFK3|C0PFK3_MAIZE Proteasome subunit alpha type OS=Zea mays OX=4577 GN=ZEAMMB73_Zm00001d027896 PE=2 SV=1
5. 仅提取第一个空格前的ID：seqkit seq pep.fa -n -i
tr|C0PFK3|C0PFK3_MAIZE通过ID过滤序列
seqkit grep -f gy.id office_cds.fa >office_cds_gy.fa
seqkit grep -n -f id.txt C1_1.fq.gz -o result.fa.gz ID中可能包含空格
碱基顺序比对工具（按碱基标颜色）- 可用于查找保守区 快捷方式
亚细胞定位预测
批量蛋白预测 5000个蛋白序列：https://busca.biocomp.unibo.it/
ngloc 原核
[亚细胞定位linu版(https://zhuanlan.zhihu.com/p/692973853)
/home/ug2065/software/WoLFPSort-master/bin/runWolfPsortSummary plant < ../del_Creinhardtii_281_v5.6.protein.fa >output001.xls直接画图
library(ggplot2)
library(tidyverse)
swolf_file <- "C:/Users/weimi/Desktop/sub.txt"
data_file <- ""
####整理整个蛋白组
df <- read.delim(swolf_file, skip = 1, header=F, check.names = F)
tidy_df <- df %>%
separate(V1, into=c("V1", "other"), sep=", ", extra = "merge", fill="right") %>%
separate(V1, into=c("V2", "loc_score"), sep=" ",extra = "merge", fill="right") %>%
mutate(loc=gsub(" .*","", loc_score)) %>% select(V2, loc)
###合并得到鉴定到的
datas <- fread(data_file)
draw <- datas %>% left_join(out[c("V2","loc")], by=c("Protein accession"="V2"))
#统计and画图 注意去重
draw_stat <- draw %>% distinct(`Protein accession`,loc) %>% dplyr::count(loc)
p <- ggplot(draw_stat, aes(reorder(loc,-n), n, fill=loc)) + geom_col() +labs(x=NULL,y=NULL)
ggsave(plot=p,"loc.png",height = 7,width=7)
ggsave(plot=p,"loc.pdf",height = 7,width=7)
fwrite(draw,"loc_anno.xls",sep="\t")
library(dplyr)
library(forcats)
df <- read.delim("~/software/WoLFPSort-master/result.txt", skip = 1, header=F, check.names = F)
out <- df %>%
separate(V1, into=c("V1", "other"), sep=", ", extra = "merge", fill="right") %>%
separate(V1, into=c("V2", "loc_score"), sep=" ",extra = "merge", fill="right") %>%
mutate(loc=gsub(" .*","", loc_score))
stat_res <- out %>% group_by(loc) %>% dplyr::summarise(n=n()) %>%
group_by(loc) %>%
dplyr::summarise(n = sum(n)) %>%
mutate(percent=paste0(" (", round(n/sum(n), 2), "%)"))
# 假设 out 是原始数据框，含 loc 列
#按数量统计
stat_res2 <- out %>%
group_by(loc) %>%
dplyr::summarise(n = n()) %>%
mutate(loc = ifelse(n < 10, "Other", loc)) %>%
group_by(loc) %>%
dplyr::summarise(n = sum(n)) %>%
mutate(percent=paste0(" (", round(n/sum(n), 2), "%)"))
#按照类别统计
stat_res2 <- out %>%
mutate(loc2 = fct_lump_n(loc, n = 2)) %>% # 只保留前5类，其余合并
dplyr::count(loc2) %>%
rename(n = n) %>%
mutate(percent=paste0(" (", round(n/sum(n), 2), "%)"))
write.table(stat_res, sep="\t", row.names = F, col.names = T, quote = F)
write.table(stat_res2, sep="\t", row.names = F, col.names = T, quote = F))
# loc n percent label ypos
# <chr> <int> <dbl> <chr> <dbl>
# 1 nucl 190 93.6 nucl (93.6%) 95
# 2 chlo 6 3 chlo (3%) 193
# 3 cyto 3 1.5 cyto (1.5%) 198.
# 4 extr 1 0.5 extr (0.5%) 200.
# 5 mito 1 0.5 mito (0.5%) 200.
# 6 pero 1 0.5 pero (0.5%) 202.
# 7 vacu 1 0.5 vacu (0.5%) 202.
out %>%
fct_lump_min(loc, n = 2)) %>% # 只保留前5类，其余合并
dplyr::count(loc2) %>%
rename(n = n) %>%
mutate(percent=paste0(" (", round(n/sum(n), 2), "%)"))
# out$loc %>% fct_lump_min(2) %>% table() #按照数目分
# out$loc %>% fct_lump_n(3) %>% table() #按照类别分
对应缩写
| chlo | 叶绿体 | chloroplast | |
| chlo_mito | 叶绿体和线粒体 | Chloroplasts;mitochondrion | |
| cysk | 细胞骨架 | cytoskeleton | |
| cyto | 细胞质 | cytoplasm | |
| cyto_nucl | 细胞质和核仁 | cytoplasm；Nucleus | |
| E.R. | 内质网 | endoplasmic reticulum | |
| E.R._plas | 内质网和质膜 | endoplasmic reticulu plasma membrane | |
| extr | 细胞外基质 | extracell | |
| golg | 高尔基体 | Golgi apparatus | |
| mito | 线粒体 | mitochondrion | |
| nucl | 细胞核 | Nucleus | |
| nucl_plas | 核仁和质膜 | Nucleus； plasma membrane | |
| pero | 过氧化酶体 | peroxisome | |
| plas | 质膜 | plasma membrane; | |
| vacu | 液泡 | vacuole | |
| 溶酶体 | Lysosome | ||
| 核内体 | Endosome | ||
| 细胞壁 | cell wall | ||
| 中心体 | centriiole | ||
| 微粒体 | microsome |
亚细胞定位下载
下载位置
Subcellular location [CC]
Gene Ontology (cellular component)|
亚细胞统计
Transclude of ../../00附件/01理论/软件与小工具/seqkit/img-20250516160925905.txt
sub.go.list=c(“GO:0009507”,
“GO:0005737”,“GO:0005856”,“GO:0005783”,“GO:0005576”,“GO:0005794”,“GO:0005739”,“GO:0016604”,“GO:0005777”,“GO:0005886”,“GO:0005774”)
term_count <- sapply(df2$location, function(go_terms) {
sapply(desc_df$go_id, function(term) grepl(term, go_terms))
})
亚细胞定位可视化
library(ggplot2)
library(ggsci)
library(tidyverse)
draw<- df
num <- nrow(df)
draw <- draw %>% arrange(count)
draw$location <- factor(draw$location, level=(draw$location))
p <- ggplot(draw, aes(reorder(location, count), count, fill=location)) + geom_col() + coord_flip()+
scale_fill_manual(values=colorRampPalette(pal_npg()(10))((num)))+
geom_text(aes(label=count), hjust=-0.5)+labs(x=NULL, fill=NULL) +
ylim(0, max(draw$count*1.2))+
guides(fill="none")+theme_classic()
p <- p+ theme(text = element_text(size=12, color="black"),
axis.title = element_text(size=12, color="black",face="bold"),
axis.text = element_text(size=12, color="black"),
legend.title = element_text(size=12, color="black",face="bold"),
panel.grid = element_blank()
)
p
ggsave(paste0(prefix, ".png"), p, width=7, height=7)
ggsave(paste0(prefix, ".pdf"), p, width=7, height=7)
###输入
df <- read.delim("C:/Users/weimi/test_count.tsv")
prefix <-
GO 统计直接画图
df=read.delim("C:/Users/weimi/123_go_count.tsv", check.names = F, header = T)
draw<- df
draw <- draw %>% arrange(count)
draw$description <- factor(draw$description, levels=draw$description)
draw$namespace <- factor(draw$namespace, levels=unique(
draw$namespace[order(draw$description, decreasing = T)]))
p <- ggplot(draw, aes(description, count, fill=namespace))+geom_col() + coord_flip() +
labs(x=NULL, fill=NULL) +geom_text(aes(label=count), hjust=-0.5) +
ylim(0,max(draw$count*1.2)) +
scale_fill_manual(values=colorRampPalette(pal_npg()(10))((11)))
p <- p + theme(text = element_text(size=12, color="black"),
axis.title = element_text(size=12, color="black",face="bold"),
axis.text = element_text(size=12, color="black"),
legend.title = element_text(size=12, color="black",face="bold"),
panel.grid = element_blank()
)