Update package version

gbertagnolli · May 13, 2024 · 344984e · 344984e
1 parent 29603c8
commit 344984e
Show file tree

Hide file tree

Showing 37 changed files with 1,377 additions and 2,769 deletions.
diff --git a/DESCRIPTION b/DESCRIPTION
@@ -1,29 +1,29 @@
 Package: diffudist
 Title: Diffusion Distance for Complex Networks
-Version: 1.0.0
+Version: 1.1
 Authors@R: c(
     person(given = "Giulia",
            family = "Bertagnolli",
            role = c("aut", "cre"),
            email = "giulia.bertagnolli@gmail.com",
            comment = c(ORCID = "0000-0001-8637-0632")),
-    person(given = "Manlio", 
-           family = "De Domenico", 
-           role = "aut", 
-           email = "manlio.dedomenico@gmail.com", 
+    person(given = "Manlio",
+           family = "De Domenico",
+           role = "aut",
+           email = "manlio.dedomenico@gmail.com",
            comment = c(ORCID = "0000-0001-5158-8594"))
            )
 URL: https://gbertagnolli.github.io/diffudist/
 BugReports: https://github.com/gbertagnolli/diffudist/issues/
 Description: Enables the evaluation of diffusion distances for complex single-layer networks.
     Given a network one can define different types of Laplacian (or transition)
-    matrices corresponding to different continuous-time random walks dynamics on the network. 
-    This package enables the evaluation of Laplacians, stochastic matrices, and the 
+    matrices corresponding to different continuous-time random walks dynamics on the network.
+    This package enables the evaluation of Laplacians, stochastic matrices, and the
     corresponding diffusion distance matrices. The metric structure induced by the network-driven
-    process is richer and more robust than the one given by shortest-paths and allows to study 
-    the geometry induced by different types of diffusion-like communication mechanisms taking 
-    place on complex networks. 
-    For more details see: De Domenico, M. (2017) <doi:10.1103/physrevlett.118.168301> and 
+    process is richer and more robust than the one given by shortest-paths and allows to study
+    the geometry induced by different types of diffusion-like communication mechanisms taking
+    place on complex networks.
+    For more details see: De Domenico, M. (2017) <doi:10.1103/physrevlett.118.168301> and
     Bertagnolli, G. and De Domenico, M. (2021) <doi:10.1103/PhysRevE.103.042301>.
 Depends: R (>= 3.5.0)
 Imports:
@@ -34,7 +34,8 @@ Imports:
     igraph,
     Matrix,
     stats,
-    Rcpp (>= 1.0.7),
+    RColorBrewer,
+    Rcpp (>= 1.0.10),
     reshape2,
     rlang,
     viridis
@@ -49,5 +50,5 @@ Suggests:
 License: GPL (>=2)
 Encoding: UTF-8
 LazyData: true
-RoxygenNote: 7.1.2
+RoxygenNote: 7.2.3
 VignetteBuilder: knitr
diff --git a/R/get_diff_prob_matrix.R b/R/get_diff_prob_matrix.R
@@ -52,6 +52,15 @@ get_laplacian <- function(g,  type = "Laplacian", weights = NULL, verbose = TRUE
       if (verbose) {
         cat("Weighted network. Using edge weights as connection strenghts.\n")
       }
+      if (!is.numeric(igraph::E(g)$weight)) {
+        cat("weight edge attribute is not numeric. Ignore weights.\n")
+        igraph::E(g)$weight <- 1
+      } else {
+        if (!all(E(g)$weight >= 0)) {
+          cat("Negative weights. Absolute value will be used.\n")
+          igraph::E(g)$weight <- abs(igraph::E(g)$weight)
+        }
+      }
     } else {
       # no (numeric) weight edge attribute
       if (verbose) {
@@ -207,7 +216,7 @@ get_diffusion_probability_matrix <- function(g, tau, type = "Normalized Laplacia
     )
   } else {
     tryCatch(
-      type <- type <- match.arg(toupper(type), types),
+      type <- match.arg(toupper(type), types),
       error = function(e) {
         cat(
           "ERROR! Wrong type of Laplacian, available types are:\n",

diff --git a/R/get_distance_matrix.R b/R/get_distance_matrix.R
@@ -5,8 +5,8 @@
 #'
 #' The diffusion distance at time \eqn{\tau} between nodes \eqn{i, j \in G}
 #' is defined as
-#' \deqn{D_{\tau}(i, j) = \vert \mathbf{p}(t|i) - \mathbf{p}(t|j) \vert_2}
-#' with \eqn{\mathbf{p}(t|i) = (e^{- \tau L})_{i\cdot} = \mathbf{e}_i e^{- \tau L}}
+#' \deqn{D_{\tau}(i, j) = \vert \mathbf{p}(t|i) - \mathbf{p}(t|j) \vert_2} with
+#' \eqn{\mathbf{p}(t|i) = (e^{- \tau L})_{i\cdot} = \mathbf{e}_i e^{- \tau L}}
 #' indicating the i-th row of the stochastic matrix \eqn{e^{- \tau L}} and
 #' representing the probability (row) vector of a random walk dynamics
 #' corresponding to the initial condition \eqn{\mathbf{e}_i}, i.e. the random
@@ -31,38 +31,40 @@
 #' Note that you can type abbreviations, e.g. "L", "N", "Q", "M" for the
 #' respective types (case is ignored). The argument match is done through
 #' \code{\link[strex]{match_arg}}.
-#' @param weights edge weights, representing the strength/intensity (not the cost!)
+#' @param weights edge weights, representing the strength/intensity (not the
+#'   cost!)
 #'   of each link. If weights is NULL (the default) and g has an edge attribute
 #'   called weight, then it will be used automatically.
-#'   If this is NA then no weights are used (even if the graph has a weight attribute).
-#' @param as_dist If the function should return a matrix or an object of class "dist" as
-#'   returned from [stats::as.dist]. Default is FALSE if the number of nodes is smaller
-#'   than 1000.
+#'   If this is NA then no weights are used (even if the graph has a weight
+#'   attribute).
+#' @param as_dist If the function should return a matrix or an object of class
+#'   "dist" as returned from [stats::as.dist]. Default is FALSE if the number
+#'   of nodes is smaller than 1000.
 #' @param verbose default TRUE
 #' @return The diffusion distance matrix \eqn{D_t}, a square numeric matrix
 #'   of the \eqn{L^2}-norm distances between posterior probability vectors, i.e.
 #'   Euclidean distances between the rows of the stochastic matrix
 #'   \eqn{P(t) = e^{-\tau L}}, where \eqn{-L = -(I - T)} is the generator of the
-#'   continuous-time random walk (Markov chain) of given \code{type} over network
-#'   \code{g}.
+#'   continuous-time random walk (Markov chain) of given \code{type} over
+#'   network \code{g}.
 #' @keywords diffusion distance
 #' @seealso \code{\link{get_diffusion_probability_matrix}}
 #' @references
 #'   De Domenico, M. (2017). Diffusion Geometry Unravels the Emergence of
 #'   Functional Clusters in Collective Phenomena. Physical Review Letters.
 #'   \doi{10.1103/PhysRevLett.118.168301}
 #'
-#'   Bertagnolli, G., & De Domenico, M. (2021). Diffusion geometry of multiplex and
-#'   interdependent systems. Physical Review E, 103(4), 042301.
+#'   Bertagnolli, G., & De Domenico, M. (2021). Diffusion geometry of multiplex
+#'   and interdependent systems. Physical Review E, 103(4), 042301.
 #'   \doi{10.1103/PhysRevE.103.042301}
 #'   \href{https://arxiv.org/abs/2006.13032}{arXiv: 2006.13032}
 #' @examples
 #' g <- igraph::sample_pa(10, directed = FALSE)
 #' dm_crw <- get_distance_matrix(g, tau = 1)
 #' dm_merw <- get_distance_matrix(g, tau = 1, type = "MERW")
 #' @export
-get_distance_matrix <- function(g, tau, type = "Normalized Laplacian", weights = NULL,
-                                as_dist = FALSE, verbose = TRUE) {
+get_distance_matrix <- function(g, tau, type = "Normalized Laplacian",
+                            weights = NULL, as_dist = FALSE, verbose = TRUE) {
   # #by default weights are considered, if present
   # if ( is.null(igraph::E(g)$weight) ) {
   #   cat("Warning: missing edge weights. Assigning 1 by default\n")
@@ -96,22 +98,23 @@ get_distance_matrix <- function(g, tau, type = "Normalized Laplacian", weights =
     DM <- stats::dist(expL)
   }
   if ((!as_dist) && (length(igraph::V(g)) < 1000)) {
-    return(as.matrix(DM))
-  } else {
-    return(DM)
+    DM <- as.matrix(DM)
   }
+  class(DM) <- c("diffudist", class(DM))
+  return(DM)
 }
 
 # get_distance_matrix <- compiler::cmpfun(getDistanceMatrixRaw)
 
 #' @describeIn get_distance_matrix Old deprecated function
-#' @usage getDistanceMatrix(g, tau, type = "Normalized Laplacian", weights = NULL,
-#'                          verbose = TRUE)
+#' @usage getDistanceMatrix(g, tau, type = "Normalized Laplacian",
+#'                          weights = NULL, as_dist = FALSE, verbose = TRUE)
 #' @export
-getDistanceMatrix <- function(g, tau, type = "Normalized Laplacian", weights = NULL,
-                              verbose = TRUE) {
+getDistanceMatrix <- function(g, tau, type = "Normalized Laplacian",
+                              weights = NULL, as_dist = FALSE, verbose = TRUE) {
   .Deprecated("get_distance_matrix")
-  return(get_distance_matrix(g, tau, type = type, weights = weights, verbose = verbose))
+  DM <- get_distance_matrix(g, tau, type = type, weights = weights, verbose = verbose)
+  return(DM)
 }
 
 #' @rdname get_distance_matrix
@@ -139,6 +142,9 @@ get_DDM <- get_distance_matrix
 #' @param Pi a transition matrix (it should be a stochastic matrix)
 #' @param tau diffusion time
 #' @param verbose default TRUE
+#' @param as_dist If the function should return a matrix or an object of class
+#'   "dist" as returned from [stats::as.dist]. Default is FALSE if the number
+#'   of nodes is smaller than 1000.
 #' @return The diffusion distance matrix \eqn{D_t}, a square numeric matrix
 #'   of the \eqn{L^2}-norm distances between posterior probability vectors, i.e.
 #'   Euclidean distances between the rows of the stochastic matrix
@@ -153,15 +159,13 @@ get_DDM <- get_distance_matrix
 #'   Functional Clusters in Collective Phenomena. Physical Review Letters.
 #'   \doi{10.1103/PhysRevLett.118.168301}
 #'
-#'   Bertagnolli, G., & De Domenico, M. (2021). Diffusion geometry of multiplex and
-#'   interdependent systems. Physical Review E, 103(4), 042301.
+#'   Bertagnolli, G., & De Domenico, M. (2021). Diffusion geometry of multiplex
+#'   and interdependent systems. Physical Review E, 103(4), 042301.
 #'   \doi{10.1103/PhysRevE.103.042301}
 #'   \href{https://arxiv.org/abs/2006.13032}{arXiv: 2006.13032}
-#' @examples
-#' g <- igraph::sample_pa(10, directed = FALSE)
-#' dm <- get_distance_matrix(g, tau = 1)
 #' @export
-get_distance_matrix_from_T <- function(Pi, tau, verbose = TRUE) {
+get_distance_matrix_from_T <- function(Pi, tau, as_dist = FALSE,
+                                       verbose = TRUE) {
   Pi <- as.matrix(Pi)
   # check square matrix
   N <- nrow(Pi)
@@ -187,11 +191,13 @@ get_distance_matrix_from_T <- function(Pi, tau, verbose = TRUE) {
     # \sqrt(\sum_i (x_i - y_i) ^ 2))
     DM <- stats::dist(expL)
   }
-  DM <- as.matrix(DM)
+  if ((!as_dist) && (length(igraph::V(g)) < 1000)) {
+    DM <- as.matrix(DM)
+  }
   # names
   colnames(DM) <- colnames(Pi)
   rownames(DM) <- colnames(Pi)
-  # class(DM) <- "DtDistMatrix"
+  class(DM) <- c("diffudist", class(DM))
   return(DM)
 }
 
@@ -206,5 +212,3 @@ get_distance_matrix_from_Pi <- get_distance_matrix_from_T
 #' @rdname get_distance_matrix_from_T
 #' @export
 get_DDM_from_Pi <- get_distance_matrix_from_T
-
-
diff --git a/R/get_mean_distance_matrix.R b/R/get_mean_distance_matrix.R
@@ -19,7 +19,9 @@
 get_mean_distance_matrix <- function(DM_list){
     #DM.list is a simple list, where each element is a distance matrix
     M <- length(DM_list)
-    return( Reduce('+', DM_list) / M )
+    DM <- Reduce('+', DM_list) / M
+    class(DM) <- c("diffudist", class(DM))
+    return(DM)
 }
 
 #' @describeIn get_mean_distance_matrix Old deprecated function

diff --git a/R/get_spectral.R b/R/get_spectral.R
@@ -2,7 +2,11 @@
 #'
 #' @description
 #' Returns the eigenvalue spectrum together with eigenvectors of a Laplacian
-#' corresponding to a network.
+#' corresponding to a network. This involves computing the eigendecomposition of
+#' a (symmetric) matrix, so it is computationally intense and may take some time.
+#' The decomposition of the normalized Laplacian \eqn{L = I - D^{-1}A} takes
+#' is computed through the decomposition of its symmetric version
+#' \eqn{L = D^{-\frac{1}{2}}AD^{-\frac{1}{2}}}. See the package vignette for details.
 #' @param g the network in the [igraph] format
 #' @param type the Laplacian type, default "Normalized Laplacian".
 #'   At the moment this is the only available option. For other types of Laplacians
@@ -32,7 +36,7 @@ get_spectral_decomp <- function(g, type = "Normalized Laplacian", verbose = FALS
   if (igraph::gsize(g) > 0) {
     # L = D - A
     L <- igraph::laplacian_matrix(g, sparse = FALSE)
-    # D^{1/2}
+    # rename D = D^{1/2}
     D <- sqrt(igraph::strength(g, mode = "out"))
     D <- diag(D)
     # D^{-1/2}
@@ -41,16 +45,17 @@ get_spectral_decomp <- function(g, type = "Normalized Laplacian", verbose = FALS
     # (symmetric) normalised Laplacian
     L <- (D_inv %*% L) %*% D_inv
     s_dec <- eigen(L, symmetric = TRUE)
-    u_L <- t(s_dec$vectors) %*% D
+    u_L <- crossprod(s_dec$vectors, D)  # t(s_dec$vectors) %*% D
     u_R <- D_inv %*% s_dec$vectors
   } else {
     stop("Edge set is empty!")
   }
-  return(list(
+  res <- list(
     "lambdas" = s_dec$values,
     "u_L" = u_L,
     "u_R" = u_R
-  ))
+  )
+  return(res)
 }
 
 #' @title Distance Matrix from Laplacian spectral decomposition
@@ -83,17 +88,19 @@ get_spectral_decomp <- function(g, type = "Normalized Laplacian", verbose = FALS
 #'   The matrix exponential is here computed using the given eigendecomposition
 #'   of the Laplacian matrix \eqn{e^{-\tau L} = Q e^{-\tau \Lambda} Q^{-1}}.
 #' @seealso \code{\link{get_spectral_decomp}}
-#' @references Bertagnolli, G., & De Domenico, M. (2021). Diffusion geometry of multiplex and
-#'   interdependent systems. Physical Review E, 103(4), 042301.
+#' @references Bertagnolli, G., & De Domenico, M. (2021). Diffusion geometry of
+#'   multiplex and interdependent systems. Physical Review E, 103(4), 042301.
 #'   \doi{10.1103/PhysRevE.103.042301}
 #'   \href{https://arxiv.org/abs/2006.13032}{arXiv: 2006.13032}
 #' @export
-get_ddm_from_eigendec <- function(tau, Q, Q_inv, lambdas,
+get_ddm_from_eigendec <- function(tau, Q, Q_inv, lambdas, as_dist = FALSE,
                                   verbose = FALSE) {
-  Nodes <- length(lambdas)
-  expL <- eigenMapMatMult(eigenMapMatMult(Q, as.matrix(diag(exp(-tau * lambdas)))), Q_inv)
+  N <- length(lambdas)
+  expL <- eigenMapMatMult(
+    eigenMapMatMult(Q, as.matrix(diag(exp(-tau * lambdas)))), Q_inv
+  )
   # expL <- Q %*% diag(exp(-tau * lambdas)) %*% Q_inv
-  if (abs(sum(expL) - Nodes) > 1e-6) {
+  if (abs(sum(expL) - N) > 1e-6) {
     stop("expL is not a stochastic matrix! Check the row sums.")
   }
   if (verbose) {
@@ -112,10 +119,12 @@ get_ddm_from_eigendec <- function(tau, Q, Q_inv, lambdas,
     # \sqrt(\sum_i (x_i - y_i) ^ 2))
     DM <- stats::dist(expL)
   }
-  DM <- as.matrix(DM)
+  if ((!as_dist) && (length(igraph::V(g)) < 1000)) {
+    DM <- as.matrix(DM)
+  }
   # names
-  colnames(DM) <- rownames(Q)
-  rownames(DM) <- colnames(DM)
+  colnames(DM) <- colnames(Pi)
+  rownames(DM) <- colnames(Pi)
+  class(DM) <- c("diffudist", class(DM))
   return(DM)
 }
-