New effects for REM: first- and second-order network coupling, balance

[sophiyaberdiyeva]

Hello! I'm writing a bachelor thesis about inetrorganizational networks, so I use goldfish package to create some tie-oriented models on time-stamped data. According to my theoretical model, I needed some additional effects for formula: (1) first- and (2) second-order network coupling as they were described by B. Uzzi with reference to W. Baker (Uzzi, B. (1996). The Sources and Consequences of Embeddedness for the Economic Performance of Organizations: The Network Effect. American Sociological Review, 61. https://doi.org/10.2307/2096399); and some (3) "balance" metric which is simply the difference between the number of interactions from i to j and vice versa. I tried to write the code using existing functions for effects as references and checked if it could be run without errors inside estimate function, and would like to leave it here to developers' consideration because of my lack of experience with Git contributions.

#  first_embeddedness-------------------------------------------------------------------
#' init stat matrix first_embeddedness using cache
#'
#' @param effectFun function with additional parameters weighted, isTwoMode,
#'   transformFun, type
#' @param network matrix n1*n2
#' @param window NULL|numeric size of the window
#' @param n1 integer nrow(network)
#' @param n2 integer ncol(network)
#'
#' @return list with named components: cache numeric vector size n2,
#'   stat matrix numeric n1*n2
#' @noRd
#'
#' @examples
#' \donttest{
#' network <- matrix(
#'   c(
#'     0, 0, 0, 1, 0, 0,
#'     0, 0, 0, 0, 0, 0,
#'     0, 2, 0, 0, 0, 3,
#'     1, 0, 0, 0, 0, 4,
#'     1, 2, 0, 0, 0, 0
#'   ),
#'   nrow = 5, ncol = 6, byrow = TRUE
#' )
#' effectFUN <- function(
    #'     weighted = TRUE, isTwoMode = TRUE, transformFun = identity) {
#'   NULL
#' }
#' init_REM_choice.first_embeddedness(effectFUN, network, NULL, 5, 6)
#' network <- matrix(
#'   c(
#'     0, 0, 0, 1, 0,
#'     0, 0, 0, 0, 0,
#'     0, 2, 0, 0, 3,
#'     1, 0, 0, 0, 4,
#'     1, 2, 0, 0, 0
#'   ),
#'   nrow = 5, ncol = 5, byrow = TRUE
#' )
#' effectFUN <- function(
    #'     weighted = TRUE, isTwoMode = FALSE, transformFun = identity,
#'     type = "ego") {
#'   NULL
#' }
#' init_REM_choice.first_embeddedness(effectFUN, network, NULL, 5, 5)
#'
#' effectFUN <- function(
    #'     weighted = TRUE, isTwoMode = FALSE, transformFun = identity,
#'     type = "alter") {
#'   NULL
#' }
#' init_REM_choice.first_embeddedness(effectFUN, network, NULL, 5, 5)
#' }
init_REM_choice.first_embeddedness <- function(effectFun, network, window, n1, n2, ...) {
  # Get arguments
  params <- formals(effectFun)
  weighted <- eval(params[["weighted"]])
  isTwoMode <- eval(params[["isTwoMode"]])
  funApply <- eval(params[["transformFun"]])
  type <- eval(params[["type"]])
  
  if (isTwoMode && type == "ego") {
    stop(dQuote("first_embeddedness"),
         "effect must not use for type 'ego' (type = 'ego') when is ",
         "a two-mode network (isTwoMode = TRUE) ",
         call. = FALSE
    )
  }
  # has window or is empty initialize empty
  if ((!is.null(window) && !is.infinite(window)) || all(network == 0)) {
    return(list(
      cache = numeric(n2),
      stat = matrix(forceAndCall(1, funApply, 0), nrow = n1, ncol = n2)
    ))
  }
  
  # Convert matrix rows to relative values (percentage by row)
  row_percentages <- t(apply(network, 1, function(row) row/sum(row)))
  
  # Calculate the sum of the squares of the elements in each line
  cache <- rowSums(row_percentages^2)
  cache[is.na(cache)] <- 0
  
  # applied transformFun instead
  stat <- forceAndCall(1, funApply, cache)
  
  # return expected n1 * n2 matrix
  # update check that type should be one of "alter" or "ego"
  byRow <- type == "alter"
  stat <- matrix(stat, nrow = n1, ncol = n2, byrow = byRow)
  #
  if (!isTwoMode) diag(stat) <- 0
  
  return(list(cache = cache, stat = stat))
}

