README.Rmd 2.0 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105
  1. ---
  2. title: "Working with model objects"
  3. output: github_document
  4. ---
  5. ## Install
  6. - Download and install R [https://www.r-project.org](https://www.r-project.org) (and RStudio [https://www.rstudio.com/products/rstudio/](https://www.rstudio.com/products/rstudio/)).
  7. - Go to R console or open scripts/README.Rmd in RStudio.
  8. - Install these packages (this needs to be done once).
  9. ```{r}
  10. if(!require(brms)){
  11. install.packages(c("brms", "here"))
  12. library(brms)
  13. library(here)
  14. }
  15. ```
  16. ## Run
  17. In R console OR in .Rmd file:
  18. - Load packages to R environment.
  19. ```{r, message=FALSE}
  20. library(brms)
  21. library(here)
  22. ```
  23. - Load the model object and print model summary.
  24. ```{r}
  25. m <- readRDS(here("models/anticons_detool.rds"))
  26. print(m)
  27. ```
  28. - Get fitted coefficients with 95% credible intervals.
  29. ```{r}
  30. posterior_summary(m)
  31. ```
  32. - Get the full posterior.
  33. ```{r}
  34. post <- posterior_samples(m)
  35. head(post)
  36. ```
  37. - What is the estimated difference in the proportion of anti-conservative p value histograms between DESeq2 and EdgeR?
  38. ```{r, eval=FALSE}
  39. posterior_deseq_edger <- inv_logit_scaled(post$b_de_tooldeseq - post$b_de_tooledger)
  40. hist(posterior_deseq_edger, breaks = 40)
  41. ```
  42. ```{r, echo=FALSE, message=FALSE}
  43. posterior_deseq_edger <- inv_logit_scaled(post$b_de_tooldeseq - post$b_de_tooledger)
  44. png(here("plots/posterior.png"))
  45. hist(posterior_deseq_edger, breaks = 40)
  46. invisible(dev.off())
  47. ```
  48. ![](plots/posterior.png)
  49. - The posterior summary for the effect size.
  50. ```{r}
  51. posterior_summary(posterior_deseq_edger)
  52. ```
  53. ```{r, echo=FALSE}
  54. ps <- posterior_summary(posterior_deseq_edger)
  55. ```
  56. The estimated effect size is somewhere between `r paste(signif(ps[1,3:4], digits = 3) * 100, collapse = " and ")` percentage points.
  57. - Get data from model object.
  58. ```{r}
  59. data <- m$data
  60. head(data)
  61. ```
  62. - Extract stan code from model object. This is the fullest model description.
  63. ```{r}
  64. stancode(m)
  65. ```
  66. ## This document
  67. This README.md was generated by running:
  68. ```{r, eval=FALSE}
  69. rmarkdown::render("scripts/README.Rmd", output_file = here::here("README.md"))
  70. ```