diff --git a/app.R b/app.R new file mode 100644 index 0000000000000000000000000000000000000000..be8b507d48761b2d04c8a6dab224ac7557184371 --- /dev/null +++ b/app.R @@ -0,0 +1,416 @@ +########################################################################### +## bootCI Shiny/R app ui.R ## +## ## +## Author Eugen Pircalabelu https://perso.uclouvain.be/eugen.pircalabelu ## +## For RShiny@UCLouvain http://sites.uclouvain.be/rshiny ## +## ## +## Licences : CC-BY for http://sites.uclouvain.be/RShiny ## +## GPL for source code on ## +## https://forge.uclouvain.be/rshiny_uclouvain/bootCI ## +########################################################################### + +library(shiny) +library(shinyjs) +library(RColorBrewer) +library(shinyWidgets) +library(shinycssloaders) +library(xtable) +library(shinyBS) +shinyUI(fluidPage( + headerPanel("Non-parametric Bootstrap"), + + sidebarLayout( + + sidebarPanel( + tags$head( + tags$style(type="text/css", "label { display: inline; }"), + tags$style(type="text/css", '.checkbox input[type="checkbox"],.radio input[type="radio"] { float: none; }') + ), + + HTML("X~Exp(μ) with μ = <br>"), + numericInput("mu", " " , min = 0,max = 100,value = 50, step=10), + + HTML("Sample size n = <br>"), + numericInput("n","",min = 10,max = 30,value = 10, step=5), + + actionButton(inputId="submit_model",label="Sample",width="150px"), + hr(), + + HTML("Number of bootstraap resamples B = <br>" ), + numericInput("B", "", value=25, min = 1, max = 50,step=1), + + HTML("Results for b-th bootstraap sample where b = <br>" ), + sliderInput("bb", "", value=1, min = 1, max = 5,step=1,animate=T), + + actionButton(inputId="refresh",label="Refresh",width="150px"), + + p(HTML("<A HREF=\"javascript:history.go(0)\">Click here to restart the experiment</A>")), + + + HTML('<hr style="border:1px solid #ccc;"/>'), + HTML('<a rel="license" href="http://creativecommons.org/licenses/by/2.0/be/"><img alt="Licence Creative Commons" style="border-width:0" + src="http://i.creativecommons.org/l/by/2.0/be/80x15.png" /></a> Ce(tte) oeuvre de <span xmlns:cc="http://creativecommons.org/ns#" + property="cc:attributionName"> <font face="Courier"> RShiny@UCLouvain </font> </span> est mise à disposition selon les termes de la <a rel="license" + href="http://creativecommons.org/licenses/by/2.0/be/">licence Creative Commons Attribution 2.0 Belgique</a>.'), + HTML('<p>Détails sur l\'utilisation de cette ressource sur <a href="http://sites.uclouvain.be/RShiny" + target="_blank"><font face="Courier"> RShiny@UCLouvain </font></a><br/> + Code source disponible sur <a href="https://forge.uclouvain.be/rshiny_uclouvain/bootCI" target="_blank">GitLab</a></p>') + + ), + + mainPanel( + plotOutput("Plot",height=800, width = "99%") + + # fluidRow(withSpinner( + # tableOutput(outputId="res_table"))) + + ) +) +) +) + + + +# Define a server for the Shiny app +function(input, output,session) { + + + # Create values which will contain the changing values + myvalues <- reactiveValues() + + values <- reactiveValues(tx_star_b=c()) +values2 <- reactiveValues(sigmax_star_b=c()) +values3 <- reactiveValues(tx_star_bMINUStheta_fn=c()) +values4 <- reactiveValues(u_star_b=c()) +values5 <- reactiveValues(xbar=c()) + + + + + # Create an empty dataframe (within values) which will contain the generated data + myvalues$df <- data.frame(Low= integer(0),Up = numeric(0), stringsAsFactors = FALSE) + + + + + + res <- reactive({ + +mu=input$mu +n=input$n +B=input$B +N=20 + +b=input$bb + + +xbar = NULL +tx_star_b=NULL +sigmax_star_b=NULL +u_star_b=NULL +tx_star_bMINUStheta_fn=NULL + + + + + +nrow = N; ncol = N;lambda=1/mu; p.col = c("deepskyblue", "red"); p.cex = c(1, 3) +set.seed(2018+as.numeric(input$submit_model)) +population <- rexp(nrow*ncol,rate=lambda) +population <- matrix(population,nrow,ncol) + +par(mar=c(0.5, 0.5, 0.5, 0.5)) +x = cbind(rep(1:ncol, nrow), gl(nrow, ncol)) +populationIndex = cbind(rep(1:ncol, nrow), gl(nrow, ncol)) + + +set.seed(2018+as.numeric(input$submit_model)) +myx=x[sample(nrow * ncol, n), ] +myxx=population[myx] + + + +set.seed(2018+as.numeric(input$submit_model)) +originalsampleIndex=x[sample(nrow * ncol, n), ] + + +for (b in 1:B){ + +set.seed(2018+b+(B*n*as.numeric(input$submit_model))) +BootSampleIndex=sample(1:n, n,replace = T) +mypcex=rep(p.cex[2],length(originalsampleIndex)) +mypcex2=seq(0,100,by=2) + +freq_vec=rep(NA,length(BootSampleIndex)) +mypcex3=rep(NA,length(BootSampleIndex)) +for (j in 1:length(unique(BootSampleIndex))){ +mylength=length(freq_vec[which(BootSampleIndex==unique(BootSampleIndex)[j])]) +freq_vec[which(BootSampleIndex==unique(BootSampleIndex)[j])]=1:mylength +mypcex3[which(BootSampleIndex==unique(BootSampleIndex)[j])]=mypcex2[1:mylength] +} + +mypcex4=mypcex+mypcex3 +mypcex4 + +BootSampleIndexMatrix=originalsampleIndex[BootSampleIndex, ] +BootSample=population[BootSampleIndexMatrix] + + # cat(file=stderr(), " ", as.numeric(guessInput5()), " ", "\n") + + orig.mean = mean(population[originalsampleIndex]) + orig.var = var(population[originalsampleIndex])*((n-1)/n) + +# cat(file=stderr(), " ", as.numeric(orig.mean), " ", "\n") + +tx_star_b <- c(tx_star_b,mean(population[BootSampleIndexMatrix])) +sigmax_star_b=c(sigmax_star_b,sqrt(var(population[BootSampleIndexMatrix])*(n-1)/n)) +tx_star_bMINUStheta_fn=c(tx_star_bMINUStheta_fn,mean(population[BootSampleIndexMatrix])-mean(population[originalsampleIndex])) +u_star_b=c(u_star_b,sqrt(n)*(mean(population[BootSampleIndexMatrix])-mean(population[originalsampleIndex]))/sqrt(var(population[BootSampleIndexMatrix])*(n-1)/n)) +xbar=c(xbar, as.numeric(mean(BootSample))) + +} + + cat(file=stderr(), " ", as.numeric(tx_star_b), " ", "\n") + + +CIboot=c(orig.mean-quantile(tx_star_bMINUStheta_fn,1-0.05/2),orig.mean-quantile(tx_star_bMINUStheta_fn,0.05/2)) +CIboot2=c(2*orig.mean-quantile(tx_star_b,1-0.05/2),2*orig.mean-quantile(tx_star_b,0.05/2)) +CItboot=c(orig.mean-sqrt(orig.var)*quantile(u_star_b,1-0.05/2)/sqrt(n),orig.mean-sqrt(orig.var)*quantile(u_star_b,0.05/2)/sqrt(n)) +CIpboot=c(orig.mean+quantile(tx_star_bMINUStheta_fn,0.05/2),orig.mean+quantile(tx_star_bMINUStheta_fn,1-0.05/2)) +CIasymp=c(orig.mean-qnorm(1-0.05/2,mean=0,sd=sqrt(orig.var))/sqrt(n),orig.mean-qnorm(0.05/2,mean=0,sd=sqrt(orig.var))/sqrt(n)) + + + return(c(as.numeric(CIboot[1]),as.numeric(CIboot[2]),as.numeric(b))) + }) + + + +newEntry <- observe({ + if(input$submit_model >= 0) { + isolate( + myvalues$df[nrow(myvalues$df) + 1,] <- c(isolate(res())[1],isolate(res())[2]) + )} + }) + + + +output$res_table <- renderTable({ + myvalues$df + }) + + +observeEvent( input$refresh, { + updateSliderInput(session, "bb", label="", value=1, min=1, max=input$B, step=1) + set.seed(2018+as.numeric(input$submit_model)) +tx_star_b=c() +sigmax_star_b=c() +tx_star_bMINUStheta_fn=c() +u_star_b=c() +xbar=c() + }) + +observeEvent( input$submit_model, { +tx_star_b=c() +sigmax_star_b=c() +tx_star_bMINUStheta_fn=c() +u_star_b=c() +xbar=c() + }) + + + + + # Fill in the spot we created for a plot + output$Plot <- renderPlot({ + + + +par(mfrow=c(4,2)) + +mu=input$mu +n=input$n +B=input$B +b=input$bb +N=20 + +xbar = NULL +tx_star_b=NULL +sigmax_star_b=NULL +u_star_b=NULL +tx_star_bMINUStheta_fn=NULL + + +nrow = N; ncol = N;lambda=1/mu; p.col = c("deepskyblue", "red"); p.cex = c(1, 3) +set.seed(2018+as.numeric(input$submit_model)) +population <- rexp(nrow*ncol,rate=lambda) +population <- matrix(population,nrow,ncol) + +par(mar=c(0.5, 0.5, 0.5, 0.5)) +x = cbind(rep(1:ncol, nrow), gl(nrow, ncol)) +populationIndex = cbind(rep(1:ncol, nrow), gl(nrow, ncol)) + +set.seed(2018+as.numeric(input$submit_model)) +myx=x[sample(nrow * ncol, n), ] +myxx=population[myx] + +mycol<-NULL +mycol<-c(mycol,"red") +mycol[]<-"gray" +mycol[1]<-"red" + +plot(x, pch = 19, col = "deepskyblue", cex = p.cex[1], axes = FALSE, ann = FALSE, xlab = "", ylab = "",bty="n") +points(myx, col = p.col[1], pch = 19, cex = p.cex[1]) +points(myx, col = p.col[2], cex = p.cex[2]) + + + +set.seed(2018+as.numeric(input$submit_model)) +originalsampleIndex=x[sample(nrow * ncol, n), ] + + + orig.mean = mean(population[originalsampleIndex]) + + +par(mar=c(0,0,0,0)) + +plot(x, pch = 19, col = "white", cex = p.cex[1], axes = FALSE, ann = FALSE, xlab = "", ylab = "",bty="n") +points(originalsampleIndex, col = p.col[1], pch = 19, cex = p.cex[1]) + +set.seed(2018+b+(B*n*as.numeric(input$submit_model))) +BootSampleIndex=sample(1:n, n,replace = T) +mypcex=rep(p.cex[2],length(originalsampleIndex)) +mypcex2=seq(0,100,by=2) + +freq_vec=rep(NA,length(BootSampleIndex)) +mypcex3=rep(NA,length(BootSampleIndex)) +for (j in 1:length(unique(BootSampleIndex))){ +mylength=length(freq_vec[which(BootSampleIndex==unique(BootSampleIndex)[j])]) +freq_vec[which(BootSampleIndex==unique(BootSampleIndex)[j])]=1:mylength +mypcex3[which(BootSampleIndex==unique(BootSampleIndex)[j])]=mypcex2[1:mylength] +} + +mypcex4=mypcex+mypcex3 +mypcex4 + +BootSampleIndexMatrix=originalsampleIndex[BootSampleIndex, ] +BootSample=population[BootSampleIndexMatrix] +points(BootSampleIndexMatrix, col = p.col[2], cex = mypcex4) + +# xb[b,]=BootSample + + +guessInput <- reactive({ + isolate({ + values$tx_star_b <- c(values$tx_star_b, as.numeric(mean(population[BootSampleIndexMatrix]))) + return(values$tx_star_b) + }) + }) + +guessInput2 <- reactive({ + isolate({ + values2$sigmax_star_b <- c(values2$sigmax_star_b, as.numeric(sqrt(var(population[BootSampleIndexMatrix])*(n-1)/n))) + return(values2$sigmax_star_b) + }) + }) + +guessInput3 <- reactive({ + isolate({ + values3$tx_star_bMINUStheta_fn <- c(values3$tx_star_bMINUStheta_fn, as.numeric(mean(population[BootSampleIndexMatrix])-mean(population[originalsampleIndex]))) + return(values3$tx_star_bMINUStheta_fn) + }) + }) + +guessInput4 <- reactive({ + isolate({ + values4$u_star_b <- c(values4$u_star_b, as.numeric(sqrt(n)*(mean(population[BootSampleIndexMatrix])-mean(population[originalsampleIndex]))/sqrt(var(population[BootSampleIndexMatrix])*(n-1)/n))) + return(values4$u_star_b) + }) + }) + + + # cat(file=stderr(), " ", as.numeric(guessInput()), " ", "\n") + # cat(file=stderr(), " ", as.numeric(guessInput2()), " ", "\n") + # cat(file=stderr(), " ", as.numeric(guessInput3()), " ", "\n") + # cat(file=stderr(), " ", as.numeric(guessInput4()), " ", "\n") + +tx_star_b <- c(tx_star_b,mean(population[BootSampleIndexMatrix])) +sigmax_star_b=c(sigmax_star_b,sqrt(var(population[BootSampleIndexMatrix])*(n-1)/n)) +tx_star_bMINUStheta_fn=c(tx_star_bMINUStheta_fn,mean(population[BootSampleIndexMatrix])-mean(population[originalsampleIndex])) +u_star_b=c(u_star_b,sqrt(n)*(mean(population[BootSampleIndexMatrix])-mean(population[originalsampleIndex]))/sqrt(var(population[BootSampleIndexMatrix])*(n-1)/n)) + +mycol2<-c(rep("gray",length(as.numeric(guessInput()))-1),"red") + +plot(0,type="n",main="",ylab='',yaxt='n',xaxt='n',bty='n') +plot(0,type="n",main="",ylab='',yaxt='n',xaxt='n',bty='n') + + +par(mar=c(2, 3, 0, 1),mgp=c(5,1,0)) +hist(as.numeric(guessInput()),freq = TRUE,xlim=c(0,mu+100),main="",ylab='',xaxt='n',bty='n',cex.axis=2) +points(as.numeric(guessInput()),rep(0,length(as.numeric(guessInput()))),cex.lab=2.2,xlim=c(0,mu+100), +xlab='',ylab='',pch=17,col=mycol2,cex=1.4, xaxt='n',bty='n') +axis(1, seq(0,mu+100,length=5), cex.axis=2, tck=-.01) +mtext("Frequency",at=-100,side=1,line =-22,cex=2,las=2) + + +points(50,0,cex.lab=2.2,xlim=c(1,3),ylim=c(0,150),xlab='',ylab='',pch=15,col="blue",cex=1.4, yaxt='n', xaxt='n',bty='n') +points(mean(population[originalsampleIndex]),0,cex.lab=2.2,xlim=c(1,3),xlab='',ylab='',pch=18,col="magenta",cex=1.8, yaxt='n', xaxt='n',bty='n') +points(mean(as.numeric(guessInput())),0,cex.lab=2.2,xlim=c(1,3),xlab='',ylab='',pch=19,col="green",cex=1.4, yaxt='n', xaxt='n',bty='n') + +CIboot=c(orig.mean-quantile(as.numeric(guessInput3()),1-0.05/2),orig.mean-quantile(as.numeric(guessInput3()),0.05/2)) + + + +plot(0,type="n",main="",ylab='',yaxt='n',xaxt='n',bty='n') +legend("center",legend=c( +as.expression(bquote(mu*"="*.(50))), +as.expression(bquote(theta*"("*F[n]*")"*"="*.(format(mean(population[originalsampleIndex]), digits = 0)))), +as.expression(bquote(T({{chi^"*("}^.(format(b, digits = 0))}^")")*"="*.(format(mean(population[BootSampleIndexMatrix]), digits = 0)))), +as.expression(bquote(T(chi^"*")*"="*.(format(mean(as.numeric(guessInput())), digits = 0)))), +paste("[",round(CIboot[1],2)," ; ",round(CIboot[2],2), "]",sep="") +),col=c("blue","magenta","red","green"),pch=c(15,18,17,19,NA),bty="n",cex=1.8) + + + +}) + + +observeEvent( input$refresh, { + updateSliderInput(session, "bb", label="", value=1, min=1, max=input$B, step=1) + set.seed(2018+as.numeric(input$submit_model)) + values$tx_star_b=c() +values2$sigmax_star_b=c() +values3$tx_star_bMINUStheta_fn=c() +values4$u_star_b=c() +values5$xbar=c() + }) + +observeEvent( input$submit_model, { + values$tx_star_b=c() +values2$sigmax_star_b=c() +values3$tx_star_bMINUStheta_fn=c() +values4$u_star_b=c() +values5$xbar=c() + }) + +observeEvent( input$n, { + values$tx_star_b=c() +values2$sigmax_star_b=c() +values3$tx_star_bMINUStheta_fn=c() +values4$u_star_b=c() +values5$xbar=c() + }) + + +observeEvent( input$B, { + updateSliderInput(session, "bb", label="", value=1, min=1, max=input$B, step=1) + values$tx_star_b=c() +values2$sigmax_star_b=c() +values3$tx_star_bMINUStheta_fn=c() +values4$u_star_b=c() +values5$xbar=c() + }) +} + + + +