Le TP précédent vous a amené à calculer des métriques globales sur un graphe. Ce TP-ci va présenter les propriétés locales.

Vous aurez à rendre un compte rendu constitué d’une suite d’intructions R et de réponses aux questions.

Préliminaires

install.packages("igraph")

Trucs et astuces

Préliminaires

Importer le paquetage dans votre environnement

library(igraph)

Changer de répertoire de travail

Modifiez cette ligne pour vous placer dans votre répertoire préféré

# setwd("~/sync/Enseignement/GRIB/igraph")

Graphes aléatoires

Plutôt que de construire des graphes manuellement ou travailler sur des données réelles, il est possible de générer des graphes aléatoires à partir d’un modèle de graphe. Nos ne détaillerons pas les différents modèles probabilistes de graphes. Sachez qu’ils sont le pendant des modèles probabilistes de séquences (cf. chaînes de Markov).

Créez les trois graphes ci-dessous et visualisez-les :

er_graph <- erdos.renyi.game(100, 2/100)
ws_graph <- watts.strogatz.game(1, 100, 4, 0.05)
ba_graph <- barabasi.game(100, directed = FALSE)
plot(er_graph, vertex.label=NA, vertex.size=3)
plot(ws_graph, vertex.label=NA, vertex.size=3)
plot(ba_graph, vertex.label=NA, vertex.size=3)

Question 1 : Qu’est-ce qu’un graphe “petit-monde” ? Un graphe “sans-échelle” ? Quel modèle génère des graphes “petits-mondes” ? Quel modèle génère des graphes “sans-échelle” ?

Question 2 : Dessinez la distribution des degrés de plusieurs réalisations de chaque modèle. D’après les définitions, est-ce ce à quoi vous vous attendez ? Pour certains modèles de graphes, la distribution des degrés doit suivre des lois spécifiques, quelles sont-elles ?

Dessins de graphes

Question 3 : Exécutez et expliquez le code suivant. Quel est le principe de l’algorithme de Kamada-Kawai ?

g <- graph.lattice( c(10,10) )
E(g)$weight <- runif(ecount(g))
E(g)$color <- "grey"
E(g)[ weight > 0.9 ]$color <- "red"
plot(g, vertex.size=2, vertex.label=NA, layout=layout.kamada.kawai, edge.width=2+3*E(g)$weight)

Question 4 : Dans g, quel est le degré des noeuds le plus fréquent ? Pour effectuer ce calcul, il est nécessaire de calculer:

Donner le code R.

Question 5 : Quel algorithme de dessin de igraph vous permettrait d’obtenir un dessin “carré” régulier de g ?

Il existe différents algorithmes de dessin, permettant de mettre en évidence certaines propriétés.

Question 6 : En vous référant à la documentation de igraph, pour chacun des trois graphes er_graph, ws_graph, ba_graph sélectionnez un autre algorithme de dessin qui apporte, selon vous, une visualisation intéressante.

Étude locale de graphe

Question 7 : Récupérez le réseau d’interactions PPI g : http://www.dil.univ-mrs.fr/~tichit/bin/tp2/NetworkDatasets/JBA.gr, puis supprimez les boucles.

Ensuite, nous pouvons calculer le maximum (ou minimum) de la betweenness, de la closeness, et de l’edge-betweenness :

maxb <- which(betweenness(g) == max(betweenness(g)))
maxc <- which(closeness(g) == max(closeness(g)))
maxeb <- which(edge.betweenness(g) == max(edge.betweenness(g)))

Les trois fonctions betweenness, closeness et edge.betweenness ci-dessus retournent un vecteur de taille |V| (|E| pour l’edge-betweenness). Elles nous permettent de trouver les sommets et arêtes relativement centraux (ou excentrés).

Question 8 : Expliquez succinctement la définition des ces trois fonctions.

