Introducción al anális de redes biológicas en R
Dr. Oscar Fontanelli
5/3/2022
En este documento veremos algunos elementos básicos de teoría de redes y cómo podemos usar el lenguaje de programación R - específicamente la biblioteca igraph - para hacer algunos análisis y visualizaciones básicas con este tipo de objetos. En particular, mostraremos algunos ejemplos usando una red de interacciones entre proteínas de la levadura. Primero vamos a cargar las bibliotecas que usaremos (si no están instaladas, pueden instalarse con la instrucción install.packages( )).
library(igraph) # Biblioteca para análisis de redes.
library(tidyverse) # Colección de bibliotecas para manipulación de datos.
library(viridis) # Biblioteca con paletas de colores.Comencemos con algunos conceptos básicos. Una red o grafo es un conjunto de puntos llamados nodos o vértices que pueden estar unidos o no a través de líneas o aristas. Estos objetos se usan para modelar sistemas de agentes que interactúan entre sí. Por ejemplo, en una red social los nodos pueden ser personas y las aristas representan relaciones de amistad; en una red neuronal los nodos son neuronas y dos nodos están conectados si existe una interacción entre ellas (por lo general una sinapsis); en una red de interacción proteína-proteína los nodos son proteínas y las aristas representan asociaciones entre ellas, las cuales pueden darse desde la perspectiva de la bioquímica, la transducción de señales, etc.
Vamos a comenzar creando una red aleatoria en R para ver algunas de sus características y propiedades. Para este ejemplo crearemos una red de acuerdo al modelo de Barabási-Albert, el cual genera redes aleatorias libres de escala mediante un mecanismo llamado conexión preferencial. Como podemos ver, tenemos una red con 15 nodos. Las aristas unen parejas de nodos, representando una interacción entre ellos. El grado de un nodo es el número de aristas adyacentes a él, es decir, el número de interacciones que tiene con otros nodos. Si dos nodos están unidos directamente por una arista decimos que estos nodos son vecinos. Un camino entre dos nodos es una secuencia de aristas que nos lleva de uno al otro. La longitud de un camino es el número de aristas que éste posee. La distancia entre dos nodos es la longitud del camino más corto que los une. El diámetro de una red es la distancia más grande entre cualquier par de nodos. El coeficiente de intermediación de un nodo es una medida del número de caminos más cortos en la red que pasan por ese nodo. La intermediación sirve para medir qué tanto un nodo actúa como un puente, uniendo diferentes regiones de la red.
# Generamos una red aleatoria de 15 nodos no dirigida.
# Para poder reproducirla, dado que tiene componentes aleatorios,
# fijamos la semilla del generador de números pseudoaleatorios.
set.seed(25032022)
g <- sample_pa(15,directed=FALSE)
# Una primera visualización de la red.
plot(g)
# Con las funciones V() y E() podemos acceder a los nodos (vertices)
# y aristas (edges) de la red.
V(g)## + 15/15 vertices, from b3f9d5a:
## [1] 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15E(g)## + 14/14 edges from b3f9d5a:
## [1] 1-- 2 2-- 3 2-- 4 2-- 5 5-- 6 2-- 7 5-- 8 2-- 9 9--10 8--11
## [11] 3--12 10--13 1--14 2--15# Hacemos una nueva visualización, personalizando algunos elementos.
set.seed(25032022)
plot(g,
vertex.color="coral", # color de los nodos
vertex.size=20, # tamaño de los nodos
edge.color="black") # color de las aristas
# Usamos funciones de igraph para calcular el grado de los nodos, su
# coeficiente de intermediación, la matriz de distancias y el diámetro
# de la red.
degree(g)## [1] 2 7 2 1 3 1 1 2 2 2 1 1 1 1 1betweenness(g)## [1] 13 80 13 0 35 0 0 13 24 13 0 0 0 0 0distances(g)## [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] [,11] [,12] [,13]
## [1,] 0 1 2 2 2 3 2 3 2 3 4 3 4
## [2,] 1 0 1 1 1 2 1 2 1 2 3 2 3
## [3,] 2 1 0 2 2 3 2 3 2 3 4 1 4
## [4,] 2 1 2 0 2 3 2 3 2 3 4 3 4
## [5,] 2 1 2 2 0 1 2 1 2 3 2 3 4
## [6,] 3 2 3 3 1 0 3 2 3 4 3 4 5
## [7,] 2 1 2 2 2 3 0 3 2 3 4 3 4
## [8,] 3 2 3 3 1 2 3 0 3 4 1 4 5
## [9,] 2 1 2 2 2 3 2 3 0 1 4 3 2
## [10,] 3 2 3 3 3 4 3 4 1 0 5 4 1
## [11,] 4 3 4 4 2 3 4 1 4 5 0 5 6
## [12,] 3 2 1 3 3 4 3 4 3 4 5 0 5
## [13,] 4 3 4 4 4 5 4 5 2 1 6 5 0
## [14,] 1 2 3 3 3 4 3 4 3 4 5 4 5
## [15,] 2 1 2 2 2 3 2 3 2 3 4 3 4
## [,14] [,15]
## [1,] 1 2
## [2,] 2 1
## [3,] 3 2
## [4,] 3 2
## [5,] 3 2
## [6,] 4 3
## [7,] 3 2
## [8,] 4 3
## [9,] 3 2
## [10,] 4 3
## [11,] 5 4
## [12,] 4 3
## [13,] 5 4
## [14,] 0 3
## [15,] 3 0diameter(g)## [1] 6En ocasiones, la arista entre dos nodos tiene una dirección, indicando una direccionalidad en la interacción. En este caso hablamos de una red dirigida. Esta direccionalidad se representa con una flecha en las aristas. En una red dirigida hay que distinguir entre tres tipos de grados: el grado de entrada, que es el número de aristas que entran hacia un cierto nodo; el grado de salida, que es el número de aristas que salen desde un cierto nodo, y el grado total, la suma de los grados de entrada y de salida. La intermediación de un nodo, la distancia entre dos nodos y el diámetro de la red dependen también de si consideramos o no la dirección de las aristas.
set.seed(25032022)
g <- sample_pa(10,directed=TRUE)
set.seed(25032022)
plot(g,
vertex.color="coral",
vertex.size=20,
edge.color="black")
# Cuando la red es dirigida tenemos que distinguir entre el grado de entrada,
# el grado de salida y el grado total.
degree(g, mode="in")## [1] 4 2 0 1 1 0 0 0 1 0degree(g, mode="out")## [1] 0 1 1 1 1 1 1 1 1 1degree(g, mode="all")## [1] 4 3 1 2 2 1 1 1 2 1# De igual modo, hay que distinguir entre el grado de intermediación con y
# sin dirección en las aristas.
betweenness(g,directed=TRUE)## [1] 0 2 0 1 1 0 0 0 1 0betweenness(g,directed=FALSE)## [1] 30 15 0 8 8 0 0 0 8 0# Las distancias entre pares de nodos dependen de si consideramos o no la
# dirección. Lo mismo ocurre con el diámetro de la red.
distances(g,mode="in")## [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10]
## [1,] 0 1 2 1 1 2 2 2 1 2
## [2,] Inf 0 1 Inf Inf Inf 1 Inf Inf Inf
## [3,] Inf Inf 0 Inf Inf Inf Inf Inf Inf Inf
## [4,] Inf Inf Inf 0 Inf Inf Inf 1 Inf Inf
## [5,] Inf Inf Inf Inf 0 1 Inf Inf Inf Inf
## [6,] Inf Inf Inf Inf Inf 0 Inf Inf Inf Inf
## [7,] Inf Inf Inf Inf Inf Inf 0 Inf Inf Inf
## [8,] Inf Inf Inf Inf Inf Inf Inf 0 Inf Inf
## [9,] Inf Inf Inf Inf Inf Inf Inf Inf 0 1
## [10,] Inf Inf Inf Inf Inf Inf Inf Inf Inf 0distances(g,mode="out")## [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10]
## [1,] 0 Inf Inf Inf Inf Inf Inf Inf Inf Inf
## [2,] 1 0 Inf Inf Inf Inf Inf Inf Inf Inf
## [3,] 2 1 0 Inf Inf Inf Inf Inf Inf Inf
## [4,] 1 Inf Inf 0 Inf Inf Inf Inf Inf Inf
## [5,] 1 Inf Inf Inf 0 Inf Inf Inf Inf Inf
## [6,] 2 Inf Inf Inf 1 0 Inf Inf Inf Inf
## [7,] 2 1 Inf Inf Inf Inf 0 Inf Inf Inf
## [8,] 2 Inf Inf 1 Inf Inf Inf 0 Inf Inf
## [9,] 1 Inf Inf Inf Inf Inf Inf Inf 0 Inf
## [10,] 2 Inf Inf Inf Inf Inf Inf Inf 1 0distances(g,mode="all")## [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10]
## [1,] 0 1 2 1 1 2 2 2 1 2
## [2,] 1 0 1 2 2 3 1 3 2 3
## [3,] 2 1 0 3 3 4 2 4 3 4
## [4,] 1 2 3 0 2 3 3 1 2 3
## [5,] 1 2 3 2 0 1 3 3 2 3
## [6,] 2 3 4 3 1 0 4 4 3 4
## [7,] 2 1 2 3 3 4 0 4 3 4
## [8,] 2 3 4 1 3 4 4 0 3 4
## [9,] 1 2 3 2 2 3 3 3 0 1
## [10,] 2 3 4 3 3 4 4 4 1 0diameter(g,directed=TRUE)## [1] 2diameter(g,directed=FALSE)## [1] 4Ahora descarguemos una red de interacción entre proteínas de la levadura. Esta red contiene 7182 interacciones entre 2361 proteínas de la levadura. Estas interacciones fueron seleccionadas por von Mering et.al [1] de entre una lista más grande de interacciones, de tal modo que se eliminaran la mayoría de los falsos positivos. Posteriormente, Bu et al. [2] utilizaron estos datos para ilustrar una metodología novedosa, basada en teoría de redes, para detectar grupos de proteínas con funciones biológicamente similares.
Estos datos vienen en un formato de lista de aristas, es decir, es una tabla indicando cada arista de la red. Para cada arista vemos los dos nodos, los cuales son proteínas, indicadas por su nombre sistemático. Primero tenemos que leer esta tabla y luego convertirla en una red, para lo cual usaremos la función graph_from_data_frame( ). De este modo tendremos una red donde los nodos son proteínas de la levadura y las aristas indican interacciones entre ellas.
# Descargamos la lista de aristas para la red de interacciones entre proteínas.
id <- "1O5NgBBewLPHGpKDFZfNmNVu0BY-Y-dDh"
yeast <- read.csv(sprintf("https://docs.google.com/uc?id=%s&export=download", id))
head(yeast)## from to
## 1 YBR236C YBR236C
## 2 YBR236C YOR151C
## 3 YBR236C YML010W
## 4 YOR151C YER095W
## 5 YOR151C YIL021W
## 6 YOR151C YGL070C# Convertimos esta lista de aristas en una red.
yeast %>%
graph_from_data_frame -> g
# Primer vistazo a la red
g## IGRAPH b4ddc00 DN-- 2361 7182 --
## + attr: name (v/c)
## + edges from b4ddc00 (vertex names):
## [1] YBR236C->YBR236C YBR236C->YOR151C YBR236C->YML010W YOR151C->YER095W
## [5] YOR151C->YIL021W YOR151C->YGL070C YOR151C->YDR404C YOR151C->YAL005C
## [9] YOR151C->YDR167W YOR151C->YBL039C YML010W->YJR017C YML010W->YBL046W
## [13] YML010W->YNL201C YML010W->YPR041W YML010W->YDR404C YNR016C->YNR016C
## [17] YNR016C->YLR386W YNR016C->YDL047W YNR016C->YJR064W YNR016C->YDL225W
## [21] YNR016C->YPR110C YNR016C->YMR106C YNR016C->YER179W YNR016C->YDR398W
## [25] YLR386W->YIL061C YLR386W->YBR143C YLR386W->YDL017W YDL047W->YDL047W
## [29] YDL047W->YJR105W YDL047W->YDL029W YDL047W->YJL026W YDL047W->YER155C
## + ... omitted several edges# Nodos de la red.
V(g)## + 2361/2361 vertices, named, from b4ddc00:
## [1] YBR236C YOR151C YML010W YNR016C YLR386W YDL047W YJR064W
## [8] YAR015W YNL220W YLR359W YOR361C YMR300C YMR235C YAL012W
## [15] YGL234W YGR061C YDR226W YJR105W YBR059C YPR180W YLR249W
## [22] YDR390C YOL086C YLR127C YNL172W YHR166C YKL022C YBR151W
## [29] YBR128C YDL185W YKL135C YPL259C YBR084W YPR010C YBL037W
## [36] YOL062C YPL195W YBR288C YGL019W YGR261C YHR174W YDR441C
## [43] YBR149W YIL062C YJR065C YNR035C YKL013C YLR370C YDL029W
## [50] YGR240C YBR234C YHR013C YGR090W YDL040C YHR023W YMR116C
## [57] YGR254W YDL137W YBR164C YDR127W YPL235W YOR136W YKL104C
## [64] YMR186W YDR427W YGR234W YGL048C YJL008C YBR249C YNR053C
## + ... omitted several verticesget.data.frame(g, what="vertices") %>% head()## name
## YBR236C YBR236C
## YOR151C YOR151C
## YML010W YML010W
## YNR016C YNR016C
## YLR386W YLR386W
## YDL047W YDL047W# Aristas de la red.
E(g)## + 7182/7182 edges from b4ddc00 (vertex names):
## [1] YBR236C->YBR236C YBR236C->YOR151C YBR236C->YML010W YOR151C->YER095W
## [5] YOR151C->YIL021W YOR151C->YGL070C YOR151C->YDR404C YOR151C->YAL005C
## [9] YOR151C->YDR167W YOR151C->YBL039C YML010W->YJR017C YML010W->YBL046W
## [13] YML010W->YNL201C YML010W->YPR041W YML010W->YDR404C YNR016C->YNR016C
## [17] YNR016C->YLR386W YNR016C->YDL047W YNR016C->YJR064W YNR016C->YDL225W
## [21] YNR016C->YPR110C YNR016C->YMR106C YNR016C->YER179W YNR016C->YDR398W
## [25] YLR386W->YIL061C YLR386W->YBR143C YLR386W->YDL017W YDL047W->YDL047W
## [29] YDL047W->YJR105W YDL047W->YDL029W YDL047W->YJL026W YDL047W->YER155C
## [33] YDL047W->YIL061C YDL047W->YDR188W YDL047W->YGL246C YDL047W->YKL139W
## [37] YDL047W->YDL171C YDL047W->YFR019W YDL047W->YPL204W YDL047W->YMR024W
## + ... omitted several edgesget.data.frame(g, what="edges") %>% head()## from to
## 1 YBR236C YBR236C
## 2 YBR236C YOR151C
## 3 YBR236C YML010W
## 4 YOR151C YER095W
## 5 YOR151C YIL021W
## 6 YOR151C YGL070Cplot(g,
vertex.size=3,
vertex.color="blue",
vertex.label=NA, # No queremos ver las etiquetas de los nodos.
edge.arrow.size=0) # No queremos ver las puntas de las flechitas.
Esta red tiene varios nodos que se conectan con ellos mismos. Estas autoconexiones no tienen mucho significado biológico, así que es buena idea suprimirlas. Podemos quitar estas auto-conexiones con la función igraph::simplify( ). Además de esto, la red viene como una red dirigida. Como en este caso las interacciones no tienen una dirección o un sentido, quitamos esta direccionalidad con la función as.undirected( ).
# Usamos la función simplify() de igraph para quitar conexiones múltiples
# y autoconexiones. Usamos la función as.unidirected() para quitar las
# direcciones de las aristas.
g <- igraph::simplify(g)
g <- as.undirected(g)
g## IGRAPH b92eaa2 UN-- 2361 6646 --
## + attr: name (v/c)
## + edges from b92eaa2 (vertex names):
## [1] YBR236C--YOR151C YBR236C--YML010W YNR016C--YLR386W YNR016C--YDL047W
## [5] YNR016C--YJR064W YLR359W--YOR361C YMR300C--YMR235C YMR300C--YAL012W
## [9] YDL047W--YJR105W YPR180W--YLR249W YPR180W--YDR390C YPR180W--YOL086C
## [13] YLR127C--YNL172W YLR127C--YHR166C YLR127C--YKL022C YBR128C--YDL185W
## [17] YKL135C--YPL259C YKL135C--YBR084W YKL135C--YPR010C YBL037W--YOL062C
## [21] YPL195W--YBR288C YPL195W--YGL019W YPL195W--YGR261C YBR288C--YGR261C
## [25] YGR261C--YHR174W YIL062C--YJR065C YIL062C--YNR035C YJR065C--YNR035C
## [29] YIL062C--YKL013C YJR065C--YKL013C YNR035C--YKL013C YOR361C--YLR370C
## + ... omitted several edgesSi tomamos de una red sólo un subcojnunto de nodos y todas las conexiones entre ellos, nos queda lo que se conoce como un subgrafo inducido. El k-núcleo o k-core de una red es el subgrafo inducido por los nodos de grado mayor o igual a k. Vamos a extraer el 4-core de esta red, es decir, el subgrafo inducido por los nodos que tienen cuatro o más vecinos. Cuando tenemos redes muy grandes, suele ser buena idea comenzar por analizar sus núcleos, pues éstos son más pequeños y contienen los nodos y las interacciones más importantes de la red.
# Extraemos el núcleo (core) de la red formado por los nodos que tienen
# grado mayor a 3.
core <- coreness(g,mode="all") # Vemos a qué core pertenecen los nodos.
core <- core[core>3] # Extraemos los nodos del core 3.
g <- induced_subgraph(g, names(core)) # Subgrafo inducido por los nodos en el vector core.
g## IGRAPH b9396a0 UN-- 849 4406 --
## + attr: name (v/c)
## + edges from b9396a0 (vertex names):
## [1] YNR016C--YLR386W YNR016C--YDL047W YNR016C--YJR064W YIL062C--YJR065C
## [5] YIL062C--YNR035C YJR065C--YNR035C YIL062C--YKL013C YJR065C--YKL013C
## [9] YNR035C--YKL013C YOR361C--YLR370C YIL062C--YLR370C YJR065C--YLR370C
## [13] YNR035C--YLR370C YKL013C--YLR370C YDL047W--YDL029W YJR064W--YDL029W
## [17] YIL062C--YDL029W YJR065C--YDL029W YNR035C--YDL029W YKL013C--YDL029W
## [21] YLR370C--YDL029W YLR370C--YGR240C YIL062C--YBR234C YJR065C--YBR234C
## [25] YNR035C--YBR234C YKL013C--YBR234C YLR370C--YBR234C YDL029W--YBR234C
## [29] YBR084W--YHR013C YGL019W--YGR090W YHR013C--YGR090W YHR013C--YDL040C
## + ... omitted several edgesPodemos agregar atributos a los nodos de una red. Vamos a agregar a cada nodo un atributo con su grado y uno más con su coeficiente de intermediación.
V(g)$degree <- degree(g)
V(g)$betweenness <- betweenness(g)
get.data.frame(g,what="vertices") %>% head() ## name degree betweenness
## YOR151C YOR151C 6 564.8803
## YML010W YML010W 4 229.7296
## YNR016C YNR016C 8 696.2778
## YLR386W YLR386W 4 148.2834
## YDL047W YDL047W 22 4522.1635
## YJR064W YJR064W 14 1992.0892Ahora podemos visualizar la red de tal modo que el tamaño de los nodos sea proporcional o que esté en función de alguno de estos atributos. Al inicio no tendremos una visualización muy clara, pero podemos ajustar el tamaño de los nodos hasta alcanzar un mejor resultado.
# Visualizamos de modo que el tamaño de cada nodo sea proporcional a
# su grado.
plot(g,
vertex.size=V(g)$degree,
vertex.color="blue",
vertex.label=NA)
# Ajustamos el tamaño de los nodos.
plot(g,
vertex.size=2*sqrt(V(g)$degree),
vertex.color="blue",
vertex.label=NA)
# Lo mismo, pero el tamaño en función de la intermediación.
plot(g,
vertex.size=log(V(g)$betweenness+1),
vertex.color="blue",
vertex.label=NA)
Nos preguntamos cuáles son las proteínas más importantes en esta red. Por un lado, podemos ver cuáles son las proteínas que interactúan con más proteínas, es decir, cuáles son los nodos con mayor grado. De igual modo, podemos ver cuáles son los nodos con mayor coeficiente de intermediación.
# Vemos cuáles son los nodos de mayor grado y mayor intermediación.
degree(g) %>% sort(decreasing=TRUE) %>% head(n=15)## YIL035C YDL213C YPR110C YPL204W YGR090W YMR049C YDR394W YLR074C YCR057C YMR106C
## 55 49 47 45 44 44 44 43 43 43
## YNL061W YKR081C YNL189W YGL137W YPR016C
## 42 41 41 41 39get.data.frame(g,what="vertices") %>%
as_tibble() %>%
arrange(-degree) %>%
head(n=15) %>%
pull(name)## [1] "YIL035C" "YDL213C" "YPR110C" "YPL204W" "YGR090W" "YMR049C" "YDR394W"
## [8] "YLR074C" "YCR057C" "YMR106C" "YNL061W" "YKR081C" "YNL189W" "YGL137W"
## [15] "YPR016C"get.data.frame(g,what="vertices") %>%
as_tibble() %>%
arrange(-betweenness) %>%
head(n=15) %>%
pull(name)## [1] "YNL189W" "YIL035C" "YPR110C" "YPL204W" "YGL137W" "YMR106C" "YDL213C"
## [8] "YBR055C" "YBR009C" "YGR090W" "YLR074C" "YER133W" "YOL139C" "YMR012W"
## [15] "YCR057C"Una comunidad en una red es un grupo de nodos que muestra un mayor número de conexiones entre ellos que hacia otros nodos de la red. Por lo general, los nodos forman comunidades cuando tienen ciertas características en común. Buscar comunidades en una red no es una tarea trivial y existen varias metodologías y algoritmos para hacerlo. Veamos un algoritmo muy conocido y muy utilizado que es el algoritmo de Louvain.
# Búsqueda de comunidades mediante el algoritmo de Louvain.
coms <- cluster_louvain(g)
coms## IGRAPH clustering multi level, groups: 13, mod: 0.52
## + groups:
## $`1`
## [1] "YOR151C" "YML010W" "YPR010C" "YPL160W" "YDL126C" "YBL058W" "YOR207C"
## [8] "YKL213C" "YOR116C" "YOR341W" "YJR017C" "YNL135C" "YPR110C" "YOR319W"
## [15] "YPR190C" "YDR037W" "YNR003C" "YMR205C" "YNL201C" "YCR079W" "YJR063W"
## [22] "YOR224C" "YBR154C" "YPR187W" "YOR210W" "YIL021W" "YGR186W" "YGL070C"
## [29] "YDR404C" "YOL005C" "YDL140C" "YKL144C" "YNL113W" "YLL034C" "YKL085W"
## [36] "YHR193C" "YGL026C" "YGL062W" "YER125W" "YOR317W" "YDR129C" "YDR353W"
## [43] "YDR343C" "YER052C" "YBR088C" "YOR335C" "YDL190C" "YPR111W" "YDR143C"
## [50] "YKR025W" "YDL150W" "YPR093C" "YGR005C" "YOL128C"
##
## + ... omitted several groups/verticescoms$membership## [1] 1 1 2 3 2 2 4 5 3 6 7 7 8 9 9 1 9 7 4 4 4 4 4 8 3
## [26] 4 2 2 2 3 2 8 8 8 8 3 8 3 8 9 9 3 10 5 6 4 4 4 4 9
## [51] 9 9 9 9 9 9 9 9 9 9 9 9 9 9 9 9 9 8 1 8 11 10 11 11 11
## [76] 11 11 11 11 11 10 11 11 11 11 11 9 2 2 9 2 9 4 9 9 9 4 7 8 8
## [101] 8 2 5 6 8 7 12 12 12 3 3 3 3 3 9 7 9 2 8 8 7 8 7 1 13
## [126] 1 6 10 2 12 8 3 8 10 3 4 3 9 9 9 9 9 9 9 9 9 4 11 11 11
## [151] 11 11 11 11 11 11 11 11 11 11 11 4 4 7 8 3 4 8 3 6 6 6 6 6 6
## [176] 6 6 6 6 6 7 7 4 1 4 9 7 7 4 6 6 9 8 8 8 8 11 11 3 2
## [201] 10 10 10 10 10 10 10 10 10 10 3 3 3 3 3 9 7 7 11 11 11 11 8 11 11
## [226] 3 3 3 7 3 11 10 10 3 3 1 1 1 5 5 8 2 7 7 3 3 2 2 2 2
## [251] 2 2 2 9 2 2 2 8 2 2 2 7 7 7 7 9 9 9 9 8 5 5 9 9 9
## [276] 7 7 1 6 1 5 8 9 9 9 12 5 7 5 5 5 5 5 5 5 11 11 11 2 11
## [301] 11 5 3 3 8 9 4 4 6 9 9 2 9 5 9 3 1 8 11 8 3 1 8 8 8
## [326] 8 8 8 8 9 7 7 9 3 1 9 3 8 5 7 7 5 8 8 8 8 2 2 2 2
## [351] 2 2 2 2 4 9 2 2 2 1 3 3 8 8 8 8 11 11 11 11 11 11 11 11 8
## [376] 12 11 11 11 11 11 9 5 10 10 4 4 7 7 7 7 7 7 7 7 7 7 7 7 7
## [401] 7 7 7 7 7 7 7 13 13 13 13 13 12 12 12 12 2 12 12 12 12 12 12 2 12
## [426] 12 12 7 9 3 3 8 3 3 3 3 7 9 3 3 12 11 3 10 3 3 3 4 7 9
## [451] 7 1 2 9 9 9 8 3 12 7 5 5 4 6 9 8 8 8 8 1 8 3 8 8 8
## [476] 8 8 8 8 8 8 8 8 2 1 8 1 9 11 11 11 11 4 4 9 6 2 8 4 4
## [501] 11 11 9 9 9 2 8 2 2 8 8 8 8 3 8 8 3 7 7 3 7 4 1 1 1
## [526] 1 1 1 1 1 1 1 1 7 1 7 4 3 4 9 9 8 8 8 8 8 8 1 8 10
## [551] 3 11 5 3 6 9 8 7 3 7 3 8 3 7 2 3 3 5 11 11 6 6 6 2 5
## [576] 7 7 7 9 9 9 3 11 10 4 11 7 7 7 7 3 2 8 8 3 3 3 4 4 8
## [601] 7 12 8 8 8 8 3 3 3 3 7 3 3 9 3 2 8 11 1 9 8 8 8 3 7
## [626] 3 3 3 8 2 3 3 1 2 3 3 1 8 9 8 3 8 3 6 6 3 3 3 9 3
## [651] 8 3 3 8 8 8 8 3 1 2 1 6 3 7 8 9 3 1 3 9 3 1 8 10 2
## [676] 3 2 10 3 3 8 3 3 11 3 3 8 7 3 8 8 8 2 8 3 8 8 9 1 1
## [701] 9 3 3 3 5 5 3 3 3 3 3 2 2 3 3 8 2 3 8 8 3 2 5 8 1
## [726] 8 3 1 8 11 1 9 1 3 9 1 8 2 2 9 3 9 8 10 8 1 2 8 3 3
## [751] 8 2 2 3 8 3 8 9 9 9 3 2 3 3 3 3 3 9 8 3 9 9 8 9 12
## [776] 10 9 9 3 9 9 8 8 3 3 9 8 8 9 11 6 3 1 3 8 7 12 8 2 11
## [801] 8 9 4 5 3 2 3 3 7 2 8 8 8 6 3 9 2 2 8 3 5 4 4 9 9
## [826] 10 9 11 2 5 5 12 3 2 5 3 11 6 11 3 3 1 1 7 1 9 1 3 1Podemos colorear los nodos de la red de acuerdo a la comunidad a la que pertenecen. De este modo, podemos ver grupos de proteínas que están fuertemenete relacionadas entre sí.
# Agregamos a los nodos un atributo que nos indique a qué comunidad
# pertenecen.
V(g)$community <- coms$membership
# Coloreamos los nodos de acuerdo a su comunidad y probamos con diferentes
# layouts.
plot(g,
vertex.size=2*sqrt(V(g)$degree),
vertex.color=V(g)$community,
vertex.label=NA,
layout=layout_with_kk)
plot(g,
vertex.size=2*sqrt(V(g)$degree),
vertex.color=V(g)$community,
vertex.label=NA,
layout=layout_with_fr)
Un layout es una disposición de los nodos en el espacio en que los graficamos. Existen varias técnicas y metodologías para colocar los nodos en el espacio. Aquí veremos la misma red graficada con nueve layouts diferentes.
# Pequeño script para ver la misma red con nueve diferentes layouts.
# El último de ellos (layout_with_gem) puede tardar un poco en correr.
layouts <- c("layout_nicely","layout_with_kk","layout_in_circle",
"layout_on_grid","layout_on_sphere","layout_randomly",
"layout_with_fr","layout_as_star","layout_with_gem")
par(mfrow=c(3,3),mar=c(1,1,1,1))
for(i in 1:length(layouts)){
LO <- layouts[i]
l <- do.call(LO,list(g))
plot(g,
vertex.size=2*sqrt(V(g)$degree),
vertex.color=V(g)$community,
vertex.label=NA,
layout=l,
main=LO)
}
En ocasiones es recomendable usar paletas de colores que se distingan muy claramente. Las paletas de colores de la colección viridis son casi siempre adecuadas, pues estas paletas están diseñadas para ser percibidas y distinguidas por personas con diferentes tipos de daltonismo.
# Seleccionamos una paleta de color de la colección viridis.
# Primero visualizamos la paleta.
pal <- turbo(max(V(g)$community))
x <- seq(1,length(pal))
y <- rep(1,length(pal))
plot(x,y,cex=6,pch=19, col=pal)
Coloreamos los nodos de la red con los colores de nuestra paleta personalizada.
# Asignamos estos colores a los nodos de la red y visualizamos.
V(g)$color <- pal[V(g)$community]
plot(g,
vertex.size=2*sqrt(V(g)$degree),
vertex.color=V(g)$color,
vertex.label=NA,
layout=layout_nicely)
Una forma en que podemos hacer que las comunidades se separen en la visualización es agregando pesos a las aristas. Si agregamos un peso grande a las aristas que unen nodos de la misma comunidad y un peso pequeño a las aristas que unen nodos de diferentes comunidades y luego de esto usamos un layout que considere el peso de las aristas, tendremos entonces una visualización donde las comunidades no quedarán revueltas y será más fácil distinguirlas.
# Primero creamos una tabla con los nodos y su comunidad.
membership.frame <- tibble(name=V(g)$name,community=V(g)$community)
# Peso de las aristas entre nodos de la misma comunidad.
edge.weight <- 50
# Lista de aristas.
edge.list <- get.edgelist(g)
# Un vector vacío para ir guardando los pesos.
weights <- c()
# Si los dos nodos de una aristas pertenecen a la comunidad,
# asignar el peso edge.weight. De lo contrario, asignar un peso de 1.
# Iterar sobre todas las aristas.
for(i in 1:nrow(edge.list)){
com.1 <- membership.frame[membership.frame$name == edge.list[i,][1],]$community
com.2 <- membership.frame[membership.frame$name == edge.list[i,][2],]$community
if(com.1 == com.2){
weight <- edge.weight
}else{
weight <- 1
}
weights <- c(weights,weight)
}
# Asignar estos pesos a las aristas.
E(g)$weight <- weights
get.data.frame(g, what="edges") %>% head(n=30)## from to weight
## 1 YNR016C YLR386W 1
## 2 YNR016C YDL047W 50
## 3 YNR016C YJR064W 50
## 4 YIL062C YJR065C 50
## 5 YIL062C YNR035C 50
## 6 YJR065C YNR035C 50
## 7 YIL062C YKL013C 50
## 8 YJR065C YKL013C 50
## 9 YNR035C YKL013C 50
## 10 YOR361C YLR370C 50
## 11 YIL062C YLR370C 50
## 12 YJR065C YLR370C 50
## 13 YNR035C YLR370C 50
## 14 YKL013C YLR370C 50
## 15 YDL047W YDL029W 1
## 16 YJR064W YDL029W 1
## 17 YIL062C YDL029W 1
## 18 YJR065C YDL029W 1
## 19 YNR035C YDL029W 1
## 20 YKL013C YDL029W 1
## 21 YLR370C YDL029W 1
## 22 YLR370C YGR240C 1
## 23 YIL062C YBR234C 50
## 24 YJR065C YBR234C 50
## 25 YNR035C YBR234C 50
## 26 YKL013C YBR234C 50
## 27 YLR370C YBR234C 50
## 28 YDL029W YBR234C 1
## 29 YBR084W YHR013C 1
## 30 YGL019W YGR090W 1# Visualizamos de nuevo.
set.seed(25032022)
plot(g,
vertex.size=2*sqrt(V(g)$degree),
vertex.color=V(g)$color,
vertex.label=NA,
layout=layout_nicely,
edge.color="gray80") # color de las aristas
Podemos también colorear las aristas de acuerdo a la comunidad de los nodos. En este ejemplo vamos a poner en gris las aristas que unen nodos de diferentes comunidades y vamos a asignar a las aristas que unen nodos de la misma comunidad el mismo color de estos nodos.
# Creamos una tabla con las aristas. Para cada aristas nos fijamos
# en sus dos nodos. Si pertenecen a la misma comunidad, le ponemos
# a la aristas el color de los nodos. Si no, la ponemos gris.
# Iteramos sobre todas las aristas. Visualizamos.
edge.list.frame <- as.data.frame(get.edgelist(g))
E(g)$color <- map2_chr(edge.list.frame$V1, edge.list.frame$V2, ~ {
ifelse(
V(g)$color[V(g)[.x]] ==
V(g)$color[V(g)[.y]],
V(g)$color[V(g)[.x]],
"gray70")
})
plot(g,
vertex.size=2*sqrt(V(g)$degree),
vertex.color=V(g)$color,
vertex.label=NA,
layout=layout_nicely,
edge.color=E(g)$color)
Por último guardamos la red en formato graphml para poder usarla en el futuro.
# Guardamos la red en formato graphml para su uso futuro.
write.graph(g,"D:/Bio/yeast_protein_interaction.graphml", format="graphml")Referencias
- Von Mering, C., Krause, R., Snel, B., Cornell, M., Oliver, S. G., Fields, S., & Bork, P. (2002). Comparative assessment of large-scale data sets of protein–protein interactions. Nature, 417(6887), 399-403.
- Bu, D., Zhao, Y., Cai, L., Xue, H., Zhu, X., Lu, H., … & Chen, R. (2003). Topological structure analysis of the protein–protein interaction network in budding yeast. Nucleic acids research, 31(9), 2443-2450.