#**SIMULATION**LES UNIVERS D'UN DÉ ET LEURS TABLEAUX DE FRÉQUENCES**ACTE1******
deacte1 <- function(de = c("aUN", "bDEUX", "cTROIS", "dQUATRE", "eCINQ", "fSIX"),
                    nbsim = 2000){
  experience <- sample(de, nbsim, replace = T)
  distfreq <- table(experience)
  print(distfreq)
  barplot(distfreq / nbsim)
}

#**SIMULATION**LES UNIVERS D'UN DÉ ET LEURS TABLEAUX DE FRÉQUENCES**ACTE2*****
deacte2 <- function(de = c("aUN", "bDEUX", "cTROIS", "dQUATRE", "eCINQ", "fSIX"),
                    proba = rep(1 / 6, 6), nbsim = 2000){
  experience <- sample(de, nbsim, replace = T, prob = proba)
  distfreq <- table(experience)
  print(distfreq)
  barplot(distfreq / nbsim)
}

#**SIMULATION**PROBLÈME HISTORIQUE D'UN GRAND DUC DE TOSCANE : 3 DÉS***ALGO 1****
toscane1 <- function(nbsim = 2000){
  de1 <- sample(seq(1, 6, 1), nbsim, replace = T)
  de2 <- sample(seq(1, 6, 1), nbsim, replace = T)
  de3 <- sample(seq(1, 6, 1), nbsim, replace = T)
  som3des <- de1 + de2 + de3
  distsom <- table(som3des)
  print(distsom)
  barplot(distsom / nbsim)
}

#**SIMULATION**PROBLÈME HISTORIQUE D'UN GRAND DUC DE TOSCANE : 3 DÉS***ALGO 2****
toscane2 <- function(nbsim = 2000){
  som3des <- rep(NA, nbsim)
  for(i in 1:nbsim){
     de1 <- sample(seq(1, 6, 1), 1)
     de2 <- sample(seq(1, 6, 1), 1)
     de3 <- sample(seq(1, 6, 1), 1)
     som3des[i] <- de1 + de2 + de3
                   }
  distsom <- table(som3des)
  print(distsom)
  barplot(distsom / nbsim)
}

#**SIMULATION**DU PROBLÈME DES "CHAINES" DE LONGUEUR 6*******
chainesdSimpl = function(L = 6, tirages = 200, nbsim = 2000){
  nblongMax <- rep(NA, nbsim)
  for(i in 1:nbsim){
      expe <- sample(c("Pile", "Face"), tirages, replace=T)
      rexpe <- rle(expe)
      nblongMax[i] <- max(rexpe$length)
                   }
  tablemax <- table(nblongMax)
  StatL <- sum(tablemax[(as.numeric(names(tablemax)) >= L)]) / nbsim
#------------- Affichage des résultats et des graphiques-----------
  cat("Une estimation de la probabilité de chaines de longueur", L, "ou plus")
  cat(" vaut :", StatL, "\n \n")
  print("Distribution* de la chaine de longueur maximale (runs)")
  print(tablemax)
  barplot(tablemax,
  xlab = paste("Longueur maximale de chaine dans les ", tirages, "tirages"),
  ylab = "Effectif")
}

#**SIMULATION**INTERVALLE DE FLUCTUATION SIMULÉ D'UNE PROPORTION*******
simulif = function(n = 100, p =.52, nsim = 1000){
vecbino <- rbinom(nsim, n, p)
# effec <- table(vecbino)
quantiles <- quantile(vecbino, probs = c(2.5,97.5)/100)
#------------- Affichage des résultats et des graphiques-----------
  print("Quantiles")
  print(quantiles)
  # barplot(effec, xlab = "Nombre de succès", ylab = "Effectif")
  hist(vecbino, xlab = "Nombre de succès", ylab = "Effectif")
  abline(v = quantiles, col = "red")
}

#*******EXPLORATION D'UN INTERVALLE DE FLUCTUATION ASYMPTOTIQUE******
pIFasy2_1 <- function(n = 200, p= .5, proba = .95){
 P <- function(n, p, proba){
      g <- qnorm(1 - (1 - proba) / 2)
      binf <- max(floor(n * p - g * sqrt(p * (1 - p) * n)), 0)
      bsup <- min(floor(n * p + g * sqrt(p * (1 - p) * n)), n)
      sum(dbinom((binf + 1):bsup, n, p))
                           }
 y <- vector(length = n)
 for(i in 1:n) y[i] <- P(i, p, proba)
#------------- Affichage des résultats et des graphiques-----------
 plot(1:n, y, cex = .3, ylim = c(.9, 1), xlab = "Taille n de l'échantillon",
      ylab = "Probabilité binomiale d'être dans l'IF asymptotique",
      main = paste("Exploration de l'intervalle de fluctuation : influence de n",
             "avec p =",p))
 abline(h = proba)
}

