Taller 2 - Análisis de Redes

suppressMessages(suppressWarnings(library(sna)))
suppressMessages(suppressWarnings(library(igraph)))
suppressMessages(suppressWarnings(library(ggraph)))
suppressMessages(suppressWarnings(library(tidygraph)))
suppressMessages(suppressWarnings(library(scales)))
suppressMessages(suppressWarnings(library(kableExtra)))
suppressMessages(suppressWarnings(library(dplyr)))
#suppressMessages(suppressWarnings(library(patchwork)))

1. El grado de un nodo en una red tanto dirigida como no dirigida se puede calcular fácilmente a partir de la matriz de adyacencia \(𝐘=[y_{i,j}]\). El out-degree \(d^{out}_i\) y el in-degree \(d^{in}_i\) del nodo \(i\) se pueden calcular respectivamente como: \[d^{out}_i=\sum_{j:j\neq i}y_{i,j} \quad \text{ y }\quad d^{in}_i=\sum_{j:j\neq i} y_{j,i}\] Muestre que si la red es no dirigida entonces \(d^{out}_i=d^{in}_i\).

   Solución Punto 1:

   En una red no dirigida, las relaciones pueden pensarse como de dos vías, es decir, una relación entre los nodos \(i\) y \(j\) implica una relación entre \(j\) e \(i\). Es por esto que la matriz de adyacencia es siempre simétrica para este tipo de redes, lo que se escribe como \[y_{i,j}=y_{j,i} \quad\text{ para todo }i\neq j\]

   Ahora bien, retomando la definición de los grados para las redes dirigidas, y usando la simetría mencionada se tiene que:

\[\begin{align*} d^{out}_i&=\sum_{j:j\neq i}y_{i,j} & \text{definición de grado de salida}\\ &= \sum_{j:j\neq i}y_{j, i} & \text{simetría de la matriz de adyacencia}\\ &= d^{in}_i\,\, \blacksquare & \text{definición del grado de entrada} \end{align*}\]

2. Tanto para redes dirigidas como no dirigidas, se define la media global de las interacciones como \[\bar{y}=\frac{1}{n(n-1)}\sum_{i,j:i\neq j}y_{i,j}\] donde \(\mathbf{Y}=[y_{i,j}]\) es la matriz de adyacencia de la red correspondiente. Tal estadístico corresponde a una descripción muy rudimentaria acerca de la plausibilidad de observar una relación entre dos nodos cualesquiera, dado que no tiene en cuenta la heterogeneidad nodal (algunos nodos son más propensos a enviar/recibir más relaciones).

   1. Muestre que para relaciones no dirigidas la media global \(\bar{y}\) es igual a la media tanto de la parte triangular superior de \(𝐘\) como de la parte triangular inferior de \(𝐘\).

      Solución:

      Primero, recordemos que la parte triangular superior de \(\mathbf{Y}\) se refiere a todos aquellos \(y_{ij}\) tal que \(i<j\). Análogamente, la parte triangular inferior se comprende de los \(y_{ij}\) tal que \(i>j\). Una propiedad importante en este caso es que como la red es no dirigida, \(\mathbf{Y}\) es simétrica, es decir, \(y_{i,j}=y_{j,i}\) para todo \(i\neq j\), lo que a su vez implica que \(\sum_{i<j}y_{i,j}=\sum_{i>j}y_{i,j}\).

      Ahora, note que se puede re-expresar la media global de las interacciones como \[\begin{align*} \bar{y}&=\frac{1}{n(n-1)}\sum_{i,j:i\neq j}y_{i,j} \\ &= \frac{1}{n(n-1)} \left[\sum_{i<j}y_{i,j}+\sum_{i>j} y_{i,j} \right] & \text{descomposición de la sumatoria}\\ &= \frac{1}{n(n-1)} \left[\sum_{i<j}y_{i,j}+\sum_{i<j} y_{i,j} \right] & \text{igualdad de las dos sumatorias}\\ &= \frac{1}{n(n-1)} \left[2\sum_{i<j}y_{i,j}\right] & \text{simplificación}\\ &=\frac{1}{n(n-1)/2} \sum_{i<j}y_{i,j} & \text{reexpresión del denominador}\\ &= \text{media de la parte triangular superior de } \mathbf{Y} \end{align*}\] La otra parte de la demostración es inmediata luego de haber visto el desarrollo anterior \[\begin{align*} \bar{y}&=\frac{1}{n(n-1)}\sum_{i,j:i\neq j}y_{i,j} = \frac{1}{n(n-1)} \left[\sum_{i<j}y_{i,j}+\sum_{i>j} y_{i,j} \right] = \frac{1}{n(n-1)} \left[\sum_{i>j}y_{i,j}+\sum_{i>j} y_{i,j} \right] \\ &= \frac{1}{n(n-1)} \left[2\sum_{i>j}y_{i,j}\right] =\frac{1}{n(n-1)/2} \sum_{i>j}y_{i,j}\\ &= \text{media de la parte triangular inferior de } \mathbf{Y} \end{align*}\]

   2. Muestre que tanto para relaciones dirigidas como no dirigidas la media global corresponde a la densidad de la red.

      Solución:

      Caso dirigido: Sea \(n\) es el número de vértices, para una red dirigida la densidad se expresa como \[\text{den}=\frac{|E|}{|V|(|V|-1)}=\frac{|E|}{n(n-1)}\] donde \(|E|\) es el número de aristas en la red.

      Dado que \(y_{i,j}=1\) si existe una relación entre \(i\) y \(j\), es natural representar el número de aristas en la red como la suma de los elementos de la matriz de adyacencia. De este modo, \[|E|=\sum_i\sum_{j:j\neq i} y_{i,j}\] lo que implica que \[\begin{align*} \text{den}&=\frac{|E|}{n(n-1)}=\frac{1}{n(n-1)}\sum_i\sum_{j:j\neq i} y_{i,j}\\ &=\frac{1}{n(n-1)}\sum_{i,j:i\neq j} y_{i,j} \\ &= \bar{y} \end{align*}\]

      Caso no dirigido: Para una red no dirigida la densidad se expresa como \[\text{den}=\frac{|E|}{|V|(|V|-1)/2}=\frac{|E|}{n(n-1)/2}\] Por el mismo razonamiento mostrado en el caso dirigido, el número de aristas corresponde, en este caso, a \[|E|=\frac{1}{2}\sum_i\sum_{j:j\neq i} y_{i,j}\]

      Así, \[\begin{align*} \text{den}&=\frac{|E|}{n(n-1)/2}=\frac{1}{n(n-1)/2}\frac{1}{2}\sum_i\sum_{j:j\neq i} y_{i,j}\\ &=\frac{2}{n(n-1)}\frac{1}{2}\sum_{i,j:i\neq j} y_{i,j} \\ &=\frac{1}{n(n-1)}\sum_{i,j:i\neq j} y_{i,j} \\ &= \bar{y} \end{align*}\]

   3. Muestre que tanto para relaciones dirigidas como no dirigidas se tiene que \((n-1)\bar{y}=\bar{d}^{out}=\bar{d}^{in}\). Es decir, el grado promedio tanto de entrada como de salida son iguales y a su vez equivalentes la densidad.

      Solución:

      En el numeral anterior se probó que \(\bar{y}=\text{den}\), y acá se usará esta igualdad.

      Caso no dirigido:

      \[\begin{align*} (n-1)\bar{y}=(n-1)\text{den}=(n-1)\frac{|E|}{n(n-1)/2}=\frac{|E|}{n/2}=\frac{2|E|}{n}=\bar{d} \end{align*}\]

      Como se vio en el primer numeral, para la red no dirigida \(\bar{d}=\bar{d}^{in}=\bar{d}^{out}\).

      Caso dirigido:

      \[\begin{align*} (n-1)\bar{y}=(n-1)\text{den}=(n-1)\frac{|E|}{n(n-1)}=\frac{|E|}{n} \end{align*}\] Antes de seguir, recordemos que \[\begin{align*} d_i^{out} &= \sum_{j:j\neq i} y_{i,j} & d_i^{in} &= \sum_{j:j\neq i} y_{j,i} & \text{definición de grado}\\ \sum_i d_i^{out} &= \sum_i \sum_{j:j\neq i} y_{i,j} & \sum_i d_i^{in} &= \sum_i \sum_{j:j\neq i} y_{j,i} & \text{sumar sobre i}\\ \sum_i d_i^{out} &= |E| & \sum_i d_i^{in} &= |E| & \text{representación del numeral (b)} \end{align*}\]

      Así, retomando la ecuación anterior nos queda \[\begin{align*} (n-1)\bar{y}&=\frac{|E|}{n} = \frac{\sum_{i}d_1^{out}}{n}= \frac{\sum_{i}d_1^{in}}{n}\\ &= \bar{d}^{out}=\bar{d}^{in} \end{align*}\]

3. Considere un grafo estrella de orden \(n\) y un grafo círculo de orden \(n\). A continuación se representan ambos grafos para \(n=9\). Estos grafos tienen aproximadamente la misma densidad, pero su estructura es muy diferente.

   1. Muestre que para el grafo círculo de orden \(n\) el grado promedio es \(\bar{d}=2\).
   2. Muestre que para el grafo estrella de orden \(n\) el grado promedio es \(\bar{d}=2\frac{n-1}{n}\) cuando \(n\rightarrow\infty\).

Solución Punto 3:

Grafo círculo: Podemos pensar a este grafo como uno en el que todo nodo tiene un nodo antes de él y otro enseguida. De este modo, es inmediato que cada nodo está conectado solo con 2 nodos, los mencionados anteriormente. Así: \[\bar{d}=\frac{1}{n}\sum_{i=1}^n d_i = \frac{1}{n} [2+2+...+2]=\frac{1}{n} 2n= 2\] Grafo estrella: Este grafo tiene una particularidad, y es que todos los nodos exteriores tienen grado 1, ya que su única conexión es con el nodo del centro; mientras que este nodo central está conectado con todos los (n-1) nodos restantes. Así: \[\begin{align*} \bar{d}&=\frac{1}{n}\sum_{i=1}^n d_i =\frac{1}{n} [1+1+...1+(n-1)]=\frac{1}{n}[(n-1)+(n-1)]=\frac{1}{n}2(n-1)\\ &= \frac{1}{n}(2n-2)=\frac{2n}{n}-\frac{2}{n}=2-\frac{2}{n} \end{align*}\] El cual, al tomar el límite cuando \(n\rightarrow \infty\) resulta en \[\lim_{n\rightarrow\infty} \bar{d} = \lim_{n\rightarrow\infty} \left(2-\frac{2}{n}\right)=2-0=2\]

4. ¿Cuáles de las siguientes secuencias son caminatas en el grafo que se presenta a continuación? ¿Cuáles senderos? ¿Cuáles circuitos? ¿Cuáles ciclos?

   1. 2−1−6−3−4.
   2. 2−1−6−3−4−1−5.
   3. 2−1−2−5−1−4.

Solución Punto 4:

Del libro de Kolaczyk (2020) tenemos las siguientes definiciones:

Una caminata de \(v_0\) a \(v_l\) es una secuencia alternante \(\{v_0, e_1, v_1,e_2, ..., v_{l-1}, e_l, v_l\}\) donde los extremos de \(e_i\) son \(\{v_{i-1}, v_i\}\). Un sendero es una caminata sin aristas repetidas. Un sendero para el cual los vertices extremos son los mismos es llamado circuito. Finalmente, un ciclo es una caminata de \(l>3\) con vertices extremos iguales, además de que todos los demás vertices deben ser distintos entre si.

Teniendo en cuenta estas definiciones, veamos cada secuencia:

g_4 <- make_graph(c( 1,2, 2,5, 5,1, 1,4, 4,3, 3,6, 6,1), directed = FALSE)

layout_fixed <- matrix(c(
  0,0,    # 1
  1,0.5,  # 2
  -1.2,0.2, # 3
  -0.7,-0.5,# 4
  1,-0.5, # 5
  -0.8,1     # 6
  ), ncol=2, byrow=TRUE)

plot_walk <- function(g, walk, title_text) {

  # aristas de la secuencia
  edge_list <- c()
  
  for(i in 1:(length(walk)-1)){
    edge_list <- c(edge_list, get_edge_ids(g, vp = c(walk[i], walk[i+1])))
  }

  edge_colors <- rep("gray80", ecount(g))
  edge_colors[edge_list] <- "red"

  plot(g,
    layout = layout_fixed,
    edge.color = edge_colors,
    edge.width = ifelse(edge_colors=="red", 3, 1),
    vertex.color = "lightblue",
    vertex.size = 30,
    vertex.label.cex = 1.2,
    main = title_text)
}

walk_a <- c(2,1,6,3,4)
walk_b <- c(2,1,6,3,4,1,5)
walk_c <- c(2,1,2,5,1,4)

par(mfrow=c(1,3), mar=c(1,1,2,1)) 
plot_walk(g_4, walk_a, "Secuencia a")
plot_walk(g_4, walk_b, "Secuencia b")
plot_walk(g_4, walk_c, "Secuencia c")

La clasificación de las secuencias se presenta en la siguiente tabla:

Secuencia	Caminata	Sendero	Circuito	Ciclo
(a)	✓	✓	✗	✗
(b)	✓	✓	✗	✗
(c)	✓	✗	✗	✗

Note que, la secuencia (a) es un sendero ya que es una caminara abierta en la que no se repiten aristas, pero no es un ciclo ya que sus vertices extremos no son iguales. A su vez, la secuencia (b) tiene la misma clasificación ya que hace el mismo recorrido de (a) añaniendo 2 aristas al final. En el caso de que (b) hubiese tenido una arista adicional de 5-2, esta secuencia sería también un circuito pero no un ciclo. Finalmente, la secuencia (c) solo es una caminata ya que repite la arista 2-1 al inicio de su configuración.

5. Considere el conjunto de datos dado en comtrade.RData (este archivo contiene una arreglo de tres dimensiones denominado comtrade), asociado con el crecimiento anual del comercio (diferencia en dólares en escala logarítmica respecto al año 2000). Este conjunto de datos involucra 30 países, 10 años desde 1996 hasta 2005, y 6 clases de productos diferentes, como se muestra a continuación:

   1. Calcule el aumento medio global \(\bar{y}\) a lo largo de los 10 años en bienes manufacturados. Para ello considere la matriz de adyacencia Y dada por:

   2. Calcule la media de todas las observaciones de cada fila de Y, es decir, calcule la media fila \(\bar{y}_{i\cdot}=\frac{1}{n-1}\sum_{j:j\neq i}y_{i, j}\) para cada país. Realice una histograma de los promedios fila \(\bar{y}_{i\cdot}\). Los promedios fila caracterizan diferentes niveles de actividad de los nodos en términos de la sociabilidad. ¿Cómo se pueden interpretar los promedios fila \(\bar{y}_{i\cdot}\)?

   3. Calcule la media de todas las observaciones de cada columna de Y, es decir, calcule la media columna \(\bar{y}_{\cdot j}=\frac{1}{n-1}\sum_{i:i\neq i}y_{i, j}\) para cada país. Realice una histograma de los promedios columna \(\bar{y}_{\cdot j}\). Los promedios columna caracterizan diferentes niveles de actividad de los nodos en términos de la popularidad. ¿Cómo se pueden interpretar los promedio columna \(\bar{y}_{\cdot j}\)?

   4. Calcule tanto la media de los promedios fila \(\bar{y}_{i\cdot}\) como la media de los promedios columna \(\bar{y}_{\cdot j}\).¿Qué se puede concluir acerca de la tendencia local en este caso?

   5. Calcule tanto la desviación estándar (DE) de los promedios fila \(\bar{y}_{i\cdot}\) como la DE de los promedios columna \(\bar{y}_{\cdot j}\). ¿Qué se puede concluir acerca de la heterogeneidad local en este caso?

   6. Calcule el coeficiente de correlación (CV) entre los promedios fila \(\bar{y}_{i\cdot}\) y los promedios columna \(\bar{y}_{\cdot j}\). Realice un dispersograma de los promedios columna \(\bar{y}_{\cdot j}\) (eje \(y\)) frente a los promedios fila \(\bar{y}_{i\cdot}\) (eje \(x\)), junto con la recta \(y=x\) como punto de referencia. ¿Qué se puede concluir?

Solución Punto 5

En primera instancia, se cargan los datos

load("comtrade.RData")
dimnames(comtrade)[c(1,3,4)]

## [[1]]
##  [1] "Australia"            "Austria"              "Brazil"              
##  [4] "Canada"               "China"                "China, Hong Kong SAR"
##  [7] "Czech Rep."           "Denmark"              "Finland"             
## [10] "France"               "Germany"              "Greece"              
## [13] "Indonesia"            "Ireland"              "Italy"               
## [16] "Japan"                "Malaysia"             "Mexico"              
## [19] "Netherlands"          "New Zealand"          "Norway"              
## [22] "Rep. of Korea"        "Singapore"            "Spain"               
## [25] "Sweden"               "Switzerland"          "Thailand"            
## [28] "Turkey"               "United Kingdom"       "USA"                 
## 
## [[2]]
## [1] "Chemicals"                                    
## [2] "Crude materials, inedible, except fuels"      
## [3] "Food and live animals"                        
## [4] "Machinery and transport equipment"            
## [5] "Manufact goods classified chiefly by material"
## [6] "Miscellaneous manufactured articles"          
## 
## [[3]]
##  [1] "1996" "1997" "1998" "1999" "2000" "2001" "2002" "2003" "2004" "2005"

