# denovo-reference
/mnt/max/a/scratch/mth/olfly/tr3/README-tr3.txt
assembly command:
$TRINITY_HOME/Trinity --seqType fq --max_memory 200G --left /mnt/max/a/scratch/mth/olfly/tr3/trimmed/HI.1557.001.Index_13.Bo_Heads_Male_R1.fastq.gz,/mnt/max/a/scratch/mth/olfly/tr3/trimmed/HI.1557.001.Index_15.Bo_Heads_Ovip_Female_R1.fastq.gz,/mnt/max/a/scratch/mth/olfly/tr3/trimmed/HI.1557.001.Index_18.Bo_ovipositor_R1.fastq.gz,/mnt/max/a/scratch/mth/olfly/tr3/trimmed/HI.1557.001.Index_6.Bo_Heads_Female_R1.fastq.gz,/mnt/max/a/scratch/mth/olfly/tr3/trimmed/HI.1557.001.Index_7.Bo_Legs_ovip_Female_R1.fastq.gz,/mnt/max/a/scratch/mth/olfly/tr3/trimmed/HI.1557.002.Index_12.SexOrgans_male_R1.fastq.gz,/mnt/max/a/scratch/mth/olfly/tr3/trimmed/HI.1557.002.Index_14.Bo_Ypin_Fem_R1.fastq.gz,/mnt/max/a/scratch/mth/olfly/tr3/trimmed/HI.1557.002.Index_16.Bo_Npin_Fem_R1.fastq.gz,/mnt/max/a/scratch/mth/olfly/tr3/trimmed/HI.1557.002.Index_4.Bo_Cpin_Fem_R1.fastq.gz,/mnt/max/a/scratch/mth/olfly/tr3/trimmed/HI.1557.002.Index_5.Testes_after_R1.fastq.gz,/mnt/max/a/scratch/mth/olfly/tr3/trimmed/HI.1557.003.Index_19.Bo_egg_R1.fastq.gz,/mnt/max/a/scratch/mth/olfly/tr3/trimmed/HI.1557.003.Index_1.Bo_larva_R1.fastq.gz,/mnt/max/a/scratch/mth/olfly/tr3/trimmed/HI.1557.003.Index_3.Bo_pupa_R1.fastq.gz,/mnt/max/a/scratch/mth/olfly/tr3/trimmed/HI.1557.003.Index_8.Bo_Thorx_Mix_R1.fastq.gz,/mnt/max/a/scratch/mth/olfly/tr3/trimmed/HI.1557.003.Index_9.Bo_adult_R1.fastq.gz,/data/images/proton/run125/R_2015_01_21_12_34_28_user_ION-125-KMathiopoulos_RNAseq_150121_KMR27-28-31-32-33.KMR27.IonXpressRNA_011.fastq,/data/images/proton/run125/R_2015_01_21_12_34_28_user_ION-125-KMathiopoulos_RNAseq_150121_KMR27-28-31-32-33.KMR28.IonXpressRNA_012.fastq,/data/images/proton/run125/R_2015_01_21_12_34_28_user_ION-125-KMathiopoulos_RNAseq_150121_KMR27-28-31-32-33.KMR31.IonXpressRNA_013.fastq,/data/images/proton/run125/R_2015_01_21_12_34_28_user_ION-125-KMathiopoulos_RNAseq_150121_KMR27-28-31-32-33.KMR32.IonXpressRNA_016.fastq,/data/images/proton/run125/R_2015_01_21_12_34_28_user_ION-125-KMathiopoulos_RNAseq_150121_KMR27-28-31-32-33.KMR33.IonXpressRNA_008.fastq,/data/images/proton/run124/R_2015_01_20_12_13_52_user_ION-124-KMathiopoulos_RNAseq_150120_KMR25-26-29-30.KMR25.IonXpressRNA_001.fastq,/data/images/proton/run124/R_2015_01_20_12_13_52_user_ION-124-KMathiopoulos_RNAseq_150120_KMR25-26-29-30.KMR26.IonXpressRNA_003.fastq,/data/images/proton/run124/R_2015_01_20_12_13_52_user_ION-124-KMathiopoulos_RNAseq_150120_KMR25-26-29-30.KMR29.IonXpressRNA_004.fastq,/data/images/proton/run124/R_2015_01_20_12_13_52_user_ION-124-KMathiopoulos_RNAseq_150120_KMR25-26-29-30.KMR30.IonXpressRNA_009.fastq --right /mnt/max/a/scratch/mth/olfly/tr3/trimmed/HI.1557.001.Index_13.Bo_Heads_Male_R2.fastq.gz,/mnt/max/a/scratch/mth/olfly/tr3/trimmed/HI.1557.001.Index_15.Bo_Heads_Ovip_Female_R2.fastq.gz,/mnt/max/a/scratch/mth/olfly/tr3/trimmed/HI.1557.001.Index_18.Bo_ovipositor_R2.fastq.gz,/mnt/max/a/scratch/mth/olfly/tr3/trimmed/HI.1557.001.Index_6.Bo_Heads_Female_R2.fastq.gz,/mnt/max/a/scratch/mth/olfly/tr3/trimmed/HI.1557.001.Index_7.Bo_Legs_ovip_Female_R2.fastq.gz,/mnt/max/a/scratch/mth/olfly/tr3/trimmed/HI.1557.002.Index_12.SexOrgans_male_R2.fastq.gz,/mnt/max/a/scratch/mth/olfly/tr3/trimmed/HI.1557.002.Index_14.Bo_Ypin_Fem_R2.fastq.gz,/mnt/max/a/scratch/mth/olfly/tr3/trimmed/HI.1557.002.Index_16.Bo_Npin_Fem_R2.fastq.gz,/mnt/max/a/scratch/mth/olfly/tr3/trimmed/HI.1557.002.Index_4.Bo_Cpin_Fem_R2.fastq.gz,/mnt/max/a/scratch/mth/olfly/tr3/trimmed/HI.1557.002.Index_5.Testes_after_R2.fastq.gz,/mnt/max/a/scratch/mth/olfly/tr3/trimmed/HI.1557.003.Index_19.Bo_egg_R2.fastq.gz,/mnt/max/a/scratch/mth/olfly/tr3/trimmed/HI.1557.003.Index_1.Bo_larva_R2.fastq.gz,/mnt/max/a/scratch/mth/olfly/tr3/trimmed/HI.1557.003.Index_3.Bo_pupa_R2.fastq.gz,/mnt/max/a/scratch/mth/olfly/tr3/trimmed/HI.1557.003.Index_8.Bo_Thorx_Mix_R2.fastq.gz,/mnt/max/a/scratch/mth/olfly/tr3/trimmed/HI.1557.003.Index_9.Bo_adult_R2.fastq.gz --CPU 20  &

