R/Generation.R

In zhibinghe/Phub: modified EM algorithm for Penalized Hub Model

#' Generate a simple pattern adjacent matrix \eqn{A}
#'
#' @param n  number of nodes
#' @param n0 number of hub nodes
#' @param pi a n-dimensional vector containing probability of appearance in non-hub group  
#' @param pf \eqn{A_{ij}} if node j in the follower set of hub i
#' @param pr \eqn{A_{ij}} if node j does not belong to the follower set of hub i

#' @return a (n0+1)*n matrix, the first row is \code{pi}
#' @export
#'
#' @examples
#' pi = rep(0.05,100)
#' A0 = GenA(100,5,0.4,0.1,pi)
#'
GenA = function(n,n0,pf,pr,pi){
  A = matrix(NA,n0,n)
  follow = vector("list",n0)
  # follower set except hub itself
  for(i in 1:n0) follow[[i]] = seq(n0+i,n,n0)
  for(i in 1:n0){
    A[i,follow[[i]]] = pf
    A[i,-follow[[i]]] = pr
  }
  diag(A) = 1
  return(rbind(pi,A))
}
#' Generate grouped data \eqn{G}
#'
#' @param A (n0+1)*n-dimensional adjacent matrix
#' @param T number of observations (groups)
#' @param rho a (n0+1) dimensional vector containing hub weight
#' @return T*n matrix
#' @export
#'
#' @examples
#' n = 100; n0 = 5; T = 1000
#' A0 = GenA(n,n0,0.4,0.1,rep(0.05,100))
#' rho = c(0.2,rep(0.8/n0,n0))
#' G0 = GenG(A0,T,rho)
#'
GenG = function(A,T,rho){
  n0=dim(A)[1]; n=dim(A)[2]
  G = matrix(0,T,n)
  for(t in 1:T){
    hub = sample(0:(n0-1),1,prob=rho)
    G[t,]=rbinom(n,1,prob=A[hub+1,])
  }
  return(G)
}