#**************SIMULATION D'UN PEIGNE D'INTERVALLES DE CONFIANCE***********
simIC <- function(n = 50, p = .3, nbsim = 100, conf = .95, nbclass = 20){
 p <- round(runif(1, .1, .9), 2)
 x <- 1:nbsim
 vectpbino <- rbinom(nbsim, n, p) / n
 b1sup <- vectpbino + 1 / sqrt(n) ; b1inf <- vectpbino - 1 / sqrt(n)
 y1max <- 1.1 ; y1min <- - .1
 dehorsic1 <- sum(b1sup > p & b1inf > p) + sum(b1sup < p & b1inf < p)
#------------- Affichage des résultats et des graphiques----------
 cat("\np extérieur à l'IC1 :", dehorsic1, "fréquence =", dehorsic1 / nbsim, "\n")
 print(p)
#------------Peigne d'intervalles de confiance-------
 par(mfrow = c(2, 1))
plot(x, vectpbino, ylim = c(y1min, y1max), cex = .5, col = "green",
  xlab = paste(nbsim, "simulations"), ylab = "proportion de succès dans l'échantillon",
  main = paste("I. C. Asymptotiques au niveau 95%, d'une proportion dans la population\n",
               "Calculés à partir de", nbsim, "échantillons simulés, de taille", n))
 points(x, b1sup, cex = .5)
 points(x, b1inf, cex = .5)
 for(i in x) lines(c(i ,i), c(b1sup[i], b1inf[i]))
#------------Histogramme d'une série simulée---------
   hist(vectpbino, breaks = seq(0, 1, 1 / nbclass), freq = F, xaxp = c(0, 1, 20),
     xlab = "Classes de proportions", ylab = "Densité", cex.axis = .7, border = "green",
     main = paste("Histogramme et boite à moustache de", nbsim, "proportions simulées"))
#--------Boite à moustaches d'une série simulée------
   boxplot(vectpbino, freq = F, horizontal = T, add = T, xaxt = "n", col = "green")
}