reczko@max:~/scratch/mth/olfly/tr3b$ $TRINITY_HOME/util/TrinityStats.pl trinity_out_dir/Trinity.fasta
Average contig: 669.78
	Total assembled bases: 203690146


PATH=/data/results/tools/align/rsem/rsem-1.2.15:/data/results/tools/samtools/samtools-0.1.19:/data/results/tools/align/bowtie2-2.2.8:/data/results/tools/align/tophat-2.1.1.Linux_x86_64:$PATH; export PATH
export TRINITY_HOME=/data/results/tools/denovo/trinity/trinityrnaseq-Trinity-v2.4.0
reczko@max:/mnt/max/a/scratch/mth/olfly/genes-paper/extra-alignments$ $TRINITY_HOME/util/align_and_estimate_abundance.pl --prep_reference --thread_count 20 --SS_lib_type F --transcripts /mnt/max/a/scratch/mth/olfly/tr3b/trinity_out_dir/Trinity.fasta-chopped-id.fa --seqType fq --single ../reads/R_2014_10_31_11_54_20_user_ION-104-KMathiopoulos_RNAseq_141031_KMR5-8.KMR5.IonXpressRNA_001.fastq  --est_method RSEM --aln_method bowtie2 --trinity_mode --output_dir tr3-KMR5 &> tr3-KMR5.edgeR-olfly.log
$TRINITY_HOME/util/align_and_estimate_abundance.pl --thread_count 20 --SS_lib_type F --transcripts /mnt/max/a/scratch/mth/olfly/tr3b/trinity_out_dir/Trinity.fasta-chopped-id.fa --seqType fq --single ../reads/R_2014_11_07_11_42_50_user_ION-108-KMathiopoulos_RNAseq_141107_KMR11-12-15-16.KMR11.IonXpressRNA_006.fastq  --est_method RSEM --aln_method bowtie2 --trinity_mode --output_dir tr3-KMR11 &> tr3-KMR11.edgeR-olfly.log
$TRINITY_HOME/util/align_and_estimate_abundance.pl --thread_count 20 --SS_lib_type F --transcripts /mnt/max/a/scratch/mth/olfly/tr3b/trinity_out_dir/Trinity.fasta-chopped-id.fa --seqType fq --single ../reads/R_2014_11_07_11_42_50_user_ION-108-KMathiopoulos_RNAseq_141107_KMR11-12-15-16.KMR12.IonXpressRNA_010.fastq  --est_method RSEM --aln_method bowtie2 --trinity_mode --output_dir tr3-KMR12 &> tr3-KMR12.edgeR-olfly.log


