load(file="TransBIG/TransBIG-expression.rdata")
clinical2=read.csv("TransBIG/TransBIG-Sample_info.csv")
load("top100.rdata")

################# ANALYZA PREZITIA ##################
# 1) chceme zistit rozdiel v preziti medzi skupinami ER+ a ER-
source("surv.plot.R") # nacitame funkciu pre zobrazenie Kaplan meierovych kriviek
library(survival)

rfs=clinical2$Time.RFS..days
rfs.event=clinical2$Event.RFS

surv.plot(survfit(Surv(rfs, rfs.event)~clinical2$ER.pos.1.yes.))

# odhadneme median prezitia
survfit(Surv(rfs, rfs.event)~clinical2$ER.pos.1.yes.)

# je rozdiel v preziti statisticky vyznamny?
survdiff(Surv(rfs, rfs.event)~clinical2$ER.pos.1.yes.)

# nie je, aj ked trend tam je
# teraz to skusme so skupinami definovanymi zhlukovanim
# musime ziskat dendrogram strom, heatmapa vnutorne aplikovala zhlukovaci algoritmus, ale my to musime spravit znovu

eukl2=dist(X[,intersect(colnames(X),top100)])
hclust2=hclust(eukl2)
plot(hclust2)
groups2=cutree(hclust2, k=2)# delime na dve skupiny
groups2

## Teraz analyza prezitia s novymi skupinami ##
surv.plot(survfit(Surv(rfs, rfs.event)~groups2))

#### Domáce úlohy:
#1. spravte podobnu analyzu s top 200 genmi a top 500 genmi
#2. potom spravte podobnu analyzu s HER2
  
# Teraz skúsime spojité modelovanie odpovede pomocou Coxovho modelu
# aplikujeme normalizáciu (median centering a IQR delenie)

# vypocítame medián a IQR expresie kazdého génu
#medians<-apply(X,2, median)
#IQRs<-apply(X,2,IQR)
#X2<-apply(X,2,FUN=function(x) x-median)

med<-median(X)
iqr<-IQR(X)
X2<-(X-med)/iqr

# má gén ESR ako spojitý vplyv?
coxph(Surv(rfs, rfs.event)~X2[,"ESR1"])


# nájdime gény, ktoré majú významný vplyv

surv.res<-c()
for (i in c(1:ncol(X2)))
{
temp<-coxph(Surv(rfs, rfs.event)~X2[,i])
a<-summary(temp)
surv.res<-rbind(surv.res, c(colnames(X2)[i],a$coefficients))
}

colnames(surv.res)<-c("gene",colnames(a$coefficients))