提取16列文件做進化樹
# 讀取文件
data <- read.csv("1vir.csv", header = TRUE)
# 檢查數(shù)據(jù)列數(shù)是否足夠
if (ncol(data) < 16) {
stop("數(shù)據(jù)列數(shù)少于 16 列")
}
# 將數(shù)據(jù)轉換為數(shù)據(jù)框（確保是 data.frame 類型）
data <- as.data.frame(data)
# 創(chuàng)建一個新的數(shù)據(jù)框只包含第一列和第十六列冠桃，同時去除包含缺失值的行
new_data <- data.frame(
col1 = data[, 1][complete.cases(data[, c(1, 16)])],
col16 = data[, 16][complete.cases(data[, c(1, 16)])]
)
# 根據(jù)第一列合并第十六列（使用 aggregate 函數(shù)）
merged_data <- aggregate(col16 ~ col1, data = new_data, FUN = function(x) paste(x, collapse = ","))
# 重命名列
colnames(merged_data) <- c("第一列", "合并后的第十六列")
# 保存結果為新的 CSV 文件
write.csv(merged_data, "merged_data.csv", row.names = FALSE)
# 讀取 CSV 文件
data <- read.csv("merged_data.csv", header = TRUE)
# 提取第一列和第二列
col1 <- data[, 1]
col2 <- data[, 2]
# 構建 FASTA 格式的字符串
fasta_str <- ""
for (i in 1:length(col1)) {
fasta_str <- paste0(fasta_str, paste0(">", col1[i], "\n", col2[i], "\n"))
}
# 保存為 FASTA 文件
write(fasta_str, file = "output.fasta")
###BiocManager::install("msa")
require(msa)
mySequenceFile <- readAAStringSet("output.fasta")
#多序列比對
myFirstAlignment <- msa(mySequenceFile)
head(mySequenceFile)
library(ggplot2)
require(seqinr)
myAlignment <- msaConvert(myFirstAlignment, type="seqinr::alignment")
d <- dist.alignment(myAlignment, "identity")
#構建NJ樹
require(ape)
tree <- nj(d)
#畫樹并輸出到PDF文件ggtree.pdf
require(ggtree)
環(huán)形進化樹
##環(huán)狀圖
p1<-ggtree(tree, layout='circular', ladderize=FALSE, size=0.8, branch.length="none",col="red")+
geom_tiplab2(hjust=-0.3)+
geom_tippoint(size=1.5,col="blue")+
geom_nodepoint(color="black", alpha=1/4, size=2) +
theme(legend.title=element_text(face="bold"), legend.position="bottom", legend.box="horizontal", legend.text=element_text(size=rel(0.5)))
p1
# 圖例位置、文字大小
###長方形圖
p2<- ggtree(tree, layout='rectangular', size=0.8, col="deepskyblue3") +
geom_tiplab(size=3, color="purple4", hjust=-0.05)+
geom_tippoint(size=1.5, color="deepskyblue3")+
geom_nodepoint(color="pink", alpha=1/4, size=5)+
theme_tree2()
ggsave(p1, file="tree_circular.pdf", width=9, height=9)
ggsave(p2, file="tree_rectangular.pdf", width=9, height=9)
