diff --git a/README.md b/README.md
index 2b60a4d..cbef854 100644
--- a/README.md
+++ b/README.md
@@ -1,24 +1,58 @@
# Maximum Entropy with Hidden Structure in R
-R support for various batch and online computations in Maximum Entropy Harmonic Grammar. Used to find MaxEnt solutions for learning problems (with or without hidden structure), generate distributions over forms in different theories, and perform associated learning tasks.
-
-## Installation
-
-1. Install R ()
-2. Clone repository (e.g. `git clone git@github.com:rstaubs/maxent-hidden-structure.git` )
-3. Open repository in R environment of choice (e.g., RStudio )
-
-## Examples
-
-Examples are shown in the `examples` directory:
-
-1. `base` - Example driver file and data files for parallel MaxEnt (with and without hidden structure)
-2. `serial` Example driver file and data files for serial MaxEnt
-
-## Notes on Optimization
-
-The gradients used in this optimization are derived here:
-
-## Citation
-
-Staubs, Robert. 2011. Harmonic Grammar in R (hgR). Software package. Amherst, MA: University of Massachusetts Amherst.
+## Harmonic Grammar in R
+This repository provides R support for various batch and online computations in Maximum Entropy Harmonic Grammar. The repository can be used to find MaxEnt solutions for learning problems (with or without hidden structure), generate distributions over forms in different theories, and perform associated learning tasks.
+
+### Installation
+1. Install R: [https://cran.rstudio.com/](https://cran.rstudio.com/)
+2. Clone the repository:
+ ```sh
+ git clone https://github.com/rstaubs/maxent-hidden-structure.git
+3. Open the repository in your preferred R environment (e.g., RStudio).
+
+### Examples
+Example files are provided in the `examples` directory:
+
+- **base**: Example driver file and data files for parallel MaxEnt (with and without hidden structure).
+- **serial**: Example driver file and data files for serial MaxEnt.
+
+### Notes on Optimization
+The gradients used in this optimization are derived from the notes available [here](https://websites.umass.edu/hgr/files/2017/07/klnotes.pdf).
+
+## Shiny App Implementation
+The Shiny app provides a user-friendly interface for implementing the HGR model. The app offers additional functionalities, such as global normalization, which is not available in the standard HGR implementation. Through the Shiny app, users can quickly upload data, generate grammars, adjust constraint weights, and view results in real-time, making the process more efficient and accessible.
+
+### Accessing the App
+You can access the Shiny app online at [MaxEnt with Hidden Structure in R](https://alingwist.shinyapps.io).
+
+### Running the App Locally
+1. Install the necessary R packages:
+ ```r
+ install.packages(c("shiny", "DT", "shinyjs", "dplyr", "htmltools"))
+2. Run the app using the following command:
+ ```r
+ shiny::runApp("shiny_app")
+
+### Sample Input Files:
+The Shiny app includes three sample input files, which can be found in the `sample_input_files` subdirectory. These files demonstrate different scenarios for using the MaxEnt model, including cases with or without hidden structures and cases of having output candidates with raw frequencies instead of normalized probabilities. these sample files can be uploaded directly into the Shiny app to see how the app handles different types of data.
+
+### User Guide
+#### Input Data Format
+- The input file should be a CSV or TXT file with the following columns:
+ - `input`: Input forms
+ - `output`: Output forms
+ - `probability`: Observed probabilities or raw frequencies
+ - Additional columns representing various constraints.
+
+#### Using the Shiny App
+1. **Upload your data file**: Use the "Choose Input File" button to upload your CSV or TXT file.
+2. **Select input type**: Choose between "Probabilities" and "Raw Frequencies".
+3. **Select normalization type**: Choose between "Within Tableau" and "Global Normalization".
+4. **Select prior type**: Choose between "L2" and "L1".
+5. **Generate Grammar**: Click "Generate Grammar" to produce the grammar.
+6. **Edit constraint weights**: If needed, adjust weights in the dynamic UI and click "Update Grammar" to recalculate probabilities and errors.
+7. **Download outputs**: Use the "Download Output" button to save the generated tableau.
+
+## Citations
+- Staubs, Robert. 2011. Harmonic Grammar in R (hgR). Software package. Amherst, MA: University of Massachusetts Amherst. https://github.com/rstaubs/maxent-hidden-structure
+- Nirheche, Ali. 2024. Shiny App for MaxEnt with Hidden Structure. Shiny Application. Amherst, MA: University of Massachusetts Amherst. https://alingwist.shinyapps.io
\ No newline at end of file
diff --git a/shiny_app/app.R b/shiny_app/app.R
new file mode 100644
index 0000000..b92f588
--- /dev/null
+++ b/shiny_app/app.R
@@ -0,0 +1,383 @@
+library(shiny)
+library(DT)
+library(shinyjs)
+library(dplyr)
+library(htmltools)
+
+# Define the UI
+ui <- fluidPage(
+ useShinyjs(), # Enables shinyjs functionalities
+ titlePanel("MaxEnt with Hidden Structure in R"),
+ sidebarLayout(
+ sidebarPanel(
+ # Input for uploading a file
+ fileInput("fileInput", "Choose Input File", accept = c("text/csv", "text/comma-separated-values,text/plain", ".csv")),
+ # Download link for sample input files
+ tags$div(style = "margin-top: -35px;", # Adjust the margin to control spacing
+ downloadLink("downloadSampleFiles", "Sample Input Files")
+ ),
+ tags$br(),
+ # Radio buttons for selecting input type
+ radioButtons("frequencyType", "Select Input Type:", choices = c("Probabilities" = "prob", "Raw Frequencies" = "raw"), selected = "prob"),
+ # Radio buttons for selecting normalization type
+ radioButtons("normalizationType", "Select Normalization:", choices = c("Within Tableau" = "within", "Global Normalization" = "global"), selected = "within"),
+ # Radio buttons for selecting prior type
+ radioButtons("priorType", "Select Prior:", choices = c("L2" = "L2", "L1" = "L1"), selected = "L2"),
+ # Button to generate grammar
+ actionButton("generateBtn", "Generate Grammar"),
+ # Button to download output
+ downloadButton("downloadOutput", "Download Output", style = "display: none;"),
+ # Placeholder for dynamic UI elements for updating weights
+ uiOutput("updateWeightsUI"),
+ # Button to load more constraints (initially hidden)
+ actionButton("loadMoreBtn", "Load More Constraints", style = "display: none;"),
+ # Button to update grammar (initially hidden)
+ actionButton("updateGrammarBtn", "Update Grammar", style = "display: none;")
+ ),
+ mainPanel(
+ # Output table to display the grammar
+ DTOutput("grammarTable"),
+ # Text output for additional information
+ textOutput("infoText")
+ )
+ )
+)
+
+# Define the server logic
+server <- function(input, output, session) {
+ # Reactive values to store various data
+ grammarData <- reactiveVal(NULL)
+ unroundedData <- reactiveVal(NULL)
+ constraints <- reactiveVal(NULL)
+ weights <- reactiveVal(NULL)
+ stored_output_df <- reactiveVal(NULL)
+ # Counter to keep track of the number of loaded constraints
+ loadedConstraints <- reactiveVal(10)
+
+ # Observe event when 'Generate Grammar' button is clicked
+ observeEvent(input$generateBtn, {
+ inFile <- input$fileInput
+ if (is.null(inFile)) {
+ showModal(modalDialog(
+ title = "Error",
+ "No file was uploaded. Please upload a file.",
+ easyClose = TRUE,
+ footer = NULL
+ ))
+ return(NULL)
+ }
+
+ fileExt <- tools::file_ext(inFile$name)
+ if (!(fileExt %in% c("txt", "csv"))) {
+ showModal(modalDialog(
+ title = "Error",
+ "The uploaded input file is not a txt or csv file. Please upload a txt or csv file.",
+ easyClose = TRUE,
+ footer = NULL
+ ))
+ return(NULL)
+ }
+
+ source("solver.R") # Load external R script for solving
+ tab <- read.table(inFile$datapath, header = TRUE, sep = "\t", check.names = FALSE)
+ hidden_structure <- "hidden" %in% colnames(tab) # Check if hidden structure is included
+
+ # Offset to determine where constraints start in the data
+ offset <- if (hidden_structure) 4 else 3
+ w <- rep(0, ncol(tab) - offset) # Initial weights for constraints
+ tab[, -(1:offset)] <- -tab[, -(1:offset)] # Negate constraints for MaxEnt calculation
+
+ # Calculate observed probabilities based on user input
+ if (input$frequencyType == "raw" && input$normalizationType == "global") {
+ observed_probability <- tab[[offset]] / sum(tab[[offset]])
+ } else if (input$frequencyType == "raw" && input$normalizationType == "within") {
+ observed_probability <- tab %>%
+ group_by(input) %>%
+ mutate(observed_probability = probability / sum(probability)) %>%
+ ungroup() %>%
+ .$observed_probability
+ } else {
+ observed_probability <- tab[[offset]]
+ if (any(observed_probability > 1)) {
+ observed_probability <- observed_probability / 100
+ }
+ }
+
+ # Create a new Tableau object
+ o.tab <- new("Tableau", data = tab, hidden = hidden_structure, theory = "maxent")
+
+ # Solve for optimal weights
+ if (input$normalizationType == "global") {
+ opt <- solve_global_norm(o.tab, w, method = "L-BFGS-B", reg = input$priorType, lower.bound = 0, upper.bound = Inf, observed_probs = observed_probability)
+ harmony <- apply(tab[, (offset + 1):ncol(tab)], 1, function(x) sum(x * opt$par))
+ expected_probability <- opt$final_probabilities
+ } else {
+ opt <- solve(o.tab, w, categorical = hidden_structure, method = "L-BFGS-B", reg = input$priorType, var = 1000000, mf = FALSE)
+ harmony <- apply(tab[, (offset + 1):ncol(tab)], 1, function(x) sum(x * opt$par))
+ }
+
+ # Round the weights for display
+ rounded_weights <- round(opt$par, 2)
+ constraints(colnames(tab)[-(1:offset)]) # Store constraint names
+ weights(rounded_weights) # Store rounded weights
+
+ # Prepare the output data frame
+ output_df <- tab[, 1:offset]
+
+ # Adjust the data frame based on frequency type
+ if (input$frequencyType == "prob") {
+ if (hidden_structure) {
+ output_df <- output_df[, -4]
+ } else {
+ output_df <- output_df[, -3]
+ }
+ }
+
+ output_df$observed_probability <- round(observed_probability, 5)
+ output_df <- cbind(output_df, tab[, (offset + 1):ncol(tab)])
+
+ output_df$harmony <- round(harmony, 5)
+ if (input$normalizationType == "global") {
+ output_df$expected_probability <- round(expected_probability, 5)
+ } else {
+ if (hidden_structure) {
+ distribution_df <- hidden.distribution(o.tab, opt$par)
+ } else {
+ distribution_df <- distribution(o.tab, opt$par)
+ }
+ output_df$expected_probability <- round(distribution_df$probability, 5)
+ }
+ expected_probability <- output_df$expected_probability
+
+ # Calculate error metrics
+ output_df$AbsoluteError <- abs(observed_probability - expected_probability)
+ output_df <- output_df %>%
+ group_by(input) %>%
+ mutate(MAE_per_tableau = mean(abs(observed_probability - expected_probability))) %>%
+ ungroup()
+ MAE_overall <- mean(abs(observed_probability - expected_probability))
+ output_df$MAE_overall <- MAE_overall
+
+ # Round error metrics for display
+ output_df$AbsoluteError <- round(output_df$AbsoluteError, 5)
+ output_df$MAE_per_tableau <- round(output_df$MAE_per_tableau, 5)
+ output_df$MAE_overall <- round(output_df$MAE_overall, 5)
+
+ # Store unrounded data for download
+ unrounded_output_df <- output_df
+ unrounded_output_df$observed_probability <- observed_probability
+ unrounded_output_df$harmony <- harmony
+ unrounded_output_df$expected_probability <- expected_probability
+ unrounded_output_df$AbsoluteError <- abs(observed_probability - expected_probability)
+ unrounded_output_df <- unrounded_output_df %>%
+ group_by(input) %>%
+ mutate(MAE_per_tableau = mean(abs(observed_probability - expected_probability))) %>%
+ ungroup()
+ unrounded_output_df$MAE_overall <- MAE_overall
+
+ # Prepare weights row for display
+ observed_prob_index <- which(names(output_df) == "observed_probability")
+ if (input$frequencyType == "raw") {
+ weight_row <- c(rep(NA, (offset+1)), rounded_weights)
+ colnames(output_df)[observed_prob_index - 1] <- "raw_frequency"
+ colnames(unrounded_output_df)[observed_prob_index - 1] <- "raw_frequency"
+ } else {
+ weight_row <- c(rep(NA, (offset)), rounded_weights)
+ }
+ weight_row <- c(weight_row, rep(NA, ncol(output_df) - length(weight_row)))
+ weight_row <- as.data.frame(t(weight_row))
+ colnames(weight_row) <- colnames(output_df)
+
+ output_df <- rbind(weight_row, output_df)
+
+ weight_row_unrounded <- c(rep(NA, (offset + 1)), opt$par)
+ weight_row_unrounded <- c(weight_row_unrounded, rep(NA, ncol(unrounded_output_df) - length(weight_row_unrounded)))
+ weight_row_unrounded <- as.data.frame(t(weight_row_unrounded))
+ colnames(weight_row_unrounded) <- colnames(unrounded_output_df)
+
+ unrounded_output_df <- rbind(weight_row_unrounded, unrounded_output_df)
+
+ # Update reactive values
+ grammarData(output_df)
+ unroundedData(unrounded_output_df)
+ stored_output_df(output_df)
+
+ shinyjs::show("downloadOutput")
+
+ # Initially load the first 10 constraints
+ loadedConstraints(10)
+
+ # Render dynamic UI for updating weights
+ updateUI <- renderDynamicConstraints()
+ output$updateWeightsUI <- renderUI(updateUI)
+ })
+
+ # Observe event when 'Load More Constraints' button is clicked
+ observeEvent(input$loadMoreBtn, {
+ currentCount <- loadedConstraints()
+ loadedConstraints(currentCount + 10)
+
+ updateUI <- renderDynamicConstraints()
+ output$updateWeightsUI <- renderUI(updateUI)
+ })
+
+ # Reactive function to render dynamic constraints
+ renderDynamicConstraints <- reactive({
+ constraints_to_load <- head(constraints(), loadedConstraints())
+
+ ui_elements <- lapply(seq_along(constraints_to_load), function(i) {
+ constraint_name <- constraints_to_load[i]
+ fluidRow(
+ column(12, numericInput(paste0("weight_", i), constraint_name, value = weights()[i]))
+ )
+ })
+
+ ui_elements <- c(
+ list(h3("Update Weights")), # Add the heading here
+ ui_elements,
+ list(actionButton("updateGrammarBtn", "Update Grammar"))
+ )
+
+ # Conditionally include the "Load More Constraints" button
+ if (length(constraints()) > loadedConstraints()) {
+ ui_elements <- c(ui_elements, list(actionButton("loadMoreBtn", "Load More Constraints")))
+ }
+
+ tagList(ui_elements)
+ })
+
+ # Observe event when 'Update Grammar' button is clicked
+ observeEvent(input$updateGrammarBtn, {
+ output_df <- stored_output_df()
+ if (is.null(output_df)) {
+ showModal(modalDialog(
+ title = "Error",
+ "No grammar has been generated. Please generate grammar first.",
+ easyClose = TRUE,
+ footer = NULL
+ ))
+ return(NULL)
+ }
+
+ hidden_structure <- "hidden" %in% colnames(output_df)
+ offset <- if (hidden_structure) 5 else 4
+
+ new_weights <- sapply(1:length(weights()), function(i) input[[paste0("weight_", i)]])
+
+ # Update weights and recalculate harmony
+ if (hidden_structure && input$frequencyType == "raw") {
+ output_df[1, (offset + 1):(offset + length(new_weights))] <- round(new_weights, 2)
+ tab <- output_df[2:nrow(output_df), c(1:(offset), (offset + 1):(offset + length(new_weights)))]
+ harmony <- rowSums(tab[, (offset + 1):ncol(tab)] * matrix(new_weights, nrow = nrow(tab), ncol = length(new_weights), byrow = TRUE))
+ } else if (!hidden_structure && input$frequencyType == "raw") {
+ output_df[1, (offset + 1):(offset + length(new_weights))] <- round(new_weights, 2)
+ tab <- output_df[2:nrow(output_df), 1:(offset + length(new_weights))]
+ harmony <- rowSums(tab[, (offset + 1):ncol(tab)] * matrix(new_weights, nrow = nrow(tab), ncol = length(new_weights), byrow = TRUE))
+ } else if (hidden_structure || input$frequencyType == "prob") {
+ output_df[1, (offset):(offset + length(new_weights) - 1)] <- round(new_weights, 2)
+ tab <- output_df[2:nrow(output_df), c(1:(offset - 1), (offset):(offset + length(new_weights) - 1))]
+ harmony <- rowSums(tab[, (offset):(offset + length(new_weights) - 1)] * matrix(new_weights, nrow = nrow(tab), ncol = length(new_weights), byrow = TRUE))
+ }
+
+ # Calculate expected probabilities
+ if (input$normalizationType == "global") {
+ scores <- exp(c(harmony))
+ total_score <- sum(scores)
+ expected_probability <- scores / total_score
+ } else if (input$normalizationType == "within") {
+ expected_probability <- numeric(length = length(harmony))
+ unique_inputs <- unique(output_df$input[2:nrow(output_df)])
+ for (input_val in unique_inputs) {
+ input_indices <- which(output_df$input[2:nrow(output_df)] == input_val) + 1
+ scores <- exp(c(harmony[input_indices - 1]))
+ total_score <- sum(scores)
+ expected_probability[input_indices - 1] <- scores / total_score
+ }
+ }
+
+ # Update the output data frame with new values
+ output_df$harmony[2:nrow(output_df)] <- round(harmony, 5)
+ output_df$expected_probability[2:nrow(output_df)] <- round(expected_probability, 5)
+ observed_probability <- output_df$observed_probability[2:nrow(output_df)]
+ new_AE <- abs(observed_probability - expected_probability)
+
+ # Calculate MAE per tableau
+ mae_per_tableau_df <- output_df[2:nrow(output_df),] %>%
+ group_by(input) %>%
+ summarise(MAE_per_tableau = mean(abs(observed_probability - expected_probability))) %>%
+ ungroup()
+
+ for (i in 2:nrow(output_df)) {
+ output_df$MAE_per_tableau[i] <- mae_per_tableau_df$MAE_per_tableau[mae_per_tableau_df$input == output_df$input[i]]
+ }
+
+ new_MAE_overall <- mean(new_AE)
+
+ for (i in 2:nrow(output_df)) {
+ output_df$AbsoluteError[i] <- round(new_AE[i - 1], 5)
+ }
+
+ for (i in 2:nrow(output_df)) {
+ output_df$MAE_overall[i] <- round(new_MAE_overall, 5)
+ }
+
+ # Update reactive values
+ grammarData(output_df)
+ stored_output_df(output_df)
+
+ unrounded_output_df <- unroundedData()
+
+ tryCatch({
+ unrounded_output_df$harmony[2:nrow(unrounded_output_df)] <- harmony
+ unrounded_output_df$expected_probability[2:nrow(unrounded_output_df)] <- expected_probability
+ unrounded_output_df$AbsoluteError[2:nrow(unrounded_output_df)] <- new_AE
+
+ for (i in 2:nrow(unrounded_output_df)) {
+ unrounded_output_df$MAE_per_tableau[i] <- mae_per_tableau_df$MAE_per_tableau[mae_per_tableau_df$input == output_df$input[i]]
+ }
+
+ unrounded_output_df$MAE_overall[2:nrow(unrounded_output_df)] <- new_MAE_overall
+ unroundedData(unrounded_output_df)
+ }, error = function(e) {
+ print(paste("Error occurred:", e$message))
+ })
+ })
+
+ # Render the grammar table
+ output$grammarTable <- renderDT({
+ data <- grammarData()
+ if (is.null(data) || nrow(data) == 0) {
+ return(NULL)
+ }
+ datatable(data, editable = TRUE, options = list(pageLength = 25), escape = FALSE) %>%
+ formatStyle(names(data), 'white-space' = 'nowrap') %>%
+ formatStyle(columns = colnames(data), escape = FALSE, formatType = "html")
+ }, server = TRUE)
+
+ # Download handler for the output tableau
+ output$downloadOutput <- downloadHandler(
+ filename = function() {
+ paste("output-tableau-", Sys.Date(), ".txt", sep = "")
+ },
+ content = function(file) {
+ data <- stored_output_df()
+ if (is.null(data)) return()
+ write.table(data, file, sep = "\t", row.names = FALSE, col.names = TRUE, quote = FALSE, na = "NA", fileEncoding = "UTF-8")
+ },
+ contentType = "text/plain"
+ )
+
+ # Download handler for the sample input files
+ output$downloadSampleFiles <- downloadHandler(
+ filename = function() {
+ "sample_input_files.zip"
+ },
+ content = function(file) {
+ file.copy("sample_input_files.zip", file)
+ },
+ contentType = "application/zip"
+ )
+}
+
+# Run the app
+shinyApp(ui = ui, server = server)
\ No newline at end of file
diff --git a/shiny_app/sample_input_files/input_with_hidden_structure.txt b/shiny_app/sample_input_files/input_with_hidden_structure.txt
new file mode 100644
index 0000000..eba40ea
--- /dev/null
+++ b/shiny_app/sample_input_files/input_with_hidden_structure.txt
@@ -0,0 +1,10 @@
+input output hidden probability Iamb Trochee *StressHi
+kata KAta (KAta) 1 1 0 0
+kata kaTA (kaTA) 0 0 1 0
+kata KAta (KA)ta 1 0 1 0
+kata kaTA ka(TA) 0 0 1 0
+pika PIka (PIka) 0 1 0 1
+pika piKA (piKA) 1 0 1 0
+pika PIka (PI)ka 0 0 1 1
+pika piKA pi(KA) 1 0 1 0
+
diff --git a/shiny_app/sample_input_files/input_with_raw_frequencies.txt b/shiny_app/sample_input_files/input_with_raw_frequencies.txt
new file mode 100644
index 0000000..bcb7808
--- /dev/null
+++ b/shiny_app/sample_input_files/input_with_raw_frequencies.txt
@@ -0,0 +1,5 @@
+input output probability Iamb Trochee *StressHi
+kata káta 35 1 0 0
+kata katá 7 0 1 0
+pika píka 40 1 0 1
+pika piká 10 0 1 0
diff --git a/shiny_app/sample_input_files/input_without_hidden_structure.txt b/shiny_app/sample_input_files/input_without_hidden_structure.txt
new file mode 100644
index 0000000..21a73d3
--- /dev/null
+++ b/shiny_app/sample_input_files/input_without_hidden_structure.txt
@@ -0,0 +1,5 @@
+input output probability Iamb Trochee *StressHi
+kata káta 1 1 0 0
+kata katá 0 0 1 0
+pika píka 1 1 0 1
+pika piká 0 0 1 0
diff --git a/shiny_app/solver.R b/shiny_app/solver.R
new file mode 100644
index 0000000..acef8d3
--- /dev/null
+++ b/shiny_app/solver.R
@@ -0,0 +1,390 @@
+# Provides objective functions and solving for MaxEnt
+# Requires: tableau
+
+# Authors: David Smith, Robert Staubs
+# Contact: rstaubs@linguist.umass.edu)
+# Last update: 3/11/2013
+
+source("tableau.R")
+
+# Solve for weights for data in tableau
+#
+# tableau Tableau object containing problem
+# weights starting weight vector
+# categorical is the target distribution categorical?
+# method optimization method to use (see ?optim)
+# reg regularization to use (L1 or L2)
+# mean mean for regularization
+# var variance for L2 regularization
+# alpha parameter for L1 regularization
+# mf is a Markedness > Faithfulness used
+# markedness boolean vector identifying Markedness constraints
+# faithfulness boolean vector identifying Faithfulness constraints
+# mf.weight weight of M > F term
+# gradient should the provided gradient be used? (FALSE generally means a numerical estimate is used)
+# lower.bound vector of lower bounds on constraint weights
+# upper.bound vector of upper bounds on constraint weights
+# Solve for weights for data in tableau
+solve <-
+ function(tableau, weights, categorical="FALSE", method="L-BFGS-B", reg="L2", mean=0, var=1, alpha=.01, mf=FALSE, markedness=0, faithfulness=0, mf.weight=1, gradient=TRUE, lower.bound=0, upper.bound=Inf, ...) {
+
+ w <- weights
+
+ if(categorical) {
+ obj <- categorical.objective(tableau)
+ }
+ else {
+ obj <- variable.objective(tableau)
+ }
+
+ if(reg == "L1") {
+ obj <- l1.prior(obj, mean=mean, alpha=alpha)
+ }
+ else if(reg == "L2") {
+ obj <- l2.prior(obj, mean=mean, var=var)
+ }
+
+ if(mf) {
+ obj <- difference.bias(obj, markedness=markedness, faithfulness=faithfulness, mf.weight=mf.weight)
+ }
+
+
+ if(gradient) {
+ grad.optim(weights, obj, lower=lower.bound, upper=upper.bound, method=method, ...)
+ }
+ else {
+ optim(weights, obj, lower=lower.bound, upper=upper.bound, method=method, ...)
+ }
+
+ }
+
+# Solve function with global normalization
+solve_global_norm <- function(tableau, weights, categorical="FALSE", method="L-BFGS-B", reg="L2", mean=0, var=1, alpha=.01, mf=FALSE, markedness=0, faithfulness=0, mf.weight=1, gradient=TRUE, lower.bound=0, upper.bound=Inf, observed_probs, ...) {
+
+ w <- weights
+
+ # Define the objective function with global normalization
+ if(categorical) {
+ obj <- global.normalization("categorical.objective", tableau, observed_probs)
+ }
+ else {
+ obj <- global.normalization("variable.objective", tableau, observed_probs)
+ }
+
+ if(reg == "L1") {
+ obj <- l1.prior(obj, mean=mean, alpha=alpha)
+ }
+
+ if(mf) {
+ obj <- difference.bias(obj, markedness=markedness, faithfulness=faithfulness, mf.weight=mf.weight)
+ }
+
+ if(gradient) {
+ result <- grad.optim(w, obj, lower=lower.bound, upper=upper.bound, method=method, ...)
+ } else {
+ result <- optim(w, obj, lower=lower.bound, upper=upper.bound, method=method, ...)
+ }
+
+ # Calculate the final probabilities with the optimized weights
+ if(tableau@hidden) {
+ violations <- as.matrix(tableau@data[, -(1:4)])
+ }
+ else {
+ violations <- as.matrix(tableau@data[, -(1:3)])
+ }
+ scores <- exp(c(violations %*% result$par))
+ total_score <- sum(scores)
+ final_probabilities <- scores / total_score
+
+ result$scores <- scores
+ result$final_probabilities <- final_probabilities
+
+ return(result)
+}
+
+global.normalization <- function(fun, tableau, observed_probabilities) {
+ fun <- match.fun(fun)
+
+ function(w) {
+ # constraint violations
+ if(tableau@hidden) {
+ violations <- as.matrix(tableau@data[, -(1:4)])
+ }
+ else {
+ violations <- as.matrix(tableau@data[, -(1:3)])
+ }
+
+ scores <- exp(c(violations %*% w))
+ total_score <- sum(scores)
+ probabilities <- scores / total_score
+
+ res <- -sum(observed_probabilities * log(probabilities))
+
+ # Calculate gradient
+ gradient <- colSums((probabilities - observed_probabilities) * violations)
+ attr(res, "gradient") <- gradient
+
+ res
+ }
+}
+
+# OBJECTIVE FUNCTIONS
+
+# Objective function for categorical problems
+categorical.objective <- function(tableau) {
+
+ if(tableau@theory == "MaxEnt") {
+ maxent.categorical.objective(tableau)
+ }
+ else {
+ warning("Objective functions not defined for theories other than Maximum Entropy.")
+ }
+}
+
+# Objective function for non-categorical problems
+variable.objective <- function(tableau) {
+
+ if(tableau@theory == "MaxEnt") {
+ maxent.variable.objective(tableau)
+ }
+ else {
+ warning("Objective functions not defined for theories other than Maximum Entropy.")
+ }
+}
+
+# Objective function for Maximum Entropy Grammar
+# (Categorical using maximum likelihood)
+maxent.categorical.objective <- function(tableau) {
+
+ # constraint violations
+ if(tableau@hidden) {
+ violations <- as.matrix(tableau@data[, -(1:4)])
+ }
+ else {
+ violations <- as.matrix(tableau@data[, -(1:3)])
+ }
+
+ obj <- function(w) {
+
+ scores <- exp(c(violations %*% w))
+
+ # maximum entropy calculation
+ res <-
+ -sum(by(data.frame(p=tableau@data$probability, scores=scores), tableau@data$input,
+ function(form) {
+ Z <- sum(form$scores)
+ #num <- form$scores[form$p == 1]
+ num <- sum(form$scores[form$p == 1])
+
+ #sum(log(num) - log(Z))
+ log(num) - log(Z)
+ }
+ ))
+ grad(res) <- -colSums(do.call(rbind,
+ by(tableau@data, tableau@data$input, function(form) {
+
+ if(tableau@hidden) {
+ features <- as.matrix(form[, -(1:4)])
+ }
+ else {
+ features <- as.matrix(form[, -(1:3)])
+ }
+
+ scores <- exp(c(features %*% w))
+ den.p <- scores / sum(scores)
+ num.p <- scores[form$probability == 1] / sum(scores[form$probability == 1])
+ den <- colSums(den.p * features)
+ num <- colSums(num.p * features[form$probability == 1, , drop=FALSE])
+ num - den
+ })))
+
+ res
+
+ }
+
+ obj
+}
+
+# Objective function for Maximum Entropy Grammar
+# (Variable using minimization of Kullback-Leibler divergence)
+maxent.variable.objective <- function(tableau) {
+
+ # constraint violations
+ if(tableau@hidden) {
+ violations <- as.matrix(tableau@data[, -(1:4)])
+ }
+ else {
+ violations <- as.matrix(tableau@data[, -(1:3)])
+ }
+
+ obj <- function(w) {
+
+ scores <- exp(c(violations %*% w))
+
+ # K-L divergence calculation
+ res <-
+ # sum over inputs
+ -sum(by(data.frame(y=tableau@data$output, q=tableau@data$probability, scores=scores), tableau@data$input,
+ function(form) {
+ Z <- sum(form$scores)
+
+ # sum over outputs
+ sum(by(form, form$y[,drop=TRUE],
+ function(output) {
+ q.target <- output$q[1]
+ num <- sum(output$scores)
+
+ p.weights <- sum(log(num) - log(Z))
+
+ if(q.target == 0.0) {
+ 0.0
+ }
+ else {
+ q.target * (-log(q.target) + p.weights)
+ }
+
+ }
+ ))
+
+ }
+ ))
+
+ grad(res) <- -colSums(do.call(rbind,
+ by(tableau@data, tableau@data$input,
+ # sum over inputs
+ function(form) {
+
+ if(tableau@hidden) {
+ features <- as.matrix(form[, -(1:4)])
+ }
+ else {
+ features <- as.matrix(form[, -(1:3)])
+ }
+
+ scores <- exp(c(features %*% w))
+ den.p <- scores / sum(scores)
+ den <- colSums(den.p * features)
+
+ colSums(do.call(rbind,
+ by(data.frame(p=form$probability, scores=scores, features=features), form$output[,drop=TRUE],
+ #sum over all output forms
+ function(output) {
+ p <- output$p[1]
+ feat <- as.matrix(output[,-(1:2)])
+
+ num.p <- output$scores / sum(output$scores)
+
+ num <- colSums(num.p * feat)
+
+ p * (num - den)
+ })
+ ))
+ }
+ )))
+
+ res
+ }
+
+ obj
+}
+
+
+# REGULARIZATION
+
+# Puts an L1 prior on an optimization
+l1.prior <- function(fun, alpha=1, mean=0, ...) {
+ fun <- match.fun(fun)
+ function(pv, ...) {
+ res <- fun(pv, ...) + sum(alpha * abs(pv - mean))
+
+ if ( !is.null(grad(res)) ) grad(res) <- grad(res) + alpha
+
+ res
+ }
+}
+
+# Puts an L2 (log-space Gaussian) prior on an optimization
+l2.prior <- function(fun, var=1, mean=0, ...) {
+ fun <- match.fun(fun)
+ function(pv) {
+
+ res <- fun(pv) + sum( (pv - mean) / sqrt(2 * var))^2
+
+
+ if ( !is.null(grad(res)) ) grad(res) <- grad(res) + ((pv - mean) / var)
+
+ res
+ }
+}
+
+# Adds a term for Markedness over Faithfulness
+# Sum(M) - Sum(F) is maximized
+# Needs a vector of booleans stating whether something
+# is Markedness or Faithfulness
+difference.bias <- function(fun, markedness, faithfulness, mf.weight=1) {
+ fun <- match.fun(fun)
+
+ function(pv, ...) {
+
+ faith <- faithfulness * pv
+ mark <- markedness * pv
+
+ diff <- faith - mark
+
+ res <- fun(pv, ...) + mf.weight * sum(diff)
+
+ if ( !is.null(grad(res)) ) grad(res) <- grad(res) + mf.weight * (faithfulness - markedness)
+
+ res
+ }
+}
+
+# HELPER FUNCTIONS
+
+# Computes an "information theoretic" logarithm
+# That is, logarithm with 0 log 0 = 0
+# Helper function for K-L divergence minimization
+safe.log <- function(x) {
+ sapply(
+ x,
+ function(y) {
+ minval <- 1E-320
+ if(is.nan(y)) {
+ 0
+ }
+ else if(y < minval) {
+ log(minval)
+ }
+ else {
+ log(y)
+ }
+ }
+ )
+}
+
+grad <- function(x) attr(x, "gradient")
+
+"grad<-" <- function(x,value) {
+ attr(x, "gradient") <- value
+ x
+}
+
+cache.grad <- function(fun) {
+ fun <- match.fun(fun)
+ saved.gradient <- numeric(0)
+ obj <- function(...) {
+ res <- fun(...)
+ saved.gradient <<- grad(res)
+ as.numeric(res)
+ }
+ grad(obj) <- function(pv) saved.gradient
+ obj
+}
+
+grad.optim <- function(par, fun, method="L-BFGS-B", ...) {
+ fun <- match.fun(fun)
+ gfun <- cache.grad(fun)
+
+ optim(par, gfun, grad(gfun), method=method, ...)
+
+}
diff --git a/shiny_app/tableau.R b/shiny_app/tableau.R
new file mode 100644
index 0000000..49fec83
--- /dev/null
+++ b/shiny_app/tableau.R
@@ -0,0 +1,639 @@
+# Provides objective functions and solving for MaxEnt
+# Requires: tableau
+
+# Authors: Robert Staubs
+# Contact: rstaubs@linguist.umass.edu)
+# Last update: 2/6/2012
+
+
+# Class for storing tableau data and theory type
+setClass("Tableau",
+ representation(
+ data="data.frame",
+ hidden="logical",
+ theory="character"
+ )
+)
+
+setGeneric("setTheory<-", function(.Object, ...) standardGeneric("setTheory<-"))
+setReplaceMethod("setTheory",
+ "Tableau",
+ function(.Object, value) {
+ th <- toupper(value)
+ if(th == "HG" || th == "HARMONIC GRAMMAR") {
+ slot(.Object, "theory") <- "HG"
+ }
+ else if(th == "ME" || th == "MAXENT" || th == "MAXIMUM ENTROPY") {
+ slot(.Object, "theory") <- "MaxEnt"
+ }
+ else if(th == "NHG" || th == "NOISYHG" || th == "NOISY HG" || th == "NOISY HARMONIC GRAMMAR") {
+ slot(.Object, "theory") <- "NoisyHG"
+ }
+ else {
+ warning(paste("Cannot find theory", value, ", using HG."))
+ slot(.Object, "theory") <- "HG"
+ }
+
+ return(.Object)
+ }
+)
+
+# Construct a tableau from data frame and theory
+setMethod("initialize",
+ "Tableau",
+ function(.Object, data, hidden=logical(0), theory=character(0)) {
+ .Object@data <- data
+
+ if(!hidden) {
+ .Object@data$probability <- normalize(data,hidden)$probability
+ }
+
+ setTheory(.Object) <- theory
+
+ .Object@hidden <- hidden
+
+ .Object
+ }
+)
+
+# Produce data frame with output distribution given weights
+setGeneric("distribution", function(.Object, ...) standardGeneric("distribution"))
+setMethod("distribution",
+ "Tableau",
+ function(.Object, w, ...) {
+
+ if(.Object@theory == "HG") {
+ overt.view(hg.distribution(.Object, w, ...), theory=.Object@theory, sum=TRUE)
+ }
+ else if (.Object@theory == "MaxEnt") {
+ maxent.distribution(.Object, w, sum=TRUE, ...)
+ }
+ else if (.Object@theory == "NoisyHG") {
+ overt.view(noisyhg.distribution(.Object, w, ...), theory=.Object@theory, sum=TRUE)
+ }
+ }
+)
+
+# Produce data frame with output and hidden structure distribution given weights
+setGeneric("hidden.distribution", function(.Object, ...) standardGeneric("hidden.distribution"))
+setMethod("hidden.distribution",
+ "Tableau",
+ function(.Object, w, ...) {
+
+ if(.Object@theory == "HG") {
+ hg.distribution(.Object, w, ...)
+ }
+ else if (.Object@theory == "MaxEnt") {
+
+ maxent.distribution(.Object, w, sum=FALSE, ...)
+ }
+ else if (.Object@theory == "NoisyHG") {
+ noisyhg.distribution(.Object, w, ...)
+ }
+ }
+)
+
+# Takes data frame of input, output, hidden, and probability
+# Yields data frame of input, output, and probability
+overt.view <- function(data, theory, sum=FALSE) {
+
+ res <- list()
+ inputs <- unique(data$input)
+ for(i in 1:length(inputs)) {
+ if(sum == FALSE) {
+ dat <- data[data$input == inputs[i],]
+ dat <- dat[!duplicated(dat$output),]
+ res <- rbind(res, dat[,c(1,2,4)])
+ }
+ else {
+ outputs <- unique(data$output[data$input == inputs[i]])
+ for(j in 1:length(outputs)) {
+ dd <- data[data$input == inputs[i] & data$output == outputs[j],]
+ if(theory == "HG") {
+ p.sum <- max(dd$probability)
+ }
+ else {
+ p.sum <- sum(dd$probability)
+ }
+ res <- rbind(res, list(input=as.character(inputs[i]), output=as.character(outputs[j]), probability=p.sum))
+ }
+ }
+ }
+ data.frame(res)
+
+}
+
+# Produce data frame with observed distribution over outputs in tableau
+setGeneric("observed", function(.Object) standardGeneric("observed"))
+setMethod("observed",
+ "Tableau",
+ function(.Object) {
+ if(.Object@hidden) {
+ overt.view(.Object@data, theory=.Object@theory, sum=FALSE)
+ }
+ else {
+ .Object@data[1:3]
+ }
+ }
+)
+
+# Produce data frame with error over outputs in tableau
+setGeneric("error", function(.Object, ...) standardGeneric("error"))
+setMethod("error",
+ "Tableau",
+ function(.Object, w, ...) {
+ expected <- distribution(.Object, w, ...)
+ observed <- observed(.Object)
+
+ data.frame(observed[1:2], probability=observed$probability - expected$probability)
+ }
+)
+
+# Sum squared error over outputs in tableau
+setGeneric("sse", function(.Object, ...) standardGeneric("sse"))
+setMethod("sse",
+ "Tableau",
+ function(.Object, w, ...) {
+ sum(error(.Object, w, ...)$probability ^ 2)
+ }
+)
+
+# Produce data frame with error over outputs in tableau
+# Takes the output distribution as categorical -- thus error is 0 or total for each input
+setGeneric("categorical.error", function(.Object, ...) standardGeneric("categorical.error"))
+setMethod("categorical.error",
+ "Tableau",
+ function(.Object, w, ...) {
+ expected <- overt.view(hg.distribution(.Object, w, ...), theory=.Object@theory, sum=TRUE)
+ observed <- observed(.Object)
+
+ data.frame(observed[1:2], probability=observed$probability - as.numeric(expected$probability))
+ }
+)
+
+# Sum squared error over outputs in tableau
+# Takes the output distribution as categorical -- thus error is 0 or total for each input
+setGeneric("categorical.sse", function(.Object, ...) standardGeneric("categorical.sse"))
+setMethod("categorical.sse",
+ "Tableau",
+ function(.Object, w, ...) {
+ sum(categorical.error(.Object, w, ...)$probability ^ 2)
+ }
+)
+
+# Produce input-output pair (or input-output-hidden triple) given weights
+setGeneric("sample.pair", function(.Object, ...) standardGeneric("sample.pair"))
+setMethod("sample.pair",
+ "Tableau",
+ function(.Object, w, distribution=numeric(0), hidden=FALSE, ...) {
+ if(length(distribution) == 0) {
+ inputs <- unique(.Object@data$input)
+ dist <- data.frame(input=inputs, probability=rep(1.0/length(inputs), length(inputs)))
+ }
+ else {
+ dist <- distribution
+ }
+
+ # sample an input according to distribution
+ sample.input <- input.sample(.Object, dist)
+
+ sample.list <- conditional.sample(.Object, w, sample.input, hidden=hidden...)
+ sample.output <- sample.list$output
+
+ if(hidden == FALSE || !.Object@hidden) {
+ list(input=as.character(sample.input), output=as.character(sample.output))
+ }
+ else {
+ list(input=as.character(sample.input), output=as.character(sample.output), hidden=as.character(sample.list$hidden))
+ }
+ }
+)
+
+# Produce a sample output given the input and weights
+setGeneric("conditional.sample", function(.Object, ...) standardGeneric("conditional.sample"))
+setMethod("conditional.sample",
+ "Tableau",
+ function(.Object, w, input, specified.output=character(0), hidden=FALSE, ...) {
+
+ if(.Object@theory == "HG") {
+ hg.conditional.sample(.Object, w, input, specified.output, hidden, ...)
+ }
+ else if (.Object@theory == "MaxEnt") {
+ maxent.conditional.sample(.Object, w, input, specified.output, hidden, ...)
+ }
+ else if (.Object@theory == "NoisyHG") {
+ noisyhg.conditional.sample(.Object, w, input, specified.output, hidden, ...)
+ }
+ }
+)
+
+
+# Produce data frame with HG output distribution given weights
+setGeneric("hg.distribution", function(.Object, ...) standardGeneric("hg.distribution"))
+setMethod("hg.distribution",
+ "Tableau",
+ function(.Object, w) {
+ data <- .Object@data
+
+ dist <- data.frame()
+
+ inputs <- unique(.Object@data$input)
+ for(i in 1:length(inputs)) {
+
+ this.input <- inputs[i]
+
+ # section of tableaux for inputs and outputs
+ if(!.Object@hidden) {
+ this.data <- data[data$input == this.input,1:2]
+ }
+ else {
+ # use one extra column for hidden structure
+ this.data <- data[data$input == this.input,1:3]
+ }
+ this.data <- cbind(this.data, probability=rep(0, length(this.data$output)))
+
+ this.sample <- hg.conditional.sample(.Object, w, this.input, ties=TRUE)
+ if(!.Object@hidden) {
+ this.data$probability[this.data$output == as.character(this.sample$output)] <- 1
+ }
+ else {
+ this.data$probability[this.data$output == as.character(this.sample$output) & this.data$hidden == as.character(this.sample$hidden)] <- 1
+ }
+
+ dist <- rbind(dist, this.data)
+ }
+
+ dist
+ }
+)
+
+# Produce data frame with MaxEnt output distribution given weights
+setGeneric("maxent.distribution", function(.Object, ...) standardGeneric("maxent.distribution"))
+setMethod("maxent.distribution",
+ "Tableau",
+ function(.Object, w, sum = TRUE, ...) {
+ data <- .Object@data
+
+ # constraint violations
+ if(.Object@hidden) {
+ violations <- as.matrix(data[, -(1:4)])
+ }
+ else {
+ violations <- as.matrix(data[, -(1:3)])
+ }
+
+ scores <- exp(c(violations %*% w))
+
+ # probabilities over inputs
+ e <- unlist(by(cbind(data, scores=scores), data$input[,drop=TRUE],
+ function(form) {
+ f.s <- form$scores
+
+ input <- form$input
+ output <- form$output
+ if(.Object@hidden) {
+ hidden <- form$hidden
+ }
+
+ uniq <- unique(form$output)
+ # sum probabilities over hidden structures
+ if(sum == TRUE && .Object@hidden) {
+ uniq <- unique(form$output)
+
+ f.s <- c()
+ for(j in 1:length(uniq)) {
+ f.s <- c(f.s, sum(form$scores[form$output == uniq[j]]))
+ }
+
+ f.s <- as.numeric(f.s)
+ }
+
+ # calculate probabilities
+ Z <- sum(f.s)
+ num <- f.s
+
+ p <- num / Z
+
+ inner <- c()
+ if(sum == TRUE || !.Object@hidden) {
+ for(i in 1:length(uniq)) {
+ inner <- c(inner, as.character(unlist(form$input[1])), as.character(unlist(uniq[i])), as.character(p[i]))
+ }
+ }
+ else {
+ for(i in 1:length(form$output)) {
+ inner <- c(inner, as.character(unlist(form$input[1])), as.character(unlist(form$output[i])), as.character(unlist(form$hidden[i])), as.character(p[i]))
+ }
+ }
+
+ unlist(inner)
+ }
+ ))
+
+ if(sum == TRUE) {
+ rez <- data.frame(matrix(e, ncol=3, byrow=TRUE))
+ colnames(rez) <- c("input", "output", "probability")
+
+ rez$probability <- as.numeric(as.character(rez$probability))
+
+ rez[order(match(
+ paste(rez$input,rez$output),
+ paste(data$input,data$output))
+ ),]
+
+ }
+ else {
+ rez <- data.frame(matrix(e, ncol=4, byrow=TRUE))
+ colnames(rez) <- c("input", "output", "hidden", "probability")
+
+ rez$probability <- as.numeric(as.character(rez$probability))
+
+ rez[order(match(
+ paste(rez$input,rez$output,rez$hidden),
+ paste(data$input,data$output,data$hidden))
+ ),]
+ }
+
+ }
+)
+
+# Produce data frame with NoisyHG output distribution given weights
+setGeneric("noisyhg.distribution", function(.Object, ...) standardGeneric("noisyhg.distribution"))
+setMethod("noisyhg.distribution",
+ "Tableau",
+ function(.Object, w, samples=1000, mean=0, sd=1, ...) {
+
+ data <- .Object@data
+
+ dist <- data.frame()
+
+ inputs <- unique(.Object@data$input)
+ iterat <- samples / length(inputs)
+
+ for(i in 1:length(inputs)) {
+
+ this.input <- inputs[i]
+
+ # section of tableaux for inputs and outputs
+ if(!.Object@hidden) {
+ this.data <- data[data$input == this.input,1:2]
+ }
+ else {
+ # use one extra column for hidden structure
+ this.data <- data[data$input == this.input,1:3]
+ }
+ this.data <- cbind(this.data, probability=rep(0, length(this.data$output)))
+
+ for(iter in 1:iterat) {
+ this.sample <- noisyhg.conditional.sample(.Object, w, this.input, mean=mean, sd=sd, ...)
+ this.output <- this.sample$output
+ if(!.Object@hidden) {
+ this.data$probability[which(this.data$output == this.output)] <- 1 + this.data$probability[which(this.data$output == this.output)]
+ }
+ else {
+ this.data$probability[this.data$output == this.output & this.data$hidden == this.sample$hidden] <- 1 + this.data$probability[which(this.data$output == this.output & this.data$hidden == this.sample$hidden)]
+ }
+ }
+
+ this.data$probability <- this.data$probability / (1.0 * iterat)
+ dist <- rbind(dist, this.data)
+ }
+
+ dist
+
+
+ }
+)
+
+# Produce an input according to sampling distribution
+setGeneric("input.sample", function(.Object, ...) standardGeneric("input.sample"))
+setMethod("input.sample",
+ "Tableau",
+ function(.Object, distribution) {
+
+ # sample an input according to distribution
+ inputs <- unique(.Object@data$input)
+
+ p <- runif(1)
+ p.sum <- 0
+ for(i in 1:length(inputs)) {
+ p.sum <- p.sum + distribution$probability[distribution$input == inputs[i]]
+ if(p.sum >= p) {
+ sample.input <- inputs[i]
+ break
+ }
+ }
+
+ sample.input
+ }
+)
+
+# Produce an output according to data distribution
+setGeneric("output.sample", function(.Object, ...) standardGeneric("output.sample"))
+setMethod("output.sample",
+ "Tableau",
+ function(.Object, input) {
+
+ obs <- observed(.Object)
+ segment <- obs[obs$input == input,]
+
+ # sample an output according to distribution
+ outputs <- unique(segment$output)
+
+ p <- runif(1)
+ p.sum <- 0
+ for(i in 1:length(outputs)) {
+ p.sum <- p.sum + segment$probability[segment$output == outputs[i]]
+ if(p.sum >= p) {
+ sample.output <- outputs[i]
+ break
+ }
+ }
+
+ sample.output
+ }
+)
+
+# Produce input-output pair from data
+setGeneric("data.sample", function(.Object, ...) standardGeneric("data.sample"))
+setMethod("data.sample",
+ "Tableau",
+ function(.Object, distribution) {
+ if(length(distribution) == 0) {
+ inputs <- unique(.Object@data$input)
+ dist <- data.frame(input=inputs, probability=rep(1.0/length(inputs), length(inputs)))
+ }
+ else {
+ dist <- distribution
+ }
+
+ # sample an input according to distribution
+ sample.input <- input.sample(.Object, dist)
+
+ # find winner for that input
+ sample.output <- output.sample(.Object, sample.input)
+
+ list(input=as.character(sample.input), output=as.character(sample.output))
+ }
+)
+
+# Produce the MaxEnt winning output given the input
+setGeneric("maxent.conditional.sample", function(.Object, ...) standardGeneric("maxent.conditional.sample"))
+setMethod("maxent.conditional.sample",
+ "Tableau",
+ function(.Object, w, sample.input, specified.output=character(0)) {
+
+ # get probability distribution over output forms
+ output.distribution <- maxent.distribution(.Object, w)
+
+ output.distribution <- output.distribution[output.distribution$input == sample.input,]
+
+ # if output specified, remove conflicting forms
+ # (used for hidden structure)
+ if(length(specified.output) != 0) {
+ output.distribution <- output.distribution[output.distribution$output == specified.output, ]
+ }
+
+ #TODO get hidden structures selected given specified.output
+ # currently don't get any hidden structures in for maxent
+
+ outputs <- unique(.Object@data$output[.Object@data$input == sample.input])
+ q <- runif(1)
+ q.sum <- 0
+ for(i in 1:length(outputs)) {
+ q.sum <- q.sum + output.distribution$probability[output.distribution$output == outputs[i]]
+ if(q.sum >= q) {
+ sample.output <- outputs[i]
+ break
+ }
+ }
+
+ as.character(sample.output)
+ }
+)
+
+# Produce the HG winning output given the input
+setGeneric("hg.conditional.sample", function(.Object, ...) standardGeneric("hg.conditional.sample"))
+setMethod("hg.conditional.sample",
+ "Tableau",
+ function(.Object, w, sample.input, specified.output=character(0), ties=FALSE) {
+
+ data <- .Object@data[.Object@data$input == sample.input,]
+
+ # if output specified, remove conflicting forms
+ # (used for hidden structure)
+ if(length(specified.output) != 0) {
+ data <- data[data$output == specified.output, ]
+ }
+
+ # constraint violations and scores
+ if(.Object@hidden) {
+ violations <- as.matrix(data[, -(1:4)])
+ }
+ else {
+ violations <- as.matrix(data[, -(1:3)])
+ }
+ scores <- c(violations %*% w)
+
+ winners <- which(max(scores) == scores)
+
+
+ if(ties == FALSE) {
+ winners <- sample(winners, 1)
+ }
+
+ if(.Object@hidden) {
+ list(output=data$output[winners], hidden=data$hidden[winners])
+ }
+ else {
+ list(output=data$output[winners])
+ }
+ }
+)
+
+# Produce the interpretive parse for the given input-output pair
+setGeneric("interpretive.parse", function(.Object, ...) standardGeneric("interpretive.parse"))
+setMethod("interpretive.parse",
+ "Tableau",
+ function(.Object, w, input, output) {
+ hg.conditional.sample(.Object, w, input, output, FALSE)
+ }
+)
+
+# Produce a sample NoisyHG output given the input
+setGeneric("noisyhg.conditional.sample", function(.Object, ...) standardGeneric("noisyhg.conditional.sample"))
+setMethod("noisyhg.conditional.sample",
+ "Tableau",
+ function(.Object, w, sample.input, specified.output=character(0), mean=0, sd=1) {
+ data <- .Object@data[.Object@data$input == sample.input,]
+
+ # if output specified, remove conflicting forms
+ # (used for hidden structure)
+ if(length(specified.output) != 0) {
+ data <- data[data$output == specified.output, ]
+ }
+
+ # constraint violations and scores
+ if(.Object@hidden) {
+ violations <- as.matrix(data[, -(1:4)])
+ }
+ else {
+ violations <- as.matrix(data[, -(1:3)])
+ }
+ scores <- c(violations %*% (w + rnorm(length(w), mean = mean, sd = sd)))
+
+ winners <- which(max(scores) == scores)
+ winners <- sample(winners, 1)
+
+ if(.Object@hidden) {
+ list(output=data$output[winners], hidden=data$hidden[winners])
+ }
+ else {
+ list(output=data$output[winners])
+ }
+
+ }
+
+)
+
+# Normalizes a dataframe with tableau data by input (to make probability distributions)
+normalize <- function(data, hidden) {
+ e <- unlist(by(data, data$input,
+ function(t) {
+ pr <- t$probability / sum(t$probability)
+ inner <- c()
+ for(i in 1:length(t$output)) {
+ if(!hidden) {
+ inner <- c(inner, as.character(unlist(t$input[1])), as.character(unlist(t$output[i])), as.character(pr[i]))
+ }
+ else {
+ inner <- c(inner, as.character(unlist(t$input[1])), as.character(unlist(t$output[i])), as.character(unlist(t$hidden[i])), as.character(pr[i]))
+ }
+ }
+ unlist(inner)
+ }
+ ))
+ if(hidden) {
+ rez <- data.frame(matrix(e, ncol=4, byrow=TRUE))
+ colnames(rez) <- c("input", "output", "hidden", "probability")
+
+ rez$probability <- as.numeric(as.character(rez$probability))
+
+ rez[order(match(
+ paste(rez$input,rez$output,rez$hidden),
+ paste(data$input,data$output,data$hidden))
+ ),]
+ }
+ else {
+ rez <- data.frame(matrix(e, ncol=3, byrow=TRUE))
+ colnames(rez) <- c("input", "output", "probability")
+
+ rez$probability <- as.numeric(as.character(rez$probability))
+
+ rez[order(match(
+ paste(rez$input,rez$output),
+ paste(data$input,data$output))
+ ),]
+ }
+}
+