#| '!! shinylive warning !!': |
#| shinylive does not work in self-contained HTML documents.
#| Please set `embed-resources: false` in your metadata.
#| standalone: true
#| viewerHeight: 1000
library(shiny)
library(tidyverse)
# Define UI for application that draws a histogram
ui <- fluidPage(
# Application title
titlePanel("Power and Sample Size"),
# Sidebar with a slider input for number of bins
sidebarLayout(
sidebarPanel(
sliderInput("n",
"Sample Size (n):",
min = 1,
max = 5000,
value = 500),
sliderInput("sd",
HTML("Standard Deviation (σ):"),
min = 1,
max = 20,
value = 10),
sliderInput("diff",
HTML("True Difference ( μ<sub>A</sub>- μ<sub>0</sub>)"),
min = 0.1,
max = 5,
value = 2),
numericInput("alpha",
HTML("α (type I error threshold)") ,
min = 0.01,
max = 0.1,
value = 0.05),
radioButtons("sides", "One vs. Two-Tailed",
choices = list("1-tailed" , "2-tailed" ), selected = "2-tailed"),
actionButton("rescale","Rescale Plot")),
# Show a plot of the generated distribution
mainPanel(
tabsetPanel(type = "tabs",
tabPanel("One-Sample\nT-test",
plotOutput("distPlot1"),
textOutput("power1"),
textOutput("print1")),
tabPanel("Two-Sample\n Difference in Means", plotOutput("distPlot"),
textOutput("power"),
textOutput("note"),
textOutput("print"))
))
)
)
server <- function(input, output) {
rvmin <- reactiveValues(lim = NULL)
rvmax <- reactiveValues(lim = NULL)
l <- reactiveValues(val = 0)
observeEvent(input$rescale, {
newmn = -5*sqrt(2*input$sd^2/input$n)
newmx = input$diff + 5*sqrt(2*input$sd^2/input$n)
rvmin$lim <- c(rvmin$lim,newmn)
rvmax$lim <- c(rvmax$lim,newmx)
l$val = l$val+1
})
output$distPlot <- renderPlot({
mn <- ifelse(l$val<1, -5*sqrt(2*10^2/500), rvmin$lim[l$val] )
mx <- ifelse(l$val<1, 2 + 5*sqrt(2*10^2/500), rvmax$lim[l$val] )
if(input$sides == "2-tailed"){
ct = sqrt(input$sd^2/input$n)*qnorm(1-(input$alpha/2))
ctl = -sqrt(input$sd^2/input$n)*qnorm(1-(input$alpha/2))}
if(input$sides =="1-tailed"){
ct = sqrt(input$sd^2/input$n)*qnorm(1-(input$alpha))
ctl =-Inf
}
funcShaded <- function(x,ct,ctl) {
y <- dnorm(x, mean = input$diff, sd = sqrt(input$sd^2/input$n))
y[x < ct] <- NA
return(y)
}
funcShadedRej <- function(x,ct,ctl) {
y <- dnorm(x, mean = 0, sd = sqrt(input$sd^2/input$n))
y[x < ct & x>ctl] <- NA
return(y)
}
p1 <- ggplot(data.frame(x =c(seq(mn,mx,by=0.00001),ct,ctl)),
aes(x = x))+
stat_function(fun = dnorm, args = list(mean = 0, sd = sqrt(input$sd^2/input$n)),
aes(color = "Null"))+
stat_function(fun = dnorm, args = list(mean = input$diff, sd = sqrt(input$sd^2/input$n)),
aes(color = "Alternative"))+
theme_bw()+
labs(color = "Hypothesis")+
geom_vline(xintercept = ct, linetype = 2, linewidth = 0.8, color = "darkcyan")+
xlab(expression(mu[A]-mu[0]))+
labs( y = "Density", fill = element_blank())+
stat_function(fun=funcShaded, geom="area", aes(fill = "Power"), alpha=0.5,args = list(ct=ct,ctl=ctl))+
stat_function(fun=funcShadedRej, geom="area", aes(fill = "Rejection Region"), alpha=0.25,args = list(ct=ct,ctl=ctl))+
theme(legend.text = element_text(size = 10))
if(input$sides == "2-tailed"){
p1 <- p1+ geom_vline(xintercept = ctl, linetype = 2, linewidth =0.8,color = "darkcyan")
}
p1
})
output$power <- renderText({
if(input$sides == "2-tailed"){
power = power.t.test(n=input$n, delta=input$diff, sd=input$sd, sig.leve=input$alpha, power=NULL,
type="two.sample", alternative="two.sided")$power%>%round(digits = 4)}
if(input$sides == "1-tailed"){
power = power.t.test(n=input$n, delta=input$diff, sd=input$sd, sig.leve=input$alpha, power=NULL,
type="two.sample", alternative="one.sided")$power%>%round(digits = 4)}
paste("Power: ", power, collapse = "")
})
output$note<-renderText( "Note: Here we assume equal sample size and equal standard deviation in both groups.")
output$distPlot1 <- renderPlot({
mn <- ifelse(l$val<1, -5*sqrt(10^2/500), rvmin$lim[l$val] )
mx <- ifelse(l$val<1, 2 + 5*sqrt(10^2/500), rvmax$lim[l$val] )
if(input$sides == "2-tailed"){
ct = sqrt(input$sd^2/input$n)*qnorm(1-(input$alpha/2))
ctl = -sqrt(input$sd^2/input$n)*qnorm(1-(input$alpha/2))}
if(input$sides =="1-tailed"){
ct = sqrt(input$sd^2/input$n)*qnorm(1-(input$alpha))
ctl =-Inf
}
funcShaded <- function(x,ct,ctl) {
y <- dnorm(x, mean = input$diff, sd = sqrt(input$sd^2/input$n))
y[x < ct] <- NA
return(y)
}
funcShadedRej <- function(x,ct,ctl) {
y <- dnorm(x, mean = 0, sd = sqrt(input$sd^2/input$n))
y[x < ct & x>ctl] <- NA
return(y)
}
p2 <- ggplot(data.frame(x = rep(seq(mn,mx,by=0.00001))),
aes(x = x))+
stat_function(fun = dnorm, args = list(mean = 0, sd = sqrt(input$sd^2/input$n)),
aes(color = "Null"))+
stat_function(fun = dnorm, args = list(mean = input$diff, sd = sqrt(input$sd^2/input$n)),
aes(color = "Alternative"))+
theme_bw()+
labs(color = "Hypothesis")+
geom_vline(xintercept = ct, linetype = 2, size = 0.8, color = "darkcyan")+
xlab(expression(mu[A]-mu[0]))+
labs( y = "Density", fill = element_blank())+
stat_function(fun=funcShaded, geom="area", aes(fill = "Power"), alpha=0.5,args = list(ct=ct,ctl=ctl))+
stat_function(fun=funcShadedRej, geom="area", aes(fill = "Rejection Region"), alpha=0.25,args = list(ct=ct,ctl=ctl))+
theme(legend.text = element_text(size = 10))
if(input$sides == "2-tailed"){
p2 <- p2 + geom_vline(xintercept = ctl, linetype = 2, size = 0.8,color = "darkcyan")
}
p2
})
output$power1 <- renderText({
if(input$sides == "2-tailed"){
power = power.t.test(n=input$n, delta=input$diff, sd=input$sd, sig.leve=input$alpha, power=NULL,
type="one.sample", alternative="two.sided")$power%>%round(digits = 4)}
if(input$sides == "1-tailed"){
power = power.t.test(n=input$n, delta=input$diff, sd=input$sd, sig.leve=input$alpha, power=NULL,
type="one.sample", alternative="one.sided")$power%>%round(digits = 4)}
paste("Power: ", power, collapse = "")
})
}
# Run the application
shinyApp(ui = ui, server = server)