DM_2_FUSTES

set.seed(77)
N = 100 # Nombre de sommets dans les graphes générés
Nb = 500; # Nombre de graphes pour la simulation

Question1

Fonction de génération d’un graphe d’Erdos-renyi (par sa matrice d’adjacence et comptage du nombre de sommets isolés) Note: On considerera que les graphes d’Erdos-renyi n’ont pas de boucle (Wikipedia) donc Gij (i=j) =0

#Generer le graphe
generer_graphe = function(P){
  matrice = matrix(data = 2, nrow = N, ncol = N)
  for(ligne in 1:N){
    for(colonne in ligne:N){
      if(ligne == colonne){
        matrice[ligne,colonne] = 0; #Pas de boucle
      }else{
        matrice[ligne,colonne] = rbinom(1,1,P);       #Bernoulli
        matrice[colonne,ligne] = matrice[ligne,colonne];  #Faire la symetrie
      }
    }
  }
  return(matrice)
}

#Compter les sommets isolés
nb_sommets_isoles = function(graphe){
  nb_isoles = 0
  for(ligne in 1:N){
    compteur = 0;
    for(colonne in 1:N){
        compteur = compteur + graphe[ligne,colonne]
    }
    if(compteur == 0){
      nb_isoles = nb_isoles + 1
    }
    else {}
  }
  return(nb_isoles)
}

Fonction qui crée Nb graphes aléatoire (probabilité p) et qui renvoie un tableau de Nb nombre de sommets isolés.

simulation = function(Nb, p) {
  tab = 1:Nb
  for (i in tab) {
    tab[i] = nb_sommets_isoles(generer_graphe(p))
  }
  return(tab)
}

Fonction qui regarde si un graphe est connexe

est_connexe = function(g) {
  if (nb_sommets_isoles(g) == 0) 
  { 
    return (TRUE)
  }
  else { return (FALSE);}
}

Fonction qui renvoie la proportion de graphe connexe dans un tirage de Nb graphe aléatoire

simulation_connexe = function(Nb, p) {
  tab = 1:Nb
  res = 0
  for (i in tab) {
    if (est_connexe(generer_graphe(p)))
    {
      res = res +1
    }
  }
  return(res/Nb)
}

c = -4;
t = simulation(Nb, ((c+log(N))/N));
hist(t);
abline(v = exp(-c), col=c("red"))

mean(t)

## [1] 54.686

exp(-c)

## [1] 54.59815

c = -2;
t = simulation(Nb, ((c+log(N))/N));
hist(t);
abline(v = exp(-c), col=c("red"))

mean(t)

## [1] 7.002

exp(-c)

## [1] 7.389056

c = 3;
t = simulation(Nb, ((c+log(N))/N));
hist(t);
abline(v = exp(-c), col=c("red"))

mean(t)

## [1] 0.034

exp(-c)

## [1] 0.04978707

On voit que \[e^{-c}\] est une bonne estimation. On retrouve en rouge sur les histogramme la valeur de \[e^{-c}\]. On remarque que le nombre de sommets isolés est centrée autour de cette valeur.

Pour N = 25

N=25
valeurs=seq(-2,4,0.5)
res=c()
obj=c()
for(i in 1:length(valeurs))
{
  c = valeurs[i]
  res[i] = mean(simulation(100, ((c+log(N))/N)))
  obj[i]=exp(-c)
}

  plot(obj,x = valeurs, type = "l", col = "red",xlab = "c",ylab = "exp(-c)", main = "Théorie")

  plot(res,x = valeurs, type = "l", col = "red",xlab = "c",ylab = "probabilté d'isolation", main = "Expérimentation")

Et pour N = 100

N=100
valeurs=seq(-3,4,0.5)
res=c()
obj=c()
for(i in 1:length(valeurs))
{
  c = valeurs[i]
  res[i] = mean(simulation(200, ((c+log(N))/N)))
  obj[i]=exp(-c)
}

  plot(obj,x = valeurs, type = "l", col = "red",xlab = "c",ylab = "exp(-c)", main = "Théorie")

  plot(res,x = valeurs, type = "l", col = "red",xlab = "c",ylab = "probabilté d'isolation", main = "Expérimentation")

On peut voir que les courbes sont tres similaire

Question2

Probabilité qu’un sommet soit isolé (n-1 car on ne compte pas les boucles) : \[P(X)=(1−p)^{n-1} =(1− \frac{c+log(n)}{n})^{n-1}\]

Calcul de l’espérance (n*p) : \[E(X)=n(1− \frac{c+log(n)}{n})^{n-1}\]

\[=(1+ \frac{−c−log(n)}{n})^{n}\] \[=e^{−c−log(n)}\] \[=e^{−c}e^{−log(\frac{1}{n})}\] \[e^{−c}∗\frac{1}{n}\] \[= \frac{e^{−c}}{n}\]

\[avec (1+\frac{a}{n]})^{n}=e^{a}\]

Sur cette base de calcul j’aimerais essayer de trouver e^c mais je trouve (e^c)/n. Je n’arrive pas à aller plus loin dans la résolution pour corriger ce petit probleme.

Question3

On va comparer l’histogramme d’une loi de poisson de parametre e^-c avec celui de notre graphe aléatoire. On prendra N sommets sur des grandes générations de 2000 essais.

c = -2
N= 100

intervalle = seq(0, 20, 2);

hist1 <- hist(simulation(2000, ((c+log(N))/N)), breaks = intervalle, plot = FALSE);
hist2 <- hist(rpois(2000, exp(-c)), breaks = intervalle, plot = FALSE);
m = rbind(hist1$counts, hist2$counts); #Aide sur internet
barplot(m, beside = TRUE, las = 2, legend = c("Graphe aléatoire", "Poisson")) #Aide sur internet

On remarque ici que les 2 histogrammes ont la meme tendance.

Question4

c = 0.5
N = 100

exp(-exp(-c));

## [1] 0.5452392

simulation_connexe(1000,((c+log(N))/N))

## [1] 0.575

c = 2
N = 100

exp(-exp(-c));

## [1] 0.873423

simulation_connexe(1000,((c+log(N))/N))

## [1] 0.902

On voit que les résultats sont tres proche, l’égalité est donc verifiée. Les résultats peuvent etre encore plus précis avec des variable (arguments) tres grande.