Generalized Akaike Information Criterion (GAIC)

Esta aplicación hace la comparación de modelos usando el Generalized Akaike Information Criterion propuesto por Akaike (1983) dado por la expresión

\[GAIC = −2 \times logLik + k \times df,\]

donde $l$ corresponde al valor de log-verosimilitud estimada para el modelo, $df$ es el número de parámetros estimados del modelo y el valor de $k$ es una penalización que el analista asigna por el exceso de parámetros. Cuando $k=2$ el $GAIC$ coincide con el $AIC$, cuando el valor de penalidad es $k = log(n)$ se obtine el Schwarz Bayesian criterion ($SBC$ or $BIC$) propuesto por Schwarz (1978) en el que $n$ es el número de observaciones del modelo.

Distribuciones disponibles

Las distribuciones que se pueden elegir en la aplicación son las que están disponibles en el paquete gamlss de R, a continuación el listado de familias y las distribuciones que pertenecen a cada una de ellas.

Reales positivos: “EXP”,“GA”,“IG”,“LNO”, “WEI3”, “BCCGo”, “exGAUS”, “GG”, “GIG”, “BCTo”, “BCPEo”.
Reales: “GU”, “RG” ,“LO”, “NET”, “TF”, “PE”, “SN1”, “SN2”, “SHASH”, “EGB2”, “JSU”, “SEP1”, “SEP2”, “SEP3”, “SEP4”,“ST1”, “ST2”, “ST3”, “ST4”, “ST5”, “GT”.
Reales entre 0 y 1: “BE”, “BEINF”, “BEINF0”, “BEINF1”, “BEOI”, “BEZI”, “GB1”.
Conteos: “PO”, “LG”, “NBI”, “NBII”, “PIG”, “DEL”, “SI”, “ZIP”, “ZAP”, “ZALG”, “ZANBI”, “ZIP2”, “ZIPIG”.
Binomiales: “BI”, “BB”, “ZIBI”, “ZIBB”, “ZABI”, “ZABB”.

Para obtener información sobre cada una de las anteriores distribuciones se recomienda consultar la tabla 4.2 y Apéndice A del Manual de gamlss.