ln -s ./tr3-KMR12/RSEM.genes.results 12
ln -s ./tr3-KMR5/RSEM.genes.results 5
ln -s ./tr3-KMR11/RSEM.genes.results 11

ln -s /data/images/proton/run124/tr3-004_v2/RSEM.genes.results 29
ln -s /data/images/proton/run124/tr3-001_v2/RSEM.genes.results 25
ln -s /data/images/proton/run124/tr3-009_v2/RSEM.genes.results 30
ln -s /data/images/proton/run124/tr3-003_v2/RSEM.genes.results 26

ln -s /data/images/proton/run125/tr3-016/RSEM.genes.results 32
ln -s /data/images/proton/run125/tr3-008/RSEM.genes.results 33
ln -s /data/images/proton/run125/tr3-013/RSEM.genes.results 31
ln -s /data/images/proton/run125/tr3-012/RSEM.genes.results 28
ln -s /data/images/proton/run125/tr3-011/RSEM.genes.results 27

ln -s 27 FAGsB1
ln -s 28 FAGsB2
ln -s 31 MAGSB1
ln -s 32 MAGSB2
ln -s 25 FAGsA1
ln -s 26 FAGsA2
ln -s 29 MAGSA1
ln -s 30 MAGSA2
ln -s 5  Testes_A
ln -s 12 Testes_B


KMR25: FAGsA1
KMR26: FAGsA2
KMR27: FAGsB1
KMR28:FAGsB2
KMR29:MAGsA1
KMR30:MAGsA2
KMR31:MAGsB1
KMR32:MAGsB2
KMR5:Testes_A
KMR12:Testes_B



$TRINITY_HOME/util/abundance_estimates_to_matrix.pl --est_method RSEM  --out_prefix rsem-gene-expressions FAGsB1  FAGsB2  MAGSB1  MAGSB2  FAGsA1  FAGsA2  MAGSA1  MAGSA2  Testes_A Testes_B


library(ggfortify)
t=read.table("rsem-gene-expressions.TMM.EXPR.matrix")
d=t(log(1+t))
pca=prcomp(d,scale.=FALSE)
pdf(file=paste("PCA_log_TMM_",i,"_vs_",j,".pdf",sep=""))
plot(pca$x[,c(i,j)],pch=".",xlab=paste("PC",i,sep=""),ylab=paste("PC",j,sep=""),main=paste("PCA PC",i," vs. PC",j,sep=""))
dev.off()
pdf(file=paste("PCA_log_TMM_",i,"_vs_",j,".pdf",sep=""))
autoplot(pca,data=d,label = TRUE, label.size = 4,shape=FALSE)
dev.off()

t=read.table("rsem-gene-expressions.TMM.EXPR.matrix")
#sel=c(1,2,3,4,5,6,7,8,9,10,11)
mydata2=t  #[,sel]
mydata=t(mydata2)
d <- dist(mydata) # euclidean distances between the rows
fit <- cmdscale(d,eig=TRUE, k=2) # k is the number of dim
pdf("mds-rsem-tmm-norm.pdf")
x <- fit$points[,1]
y <- fit$points[,2]
plot(x, y, xlab="Coordinate 1", ylab="Coordinate 2", 
  main="Multidimensional Scaling", type="n")
