R Markdown

Pour ce devoir j’ai réflechi avec Rémy Palomo

1.1

set.seed(3)
a = 10
n = 15
k = 100
m = replicate(k, max(as.integer(runif(n,0,a+1))))
mean(m)
## [1] 9.69
a = 10
n = 10
k = 100
m = replicate(k, max(as.integer(runif(n,0,a+1))))
mean(m)
## [1] 9.39
a = 100
n = 10
k = 100
m = replicate(k, max(as.integer(runif(n,0,a+1))))
mean(m)
## [1] 91.8
a = 10
n = 10
k = 100
m = replicate(k, max(as.integer(runif(n,0,a+1))) + a/n)
hist(m)

1.2

set.seed(1)
a = 10
n = 15
k = 100
m = replicate(k, (2/n)*sum(as.integer(runif(n,0,a+1))))
mean(m)
## [1] 9.836
a = 10
n = 10
k = 100
m = replicate(k, (2/n)*sum(as.integer(runif(n,0,a+1))))
mean(m)
## [1] 9.852
a = 100
n = 1000
k = 100
m = replicate(k, (2/n)*sum(as.integer(runif(n,0,a+1))))
mean(m)
## [1] 99.92454

On observe que les M’ obtenus sont plus proches de A pour les mêmes tailles d’échantillons

Démonstration: \(M' = \frac{2}{n} * \sum_{k=0}^{a}X_i\)

or, \(E[\sum_{i=0}^{n}X_i]=\sum_{i=0}^{n}E[X_i]\) et \(E[X_i]= \frac{1}{a+1}*\sum_{k=0}^{a}k = \frac{1}{a+1}*\frac{a*(a+1)}{2} = \frac{a}{2}\)

Donc \(E[M'] = n*\frac{a}{2}*\frac{2}{n}=a\)

On calcule la variance : \(E[M']=a\) donc \(E[M']^2=a^2\) et \(E[M'^2]=\frac{2}{n}*\sum_{i=0}^{n}E[X_i^2]\)

Or, \(E[X_i^2]=\frac{1}{a+1}*\sum_{k=0}^{a}k^2\) \(=\frac{1}{a+1}*\frac{a(a+1)(2a+1)}{6}\) \(=\frac{a(2a+1)}{6}\)

Donc \(E[M'^2]=\frac{2}{n}*n*\frac{a (2a+1)}{6}\) \(=\frac{a(2a+1)}{3}\)

Enfin \(Var(M')= E[M'^2] - E[M']^2\) \(=\frac{a(2a+1)}{3}-a^2\) \(=\frac{2a^2 +a}{3} - \frac{3a^2}{3}\) \(=\frac{a-a^2}{3}\) \(=\frac{a(1-a)}{3}\)

Cet estimateur semble donc meilleur qu’avec M.

2.3

set.seed(3)
res <- NULL
n=10
for (i in 1:1000){
  a= runif(1,0,1)
  M= max(runif(n,0,a))
  r = 1.1 *M
  if (r < a){
    res[i] <- r-M
  } else {
    res[i] <- 0
  }
}
mean(res)
## [1] 0.01631215

C’est une stratégie peu payante puisque l’ésperance du gain est assez faible

On le montre grâce à : \(E[M]=A - \frac{A}{n}= A - \frac{A}{10} = \frac{9}{10}A\) et \(E[r(M)]= 1.1 * E[M] = 0.99*A\) donc \(E[G]= E[ r(M)-M ] = E[ r(M)] - E[M] = 0.99A - 0.9A = 0.09A\)

2.4

Ce code permet d’estimer \(P[A=a,M=m]\) pour les valeurs voulues de a et m :

res <- NULL
n=3
for (i in 1:100){
  A= as.integer(runif(1,0,10+1))/10
  M= max(replicate(n,as.integer(runif(1,0,A*10+1))/10))
  a= 0.6
  m= 0.5
  if (A==a && M==m){
    res[i] <- 1
  } else {
    res[i] <- 0
  }
}
mean(res)
## [1] 0.04

Avec ce code, on trace une estimation de la probabilité de gagner en fonction de la valeur de A lorsque M vaut 0.5 :

set.seed(10)
resf <- NULL
for (j in 5:10){
  res <- NULL
  n=3
  for (i in 1:100){
    A= as.integer(runif(1,0,10+1))/10
    M= max(replicate(n,as.integer(runif(1,0,A*10+1))/10))
    a= 0.1*j
    m= 0.5
    if (A==a && M==m){
      res[i] <- 1
    } else {
      res[i] <- 0
    }
  }
  resf[j] <- mean(res)
}
plot(resf)

On observe que la meilleure valeur pour A dans ce cas est 0.5

2.5

C’est une meilleure stratégie que le facteur car elle permet d’approcher plus rapidement la valeur de A.

set.seed(1)
res <- NULL
  n= 2
  alpha= 0.7 
  for (i in 1:100){
    a= runif(1,0,1)
    M= max(runif(n,0,a))
    r = M^alpha
    if (r < a){
      res[i] <- r-M
    } else {
      res[i] <- 0
    }
  }
"pour alpha = 0.7"
## [1] "pour alpha = 0.7"
mean(res)
## [1] 0.05974903
res <- NULL
  n= 2
  alpha= 0.5
  for (i in 1:100){
    a= runif(1,0,1)
    M= max(runif(n,0,a))
    r = M^alpha
    if (r < a){
      res[i] <- r-M
    } else {
      res[i] <- 0
    }
  }
"pour alpha = 0.5"
## [1] "pour alpha = 0.5"
mean(res)
## [1] 0.06802448
res <- NULL
  n= 2
  alpha= 0.3
  for (i in 1:100){
    a= runif(1,0,1)
    M= max(runif(n,0,a))
    r = M^alpha
    if (r < a){
      res[i] <- r-M
    } else {
      res[i] <- 0
    }
  }
"pour alpha = 0.3"
## [1] "pour alpha = 0.3"
mean(res)
## [1] 0.05468087

La meilleure valeur pour alpha parmi ces 3 nombres est 0.5 Pour essayer de trouver mieux, on peut représenter une estimation de l’ésperance du gain G en fonction d’alpha

set.seed(1)
resf<- NULL
for (j in 1:10){
  res <- NULL
  n= 2
  alpha= j/10 
  for (i in 1:100){
    a= runif(1,0,1)
    M= max(runif(n,0,a))
    r = M^alpha
    if (r < a){
      res[i] <- r-M
    } else {
      res[i] <- 0
    }
  }
  resf[j]<-mean(res)
}
plot(resf)

On observe que le gain varie beaucoup mais que le pic se trouve souvent entre 0.4 et 0.6

set.seed(1)
res <- NULL
  n= 10
  alpha= 0.7 
  for (i in 1:100){
    a= runif(1,0,1)
    M= max(runif(n,0,a))
    r = M^alpha
    if (r < a){
      res[i] <- r-M
    } else {
      res[i] <- 0
    }
  }
"pour alpha = 0.7"
## [1] "pour alpha = 0.7"
mean(res)
## [1] 0.005931772
res <- NULL
  n= 10
  alpha= 0.5
  for (i in 1:100){
    a= runif(1,0,1)
    M= max(runif(n,0,a))
    r = M^alpha
    if (r < a){
      res[i] <- r-M
    } else {
      res[i] <- 0
    }
  }
"pour alpha = 0.5"
## [1] "pour alpha = 0.5"
mean(res)
## [1] 0.01075223
res <- NULL
  n= 10
  alpha= 0.3
  for (i in 1:100){
    a= runif(1,0,1)
    M= max(runif(n,0,a))
    r = M^alpha
    if (r < a){
      res[i] <- r-M
    } else {
      res[i] <- 0
    }
  }
"pour alpha = 0.3"
## [1] "pour alpha = 0.3"
mean(res)
## [1] 0.005848023

Pour n=10, c’est pour alpha=0.7 que le gain est le plus grand

set.seed(8)
resf<- NULL
for (j in 1:10){
  res <- NULL
  n= 10
  alpha= j/10 
  for (i in 1:100){
    a= runif(1,0,1)
    M= max(runif(n,0,a))
    r = M^alpha
    if (r < a){
      res[i] <- r-M
    } else {
      res[i] <- 0
    }
  }
  resf[j]<-mean(res)
}
plot(resf)

Pour n=10, le pic se trouve entre 0.7 et 0.8, donc la meilleure valeur pour alpha se trouve aux alentours de 0.75

En définitif, la méthode du facteur est plus efficace lorsque le nombre n de tirages est grand tandis que celle de la puissance a une meilleure efficacité lorsque n est petit (si on choisit bien alpha).

set.seed(4)
resf<- NULL
resf2<- NULL
for (j in 1:10){
  res <- NULL
  res2 <- NULL
  n= j
  alpha= 0.6 
  for (i in 1:100){
    a= runif(1,0,1)
    M= max(runif(n,0,a))
    r = M^alpha
    r2 = 1.1*M
    if (r < a){
      res[i] <- r-M
    } else {
      res[i] <- 0
    }
    if (r2 < a){
      res2[i] <- r2-M
    } else {
      res2[i] <- 0
    }
  }
  resf[j]<-mean(res)
  resf2[j]<-mean(res2)
}
i=1:10
plot(resf,type ="p",col="blue",xlab ="nb de tirages",ylab = "E(G)")
par(new = TRUE)
plot(resf2,type ="p",col="red",xaxt="n",yaxt="n",xlab ="",ylab = "")
legend("topright", c("^alpha", "* 1.1"),
       col = c("blue", "red"), pch = c(1, 1))

```