script: /data/user2/work/80/statistics~sleuth "$scriptdir"/statistics~scatter_with_variance c2997108/centos6:3 c2997108/centos7:1-trinity_2.8.5-kallisto_0.46.0 c2997108/centos7:1-trinity_2.8.5-kallisto_0.46.0-blast_2.9.0-trinotate-3.1.1-R_3-kegg_2 c2997108/centos7:1-trinity_2.8.5-kallisto_0.46.0-blast_2.9.0-trinotate-3.1.1-R_4-kegg_4 centos:centos6 using docker + set -o pipefail + script=run-sleuth.0h_A_1.fq.gz.1week_A2_1.fq.gz.isoform.R + inputsamplematrix=sample.input.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt + inputmatrix=kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt + inputgenematrix=kallisto.gene.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt + inputgenetpmmatrix=kallisto.gene.TPM.not_cross_norm + inputt2g=Trinity.fasta.t2g.txt + g1=0h_A_1.fq.gz + g2=1week_A2_1.fq.gz + '[' 0h_A_1.fq.gz = '' ']' + '[' 1week_A2_1.fq.gz = '' ']' + p=0.05 ++ basename kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt + outputfile=sleuth.kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt ++ basename sample.input.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt + samplex=sleuth.sample.input.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.x ++ basename sample.input.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt + sampley=sleuth.sample.input.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.y ++ basename kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt + outputfilexup=result.sleuth.kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.0h_A_1.fq.gz.up.1week_A2_1.fq.gz.down.txt ++ basename kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt + outputfileyup=result.sleuth.kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.0h_A_1.fq.gz.down.1week_A2_1.fq.gz.up.txt ++ basename kallisto.gene.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt + outputgenefile=sleuth.kallisto.gene.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt ++ basename kallisto.gene.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt + outputgenefilexup=result.sleuth.kallisto.gene.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.0h_A_1.fq.gz.up.1week_A2_1.fq.gz.down.txt ++ basename kallisto.gene.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt + outputgenefileyup=result.sleuth.kallisto.gene.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.0h_A_1.fq.gz.down.1week_A2_1.fq.gz.up.txt + '[' Trinotate.xls3.isoform.cnt = '' ']' + annotationfile=Trinotate.xls3.isoform.cnt + '[' Trinotate.xls3.gene.cnt = '' ']' + annotationgenefile=Trinotate.xls3.gene.cnt + '[' kallisto.output = '' ']' + ln -s kallisto.output input.sleuth.extract.0h_A_1.fq.gz.1week_A2_1.fq.gz + cat + docker run -v /data/user2/work/80:/data/user2/work/80 -w /data/user2/work/80 -u root -i --rm c2997108/centos7:1-trinity_2.8.5-kallisto_0.46.0-blast_2.9.0-trinotate-3.1.1-R_4-kegg_4 R --vanilla R version 3.5.2 (2018-12-20) -- "Eggshell Igloo" Copyright (C) 2018 The R Foundation for Statistical Computing Platform: x86_64-redhat-linux-gnu (64-bit) R is free software and comes with ABSOLUTELY NO WARRANTY. You are welcome to redistribute it under certain conditions. Type 'license()' or 'licence()' for distribution details. R is a collaborative project with many contributors. Type 'contributors()' for more information and 'citation()' on how to cite R or R packages in publications. Type 'demo()' for some demos, 'help()' for on-line help, or 'help.start()' for an HTML browser interface to help. Type 'q()' to quit R. > library("sleuth") > s2c <- read.table("sample.input.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt", header = TRUE, stringsAsFactors=FALSE) > s2c <- dplyr::select(s2c, sample = id, condition) > sample_id <- s2c[,1] > kal_dirs <- file.path("input.sleuth.extract.0h_A_1.fq.gz.1week_A2_1.fq.gz", sample_id) > s2c <- dplyr::mutate(s2c, path = kal_dirs) > t2g <- read.table("Trinity.fasta.t2g.txt", header = TRUE, stringsAsFactors=FALSE) > so <- sleuth_prep(s2c, ~ condition, target_mapping = t2g, extra_bootstrap_summary = TRUE) reading in kallisto results dropping unused factor levels .. normalizing est_counts 10479 targets passed the filter normalizing tpm merging in metadata summarizing bootstraps .. > so <- sleuth_fit(so) fitting measurement error models shrinkage estimation Error in if (sum(valid) == 0) { : missing value where TRUE/FALSE needed Calls: sleuth_fit ... sliding_window_grouping -> do -> do.grouped_df -> eval_tidy -> get_quantile Execution halted + true + '[' '!' -e sleuth.kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt ']' + touch sleuth.kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt + '[' '!' -e sleuth.kallisto.gene.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt ']' + touch sleuth.kallisto.gene.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt + tail -n+2 sleuth.kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt + cut -f 1 + awk '-F\t' -v g1=0h_A_1.fq.gz '$2==g1{print $1}' sample.input.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt + awk '-F\t' -v g2=1week_A2_1.fq.gz '$2==g2{print $1}' sample.input.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt + bash /data/user2/work/80/statistics~scatter_with_variance -x 0h_A_1.fq.gz -y 1week_A2_1.fq.gz kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt sleuth.kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.id sleuth.sample.input.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.x sleuth.sample.input.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.y c2997108/centos6:3 c2997108/centos7:1-trinity_2.8.5-kallisto_0.46.0 c2997108/centos7:1-trinity_2.8.5-kallisto_0.46.0-blast_2.9.0-trinotate-3.1.1-R_3-kegg_2 c2997108/centos7:1-trinity_2.8.5-kallisto_0.46.0-blast_2.9.0-trinotate-3.1.1-R_4-kegg_4 centos:centos6 using docker + set -o pipefail + countmat=kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt + sigmat=sleuth.kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.id + listx=sleuth.sample.input.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.x + listy=sleuth.sample.input.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.y + label_x=0h_A_1.fq.gz + label_y=1week_A2_1.fq.gz + p=0.05 ++ basename kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt ++ basename sleuth.kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.id + output=kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.sleuth.kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.id + cat kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt + awk '-F\t' 'FILENAME==ARGV[1]{if(FNR>1){for(i=2;i<=NF;i++){a[i]+=$i}}} FILENAME==ARGV[2]{if(FNR==1){print $0}else{ORS=""; print $1;for(i=2;i<=NF;i++){print "\t"$i/a[i]*1000*1000}; print "\n"}}' /dev/stdin kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt + cat sleuth.sample.input.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.x + awk '-F\t' 'FILENAME==ARGV[1]{a[$1]=1} FILENAME==ARGV[2]{if(FNR==1){for(i=1;i<=NF;i++){if(a[$i]==1){b[i]=1}}}; ORS=""; print $1; for(i=2;i<=NF;i++){if(b[i]==1){print "\t"$i}}; print "\n"}' /dev/stdin kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.sleuth.kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.id.cpm + cat sleuth.sample.input.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.y + awk '-F\t' 'FILENAME==ARGV[1]{a[$1]=1} FILENAME==ARGV[2]{if(FNR==1){for(i=1;i<=NF;i++){if(a[$i]==1){b[i]=1}}}; ORS=""; print $1; for(i=2;i<=NF;i++){if(b[i]==1){print "\t"$i}}; print "\n"}' /dev/stdin kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.sleuth.kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.id.cpm + cat sleuth.kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.id + awk '-F\t' 'FILENAME==ARGV[1]{a[$1]=1} FILENAME==ARGV[2] && (FNR==1 || a[$1]==1){print $0}' /dev/stdin kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.sleuth.kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.id.cpm_x + cat sleuth.kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.id + awk '-F\t' 'FILENAME==ARGV[1]{a[$1]=1} FILENAME==ARGV[2] && (FNR==1 || a[$1]==1){print $0}' /dev/stdin kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.sleuth.kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.id.cpm_y + cat + docker run -v /data/user2/work/80:/data/user2/work/80 -w /data/user2/work/80 -u root -i --rm c2997108/centos7:1-trinity_2.8.5-kallisto_0.46.0-blast_2.9.0-trinotate-3.1.1-R_4-kegg_4 R --vanilla R version 3.5.2 (2018-12-20) -- "Eggshell Igloo" Copyright (C) 2018 The R Foundation for Statistical Computing Platform: x86_64-redhat-linux-gnu (64-bit) R is free software and comes with ABSOLUTELY NO WARRANTY. You are welcome to redistribute it under certain conditions. Type 'license()' or 'licence()' for distribution details. R is a collaborative project with many contributors. Type 'contributors()' for more information and 'citation()' on how to cite R or R packages in publications. Type 'demo()' for some demos, 'help()' for on-line help, or 'help.start()' for an HTML browser interface to help. Type 'q()' to quit R. > > label_x="0h_A_1.fq.gz" > label_y="1week_A2_1.fq.gz" > file_data_x="kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.sleuth.kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.id.cpm_x" > file_data_y="kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.sleuth.kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.id.cpm_y" > file_sig_data_x="kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.sleuth.kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.id.cpm_x_sig" > file_sig_data_y="kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.sleuth.kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.id.cpm_y_sig" > data_x=round(read.table(file_data_x,sep="\t",header=T,row.names=1)) > data_y=round(read.table(file_data_y,sep="\t",header=T,row.names=1)) > sig_data_x=read.table(file_sig_data_x,sep="\t",header=T,row.names=1) > if(dim(sig_data_x)[1]!=0){sig_data_x=round(sig_data_x)} > sig_data_y=read.table(file_sig_data_y,sep="\t",header=T,row.names=1) > if(dim(sig_data_y)[1]!=0){sig_data_y=round(sig_data_y)} > mean_data_x=apply(data_x,1,mean) > mean_data_y=apply(data_y,1,mean) > if(dim(data_x)[2]>1){ + dx=apply(data_x,1,function(i){return(abs(qt(0.05/2, length(i)-1)) * sd(i) / sqrt(length(i)))}) + }else{ + dx=0 + } > if(dim(data_y)[2]>1){ + dy=apply(data_y,1,function(i){return(abs(qt(0.05/2, length(i)-1)) * sd(i) / sqrt(length(i)))}) + }else{ + dy=0 + } > vx=dx/mean_data_x > vy=dy/mean_data_y > vx[is.na(vx)]=0 > vy[is.na(vy)]=0 > #total_var_data=sqrt(dx^2+dy^2)/sqrt(mean_data_x^2+mean_data_y^2) > #total_var_data[is.na(total_var_data)]=0 > total_var_data=sqrt(vx^2+vy^2) > mean_sig_data_x=apply(sig_data_x,1,mean) > mean_sig_data_y=apply(sig_data_y,1,mean) > min_x=min(mean_data_x[mean_data_x>0])/10 > min_y=min(mean_data_y[mean_data_y>0])/10 > plot_mean_data_x=mean_data_x+min_x > plot_mean_data_y=mean_data_y+min_y > plot_mean_sig_data_x=mean_sig_data_x+min_x > plot_mean_sig_data_y=mean_sig_data_y+min_y > max_x=max(plot_mean_data_x) > max_y=max(plot_mean_data_y) > library(Cairo) > mainplotfunc=function(){ + if(max(dx,dy)>0){ + start_class=quantile(total_var_data[total_var_data>0],c(0,1,25,50,75,90,100)/100)[2] + end_class=quantile(total_var_data[total_var_data>0],c(0,1,25,50,75,90,100)/100)[6] + + CairoPNG("kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.sleuth.kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.id.scatter.png",width=800,height=600) + num_class=22 + library(RColorBrewer) + #mycol=colorRampPalette(c(rgb(0.5,0,0.2,0.4),rgb(1,0.8,0.7,1)), alpha=TRUE)(num_class) + colPal3 <- colorRampPalette(brewer.pal(9,"YlOrRd")[3:9]) + mycol=colPal3(num_class) + + layout(matrix(c(1,1,1,1,2),1,5)) + par(ps = 18) + k=0 + i=end_class + n=(total_var_data>=i) + xt=plot_mean_data_x[n] + yt=plot_mean_data_y[n] + plot(xt,yt,log="xy",xlim=c(min_x,max_x),ylim=c(min_y,max_y),cex=2*(k+1)^3/num_class^3,pch=21,col=adjustcolor(mycol[k+1],0),bg=adjustcolor(mycol[k+1],1-0.4*k/num_class),xlab=label_x,ylab=label_y) + + for(k in 1:(num_class-2)){ + i=end_class+(start_class-end_class)*k/(num_class-2) + n=(total_var_data>=i & total_var_data=i & total_var_data=",i),pos=4) + } + k=num_class-1 + i=min(total_var_data) + par(new=T) + plot(c(0),c(k),xlim=c(0,0),ylim=c(-1,num_class),col=adjustcolor(mycol[k+1],0),bg=adjustcolor(mycol[k+1],1-0.4*k/num_class),pch=21,cex=2*(k+1)^3/num_class^3,bty="n",axes=F,xlab="",ylab="") + text(c(0.05),c(k),paste(">=",i),pos=4) + text(0,num_class+0.2,"95%conf/mean") + par(new=T) + plot(c(-0.3),c(-1),xlim=c(0,0),ylim=c(-1,num_class),col="green",cex=1,bty="n",axes=F,xlab="",ylab="") + text(c(-0.3),c(-1),paste("significant"),pos=4) + + dev.off() + } + } > subplotfunc=function(){ + CairoPNG("kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.sleuth.kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.id.scatter.normal.png",width=600,height=600) + plot(plot_mean_data_x,plot_mean_data_y,log="xy",xlim=c(min_x,max_x),ylim=c(min_y,max_y),cex=0.1,pch=21,col=adjustcolor("#000000",0.5),xlab=label_x,ylab=label_y) + par(new=T) + plot(plot_mean_sig_data_x,plot_mean_sig_data_y,log="xy",xlim=c(min_x,max_x),ylim=c(min_y,max_y),col="green",cex=2,axes=F,xlab="",ylab="") + dev.off() + } > > wflag=0 > while(wflag<=5){ + restry=NULL + restry=try(mainplotfunc(), silent = FALSE) + if(class(restry)!="try-error"){break} + try(dev.off()) + Sys.sleep(10+runif(1, min=0,max=100)/10) + wflag=wflag+1 + } > wflag=0 > while(wflag<=5){ + restry=NULL + restry=try(subplotfunc(), silent = FALSE) + if(class(restry)!="try-error"){break} + try(dev.off()) + Sys.sleep(10+runif(1, min=0,max=100)/10) + wflag=wflag+1 + } > > > + rm -f kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.sleuth.kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.id.cpm kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.sleuth.kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.id.cpm_x kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.sleuth.kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.id.cpm_y kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.sleuth.kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.id.cpm_x_sig kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.sleuth.kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.id.cpm_y_sig + post_processing + '[' 3 = 1 ']' + exit + awk '-F\t' 'NR==1{print "id\t"$0} NR>1{print $0}' sleuth.kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.tpm + awk '-F\t' 'FILENAME==ARGV[1]{x[$1]=1} FILENAME==ARGV[2]{y[$1]=1} FILENAME==ARGV[3]{if(FNR==1){for(i=2;i<=NF;i++){if(x[$i]==1){t[i]=1; xn++}else if(y[$i]==1){t[i]=2; yn++}else{t[i]=3}}}else{delete ave; for(i=2;i<=NF;i++){ave[t[i]]+=$i}; if(ave[1]/xn >= ave[2]/yn){print $1"\tx"}else{print $1"\ty"}}} ' sleuth.sample.input.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.x sleuth.sample.input.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.y /dev/stdin awk: 致命的: ファイル `sleuth.kallisto.isoform.counts.matrix.0h_A_1.fq.gz.1week_A2_1.fq.gz.txt.tpm' を読み込み用に開けません (そのようなファイルやディレクトリはありません) ++ onerror 111 ++ status=2 ++ script=/data/user2/work/80/statistics~sleuth ++ line=111 ++ shift ++ set +x ------------------------------------------------------------ Error occured on /data/user2/work/80/statistics~sleuth [Line 111]: Status 2 PID: 7020 User: user2 Current directory: /data/user2/work/80 Command line: /data/user2/work/80/statistics~sleuth ------------------------------------------------------------