text(x, y, labels = row.names(mydata), cex=.7)
dev.off()






$TRINITY_HOME/util/abundance_estimates_to_matrix.pl --est_method RSEM  --out_prefix genes 5 11 12 25 26 27 28 29 30 31 32 33
* Outputting combined matrix.

/mnt/fix/c/solid.data_results.recovered.data/tools/denovo/trinity/trinityrnaseq-Trinity-v2.4.0/util/support_scripts/run_TMM_scale_matrix.pl --matrix genes.TPM.not_cross_norm > genes.TMM.EXPR.matrixCMD: R --vanilla -q < __tmp_runTMM.R 1>&2 
> library(edgeR)
Loading required package: limma
> 
> rnaseqMatrix = read.table("genes.TPM.not_cross_norm", header=T, row.names=1, com='', check.names=F)
> rnaseqMatrix = as.matrix(rnaseqMatrix)
> rnaseqMatrix = round(rnaseqMatrix)
> exp_study = DGEList(counts=rnaseqMatrix, group=factor(colnames(rnaseqMatrix)))
> exp_study = calcNormFactors(exp_study)
> exp_study$samples$eff.lib.size = exp_study$samples$lib.size * exp_study$samples$norm.factors
> write.table(exp_study$samples, file="genes.TPM.not_cross_norm.TMM_info.txt", quote=F, sep="\t", row.names=F)
> 
Done.


t=read.table("genes.counts.matrix")
sel=c(1,2,3,4,5,6,7,8,9,10,11)
mydata2=t[,sel]
mydata=t(mydata2)
d <- dist(mydata) # euclidean distances between the rows
fit <- cmdscale(d,eig=TRUE, k=2) # k is the number of dim
pdf("mds-all.pdf")
x <- fit$points[,1]
y <- fit$points[,2]
plot(x, y, xlab="Coordinate 1", ylab="Coordinate 2", 
  main="Multidimensional Scaling", type="n")
text(x, y, labels = row.names(mydata), cex=.7)
dev.off()

t=read.table("genes.TMM.EXPR.matrix")
sel=c(1,2,3,4,5,6,7,8,9,10,11)
mydata2=t[,sel]
mydata=t(mydata2)
d <- dist(mydata) # euclidean distances between the rows
fit <- cmdscale(d,eig=TRUE, k=2) # k is the number of dim
pdf("mds-all-norm.pdf")
x <- fit$points[,1]
y <- fit$points[,2]
plot(x, y, xlab="Coordinate 1", ylab="Coordinate 2", 
  main="Multidimensional Scaling", type="n")
text(x, y, labels = row.names(mydata), cex=.7)
dev.off()

t=read.table("genes.counts.matrix")
sel=c(1,2,4,5,6,7,8,9,10,11)
mydata2=t[,sel]
mydata=t(mydata2)
d <- dist(mydata) # euclidean distances between the rows
fit <- cmdscale(d,eig=TRUE, k=2) # k is the number of dim
pdf("mds-all-wo12.pdf")
x <- fit$points[,1]
y <- fit$points[,2]
plot(x, y, xlab="Coordinate 1", ylab="Coordinate 2", 
  main="Multidimensional Scaling", type="n")
text(x, y, labels = row.names(mydata), cex=.7)
dev.off()

t=read.table("genes.TMM.EXPR.matrix")
sel=c(1,2,4,5,6,7,8,9,10,11)
mydata2=t[,sel]
mydata=t(mydata2)
d <- dist(mydata) # euclidean distances between the rows
fit <- cmdscale(d,eig=TRUE, k=2) # k is the number of dim
pdf("mds-all-wo12-norm.pdf")
x <- fit$points[,1]
y <- fit$points[,2]
plot(x, y, xlab="Coordinate 1", ylab="Coordinate 2", 
  main="Multidimensional Scaling", type="n")
text(x, y, labels = row.names(mydata), cex=.7)
dev.off()