Ese arreglo de datos 4-dimensional se va a tratar solamente como una matriz (2-dimensional), al promediar todos los años y quedarse solamente con la información de la quinta y sexta clase de productos diferentes.

Y <- apply(X = comtrade[,,c(5,6),], MARGIN = c(1,2), FUN = mean)
# Se están considerando solamente los productos 5 y 6
# Se están promediando todos los años

La matriz de adyacencia nos muestra que esta es una red dirigida ponderada de tipo simple.

Primero, se calcula el aumento medio global dado por la formula \(\bar{y}=\frac{1}{n(n-1)}\sum_{i,j:i\neq j}y_{i,j}\)

n <- nrow(Y)
y_bar <- sum(Y, na.rm = TRUE) / (n*(n-1)); y_bar

## [1] 0.03778362

El aumento medio global es de 0.0377, este valor nos dice que en promedio los países involucrados presentaron un crecimiento leve de comercio en el sector manufacturero durante los años 1996 a 2005, ya que este es positivo.

Ahora, calculamos la media por filas de todas las observaciones \(\bar{y}_{i.}=\frac{1}{n-1}\sum_{j:j\neq i} y_{i,j}\)

y_bar_i. <- (apply(Y, 1, sum, na.rm = TRUE))/(n-1)
hist(y_bar_i., main = "Distribucion del crecimiento promedio de exportaciones \n del pais i hacia los demas",
     xlab = "Crecimiento promedio", ylab = "Frecuencia") # main = "Distribución de la sociabilidad de los países

Como se indica en el enunciado, los promedios fila sirven para caracterizar la sociabilidad de los nodos, lo que en este contexto es el nivel de actividad promedio de cada país como exportador. El comportamiento asimétrico a derecha del histograma deja ver que la tendencia es que los países no sean tan sociables, es decir, que presenten un crecimiento moderado. Aunque hay algunos que sobresalen de la distribución.

Particularmente, los países con un promedio por fila más grande son China (0.14), Turquía (0.10) y la República Checa (0.09). Esto se traduce en que estos son los países más sociables comercialmente hablando en el periodo estudiado.

A continuación, calculamos la media por columnas de todas las observaciones \(\bar{y}_{.j}=\frac{1}{n-1}\sum_{i:i\neq i} y_{i,j}\)

y_bar_.j <- (apply(Y, 2, sum, na.rm = TRUE))/(n-1)
hist(y_bar_.j, main = "Distribucion del crecimiento promedio de importaciones \n del pais j",
     xlab = "Crecimiento promedio", ylab = "Frecuencia") # main = "Distribución de la popularidad de los países

Los promedios columna ayudan a caracterizar la popularidad de los nodos, es decir, el nivel de actividad promedio que tienen los países como importadores. En este caso, si vemos países que presentaron disminución en sus importaciones, y la tendencia sigue siendo una popularidad leve (0 a 0.05).

No obstante, hay un grupo de países cuya popularidad fue mayor y estos son China (0.145), México (0.142) y Turquía (0.11). Note que China y Turquía podrían tomarse como los países más activos económicamente hablando ya que fueron tanto sociables como populares.

Como complemento al análisis, se calculan las medias y deasviaciones estándar de los promedios calculados por fila y por columna.

mean(y_bar_i.)

## [1] 0.03778362

mean(y_bar_.j)

## [1] 0.03778362

La igualdad de las medias refleja que en promedio, el sistema comercial está balanceado aunque las importaciones y exportaciones de cada país no sean las mismas. Se puede afirmar que estamos frente a una tendencia local positiva ya que se presentó un aumento (leve) medio global del comercio en bienes manufacturados

sd(y_bar_i.)

## [1] 0.03019967

sd(y_bar_.j)

## [1] 0.04101555

Por otro lado, las desviaciones estándar nos dicen que tan variables fueron estos comportamientos. Estos valores no son pequeños respecto de la media que se tiene. Por ejemplo, para la media por fila: \(CV: 0.0302/0.038=79,5\%\) lo que significa que la actividad comercial de los países fue muy dispersa.

Que la desviación estándar de los promedios columna sea mayor al de los promedios fila revela que hay más heterogeneidad en la popularidad de los nodos que en la sociabilidad.

# Coeficiente de correlación
cor(y_bar_i., y_bar_.j)

## [1] 0.7002526

Este coeficiente de correlación superior al 70% muestra que la sociabilidad y popularidad están fuertemente correlacionadas linealmente en sentido positivo.

plot(x = y_bar_i., y = y_bar_.j, xlab = "Promedios fila (exportaciones)", ylab = "Promedios columna (importaciones)", 
     pch = 19, main = "Dispersograma de los promedios")

abline(a = 0, b = 1, col = "red", lwd = 2)

# Resaltar el punto más lejano
points(y_bar_i.[18], y_bar_.j[18],
       col = "blue", pch = 19, cex = 1.5)
# Etiquetar el punto
text(y_bar_i.[18], y_bar_.j[18],
     labels = names(y_bar_i.)[18], pos = 4, col = "blue")

El gráfico de dispersión muestra una asociación positiva entre ambos promedios, lo cual va en concordancia con el valor del coeficiente de correlación calculado anteriormente. Así mismo, se observa que la nube de puntos es más densa en la parte inferior izquierda del gráfico ya que priman la sociabilidad y popularidad baja, aunque hay excepciones de países con una alta movilidad en el aspecto económico.

Finalmente, se observa un punto fuera de la tendencia general: México. Este país es el único que presenta un bajo nivel de exportaciones pero un alto nivel de importaciones.

6. Considere el conjunto de datos dado en conflict.RData recopilado por Mike Ward y Xun Cao del departamento de Ciencias Políticas de la Universidad de Washington, asociado con datos de conflictos entre países en los años 90. El archivo conflict.RData contiene una lista con tres arreglos, X, Y, y D. X tiene tres campos: population (población en millones), gdp (PIB en millones de dolares) polity (puntuación política, un índice de democracia). Y hace referencia a una matriz \(𝐘=[y_{i,j}]\) en la que \(y_{i,j}\) representa el número de conflictos iniciados por el país \(i\) hacia el país \(j\). Finalmente, D es un arreglo de tres dimensiones dimensiones cuya tercera dimensión contiene indices entre cada par de países asociados con: comercio (dimensión 1), importaciones (dimensión 2), organizaciones intergubernamentales (dimensión 3), y distancia geográfica (dimensión 4).

   1. Hacer una visualización decorada de la red de conflictos teniendo en cuenta diferentes diseños.

   2. Calcule e interprete la media global.

   3. Obtenga y grafique la distribución del out-degree y del in-degree. Calcule e interprete la media y la desviación estándar de esta distribución.

   4. Calcule el coeficiente de correlación entre los valores del out-degree y el in-degree. Realice un dispersograma de los grados de entrada (eje \(y\)) frente a los grados de salida (eje \(x\)), junto con la recta \(y=x\) como punto de referencia. ¿Qué se puede concluir?

   5. Identifique los países mas activos.

Solución Punto 6

load("conflict.RData")
g6_conflict <- graph_from_adjacency_matrix(dat$Y, mode = "directed", weighted = TRUE, diag = FALSE)
vertex_attr(g6_conflict) <- as.data.frame(dat$X)

Como Ward & Cao (2007) mencionan en su artículo, las variables PIB (GDP) y población son de tipo continuo, esta última de carácter positivo, y la variable Democracia y Autocracia (polity) es un puntaje que se mueve en el intervalo (-10, 10). Con fines ilustrativos, esta última variable se va a categorizar como: Democracia si el puntaje es positivo y Autocracia si este es negativo.

# Categorización de la variable 'polity'
V(g6_conflict)$polity <- ifelse(V(g6_conflict)$polity>0, "Democracia", "Autocracia")
V(g6_conflict)$name <- rownames(dat$Y)

# Establecer entre quienes se da la relación
edges <- ends(g6_conflict, E(g6_conflict))
pol_from <- V(g6_conflict)$polity[match(edges[,1], V(g6_conflict)$name)]
pol_to <- V(g6_conflict)$polity[match(edges[,2], V(g6_conflict)$name)]
E(g6_conflict)$tipo_regimen <- paste(pol_from,pol_to,sep = "-")
E(g6_conflict)$tipo_regimen <- ifelse(pol_from == "Autocracia" & pol_to == "Autocracia","Autocracia-Autocracia",
                                      ifelse(pol_from == "Democracia" & pol_to == "Democracia","Democracia-Democracia",
                                             "Mixta")
                                      )

Adicionalmente, recordando los análisis del Taller 1, hay dos actores especialmente relevantes para esta red: Iraq y Estados Unidos. Su importancia se debe a que son los dos países que más conflictos han recibido, pero al mismo tiempo hacen parte de los más conflictivos. Es por esto que en las visualizaciones los vamos a resaltar.

Nota: Se usa la variable PIB (GDP) como la que caracteriza el tamaño de los vértices, no se muestran visualizaciones en las que la variable Población cumpla este rol dado que esta última variable generó diseños muy parecidas a la primera, y se considera que no aporta nada adicional al análisis.

g6_tbl <- as_tbl_graph(g6_conflict)

set.seed(9)
set_graph_style(plot_margin = margin(1,1,1,1))

# Layout
layout6 <- create_layout(g6_tbl, layout = 'igraph', algorithm = "star", center = 58)

ggraph(layout6) +
  geom_edge_link(aes(width = weight,
                     color = tipo_regimen),
                 alpha = 0.5,
                 arrow = arrow(length = unit(2, "mm")),
                 end_cap = circle(2, "mm")) +
  geom_node_point(aes(size = log(gdp),
                      color = polity)) +
  geom_node_label(
    aes(label = ifelse(name %in% c("IRQ","USA"), name, NA),
        fill = polity),
    linewidth = 0,
    size = 4,
    show.legend = FALSE
  ) +
  scale_edge_width(range = c(0.5, 4), guide = "none") +
  scale_edge_color_manual(
  name = "Tipo de relación",
  values = c(
    "Autocracia-Autocracia" = "dodgerblue4",  # azul coherente
    "Democracia-Democracia" = "darkgreen",    # verde coherente
    "Mixta" = "mediumpurple3"                 # mezcla visual azul+verde
  )
)+
  scale_size(range = c(1, 7), guide = "none") +
  scale_color_manual(
    name = "Regimen politico",
    values = c("Democracia" = "green",
               "Autocracia" = "dodgerblue3")
  ) +
  scale_fill_manual(
    values = c("Democracia" = "green",
               "Autocracia" = "dodgerblue3"
    ),
    guide = "none"
  ) +
  theme_void() +
  theme(legend.position = "right") +
  labs(
    title = "Red para los conflictos entre paises - Layout = star",
    subtitle = "El tamano de los vertices representa el log(PIB)"
  )

set.seed(9)
set_graph_style(plot_margin = margin(1,1,1,1))

# Layout
layout6 <- create_layout(g6_tbl, layout = 'igraph', algorithm = "nicely")

ggraph(layout6) +
  geom_edge_link(aes(width = weight,
                     color = tipo_regimen),
                 alpha = 0.5,
                 arrow = arrow(length = unit(2, "mm")),
                 end_cap = circle(2, "mm")) +
  geom_node_point(aes(size = log(gdp),
                      color = polity)) +
  geom_node_label(
    aes(label = ifelse(name %in% c("IRQ","USA"), name, NA),
        fill = polity),
    linewidth = 0,
    size = 4,
    show.legend = FALSE
  ) +
  scale_edge_width(range = c(0.5, 4), guide = "none") +
  scale_edge_color_manual(
  name = "Tipo de relación",
  values = c(
    "Autocracia-Autocracia" = "dodgerblue4",  # azul coherente
    "Democracia-Democracia" = "darkgreen",    # verde coherente
    "Mixta" = "mediumpurple3"                 # mezcla visual azul+verde
  )
)+
  scale_size(range = c(1, 7), guide = "none") +
  scale_color_manual(
    name = "Regimen politico",
    values = c("Democracia" = "green",
               "Autocracia" = "dodgerblue3")
  ) +
  scale_fill_manual(
    values = c("Democracia" = "green",
               "Autocracia" = "dodgerblue3"
    ),
    guide = "none"
  ) +
  theme_void() +
  theme(legend.position = "right") +
  labs(
    title = "Red para los conflictos entre paises - Layout = nicely",
    subtitle = "El tamano de los vertices representa el log(PIB)"
  )

set.seed(123)
set_graph_style(plot_margin = margin(1,1,1,1))

# Layout
layout6 <- create_layout(g6_tbl, layout = 'igraph', algorithm = "dh")

ggraph(layout6) +
  geom_edge_link(aes(width = weight,
                     color = tipo_regimen),
                 alpha = 0.5,
                 arrow = arrow(length = unit(2, "mm")),
                 end_cap = circle(2, "mm")) +
  geom_node_point(aes(size = log(gdp),
                      color = polity)) +
  geom_node_label(
    aes(label = ifelse(name %in% c("IRQ","USA"), name, NA),
        fill = polity),
    linewidth = 0,
    size = 4,
    show.legend = FALSE
  ) +
  scale_edge_width(range = c(0.5, 4), guide = "none") +
  scale_edge_color_manual(
  name = "Tipo de relación",
  values = c(
    "Autocracia-Autocracia" = "dodgerblue4",  # azul coherente
    "Democracia-Democracia" = "darkgreen",    # verde coherente
    "Mixta" = "mediumpurple3"                 # mezcla visual azul+verde
  )
)+
  scale_size(range = c(1, 7), guide = "none") +
  scale_color_manual(
    name = "Regimen politico",
    values = c("Democracia" = "green",
               "Autocracia" = "dodgerblue3")
  ) +
  scale_fill_manual(
    values = c("Democracia" = "green",
               "Autocracia" = "dodgerblue3"
    ),
    guide = "none"
  ) +
  theme_void() +
  theme(legend.position = "right") +
  labs(
    title = "Red para los conflictos entre paises - Layout = dh",
    subtitle = "El tamano de los vertices representa el log(PIB)"
  )

set.seed(123)
set_graph_style(plot_margin = margin(1,1,1,1))

# Layout
layout6 <- create_layout(g6_tbl, layout = 'igraph', algorithm = "graphopt")

ggraph(layout6) +
  geom_edge_link(aes(width = weight,
                     color = tipo_regimen),
                 alpha = 0.5,
                 arrow = arrow(length = unit(2, "mm")),
                 end_cap = circle(2, "mm")) +
  geom_node_point(aes(size = log(gdp),
                      color = polity)) +
  geom_node_label(
    aes(label = ifelse(name %in% c("IRQ","USA"), name, NA),
        fill = polity),
    linewidth = 0,
    size = 4,
    show.legend = FALSE
  ) +
  scale_edge_width(range = c(0.5, 4), guide = "none") +
  scale_edge_color_manual(
  name = "Tipo de relación",
  values = c(
    "Autocracia-Autocracia" = "dodgerblue4",  # azul coherente
    "Democracia-Democracia" = "darkgreen",    # verde coherente
    "Mixta" = "mediumpurple3"                 # mezcla visual azul+verde
  )
)+
  scale_size(range = c(1, 7), guide = "none") +
  scale_color_manual(
    name = "Regimen politico",
    values = c("Democracia" = "green",
               "Autocracia" = "dodgerblue3")
  ) +
  scale_fill_manual(
    values = c("Democracia" = "green",
               "Autocracia" = "dodgerblue3"
    ),
    guide = "none"
  ) +
  theme_void() +
  theme(legend.position = "right") +
  labs(
    title = "Red para los conflictos entre paises - Layout = graphopt",
    subtitle = "El tamano de los vertices representa el log(PIB)"
  )

set.seed(1)
set_graph_style(plot_margin = margin(1,1,1,1))

# Layout
layout6 <- create_layout(g6_tbl, layout = 'igraph', algorithm = "fr")

ggraph(layout6) +
  geom_edge_link(aes(width = weight,
                     color = tipo_regimen),
                 alpha = 0.5,
                 arrow = arrow(length = unit(2, "mm")),
                 end_cap = circle(2, "mm")) +
  geom_node_point(aes(size = log(gdp),
                      color = polity)) +
  geom_node_label(
    aes(label = ifelse(name %in% c("IRQ","USA"), name, NA),
        fill = polity),
    linewidth = 0,
    size = 4,
    show.legend = FALSE
  ) +
  scale_edge_width(range = c(0.5, 4), guide = "none") +
  scale_edge_color_manual(
  name = "Tipo de relación",
  values = c(
    "Autocracia-Autocracia" = "dodgerblue4",  # azul coherente
    "Democracia-Democracia" = "darkgreen",    # verde coherente
    "Mixta" = "mediumpurple3"                 # mezcla visual azul+verde
  )
)+
  scale_size(range = c(1, 7), guide = "none") +
  scale_color_manual(
    name = "Regimen politico",
    values = c("Democracia" = "green",
               "Autocracia" = "dodgerblue3")
  ) +
  scale_fill_manual(
    values = c("Democracia" = "green",
               "Autocracia" = "dodgerblue3"
    ),
    guide = "none"
  ) +
  theme_void() +
  theme(legend.position = "right") +
  labs(
    title = "Red para los conflictos entre paises - Layout = fr",
    subtitle = "El tamano de los vertices representa el log(PIB)"
  )

set.seed(123)
set_graph_style(plot_margin = margin(1,1,1,1))

# Layout
layout6 <- create_layout(g6_tbl, layout = 'igraph', algorithm = "lgl")