#' update stat first_embeddedness using cache
#'
#' @param network matrix n1*n2
#' @param sender integer
#' @param receiver integer
#' @param replace numeric
#' @param cache numeric vector size n2
#' @param n1 integer nrow(network)
#' @param n2 integer ncol(network)
#' @param isTwoMode logical
#' @param weighted logical
#' @param transformFun function to apply to the stat
#' @param type character should be 'alter' or 'ego' type
#'
#' @return list:
#'   cache numeric vector size n2,
#'   changes NULL || array cbind(node1 = x, node2 = y, replace = z) stat updates
#' @noRd
#'
#' @examples
#' \donttest{
#' network <- matrix(
#'   c(
#'     0, 0, 0, 1, 0,
#'     0, 0, 0, 0, 0,
#'     0, 2, 0, 0, 3,
#'     1, 0, 0, 0, 4,
#'     1, 2, 0, 0, 0
#'   ),
#'   nrow = 5, ncol = 5, byrow = TRUE
#' )
#' cache <- c(2, 7, 0, 1, 0, 7)
#' update_REM_choice_first_embeddedness(
#'   network,
#'   1, 2, 3,
#'   cache, 5, 6,
#'   isTwoMode = TRUE, weighted = TRUE, transformFun = sqrt, type = "ego"
#' )
#'
#' update_REM_choice_first_embeddedness(
#'   network,
#'   1, 2, 3,
#'   cache, 5, 6,
#'   isTwoMode = TRUE, weighted = TRUE, transformFun = sqrt, type = "alter"
#' )
#' }
update_REM_choice_first_embeddedness <- function(
    network,
    sender, receiver, replace, cache,
    n1, n2, isTwoMode = FALSE,
    weighted = FALSE, transformFun = identity,
    type = c("alter", "ego")) {
  type <- match.arg(type)
  # init res: changes NULL if not changes
  res <- list(cache = NULL, changes = NULL)
  
  # Get old value
  oldValue <- network[sender, receiver]
  
  # change weighted
  if (!weighted) {
    oldValue <- sign(oldValue)
    replace <- sign(replace)
  }
  
  # Check if old value has changed
  if (is.na(oldValue) && is.na(replace)) {
    return(res)
  } else if (!is.na(oldValue) && !is.na(replace) && oldValue == replace) {
    return(res)
  }
  
  if (is.na(oldValue)) oldValue <- 0
  if (is.na(replace)) replace <- 0
  # update cache for sender, event/tie
  
  new_network <- network
  new_network[sender, receiver]<- replace
  
  # Convert matrix rows to relative values (percentage by row)
  row_percentages <- t(apply(new_network, 1, function(row) row/sum(row)))
  
  # Calculate the sum of the squares of the elements in each line
  row_sums <- rowSums(row_percentages^2)
  row_sums[is.na(row_sums)] <- 0
  
  cache[sender] <-row_sums[sender]
    
  if (type == "alter") {
    # when the cache had change
    if (isTwoMode) {
      other <- seq_len(n1)
    } else {
      other <- setdiff(seq_len(n1), sender)
    }
    
    changes <- cbind(
      node1 = other,
      node2 = sender,
      replace = forceAndCall(1, transformFun, cache[sender])
    )
  } else if (type == "ego") {
    # when the cache had change
    if (isTwoMode) {
      other <- seq_len(n2)
    } else {
      other <- setdiff(seq_len(n2), sender)
    }
    
    changes <- cbind(
      node1 = sender,
      node2 = other,
      replace = forceAndCall(1, transformFun, cache[sender])
    )
  }
  return(list(cache = cache, changes = changes))
}

#second_embeddedness-------------------------------------------------------------------
#' init stat matrix second_embeddedness using cache
#'
#' @param effectFun function with additional parameters weighted, isTwoMode,
#'   transformFun, type
#' @param network matrix n1*n2
#' @param window NULL|numeric size of the window
#' @param n1 integer nrow(network)
#' @param n2 integer ncol(network)
#'
#' @return list with named components: cache numeric vector size n2,
#'   stat matrix numeric n1*n2
#' @noRd
#'
#' @examples
#' \donttest{
#' network <- matrix(
#'   c(
#'     0, 0, 0, 1, 0, 0,
#'     0, 0, 0, 0, 0, 0,
#'     0, 2, 0, 0, 0, 3,
#'     1, 0, 0, 0, 0, 4,
#'     1, 2, 0, 0, 0, 0
#'   ),
#'   nrow = 5, ncol = 6, byrow = TRUE
#' )
#' effectFUN <- function(
    #'     weighted = TRUE, isTwoMode = TRUE, transformFun = identity) {
