1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
| f_getSA_1 <- function(subjectID){
res = list()
res$subjectID <- subjectID
res$cases <- getCaseLists(mycgds, subjectID)
res$GP <- getGeneticProfiles(mycgds,subjectID)
View(res$GP)
View(res$cases)
res
}
f_getSA_5 <- function(lc_a, lc_sGP, lc_cg,lc_cc=1){
lc_a$cc <- lc_a$cases[lc_cc,1]
lc_a$sGP <- lc_a$GP[lc_sGP,1]
lc_a$cg <- unlist(lc_cg)
lc_a
}
f_getSA_2 <- function(lc_a, lc_c=c('OS_MONTHS','OS_STATUS')){
lc_a$cd <- getClinicalData(mycgds, lc_a$cc)
View(head(lc_a$cd))
View(dim(lc_a$cd))
View(colnames(lc_a$cd))
lc_a$PD <- getProfileData(mycgds,lc_a$cg,lc_a$sGP,lc_a$cc)
View(head(lc_a$PD))
View(dim(lc_a$PD))
if(all(lc_c %in% colnames(lc_a$cd))){
lc_a$os <- lc_a$cd[,lc_c]
lc_a$os <- subset(lc_a$os, lc_a$os[[lc_c[1]]]>0)
print(unique(lc_a$os[[2]]))
}else{
print(lc_c %in% colnames(lc_a$cd))
}
lc_a
}
f_getSA_4 <- function(lc_b){
p.val = 1 - pchisq(lc_b$chisq, length(lc_b$n) - 1)
p.val
}
f_getSA_3 <- function(lc_a, lc_g, lc_ST = '1:DECEASED', lc_c=c('OS_MONTHS','OS_STATUS'), lc_force = F, lc_d=F){
dat2 <- subset(lc_a$PD, !is.na(lc_a$PD[[lc_g]]))
common <- intersect(rownames(lc_a$os), rownames(dat2))
lc_a$dat <- lc_a$os[common,]
lc_a$dat2 <- dat2[common,]
res <- lc_a$dat2[[lc_g]]
if(lc_d){
res <- data.frame(res)
}else{
res <- data.frame(ifelse(res > median(res),'high','low'))
}
if(length(unique(res[[1]]))<2){
lc_a$pv <- 1
return(lc_a)
}
colnames(res) <- lc_g
nmsl <<- res
nm_my.surv <<- Surv(lc_a$dat[[lc_c[1]]],lc_a$dat[[lc_c[2]]]==lc_ST)
nmb <<- paste("nm_my.surv~", lc_g)
lc_b <- survdiff(nm_my.surv~., data = res)
lc_a$pv <- f_getSA_4(lc_b)
if(lc_a$pv < 0.05 lc_force){
kmfit1 <- survfit(formula(nmb), data = nmsl)
lc_a$pr <- ggsurvplot(kmfit1,conf.int =F, pval = T,risk.table =T, ncensor.plot = TRUE)
}
lc_a
}
library(fdrtool)
f_getSA_6 <- function(lc_a, ...){
res1 = NULL
res2 = list()
for(hhh in colnames(lc_a$PD)){
lc_a <- f_getSA_3(lc_a, hhh, ...)
res1 <- rbind(res1, c(hhh, lc_a$pv))
if(lc_a$pv < 0.05){
res2[[hhh]] <- lc_a$pr
}
}
if(length(res1) > 0){
colnames(res1) <- c('gene','p')
lc_a$srl <- list(pv=data.frame(res1), pr=res2)
lc_a$srl$pv$p <- as.numeric(lc_a$srl$pv$p)
lc_a$srl$fdr <- fdrtool(lc_a$srl$pv$p, statistic="pvalue")
lc_a$srl$pv$q <- lc_a$srl$fdr$qval
lc_a$srl$s <- subset(lc_a$srl$pv, p <0.05 & q <0.05)
}else{
print('error: NULL lc_a$srl')
}
lc_a
}
|
