R/graph_eigen.R

In baruuum/btoolbox: random toolbox

#' Eigen decomposition of matrices associated with graphs
#' 
#' @param x adjacency matrix (or affiliation matrix) of a graph; must be a \code{Matrix::sparseMatrix} and symmetric
#' @param k number of eigenparis to extract (largest algebraic values for adjacency matrix and smallest algebraic values for Laplacian matrix)
#' @param type string that is either \code{"adj", "adj_sym", "adj_left", "lap", "lap_sym", "lap_left"}, standing for, respectively, the adjacency matrix (\code{x}), the normalized (symmetric) adjacency matrix, the (left) normalized adjacency matrix, the Laplacian, the normalized (symmetric) Laplacian, or the (left) normalized Laplacian.
#' @param check string of ether "all" (default), "isolates," "symmetric," or "none." Checks whether there are dangling nodes or whether \code{x} is asymmetric (or both/none), and throws error if so
#' @param set_diag_zero if \code{TRUE}, the diagonals of \code{x} are set to zero before computing decomposition
#' @details Let \eqn{A} be the adjacency matrix associated with graph \eqn{G}. The  Laplacian is given as \deqn{L = D - A}, where \eqn{D} is a diagonal matrix containing the degrees. Let \eqn{I} be the identity matrix of same dimensions as \eqn{A}. Then the different versions of normalized adjacency and Laplacian matrices are defined as follows: 
#' 
#' Normalized (symmetric) adjacency matrix: \deqn{A_{s} = D^{-1/2}AD^{-1/2}}
#' (Left) normalized adjacency matrix: \deqn{A_{l} = D^{-1}A}
#' Normalized (symmetric) Laplacian: \deqn{L_{s} = I - A_{s}}
#' (Left) normalized Laplacian: \deqn{L_{l} = I - A_{l}}
#' @return returns the \code{k} largest/smallest eigenpairs of the requested matrix. Eigenvalues, and their corresponding eigenvectors, are ordered in decreasing order for adjacency matrices and in increasing order for Laplacian matrices.
#' @export
sp_eigen_adj = function(
    x, 
    k, 
    type = c("adj", "adj_sym", "adj_left", "lap", "lap_sym", "lap_left"),
    check = c("all", "symmetric", "isolates", "none"),
    set_diag_zero = FALSE
){
    
    type = match.arg(type, several.ok = FALSE)
    check = match.arg(check, several.ok = FALSE)
    
    check_iso = check %in% c("all", "isolates")
    check_sym = check %in% c("all", "symmetric")
    
    if (k %% 1 != 0) stop("k has to be in integer")
    
    #NOTE:
    # notice that the eigenvectors of the normalized adjacency matrix
    # and the normalized Laplacian are the same and only the eigenvalues differ
    e_pair = switch(
        type,
        adj      = .adj_eigen(x, k, check_iso, check_sym, set_diag_zero),
        adj_sym = .adj_eigen_norm(x, k, check_iso, check_sym, set_diag_zero),
        adj_left = .adj_eigen_left_norm(x, k, check_iso, check_sym, set_diag_zero),
        lap      = .lap_eigen(x, k, check_iso, check_sym, set_diag_zero),
        lap_sym = .adj_eigen_norm(x, k, check_iso, check_sym, set_diag_zero),
        lap_left = .adj_eigen_left_norm(x, k, check_iso, check_sym, set_diag_zero)
    )
    
    n = length(e_pair$e_val)
    list(
        values = if (type %in% c("adj", "adj_sym", "adj_left")) {
            
            as.numeric(e_pair$e_val)[n:1]
            
        } else if (type %in% c("lap_sym", "lap_left")) {
            
            1 - as.numeric(e_pair$e_val)[n:1]
            
        } else {
            
            as.numeric(e_pair$e_val)
            
        },
        vectors = if (type == "lap") {
            
            e_pair$e_vec[, 1:n, drop = FALSE]
            
        } else {
            
            e_pair$e_vec[, n:1, drop = FALSE]
            
        }
    )
    
}