This is an R Markdown Notebook. When you execute code within the notebook, the results appear beneath the code.

Try executing this chunk by clicking the Run button within the chunk or by placing your cursor inside it and pressing Ctrl+Shift+Enter.

library(tidyverse)
── Attaching packages ──────────────────────────────────────────────────────────────────────────────────────────────── tidyverse 1.2.1 ──
✔ ggplot2 3.1.0       ✔ purrr   0.3.1  
✔ tibble  2.0.1       ✔ dplyr   0.8.0.1
✔ tidyr   0.8.3       ✔ stringr 1.4.0  
✔ readr   1.3.1       ✔ forcats 0.4.0  
── Conflicts ─────────────────────────────────────────────────────────────────────────────────────────────────── tidyverse_conflicts() ──
✖ dplyr::filter() masks stats::filter()
✖ dplyr::lag()    masks stats::lag()
adjust_p <- function(p, S) {
  padj <-((4*S-1)*p-3)/(4*(S-1))
  return(padj)
}
get_p <- function(p, S) {
  padj <- p*((S-1)/(S - 0.25)) +3/(4*S-1)
  return(padj)
}
S<- c(8)
p<- .97
p_t <- get_p(p,S)
dt <- tibble(r = 1:100) %>% mutate(mu_p =1-p^r,
                                     mu_pt=1-p_t^r )
dt_s <- tibble(S = 2:20) %>% mutate(p_97 =get_p(.97, S),
                                    p_98 =get_p(.98, S),
                                  p_99 =get_p(.99, S))
ggplot(dt_s) +
  geom_point(aes(S, p_97, colour="0.97"))+
   geom_point(aes(S, p_98,colour="0.98"))+
  geom_point(aes(S,p_99,colour="0.99"))

ggplot(dt) +
  geom_point(aes(r,mu_p),colour="blue")+
  geom_point(aes(r,mu_pt),colour="red")

LS0tCnRpdGxlOiAiUiBOb3RlYm9vayIKb3V0cHV0OiBodG1sX25vdGVib29rCi0tLQoKVGhpcyBpcyBhbiBbUiBNYXJrZG93bl0oaHR0cDovL3JtYXJrZG93bi5yc3R1ZGlvLmNvbSkgTm90ZWJvb2suIFdoZW4geW91IGV4ZWN1dGUgY29kZSB3aXRoaW4gdGhlIG5vdGVib29rLCB0aGUgcmVzdWx0cyBhcHBlYXIgYmVuZWF0aCB0aGUgY29kZS4gCgpUcnkgZXhlY3V0aW5nIHRoaXMgY2h1bmsgYnkgY2xpY2tpbmcgdGhlICpSdW4qIGJ1dHRvbiB3aXRoaW4gdGhlIGNodW5rIG9yIGJ5IHBsYWNpbmcgeW91ciBjdXJzb3IgaW5zaWRlIGl0IGFuZCBwcmVzc2luZyAqQ3RybCtTaGlmdCtFbnRlciouIAoKCmBgYHtyfQpsaWJyYXJ5KHRpZHl2ZXJzZSkKYGBgCgpgYGB7cn0KYWRqdXN0X3AgPC0gZnVuY3Rpb24ocCwgUykgewogIHBhZGogPC0oKDQqUy0xKSpwLTMpLyg0KihTLTEpKQogIHJldHVybihwYWRqKQp9CmdldF9wIDwtIGZ1bmN0aW9uKHBfdCwgUykgewogIGEgPC0gMy8oNCpTLTEpCiAgcCA8LSBwX3QqKDEtYSkgKyBhCiAgcmV0dXJuKHApCn0KClM8LSBjKDgpCnBfdDwtIApwPC0gZ2V0X3AocF90LFMpCmBgYAoKYGBge3J9CmR0IDwtIHRpYmJsZShyID0gMToxMDApICU+JSBtdXRhdGUobXVfcCA9MS1wXnIsCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICBtdV9wdD0xLXBfdF5yICkKCmR0X3MgPC0gdGliYmxlKFMgPSAyOjIwKSAlPiUgbXV0YXRlKHBfOTcgPWdldF9wKC45NywgUyksCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIHBfOTggPWdldF9wKC45OCwgUyksCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICBwXzk5ID1nZXRfcCguOTksIFMpKQoKYGBgCgpgYGB7cn0KZ2dwbG90KGR0X3MpICsKICBnZW9tX3BvaW50KGFlcyhTLCBwXzk3LCBjb2xvdXI9IjAuOTciKSkrCiAgIGdlb21fcG9pbnQoYWVzKFMsIHBfOTgsY29sb3VyPSIwLjk4IikpKwogIGdlb21fcG9pbnQoYWVzKFMscF85OSxjb2xvdXI9IjAuOTkiKSkKYGBgCgoKYGBge3J9CmdncGxvdChkdCkgKwogIGdlb21fcG9pbnQoYWVzKHIsbXVfcCksY29sb3VyPSJibHVlIikrCiAgZ2VvbV9wb2ludChhZXMocixtdV9wdCksY29sb3VyPSJyZWQiKQpgYGAKCgo=