Skip to content

Instantly share code, notes, and snippets.

@jrnold
Last active December 26, 2015 08:49
  • Star 0 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
Star You must be signed in to star a gist
Save jrnold/7124461 to your computer and use it in GitHub Desktop.
Confidence interval simulation
library("shiny")
library("ggplot2")
library("plyr")
norm_mean <- 0
norm_sd <- 1
##' Take sample from normal dist and calculate confidence interval for normal
smpl_mean_sample_ci <- function(n, p=0.95, pop_mean=0, pop_sd=1) {
smpl <- rnorm(n, norm_mean, norm_sd)
smpl_sd <- sd(smpl)
smpl_mean <- mean(smpl)
tailprob <- (1 - p) / 2
q <- -qt(tailprob, df=(n - 1), lower.tail=TRUE)
se <- smpl_sd / sqrt(n)
ci <- data.frame(ci_lb = smpl_mean - q * se,
ci_ub = smpl_mean + q * se)
ci
}
shinyServer(function(input, output) {
sample_ci <- reactive({
input$draw
mutate(rdply(input$nsamples, smpl_mean_sample_ci(input$smplsize,
input$confidence / 100)),
n = seq_len(input$nsamples),
contains_mean = ((0 > ci_lb) & (0 < ci_ub)))
})
output$plot <- renderPlot({
print(ggplot(sample_ci(), aes(x = n, ymin = ci_lb, ymax = ci_ub,
colour = contains_mean))
+ geom_hline(xintercept = 0, colour="blue")
+ geom_linerange()
+ coord_flip()
+ scale_x_continuous("samples")
+ scale_y_continuous(sprintf("%d%% CI", input$confidence))
+ scale_colour_manual(values = c("red", "black"))
+ theme(legend.position = "none",
axis.text.y = element_blank(),
axis.ticks.y = element_blank()))
})
output$npct <- renderText({
c("Percent of CI containing the population mean:", round(mean(sample_ci()$contains_mean) * 100, 2), "%")
})
})
library("shiny")
shinyUI(pageWithSidebar(
headerPanel("Confidence Intervals"),
sidebarPanel(
actionButton("draw", "Draw samples"),
numericInput("confidence", "Confidence (%):", 95, min = 1, max = 99,
step=1),
numericInput("smplsize", "Sample sizes:", 100, min = 2, max = 1000,
step=1),
numericInput("nsamples", "Number of samples:", 100, min = 1, max = 500,
step=1)
),
mainPanel(
textOutput("npct"),
plotOutput("plot"),
p("The population has a mean of 0 and a standard deviation of 1.",
"Intervals including the population mean are colored red, those not including the population mean are colored black")
)))
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment