# # Sitzung des Stochastik-Praktikums, WS 2013/2014, JGU Mainz # 11.11.2013 ## "Casino R" : ## Kombinatorische Wahrscheinlichkeiten per Simulation ## einschaetzen, ## Gluecksspiele ## ########################################### ## Kombinatorik: Ziehen mit Zuruecklegen dinge <- c(0:9,'a','b','c','d','e','f') k <- 5 ## 'von Hand': res <- character(k) for (i in 1:k) { res[i] <- dinge[ceiling(runif(1,min=0, max=length(dinge)))] } res ## (funktional) dasselbe mit einem R-Befehl: sample(dinge, k, replace=TRUE) ## Siehe auch: ?sample ## ########################################### ## Kombinatorik: Ziehen ohne Zuruecklegen dinge <- c(0:9,'a','b','c','d','e','f') k <- 5 ## 'von Hand': res <- character(0) dinge.kopie <- dinge for (i in 1:k) { w <- ceiling(runif(1,min=0, max=length(dinge.kopie))) res[i]<-dinge.kopie[w] dinge.kopie <- dinge.kopie[-w] ## (laesst das w-te Element aus) } res ## (funktional) dasselbe mit einem R-Befehl: sample(dinge, k, replace=FALSE) ## ########################################### ## Beispiel: Zufaellige Permutationen n <- 10 perm <- sample(1:n, n, replace=FALSE) ## ergibt eine zufaellige Permutation von 1,...,n perm zyklusder1 <- function(p) { ## bestimme den Zyklus der 1 in der Permutation p zykl <- numeric(0) b <- 1 while (TRUE) { zykl <- c(zykl, b) b <- p[b] if (b==1) break } zykl } zyklusder1(perm) enthaeltderzyklusder1die2 <- function(p) { ## pruefe, ob der Zyklus der 1 in der Permutation p die 2 enthaelt b <- 1 while (TRUE) { b <- p[b] if (b==2) return(TRUE) if (b==1) break } FALSE } enthaeltderzyklusder1die2(perm) ## Wie oft kommt es vor, dass der Zyklus der 1 in ## einer zufaelligen Permutation die 2 enthaelt? rep <- 2000 res <- logical(rep) ## Simuliere rep mal: for (i in 1:rep) { perm <- sample(1:n, n, replace=FALSE) res[i] <- enthaeltderzyklusder1die2(perm) } sum(res)/rep ## Anteil der "Erfolge" ## dasselbe mit einem R-Befehl: res <- replicate(rep, enthaeltderzyklusder1die2(sample(1:n, n, replace=FALSE))) sum(res)/rep ## Siehe auch: ?replicate ## auch: ?sapply ## Fragen beim "Kniffel"/"Yahtzee": z.B. ## Wie oft wirft man einen Dreier-Pasch (auf ein Mal)? anz.seiten <- 6 ww <- sample(1:anz.seiten,5,replace=TRUE) ## 5-facher Wuerfelwurf ww hs <- function(x, anzseiten=6) { ## das "Haeufigkeitsspektrum" der Augenzahlen sapply(1:anzseiten, function(i) sum(x==i)) } hs(ww) ifelse(3 %in% hs(ww), 1, 0) ## Pruefen, ob eine Augenzahl 3 mal vorkommt ## (siehe auch ?%in%) res <- replicate(10000, ifelse(3 %in% hs(sample(1:anz.seiten,5,replace=TRUE)), 1, 0)) sum(res) mean(res) ## #################################################### # Simuliere einen "Irrfahrtspfad"=Gewinnprozess im Muenzwurfspiel: # In jeder Runde wird Münze mit Erfolgsw'keit p geworfen # bei Ergebnis "1" erhält Spieler A einen Euro von Spieler B, # bei Ergebnis "0" zahlt Spieler A einen Euro an Spieler B. n <- 20 # Anz. Schritte p <- 0.5 # Erfolgsw'keit der Münze k <- 0 # Startwert wuerfe <- (2*rbinom(n, size=1, prob=p)-1) plot(wuerfe) pfad <- k+cumsum(wuerfe) # pfad[i]=Gewinn von Spieler A nach i Runden plot(pfad, xlab="Schritte", type="b", main="Ein Spielverlauf") # füge zur besseren Sichtbarkeit Linien hinzu # Nehmen wir an, A und B spielen, bis einer 10 Euro verloren hat # (und der aktuelle Gewinn von A sei k) MW.Spiel <- function(k) { gewinn <- k while(gewinn > -10 && gewinn <10) gewinn <- gewinn+(2*rbinom(1,1,p)-1) if (gewinn==10) 1 else 0 } MW.Spiel(0) replicate(10, MW.Spiel(0)) sum(replicate(1000, MW.Spiel(0)))/1000 # Wie sieht die Gewinnwahrscheinlichkeit von A als Funktion des # Startwerts (="Vorsprung") aus? startwerte <- (-9):9 ewsA <- numeric(19) versuche <- 1000 for (i in 1:19) { ewsA[i]<-sum(replicate(versuche, MW.Spiel(startwerte[i])))/versuche cat("Startwert",startwerte[i],"geschaetzte Gewinnw'keit für A:",ewsA[i],"\n") } ewsA plot(startwerte, ewsA, xlab="Startwert", ylab="geschaetzte Gewinnw'keit für A") abline(0.5,1/20,lty=2) # Gerade 0.5+(1/20)*x einzeichnen # Analog kann man die Frage für p != 1/2 angehen # (dann sieht die exakte Antwort allerdings anders aus). ## #################################################### ## Black Jack mit R spielen/simulieren ## Spass mit unicode-Symbolen: spielkartenfarben2 <- c('\u2662', '\u2661', '\u2660', '\u2663') ## nur schwarze Symbole: spielkartenfarben <- c('\u2666', '\u2665', '\u2660', '\u2663') spielkartenwerte <- c('2','3','4','5','6','7','8','9','Z','B','D','K','A') ## Einen Satz Spielkarten erzeugen: deck <- character(52) for (i in 1:13) { for (j in 1:4) { deck[(j-1)*13+i] <- paste(spielkartenwerte[i], spielkartenfarben[j], sep='') } } deck wert <- function(h) { ## Berechne den Wert des Kartensatzes h anz.asse <- 0 wert <- 0 ## Addiere Kartenwerte (und zaehle mit, wieviele Asse) for (k in h) { kw <- substr(k,1,1) ## eigentlich interessiert uns nur das erste Zeichen if (kw %in% c('Z','B','D','K')) { ## %in% prueft Mitgliedschaft, vgl. ?%in% wert <- wert+10 } else if (kw=='A') { wert <- wert+1 anz.asse <- anz.asse+1 } else if (kw %in% c('2','3','4','5','6','7','8','9')) { wert <- wert + as.numeric(kw) } else { ## fange ab, falls eine "inexistente" Karte gegeben return(NA) } } ## Asse koennen als 11 zaehlen, solange der Gesamtwert <= 21 ist ## (wir haben sie oben aber nur als 1 gewertet) if (anz.asse>0 && wert <=11) { wert <- wert+10 } wert } ## Beispiel: wert(c('A','6','2')) ## (ergibt 19) ## Wir verwenden (globale) Variable fuer folgende Spielparameter: anz.decks.im.stapel <- 2 ## Wieviele Kartensaetze am Anfang im Stapel neu.mischen.grenze <- 20 ## wenn der Stapel unter diese Groesse faellt, ## wird neu gemischt strategie0 <- function(offen, geber, spieler) { ## Eine "Strategie" fuer Black Jack ## Eingaben (jeweils Vektor vom Typ character): ## offen = Karten, die bereits gefallen sind ## geber = Karten, die der Geber in der aktuellen Runde bereits hat ## spieler = Karten, die der Spieler in der aktuellen Runde bereits hat ## Ausgabe: eine ganze Zahl, Interpretation: ## Wenn length(geber)==0 (und dann sollte auch length(spieler)==0) sein: ## Dann beginnt eine neue Runde, die Ausgabe ist der Wetteinsatz fuer ## diese Runde. ## Wenn length(geber)>0 (und dann sollte auch length(spieler)>0) sein: ## Ausgabe 1 bedeutet: eine weitere Karte nehmen ## Ausgabe 0 bedeutet: keine Karte mehr nehmen ## Eine simple Strategie if (length(geber)==0) return(1) # wir setzen immer 1 ## Wir ahmen den Geber nach: ## nehme stets eine neue Karte, wenn der aktuelle Wert unter 17 ist if (wert(spieler)<17) return(1) else return(0) } strategie.interaktiv <- function(offen, geber, spieler) { ## Wir fragen den Benutzer: if (length(geber)==0) { cat('Ihr Einsatz (0=Spielende)') a <- scan(file='', what=integer(), n=1, quiet=TRUE) return(a) } else { cat('Der Geber hat:',geber,'\n','Sie haben:',spieler,'\n') cat('Eine weitere Karte nehmen? (j/n)') a <- scan(file='', what=character(), n=1, quiet=TRUE) if (a=='j') return(1) else return(0) } } strategie.zaehle <- function(offen, geber, spieler) { ## Variiere Strategie, je nachdem welche Karten schon gefallen sind: ## Zunaechst: Zaehle (geeignet) mit: karten <- c('2','3','4','5','6','7','8','9','Z','B','D','K','A') zaehlwerte <- c(1,1,1,1,1,0,0,0,-1,-1,-1,-1,-1) ## 2-6 geben +1, Z-A geben -1 zaehlsumme <- 0 for (k in offen) { zaehlsumme <- zaehlsumme+zaehlwerte[match(substr(k,1,1),karten)] } zaehlsumme <- zaehlsumme/ceiling((52*anz.decks.im.stapel-length(offen))/52) ## Faustregel: Dividiere zaehlsumme durch Anzahl noch im Stapel ## befindlicher ganzer Kartensaetze ## Variiere den Einsatz, je nachdem welche Karten schon gefallen sind: if (length(geber)==0) { ## Einsatz in Abhaengigkeit des Wertes von zaehlsumme if (zaehlsumme>3) { return(6) } else if (zaehlsumme>2) { return(4) } else if (zaehlsumme>1) { return(2) } else { return(0.1) } } ## Im eigentlichen Spiel: ## nehme eine neue Karte, wenn der aktuelle Wert unter limit ist, ## wobei limit je nach Wert von zaehlsumme variiert if(zaehlsumme <= -3) {limit <- 18} else if(zaehlsumme <= -2) {limit = 17} else if(zaehlsumme >= 3) {limit = 14} else if(zaehlsumme >= 2) {limit = 15} else {limit = 16} ## limit <- 17 ## oder sehr simpel ... if (wert(spieler)=52*anz.decks.im.stapel-neu.mischen.grenze) { ## und erzeugen ggfs. einen neu gemischten Kartenstapel: stapel <- rep(deck, times=anz.decks.im.stapel) ## und mischen: stapel <- sample(stapel, size=52*anz.decks.im.stapel, replace=FALSE) si <- 1 ## Index der aktuell obersten Karte nun = 1 } ## Es geht los ... if (si==1) { offen <- character(0) } else { offen <- stapel[1:(si-1)] } ## offen enthaelt nun die bereits ausgespielten Karten einsatz <- strat(offen,character(0),character(0)) ## Abfangen, ob (interaktiver) Spieler beendet: if(einsatz <= 0) break ## Spieler erhaelt zunaechst zwei Karten und Geber erhaelt eine Karte spi <- stapel[si:(si+1)] geb <- stapel[si+2] si <- si+3 ## Der Spieler koennte direkt einen Black Jack haben, ## dann gibt es (fuer ihn) nichts weiter zu entscheiden, ## sonst geht es los: while(wert(spi)<21) { if (strat(offen, geb, spi)==0) { ## dann bricht der Spieler ab break } ## sonst bekommt der Spieler eine weitere Karte: spi <-c(spi, stapel[si]) si <- si+1 } ## Wenn der Spieler bis jetzt nicht verloren hat, ist nun der ## Geber an der Reihe: if (wert(spi)<=21) { while (wert(geb)<17) { ## der Geber erhaelt solange Karten, bis er Kartenwert 17 erreicht geb <- c(geb,stapel[si]); si <- si+1 } ## Wir pruefen, wer nun gewonnen hat: ## (und fangen zunaechst alle "Black Jack"-Moeglichkeiten ab if (wert(geb)==21 && length(geb)==2 && wert(spi)==21 && length(spi)==2) { resultat[runde] <- 0 if (interaktiv) cat('Unentschieden: Geber', geb,' und Spieler', spi,'haben beide Black Jack.\n') } else if (wert(geb)==21 && length(geb)==2) { resultat[runde] <- -einsatz if (interaktiv) cat('Verloren: Geber', geb,' hat Black Jack und Spieler', spi,'nicht.\n') } else if (wert(spi)==21 && length(spi)==2) { resultat[runde] <- einsatz ## auch 1.5*einsatz if (interaktiv) cat('Gewonnen: Spieler', spi,' hat Black Jack und Geber', geb,'nicht.\n') } else { ## sonst sind wir im "Brot&Butter-Fall": niemand hat Black Jack if (wert(geb)>21) { resultat[runde] <- einsatz if (interaktiv) cat('Gewonnen: Spieler', spi,'hat',wert(spi),'<=21','und Geber', geb,'nicht.\n') } else if (wert(spi)>wert(geb)) { resultat[runde] <- einsatz if (interaktiv) cat('Gewonnen: Spieler', spi,'hat',wert(spi),'und Geber', geb,'hat',wert(geb),'\n') } else if (wert(spi) 21, d.h. verloren resultat[runde] <- -einsatz if (interaktiv) cat('Verloren: Spieler',spi,'hat',wert(spi),'>21.\n') } if (interaktiv) cat('Aktuell kumulierter Gewinn:',sum(resultat[1:runde]),'\n') runde <- runde+1 } ## Rueckgabe: Vektor der Resultate der einzelnen Runden resultat[1:(runde-1)] } ## Probieren wir's aus: anz.runden <- 1000 ## Interaktiv (R als Computerspiel) rstri <- spiele.BlackJack(strategie.interaktiv, anz.runden, interaktiv=TRUE) sum(rstri) plot(cumsum(rstri), type='l', xlab='Spiele', ylab='Kumulierter Gewinn') ## Simple Strategie rstr0 <- spiele.BlackJack(strategie0, anz.runden) sum(rstr0) plot(cumsum(rstr0), type='l', xlab='Spiele', ylab='Kumulierter Gewinn') ## Strategie mit Zaehlung rstr1 <- spiele.BlackJack(strategie.zaehle, anz.runden) sum(rstr1) plot(cumsum(rstr1), type='l', xlab='Spiele', ylab='Kumulierter Gewinn') ## anz.runden <- 10000 rstr1 <- spiele.BlackJack(strategie.zaehle, anz.runden) sum(rstr1) #### ## riesige Stichprobe (Achtung: laeuft u.U. _sehr_ lange) anz.runden <- 1e6 rstr0 <- spiele.BlackJack(strategie0, anz.runden) sum(rstr0) mean(rstr0) sd(rstr0) rstr1 <- spiele.BlackJack(strategie.zaehle, anz.runden) sum(rstr1) mean(rstr1) sd(rstr1)