ggraph(layout6) +
  geom_edge_link(aes(width = weight,
                     color = tipo_regimen),
                 alpha = 0.5,
                 arrow = arrow(length = unit(2, "mm")),
                 end_cap = circle(2, "mm")) +
  geom_node_point(aes(size = log(gdp),
                      color = polity)) +
  geom_node_label(
    aes(label = ifelse(name %in% c("IRQ","USA"), name, NA),
        fill = polity),
    linewidth = 0,
    size = 4,
    show.legend = FALSE
  ) +
  scale_edge_width(range = c(0.5, 4), guide = "none") +
  scale_edge_color_manual(
  name = "Tipo de relación",
  values = c(
    "Autocracia-Autocracia" = "dodgerblue4",  # azul coherente
    "Democracia-Democracia" = "darkgreen",    # verde coherente
    "Mixta" = "mediumpurple3"                 # mezcla visual azul+verde
  )
)+
  scale_size(range = c(1, 7), guide = "none") +
  scale_color_manual(
    name = "Regimen politico",
    values = c("Democracia" = "green",
               "Autocracia" = "dodgerblue3")
  ) +
  scale_fill_manual(
    values = c("Democracia" = "green",
               "Autocracia" = "dodgerblue3"
    ),
    guide = "none"
  ) +
  theme_void() +
  theme(legend.position = "right") +
  labs(
    title = "Red para los conflictos entre paises - Layout = lgl",
    subtitle = "El tamano de los vertices representa el log(PIB)"
  )

Se calcula ahora la media global \(\bar{y}\)

n <- nrow(dat$Y)
(y_bar6 <- sum(dat$Y, na.rm = TRUE)/(n*(n-1)) )

## [1] 0.01818724

La media global, al coincidir con la densidad de la red, se puede interpretar como que tan cerca está la red de ser un ‘clique’. Para la presente red, un valor de 0.018 evidencia que la red está lejos de considerarse un grafo completo, esto es, que es poco densa.

A continuación, se grafican la distribución del grado de salida y la del grado de entrada

d6_out <- igraph::degree(graph = g6_conflict, mode = "out")
d6_in <- igraph::degree(graph = g6_conflict, mode = "in")

par(mfrow=c(1, 2))
plot(x = NA, 
     y = NA, 
     type = "n", 
     xlim = c(0,30), 
     ylim = c(0,0.4), 
     xlab = "Grado", 
     ylab = "Densidad", 
     main = "Distribucion del grado de salida")
hist(d6_out, freq = F, col = "lightskyblue", border = "royalblue", add = T, breaks = seq(0, 30, 2))

plot(x = NA, 
     y = NA, 
     type = "n", 
     xlim = c(0,30), 
     ylim = c(0,0.4), 
     xlab = "Grado", 
     ylab = "Densidad", 
     main = "Distribucion del grado de entrada")
hist(d6_in, freq = F, col = "lightskyblue", border = "royalblue", add = T, breaks = seq(0, 16, 2))

La distribución de los grados sigue un comportamiento asimétrico con cola a derecha. Concluimos que predomina en la red la heterogeneidad de la conectividad, con la mayoría de nodos poco conectados y unos muy pocos bastante enlazados.

# Media y desviación estándar de la distribución
mean(d6_in)

## [1] 1.561538

sd(d6_in)

## [1] 1.984451

mean(d6_out)

## [1] 1.561538

sd(d6_out)

## [1] 3.589398

Las medias de los grados de salida y entrada de la red son iguales, resultado acorde a la teoría que nos dice \(\bar{d}^{in}=\bar{d}^{out}\). El valor de 1.561 indica que en promedio, los actores de la red están conectados con 1 a 2 nodos. La diferencia entre ambas distribuciones radica en su variabilidad. Los grados de salida, es decir, el número de conflictos que inicia un país, son mucho más dispersos que los de entrada; Iraq y Jordania son los países cuyos conflictos iniciados superan los 25. Mientras que los grados de entrada toman su valor máximo de 15 con Iraq.

cor(d6_in, d6_out)

## [1] 0.6040145

La correlación muestra una asociación lineal moderada-alta positiva entre los grados de salida y los grados de entrada. Esto nos mostraría que hay cierto equilibrio entre los conflictos que genera un país y los que recibe, sabiendo que hay casos extremos que no siguen esta tendencia.

plot(x = d6_out, y = d6_in, xlab = "Grados de salida", ylab = "Grados de entrada",
  pch = 19, main = "Dispersograma de los grados de entrada y salida")
abline(a = 0, b = 1, col = "red", lwd = 2)

points(d6_out["IRQ"], d6_in["IRQ"], col = "blue", pch = 19, cex = 1.5)
points(d6_out["JOR"], d6_in["JOR"], col = "blue", pch = 19, cex = 1.5)

text(d6_out["IRQ"], d6_in["IRQ"], labels = names(d6_out)[58], pos = 2, col = "blue")
text(d6_out["JOR"], d6_in["JOR"], labels = names(d6_out)[62], pos = 2, col = "blue")

El gráfico de dispersión deja ver una asociación positiva entre ambos grados de la red, aunque los puntos no se acercan del todo a la recta \(y=x\). Se observan dos puntos alejados de la nube de puntos dado su comportamiento atípico en el eje x: Iraq que siempre está metido en conflictos (ya sea que los inicie o no), y Jordania, país bastante conflictivo. Podemos concluir que, en general hay una asimetría estructural en la red debido a que hay países muy activos en temas conflictos y otros muy pasivos.

ord_grados_in6 <- sort(d6_in, decreasing = TRUE)
ord_grados_out6 <- sort(d6_out, decreasing = TRUE)

head(ord_grados_in6, 5)

## IRQ USA HAI TUR JPN 
##  15   8   7   6   5

head(ord_grados_out6, 5)

## IRQ JOR USA UGA CHN 
##  27  26  11   7   6

Para identificar a los países más activos formalmente se calculan los grados de salida y entrada de los vértices. Los países que han sido foco de más conflictos son Iraq, Estados Unidos, Haiti, Turquia y Japon, con 15, 8, 7, 6 y 5 conflictos recibidos. Mientras tanto, los países que más conflictos han iniciado son Iraq, Jordania, Estados Unidos, Uganda y China, con 27, 26, 11, 7 y 6 conflictos.

7. Para todos los vértices de los cuatro grafos que se presentan a continuación, calcular el grado y las medidas de centralidad. Para cada grafo completar e interpretar la siguiente tabla. Interpretar los resultados. ¿Coincidencia que Estados Unidos sea otro país inmerso en conflictos de mucho tipo? Ahora sabemos que no.

g7_1 <- make_star(n = 5, mode = "undirected")
g7_2 <- make_tree(n = 5, children = 2, mode = "undirected")
g7_3 <- graph_from_literal(1-2-3-4-5)
g7_4 <- make_ring(n = 5, directed = FALSE)

par(mfrow = c(2, 2), mar = c(1, 3, 1, 3)  
)
plot(g7_1, layout = layout_with_fr, main = "Primer grafo")
plot(g7_2, layout = layout_with_fr, main = "Segundo grafo")
plot(g7_3, layout = layout_with_fr, main = "Tercer grafo")
plot(g7_4, main = "Cuarto grafo")

# Lista de todos los grafos para automatizar el proceso
grafos7 <- list(g7_1 = g7_1, g7_2 = g7_2, g7_3 = g7_3, g7_4 = g7_4)

medidas <- lapply(grafos7, function(g) {
  list(
    grado = igraph::degree(g),
    cercania = igraph::closeness(g, normalized = TRUE),
    intermediacion = igraph::betweenness(g, normalized = TRUE),
    eigen = igraph::eigen_centrality(g)$vector
  )
})

resumen <- lapply(medidas, function(med) {
  lapply(med, function(x) {
    c(mean = mean(x),
      sd = sd(x) )
  })
})

tabla <- data.frame(
  Grafo = c("Primer", "Primer", "Segundo", "Segundo", "Tercero", "Tercero", "Cuarto", "Cuarto"),
  Medida = c("Media", "DE", "Media", "DE", "Media", "DE", "Media", "DE"),
  Grado = c(1.6, 1.342, 1.6, 0.894, 1.6, 0.547, 2, 0),
  Cercania = c(0.657, 0.192, 0.582, 0.147, 0.522, 0.118, 0.667, 0),
  Intermediacion = c(0.2, 0.447, 0.266, 0.384, 0.333, 0.312, 0.166, 0),
  Propia = c(0.6, 0.224, 0.652, 0.232, 0.746, 0.231, 1, 0)
)

kable(tabla, align = "c") %>%
  collapse_rows(columns = 1) %>%
  kable_styling(
   # full_width = FALSE,
    bootstrap_options = c("bordered")
  )

Grafo	Medida	Grado	Cercania	Intermediacion	Propia
Primer	Media	1.600	0.657	0.200	0.600
	DE	1.342	0.192	0.447	0.224
Segundo	Media	1.600	0.582	0.266	0.652
	DE	0.894	0.147	0.384	0.232
Tercero	Media	1.600	0.522	0.333	0.746
	DE	0.547	0.118	0.312	0.231
Cuarto	Media	2.000	0.667	0.166	1.000
	DE	0.000	0.000	0.000	0.000

Primero, hablemos del grado. Todos los grafos considerados tienen el mismo grado medio de 1.6 a excepción del último que tiene un valor de 2, debido a que este último es un ciclo, por lo que todos los nodos tienen grado 2 (se probó en un punto anterior). La diferencia en este criterio es la desviación estándar, en la que el primer grafo tipo estrella presenta grados más dispersos.

Se construyeron las centralidad normalizadas para facilitar su interpretación.

La centralidad por cercanía más alta la obtiene el cuarto grafo, aspecto natural teniendo en cuenta que esta red es la que presenta una menor distancia geodésica promedio. Cabe resaltar que para todos los grafos se tiene que en promedio sus nodos son importantes gracias a que tienen varios nodos cerca (todos los valores superaron el 0.5).

La centralidad media por intermediación toma valores más pequeños, el más alto de 0.33 alcanzado por el tercer grafo, el cual tiene la menor dispersión aparte del cuarto grafo. Vemos que la mayor dispersión la obtiene el primer grafo porque solo uno de sus nodos está entre algún par de vértices; de hecho, está entre todos los pares.

Finalmente, teniendo en cuenta que se vio que al rededor de todas las redes en promedio los nodos se consideran centrales, la centralidad propia media también toma valores grandes. La red con mayor nivel de importancia promedio de sus nodos bajo este criterio fue la cuarta.

En términos generales: las medidas en promedio son iguales para el grafo cuatro dada su estructura circular y el tercer grafo presentó mayor dispersión en todas sus medidas

8. Considere los datos relacionales acerca de los conflictos internacionales del archivo conflict.RData después de simetrizarla débilmente y remover los nodos aislados:

   1. Hacer una visualización decorada de la red.
   2. Caracterizar local y estructuralmente la red, en términos de la distancia, la centralidad, la cohesión, la conectivaidad, y el agrupamento. Utilizar todas las métricas disponbles.
   3. Interpretar los resultados.

Solución Punto 8:

La simetrización de la matriz permite que podamos contar ahora con una red no dirigida.

Nota: La simetrización también hace que ahora la red se considere no ponderada, por eso aunque se diga que es ponderada, los pesos son todos iguales a 1.

load("conflict.RData")
# Simetrización débil de la matriz de adyacencia
Y8 <- symmetrize(dat$Y, rule = "weak")
dimnames(Y8) <- dimnames(dat$Y)

g8_conflict <- graph_from_adjacency_matrix(Y8, mode = "undirected", weighted = TRUE, diag = FALSE)
vertex_attr(g8_conflict) <- as.data.frame(dat$X)
V(g8_conflict)$name <- rownames(dat$X)

g8_conflict <- delete_vertices(g8_conflict, which(igraph::degree(g8_conflict)==0))
#V(g8_conflict)$name <- rownames(dat$Y)

vcount(g6_conflict); vcount(g8_conflict)

## [1] 130

## [1] 91

De los 130 actores que se tenían inicialmente, se va a trabajar con 91 que no son aislados.

ecount(g6_conflict); ecount(g8_conflict)

## [1] 203

## [1] 160

En términos de las aristas, se remueven un total de 43.

Al igual que en el punto (6) se categoriza la variable ‘polity’ para facilitar los análisis

# Categorización de la variable 'polity'
V(g8_conflict)$polity <- ifelse(V(g8_conflict)$polity>0, "Democracia", "Autocracia")

# Identificación del tipod e relación entre países
edges <- ends(g8_conflict, E(g8_conflict))
pol_from <- V(g8_conflict)$polity[match(edges[,1], V(g8_conflict)$name)]
pol_to <- V(g8_conflict)$polity[match(edges[,2], V(g8_conflict)$name)]
E(g8_conflict)$tipo_regimen <- ifelse(pol_from == "Autocracia" & pol_to == "Autocracia","Autocracia-Autocracia",
                                      ifelse(pol_from == "Democracia" & pol_to == "Democracia","Democracia-Democracia",
                                             "Mixta")
                                      )

g8_tbl <- as_tbl_graph(g8_conflict)
set_graph_style(plot_margin = margin(1,1,1,1))

# Layout
set.seed(1)
layout8 <- create_layout(g8_tbl, layout = 'igraph', algorithm = "fr")

ggraph(layout8) +
  geom_edge_link(aes(width = weight,
                     color = tipo_regimen),
                 alpha = 0.5,
                 end_cap = circle(2, "mm")) +
  geom_node_point(aes(size = log(gdp),
                      color = polity)) +
  geom_node_label(
    aes(label = ifelse(name %in% c("IRQ","USA"), name, NA),
        fill = polity),
    linewidth = 0,
    size = 4,
    show.legend = FALSE
  ) +
  scale_edge_width(range = c(1, 1), guide = "none") +
  scale_edge_color_manual(
  name = "Tipo de relación",
  values = c(
    "Autocracia-Autocracia" = "dodgerblue4",  # azul coherente
    "Democracia-Democracia" = "darkgreen",    # verde coherente
    "Mixta" = "mediumpurple1"                 # mezcla visual azul+verde
  )
)+
  scale_size(range = c(1, 7), guide = "none") +
  scale_color_manual(
    name = "Regimen politico",
    values = c("Democracia" = "green",
               "Autocracia" = "dodgerblue3")
  ) +
  scale_fill_manual(
    values = c("Democracia" = "green",
               "Autocracia" = "dodgerblue3"
    ),
    guide = "none"
  ) +
  theme_void() +
  theme(legend.position = "right") +
  labs(
    title = "Red para los conflictos entre paises",
    subtitle = "El tamano de los vertices representa el log(PIB)"
  )

Note que hay heterogeneidad entre los tipos de relación identificados por los colores de los enlaces. Es decir, salvo ciertos casos que se ven en el centro de la red, parece que los conflictos no están dirigidos necesariamente hacia países de un régimen político distinto.

Caracterización en términos de distancia: Se calculan el diámetro de la red y la distancia geodésica promedio. Adicionalmente se presenta la distribución de distancias geodésicas de la red.

diameter(g8_conflict)

## [1] 9

(d <- get_diameter(g8_conflict))

## + 10/91 vertices, named, from b11620d:
##  [1] BOT NAM UGA SUD USA JOR BNG MYA THI CAM

mean_distance(g8_conflict)

## [1] 3.65248

# Visualización
caminos <- distance_table(g8_conflict)$res
names(caminos) <- 1:length(caminos)
barplot(
  prop.table(caminos), 
  xlab = "Distancia geodésica", 
  ylab = "F. Relativa", 
  border = "grey", 
  col = "grey", 
  main = "Distribución de distancias geodésicas"
)

La distancia geodésica entre dos vértices es la longitud del camino más corto que los conecta. Un diámetro de la red igual a 9, como valor máximo de estas distancias, nos dice que los dos países más alejados necesitan 9 aristas para conectarse por el camino más corto. Para una red bastante dispersa que solo cuenta con 91 actores este valor es acorde.

La distancia geodésica promedio es 3.6, indicando que en promedio, dos países se conectan mediante aproximadamente tres o cuatro intermediarios. Aunque el diámetro de la red es 9, este valor refleja solo algunos casos extremos, mientras que la mayoría de los países se encuentran relativamente cercanos entre sí dentro de la red; afirmación que se soporta con el gráfico de la distribución de estas distancias, el cual presenta asimetría a derecha.

Caracterización en términos de centralidad: Se calcula el top 5 de países más centrales según los criterios vistos en clase.

# top 5
cc <- igraph::closeness(graph = g8_conflict, normalized = T)
head(sort(cc, decreasing = T), n = 5)

## BUI LES MAL MZM SAF 
##   1   1   1   1   1

# top 5
bc <- igraph::betweenness(graph = g8_conflict, normalized = T)
head(sort(bc, decreasing = T), n = 5)

##        IRQ        USA        JOR        VEN        SIE 
## 0.29460000 0.29159005 0.23590191 0.13383271 0.09171492

# top 5
ec <- igraph::eigen_centrality(graph = g8_conflict)$vector
head(sort(ec, decreasing = T), n = 5)

##       IRQ       JOR       USA       CAN       SYR 
## 1.0000000 0.9522752 0.4660763 0.3671518 0.3570798

En términos de la cercanía, los países Burundi, Lesotho, Malasia, Mozambique y Sudáfrica son 100% importantes dada su cercanía con muchos otros nodos. Esto es, son países que están en medio de actividades conflictivas.

Los países Iraq, Estados Unidos, Jordania, Unión Soviética y Sierra Leona son importantes cuando el criterio es la presencia de estos entre los demás pares de países. Esto se traduce a que son países que pueden estar actuando como intermediarios en la estructura de conflictos, conectando conflictos distintos entre si. Note que estos valores no son tan altos como en el caso anterior.

