Patricia=read.table('Comparison_array_seq.dat')
DRIP_low=read.csv('Patricized_DRIP_Lower_in_KD.tsv',sep='\t')
DRIP_high=read.csv('Patricized_DRIP_Higher_in_KD.tsv',sep='\t')

low_array=c()
low_seq=c()
medium_array=c()
medium_seq=c()
high_array=c()
high_seq=c()

for (gene in Patricia$V1)
{
class=1
individual_genes=unique(strsplit(as.character(gene),'/',fixed=T)[[1]])
for (lo in DRIP_low$Gene.name)
{
individual=unique(strsplit(as.character(lo),' ',fixed=T)[[1]])
for (low_gene in individual)
{
if (is.element(low_gene,individual_genes))
{
class=0
}
}
}
for (hi in DRIP_high$Gene.name)
{
individual=unique(strsplit(as.character(hi),' ',fixed=T)[[1]])
for (high_gene in individual)
{
if (is.element(high_gene,individual_genes))
{
class=2
}
}
}
if (class==0)
{
low_array=append(low_array,Patricia$V2[Patricia$V1==gene])
low_seq=append(low_seq,Patricia$V4[Patricia$V1==gene])
}
if (class==1)
{
medium_array=append(medium_array,Patricia$V2[Patricia$V1==gene])
medium_seq=append(medium_seq,Patricia$V4[Patricia$V1==gene])
}
if (class==2)
{
high_array=append(high_array,Patricia$V2[Patricia$V1==gene])
high_seq=append(high_seq,Patricia$V4[Patricia$V1==gene])
}
}

pdf('DRIP_vs_microarray_transcriptomics.pdf',width=6,height=6)
par(mar=c(5,6,4,2))
p=kruskal.test(c(low_array,medium_array,high_array)~as.factor(c(rep(0,times=length(low_array)),rep(1,times=length(medium_array)),rep(2,times=length(high_array)))))$p.value
boxplot(low_array,medium_array,high_array,axes=F,ylab='Expression change upon SETX k.d. (log2(FC))\naccording to microarray',main=paste('Kruskal-Wallis p-value=',signif(p,digits=3),sep=''))
axis(1,labels=paste(c('DRIP lower\nin SETX k.d.\n(n=','DRIP unaffected\n(n=','DRIP higher\nin SETX k.d.\n(n='),c(length(low_array),length(medium_array),length(high_array)),c(')',')',')'),sep=''),at=c(1:3),padj=0.75)
axis(2)
dev.off()
pdf('DRIP_vs_RNA-Seq_transcriptomics.pdf',width=6,height=6)
par(mar=c(5,6,4,2))
p=kruskal.test(c(low_seq,medium_seq,high_seq)~as.factor(c(rep(0,times=length(low_seq)),rep(1,times=length(medium_seq)),rep(2,times=length(high_seq)))))$p.value
boxplot(low_array,medium_array,high_array,axes=F,ylab='Expression change upon SETX k.d. (log2(FC))\naccording to RNA-Seq',main=paste('Kruskal-Wallis p-value=',signif(p,digits=3),sep=''))
axis(1,labels=paste(c('DRIP lower\nin SETX k.d.\n(n=','DRIP unaffected\n(n=','DRIP higher\nin SETX k.d.\n(n='),c(length(low_seq),length(medium_seq),length(high_seq)),c(')',')',')'),sep=''),at=c(1:3),padj=0.75)
axis(2)
dev.off()