#'   NULL
#' }
#' init_REM_choice.second_embeddedness(effectFUN, network, NULL, 5, 6)
#' network <- matrix(
#'   c(
#'     0, 0, 0, 1, 0,
#'     0, 0, 0, 0, 0,
#'     0, 2, 0, 0, 3,
#'     1, 0, 0, 0, 4,
#'     1, 2, 0, 0, 0
#'   ),
#'   nrow = 5, ncol = 5, byrow = TRUE
#' )
#' effectFUN <- function(
    #'     weighted = TRUE, isTwoMode = FALSE, transformFun = identity,
#'     type = "ego") {
#'   NULL
#' }
#' init_REM_choice.second_embeddedness(effectFUN, network, NULL, 5, 5)
#'
#' effectFUN <- function(
    #'     weighted = TRUE, isTwoMode = FALSE, transformFun = identity,
#'     type = "alter") {
#'   NULL
#' }
#' init_REM_choice.second_embeddedness(effectFUN, network, NULL, 5, 5)
#' }
init_REM_choice.second_embeddedness <- function(effectFun, network, window, n1, n2, ...) {
  # Get arguments
  params <- formals(effectFun)
  weighted <- eval(params[["weighted"]])
  isTwoMode <- eval(params[["isTwoMode"]])
  funApply <- eval(params[["transformFun"]])
  type <- eval(params[["type"]])
  
  if (isTwoMode && type == "ego") {
    stop(dQuote("second_embeddedness"),
         "effect must not use for type 'ego' (type = 'ego') when is ",
         "a two-mode network (isTwoMode = TRUE) ",
         call. = FALSE
    )
  }
  # has window or is empty initialize empty
  if ((!is.null(window) && !is.infinite(window)) || all(network == 0)) {
    return(list(
      cache = numeric(n2),
      stat = matrix(forceAndCall(1, funApply, 0), nrow = n1, ncol = n2)
    ))
  }
  
  # Convert matrix rows to relative values (percentage by row)
  row_percentages <- t(apply(network, 1, function(row) row/sum(row)))
  
  # Calculate the sum of the squares of the elements in each line
  rowsums <- rowSums(row_percentages^2)
  rowsums[is.na(rowsums)] <- 0
  
  
  num_vertices <- nrow(network)  # Number of vertices in the graph
  metrics <- numeric(num_vertices)  # Initialize the vector to store metrics
  
  for (i in 1:num_vertices) {
    neighbors <- which(m[,i ] != 0)  # Vertices from which edges go to vertex i
    
    # Sum the values for all connected vertices
    sum_Q <- sum(rowsums[neighbors])
    
    # Calculate the metric for vertex i
    metric <- sum_Q / length(neighbors)
    
    # Save the metric value
    metrics[i] <- metric
  }
  metrics[is.na(metrics)] <- 0
  cache <- metrics
  
  
  # applied transformFun instead
  stat <- forceAndCall(1, funApply, cache)
  
  # return expected n1 * n2 matrix
  # update check that type should be one of "alter" or "ego"
  byRow <- type == "alter"
  stat <- matrix(stat, nrow = n1, ncol = n2, byrow = byRow)
  #
  if (!isTwoMode) diag(stat) <- 0
  
  return(list(cache = cache, stat = stat))
}