Por último, Iraq, Jordania, Estados Unidos, Canadá y Siria son los países más centrales gracias a que cuentan con vecinos igualmente centrales. Es lógico, dados los análisis previos, que estos estén involucrados en conflictos con países muy activos, porque de por sí han tenido conflictos entre ellos.

En términos generales note que Iraq, Jordania y Estados Unidos son los tres actores que se consideran centrales bajo más de un criterio.

Caracterización en términos de cohesión: En este apartado, se encuentran los ‘cliques’ más grandes de la red junto con su frecuencia y se llevan a cabo los censos de estos diádicos y triádicos. Además de calcular la densidad y transitividad de la red.

# frecuencias de clanes
table(sapply(X = cliques(graph = g8_conflict), FUN = length)) # de tamaño 4, hay 9 cliques

## 
##   1   2   3   4 
##  91 160  67   9

largest_cliques(graph = g8_conflict)

## [[1]]
## + 4/91 vertices, named, from b11620d:
## [1] TUR IRQ JOR GRC
## 
## [[2]]
## + 4/91 vertices, named, from b11620d:
## [1] TUR IRQ JOR SYR
## 
## [[3]]
## + 4/91 vertices, named, from b11620d:
## [1] UGA RWA DRC ANG
## 
## [[4]]
## + 4/91 vertices, named, from b11620d:
## [1] PRK CHN ROK JPN
## 
## [[5]]
## + 4/91 vertices, named, from b11620d:
## [1] QAT IRQ SAU JOR
## 
## [[6]]
## + 4/91 vertices, named, from b11620d:
## [1] USA IRQ JOR CAN
## 
## [[7]]
## + 4/91 vertices, named, from b11620d:
## [1] USA IRQ JOR SYR
## 
## [[8]]
## + 4/91 vertices, named, from b11620d:
## [1] SPN CAN JOR IRQ
## 
## [[9]]
## + 4/91 vertices, named, from b11620d:
## [1] SYR IRQ JOR ISR

# número clan
clique_num(graph = g8_conflict)

## [1] 4

Hay 9 clanes con el tamaño más grande de la red de conflictos, los cuales se comprenden de 4 actores. Note que Iraq y Jordania hacen parte de 7 de estos nueve clanes, es decir que son países que han tenido varios enfrentamientos en los que todos los países involucrados se atacan unos con otros.

El número clan, por esta razón, es solamente de 4, indicando que los conflictos tienden a fomentarse conjuntamente en grupos pequeños de países, no en grandes grupos totalmente conectados. Lógicamente, sería muy extraño tener un conflicto con 20 países involucrados, en el que todos se atacaran.

dyad_census(g8_conflict)

## $mut
## [1] 160
## 
## $asym
## [1] 0
## 
## $null
## [1] 3935

triad_labels <- c("003",  "012",  "102",  "021D", "021U", "021C",
                  "111D", "111U", "030T", "030C", "201",
                  "120D", "120U", "120C", "210",  "300")
censo3 <- triad_census(g8_conflict)
names(censo3) <- triad_labels

Del censo de díadas obtenemos que hay 160 parejas, lo cual coincide con el número de cliques de 2 nodos visto arriba.

# número de redes nulas
B <- 1000

# grados de salida y entrada del grafo observado
deg8 <- igraph::degree(g8_conflict, mode = "all")

# matriz para guardar los censos triádicos bajo el modelo nulo
null_counts <- matrix(
  NA_real_,
  nrow = B,
  ncol = length(triad_labels),
  dimnames = list(NULL, triad_labels)
)

# redes nulas dirigidas con la misma secuencia de grados
set.seed(123)
for (b in seq_len(B)) {
  g_null <- igraph::sample_degseq(
    out.deg = deg8, # Solo este porque la red es no dirigida
    method  = "edge.switching.simple"
  )
  
  null_counts[b, ] <- igraph::triad_census(g_null)
}

# media y desviación estándar bajo el modelo nulo
mu_null <- colMeans(null_counts)
sd_null <- apply(X = null_counts, MARGIN = 2, FUN = sd)

# puntuaciones Z
# puntuaciones Z
z_scores <- ifelse(
  sd_null > 0,
  (censo3 - mu_null) / sd_null,
  NA_real_
)

# versión normalizada
z_scores_norm <- z_scores / sqrt(sum(z_scores^2, na.rm = TRUE))

# tabla de resultados
triad_z <- data.frame(
  triad   = triad_labels,
  obs     = as.numeric(censo3),
  mean    = mu_null,
  sd      = sd_null,
  z       = z_scores,
  z_norm  = z_scores_norm,
  row.names = NULL
)

# tabla
triad_z[c(1, 3, 11, 16),]

##    triad    obs       mean        sd         z z_norm
## 1    003 108408 108438.359  5.918543 -5.129472   -0.5
## 3    102  11981  11889.923 17.755628  5.129472    0.5
## 11   201   1029   1120.077 17.755628 -5.129472   -0.5
## 16   300     67     36.641  5.918543  5.129472    0.5

El censo de tríadas nos muestra que las 108408 tríadas nulas y las 1029 tríadas de tipo 201 aparecieron menos de lo esperado respecto de un modelo nulo, mientras que las triadas de tipo 102 y las tríadas completas se observaron más de lo esperado en la red de conflictos, lo que significa que podemos considerar que en la red hay un número significativo de estas dos últimas configuraciones.

# Densidad
edge_density(graph = g8_conflict)

## [1] 0.03907204

# Puedo hacer la densidad para los subgrafos inducidos por Democracia y Autocracia
g8_demo <- induced_subgraph(g8_conflict, vids = V(g8_conflict)[polity == "Democracia"])
g8_auto <- induced_subgraph(g8_conflict, vids = V(g8_conflict)[polity == "Autocracia"])

edge_density(graph = g8_demo)

## [1] 0.02806122

edge_density(graph = g8_auto)

## [1] 0.06852497

La densidad total de la red de conflictos entre países es de 0.039, valor que rectifica el hecho de que esta red es desconectada y está bastante lejos de volverse un clan.

Se decidió además calcular la densidad para los subgrafos inducidos por la variable nodal ‘polity’. Note que el subgrafo generado por los países democráticos es menos denso que el grafo total, lo que se interpreta como que entre países de este tipo, hay un poco menos de relaciones conflictivas. Por otro lado, el grafo generado por los países autocráticos si es levemente más denso que el general.

# Transitividad
transitivity(graph = g8_conflict, type = "global")

## [1] 0.1634146

# intransitividad local
hist(transitivity(g8_conflict, type = "local"), xlab = "Transitividad local", ylab = "Frecuencia",
     main = "Distribución de la transitividad local")

La transitividad global toma un valor de 0.163. Esta se puede interpretar como que las triplas tienen un 16,3% de probabilidad de formar triángulos, mostrando un agrupamiento global relativamente bajo.

En términos de la transitividad local vemos que hay dos grupos de nodos: unos con nivel de agrupamiento bajo y otros con un agrupamiento alto. Este último grupo, hay que resaltar, es menos frecuente que el primero.

Caracterización en términos de conectividad: Se estudia la conectividad por vértices y aristas, la componente gigante y sus características de conectividad.

# componente gigante
is_connected(g8_conflict)

## [1] FALSE

vertex_connectivity(g8_conflict)

## [1] 0

edge_connectivity(g8_conflict)

## [1] 0

Como podíamos intuir, esta red no está conectada y sus conectividades nodal y de aristas toman el valor de cero.

Para estudiar mejor la red, procedemos a encontrar la componente gigante de esta.

# componentes
componentes <- decompose(g8_conflict)
length(componentes)

## [1] 5

table(sapply(X = componentes, FUN = vcount)) # 83 es el tamaño de la componente gigante

## 
##  2 83 
##  4  1

g8_cg <- decompose(g8_conflict)[[1]]
# conectividad nodal
vertex_connectivity(g8_cg)

## [1] 1

edge_connectivity(g8_cg)

## [1] 1

articulation_points(g8_cg)

## + 22/83 vertices, named, from b4e7e21:
##  [1] USA GHA SIE NIG NIR UGA NAM DRC TUR GRC NIC COL VEN CHN THI MYA BNG IND SAU
## [20] SEN UKG HAI

length(articulation_points(g8_cg))/vcount(g8_cg)

## [1] 0.2650602

La descomposición de la red muestra que hay 5 componentes: cuatro de ellas de tamaño 2 y una con tamaño de 83 nodos, la cual es nuestra componente gigante que representa el 91.2% de los nodos de la red no dirigida original.

La pregunta ahora es: ¿Es esta componente sensible a remover nodos de ella? La respuesta es si, ya que la conectividad nodal y de aristas es de 1, es decir, que basta con eliminar un nodo o una arista de la componente gigante para que esta se desconecte. Además, los puntos de articulación son el 26.5% del total de nodos (22 de 83). Algunos de los países que hacen parte de este grupo son: Estados Unidos, Nigeria, Colombia y China.

Caracterización en términos del agrupamiento: Por último, llevaremos a cabo la detección de comunidades de la red para encontrar vértices homogéneos.

Se van a correr varios algoritmos de agrupamiento, tanto jerárquico como de partición espectral, y se van a comparar en términos del que alcance la mayor modularidad.

# algoritmos
clust8_fast_greedy <- cluster_fast_greedy(g8_conflict)
clust8_leading_eigen <- cluster_leading_eigen(g8_conflict)
clust8_walktrap <- cluster_walktrap(g8_conflict)
clust8_louvain <- cluster_louvain(g8_conflict)
clust8_label_prop <- cluster_label_prop(g8_conflict)
#clust8_spinglass <- cluster_spinglass(g8_conflict)
clust8_optimal <- cluster_optimal(g8_conflict)
clust8_infomap <- cluster_infomap(g8_conflict)

# gráficos

igraph_options(vertex.size = 10, vertex.frame.color = "black")
par(mfrow = c(3, 3), mar = c(0, 0, 2, 0))

# extraer coordenadas del layout de ggraph
layout8_mat <- as.matrix(layout8[, c("x", "y")])
 
plot(g8_conflict, vertex.label = NA, layout = layout8_mat, vertex.color = clust8_fast_greedy$membership, main = paste0("fast greedy: ", "Mod = ", round(modularity(clust8_fast_greedy), 4)))
plot(g8_conflict, vertex.label = NA, layout = layout8_mat, vertex.color = clust8_leading_eigen$membership, main = paste0("leading eigen: ", "Mod = ", round(modularity(clust8_leading_eigen), 4)))
plot(g8_conflict, vertex.label = NA, layout = layout8_mat, vertex.color = clust8_walktrap$membership, main = paste0("walktrap: ", "Mod = ", round(modularity(clust8_walktrap), 4)))
plot(g8_conflict, vertex.label = NA, layout = layout8_mat, vertex.color = clust8_louvain$membership, main = paste0("louvain: ", "Mod = ", round(modularity(clust8_louvain), 4)))
plot(g8_conflict, vertex.label = NA, layout = layout8_mat, vertex.color = clust8_label_prop$membership, main = paste0("label prop: ", "Mod = ", round(modularity(clust8_label_prop), 4)))
plot(g8_conflict, vertex.label = NA, layout = layout8_mat, vertex.color = clust8_optimal$membership, main = paste0("optimal: ", "Mod = ", round(modularity(clust8_optimal), 4)))
plot(g8_conflict, vertex.label = NA, layout = layout8_mat, vertex.color = clust8_infomap$membership, main = paste0("infomap: ", "Mod = ", round(modularity(clust8_infomap), 4)))

Note que para la red de conflictos maximiza su modularidad bajo al algoritmo de agrupamiento ‘optimal’. Este calcula la estructura de comunidades óptimas de un grafo maximizando la medida de modularidad sobre todas las posibles particiones. Ahora, vamos a analizar esta clasificación más a fondo

# agregar comunidades como atributo de nodo
g8_tbl <- as_tbl_graph(g8_conflict)
g8_tbl <- g8_tbl %>%
  mutate(comm_optimal = factor(clust8_optimal$membership))

set_graph_style(plot_margin = margin(1,1,1,1))

set.seed(1)
layout8 <- create_layout(g8_tbl, layout = "igraph", algorithm = "fr")

ggraph(layout8) +
  geom_edge_link(aes(width = weight),
    alpha = 0.4,
    colour = "gray62",
    end_cap = circle(2, "mm")
  ) +
  geom_node_point(
    aes(color = comm_optimal,
        size = log(gdp) )
  ) +
  scale_edge_width(range = c(0, 2.5), guide = "none") +
  scale_size(range = c(1, 7), guide = "none") +
  theme_void() +
  labs(
    title = "Agrupamiento de la red de conflictos mediante el algoritmo cluster_optimal",
    subtitle = "El tamaño de los vertices representa el log(PIB)\n Modularidad = 0.5664",
    color = "Comunidad"
  )

is_hierarchical(clust8_optimal)

## [1] FALSE

sizes(clust8_optimal)

## Community sizes
##  1  2  3  4  5  6  7  8  9 10 11 12 
##  7 26 13  9  9  8  2  8  3  2  2  2

El algoritmo está clasificando a los actores en 12 clusters. El cluster más grande lo conforman 26 países, los cuales según la gráfica se ubican hacia la parte más concurrida de la red. El cluster que le sigue en tamaño tiene la mitad de nodos del anterior y lo podemos ubicar en la parte inferior derecha de la red.

Que hayan tantos clusters nos habla de una heterogeneidad de los datos, aspecto ya estudiado en los análisis previos.

pol <- as.numeric(as.factor(vertex_attr(graph = g8_conflict, name = "polity")))
table(pol)

## pol
##  1  2 
## 42 49

compare(comm1 = pol, comm2 = clust8_optimal$membership, method = "rand")

## [1] 0.5203907

Véase además que, comparando esta partición con la inducida por la política del país (Democracia o Autocracia) encontramos que ambas agrupaciones condicionen en aproximadamente la mitad de par de nodos.

9. Sintetizar y replicar los Capítulos 6, 8, y 9 de Luke, D. A. (2015).

CAPÍTULO 6: GRÁFICOS DE RED AVANZADOS

En este capítulo se exploran técnicas especializadas para graficar redes.

Gráficos interactivos:

Mucha veces las visualizaciones estáticas clásicas de R no son suficientes para el usuario, por eso se ve la necesidad de generar representaciones interactivas.

El paquete igraph cuenta con la función tkplot() que permite hacer gráficos interactivos de redes simples mediante una ventana gráfica emergente.

#library(devtools)
#install_github("DougLuke/UserNetR")
suppressMessages(suppressWarnings(library(statnet))) 
suppressMessages(suppressWarnings(library(intergraph)))
suppressMessages(suppressWarnings(library(UserNetR)))

data(Bali): Esta red muestra las interacciones entre el grupo terrorista Jemaah Islamiyah que llevó a cabo los atentados en Bali en 2002. Se cuenta con 17 vértices (miembros del grupo terrorista) y 126 aristas (contactos entre los actores).

data(Bali)
iBali <- asIgraph(Bali)
Coord <- tkplot(iBali, vertex.size=3,
vertex.label=V(iBali)$role,
vertex.color="darkgreen")
# Edit plot in Tk graphics window before
# running next two commands.
MCoords <- tkplot.getcoords(Coord)

## Warning: `tkplot.getcoords()` was deprecated in igraph 2.0.0.
## ℹ Please use `tk_coords()` instead.
## This warning is displayed once per session.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.

plot(iBali, layout=MCoords, vertex.size=5,
vertex.label=NA, vertex.color="lightblue")

El paquete networkD3, por su lado, permite construir gráficos interactivos que se desplieguen en páginas HTML propiamente (el tipo de visualización apropiado para este Taller). La forma más sencilla de hacerlo es mediante simpleNetwork().

library(networkD3)

## Warning: package 'networkD3' was built under R version 4.5.3

src <- c("A","A","B","B","C","E")
target <- c("B","C","C","D","B","C")
net_edge <- data.frame(src, target)
simpleNetwork(net_edge)

Si se quiere algo más sofisticado se puede optar por la función forceNetwork(). A esta función se le deben dar de argumentos dos dataframes: uno con la lista de aristas y otro con los ID de los nodos y sus propiedades. La siguiente visualización nos permite visualizar la red coloreada de acuerdo a los roles de los actores.

# Este código no se corre porque estaba haciendo que se me dañara el R :b

iBali_edge <- get.edgelist(iBali)
iBali_edge <- iBali_edge - 1
iBali_edge <- data.frame(iBali_edge)
iBali_nodes <- data.frame(NodeID=as.numeric(V(iBali)-1),
Group=V(iBali)$role,
Nodesize=(degree(iBali)))
forceNetwork(Links = iBali_edge, Nodes = iBali_nodes,
Source = "X1", Target = "X2",
NodeID = "NodeID",Nodesize = "Nodesize",
radiusCalculation="Math.sqrt(d.nodesize)*3",
Group = "Group", opacity = 0.8,
legend=TRUE)

Otro paquete para efectuar estos gráficos es visNetwork. Este último es bastante amplio en lo que el usuario puede controlar de la apariencia de la red.

# Este código no se corre porque estaba haciendo que se me dañara el R :b

library(visNetwork)
iBali_edge <- igraph::as_edgelist(iBali)
iBali_edge <- data.frame(from = iBali_edge[,1], to = iBali_edge[,2])
iBali_nodes <- data.frame(id = as.numeric(V(iBali)))
visNetwork(iBali_nodes, iBali_edge, width = "100%")