#????????????????????????????
# ***FONCTION ***SIMULATIONS NUMÉRIQUES LOI DES GRANDS NOMBRES  lgn0*****
# Illustration graphique de la loi des grands nombres :
# Exemple du jet d'un dé (variable de bernoulli), occurence du CINQ
# Lorsque n augmente, on observe les suites de distributions
# Les fréquences tendent vers une valeur limite : la probabilité
# Les écarts à cette valeur limite sont de plus en plus faibles
lgn0 <- function(nbsim = 30, p = .3, n7 = c(10, 20, 50, 100, 500, 1000, 2000)){
  de <- seq(1, 6, 1)
  nechant <- rep(n7, each = nbsim)
  simfreqar <- rep(NA, nbsim * length(n7))
  for(
  simfreqar <- c(sample(de, n7[1], replace = T) / n7[1],
               rbinom(nbsim, n7[2], p) / n7[2], rbinom(nbsim, n7[3], p) / n7[3],
               rbinom(nbsim, n7[4], p) / n7[4], rbinom(nbsim, n7[5], p) / n7[5],
               rbinom(nbsim, n7[6], p) / n7[6], rbinom(nbsim, n7[7], p) / n7[7])
#--------------- Affichage des résultats et des graphiques------------
  graphics.off()
  plot(nechant, simfreqar, xaxp = c(0, max(n7), 10),
    main = "Distributions des fréquences des succès")
  dev.new()
  plot(as.factor(nechant), simfreqar,
    xlab = "Échantillons par tailles",
    ylab = "Fréquence des succès",
    main = "Résumé des distributions")
}
#????????????????????????????

# ***FONCTION ***SIMULATIONS NUMÉRIQUES LOI DES GRANDS NOMBRES **lgn1***
# Illustration graphique de la loi des grands nombres :
# Exemple de variables binomiales
# Lorsque n augmente, on observe les suites de distributions
# Les fréquences tendent vers une valeur limite : la probabilité
# Les écarts à cette valeur limite sont de plus en plus faibles
lgn <- function(nbsim = 20, p = .3, n7 = c(10, 20, 50, 100, 500, 1000, 2000)){
  nechant <- rep(n7, each = nbsim)
  simfreqar <- c(rbinom(nbsim, n7[1], p) / n7[1],
               rbinom(nbsim, n7[2], p) / n7[2], rbinom(nbsim, n7[3], p) / n7[3],
               rbinom(nbsim, n7[4], p) / n7[4], rbinom(nbsim, n7[5], p) / n7[5],
               rbinom(nbsim, n7[6], p) / n7[6], rbinom(nbsim, n7[7], p) / n7[7])
#--------------- Affichage des résultats et des graphiques------------
  graphics.off()
  plot(nechant, simfreqar, xaxp = c(0, max(n7), 10),
    main = "Distributions des fréquences des succès")
  dev.new()
  plot(as.factor(nechant), simfreqar,
    xlab = "Échantillons par tailles",
    ylab = "Fréquence des succès",
    main = "Résumé des distributions")
}

# ******FONCTION **SIMULATION** PILE FACE DISTRIBUTION DE PROBA *****
pileface <- function(nbsim = 2000){
 resultats <- vector(length = 3)
 names(resultats) <- c("DeuxPiles", "DeuxFaces", "Autre")
 for(i in 1:nbsim){
    pieceA <- sample(c("Pile", "Face"), 1)
    pieceB <- sample(c("Pile", "Face"), 1)
    if(pieceA == "Pile" & pieceB == "Pile") {
            resultats[1] <- resultats[1] + 1} else {
            if(pieceA == "Face" & pieceB == "Face")  {
                  resultats[2] <- resultats[2] + 1   } else {
                  resultats[3] <- resultats[3] + 1          }
                                                   }
                  }
#--------------- Affichage des résultats et des graphiques------------
 print(resultats)
 print(resultats / nbsim)
 barplot(resultats / nbsim)
}

# ******FONCTION **SIMULATION**CROIX PILE DE D'ALEMBERT*****
croixpile <- function(nbsim = 2000){
 resultats <- vector(length = 3)
 names(resultats) <- c("GagnéCoup1", "GagnéCoup2", "Perdu")
 for(i in 1:nbsim){
    coup1 <- sample(c("Pile", "Croix"), 1)
    if(coup1 == "Croix") {resultats[1] <- resultats[1] + 1} else {
       coup2 <- sample(c("Pile", "Croix"), 1)
       if(coup2 == "Croix") {resultats[2] <- resultats[2] + 1} else {
          resultats[3] <- resultats[3] + 1                          }
                                                                 }
                  }
#--------------- Affichage des résultats et des graphiques------------
 print(resultats)
 print(resultats / nbsim)
 barplot(resultats / nbsim)
}

# **FONCTION **SIMULATION**RANG DU PREMIER CROIX LORS DE 2 JETS***
premiercroix <- function(nbsim = 2000){
 resultats <- vector(length = 3)
 names(resultats) <- c("CroixRang1", "CroixRang2", "PasDeCroix")
 for(i in 1:nbsim){
    jet1 <- sample(c("Pile", "Croix"), 1)
    jet2 <- sample(c("Pile", "Croix"), 1)
    if(jet1 == "Croix") {resultats[1] <- resultats[1] + 1}
    if(jet1 != "Croix" & jet2 == "Croix") {resultats[2] <- resultats[2] + 1}
    if(jet1 != "Croix" & jet2 != "Croix") {resultats[3] <- resultats[3] + 1}
                  }
#--------------- Affichage des résultats et des graphiques------------
 print(resultats)
 print(resultats / nbsim)
 barplot(resultats / nbsim)
}

#**FONCTION **SIMULATION**TIRAGE SANS REMISE D'UN ÉCHANTILLON DE n BOULES ****
#***DANS UNE URNE CONTENANT R BOULES ROUGES ET Q BOULES D'AUTRES COULEURS ****
#***LA VARIABLE ÉTUDIÉE EST LE NOMBRE K DE BOULES ROUGES DANS L'ÉCHANTILLON ****
#*** OUVRIR DEUX FENÊTRES GRAPHIQUES AVEC dev.new() dev.new()
sansremise2 <- function(R = 6, Q = 14, n = 5, nbsim = 10000){
  M <- R + Q
  urne <- rep(c("rouge", "autre"), c(R, Q))
  vecnbrouges <- rep(NA, nbsim)
  for(i in 1:nbsim){
     echant <- sample(urne, n, replace = F)
     nbrouges <- sum(echant == "rouge")
     vecnbrouges[i] <- nbrouges
                   }
  distribnbrouges <- table(vecnbrouges) / nbsim
#--------------- Affichage des résultats et des graphiques------------
  cat("\nProportion de rouges dans l'urne :", R / M, "\n\n")
  print(distribnbrouges)
  dev.set(2)
  barplot(distribnbrouges, main = "TIRAGES SANS REMISE")
}

#**FONCTION **SIMULATION**TIRAGE AVEC REMISE D'UN ÉCHANTILLON DE n BOULES ****
#***DANS UNE URNE CONTENANT R BOULES ROUGES ET Q BOULES D'AUTRES COULEURS ****
#***LA VARIABLE ÉTUDIÉE EST LE NOMBRE K DE BOULES ROUGES DANS L'ÉCHANTILLON ****
avecremise3 <- function(R = 6, Q = 14, n = 5, nbsim = 10000){
  M <- R + Q
  urne <- rep(c("rouge", "autre"), c(R, Q))
  vecnbrouges <- rep(NA, nbsim)
  for(i in 1:nbsim){
     echant <- sample(urne, n, replace = T)
     nbrouges <- sum(echant == "rouge")
     vecnbrouges[i] <- nbrouges
                   }
  distribnbrouges <- table(vecnbrouges) / nbsim
#--------------- Affichage des résultats et des graphiques------------
  cat("\nProportion de rouges dans l'urne :", R / M, "\n\n")
  print(distribnbrouges)
  dev.set(3)
  barplot(distribnbrouges, main = "TIRAGES AVEC REMISE")
}