En ce qui concerne l’edge-betweenness, nous obtenons le numéro de l’arête de plus haute betweenness. Il serait intéressant de trouver quels sont ses deux sommets incidents.

Question 9 : Quelle est la fonction R qui retourne le couple de sommets incidents à une arête donnée ?

Trouver le centre du graphe

Question 10 : Comment calculer l’excentricité des sommets du graphe ? Comment retourner la liste (le vecteur) des sommets appartenant au centre du graphe ? Puis ceux d’exentricité maximale ?

Question 11 : Comment calculer le sous-graphe formé du centre du graphe ?

K-core

Il est également possible de calculer les K-core, et de ne conserver que l’ensemble des commets les plus centraux.

Appliquez la fonction graph.coreness() à un graphe aléatoire de type Erdős-Rényi.

Question 12 : Donnez le code R permettant de supprimer tous les sommets n’appartenant pas au core d’indice maximal.

Déconnecter le graphe

Un graphe de trop grande taille est relativement difficile à traiter (problème de visualisation/navigation et problème de temps de calcul).

Nous allons produire un graphe g3 dans lequel l’arête d’edge-betweenness maximale a été elevée. Nous itérons le processus jusqu’à ce que notre graphe soit déconnecté. Testez. Utilisons la fonction is.connected() à chaque étape pour vérifier si notre graphe est toujours connecté ou non.

g3 <- g
i <- 0                        # pour afficher le nombre de tours
while(is.connected(g3)) {
  i <- i+1                    # pour afficher le nombre de tours
  message(cbind("round ", i)) # pour afficher le nombre de tours
  g3 <- delete_edges(g3, which(edge.betweenness(g3) == max(edge.betweenness(g3))))
}

Ce processus nous permet de retirer un minimum d’arêtes très centrales, afin de déconnecter le graphe en deux parties de taille similaire.

Pour savoir quels sont les sommets situés dans chaque partie :

clusters(g3)

Question 13 : Afin de pouvoir travailler sur chaque partie du graphe, trouvez une fonction qui décompose notre graphe g3 (formé de deux composantes connexes) en deux graphes distincts.

Visualisation de communautés

Exécutez le code suivant :

com <- spinglass.community(g, spins=3)
# V(g)$color <- com$membership+1
g <- set.graph.attribute(g, "layout", layout.kamada.kawai(g))
plot(g, vertex.label=NA, vertex.size=5,vertex.color=com$membership+1)

Enrichissement des communautés

Question 13bis : Exécutez le code ci-dessous. Documentez-vous sur le package gProfileR et expliquez le code ci-dessous. Pour chaque cluster, indiquez les 2 ou 3 différents processus biologiques dans lesquels un grand nombre de gènes du cluster sont impliqués.

# install.packages("gProfileR")
library("gProfileR")
genes_c1 <- V(g)$name[which(com$membership == 1)]
results_c1 <- gprofiler(genes_c1, organism = "dmelanogaster",
    hier_filtering="strong",src_filter=c("GO:BP","GO:CC","REAC"),
    correction_method="fdr", max_p_value =0.001)
genes_c2 <- V(g)$name[which(com$membership == 2)]
results_c2 <- gprofiler(genes_c2, organism = "dmelanogaster",
    hier_filtering="strong",src_filter=c("GO:BP","GO:CC","REAC"),
    correction_method="fdr", max_p_value =0.001)
genes_c3 <- V(g)$name[which(com$membership == 3)]
results_c3 <- gprofiler(genes_c3, organism = "dmelanogaster",
    hier_filtering="strong",src_filter=c("GO:BP","GO:CC","REAC"),
    correction_method="fdr", max_p_value =0.001)

Question 14 : Lisez http://stackoverflow.com/questions/9471906/what-are-the-differences-between-community-detection-algorithms-in-igraph et expliquez au mieux les 3 méthodes de détection des communautés qui vous plaisent le plus.

Interaction protéine-protéine, Co-expression, Complexes protéiques