iBali_nodes$group <- V(iBali)$role
iBali_nodes$value <- degree(iBali)
net <- visNetwork(iBali_nodes, iBali_edge, width = "100%",legend=TRUE)
visOptions(net,highlightNearest = TRUE)

Diagramas especializados

Diagramas de Arco: Estos son útiles cuando el posicionamiento de los nodos es menos importante que el patrón de las aristas.

#library(devtools)
#install_github("gastonstat/arcdiagram")
library(arcdiagram)

data(Simpsons)
iSimp <- asIgraph(Simpsons)
simp_edge <- as_edgelist(iSimp)
arcplot(simp_edge)

Este diagrama se puede organizar para resaltar características de la red. Se organizan en grupos identificados por sus colores, por ejemplo.

s_grp <- V(iSimp)$group
s_col = c("#a6611a", "#dfc27d","#80cdc1","#018571")
cols = s_col[s_grp]
node_deg <- igraph::degree(iSimp)
arcplot(simp_edge, lwd.arcs=2, cex.nodes=node_deg/2, labels=V(iSimp)$vertex.names,
        col.labels="darkgreen",font=1, pch.nodes=21,line=1,col.nodes = cols, bg.nodes = cols, show.nodes = TRUE)

Diagramas de Acordes (Chord): Los diagramas de acordes son un tipo especializado de gráfico que utiliza una disposición circular para mostrar las interrelaciones entre los nodos de una red. Son de especial interés cuando las redes son ponderadas o dirigidas.

La base FIFA_Nether es una red de del equipo de fútbol de la Copa Mundial de los Países Bajos 2010. Esta es de tipo dirigida, que se puede ver como las personas que hacen los pases y los que los reciben.

library(circlize)
data(FIFA_Nether)
FIFAm <- as.sociomatrix(FIFA_Nether,attrname='passes') 
names <- c("GK1","DF3","DF4","DF5","MF6", "FW7","FW9","MF10","FW11","DF2","MF8")
rownames(FIFAm) = names
colnames(FIFAm) = names

# No considerar los pases menores a 10
FIFAm[FIFAm < 10] <- 0
FIFAm

##      GK1 DF3 DF4 DF5 MF6 FW7 FW9 MF10 FW11 DF2 MF8
## GK1    0  42  67  21   0  27   0    0    0  17   0
## DF3   30   0  44  14  42  15   0    0   10  36  29
## DF4   38  43   0  57  18  11   0   21    0   0  28
## DF5    0  14  47   0  11  50  20   40    0   0  42
## MF6    0  28  25  10   0  41  28   37   14  34  21
## FW7    0  12   0  21  21   0  15   33    0  25  18
## FW9    0   0   0   0   0  12   0   31   16   0   0
## MF10   0  11  11  22  43  29  20    0   28  13  21
## FW11   0   0   0   0   0   0  11   15    0  21  12
## DF2   29  38   0   0  45  38  10   18   26   0  15
## MF8   12  25  26  38  23  13  12   32   11  24   0

chordDiagram(FIFAm)

Se quiere ahora que los jugadores que juegan en la misma posición compartan color

grid.col <- c("#AA3939",rep("#AA6C39",4),rep("#2D882D",3),rep("#226666",3))
chordDiagram(FIFAm,directional = TRUE,grid.col = grid.col,
             order=c("GK1","DF2","DF3","DF4","DF5","MF6","MF8","MF10","FW7",
                     "FW9","FW11"))

Mapas de calor (Heatmaps): Gráficos especialmente valorados para redes ponderadas.

data(FIFA_Nether)
FIFAm <- as.sociomatrix(FIFA_Nether,attrname='passes')
colnames(FIFAm) <- c("GK1","DF3","DF4","DF5", "MF6","FW7","FW9","MF10",
                     "FW11","DF2","MF8")
rownames(FIFAm) <- c("GK1","DF3","DF4","DF5", "MF6","FW7","FW9","MF10",
                     "FW11","DF2","MF8")

palf <- colorRampPalette(c("#669999", "#003333"))
# Entre más oscutro el color, más pases realizados
heatmap(FIFAm[,11:1],Rowv = NA,Colv = NA,col = palf(60), scale="none", margins=c(11,11) )

Gráficos para redes usando otros paquetes de R

Aunque ggplot2 no está diseñado para manejar todos los requisitos de un paquete completo de visualización de red, algunas de sus capacidades gráficas avanzadas se pueden utilizar para crear rutinas especializadas de trazado de red.

La siguiente función edgeMaker se usa para crear gráficos atractivos y funcionales de redes dirigidas usando bordes curvos de intensidad

# Función que crea aristas curvas entre díadas conectadas
edgeMaker <- function(whichRow,len=100, curved = TRUE){
  fromC <- layoutCoordinates[adjacencyList[whichRow,1],]
  toC <- layoutCoordinates[adjacencyList[whichRow,2],]
  graphCenter <- colMeans(layoutCoordinates)
  bezierMid <- c(fromC[1], toC[2])
  distance1 <- sum((graphCenter - bezierMid)^2)
  if(distance1 < sum((graphCenter - c(toC[1], fromC[2]))^2)){
    bezierMid <- c(toC[1], fromC[2])
    }
  bezierMid <- (fromC + toC + bezierMid) / 3
  if(curved == FALSE){bezierMid <- (fromC + toC) / 2}
  edge <- data.frame(bezier(c(fromC[1], bezierMid[1], toC[1]), c(fromC[2], bezierMid[2],
                                                                 toC[2]), evaluation = len))
  edge$Sequence <- 1:len
  edge$Group <- paste(adjacencyList[whichRow, 1:2], collapse = ">")
  return(edge)
  }

library(sna)
library(ggplot2)
library(Hmisc)

data(FIFA_Nether)
fifa <- FIFA_Nether
# Los datos deben estar en formato de matriz de aristas
fifa.edge <- as.edgelist.sna(fifa,attrname='passes')
fifa.edge <- data.frame(fifa.edge)
names(fifa.edge)[3] <- "value"
# No se tienen en cuenta los pases menores a 10
fifa.edge <- fifa.edge[fifa.edge$value > 9,]
adjacencyList <- fifa.edge

layoutCoordinates <- gplot(network(fifa.edge)) # función de sna

allEdges <- lapply(1:nrow(fifa.edge), edgeMaker, len = 500, curved = TRUE)
allEdges <- do.call(rbind, allEdges)

new_theme_empty <- theme_bw()
new_theme_empty$line <- element_blank()
new_theme_empty$rect <- element_blank()
new_theme_empty$strip.text <- element_blank()
new_theme_empty$axis.text <- element_blank()
new_theme_empty$plot.title <- element_blank()
new_theme_empty$axis.title <- element_blank()
new_theme_empty$plot.margin <- structure(c(0,0,-1,-1), unit = "lines",
                                         valid.unit = 3L, class = "unit")

zp1 <- ggplot(allEdges)
zp1 <- zp1 + geom_path(aes(x = x, y = y, group = Group, colour=Sequence, size=-Sequence))

## Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
## ℹ Please use `linewidth` instead.
## This warning is displayed once per session.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.

zp1 <- zp1 + geom_point(data = data.frame(layoutCoordinates),
                        aes(x = x, y = y),
                        size = 4, pch = 21,
                        colour = "black", fill = "gray")
zp1 <- zp1 + scale_colour_gradient(low = gray(0),
                                   high = gray(9/10),
                                   guide = "none")
zp1 <- zp1 + scale_size(range = c(1/10, 1.5),
                        guide = "none")
zp1 <- zp1 + new_theme_empty
print(zp1)

CAPÍTULO 8: SUBGRUPOS

El objetivo de este capítulo es cubrir varias técnicas presentes en R para identificar y examinar subgrupos que redes sociales grandes puedan tener.

Cohesión social: Los subgrupos cohesivos son conjuntos de actores que están unidos por vínculos frecuentes, fuertes y directos. ¿Como podemos identificarlos?

• Cliques: Un clan es un subgrafo máximo completo; es decir, un subconjunto de nodos que tienen todas las aristas posibles entre ellos. El número clique es el tamaño del clan más grande. Los cliques, sin embargo, presentan dos desventajas: son un subgrupo cohesivo muy conservativo y no tan comunes en redes grandes.

• k-Cores: Un núcleo-k es un subgrafo máximo donde cada vértice está conectado a al menos k otros vértices en el subgrafo. Además, están anidados, no se superponen y son fáciles de identificar. Su análisis empieza identificando todo el núcleo-l, para luego examinar visualmente sus estructuras.

  Un coreness alto implica que el nodo está integrado en el núcleo de la red, mientras que un coreness bajo nos dice que el nodos es más periférico.

Ejemplo de juguete: En el siguiente ejemplo se ven dos cliques: A, B, C, D y E, F, G

clqexmp <- graph_from_literal(A:B:C:D--A:B:C:D,D-E,E-F-G-E)
clique_num(clqexmp)

## [1] 4

cliques(clqexmp, min = 3)

## [[1]]
## + 3/7 vertices, named, from ed60ad7:
## [1] B C D
## 
## [[2]]
## + 3/7 vertices, named, from ed60ad7:
## [1] E F G
## 
## [[3]]
## + 3/7 vertices, named, from ed60ad7:
## [1] A B C
## 
## [[4]]
## + 4/7 vertices, named, from ed60ad7:
## [1] A B C D
## 
## [[5]]
## + 3/7 vertices, named, from ed60ad7:
## [1] A B D
## 
## [[6]]
## + 3/7 vertices, named, from ed60ad7:
## [1] A C D

max_cliques(clqexmp, min = 3)

## [[1]]
## + 3/7 vertices, named, from ed60ad7:
## [1] E F G
## 
## [[2]]
## + 4/7 vertices, named, from ed60ad7:
## [1] A B D C

largest_cliques(clqexmp)

## [[1]]
## + 4/7 vertices, named, from ed60ad7:
## [1] A B D C

Ejemplo: Red DHHS:

data(DHHS)
iDHHS <- asIgraph(DHHS) # convertirlo a objeto igraph
iDHHS <- subgraph_from_edges(iDHHS,E(iDHHS)[collab > 2]) # subgrafo

La base DHHS muestra los vínculos de colaboración entre expertos en control del tabaco que trabajan en varios institutos y agencias dentro del Departamento de Salud y Servicios Humanos en 2005. Identificamos la estructura k-core en la red como

core <- coreness(iDHHS)
table(core)

## core
##  1  2  3  4  5  6 
##  7  6  2  5  2 26

maxCoreness <- max(core)
maxCoreness

## [1] 6

Para mejorar la interpretación se etiqueta a los nodos con su pertenencia a los distintos núcleos obtenidos. Naturalmente, los ndoos más conectados (coreness=6) están en el centro de la red, y los que tienen coreness=1 están en las periferias.

set.seed(1)
colors <- rainbow(maxCoreness)
op <- par(mar = rep(0, 4))
plot(iDHHS,vertex.label=core, vertex.color=colors[core])

par(op)

V(iDHHS)$name <- core
V(iDHHS)$color <- colors[core]
iDHHS1_6 <- iDHHS
iDHHS2_6 <- induced_subgraph(iDHHS, vids=which(core > 1))
iDHHS3_6 <- induced_subgraph(iDHHS, vids=which(core > 2))
iDHHS4_6 <- induced_subgraph(iDHHS, vids=which(core > 3))
iDHHS5_6 <- induced_subgraph(iDHHS, vids=which(core > 4))
iDHHS6_6 <- induced_subgraph(iDHHS, vids=which(core > 5))

lay <- layout_with_fr(iDHHS)
op <- par(mfrow=c(3,2),mar = c(0,0,2,0))
plot(iDHHS1_6,layout=lay,main="All k-cores")
plot(iDHHS2_6,layout=lay[which(core > 1),], main="k-cores 2-6")
plot(iDHHS3_6,layout=lay[which(core > 2),], main="k-cores 3-6")
plot(iDHHS4_6,layout=lay[which(core > 3),], main="k-cores 4-6")
plot(iDHHS5_6,layout=lay[which(core > 4),], main="k-cores 5-6")
plot(iDHHS6_6,layout=lay[which(core > 5),], main="k-cores 6-6")

par(op)

Detección de comunidades: Un subgrupo en una red es un conjunto de nodos que tiene un número grande de vínculos internos, y también pocos vínculos del grupo con otras partes de la red.

• Modularidad: Es una medida de la estructura de la red, específicamente el grado en que los nodos exhiben agrupación donde hay mayor densidad dentro de los clusters y menor densidad entre ellos. Esta se define como la fracción de lazos que caen dentro de los grupos dados menos la fracción esperada si los lazos se distribuyeran al azar. Algoritmos tratan de maximizarla para devolver una clasificación de los nodos que mejor se desempeñe.

• Algoritmos para la detección de comunidades: se cuenta con 9 algoritmos para la detección de comunidades en igraph, cada una con el tipo de red al que van dirigidos.

Ejemplo de juguete:

g1 <- graph_from_literal(A-B-C-A,D-E-F-D,G-H-I-G,A-D-G-A)
V(g1)$grp_good <- c(1,1,1,2,2,2,3,3,3)
V(g1)$grp_bad <- c(1,2,3,2,3,1,3,1,2)

op <- par(mfrow=c(1,2))
plot(g1,vertex.color=(V(g1)$grp_good), vertex.size=20, main="Good Grouping")
plot(g1,vertex.color=(V(g1)$grp_bad), vertex.size=20, main="Bad Grouping")

par(op)

Este análisis se corrobora al calcular la modularidad:

modularity(g1,V(g1)$grp_good)

## [1] 0.4166667

modularity(g1,V(g1)$grp_bad)

## [1] -0.3333333

Ejemplo: Red DHHS

Aplicado el concepto de modularidad a la red DHHS se obtiene que la variable ‘agency’ explica un poco del agrupamiento presente en la red.

iDHHS <- asIgraph(DHHS)
modularity(iDHHS,(V(iDHHS)$agency+1))

## [1] 0.1402264

Ejemplo: Red Moreno

Se analizan las modularidades de los datos Moreno y Facebook con el fin de compararlas. En la primera base se usa la clasificación dada por el género, mientras que en la segunda base este agrupamiento lo da la variable grupo. Ambas redes tienen mayor modularidad que DHHS

# 1ros datos
data(Moreno)
iMoreno <- asIgraph(Moreno)
table(V(iMoreno)$gender)

## 
##  1  2 
## 16 17

# modularidad
modularity(iMoreno,V(iMoreno)$gender)

## [1] 0.4761342

# 2dos datos
data(Facebook)
levels(factor(V(Facebook)$group))

## This graph was created by an old(er) igraph version.
## ℹ Call `igraph::upgrade_graph()` on it to use with the current igraph version.
## For now we convert it on the fly...

## [1] "B" "C" "F" "G" "H" "M" "S" "W"

grp_num <- as.numeric(factor(V(Facebook)$group))

# modularidad 
modularity(Facebook,grp_num)

## [1] 0.6145798

Para los datos de Moreno la modularidad mediante el algoritmo es bastante alta

cw <- cluster_walktrap(iMoreno)
membership(cw)

##  [1] 1 1 1 1 1 1 1 1 3 3 3 5 5 5 5 1 3 2 2 2 4 4 4 2 2 2 2 2 2 2 2 6 6

modularity(cw)

## [1] 0.6181474

plot(cw, iMoreno)

Ejemplo: Red Bali

Una práctica común es usar más de un algoritmo de detección y comparar los resultados

data(Bali)
iBali <- asIgraph(Bali)

cw <- cluster_walktrap(iBali)
modularity(cw)

## [1] 0.2830688

membership(cw)

##  [1] 2 1 2 1 2 2 1 2 2 3 3 3 3 3 2 2 2

ceb <- cluster_edge_betweenness(iBali)
modularity(ceb)

## [1] 0.2387251

membership(ceb)

##  [1] 1 1 1 1 1 1 1 1 1 2 2 2 2 2 1 1 1

cs <- cluster_spinglass(iBali)
modularity(cs)

## [1] 0.2966742

membership(cs)

##  [1] 2 3 2 1 2 2 3 2 2 1 1 1 1 1 2 2 2

cfg <- cluster_fast_greedy(iBali)
modularity(cfg)

## [1] 0.2629126

membership(cfg)

##  [1] 2 2 1 2 1 2 2 1 1 3 3 3 3 3 1 1 1

clp <- cluster_label_prop(iBali)
modularity(clp)

## [1] 0.2387251

membership(clp)

##  [1] 1 1 1 1 1 1 1 1 1 2 2 2 2 2 1 1 1

cle <- cluster_leading_eigen(iBali)
modularity(cle)

## [1] 0.2750063

membership(cle)

##  [1] 1 1 1 2 1 1 2 1 1 2 2 2 2 2 1 1 1

cl <- cluster_louvain(iBali)
modularity(cl)

## [1] 0.2966742

membership(cl)

##  [1] 1 2 1 3 1 1 2 1 1 3 3 3 3 3 1 1 1

co <- cluster_optimal(iBali)
modularity(co)

## [1] 0.2966742

membership(co)

##  [1] 1 2 1 3 1 1 2 1 1 3 3 3 3 3 1 1 1

Los resultados también se pueden comparar utilizando alguna métrica como el Indice ajustado de Rand

compare(as.numeric(factor(V(iBali)$role)),cw,method="adjusted.rand")

## [1] 0.3504908

compare(cw,ceb,method="adjusted.rand")

## [1] 0.6155779

compare(cw,cs,method="adjusted.rand")

## [1] 0.8898148

