library(shiny)
library(mosaic)
library(lattice)
library(latticeExtra)
shinyServer(function(input, output) {
FindingT <- function(Lamb, MIN, MAX){
N <- rpois(n = 1, lambda = Lamb)
xSubi <- runif(n = N, min = MIN, max = MAX)
sumXis <- sum(xSubi)
}
tableData <- reactive({
Simulation <- do(input$numsim)*FindingT(input$Lamb, input$MIN, input$MAX)
ExpectationT <- mean(Simulation$FindingT)
ExpectationN <- input$Lamb
ExpectationXi <- .5*(input$MIN + input$MAX)
Theoretical_ExpectationT <- ExpectationN*ExpectationXi
ForConfindence <- Simulation - Theoretical_ExpectationT
confidence <- quantile(ForConfindence$FindingT, c(0.05, 0.95))
confidence2 <- quantile(Simulation$FindingT, c(0.05, 0.95))
c1 <- confidence2[1]
c2 <- confidence2[2]
Simulated_VarianceT <- var(Simulation$FindingT)
VarianceXi <- (1/12)*(input$MAX - input$MIN)^2
VarianceN <- input$Lamb
Theoretical_VarianceT <- (ExpectationXi^2)*VarianceN + ExpectationN*VarianceXi
list(Theoretical_ExpectationT=Theoretical_ExpectationT,
ExpectationT=ExpectationT,
Theoretical_VarianceT=Theoretical_VarianceT,
Simulated_VarianceT=Simulated_VarianceT,
ExpectationN=ExpectationN,
VarianceN=VarianceN,
ExpectationXi=ExpectationXi,
VarianceXi=VarianceXi,
Simulation=Simulation,
Theoretical_ExpectationT=Theoretical_ExpectationT,
c1=c1,
c2=c2)
})
output$d1 <- renderTable({
data.frame(Theoretical = c("Expectation of N", "Variance of N",
"Expectation of Xi", "Variance of Xi"),
Values = c(tableData()$ExpectationN, tableData()$VarianceN,
tableData()$ExpectationXi, tableData()$VarianceXi))
})
output$d2 <- renderTable({
data.frame(ValueofT = c("Theoretical", "Simulated"),
Expectation = c(tableData()$Theoretical_ExpectationT, tableData()$ExpectationT),
Variance = c(tableData()$Theoretical_VarianceT, tableData()$Simulated_VarianceT),
StandardDeviation = c(sqrt(tableData()$Theoretical_VarianceT), sqrt(tableData()$Simulated_VarianceT)))
})
output$debug <- renderText({
c("expectation:", tableData()$Theoretical_ExpectationT,
"lower bound:", round(tableData()$c1, 3), "upper bound:",
round(tableData()$c2, 3))
})
graphing <- function(points, conf, ex){
cat("conf=", conf, "\n")
a <- histogram(~points, xlab="E[T]", main="Distribution of Simulated E[X]", col="slategray2")
b <- latticeExtra :: layer(panel.abline(v=conf, col="red", lwd=3), data=list(conf=conf))
c <- latticeExtra :: layer(panel.abline(v=ex, col="blue", lwd=3), data=list(ex=ex))
d <- a + b + c
return(d)
}
output$plotUnif <- renderPlot({
datas <- tableData()$Simulation$FindingT
lower <- tableData()$c1
upper <- tableData()$c2
conf <<- c(lower, upper)
theoretical <- tableData()$Theoretical_ExpectationT
graphing(points = datas, conf= c(lower, upper), ex= theoretical)
})
output$TEu1 <- renderUI({
withMathJax(
helpText("$$= \\left(\\frac{\\theta_{min}+\\theta_{max}}{2}\\right) \\lambda $$")
)})
output$TEu2 <- renderText({
c("= ", tableData()$Theoretical_ExpectationT)
})
output$summation <- renderUI({
withMathJax(
helpText('$$T = \\sum_{i=1}^{N} X_i$$')
)
})
output$expect1 <- renderUI({
withMathJax(
helpText('$$E[T|N] = NE[Xi]$$')
)
})
output$expect2 <- renderUI({
withMathJax(
helpText("$$E(T) = E[NE(X)] = E(N)E(X)$$"),
helpText("$$E(T) = E[NE[S]] = E[N]E(S)$$")
)
})
output$var1 <- renderUI({
withMathJax(
helpText("$$Var(T|N) = N Var(X)$$")
)
})
output$var2 <- renderUI({
withMathJax(
helpText("$$E[Var(T|N)] = E[N]Var(X)$$")
)
})
output$var3 <- renderUI({
withMathJax(
helpText("$$E(T|N)=NE[X]$$")
)
})
output$var4 <- renderUI({
withMathJax(
helpText("$$Var(E[T|N]) = E[X]^2 Var(N)$$")
)
})
output$var5 <- renderUI({
withMathJax(
helpText("$$Var(T) = E[X]^2Var(N) + E[N]Var(X)$$")
)
})
output$adams <- renderUI({
withMathJax(
helpText("$$E[T] = E[E[T|N]]$$")
)
})
output$eves <- renderUI({
withMathJax(
helpText("$$Var(T) = E[Var(T|N)] + Var(E[T|N])$$")
)
})
output$TVu1 <- renderUI({
withMathJax(
helpText(c("$$Var(T) = [E[X]]^2Var(N) + E[N]Var(X)$$")))
})
output$TVu2 <- renderUI({
withMathJax(
helpText(c("$$= \\left(\\frac{(\\theta_{min} + \\theta_{max})}{2}\\right)^2 \\lambda +
\\lambda \\left( \\frac{(\\theta_{min} + \\theta_{max})^2}{12} \\right)^2$$")))
})
output$equations5LATEX <- renderUI({
withMathJax(
helpText(c("$$= \\left(\\frac{(", input$MIN, "+", input$MAX, ")}{2}\\right)^2 \\lambda +
\\lambda \\left( \\frac{(\\theta_{min} + \\theta_{max})^2}{12} \\right)^2$$")))
})
output$TVu3<- renderText({c("= (.5(", input$MIN, " + ", input$MAX, "))^2 x ",
input$lamb, " + ", input$lamb, " x (1/12(", input$MIN, " + ",
input$MAX, ")^2)^2")})
output$TVu4 <- renderText({c("= ", round(tableData()$Theoretical_VarianceT,2))})
output$BelowPlot <- renderText("As shown above, the Theoretical E[T] (shown in blue), the total amount
spent at the store in one day, is within the bounds of a 95% confidence
interval of the Simulated E[T] (shown in red).")
output$BelowPlot2 <- renderText("As shown above, the Theoretical E[T] (shown in blue), the total amount
spent at the store in one day, is within the bounds of a 95% confidence
interval of the Simulated E[T] (shown in red).")
output$BelowPlot3 <- renderText("As shown above, the Theoretical E[T] (shown in blue), the total amount
spent at the store in one day, is within the bounds of a 95% confidence
interval of the Simulated E[T] (shown in red).")
FindingTNorm <- function(LambN, mean1, sd1){
Lamb <- input$LambN
mean1 <- input$mean1
sd1 <- input$sd1
N <- rpois(n = 1, lambda = LambN)
xSubi <- rnorm(n = N,mean = mean1, sd = sd1)
NormXsum <- sum(xSubi)
}
tableDataNorm <- reactive({
Simu <- do(input$numsimN)*FindingTNorm(input$LambN, input$mean1, input$sd1)
ExpectT <- mean(Simu$FindingTNorm)
ExpectN <- input$LambN
ExpectXi <- input$mean1
TExpectationT <- ExpectN*ExpectXi
ForConfi <- Simu - TExpectationT
confi<- quantile(ForConfi$FindingTNorm, c(0.05, 0.95))
confi2 <- quantile(Simu$FindingTNorm, c(0.05, 0.95))
c1N <- confi2[1]
c2N <- confi2[2]
SVarianceT <- var(Simu$FindingTNorm)
VarXi <- (input$sd1)^2
VarN <- input$LambN
TVarianceT <- (ExpectXi^2)*VarN + ExpectN*VarXi
list(TExpectationT=TExpectationT,
ExpectT=ExpectT,
TVarianceT=TVarianceT,
SVarianceT=SVarianceT,
ExpectN=ExpectN,
VarN=VarN,
ExpectXi=ExpectXi,
VarXi=VarXi,
Simu=Simu,
TExpectationT=TExpectationT,
c1N=c1N,
c2N=c2N)
})
output$d4 <- renderTable({
data.frame(ValuesforT = c("Theoretical", "Simulated"),
Expectation = c(tableDataNorm()$TExpectationT, tableDataNorm()$ExpectT),
Variance = c(tableDataNorm()$TVarianceT, tableDataNorm()$SVarianceT),
StandardDeviation = c(sqrt(tableDataNorm()$TVarianceT), sqrt(tableDataNorm()$SVarianceT)))
})
output$d3 <- renderTable({
data.frame(Theoretical = c("Expectation of N", "Variance of N",
"Expectation of Xi", "Variance of Xi"),
Values = c(tableDataNorm()$ExpectN, tableDataNorm()$VarN,
tableDataNorm()$ExpectXi, tableDataNorm()$VarXi))
})
output$plotNorm <- renderPlot({
datas <- tableDataNorm()$Simu$FindingTNorm
lower1 <- tableDataNorm()$c1N
upper1 <- tableDataNorm()$c2N
theory <- tableDataNorm()$TExpectationT
graphing(points = datas, conf= c(lower1, upper1), ex= theory)
})
output$debug2 <- renderText({
c("expectation:", tableDataNorm()$TExpectationT, "lower bound:",
round(tableDataNorm()$c1N, 3), "upper bound:", round(tableDataNorm$c2N, 3))
})
output$TEn1 <- renderUI({
withMathJax(
helpText("$$= \\mu\\lambda $$")
)})
output$TEn2 <- renderText({
c("= ", tableDataNorm()$TExpectationT)
})
output$TVn1 <- renderUI({
withMathJax(
helpText(c("$$Var(T) = E[X]^2Var(N) + E[N]Var(X)$$")))
})
output$TVn2 <- renderUI({
withMathJax(
helpText(c("$$= \\mu^2\\lambda + \\lambda\\sigma^2$$")))
})
output$TVn3<- renderText({c("= (.5(", input$MIN, " + ", input$MAX, "))^2 x ",
input$lamb, " + ", input$lamb, " x (1/12(", input$MIN, " + ",
input$MAX, ")^2)^2")})
output$TVn4 <- renderText({c("= ", round(tableDataNorm()$TVarianceT,2))})
FindingTGam <-
function(LambG, shape1, rate1){
LambG <- input$LambG
shape1 <- input$shape1
rate1 <- input$rate1
NG <- rpois(n = 1, lambda = LambG)
xSubiG <- rgamma(n = NG, rate = rate1, shape = shape1)
GamXsum <- sum(xSubiG)
return(GamXsum)
}
tableDataGamma <- reactive({
SimuGam <- do(input$numsimG)*FindingTGam(input$LambG, input$shape1, input$rate1)
ExT <- mean(SimuGam$FindingTGam)
Mid <- median(SimuGam$FindingTGam)
ExN <- input$LambG
ExXi <- input$shape1 / input$rate1
TExT <- ExN*ExXi
GamConf <- SimuGam - TExT
Gconf2 <- quantile(SimuGam$FindingTGam, c(0.05, 0.95))
c1G <- Gconf2[1]
c2G <- Gconf2[2]
SimVarG <- var(SimuGam$FindingTGam)
VarGXi <- input$shape1/(input$rate1)^2
VarGN <- input$LambG
GamTVarianceT <- (ExXi^2)*VarGN + ExN*VarGXi
list(TExT=TExT,
ExT=ExT,
GamTVarianceT=GamTVarianceT,
SimVarG=SimVarG,
ExN=ExN,
VarGN=VarGN,
ExXi=ExXi,
VarGXi=VarGXi,
SimuGam=SimuGam,
TExT=TExT,
c1G=c1G,
c2G=c2G)
})
output$tableG1 <- renderTable({
data.frame(Theoretical = c("Expectation of N", "Variance of N",
"Expectation of Xi", "Variance of Xi"),
Values = c(tableDataGamma()$ExN, tableDataGamma()$VarGN,
tableDataGamma()$ExXi, tableDataGamma()$VarGXi))
})
output$tableG2 <- renderTable({
data.frame(ValueofT = c("Theoretical", "Simulated"),
Expectation = c(tableDataGamma()$TExT, tableDataGamma()$ExT),
Variance = c(tableDataGamma()$GamTVarianceT, tableDataGamma()$SimVarG),
StandardDeviation = c(sqrt(tableDataGamma()$GamTVarianceT), sqrt(tableDataGamma()$SimVarG)))
}
)
output$plotGam <- renderPlot({
datas <- tableDataGamma()$SimuGam$FindingTGam
lower1 <- tableDataGamma()$c1G
upper1 <- tableDataGamma()$c2G
theory <- tableDataGamma()$TExT
graphing(points = datas, conf= c(lower1, upper1), ex= theory)
})
output$TEg1 <- renderUI({
withMathJax(
helpText("$$=\\left(\\frac{\\theta_{shape}}{\\theta_{rate}}\\right) \\lambda $$")
)})
output$TEg2 <- renderText({
c("= ", tableDataGamma()$TExT)
})
output$TVg1 <- renderUI({
withMathJax(
helpText(c("$$Var(T) = [E[X]]^2Var(N) + E[N]Var(X)$$")))
})
output$TVg2 <- renderUI({helpText(c("$$= \\left(\\frac{\\theta_{shape}}{\\theta_{rate}}\\right)^2 \\lambda +
\\lambda \\left( \\frac{\\theta_{shape}}{\\theta_{rate}^2} \\right)$$"))
})
output$TVg3<- renderText({c("= (.5(", input$MIN, " + ", input$MAX, "))^2 x ",
input$lamb, " + ", input$lamb, " x (1/12(", input$MIN, " + ",
input$MAX, ")^2)^2")})
output$TVg4 <- renderText({c("= ", round(tableDataGamma()$GamTVarianceT,2))})
})
library(shiny)
library(ggplot2)
shinyUI(navbarPage(
title = "Adam's and Eve's laws",
tabPanel("An Example",
sidebarLayout(position ="right",
sidebarPanel(
withMathJax(),
h1("Adam's Law"),
p("Adam's Law or the Law of Total Expectation states that when
given the coniditonal expectation of a random variable T which
is conditioned on N, you can find the expected value of unconditional
T with the following equation:"),
uiOutput('adams'),
hr(),
h1("Eve's Law"),
p("Eve's Law (EVVE's Law) or the Law of Total Variance is used to
find the variance of T when it is conditional on N,
it states that:"),
uiOutput('eves'),
hr(),
h4("Recall for this example, "),
p(strong("T"), " = total amount spent at the store"),
p(strong("N"), " = total number of customers that day"),
p(strong("Xi"), " = amount spent by the ith customer")
),
mainPanel(
withMathJax(),
p("The following example is drawn from examples D and E in Section 4 in Chapter 4 of the Third Edition of", em(
"Mathematical Statistics and Data Analysis"), "by John A. Rice."),
h1("An Example"),
p("Say we want to estimate
the total amount of money spent by all of the
customers that enter a store in one day. We'll call this E(T). This amount is conditional on
the amount of customers which enter the store on that day.
This is N. Let's say that N is a random Poisson variable
with a rate parameter \\(\\lambda\\). Xi is the expected amount of money each of the
N customers will spend. Each customer's spending is independent."
),
uiOutput('summation'),
hr(),
h3("Finding the Expectation"),
p("Let's say that each customer is expected to spend the same amount
E(X). Additionally, each customer's spending has the variance, Var(X).
Because all Xis are independent, we know that "),
uiOutput('expect1'),
p("Using", strong("Adam's law"),"and the fact above, we can conclude that"),
uiOutput('expect2'),
hr(),
h3("Finding the Variance"),
p("Because we know that"),
uiOutput('var1'),
p("we also know that"),
uiOutput("var2"),
p("With the knowledge that"),
uiOutput('var3'),
p("we can deduce"),
uiOutput('var4'),
p("Using", strong("Eve's Law"),"we can finally conclude that"),
uiOutput('var5'),
p("This app is designed to illustrate that these two conclusions are true with some
simple simulations. The panels on the top of this page give different options for
distributions of S, the amount each customer spends.")
)
)),
tabPanel('Uniform',
sidebarPanel(
withMathJax(),
p("This example models Xi (the amount each customer spends) from a", strong("uniform distribution."), "
Recall that N (the number of customers entering the store) comes from a",
strong("Poisson distribution."), " Enter parameter values below to define these distributions."),
numericInput("MIN", "\\(\\theta_{min}\\): The Minimum Amount Spent at Store by Each Customer", "0"),
numericInput("MAX", "\\(\\theta_{max}\\): The Maximum Amount Spent at Store by Each Customer", "100"),
numericInput("Lamb", "\\(\\lambda\\): The Rate to Determine How Many Customers Enter the Store", 5),
numericInput("numsim", "Number of simulations", 1000),
hr(),
h5("Theoretical Expectation of T"),
p("$$E(T)$$"),
p("$$= E[E(T|N)]$$"),
p("$$= E(X) E(N)$$"),
uiOutput("TEu1"),
textOutput("TEu2"),
h5("Theoretical Variance of T"),
uiOutput("TVu1"),
uiOutput("TVu2"),
textOutput("TVu4")
),
mainPanel(
plotOutput('plotUnif'),
textOutput('BelowPlot'),
hr(),
h4("Theoretical Values"),
tableOutput('d1'),
h4("Comparison to Simulated Values"),
tableOutput('d2'))),
tabPanel('Normal',
sidebarPanel(
p("This example models Xi (the amount each customer spends) from a", strong("normal distribution."), "
Recall that N (the number of customers entering the store) comes from a",
strong("Poisson distribution."), " Enter parameter values below to define these distributions."),
numericInput("mean1", "\\(\\mu\\): The Mean Amount Spent at Store by Each Customer", "50"),
numericInput("sd1", "\\(\\sigma\\): The Standard Deviation Amount Spent at Store by Each Customer", "5"),
numericInput("LambN", "\\(\\lambda\\): The Rate to Determine How Many Customers Enter the Store", 5),
numericInput("numsimN", "Number of simulations", 1000),
hr(),
h5("Theoretical Expectation of T"),
p("$$E(T)$$"),
p("$$= E[E(T|N)]$$"),
p("$$= E(X) E(N)$$"),
uiOutput("TEn1"),
textOutput("TEn2"),
h5("Theoretical Variance of T"),
uiOutput("TVn1"),
uiOutput("TVn2"),
textOutput("TVn4")
),
mainPanel(
plotOutput('plotNorm'),
textOutput('BelowPlot2'),
hr(),
h4("Theoretical Values"),
tableOutput('d3'),
h4("Comparison to Simulated Values"),
tableOutput('d4'))),
tabPanel('Gamma',
sidebarPanel(
p("This example models Xi (the amount each customer spends) from a", strong("gamma distribution."), "
Recall that N (the number of customers entering the store) comes from a",
strong("Poisson distribution."), " Enter parameter values below to define these distributions."),
numericInput("shape1", "\\(\\theta_{shape}\\): Shape Parameter for theAmount Spent at Store by Each Customer", "500"),
numericInput("rate1", "\\(\\theta_{rate}\\): Rate Parameter for the Amount Spent at Store by Each Customer", "10"),
numericInput("LambG", "\\(\\lambda\\): The Rate to Determine How Many Customers Enter the Store", 5),
numericInput("numsimG", "Number of simulations", 1000),
hr(),
h5("Theoretical Expectation of T"),
p("$$E(T)$$"),
p("$$= E[E(T|N)]$$"),
p("$$= E(X) E(N)$$"),
uiOutput("TEg1"),
textOutput("TEg2"),
h5("Theoretical Variance of T"),
uiOutput("TVg1"),
uiOutput("TVg2"),
textOutput("TVg4")
),
mainPanel(
plotOutput('plotGam'),
textOutput('BelowPlot3'),
hr(),
h4("Theoretical Values"),
tableOutput('tableG1'),
h4("Comparison to Simulated Values"),
tableOutput('tableG2')
)),
tabPanel('About',
p("This App was created by Coco Kusiak to create a visualization
of Adam's and Eve's laws designed to help students learning probability and its applications.
Further information about these laws can be found at:"),
wellPanel(
helpText( a("Adam's Law", href="https://en.wikipedia.org/wiki/Law_of_total_expectation", target="_blank"), br(),
a("Eve's Law", href="https://en.wikipedia.org/wiki/Law_of_total_variance", target="_blank")
)
),
p("A special thanks to Professor Nicholas Horton and Professor Susan Wang."),
p("Other resources include", em("Mathematical Statisics and Data Analysis"),
"by John A. Rice and ", em("Introduction to Probability"), "by Joseph K.
Blitzstein and Jessica Hwang."))
))