kraken获得丰度、lefese分析

for i in 141 97 83
do
 kraken2 --db /mnt/database/kraken/minidb/minikraken2_v1_8GB --report $i.RNA --use-mpa-style --threads 32 --paired /data/72sample/sub72/RNA/$i.1.fq /data/72sample/sub72/RNA/$i.2.fq > $i.RNA.sam && rm -f $i.RNA.sam
done

#DNAkraken
setwd("/mnt/mzy/dairycow/72samples/kraken")
DNAsample<-c("123","98","126","4","128","32","134","8","130","29","131","42","119","19","139","24","125","37","144","14","117","5","142","13","127","22","129","48","138","47","140","20","132","18","149","73","120","105","148","44","122","30","124","26","143","34","133","12","145","21","152","10","118","102","147","80","121","1","151","35","135","38","136","2","150","17","146","55","141","11","137","97")
system("awk -F '\t' '{print $1}' *DNA|sort|uniq >uniq.names")
library(tidyverse)
read.table("uniq.names",header=F,sep="\t",quote="")->a
unique(a)->a
for (i in DNAsample){read.table(paste(i,".DNA",sep=""),header=F,sep="\t",quote="",fill=T)->b; colnames(b)[2]<-i;left_join(a,b,by="V1")->a;a[!duplicated(a[,1]),]->a;print(i);print(date())}
a[is.na(a)]<-0
write.table(a,file="krakenDNA.csv",row.names=F,col.names=T,quote=F,sep=",")

#RNAkraken
setwd("/mnt/mzy/dairycow/72samples/kraken")
RNAsample<-c("123","98","126","107","128","32","134","8","130","29","131","42","119","108","139","24","125","37","144","83","117","5","142","74","127","22","129","81","138","99","140","20","132","82","149","73","120","105","148","44","122","101","124","100","143","76","133","12","145","79","152","77","118","102","147","80","121","106","151","104","135","38","136","2","150","78","146","55","141","11","137","97")
system("awk -F '\t' '{print $1}' *RNA|sort|uniq >uniq.names")
library(tidyverse)
read.table("uniq.names",header=F,sep="\t",quote="")->a
unique(a)->a
for (i in RNAsample){read.table(paste(i,".RNA",sep=""),header=F,sep="\t",quote="",fill=T)->b; colnames(b)[2]<-i;left_join(a,b,by="V1")->a;a[!duplicated(a[,1]),]->a;print(i);print(date())}
a[is.na(a)]<-0
write.table(a,file="krakenRNA.csv",row.names=F,col.names=T,quote=F,sep=",")

  • 0
    点赞
  • 1
    收藏
    觉得还不错? 一键收藏
  • 0
    评论

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值