compare(cw,cfg,method="adjusted.rand")

## [1] 0.6691802

op <- par(mfrow=c(3,2),mar=c(0,0,2,0))
plot(ceb, iBali,vertex.label=V(iBali)$role, main="Edge Betweenness")
plot(cfg, iBali,vertex.label=V(iBali)$role, main="Fastgreedy")
plot(clp, iBali,vertex.label=V(iBali)$role, main="Label Propagation")
plot(cle, iBali,vertex.label=V(iBali)$role, main="Leading Eigenvector")
plot(cs, iBali,vertex.label=V(iBali)$role, main="Spinglass")
plot(cw, iBali,vertex.label=V(iBali)$role, main="Walktrap")

par(op)

CAPÍTULO 9: REDES DE AFILIACIÓN

Una red de afiliación es aquella en la que los miembros están afiliados entre sí sobre la base de la coparticipación en algún tipo de evento. Por ejemplo, se puede pensar que los estudiantes que pertenecen a la misma clase están conectados entre sí, aunque tal vez no sepamos si tienen vínculos sociales directos.

• Afiliaciones como redes modo-2: La siguiente es una matriz de incidencia que nos dice como los actores (estudiantes) pertenecen a ciertos grupos (clases). Esta matriz tiene dimensión \(n\times g\)

C1 <- c(1,1,1,0,0,0)
C2 <- c(0,1,1,1,0,0)
C3 <- c(0,0,1,1,1,0)
C4 <- c(0,0,0,0,1,1)
aff.df <- data.frame(C1,C2,C3,C4)
row.names(aff.df) <- c("S1","S2","S3","S4","S5","S6")
aff.df

##    C1 C2 C3 C4
## S1  1  0  0  0
## S2  1  1  0  0
## S3  1  1  1  0
## S4  0  1  1  0
## S5  0  0  1  1
## S6  0  0  0  1

• Grafos bipartitos: En las redes de afiliación siempre hay dos tipos de nodos: los actores y los grupos a los que pertenecen. Entonces, las aristas conectan actores a sus grupos, razón por la cual no hay aristas directas entre actores o grupos

bn <- graph_from_biadjacency_matrix(as.matrix(aff.df))
plt.x <- c(rep(2,6),rep(4,4))
plt.y <- c(7:2,6:3)
lay <- as.matrix(cbind(plt.x,plt.y))
shapes <- c("circle","square")
colors <- c("blue","red")
plot(bn,vertex.color=colors[V(bn)$type+1], vertex.shape=shapes[V(bn)$type+1], vertex.size=10,vertex.label.degree=-pi/2,
     vertex.label.dist=1.2,vertex.label.cex=0.9, layout=lay)

Conceptos básicos de las redes de afiliación

• Red de afiliación a partir de matrices de incidencia: La función graph_from_biadjacency_matrix() toma una matriz y la transforma en una red de afiliación leyendo las filas como actores y las columnas como grupos.

• Red de afiliación a partir de la lista de aristas: Las listas de aristas se pueden transformar en redes de afiliación siempre bque los nodos de un tipo solo estén conectados a nodos del otro tipo.

• Graficar las redes de afiliación: Es usual designar diferentes colores y formas a los vértices de acuerdo a su afiliación, esto facilita la interpretación.

Ejemplos de juguete:

# Red desde la matriz de incidencia
(bn <- graph_from_biadjacency_matrix(as.matrix(aff.df)) )

## IGRAPH ef874df UN-B 10 11 -- 
## + attr: type (v/l), name (v/c)
## + edges from ef874df (vertex names):
##  [1] S1--C1 S2--C1 S2--C2 S3--C1 S3--C2 S3--C3 S4--C2 S4--C3 S5--C3 S5--C4
## [11] S6--C4

V(bn)$type

##  [1] FALSE FALSE FALSE FALSE FALSE FALSE  TRUE  TRUE  TRUE  TRUE

V(bn)$name

##  [1] "S1" "S2" "S3" "S4" "S5" "S6" "C1" "C2" "C3" "C4"

# Red desde la matriz de lista de aristas
el.df <- data.frame(rbind(c("S1","C1"), c("S2","C1"), c("S2","C2"), c("S3","C1"),
                          c("S3","C2"), c("S3","C3"), c("S4","C2"), c("S4","C3"),
                          c("S5","C3"), c("S5","C4"), c("S6","C4")))
el.df

##    X1 X2
## 1  S1 C1
## 2  S2 C1
## 3  S2 C2
## 4  S3 C1
## 5  S3 C2
## 6  S3 C3
## 7  S4 C2
## 8  S4 C3
## 9  S5 C3
## 10 S5 C4
## 11 S6 C4

bn2 <- graph_from_data_frame(el.df,directed=FALSE)
bn2

## IGRAPH ef91e7c UN-- 10 11 -- 
## + attr: name (v/c)
## + edges from ef91e7c (vertex names):
##  [1] S1--C1 S2--C1 S2--C2 S3--C1 S3--C2 S3--C3 S4--C2 S4--C3 S5--C3 S5--C4
## [11] S6--C4

# Establecemos esto para que R sepa que la red es bipartita
V(bn2)$type <- V(bn2)$name %in% el.df[,1]
bn2

## IGRAPH ef91e7c UN-B 10 11 -- 
## + attr: name (v/c), type (v/l)
## + edges from ef91e7c (vertex names):
##  [1] S1--C1 S2--C1 S2--C2 S3--C1 S3--C2 S3--C3 S4--C2 S4--C3 S5--C3 S5--C4
## [11] S6--C4

# Ejemplo de gráfico de una red de afiliación
shapes <- c("circle","square")
colors <- c("blue","red")
plot(bn,vertex.color=colors[V(bn)$type+1], vertex.shape=shapes[V(bn)$type+1], vertex.size=10,vertex.label.degree=-pi/2,
     vertex.label.dist=1.2,vertex.label.cex=0.9)

• Proyecciones: Examinar ambos tipos de nodos en una red de afiliación es a menudo el primer paso para estudiar este tipo de red. Sin embargo, también es útil examinar las conexiones directas entre los nodos de un tipo a la vez (clases y estudiantes, en este caso). Esto se puede hacer extrayendo y visualizando las proyecciones de modo único de la red de afiliación de dos modos.

Ejemplo de juguete:

# Proyecciónd e la red bipartita
bn.pr <- bipartite_projection(bn)
bn.pr

## $proj1
## IGRAPH efb619d UNW- 6 8 -- 
## + attr: name (v/c), weight (e/n)
## + edges from efb619d (vertex names):
## [1] S1--S2 S1--S3 S2--S3 S2--S4 S3--S4 S3--S5 S4--S5 S5--S6
## 
## $proj2
## IGRAPH efb61e9 UNW- 4 4 -- 
## + attr: name (v/c), weight (e/n)
## + edges from efb61e9 (vertex names):
## [1] C1--C2 C1--C3 C2--C3 C3--C4

graph.density(bn.pr$proj1)

## Warning: `graph.density()` was deprecated in igraph 2.0.0.
## ℹ Please use `edge_density()` instead.
## This warning is displayed once per session.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.

## [1] 0.5333333

bn.student <- bn.pr$proj1
bn.class <- bn.pr$proj2
graph.density(bn.student)

## [1] 0.5333333

shapes <- c("circle","square")
colors <- c("blue","red")
op <- par(mfrow=c(1,2))
plot(bn.student,vertex.color="blue", vertex.shape="circle",main="Students", edge.width=E(bn.student)$weight*2,
     vertex.size=15,vertex.label.degree=-pi/2, vertex.label.dist=1.2,vertex.label.cex=1)
plot(bn.class,vertex.color="red", vertex.shape="square",main="Classes", edge.width=E(bn.student)$weight*2,
     vertex.size=15,vertex.label.degree=-pi/2, vertex.label.dist=1.2,vertex.label.cex=1)

par(op)

Ejemplo: Actores de Hollywood como una red de afiliación

Los actores de Hollywood son un buen ejemplo de una red de afiliación, los actores están conectados entre sí a través de las películas en las que aparecen juntos. El conjunto de datos hwd es un objeto gráfico bipartito igraph. Los datos contienen las 10 películas más populares para cada año entre 199 hasta 2014, y los primeros 10 actores enlistados en cada película.

data(hwd)
h1 <- hwd
h1

## This graph was created by an old(er) igraph version.
## ℹ Call `igraph::upgrade_graph()` on it to use with the current igraph version.
## For now we convert it on the fly...

## IGRAPH 9cdab39 UN-B 1365 1600 -- 
## + attr: name (v/c), type (v/l), year (v/n), IMDBrating (v/n),
## | MPAArating (v/c)
## + edges from 9cdab39 (vertex names):
##  [1] Inception          --Leonardo DiCaprio   
##  [2] Inception          --Joseph Gordon-Levitt
##  [3] Inception          --Ellen Page          
##  [4] Inception          --Tom Hardy           
##  [5] Inception          --Ken Watanabe        
##  [6] Inception          --Dileep Rao          
##  [7] Inception          --Cillian Murphy      
## + ... omitted several edges

V(h1)$name[1:10]

##  [1] "Inception"                                   
##  [2] "Alice in Wonderland"                         
##  [3] "Kick-Ass"                                    
##  [4] "Toy Story 3"                                 
##  [5] "How to Train Your Dragon"                    
##  [6] "Despicable Me"                               
##  [7] "Scott Pilgrim vs. the World"                 
##  [8] "Hot Tub Time Machine"                        
##  [9] "Harry Potter and the Deathly Hallows: Part 1"
## [10] "Tangled"

V(h1)$type[1:10]

##  [1] TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE TRUE

Esta red cuenta con 1365 nodos (160 películas y 1205 actores) y 1600 aristas. Como la red es tan masiva, graficamos solo un pedazo de ella

par(mar=c(0, 0, 2, 0))
#set.seed(2453562)
V(h1)$shape <- ifelse(V(h1)$type==TRUE, "square","circle")
V(h1)$shape[1:10]

##  [1] "square" "square" "square" "square" "square" "square" "square" "square"
##  [9] "square" "square"

V(h1)$color <- ifelse(V(h1)$type==TRUE, "red","lightblue")

h2 <- subgraph.edges(h1, E(h1)[.inc(V(h1)[name %in% c("The Wolf of Wall Street", "Gangs of New York", "The Departed")])])

## Warning: `subgraph.edges()` was deprecated in igraph 2.1.0.
## ℹ Please use `subgraph_from_edges()` instead.
## This warning is displayed once per session.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.

plot(h2, layout = layout_with_kk, main = "Red de afiliacion de Scorsese")

Para esta red, el grado del nodo de los actores es informativo (todas las pelícualas tienen grado=10, así que este no es relevante)

table(igraph::degree(h1,v=V(h1)[type==FALSE]))

## 
##   1   2   3   4   5   6   7   8 
## 955 165  47  23  11   2   1   1

mean(igraph::degree(h1,v=V(h1)[type==FALSE]))

## [1] 1.327801

# Esto deja ver que la mayoria de actores solo aparecen en una película

Los actores más ocupados de la última década se identifican así:

V(h1)$deg <- igraph::degree(h1)
V(h1)[type==FALSE & deg > 4]$name

##  [1] "Leonardo DiCaprio" "Emma Watson"       "Richard Griffiths"
##  [4] "Harry Melling"     "Daniel Radcliffe"  "Rupert Grint"     
##  [7] "James Franco"      "Ian McKellen"      "Martin Freeman"   
## [10] "Bradley Cooper"    "Christian Bale"    "Samuel L. Jackson"
## [13] "Natalie Portman"   "Brad Pitt"         "Liam Neeson"

busy_actor <- data.frame(cbind(Actor = V(h1)[type==FALSE & deg > 4]$name,
                               Movies = V(h1)[type==FALSE & deg > 4]$deg))
busy_actor[order(busy_actor$Movies,decreasing=TRUE),]

##                Actor Movies
## 5   Daniel Radcliffe      8
## 11    Christian Bale      7
## 1  Leonardo DiCaprio      6
## 2        Emma Watson      6
## 3  Richard Griffiths      5
## 4      Harry Melling      5
## 6       Rupert Grint      5
## 7       James Franco      5
## 8       Ian McKellen      5
## 9     Martin Freeman      5
## 10    Bradley Cooper      5
## 12 Samuel L. Jackson      5
## 13   Natalie Portman      5
## 14         Brad Pitt      5
## 15       Liam Neeson      5

Se puede acceder a la popularidad de los actores basada en la popularidad de las películas en las que aparecen, esta última medida con la variable IMDBrating.

for (i in 161:1365) {
V(h1)[i]$totrating <- sum(V(h1)[.nei(i)]$IMDBrating)
}
max(V(h1)$totrating,na.rm=TRUE)

## [1] 60.9

pop_actor <- data.frame(cbind(Actor = V(h1)[type==FALSE & totrating > 40]$name,
                              Popularity = V(h1)[type==FALSE &
                                                   totrating > 40]$totrating))
pop_actor[order(pop_actor$Popularity,decreasing=TRUE),]

##               Actor Popularity
## 3  Daniel Radcliffe       60.9
## 4    Christian Bale       55.5
## 1 Leonardo DiCaprio       49.6
## 2       Emma Watson         45
## 5         Brad Pitt       40.5

Se quiere analizar ahora si los actores más ocupados son protagonistas en más películas populares en promedio.

for (i in 161:1365) {
  # claular la media de la puntuación
V(h1)[i]$avgrating <- mean(V(h1)[.nei(i)]$IMDBrating)}
num <- V(h1)[type==FALSE]$deg
avgpop <- V(h1)[type==FALSE]$avgrating
summary(lm(avgpop ~ num))

## 
## Call:
## lm(formula = avgpop ~ num)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -3.9858 -0.4330  0.1977  0.6170  1.6142 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  7.33868    0.05440 134.911   <2e-16 ***
## num          0.04714    0.03527   1.337    0.182    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.9605 on 1203 degrees of freedom
## Multiple R-squared:  0.001483,   Adjusted R-squared:  0.0006528 
## F-statistic: 1.786 on 1 and 1203 DF,  p-value: 0.1816

# Relación entre número de películas y la calificación promedio de ellas
scatter.smooth(num,avgpop,col="lightblue", ylim=c(2,10),span=.8, xlab="Number of Movies", ylab="Avg. Popularity")

Proyecciones de actores y películas

h1.pr <- bipartite_projection(h1)
h1.act <- h1.pr$proj1
h1.mov <- h1.pr$proj2
h1.act

## IGRAPH f415f30 UNW- 1205 6903 -- 
## + attr: name (v/c), year (v/n), IMDBrating (v/n), MPAArating (v/c),
## | shape (v/c), color (v/c), deg (v/n), totrating (v/n), avgrating
## | (v/n), weight (e/n)
## + edges from f415f30 (vertex names):
##  [1] Leonardo DiCaprio--Joseph Gordon-Levitt
##  [2] Leonardo DiCaprio--Ellen Page          
##  [3] Leonardo DiCaprio--Tom Hardy           
##  [4] Leonardo DiCaprio--Ken Watanabe        
##  [5] Leonardo DiCaprio--Dileep Rao          
##  [6] Leonardo DiCaprio--Cillian Murphy      
## + ... omitted several edges

h1.mov

## IGRAPH f415f82 UNW- 160 472 -- 
## + attr: name (v/c), year (v/n), IMDBrating (v/n), MPAArating (v/c),
## | shape (v/c), color (v/c), deg (v/n), totrating (v/n), avgrating
## | (v/n), weight (e/n)
## + edges from f415f82 (vertex names):
##  [1] Inception--The Wolf of Wall Street    Inception--Django Unchained          
##  [3] Inception--The Departed               Inception--Gangs of New York         
##  [5] Inception--Catch Me If You Can        Inception--The Dark Knight Rises     
##  [7] Inception--10 Things I Hate About You Inception--Batman Begins             
##  [9] Inception--The Dark Knight            Inception--Training Day              
## [11] Inception--Big Fish                  
## + ... omitted several edges

par(mar = c(0, 0, 2, 0))
set.seed(1)
plot(h1.mov, vertex.color="red", vertex.shape="circle", vertex.size=(V(h1.mov)$IMDBrating)-3, vertex.label=NA,
     main = "Red de afiliación de las películas")

graph.density(h1.mov)

## [1] 0.03710692

count_components(h1.mov)

## [1] 12

components(h1.mov)$csize

##  [1] 148   1   1   1   1   1   1   2   1   1   1   1

table(E(h1.mov)$weight)

## 
##   1   2   3   4   5   6   7  10 
## 411  21  12  16   6   1   2   3

par(mar = c(0, 0, 2, 0))
set.seed(1)
h2.mov <- induced_subgraph(h1.mov, vids=components(h1.mov)$membership==1)
plot(h2.mov,vertex.color="red", edge.width=sqrt(E(h1.mov)$weight), vertex.shape="circle",
     vertex.size=(V(h2.mov)$IMDBrating)-3, vertex.label=NA, main = "Componente mas grande de la red de peliculas")

table(coreness(h2.mov))

## 
##  1  2  3  4  5  6  7 
## 11  5 23 65 29  7  8

par(mar = c(0, 0, 2, 0))
set.seed(1)

# Grafo inducido por los nodos que están más conectados en la red
  # coreness más alto
h3.mov <- induced_subgraph(h2.mov, vids = coreness(h2.mov)>4)
plot(h3.mov,vertex.color="red", vertex.shape="circle", edge.width=sqrt(E(h1.mov)$weight),
     vertex.label.cex=0.7,vertex.label.color="darkgreen", vertex.label.dist=0.3, 
     vertex.size=(V(h3.mov)$IMDBrating)-3, 
     main = "Grafo inducido por los nodos que están más conectados \n en la red de películas")

