PS_ROUQUIER_DM_3
library(dplyr)## Warning: package 'dplyr' was built under R version 3.3.3##
## Attaching package: 'dplyr'## The following objects are masked from 'package:stats':
##
## filter, lag## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, unionlibrary(ggplot2)## Warning: package 'ggplot2' was built under R version 3.3.3library(reshape2)## Warning: package 'reshape2' was built under R version 3.3.3library(Rmisc) ## Warning: package 'Rmisc' was built under R version 3.3.3## Loading required package: lattice## Loading required package: plyr## Warning: package 'plyr' was built under R version 3.3.3## -------------------------------------------------------------------------## You have loaded plyr after dplyr - this is likely to cause problems.
## If you need functions from both plyr and dplyr, please load plyr first, then dplyr:
## library(plyr); library(dplyr)## -------------------------------------------------------------------------##
## Attaching package: 'plyr'## The following objects are masked from 'package:dplyr':
##
## arrange, count, desc, failwith, id, mutate, rename, summarise,
## summarizeset.seed(4200)Code pour Marche Aléatoire:
rw_line<-function(N,ITER,disp){
t1 = Sys.time()
A = matrix(data = rep(0,N*N),nrow=N)
for(i in 2:N-1){
A[i,i-1] = 1/2
A[i,i+1] = 1/2
}
A[1,2] = 1
A[N,N-1] = 1
B = c(rep(0,N))
current = sample(1:N,1)
for(k in 1:ITER){
last = current
if(max(A[current,])!= 0){ #dangling node? Si oui, on ne tire pas un nombre car proba = 0
current = sample(c(1:N), 1, prob=A[current,],replace = TRUE)
} else {
current = sample(c(1:N),1)
}
if(runif(1)>0.85){ #Permet de se sortir des buckets. 1/6 du temp on se TP à un rrandom quoi qu'il arrive
current = sample(c(1:N),1)
}
B[current] = B[current] + 1
}
C = B/ITER
t2 = difftime(Sys.time(),t1)
if(disp){
print(t2)
barplot(C,ylim=c(0.0,max(C)*1.2),names.arg = c(1:N))
}
return(t2)
}
rw_ring<-function(N,ITER,disp){
t1 = Sys.time()
A = matrix(data = rep(0,N*N),nrow=N)
for(i in 2:N-1){
A[i,i-1] = 1/2
A[i,i+1] = 1/2
}
A[1,N] = 1/2
A[1,2] = 1/2
A[N,1] = 1/2
A[N,N-1] = 1/2
B = c(rep(0,N))
current = sample(1:N,1)
for(k in 1:ITER){
last = current
if(max(A[current,])!= 0){ #dangling node? Si oui, on ne tire pas un nombre car proba = 0
current = sample(c(1:N), 1, prob=A[current,],replace = TRUE)
} else {
current = sample(c(1:N),1)
}
if(runif(1)>0.85){ #Permet de se sortir des buckets. 1/6 du temp on se TP à un random quoi qu'il arrive
current = sample(c(1:N),1)
}
B[current] = B[current] + 1
}
C = B/ITER
t2 = difftime(Sys.time(),t1)
if(disp){
print(t2)
barplot(C,ylim=c(0.0,max(C)*1.2),names.arg = c(1:N))
}
return(t2)
}
rw_chupa<-function(N_head,N_tail,ITER,disp){
t1 = Sys.time()
A = matrix(data = rep(0,(N_head+N_tail)*(N_head+N_tail)),nrow=(N_head+N_tail))
for(i in 2:N_head-1){
A[i,i-1] = 1/2
A[i,i+1] = 1/2
}
A[1,N_head] = 1/3
A[1,2] = 1/3
A[1,N_head+1] = 1/3
A[N_head,1] = 1/2
A[N_head,N_head-1] = 1/2
if(N_tail>1){
for(r in (N_head+2):(N_head+N_tail-1)){
A[r,r-1] = 1/2
A[r,r+1] = 1/2
}
A[N_head+1,1] = 1/2
A[N_head+1,N_head+2] = 1/2
A[N_head+N_tail,N_head+N_tail-1] = 1
}else{
A[N_head+1,1] = 1
}
B = c(rep(0,N_head+N_tail))
current = sample(1:(N_head+N_tail),1)
for(k in 1:ITER){
last = current
if(max(A[current,])!= 0){ #dangling node? Si oui, on ne tire pas un nombre car proba = 0
current = sample(c(1:(N_head+N_tail)), 1, prob=A[current,],replace = TRUE)
} else {
current = sample(c(1:(N_head+N_tail)),1)
}
if(runif(1)>0.85){ #Permet de se sortir des buckets. 1/6 du temp on se TP à un rrandom quoi qu'il arrive
current = sample(c(1:(N_head+N_tail)),1)
}
B[current] = B[current] + 1
}
C = B/ITER
t2 = difftime(Sys.time(),t1)
if(disp){
print(t2)
barplot(C,ylim=c(0.0,max(C)*1.2),names.arg = c(1:(N_head+N_tail)))
}
return(t2)
}Code pour itération de la fonction de transition:
ft_line<-function(N,ITER,disp){
t1 = Sys.time()
A = matrix(data = rep(0,N*N),nrow=N)
for(i in 2:N-1){
A[i,i-1] = 1/2
A[i,i+1] = 1/2
}
A[1,2] = 1
A[N,N-1] = 1
B = c(rep(1/N,N))
for(k in 1:ITER){
B = B%*%A
}
t2 = difftime(Sys.time(),t1)
if(disp){
print(t2)
barplot(B,ylim=c(0.0,max(B)*1.2),names.arg = c(1:N))
}
return(t2)
}
ft_ring<-function(N,ITER,disp){
t1 = Sys.time()
A = matrix(data = rep(0,N*N),nrow=N)
for(i in 2:N-1){
A[i,i-1] = 1/2
A[i,i+1] = 1/2
}
A[1,N] = 1/2
A[1,2] = 1/2
A[N,1] = 1/2
A[N,N-1] = 1/2
B = c(rep(1/N,N))
for(k in 1:ITER){
B = B%*%A
}
t2 = difftime(Sys.time(),t1)
if(disp){
print(t2)
barplot(B,ylim=c(0.0,max(B)*1.2),names.arg = c(1:N))
}
return(t2)
}
ft_chupa<-function(N_head,N_tail,ITER,disp){
t1 = Sys.time()
A = matrix(data = rep(0,(N_head+N_tail)*(N_head+N_tail)),nrow=(N_head+N_tail))
for(i in 2:N_head-1){
A[i,i-1] = 1/2
A[i,i+1] = 1/2
}
A[1,N_head] = 1/3
A[1,2] = 1/3
A[1,N_head+1] = 1/3
A[N_head,1] = 1/2
A[N_head,N_head-1] = 1/2
if(N_tail>1){
for(r in (N_head+2):(N_head+N_tail-1)){
A[r,r-1] = 1/2
A[r,r+1] = 1/2
}
A[N_head+1,1] = 1/2
A[N_head+1,N_head+2] = 1/2
A[N_head+N_tail,N_head+N_tail-1] = 1
}else{
A[N_head+1,1] = 1
}
B = c(rep(1/(N_head+N_tail),N_head+N_tail))
for(k in 1:ITER){
B = B%*%A
}
t2 = difftime(Sys.time(),t1)
if(disp){
print(t2)
barplot(B,ylim=c(0.0,max(B)*1.2),names.arg = c(1:(N_head+N_tail)))
}
return(t2)
}Code pour calcul théorique avec valeurs propres:
vp_line<-function(N,ITER,disp){
t1 = Sys.time()
A = matrix(data = rep(0,N*N),nrow=N)
for(i in 2:N-1){
A[i,i-1] = 1/2
A[i,i+1] = 1/2
}
A[1,2] = 1
A[N,N-1] = 1
B = DTMCPack::statdistr(A)
t2 = difftime(Sys.time(),t1)
if(disp){
print(t2)
barplot(B,width = 0.33,ylim=c(0.0,max(B)*1.2),names.arg = c(1:N))
}
return(t2)
}
vp_ring<-function(N,ITER,disp){
t1 = Sys.time()
A = matrix(data = rep(0,N*N),nrow=N)
for(i in 2:N-1){
A[i,i-1] = 1/2
A[i,i+1] = 1/2
}
A[1,N] = 1/2
A[1,2] = 1/2
A[N,1] = 1/2
A[N,N-1] = 1/2
B = c(rep(1/N,N))
B = DTMCPack::statdistr(A)
t2 = difftime(Sys.time(),t1)
if(disp){
print(t2)
barplot(B,width = 0.33,ylim=c(0.0,max(B)*1.2),names.arg = c(1:N))
}
return(t2)
}
vp_chupa<-function(N_head,N_tail,ITER,disp){
t1 = Sys.time()
A = matrix(data = rep(0,(N_head+N_tail)*(N_head+N_tail)),nrow=(N_head+N_tail))
for(i in 2:N_head-1){
A[i,i-1] = 1/2
A[i,i+1] = 1/2
}
A[1,N_head] = 1/3
A[1,2] = 1/3
A[1,N_head+1] = 1/3
A[N_head,1] = 1/2
A[N_head,N_head-1] = 1/2
if(N_tail>1){
for(r in (N_head+2):(N_head+N_tail-1)){
A[r,r-1] = 1/2
A[r,r+1] = 1/2
}
A[N_head+1,1] = 1/2
A[N_head+1,N_head+2] = 1/2
A[N_head+N_tail,N_head+N_tail-1] = 1
}else{
A[N_head+1,1] = 1
}
B = DTMCPack::statdistr(A)
t2 = difftime(Sys.time(),t1)
if(disp){
print(t2)
barplot(B,width = 1,ylim=c(0.0,max(B)*1.2),names.arg = c(1:(N_head+N_tail)))
}
return(t2)
}Fonction de test:
test_line <- function(N,ITER){
ITER2 = 100
ratio = ITER/ITER2
ITER = ITER2
temp = c(1:ITER)
temp = temp*ratio
line= rep(0,ITER)
for(m in 1:ITER2){
ITER = m*ratio
line[m] = rw_line(N,ITER,0)
}
Implementation=c(rep("line_rw",ITER))
rw_df_line = data.frame(temp,Implementation,line)
for(m in 1:ITER2){
ITER = m*ratio
line[m] = ft_line(N,ITER,0)
}
Implementation=c(rep("line_ft",ITER))
ft_df_line = data.frame(temp,Implementation,line)
Implementation=c(rep("line_vp",ITER))
vp_temp = vp_line(N,ITER,0)
line = rep(vp_temp,ITER)
vp_df_line = data.frame(temp,Implementation,line)
line_df=rbind.data.frame(rw_df_line,ft_df_line,vp_df_line)
ggplot(line_df, aes(x=temp, y=line, color=Implementation)) + geom_point(alpha=.5) + geom_smooth(method="lm",formula=y~x*I(log(x)))
}
test_ring <- function(N,ITER){
ITER2 = 100
ratio = ITER/ITER2
ITER=ITER2
temp = c(1:ITER)*ratio
ring= rep(0,ITER)
for(m in 1:ITER2){
ITER = m*ratio
ring[m] = rw_ring(N,ITER,0)
}
Implementation=c(rep("ring_rw",ITER))
rw_df_ring = data.frame(temp,Implementation,ring)
for(m in 1:ITER2){
ITER = m*ratio
ring[m] = ft_ring(N,ITER,0)
}
Implementation=c(rep("ring_ft",ITER))
ft_df_ring = data.frame(temp,Implementation,ring)
Implementation=c(rep("ring_vp",ITER))
vp_temp = vp_ring(N,ITER,0)
ring = rep(vp_temp,ITER)
vp_df_ring = data.frame(temp,Implementation,ring)
ring_df=rbind.data.frame(rw_df_ring,ft_df_ring,vp_df_ring)
ggplot(ring_df, aes(x=temp, y=ring, color=Implementation)) + geom_point(alpha=.5) + geom_smooth(method="lm",formula=y~x*I(log(x)))
}
test_chupa <- function(N_head,N_tail,ITER){
ITER2 = 100
ratio = ITER/ITER2
ITER=ITER2
temp = c(1:ITER)*ratio
chupa= rep(0,ITER)
for(m in 1:ITER2){
ITER = m*ratio
chupa[m] = rw_chupa(N_head,N_tail,ITER,0)
}
Implementation=c(rep("chupa_rw",ITER))
rw_df_chupa = data.frame(temp,Implementation,chupa)
for(m in 1:ITER2){
ITER = m*ratio
chupa[m] = ft_chupa(N_head,N_tail,ITER,0)
}
Implementation=c(rep("chupa_ft",ITER))
ft_df_chupa = data.frame(temp,Implementation,chupa)
Implementation=c(rep("chupa_vp",ITER))
vp_temp = vp_chupa(N_head,N_tail,ITER,0)
chupa = rep(vp_temp,ITER)
vp_df_chupa = data.frame(temp,Implementation,chupa)
chupa_df=rbind.data.frame(rw_df_chupa,ft_df_chupa,vp_df_chupa)
ggplot(chupa_df, aes(x=temp, y=chupa, color=Implementation)) + geom_point(alpha=.5) + geom_smooth(method="lm",formula=y~x*I(log(x)))
}Tests:
Avant toute chose, j’aimerais définir la série de tests suivante
comme incomplète. Elle permet bien de se donner une idée du temps
nécessaire pour chaque valeur de paramètres, mais il manque une
information importante: Elle ne renseigne pas sur la précision des
calculs
Les calculs effectués par vp sont parfaits. Il faudrait
effectuer un plot avec les 3 dataframe. Les valeurs de vp agiraient
comme référence, et il faudrait regarder la différence de chacune des
deux autres courbes avec ces valeurs. J’aurai du le faire, voila.
De plus, les tests temporels sont réalisés en regardant la différence
entre deux appels de sys.time(). Ils y a donc de gros faCteurs modulant
les tests, et les rendants décorellés les uns des autres. Il aurait été
mieux de séparer totalement les tests, en rénitialisant l’environnement à
chaque fois. De plus, le calcul des vp est toujours effectué une seule
fois dans les tests, et est donc trop exposé. Il aurait été bien mieux
de faire une moyenne
Chaque test en réalise 100, la variable
ITER partant de 1 (les résultats à ITER très petits sont d’ailleurs
potentiellement trop incorrects pour être pris en compte) jusqu’a l’ITER
donné en paramètre, avec un pas de ITER/100 à chaque fois. Cela permet
de se donner une bonne approximation de l’évolution des courbes sans
avoir d’explosions calculatoires à ITER grand
En ordonnée il y a le temps en secondes, et en abcsisse la valeur d’ITER
Evidemment, plus les points sont bas, mieux c’est. Egalement, seul le nombre de noeuds impacte les calculs de vp.
Test 1: N très petit (10), ITER moyen(1000)
test_line(10,1000)
test_ring(10,1000)
test_chupa(5,5,1000)
Test 1bis: N très petit (10), ITER très grand(100.000)
test_line(10,100000)
test_ring(10,100000)
test_chupa(5,5,100000)
On voit de suite que lorsque N reste petit, le calcul par vp est
toujour plus efficace. Vu les capacités de ma machine, il ne m’est pas
possible de tester jusqu’a quelle valeur de N cela resterai plus
efficace. Mais pour dess graphes comme celui d’internet, ce calcul sur
un graphe de ce type est impossible à calculer. En effet on passe par
l’inverse d’une matrice, et sur des graphes aussi gros, cela crée une
explosion calculatoire enorme.
la opremière valeur de vp_line est étrange d’ailleurs, sûrement un soucis lors du calcul.
On voit d’ailleurs que la forme du graphe ne change pas les temps de calculs, ce qui est normal.
lorsque ITER est grand par rapport au nombre de noeuds, on voit que rw est toujours plus efficace que ft
Test 2: N moyen (100), ITER moyen (1.000)
test_line(100,1000)
test_ring(100,1000)
test_chupa(50,50,1000)
Ce test met en excergue que lrosque N devient grand, il est devient plus efficace de faire rw de ft
Test 3: N grand (1.000), ITER moyen(1.000)
test_line(1000,1000)
test_ring(1000,1000)
test_chupa(500,500,1000)
celui ci permet une mise en évidence que quand N deviendra très grand,
les autres tests seront plus efficaces. Attention tout de même à ne pas
se faire feinter! On voit ici qu’il est ici plus efficace de faire rw et
ft quand ITER est plius petit. OR, quand ITER est trop petit par
rapport à la valeur de N, les résultats ne sont pas assez précis!
Réaliser à chaque fois les tests dont j’ai parlé préalablement auraient
permi de bien le voir.
Par intuition, pour rw, les valeurs qui
tendent le plus rapidement à être correctes doivent être celles de
l’anneau, étant donné qu’il n’y a aucune particularité.
Pour la ligne, il faut être sûr que le jeton arrive de façon significative sur les bords
Même choses pour la chupa qui doit être la forme la plus difficile à
estimer pour rw, il y a des possibilité de rester bloquer dans une ou
l’autre des parties de la chupa.
Pour ft, la forme du graphe n’influe pas sur la valeur d’ITER minimum à avoir pour avoir des résultat satifaisemment précis.
Pour vp, les valeurs sont toujours exactes.
Test 4: Chupa uniquement: Variation du rapport size_head/size_tail (10/90, 50/50, 90/10), N très petit, ITER moyen
test_chupa(1,9,1000)
test_chupa(5,5,1000)
test_chupa(9,1,1000)
Test 4b: Chupa uniquement: Variation du rapport size_head/size_tail (10/90, 50/50, 90/10), N très petit, ITER trè grand
test_chupa(1,9,100000)
test_chupa(5,5,100000)
test_chupa(9,1,100000)
Test 4c: Chupa uniquement: Variation du rapport size_head/size_tail (10/90, 50/50, 90/10), N moyen, ITER moyen
test_chupa(10,90,1000)
test_chupa(50,50,1000)
test_chupa(90,10,1000)
Test 4d: Chupa uniquement: Variation du rapport size_head/size_tail (10/90, 50/50, 90/10), N grand, ITER moyen
test_chupa(100,900,1000)
test_chupa(500,500,1000)
test_chupa(900,100,1000)
Les changements de valeurs de chupa n’influencent pas les résultats temporels.