利用shell对物种丰度表进行处理：
1.原始丰度表如下：
image.png
每列由空格隔开，物种由分号隔开
处理如下：（注：处理完需要在物种分类手动加入界门纲目科属种）
#因为物种分类也有空格，因此需要和丰度分开处理，首先去掉最前面的空格，然后以分号分隔取丰度这一边，然后将空格转换为制表符
cat tax.xls | sed 's/^\ //g' | awk -F ";" '{print $1}' | sed 's/\ /\t/g' > tax11.txt
#然后取物种分类这边，将分号替换为制表符
cat tax.xls | sed 's/^\ //g' | awk -F " " '{for (i=37;i<=NF;i++)printf("%s ", $i);print ""}' | sed 's/\;/\t/g' > tax_z.txt
处理完物种分类如下：
image.png
然后将两个文件导入R语言：
ab <- read.delim("tax11.txt",sep = "\t",header = T,row.names= NULL)
species <- read.delim("tax_z.txt",sep = "\t",header = T,row.names= NULL)
species_ab <- cbind(ab,species)
species_ab <- species_ab[, ！colnames(species_ab) == "Taxonomy"]
#按照不同的界提取子表格
k__Archaea <- subset(species_ab,species_ab$Kingdom == 'k__Archaea')
k__Bacteria <- subset(species_ab,species_ab$Kingdom == 'k__Bacteria')
k__Eukaryota <- subset(species_ab,species_ab$Kingdom == 'k__Eukaryota')
k__unidentified <- subset(species_ab,species_ab$Kingdom == 'k__unidentified')
k__Viruses <- subset(species_ab,species_ab$Kingdom == 'k__Viruses')
然后就能对这些界的物种进行进一步的处理，比如，想得到细菌界所有纲的累积丰度：
之前在diamond比对以及结果数据中写过一个叠加小代码。
然而，发现这代码运行在R中运行非常慢，不适合大数据处理。
今天换一个思路写下代码，一般几百Mb大小的丰度文件，也就一分钟不到就运行完了。
附上代码：（注：文件处理格式，前几列为样本丰度，后几列是物种分类）
#参数data指数据需要处理的数据框，lie指需要对哪个分类等级进行叠加，sample指样本数目。
diejia <- function(data = data.frame(),lie = "Class" ,sample = 1){
list1 <- list()
x2 <- data.frame()
sum_ab <- data.frame()
xx <- data.frame()
xxx <- data.frame()
engame <- data.frame()
z <- levels(factor(data[,lie])) #提取指定列的不同的分类名
for(i in z){
x2 <- subset(data,data[,lie] %in% i)
list1[[i]] <- x2
} #第一个循环将相同名字的所有行形成新的数据框，并加入到列表中
for (o in 1:length(z)) {
specie <- list1[[o]][1,(sample+1):ncol(list1[[o]])]
xxx <- list1[[o]][,1:sample]
sum_ab <- colSums(xxx)
sum_ab <- as.data.frame(t(as.data.frame(sum_ab)))
xx <- cbind(sum_ab,specie)
engame <- rbind(engame,xx)
} #第二个循环，从列表中取每一个数据框，对他们求每列的和，也就是该物种在样本中的所有丰度。然后再将这个丰度与每一个数据框所处理 的丰度合并成一个新的数据框。
return(engame)
}
#写完代码就可以运行了
Bacteria_class <- diejia(data = k__Bacteria , lie = "Class" ,sample = 36)
最后得到每个纲的丰度表格
image.png