10. Para todas las visualizaciones decoradas de este taller, realizar visualizaciones alternativas usando heat maps y hive plots.

Solución Punto 10

# funciones
get_adjacency_ordered <- function(xi, A) {
  xi2 <- xi[order(xi)]
  indices <- order(xi)
  d <- NULL
  for (i in 1:(length(xi)-1)) if (xi2[i] != xi2[i+1]) d <- c(d, i)
  list(A = A[indices,indices], d = d)
}

heat.plot0 <- function (mat, show.grid = FALSE, cex.axis, tick, labs, col.axis, ...) { 
  JJ <- dim(mat)[1]
  colorscale <- colorRampPalette(c("white", "lightblue", "blue", "darkblue"))(50)
  #colorscale <- colorRampPalette(c("white", "orange", "red", "darkred"))(200)
  if(missing(labs))     labs <- 1:JJ
  if(missing(col.axis)) col.axis <- rep("black", JJ)
  if(missing(cex.axis)) cex.axis <- 0.5
  if(missing(tick))     tick <- TRUE
  # adjacency matrix
  image(seq(1, JJ), seq(1, JJ), mat, axes = FALSE, xlab = "", ylab = "", col = colorscale, ...)
  for (j in 1:JJ) {
    axis(1, at = j, labels = labs[j], las = 2, cex.axis = cex.axis, tick, col.axis = col.axis[j], col.ticks = col.axis[j])
    axis(2, at = j, labels = labs[j], las = 2, cex.axis = cex.axis, tick, col.axis = col.axis[j], col.ticks = col.axis[j])
  }
  box()
  if (show.grid) grid(nx = JJ, ny = JJ)
}

Las visualizaciones decoradas presentes en este taller corresponden a la base de datos conflict.RData, tratada en los puntos (6) y (8).

Red de conflictos dirigida y ponderada: La primera matriz de adyacencia a considerar dat$Y es la resultante de la red conflict.RData. Esta es dirigida, ponderada y cuenta con 130 actores.

Veamos su heatmap.

Y_6 <- dat$Y
heatmap(Y_6,Rowv = NA, Colv = NA, scale = "none",
  #col = rev(heat.colors(100)),
  margins = c(5,5), main = "Heatmap de la Red de conflictos dirigida")

heat.plot0(mat = Y_6, labs = rownames(Y_6), main = "Heatmap para la red de conflictos dirigida")

El heatmap, que en este caso es bastante ‘pobre’ en temas de dinámicas de color va en concordancia con una red bastante dispersa en la que sobresalen países notorios como Iraq y Jordania.

Ahora, para graficar el hiveplot es importante establecer primero que van a representar los ejes de la colmena. Vamos a considerar dos escenarios:

1. Los ejes están dados por la clasificación dada por coreness
2. Los ejes están dados por una categorización del PIB

En ambos casos, el color de las aristas estará dado por el tipo de relación entre dos países y su grosor por el peso de la red (número de conflictos).

library(HiveR)
V(g6_conflict)$core <- coreness(g6_conflict)
V(g6_conflict)$core_cat <- cut(V(g6_conflict)$core, breaks = 3,
                               labels = c("Periferia", "Intermedio", "Núcleo"))
V(g6_conflict)$gdp_cat <- cut(log(V(g6_conflict)$gdp), breaks = 3,
                              labels = c("Bajo", "Medio", "Alto"))

ggraph(g6_conflict, 'hive', axis = core_cat, sort.by = 'degree') + 
  geom_edge_hive(aes(colour = tipo_regimen,
                     width = weight), 
                 alpha = 0.7) +
  geom_axis_hive(aes(colour = core_cat), 
                 size = 3, 
                 label = FALSE) + 
  scale_edge_width(range = c(0.5, 4), 
                   guide = "none") +
  scale_edge_color_manual(name = "Tipo de relación",
                          values = c("Autocracia-Autocracia" = "dodgerblue4", 
                                     "Democracia-Democracia" = "darkgreen",    
                                     "Mixta" = "mediumpurple1") 
                          ) +
  scale_color_discrete(name = "Coreness")+
  coord_fixed() + 
  theme_gray() +
  theme(axis.text = element_blank(), axis.title = element_blank()
) +
  labs(
    title = "Hiveplot para la red dirigida ponderada de conflictos",
    subtitle = "Los ejes están dados por el coreness de la red"
  )

ggraph(g6_conflict, 'hive', axis = gdp_cat, sort.by = 'degree') + 
  geom_edge_hive(aes(colour = tipo_regimen,
                     width = weight), 
                 alpha = 0.7) +
  geom_axis_hive(aes(colour = gdp_cat), 
                 size = 3, 
                 label = FALSE) + 
  scale_edge_width(range = c(0.5, 4), 
                   guide = "none") +
  scale_edge_color_manual(name = "Tipo de relación",
                          values = c("Autocracia-Autocracia" = "dodgerblue4", 
                                     "Democracia-Democracia" = "darkgreen",    
                                     "Mixta" = "mediumpurple1") 
                          ) +
  scale_color_discrete(name = "PIB")+
  coord_fixed() + 
  theme_gray() +
  theme(axis.text = element_blank(), axis.title = element_blank()
) +
  labs(
    title = "Hiveplot para la red dirigida ponderada de conflictos",
    subtitle = "Los ejes están dados por el PIB de los países"
  )

Es interesante ver como nodos que están en la periferia están conectados mayormente con nodos del núcleo de la red a través de relaciones Autocracia-Autocracia, así como ver la alta interacción que se presenta en la parte central e intermedia de la red.

En segundo lugar, el Hiveplot cuyos ejes son escalas del PIB de los países revela una tendencia. Entre países con alto y medio nivel de PIB son más frecuentes las relaciones de tipo mixto: Democracia-Autocracia o viceversa, mientras que los países con bajo producto interno bruto parecen tener má conflictos de tipo Autocracia-Autocracia con países de nivel medio de PIB. De igual manera, los conflictos entre extremos (países con PIB bajo y alto) se presentan menos.

Red de conflictos no dirigida y sin datos aislados: En este caso, la matriz de adyacencia de la red conflict.RData se transformó de forma que las relaciones ahora son no dirigidas y no se consideran los puntos aislados, resultando en una red de 91 actores.

Y_8 <- Y8
heatmap(Y_8,Rowv = NA, Colv = NA, scale = "none",
  #col = rev(heat.colors(100)),
  margins = c(5,5), main = "Heatmap de la Red de conflictos dirigida")

heat.plot0(mat = Y_8, labs = rownames(Y_8))

El hecho de que en este segundo ejercicio no se tengan puntos aislados se nota en el heatmap, se ve mayor presencia de conflictos entre países. Esta representación hace aún más evidente el papel importante que juegan países como Iraq, Jordania y Estados Unidos.

Adicionalmente, esta red fue caracterizada en términos de un agrupamiento que maximizó la modularidad, véase así mismo el Diagrama heatmap resultante al agrupar a los nodos por su respectivo cluster.

# asignaciones 
xi <- clust8_optimal$membership
# asignaciones ordenadas 
xi2 <- xi[order(xi)]
# matriz de adyacencia original
Y <- as_adjacency_matrix(graph = g8_conflict, sparse = F)
# matriz de adyacencia ordenada y lineas divisorias de acuerdo con las comunidades
tmp <- get_adjacency_ordered(xi = xi, A = Y)
A <- tmp$A
d <- tmp$d


# visualización
#par(mfrow = c(1,2))
n_comm <- length(unique(clust8_optimal$membership))
colores_ggraph <- scales::hue_pal()(n_comm)
colores_nodos <- colores_ggraph[clust8_optimal$membership]

# mismo orden que la matriz reordenada
indices <- order(xi)

labs_ord <- rownames(Y)[indices]

heat.plot0(mat = A, col.axis = colores_nodos[indices],labs = labs_ord)
#heat.plot0(mat = A, col.axis = colores_nodos[indices], labs = rownames(Y))
abline(v = d + .5, h = d + .5)

Este mapa de color deja ver que efectivamente la clasificación elegida maximiza los conflictos entre países del mismo cluster y minimiza conflictos entre actores de diferentes grupos. De forma descriptiva note que, Colombia está en el mismo grupo que Venezuela, Costa Rica y Honduras :b.

Ahora, veamos el hiveplot para esta red, usando como ejes las dos siguientes opciones:

1. La categorización del PIB de los países
2. Los clusters obtenidos por el algoritmo de agrupamiento

V(g8_conflict)$gdp_cat <- cut(log(V(g8_conflict)$gdp), breaks = 3,
                              labels = c("Bajo", "Medio", "Alto"))
V(g8_conflict)$cluster <- as.factor(clust8_optimal$membership)

ggraph(g8_conflict, 'hive', axis = gdp_cat, sort.by = 'degree') + 
  geom_edge_hive(aes(colour = tipo_regimen,
                     width = weight), 
                 alpha = 0.7) +
  geom_axis_hive(aes(colour = gdp_cat), 
                 size = 3, 
                 label = FALSE) + 
  scale_edge_width(range = c(1, 1), 
                   guide = "none") +
  scale_edge_color_manual(name = "Tipo de relación",
                          values = c("Autocracia-Autocracia" = "dodgerblue4", 
                                     "Democracia-Democracia" = "darkgreen",    
                                     "Mixta" = "mediumpurple1") 
                          ) +
  scale_color_discrete(name = "Categorización del PIB")+
  coord_fixed() + 
  theme_gray() +
  theme(axis.text = element_blank(), axis.title = element_blank()
) +
  labs(
    title = "Hiveplot para la red dirigida ponderada de conflictos",
    subtitle = "Los ejes están dados por el PIB de los países"
  )

ggraph(g8_conflict, 'hive', axis = cluster, sort.by = 'degree') + 
  geom_edge_hive(aes(colour = tipo_regimen,
                     width = weight), 
                 alpha = 0.7) +
  geom_axis_hive(aes(colour = cluster), 
                 size = 3, 
                 label = FALSE) + 
  scale_edge_width(range = c(1, 1), 
                   guide = "none") +
  scale_edge_color_manual(name = "Tipo de relación",
                          values = c("Autocracia-Autocracia" = "dodgerblue4", 
                                     "Democracia-Democracia" = "darkgreen",    
                                     "Mixta" = "mediumpurple1") 
                          ) +
  scale_color_discrete(name = "Cluster")+
  coord_fixed() + 
  theme_gray() +
  theme(axis.text = element_blank(), axis.title = element_blank()
) +
  labs(
    title = "Hiveplot para la red dirigida ponderada de conflictos",
    subtitle = "Los ejes están dados por el agrupamiento hecho con cluster_optimal"
  ) + guides(
  colour = guide_legend(ncol = 2),
  edge_colour = guide_legend(ncol = 1)
)

El hiveplot en el que los ejes son los agrupamientos creados deja ver que las pocas relaciones que se dan entre actores de distinto cluster, son en su mayoría del tipo Mixto.

11. Además de la ley de potencias para la distribución del grado, existen alternativas comunes como la distribución log normal, la exponencial, la Poisson y la ley de potencias con corte exponencial. Caracterizar cada una de estas opci ones de manera detallada y concisa, describiendo su forma y la intuición que representa en términos de heterogeneidad y presencia de nodos altamente conectados. Explicar para qué sirve modelar la distribución del grado, en qué contextos resulta apropiado usar cada familia, y cuáles son sus principales ventajas y desventajas desde el punto de vista interpretativo y empírico. Seleccionar una de estas distribuciones y ejemplificar su uso con un conjunto de datos real de su elección, incluyendo una evaluación del ajuste y una interpretación de los resultados.

Solución Punto 11

¿Por qué modelar la distribución del grado?

Caracterizar la distribución del grado con una una familia de distribuciones sirve para realizar un resumen interpretable que captura cuánta heterogeneidad hay en las conexiones. En términos prácticos permite:

• Identificar si existen hubs (nodos muy conectados)

• Evaluzar la heterogeneidad estructural

• Comparar redes

• Inferir mecanismos de formación de enlaces

Por ejemplo, en redes sociales no es común que hayan nodos muy influyentes , y en redes de conflicto hay actores dominantes.

Con esto presente, surge una preguna clave_ ¿La conectividad está distribuida uniformemente o concentrada en pocos nodos?

MODELOS PARA LA DISTRIBUCIÓN DEL GRADO

Recordemos que la distribución del grado de G es la colección de frecuencias relativas \(f_0, f_1, ...\), donde \[f_d = \frac{|\{v\in V:d_v=d\}|}{|V|},\] es decir, la fracción de vertices que tienen grado \(d_v=d\)

Distribución log-normal \[f_d = \propto\frac{1}{d}\exp\left(-\frac{(\ln d-\mu)^2}{2\sigma^2}\right)\] lo que en escala \(\log\) corresponde a \[\log(f_d)\sim N(\mu, \sigma^2).\] Esta produce una distribución asimétria en la que muchos nodos tienen bajo grado y hay una cola pesada a derecha (no ta pesada como la ley de potencias) con nodos que tienen grados altos. Es común en sistemas complejos no lineales en las que el crecimiento del grado es acumulativo multiplicativo, como en redes biológicas.

Veamos que por su forma representa una alta heterogeneidad en la conectividad y presencia de hubs; hasta el punto en el que muchos casos que parecían ser power-law en realidad se ajustaban mejor con una log-normal. Surge entonces la pregunta: ¿Como saber cual de estos es el más apropiado para mis datos? Graficando la distribución empírica del grado en escala \(\log\).

Distribución exponencial

\[f_d \propto e^{-\lambda d}\]

La forma de esta distribución es de decaimiento corto y una cola corta pero asimétrica.

El modelo ayuda a representar una heterogeneidad moderada gracias a su cola corta. En esta si hay algunos nodos más conectados que otros, sin tener hubs grandes. En comparación, este modelo es más homogeneo que una ley de potencias pero más heterogeneo que un Poisson.

Se usa cuando hay indicios de desigualdad en la conectividad pero restricciones fuertes en los actores con los que los nodos pueden estar conectados.

Distribución Poisson \[f_d = \frac{e^{-\lambda}\lambda^d}{d!}\] donde \(\alpha\) es un valor real positivo correspondiente al grado medio y la varianza de los grados.

Esta es una distribución que alcanza su punto más alto alrededor del promedio \(\lambda\). Debido a esto, todos los vértices en el grafo tienen valores de grado muy similares, haciéndolo el modelo clásico usado para generar redes aleatorias, pero al mismo tiempo, es un modelo insatisfactorio en redes reales donde la conectividad no se comporta de forma tan homogénea.

Como hay homogeneidad en los grados, los nodos muy conectados (hubs) son extremadamente raros.

Luego, como se mencionó anteriormente, este modelo es ideal para formar enlaces al azar en redes donde no hay preferencias estructurales.

Ley de potencias con corte exponencial \(\beta\)

\[f_d = c d^{-\alpha}\exp\left(-\frac{d}{\beta}\right)\] donde \(c\) es la constante de normalziación, \(\alpha\) es el exponente de la ley de potencias y \(\beta\) es el corte.

Conocida también como ley de potencias truncada, esta se comporta como una ley de potencias pero con un límite en el grado “probablemente” más alto. Por ejemplo, las probabilidad de encontrar un nodo con grado \(d=100\) puede ser igual a \(0.0001\) bajo el modelo power-law; mientras que con unn corte en \(\beta=10\), esta probabilidad deciende a \(0.00000001\).

Es correcto usarla cuando hay presencia de hubs con restricciones reales a la cantidad de conexciones que pueden tener los nodos. Por ejemplo, redes tecnológicas o redes sociales con límites cognitivos.

library(knitr)

## Warning: package 'knitr' was built under R version 4.5.3

library(kableExtra)

tabla <- data.frame(
  Criterio = c("Ventajas","Desventajas"),
  Poisson = c(
    "Fácil interpretación Simplicidad matemática Homogeneidad",
    "Poco realista Subestima la presencia de nodos muy conectados Decaimiento rápido"
  ),
  Exponencial = c(
    "Buen modelo intermedio Modelado de redes finitas",
    "No reproduce colas pesadas reales Subestima la presencia de nodos muy conectados"
  ),
  Log_normal = c(
    "Flexibilidad Ajuste empírico frecuente",
    "Complejidad de interpretación.      No es libre de escala"
  ),
  Powerlaw_truncada = c(
    "Probabilidades más realistas   Modela limitaciones",
    "Dificultad de estimación  Sensibilidad a cambios en los datos"
  )
)

kable(tabla, 
      escape = FALSE,
      booktabs = FALSE,   # ← permite líneas completas
      align = "c") %>%    # ← centra contenido
  kable_styling(full_width = TRUE) %>%
  column_spec(1, width = "2.5cm", border_right = TRUE) %>%
  column_spec(2:5, width = "3.5cm", border_right = TRUE) %>%
  row_spec(nrow(tabla), extra_css = "border-bottom: 1px solid;")

Criterio	Poisson	Exponencial	Log_normal	Powerlaw_truncada
Ventajas	Fácil interpretación Simplicidad matemática Homogeneidad	Buen modelo intermedio Modelado de redes finitas	Flexibilidad Ajuste empírico frecuente	Probabilidades más realistas Modela limitaciones
Desventajas	Poco realista Subestima la presencia de nodos muy conectados Decaimiento rápido	No reproduce colas pesadas reales Subestima la presencia de nodos muy conectados	Complejidad de interpretación. No es libre de escala	Dificultad de estimación Sensibilidad a cambios en los datos