#' update stat second_embeddedness using cache
#'
#' @param network matrix n1*n2
#' @param sender integer
#' @param receiver integer
#' @param replace numeric
#' @param cache numeric vector size n2
#' @param n1 integer nrow(network)
#' @param n2 integer ncol(network)
#' @param isTwoMode logical
#' @param weighted logical
#' @param transformFun function to apply to the stat
#' @param type character should be 'alter' or 'ego' type
#'
#' @return list:
#'   cache numeric vector size n2,
#'   changes NULL || array cbind(node1 = x, node2 = y, replace = z) stat updates
#' @noRd
#'
#' @examples
#' \donttest{
#' network <- matrix(
#'   c(
#'     0, 0, 0, 1, 0, 0,
#'     0, 0, 0, 0, 0, 0,
#'     0, 2, 0, 0, 0, 3,
#'     1, 0, 0, 0, 0, 4,
#'     1, 2, 0, 0, 0, 0
#'   ),
#'   nrow = 5, ncol = 6, byrow = TRUE
#' )
#' cache <- c(2, 7, 0, 1, 0, 7)
#' update_REM_choice_second_embeddedness(
#'   network,
#'   1, 2, 3,
#'   cache, 5, 6,
#'   isTwoMode = TRUE, weighted = TRUE, transformFun = sqrt, type = "ego"
#' )
#'
#' update_REM_choice_second_embeddedness(
#'   network,
#'   1, 2, 3,
#'   cache, 5, 6,
#'   isTwoMode = TRUE, weighted = TRUE, transformFun = sqrt, type = "alter"
#' )
#' }
update_REM_choice_second_embeddedness <- function(
    network,
    sender, receiver, replace, cache,
    n1, n2, isTwoMode = FALSE,
    weighted = FALSE, transformFun = identity,
    type = c("alter", "ego")) {
  type <- match.arg(type)
  # init res: changes NULL if not changes
  res <- list(cache = NULL, changes = NULL)
  
  # Get old value
  oldValue <- network[sender, receiver]
  
  # change weighted
  if (!weighted) {
    oldValue <- sign(oldValue)
    replace <- sign(replace)
  }
  
  # Check if old value has changed
  if (is.na(oldValue) && is.na(replace)) {
    return(res)
  } else if (!is.na(oldValue) && !is.na(replace) && oldValue == replace) {
    return(res)
  }
  
  if (is.na(oldValue)) oldValue <- 0
  if (is.na(replace)) replace <- 0
  # update cache for sender, event/tie
  
  new_network <- network
  new_network[sender, receiver]<-replace
  
  # Convert matrix rows to relative values (percentage by row)
  row_percentages <- t(apply(new_network, 1, function(row) row/sum(row)))
  
  # Calculate the sum of the squares of the elements in each line
  rowsums <- rowSums(row_percentages^2)
  rowsums[is.na(rowsums)] <- 0
  

  num_vertices <- nrow(new_network)  # Number of vertices in the graph
  metrics <- numeric(num_vertices)  # Initialize the vector to store metrics
  
  for (i in 1:num_vertices) {
    neighbors <- which(new_network[,i ] != 0)  # Vertices from which edges go to vertex i 
    
    # Sum the values for all connected vertices
    sum_Q <- sum(rowsums[neighbors])
    
    # Calculate the metric for vertex i
    metric <- sum_Q / length(neighbors)
    
    # Save the metric value
    metrics[i] <- metric
  }
  metrics[is.na(metrics)] <- 0
  cache <- metrics
  
  
  if (type == "alter") {
    # when the cache had change
    if (isTwoMode) {
      other <- seq_len(n1)
    } else {
      other <- setdiff(seq_len(n1), sender)
    }
    
    changes <- cbind(
      node1 = other,
      node2 = sender,
      replace = forceAndCall(1, transformFun, cache[sender])
    )
  } else if (type == "ego") {
    # when the cache had change
    if (isTwoMode) {
      other <- seq_len(n2)
    } else {
      other <- setdiff(seq_len(n2), sender)
    }
    
    changes <- cbind(
      node1 = sender,
      node2 = other,
      replace = forceAndCall(1, transformFun, cache[sender])
    )
  }
  return(list(cache = cache, changes = changes))
}

# balance ---------------------------------------------------------------------
#' init stat matrix balance
#' 
#' @param effectFun function with additional parameters weighted, transformFun
#' @param network matrix n1*n2
#' @param window NULL|numeric size of the window
#' @param n1 integer nrow(network)
#' @param n2 integer ncol(network)
#'
#' @return list: stat matrix numeric n1*n2
#' @noRd
#'
#' @examples
#' \donttest{
#' network <- matrix(
#'   c(
#'     0, 0, 0, 1, 0, 0,
#'     0, 0, 0, 0, 0, 0,
#'     0, 2, 0, 0, 0, 3,
#'     1, 0, 0, 0, 0, 4,
#'     1, 2, 0, 0, 0, 0
#'   ),
#'   nrow = 5, ncol = 6, byrow = TRUE
#' )
#' effectFUN <- function(weighted = TRUE, transformFun = identity) {
#'   NULL
#' }
#' init_DyNAM_choice.balance(effectFun,network)
#' }
init_DyNAM_choice.balance <- function(effectFun, network, n1, n2, ...) {
  # get arguments
  params <- formals(effectFun)
  weighted <- eval(params[["weighted"]])
  funApply <- eval(params[["transformFun"]])
  
  diff_network <- network - t(network)
  return(list(stat = unname(diff_network)))
}