Question 15 : Cette question est la suite de la question 15 du TP précédent. Utilisez au mieux votre script pour mesurer et classifier l’ensemble des réseaux PPI fournis dans la section Data ci-dessous ainsi que le réseau CORUM. Vous devez donc regrouper les réseaux dans certaines classes selon leur “ressemblance” et expliquer le pourquoi de votre démarche.

Des exemples de rendu sont affichés ci-dessous. Vous pouvez par exemple utiliser les fonctions dist, hclust, etc. pour faire du clustering sur votre dataframe contenant les résultats. Une seconde possibilité est de faire une analyse en composantes principales, par exemple grâce au package FactoMineR, la fonction PCA. Évitez de travailler sur la maille si vous faites une ACP (rappel: une ACP sélectionne les composantes qui expliquent le mieux la variance des données. I.e. les variables dont la variance est la plus grande dans l’échantillon. La maille de gros graphes orientés tend vers 2, celle de graphes non orientés tend vers 3. La variance de cette variable au sein de vos échantillons est donc souvent nulle, ce que l’ACP n’aime pas du tout du tout).

Data

##                         E    V Deg Density Diam Path.Len Clust Core
## droso158.gr           249  149  10  0.0226   12      5.3 0.073    4
## human1.gr             360  228  25  0.0139   15      5.1 0.062    3
## Pathway.gr         166761 8839 907  0.0043   13      3.7 0.350  110
## small_umbilical.gr    212   93  17  0.0496    8      3.8 0.168    4
## Complexes.gr        36762 2528 275  0.0115   14      4.3 0.788  140
## droso205.gr           343  196  14  0.0179   13      5.6 0.188    4
## JBA.gr                346  200  14  0.0174   16      6.5 0.191    4
## umbilical.gr         2002 1126  54  0.0032   13      4.9 0.046    4

Question 16: Suite de la question 16 du TP précédent. Récupérez les deux scripts PPI_network.Rmd et co-expression_network.Rmd. Placez-les dans le même répertoire que CORUM_complexes.Rmd et Utils.R. Exécutez PPI_network.Rmd. Vous pouvez également tester co-expression_network.Rmd mais attention, il vous faut 2GB de mémoire vive disponible car le réseau produit est très grand. Expliquez, de la même manière que pour CORUM_complexes.Rmd la partie centrale de ces deux scripts. Insistez en particulier sur la méthodologie permettant de construire un graphe PPI à partir de données d’expression.

Co-expression

Question 17: Cette question est longue en temps d’exécution. Corrélation de Spearman : combien d’arêtes obtient-on si on passe le seuil à 0.75, 0.8, 0.85, 0.90, 0.95 ? Donnez le code R. Générez éventuellement un graphique.

PPI

Question 18: Parmi les bases de données listées dans la variable DB, lesquelles contiennent la plus grande quantité de données PPI ? Donner le code R produisant des graphes distincts et leur taille. Générez éventuellement un graphique.

Corum, PPI, Co-expression

Question 19: Les complexes protéiques sont des protéines interagissant ensemble. On peut donc à juste titre supposer qu’un grand nombre d’arêtes sont communes entre le graphe PPI et le graphe CORUM. Est-ce vraiment le cas. Calculez, et éventuellement plottez le nombre (ou la proportion) d’arêtes communes et distinctes.

Question 20: Dans le même ordre d’idées, si la vocation d’un ensemble de protéines est d’interagir ensemble, ce serait une bonne idée qu’elles soient exprimées au même moment. Donc, même question que ci-dessus, avec les trois réseaux. Les diagrammes de Venn permettent une représentation assez sympathique de ce type de résultat. Les packages VennDiagram et gplots contiennent tous les deux des fonctions permettant de générer de tels diagrammes.

Produisez un diagramme de Venn. Les résultats vous paraissent-ils en accord avec les hypothèses avancées ?