APLICACIÓN:

La base de ICEWS es una de datos que se comprende de interacciones entre actores socio políticos de países alrededor del mundo. Los actores, pueden interactuar de 4 formas: Cometiendo un acto material negativo contra otro, cometiendo un acto material positivo hacia otro, dando un mensaje verbal negativo a otro y dando un mensaje verbal positivo a otro. Para el presente ejercicio, consideramos la relación de cometer un acto material positivo.

Así mismo, los datos se tienen para cada día del año 2014. Acá consideraremos la suma de las interacciones entre países en este periodo, resultando en una red podnerada y dirigida.

load("icews_2014.RData")
g11 <- graph_from_adjacency_matrix(Y_material_bueno, mode = "directed", weighted = TRUE, diag = FALSE)
vcount(g11)

## [1] 233

ecount(g11)

## [1] 1612

#g11_sub <- induced_subgraph(g11, vids = V(g11)[igraph::degree(g11, mode = "out") > 2])
g11_sub <- delete_vertices(g11, which(igraph::degree(g11, mode = "total") == 0) )

comp <- igraph::components(g11_sub)
giant_comp <- which.max(comp$csize)
g11_sub <- induced_subgraph(g11_sub,vids = V(g11_sub)[comp$membership == giant_comp])
vcount(g11_sub)

## [1] 190

ecount(g11_sub)

## [1] 1609

Para facilidad de análisis se va a trabajar con la componente convexa más grande de la red, la cual cuenta con 190 nodos y 1609 aristas.

library(networkD3)
library(htmltools)

V(g11_sub)$out_degree <- igraph::degree(g11_sub, mode = "out")

# Extraer aristas
edges <- igraph::as_data_frame(g11_sub, what = "edges")

# Extraer nodos
nodes <- data.frame(
  name = V(g11_sub)$name,
  out_degree = V(g11_sub)$out_degree
)
nodes$group <- 1

edges$source <- match(edges$from, nodes$name) - 1
edges$target <- match(edges$to, nodes$name) - 1

p <- forceNetwork(
  Links = edges,
  Nodes = nodes,
  Source = "source",
  Target = "target",
  NodeID = "name",
  Group = "group",
  # Tamaño proporcional al grado de salida
  Nodesize = "out_degree",
  Value = "weight",

  opacity = 0.9,
  zoom = TRUE,
  opacityNoHover = 0.2,
  linkColour = "rgba(150,150,150,0.3)",
  linkDistance = 50,
  charge = -120,
  fontSize = 10 
)
tagList(
  tags$h2("Red ICEWS de acciones materiales positivas en 2014"),
  p
)

Red ICEWS de acciones materiales positivas en 2014

El gráfico de la red ICEWS nos muestra que esta es una red a prori densa y con una conectividad variada. Para analizar esto a fondo, veamos la distribución de los grados de salida, de entrada y el total.

d11_out <- igraph::degree(graph = g11_sub, mode = "out")
d11_in <- igraph::degree(graph = g11_sub, mode = "in")
d11_total <- igraph::degree(graph = g11_sub, mode = "total")

ord_grados_in11 <- sort(d11_in, decreasing = TRUE)
ord_grados_out11 <- sort(d11_out, decreasing = TRUE)

head(ord_grados_in11, 5)

##      United States Russian Federation              China               Iraq 
##                 84                 48                 42                 35 
##             France 
##                 34

head(ord_grados_out11, 5)

##      United States              China     United Kingdom              Japan 
##                101                 60                 52                 46 
## Russian Federation 
##                 45

Los países que más ayudas materiales han brindado son Estados Unidos, Rusia y China. Del otro lado, los países que más ayuda material han recibido son Estados Unidos, China y el Reino Unido.

par(mfrow=c(1, 3))

plot(x = NA, y = NA, type = "n", xlim = c(0,110), ylim = c(0,0.12), 
     xlab = "Grado", ylab = "Densidad", main = "Distribucion del grado de salida")
hist(d11_out, breaks = seq(0, 200, by = 5), freq = FALSE, col = "lightskyblue", border = "royalblue",add = TRUE)

plot(x = NA, y = NA, type = "n", xlim = c(0,100), ylim = c(0,0.12), 
     xlab = "Grado", ylab = "Densidad", main = "Distribucion del grado de entrada")
hist(d11_in, breaks = seq(0, 200, by = 5),   freq = FALSE, col = "lightskyblue", border = "royalblue", add = TRUE)

plot(x = NA, y = NA,  type = "n", xlim = c(0,100), ylim = c(0,0.12), 
     xlab = "Grado", ylab = "Densidad", main = "Distribucion del grado total")

hist(d11_total,breaks = seq(0, 200, by = 5), freq = FALSE,col = "lightskyblue",border = "royalblue",add = TRUE)

El diagrama de dispersión soporta la hipótesis de que estamos frente a una red en la que priman los actores de bajo grado. Adicionalmente hay una clara asimetría a derecha para todos los tipos que se extiende hasta países con grados superiores a 100. Podemos afirmar que este es un caso de alta heterogeneidad en la conectividad de la red.

Dados los hallazgos ya mencionados, se elige una log-normal como modelo a ajustar en la distribución de los grados, la cual se va a comparar con la ya conocida ley de potencias. El ajuste se hará para cada tipo de grado por aparte, ya que se ven ligeros comportamientos distintos entre ellos y se quieren cubrir todos los escenarios.

library(poweRlaw)
ajustar <- function(d, nombre) {
  d <- d[d > 0]

  # --- Power-law ---
  m_pl <- displ$new(d)
  est_xmin <- estimate_xmin(m_pl)
  m_pl$setXmin(est_xmin)
  m_pl$setPars(estimate_pars(m_pl))

  # --- Log-normal ---
  m_ln <- dislnorm$new(d)
  m_ln$setXmin(m_pl$getXmin())
  m_ln$setPars(estimate_pars(m_ln))

  # --- Comparación Vuong ---
  comp <- compare_distributions(m_pl, m_ln)
  mejor <- dplyr::case_when(comp$p_two_sided >= 0.05 ~ "Indistinguibles",
                            comp$test_statistic > 0  ~ "PL",
                            TRUE                     ~ "Log-normal")
  params <- data.frame(
    Grado      = nombre,
    PL_alpha   = round(m_pl$pars, 3),
    LN_meanlog = round(m_ln$pars[1], 3),
    LN_sdlog   = round(m_ln$pars[2], 3),
    xmin       = m_pl$getXmin(),
    PL_loglik  = round(dist_ll(m_pl), 2),
    LN_loglik  = round(dist_ll(m_ln), 2),
    Vuong_Z    = round(comp$test_statistic, 3),
    Vuong_pval = round(comp$p_two_sided, 4),
    Conclusion = mejor
  )

  list( params = params, m_pl = m_pl, m_ln  = m_ln
  )
}
res_out   <- ajustar(d11_out,   "OUT")
res_in    <- ajustar(d11_in,    "IN")
res_total <- ajustar(d11_total, "TOTAL")

tabla <- bind_rows(res_out$params, res_in$params, res_total$params)
kable(tabla,
      align = "c",
      row.names = FALSE,
      booktabs = TRUE) %>%
  kable_styling(full_width = TRUE)

Grado	PL_alpha	LN_meanlog	LN_sdlog	xmin	PL_loglik	LN_loglik	Vuong_Z	Vuong_pval	Conclusion
OUT	2.606	1.083	1.182	10	-187.18	-186.33	-0.830	0.4064	Indistinguibles
IN	3.011	1.869	0.864	13	-139.74	-138.90	-0.739	0.4596	Indistinguibles
TOTAL	3.009	2.034	1.026	27	-150.88	-150.42	-0.570	0.5687	Indistinguibles

graficar_modelos <- function(res, titulo){
  m_pl <- res$m_pl
  m_ln <- res$m_ln

  plot(m_pl, main = titulo, xlab = "Grado d", ylab = "P(F ≥ d)")
  lines(m_pl, col = "red", lwd = 2)
  lines(m_ln, col = "darkgreen", lwd = 2)
  legend("bottomleft", legend = c("Datos", "Power-law", "Log-normal"), col = c("black", "red", "darkgreen"),
         lwd = c(NA,2,2), pch = c(1,NA,NA), bty = "n")
}

par(mfrow=c(1,3))
graficar_modelos(res_out, "CCDF en escala log-log \n grado de salida")
graficar_modelos(res_in, "CCDF en escala log-log \n grado de entrada")
graficar_modelos(res_total, "CCDF en escala log-log \n grado total")

Note que tanto la comparación de ajuste de los modelos como los gráficos de las CCDF en escala log-log muestran que ambos modelos, el de ley de potencias y el log-normal, no presentan mayores diferencias entre sí.

Ahora, desde el grado mínimo en el que se ajustan las modelos, ambos se ajustan correctamente a los datos puesto que se observa que los datos siguen una línea recta en la escala log-log a partir del \(d\) mínimo estimado. ¿Por qué se estima esto? Antes de ese grado mínimo la distribución del grado no suele seguir estos modelos, porque tiene un comportamiento más irregular.

A pesar de que el ajuste de ambos modelos es satisfactorio, note que hay un actor ubicado en la parte inferior derecha del gráfico que puede estar ‘jalando’ el ajuste de la ley de potencias hacia una pendiente menos pronunciada.

Al final, los grados de salida, de entrada y totales no presentaron diferencias notorias en sus modelos ajustados, los \(\alpha\) de power-law y los parámetros de la log-normal resultaron bastante parecidos. Esto justifica que en posibles análisis de la distribución del grado, pueda analizarse solamente uno de estos.

12. Además del grado, la centralidad por cercanía, la centralidad por intermediación y la centralidad propia, existen otras medidas relevantes para cuantificar la importancia de los nodos en una red. En el caso de redes no dirigidas, caracterizar de manera detallada y concisa las medidas PageRank, Katz centrality, Harmonic centrality y VoteRank, incluyendo la intuición estadística que las motiva y el tipo de estructura que cada una prioriza. Luego, con base en el conjunto de datos real utilizado en el numeral anterior, calcular estas medidas para todos los nodos, presentar los principales resultados, comparar los rankings obtenidos e interpretar las diferencias observadas en términos de la estructura de la red y del rol de los nodos más destacados.

Solución Punto 12

PageRank

El método PageRank se crea para medir la importancia relativa de las páginas web. Este calcula un ranking para cada página basado en el grafo de la web. El algoritmo mide la importancia de cada nodo dentro del grafo, basado en el número de relaciones entrantes que este tiene y la importancia de los nodos fuente correspondientes. La suposición subyacente es que una página es tan importante como las páginas que enlazan a ella.

En su forma más simple, el rango de una página se define como la suma de los rangos de sus páginas de origen divididos por su número de enlaces salientes. Sea \(PR(i)\) el PageRank del nodo \(i\) y \(d_i^{out}\) su grado de salida, el PageRank del nodo \(j\) es: \[PR(j)=\sum_{i\rightarrow j}\frac{PR(i)}{d_i^{out}}\] donde \(i\) hace parte de los nodos de entrada de \(j\). Esta versión asume que cada enlace transfiere una parte equitativa de la autoridad de la página de origen.

Para solucionar problemas de páginas sin enlaces de salida, se introduce un factor \(c\) de amortiguamiento (usualmente establecido en 0.85), el cual representa la probabilidad de que un usuario siga haciendo clic en los enlaces citados en lugar de saltar a una página aleatoria.

Así, la formula completa resulta: \[PR(j)=\frac{1-c}{N}+c\sum_{i\in B_j}\frac{PR(i)}{d_i^{out}}\] donde \(N\) es el número total de páginas (nodos) en la red, \(B_j\) es el conjunto de páginas que enlazan a \(j\) (sus nodos de entrada).

Katz centrality

La centralidad de Katz calcula la centralidad de un nodo basándose en la centralidad de sus vecinos. Es una generalización de la centralidad propia para grafos dirigidos.

Para un nodo \(i\) esta se escribe como \[c_k(i)=\alpha\sum_j y_{i,j}c_k(j)+\beta\] donde \(\beta\) es un parámetro que controla la centralidad inicial y \(\alpha<1/\lambda_{max}\) con \(\lambda\) como valor propio de la matriz de adyacencia.

Esta medida calcula la influencia relativa de un nodo dentro de una red midiendo el número de vecinos inmediatos (nodos de primer grado) y también todos los demás nodos de la red que se conectan al nodo en cuestión a través de estos vecinos inmediatos. Sin embargo, las conexiones realizadas con vecinos distantes se penalizan mediante el factor de atenuación \(\alpha\).

Harmonic centrality

Bajo la centralidad armónica un nodo es importante cuando puede alcanzar a muchos otros nodos con distancias cortas, incluso si algunos nodos están desconectados de él. \[H(i)=\sum_{i\neq j}\frac{1}{d(i,j)}\] donde \(d(i,j)\) es la distancia más corta entre los nodos.

Esta centralidad es similar a la centralidad por cercanía, con la diferencia de que en esta no se calcula el inverso de la suma de las distancias geodésicas, sino la suma del inverso de estas. El algoritmo ayuda a responder la pregunta: ¿Que tan accesibles son los demás desde mi?

VoteRank

Este método busca identificar difusores influyentes en redes complejas de acuerdo a las puntuaciones de votación de sus nodos vecinos.

El funcionamiento del algoritmo se basa en que cada nodo \(u\) posee un par de valores (\(s_u, va_u\)): una puntuación de votación (\(s_u\)), que es el número de votos recibidos de sus vecinos, y una capacidad de votación (\(va_u\)), que es la fuerza del voto que él puede dar a sus vecinos.

Pasos del algoritmo:

1. Inicialización: Las tuplas para todos los nodos se establecen en (0,1)

2. Votación: Cada nodo vota por sus vecinos. La puntuación de votación de un nodo (\(s_u\)) se calcula como la suma de las capacidades de votación (\(va\)) de todos sus vecinos inmediatos.

3. Selección: Se elige el nodo con la puntuación de votación más alta (\(v_{max}\)) como un difusor influyente. Si un nodo es seleccionado, deja de participar en las siguientes elecciones.

4. Actualización: Se debilita la capacidad de votación de los vecinos del nodo recién elegido.

5. Se repiten los pasos 2 a 4 hasta que se haya seleccionado el número deseado (\(r\)) de difusores.

Al reducir la capacidad de votación de los vecinos de un nodo elegido, el algoritmo disminuye la probabilidad de que otros nodos en la misma vecindad sean seleccionados en las rondas siguientes, asegurando que los difusores estén dispersos por toda la red, lo cual es más efectivo para la propagación de información.

APLICACIÓN

Con base en la base de datos ICEWS del numeral anterior se van a calcular y comparar las medidas explicadas

devtools::install_github("viliamzigo/voterank")
library(voterank)

# Centralidad PageRank
pr <- igraph::page_rank(graph = g11_sub)$vector
head(sort(pr, decreasing = T), n = 5)

##      United States              Syria               Iraq Russian Federation 
##         0.07665626         0.05184290         0.04102692         0.04017149 
##             Israel 
##         0.03658062

# Centralidad Armónica
ha <- igraph::harmonic_centrality(graph = g11_sub, normalized = TRUE)
head(sort(ha, decreasing = T), n = 5)

##  United States United Kingdom          Japan        Ukraine           Iran 
##      0.4656966      0.4437390      0.4385362      0.4333333      0.4258377

# VoteRank
vr <- voterank(g = g11_sub, r = 5)
vr

## + 5/190 vertices, named, from fbec841:
## [1] United States      China              United Kingdom     Russian Federation
## [5] France

# Centralidad de Katz
ka <- igraph::alpha_centrality(graph = g11_sub)
head(sort(ka, decreasing = T), n = 5) # valores altos= nodo conectado a nodos importantes

##   Moldova, Republic of                Burundi                Hungary 
##              18.873903              17.791625              16.927165 
## Bosnia and Herzegovina                Albania 
##               8.597969               7.864083

Recordemos que la red consta de 190 países los cuales se relacionan al efectuar acciones positivas hacia otros países.

Bajo la centralidad PageRank los países más importantes resultaron ser: Estados Unidos, Siria, Iraq, Rusia e Israel. Dado que este método nos dice que un nodo es tan importante como los nodos que se enlazan con él, vemos que estos países son considerados centrales porque muchos países los han ayudado materialmente.

En armonía con PageRank, la centralidad por armonía establece a Estados Unidos como el actor más importante de la red, seguido de Reino Unido y Japón. Este ranking es lógico al observar a estos actores en la visualización de la red. Sus conexiones los posicionan en una parte central del grafo, evidenciando que alcanzan a muchos otros países sin recorrer largas distancias, por lo que su accesibilidad estructural es dominante.

Como tercera medida, el VoteRank eligió como difusores influyentes a Estados Unidos, China, Reino Unido, Rusia y Francia. Su clasificación se interpreta como que son estos los actores ‘foco’ en la red, en el sentido de que es a través de ellos que se llevan a cabo más acciones materiales positivas entre países. Además, el top 3 de difusores son a su vez los nodos con mayor grado de salida en el grafo (los que más ayudas han promovido).

Por último, la centralidad Katz nos muestra un panorama totalmente distinto: República de Moldova, Burundi, Hungria, Bosnia Herzegovina y Albania son los países más importantes bajo este criterio. Según la definición, son estos los países que están “en la mitad de la acción”, ya que se rodean de nodos centrales, sin ser los países que más ayuadn ni los que más ayudas reciben.