#' update stat balance using cache
#'
#' @param network matrix n1*n2
#' @param sender integer
#' @param receiver integer
#' @param replace numeric
#' @param weighted logical
#' @param transformFun function to apply to the stat
#'
#' @return list:
#'   changes NULL || array cbind(node1 = x, node2 = y, replace = z) stat updates
#' @noRd
#' @aliases balance
#'
#' @examples
#' \donttest{
#' network <- matrix(
#'   c(
#'     0, 0, 0, 1, 0, 0,
#'     0, 0, 0, 0, 0, 0,
#'     0, 2, 0, 0, 0, 3,
#'     1, 0, 0, 0, 0, 4,
#'     1, 2, 0, 0, 0, 0
#'   ),
#'   nrow = 5, ncol = 6, byrow = TRUE
#' )
#' update_DyNAM_choice_balance(network,
#'   1, 2, 3,
#'   weighted = TRUE, transformFun = sqrt
#' )
#' }
update_DyNAM_choice_balance <- function(
    network,
    sender, receiver, replace,
    weighted = FALSE, transformFun = identity) {
  # No change check, irrelevant for two-mode network
  # if(sender == receiver) return(NULL)
  
  # init res
  res <- list(changes = NULL)
  
  # Get old value
  oldValue <- network[sender, receiver]
  
  # change for weighted effect
  if (!weighted) {
    oldValue <- sign(oldValue)
    replace <- sign(replace)
  }
  
  # Check if old value has changed
  if (is.na(oldValue) && is.na(replace)) {
    return(res)
  } else if (!is.na(oldValue) && !is.na(replace) && oldValue == replace) {
    return(res)
  }
  
  if (is.na(oldValue)) oldValue <- 0
  if (is.na(replace)) replace <- 0
  
  new_network <- network
  new_network[sender, receiver]<-replace
  diff_network <- new_network - t(new_network)
  
  # change stat
  res$changes <- cbind(
    node1 = sender,
    node2 = receiver,
    replace = if (!weighted) {
      1 * (replace > 0)
    } else {
      diff_network[sender, receiver]
    }
  )
  
  return(res)
}

# balance ---------------------------------------------------------------------
#' init stat matrix balance
#'
#' @param effectFun function with additional parameters weighted, transformFun
#' @param network matrix n1*n2
#' @param window NULL|numeric size of the window
#' @param n1 integer nrow(network)
#' @param n2 integer ncol(network)
#'
#' @return list: stat matrix numeric n1*n2
#' @noRd
#'
#' @examples
#' \donttest{
#' network <- matrix(
#'   c(
#'     0, 0, 0, 1, 0, 0,
#'     0, 0, 0, 0, 0, 0,
#'     0, 2, 0, 0, 0, 3,
#'     1, 0, 0, 0, 0, 4,
#'     1, 2, 0, 0, 0, 0
#'   ),
#'   nrow = 5, ncol = 6, byrow = TRUE
#' )
#' effectFUN <- function(weighted = TRUE, transformFun = identity) {
#'   NULL
#' }
#' init_REM_choice.balance(effectFUN, network)
#' }
init_REM_choice.balance <- function(effectFun, network, window, n1, n2, ...) {
  init_DyNAM_choice.balance(
    effectFun = effectFun, network = network,
    window = window, n1 = n1, n2 = n2, ...
  )
}



#' update stat balance using cache
#'
#' @param network matrix n1*n2
#' @param sender integer
#' @param receiver integer
#' @param replace numeric
#' @param weighted logical
#' @param transformFun function to apply to the stat
#'
#' @return list:
#'   changes NULL || array cbind(node1 = x, node2 = y, replace = z) stat updates
#' @noRd
#' @aliases balance
#'
#' @examples
#' \donttest{
#' network <- matrix(
#'   c(
#'     0, 0, 0, 1, 0, 0,
#'     0, 0, 0, 0, 0, 0,
#'     0, 2, 0, 0, 0, 3,
#'     1, 0, 0, 0, 0, 4,
#'     1, 2, 0, 0, 0, 0
#'   ),
#'   nrow = 5, ncol = 6, byrow = TRUE
#' )
#' update_RAM_choice_balance(network,
#'   1, 2, 3,
#'   weighted = TRUE, transformFun = sqrt
#' )
#' }
update_REM_choice_balance <- function(
    network,
    sender, receiver, replace,
    weighted = FALSE, transformFun = identity) {
  update_DyNAM_choice_balance(
    network = network,
    sender = sender, receiver = receiver, replace = replace,
    weighted = weighted, transformFun = transformFun
  )
}