R/CIBgraph.R

In alastair-JL/StochasticCIB: Stochastic CIB tools

#' CIBGraph
#'
#' Take your CIB transition matrix for a particular transition rule, and produce the corresponding graph.
#' @keywords CIB
#' @export
#' @param Transitions A stochastic transition matrix (rows are presumed to add to 1). Can be generated by any of the function referenced in \code{\link{CIBTransitionCalculators}}
#' @param weight How much weight do you wish to give each node; most often proportional to your forecast probability. If no value is given, all nodes will be the same size.
#' @param target All nodes are going to be colour coded. Based on which factor do you want to colour code them? This value takes on integer values.
#' @param cutoffMultiplier For the sake of reasonable visualization, CIBGraph will, by default, ignore all edges with probability less than 10^-4. If you wish to include all edges, set this value to 0. If you wish to neglect more edges, set this cut off higher, for example 10^-2.
#' @param layout CIBgraph uses the package igraph internally in order to draw your markov chain. By default we use "layout_with_dh" to layout the nodes. If you'd prefer a different layout, select one here. See the R \link{igraph} manual for details.
#' @return NULL. Function generates a picture, but has not return type.
#' @author Alastair Jamieson Lane. <aja107@@math.ubc.ca>
#' @examples
#' data(ExampleCIBdata)
#' boltzTrans<-LocalBoltzmann(ExampleCIBdata)
#' forecast<- CIBforecast(boltzTrans)
#' CIBgraph(boltzTrans,forecast[[1]])
#' CIBgraph(boltzTrans,forecast[[1]],2)
#' CIBgraph(boltzTrans,forecast[[1]],2,0)
#' CIBgraph(boltzTrans,forecast[[1]],2,10^-4,layout.grid)


CIBgraph<-function(Transitions,weight=NULL,target=1,cutoffMultiplier=10^-4,layout=igraph::layout_with_dh){
  network <- igraph::graph_from_adjacency_matrix(Transitions , mode='directed', diag=T,weighted=T)
  
  if(nrow(Transitions)!=ncol(Transitions)){
    error("Transition matrix must be square. What are you doing?")
  }
  
  if(is.null(weight)|NROW(weight)!=NROW(Transitions)){
    weight=rep(15,nrow(Transitions));
  }else{
    weight<-log(weight+10^-10);
    weight<- -sqrt(max(weight)-weight)
    weight<- 16*weight/max(abs(weight))+26
    weight<-as.numeric(t(weight))
  }
  
  names=colnames(Transitions)
  carac=c(substr(names,target,target))
  
  nodes <- data.frame(names,carac,weight,stringsAsFactors =F)
  
  linkRow= c()
  linkCol= c()
  linkThk= c()
  linkHex=c()

  grays<- gray.colors(50, start = 0.75, end = 0.25)
  
  
      cutOff<-sum(Transitions)/ncol(Transitions)*cutoffMultiplier
      for(iii in 1:nrow(Transitions)){
        for(jjj in 1:nrow(Transitions)){
          if(Transitions[iii,jjj]>cutOff){
            linkRow<-c(linkRow,nodes$names[iii])
            linkCol<-c(linkCol,nodes$names[jjj])
            linkThk<-c(linkThk,3-2.6*exp(-4*Transitions[iii,jjj]/max(Transitions)))
            linkHex<-c(linkHex,grays[floor(49.9*(Transitions[iii,jjj]/max(Transitions)))+1])
          }
        } 
      }
    
  links<- data.frame(linkRow,linkCol,linkThk,stringsAsFactors =F)
  network <- igraph::graph_from_data_frame(d=links, vertices=nodes, directed=T) 
  
  igraph::layout_nicely(network, dim = 2)
  
  coul  <- c("#88CC44","#4488CC","#CC4488","#CC8844","#44CC88","#8844CC","#22EE22","#2222EE","#EE2222")
  # Create a vector of color
  my_color <- coul[as.numeric(as.factor(igraph::V(network)$carac))]

  plot(network,vertex.size=weight, vertex.color=my_color,
       vertex.label.cex=weight/27, vertex.label.family="Ariel",
       edge.width=linkThk, edge.curved=0.35,    edge.arrow.size=0.8,                            # Arrow size, defaults to 1
       edge.arrow.width=0.8, edge.color=linkHex,loop.angle=0.85,layout=layout )
}