R/path_compare.R

path_compare<-function(list){
  profile1=list[[1]]
  profile2=list[[2]]
  pathway=colnames(profile1)
  data<-list(model_sham=profile1,quan_model=profile2)
  type<-colnames(data$quan_model)
  mcoin<-mcia(data)
  axis<-mcoin$mcoa$Tli
  # plot(mcoin,df.color=1:2, phenovec=type,sample.lab = FALSE,sample.legend = TRUE)

  distance<-function(co1,co2){
    matrix<-rbind(co1,co2)
    dis<-dist(matrix)
    return(dis)
  }

  ####### generate pathway-gene profile

  random_profile<-function(real_profile){
    random_profile<-matrix(nrow=nrow(real_profile),ncol=ncol(real_profile))
    colnames(random_profile)=colnames(real_profile)
    rownames(random_profile)=rownames(real_profile)
    exvector<-c()
    for(i in 1:nrow(real_profile)){
      exvector<-c(exvector,as.numeric(as.character(real_profile[i,])))
    }
    ##### sample gene expression
    for(i in 1:nrow(real_profile)){
      for(j in 1:ncol(real_profile)){
        random_profile[i,j]=sample(exvector,1,replace=TRUE)
      }
    }
    del_index=c()
    for(i in 1:nrow(real_profile)){
      if(sum(random_profile[i,])==0){del_index<-c(del_index,i)}
    }
    random_profile=random_profile[-del_index,]
    return(random_profile)
  }
  #############


  ############ project by mcia
  mcia_profile<-function(profile1_random,profile2_random){
    data_random<-list(x=profile1_random,y=profile2_random)
    type_random<-colnames(data_random$x)
    mcoin_random<-mcia(data_random)
    axis_random<-mcoin_random$mcoa$Tli
    #####  compute distance between pathways
    dis<-c()
    for(i in 1:ncol(profile1_random)){
      strcon1=paste(pathway[i],paste('.','x',sep=''),sep='')
      strcon2=paste(pathway[i],paste('.','y',sep=''),sep='')
      index1=c(rownames(axis_random)==strcon1)
      coordinate_con1<-axis_random[index1,]
      index2=c(rownames(axis_random)==strcon2)
      coordinate_con2<-axis_random[index2,]
      dis_con1_con2<-distance(coordinate_con1,coordinate_con2)
      dis<-c(dis,dis_con1_con2)
    }
    return(dis)
  }

  ####### function random 999 times
  dis_all_profile<-data.frame()
  for(i in 1:999){
    profile1_random<-random_profile(profile1)
    profile2_random<-random_profile(profile2)
    dis_quan<-mcia_profile(pathway,profile1_random,profile2_random)
    dis_all_profile<-rbind(dis_all_profile,dis_quan)
  }

  #########  add real condition distance
  addrealdis<-function(dataframe){
    dis_real<-mcia_profile(profile1,profile2)
    thousand_dis<-rbind(dataframe,dis_real)
    colnames(thousand_dis)=colnames(profile1)
    return(thousand_dis)
  }
  final_dis<-addrealdis(dis_all_profile)

  compute_sigpathway<-function(dataframe){
    pvalue<-c()
    for(i in 1:ncol(dataframe)){
      sig=nrow(dataframe[dataframe[,i]>dataframe[1000,i],])
      pvalue<-c(pvalue,sig/1000)
    }
    return(pvalue)
  }
  getsigpathway<-function(matrix_dis){
    pvalue<-compute_sigpathway(matrix_dis)
    sig<-pvalue<0.05
    sig_pathway<-rbind(pvalue,sig)
    colnames(sig_pathway)<-colnames(profile1)
    rownames(sig_pathway)<-c('pvalue','significant')
    return(sig_pathway)

  }
  sig_quan<-getsigpathway(final_dis)

  return(sig_quan)
}
github-gs/Quantitative-pathway-analysis documentation built on May 13, 2019, 12:18 p.m.