library(shiny) require(gamlss) shinyServer(function(input, output, session){ # Function to plot the four histograms four.hist <- function(k, f, p) { par(cex.main=0.95) inFile <- input$file1 if(is.null(inFile)) dt <- cars else dt <- read.csv(inFile$datapath, header=input$header, sep=input$sep) mod <- fitDist(dt[, p], type=f, k=k) par(mfrow=c(2, 2), bg='gray98') for(i in 1:4){ denst <- density(dt[, p]) res <- histDist(dt[, p], family=names(mod$fits)[i], main='', ylab='Density', xlab=p, las=1, line.wd=3, line.ty=1, line.col='dodgerblue2', ylim=c(0, (2 * max(denst$y)))) gaic <- round(-2 * logLik(res) + k * length(res$parameters), 2) title(main=paste(i, ')', names(mod$fits)[i], 'distribution with GAIC=', gaic), col.main='blue4') param <- c('mu', 'sigma', 'nu', 'tau') np <- length(res$parameters) fun1 <- function(x) eval(parse(text=x)) hat.param <- sapply(as.list(paste('res$', param[1:np], sep='')), fun1) hat.param <- round(hat.param, digits=2) txt <- paste('hat(', param[1:np], ')==', hat.param, sep='') txt <- paste(txt, collapse=', ') legend('topright', bty='n', legend=eval(parse(text=paste('expression(', txt, ')')))) } } four.hist.qqplot <- function(k, f, p) { par(cex.main=0.95) inFile <- input$file1 if(is.null(inFile)) dt <- cars else dt <- read.csv(inFile$datapath, header=input$header, sep=input$sep) mod <- fitDist(dt[, p], type=f, k=k) par(mfrow=c(4, 2), bg='gray98') for(i in 1:4){ denst <- density(dt[, p]) res <- histDist(dt[, p], family=names(mod$fits)[i], main='', ylab='Density', xlab=p, las=1, line.wd=3, line.ty=1, line.col='dodgerblue2', ylim=c(0, (2 * max(denst$y)))) gaic <- round(-2 * logLik(res) + k * length(res$parameters), 2) title(main=paste(i, ')', names(mod$fits)[i], 'distribution with GAIC=', gaic), col.main='blue4') param <- c('mu', 'sigma', 'nu', 'tau') np <- length(res$parameters) fun1 <- function(x) eval(parse(text=x)) hat.param <- sapply(as.list(paste('res$', param[1:np], sep='')), fun1) hat.param <- round(hat.param, digits=2) txt <- paste('hat(', param[1:np], ')==', hat.param, sep='') txt <- paste(txt, collapse=', ') legend('topright', bty='n', legend=eval(parse(text=paste('expression(', txt, ')')))) #qqnorm(runif(10)) n <- dim(dt)[1] a <- paste('q', names(mod$fits)[i], sep='') a <- paste(a, '(ppoints(n),', sep='') a <- paste(a, paste(hat.param, collapse=', '), ')') par(bty='n', col.axis='blue4', col.lab='blue4') qqplot(y=dt[, p], x=eval(parse(text=a)), xlab=a, ylab=p, col='blue4') axis(side=1, col='blue4') axis(side=2, col='blue4') } } observe({ inFile <- input$file1 if(is.null(inFile)) dt <- cars else dt <- read.csv(inFile$datapath, header=input$header, sep=input$sep) updateSelectInput(session, "response", choices = names(dt)) updateSelectInput(session, "familia") }) output$data_table <- renderTable({ inFile <- input$file1 if(is.null(inFile)) dt <- cars else dt <- read.csv(inFile$datapath, header=input$header, sep=input$sep) dt }) output$distPlot1 <- renderPlot({ four.hist(input$k, input$familia, input$response) }) output$distPlot2 <- renderPlot({ four.hist.qqplot(input$k, input$familia, input$response) }) output$downplot1 <- downloadHandler( filename = function() { paste("four_dist", input$type_plot, sep=".") }, # content is a function with argument file. content writes the plot to the device content = function(file) { if(input$type_plot == "png") png(file) # open the png device else pdf(file) # open the pdf device four.hist(input$k, input$familia, input$response) # draw the plot dev.off() # turn the device off } ) output$downplot2 <- downloadHandler( filename = function() { paste("four_dist", input$type_plot, sep=".") }, # content is a function with argument file. content writes the plot to the device content = function(file) { if(input$type_plot == "png") png(file) # open the png device else pdf(file) # open the pdf device four.hist.qqplot(input$k, input$familia, input$response) # draw the plot dev.off() # turn the device off } ) output$markdown <- renderUI({ HTML(markdown::markdownToHTML(knit('Teoria.Rmd', quiet = TRUE))) }) })

library(shiny) library(knitr) shinyUI(pageWithSidebar( headerPanel("Exploring the distribution that fits better for a variable"), # Here the sidebarpanel sidebarPanel( h6("This app can be used to identify the best four distributions that fit better for a single variable. The user can provide the dataset using the buttons below."), h6("As an example the app uses the speed variable from the cars dataset."), h4("To use the app follow the next steps:"), h6("1. Upload the file,"), h6("2. Select the separator character,"), h6("3. Select the variable."), fileInput(inputId="file1", label="Use the next button to load your file.", accept = c( "text/csv", "text/comma-separated-values", "text/tab-separated-values", "text/plain", ".csv", ".tsv" )), checkboxInput("header", label="Does your file contain variable names in the first line?", value=TRUE), selectInput(inputId="sep", label="What is the separator character in your file?", choices=list("tab"="\t", ","=",", ";"=";", "space"=" "), selected=";"), selectInput("response", label="Select the variable.", choices=""), selectInput("familia", "Select the family for the response variable:", choice=list("Real plus" = "realplus", "Real" = "realAll", "Real between entre 0 y 1" = "real0to1", "Counts" = "counts", "Binomial" = "binom")), sliderInput("k", "Choose the number \$ k \$ to penalize for the excess number of parameters in \$GAIC=-2 \\times logLik + k \\times df\$. The value \$ df \$ corresponds to the number of distribution parameters.", min=2, max=20, value=2, step=1, animate=TRUE), radioButtons(inputId = "type_plot", label = "Select the file type to download the plot.", choices = list("pdf", "png")), br(), p("This app was created by Semillero de R at Universidad Nacional de Colombia:"), tags$a(href="https://srunal.github.io", "https://srunal.github.io") ), # Here the panels mainPanel( tabsetPanel(type ="pills", tabPanel("Selected distributions I", h4("A continuación se presenta el histograma para la variable de interés con la curva de densidad para cada distribución."), plotOutput("distPlot1", width="700px", height="600px"), downloadButton(outputId = "downplot1", label = "Download the plot") ), tabPanel("Selected distributions II", h4("A continuación se presenta el histograma para la variable de interés con la curva de densidad para cada distribución y el qqplot asociado."), plotOutput("distPlot2", width="700px", height="800px"), downloadButton(outputId = "downplot2", label = "Download the plot") ), tabPanel("Data", uiOutput("data_table")), tabPanel("Goodness of fit test"), tabPanel("Theory", includeHTML("include.html")) ) ) ) )