diff --git a/Parses.R b/Parses.R index 71fefa2..1c1dc42 100644 --- a/Parses.R +++ b/Parses.R @@ -40,7 +40,7 @@ buildExpr <- function(pressStatus) { } expr<-substr(expr, 1, nchar(expr)-2) expr<-paste0(expr, ')') - + return(expr) } @@ -55,14 +55,14 @@ parseScenario <- function(press, prefix = 'p') { cat('Duplicated pressure node names found') print(pressNodes[duplicated(pressNames)]) } - + return(list( - timeSeq=press, - nodes=data.frame(name = pressNames, - code=paste0(prefix, seq(1:length(pressNames))), - growth = coefs[,'growth'], - confidence=coefs[,'confidence'], - layer=coefs[,'layer'], + timeSeq=press, + nodes=data.frame(name = pressNames, + code=paste0(prefix, seq(1:length(pressNames))), + growth = coefs[,'growth'], + confidence=coefs[,'confidence'], + layer=coefs[,'layer'], stringsAsFactors = FALSE), edges=data.frame(input=NULL, output=NULL, impact=NULL) )) @@ -75,7 +75,7 @@ getInitial <- function(string, letter) { split <- function(cell) { params <- unlist(strsplit(cell, ',')) values <- rep(0, length(states)) - + for (n in 1:length(params)) { kvp <- unlist(strsplit(params[n], '=')) ref <- match(getInitial(trimws(kvp[1])), getInitial(states)) @@ -106,11 +106,11 @@ buildGraph <- function(model, desc) { #model contains the following # node table, edge table - + #descriptor (desc) contains: #inputCode - the top layer of the model #outputCodes - all subsequent layers to be included in the model - + inputNodes <- model$nodes$code[which(startsWith(model$nodes$code, desc$inputCode))] inputText <- paste0("[", inputNodes, "]", collapse ="") @@ -128,19 +128,19 @@ buildGraph <- function(model, desc) { edges <- paste0(edges, paste0("[", outNodes[idx], "|", substr(inputsStr, start=1, stop=(nchar(inputsStr)-1)), "]")) #Make the coefficient of the distribution - coefVal <- setNames(c(model$nodes$growth[nodeRef], model$edges$values[rows]), + coefVal <- setNames(c(model$nodes$growth[nodeRef], model$edges$values[rows]), c("(Intercept)", model$edges$input[rows]) ) #str(coefVal) - outDist[[idx]] <- list(coef = coefVal, + outDist[[idx]] <- list(coef = coefVal, sd = model$nodes$confidence[nodeRef]) } - + print('about to build network') print(paste0(inputText, edges)) net <- model2network(paste0(inputText, edges), debug=TRUE) - + print('network build successful') inDist <- vector(mode="list", length=length(inputNodes)) @@ -153,19 +153,19 @@ buildGraph <- function(model, desc) { allDists = as.list(setNames(c(inDist, outDist), c(inputNodes, outNodes))) cfit = custom.fit(net, allDists) - + cat('about to calculate sample distributions') print(outNodes) - + sampleDists <- cpdist(cfit, nodes = outNodes, evidence = TRUE, n = 10000, method = "lw") summDists <- summary(sampleDists) #stdDev <- sd(sampleDists) - + print('sample distribution build successful') - + model$edges$input <- model$nodes$name[match(model$edges$input, model$nodes$code)] model$edges$output <- model$nodes$name[match(model$edges$output, model$nodes$code)] - + return( list( nodes = model$nodes, @@ -180,10 +180,10 @@ buildGraph <- function(model, desc) { getValidNodes <- function(mapping, prevOutputs, prefix) { - + #Find row id for input nodes, internal and published inputNodes <- mapping[2:nrow(mapping),1] - + #check that all input nodes are in the previous table inputNodes <- delNA(mapping[mapping[,"Node.Type"] == 'input', "Nodes"]) if (length(inputNodes)>0) { @@ -192,27 +192,27 @@ getValidNodes <- function(mapping, prevOutputs, prefix) { print(inputNodes[!inputNodes %in% prevOutputs$name]) } } else print('Invalid sheet - table must have at least one input row containing names from previous table') - + #Check the row headings concur with previous names validInputs <- delNA(inputNodes[which(unique(inputNodes) %in% prevOutputs$name)]) if (length(validInputs)==0) print('Invalid sheet - table must have at least one input row containing names from previous table') - - + + inputInts <- delNA(inputNodes[mapping$Node.Type!='link']) if (sum(duplicated(inputInts))>0) { cat('Duplicated input node names found') print(inputNodes[duplicated(inputNodes)]) } - + outNodes <- delNA(colnames(mapping)[FIRST_NODE_COL:ncol(mapping)]) if (sum(duplicated(outNodes))>0) { cat('Duplicated output node names found') print(outNodes[duplicated(outNodes)]) } - - + + #check that all internal nodes are in the columns intNodes <- delNA(mapping[mapping[,"Node.Type"] == 'internal', "Nodes"]) if (length(intNodes)>0) { @@ -221,17 +221,17 @@ getValidNodes <- function(mapping, prevOutputs, prefix) { print(intNodes[!intNodes %in% outNodes]) } } - + coefs <- matrix(data=NA, nrow=length(outNodes), ncol=3, dimnames=list(NULL, c('growth', 'confidence', 'layer'))) for (idx in 1:length(outNodes)) { col <- match(outNodes[idx], colnames(mapping)) coefs[idx,] <- as.numeric(split(mapping[1, col]))[match(c('growth', 'confidence', 'layer'), states)] } - + print(coefs) return(data.frame( - code=c(prevOutputs$code, paste0(prefix, seq(1:length(outNodes)))), + code=c(prevOutputs$code, paste0(prefix, seq(1:length(outNodes)))), name=c(prevOutputs$name, outNodes), growth=c(prevOutputs$growth, coefs[,"growth"]), confidence=c(prevOutputs$confidence, coefs[,"confidence"]), @@ -248,17 +248,17 @@ getValidEdges <- function(mapping, nodeDF, prevEdge=NULL, prefix) { str(nodeDF) edgeCols <- c('inputNode', 'outputNode', 'impact') edgeM <- matrix(data=NA, nrow=0, ncol=length(edgeCols), dimnames=list(NULL, edgeCols)) - + #to start let just get the statements and print them out.... for (col in FIRST_NODE_COL:ncol(mapping)) { count=0 - + for (row in 2:nrow(mapping)) { if (!is.na(mapping[row, col])) { - edgeM <- rbind(edgeM, - c(getCode(mapping[row, 1], nodeDF), - getCode(colnames(mapping)[col], nodeDF), + edgeM <- rbind(edgeM, + c(getCode(mapping[row, 1], nodeDF), + getCode(colnames(mapping)[col], nodeDF), split(mapping[row,col])[match('impact', states)] ) ) @@ -285,10 +285,10 @@ getValidEdges <- function(mapping, nodeDF, prevEdge=NULL, prefix) { } parseMapping <- function(mapping, prevOutputs, prefix) { - + mapping <- mapping[,-1] mapping[,1] <- cleanTitles(mapping[,1]) - + nodeDF <- getValidNodes(mapping, prevOutputs$nodes, prefix) edgeDF <- getValidEdges(mapping, nodeDF, prevEdge=prevOutputs$edges, prefix) @@ -302,32 +302,32 @@ parseMapping <- function(mapping, prevOutputs, prefix) { parseSheet <- function(fName) { #get sheet names - + print(paste('starting sheet load', fName)) - + if (file.exists(fName)) { names <- openxlsx::getSheetNames(fName) - + if (length(names)>0) { - + sheets <- sort(delNA(match(names, mappings))) - + cat('starting sheet parse') print(sheets) - + if (sum(sheets==refs)==length(refs)) { #read all mapping tables scenario <- parseScenario(readXL(fName,mappings[1], startRow=1), prefix='p') p_ba <- parseMapping(readXL(fName,mappings[2], startRow=1), scenario, prefix='ba') p_op <- parseMapping(readXL(fName,mappings[3], startRow=1), p_ba, prefix='op') p_es <- parseMapping(readXL(fName,mappings[4], startRow=1), p_op, prefix='es') - + #print('building graphs') #p_baNet <- buildGraph(p_ba, desc=list(inputCode='p', outputCodes='ba')) #p_opNet <- buildGraph(p_op, desc=list(inputCode='p', outputCodes=c('ba', 'op'))) #p_esNet <- buildGraph(p_es, desc=list(inputCode='p', outputCodes=c('ba', 'op', 'es'))) - + print('sheet load completed') return( #list( diff --git a/app.R b/app.R index 1698101..a737005 100644 --- a/app.R +++ b/app.R @@ -42,7 +42,7 @@ legends <- c('Pressures', addResourcePath("js", "./www/js") ui<-dashboardPage( - dashboardHeader(title = "JNCC MESO online", + dashboardHeader(title = "JNCC MESO online", tags$li( id = "dropdownHelp", class = "dropdown", @@ -98,8 +98,8 @@ ui<-dashboardPage( selectInput("modelSelect", "Select MESO model", choices=c(""), selected=NULL, multiple=FALSE), downloadButton("download", "", icon=icon("download")), uiOutput("pressureList") - #selectInput("layerSelect", "Select Transition", - # choices=transitions, + #selectInput("layerSelect", "Select Transition", + # choices=transitions, # selected=NULL, multiple=FALSE) ) ), @@ -108,7 +108,7 @@ ui<-dashboardPage( tabItem(tabName = "1", h2('Impact Distribution'), fluidRow( column( - width=6, + width=6, h4('Effect on bio-assemblage') ), column( @@ -164,7 +164,7 @@ ui<-dashboardPage( # fluidPage( # google_mapOutput(outputId = "map", width = "100%", height = "750px") # ) - #), + #), tabItem(tabName = "3",h4("Ingestion"), fluidPage( p("Select a spreadsheet from your network for input into the JNCC Bayesian Network Analyser:"), @@ -179,32 +179,32 @@ ui<-dashboardPage( server <- function(input, output, session) { #SERVER Constants - + print('Loading data') - + #set_key("AIzaSyAw8_btgGN1drf8qhCxNcotP6r11qEXA_M") dataStorage <- 'data/' - + models<-NULL pressures <- NULL - + .loadStatus <- reactiveValues( valid = c(p=FALSE, ba=FALSE, op=FALSE, es=FALSE), msgs = NULL ) - + .likelihoods <-reactiveValues( p_ba = NULL, ba_os = NULL, os_es = NULL, p_es = NULL ) - + setPressures <- function(newPressures) { pressures <<- newPressures } - + .resistanceScores <- c( ins= -0.01, hr = -0.2, @@ -214,15 +214,15 @@ server <- function(input, output, session) { ssgr = 0, pressSD = 0.5 ) - + .selections <- reactiveValues( - model=1, - bbnImpact=1, + model=1, + bbnImpact=1, bbnNames=FALSE, bbnEdges=FALSE, pressStatus=NULL ) - + getImpact <- function(v) { print(v) if ((v == "INS") || (v == "IV")) return(.resistanceScores[1]) @@ -232,20 +232,20 @@ server <- function(input, output, session) { if (v == "NR") return(.resistanceScores[5]) as.numeric(v) } - + getAvailableModels <- function() { fileList <- list.files(dataStorage, pattern='.xlsx') - + modelList <- list() cnt<-1 - + for (idx in 1:length(fileList)) { print(paste('attempting to load', paste0(dataStorage, fileList[idx]))) - + tmp <- parser$parseSheet(paste0(dataStorage, fileList[idx])) print(tmp) tmp$edges$values <- sapply(tmp$edges$impact, getImpact) - + if (!is.null(tmp)) { modelList[[cnt]] <- tmp models <<- c(models, substr(fileList[idx], 1, (nchar(fileList[idx])-5))) @@ -257,36 +257,36 @@ server <- function(input, output, session) { updateSelectInput(session, "modelSelect", choices=models) return(modelList) } - + #parse on load sheets in the input sheet folder - replace with R Data modelList <- getAvailableModels() - - + + calcLikelihood <- function(layer, pressStatus) { - + isolate({ #if (layer==1) layerStr='ba' else if (layer==2) layerStr='op' else if (layer==3) layerStr='es' - - + + #layerRange <- which(startsWith(thisModel$nodes$code, layerStr)) - + #nodeCodes <-thisModel$nodes$code[layerRange] #nodeNames <- thisModel$nodes$name[layerRange] - + thisModel <- modelList[[.selections$model]] - + modelList[[.selections$model]]$edges$values <<- sapply(thisModel$edges$impact, getImpact) modelList[[.selections$model]]$nodes$growth <<- .resistanceScores['ssgr'] modelList[[.selections$model]]$nodes$confidence <<- .resistanceScores['pressSD'] - + thisModel <- modelList[[.selections$model]] - + MEANPOS=1 MEANNEG=0 - + expr <- "list(" for (p in 1:nrow(pressStatus)) { if (pressStatus$status[p] == 'On') { @@ -294,12 +294,12 @@ server <- function(input, output, session) { } else { threshold = MEANNEG } - + expr <- paste0(expr, "\"", pressStatus$code[p], "\"=", threshold, ", ") } expr <-substr(expr, 1, nchar(expr)-2) expr<-paste0(expr, ')') - + thisNet <- parser$buildGraph(thisModel, desc=list(inputCode='p', outputCodes=c('ba', 'op', 'es'))) sampleDists <- cpdist( @@ -311,16 +311,16 @@ server <- function(input, output, session) { debug=TRUE ) }) - + print(sampleDists) - + #displayCols <- match(nodeCodes, colnames(sampleDists)) sampleDists <- sampleDists[,match(thisModel$nodes$code, colnames(sampleDists))] means <- apply(sampleDists, 2, mean) stdDev <- apply(sampleDists, 2, sd) - + print(paste('Building likelihoods from model, sample dists', length(thisModel$nodes$name), length(sampleDists))) - + return(data.frame( name = thisModel$nodes$name, code = thisModel$nodes$code, @@ -337,12 +337,12 @@ server <- function(input, output, session) { stringsAsFactors=FALSE )) } - - + + observeEvent(input$modelSelect, { .selections$model <<- match(input$modelSelect, models) }) - + observeEvent(reactiveValuesToList(input), { isolate(myList <- reactiveValuesToList(input)) matches <- match(pressures$code, names(myList)) @@ -350,9 +350,9 @@ server <- function(input, output, session) { if (length(matches)>0) { status <-NULL for (n in 1:length(matches)) status[n] = myList[[matches[n]]] - + newStatus <- data.frame(code=pressures$code, status=status, stringsAsFactors = FALSE) - + if (!identical(newStatus, .selections$pressStatus)) { print('Running calc') #.likelihoods$p_ba <<- calcLikelihood(1, newStatus) @@ -362,42 +362,42 @@ server <- function(input, output, session) { #write.xlsx(.likelihoods$p_es, 'tmp.xlsx') .selections$pressStatus <<- newStatus } - + } }) - + makeRadioButtons <- function(row) { radioButtons(row['code'], row['name'], choices=c('Off', 'On'), selected='Off', inline=TRUE) } - + output$pressureList <- renderUI({ #isolate({ if (!is.null(modelList[[.selections$model]]$nodes)) { pressCodes <- which(startsWith(modelList[[.selections$model]]$nodes$code, 'p')) - pressures <- data.frame(code = modelList[[.selections$model]]$nodes$code[pressCodes], + pressures <- data.frame(code = modelList[[.selections$model]]$nodes$code[pressCodes], name = modelList[[.selections$model]]$nodes$name[pressCodes], stringsAsFactors=FALSE) setPressures(pressures) btnList <- apply(pressures, 1, makeRadioButtons) } }) - + observeEvent(input$bbnImpactSelect, { - #filter nodes and edges to + #filter nodes and edges to .selections$bbnImpact <- thresholds[match(input$bbnImpactSelect, impacts)] print(paste("Setting bbn impact", .selections$bbnImpact)) }) - + observeEvent(input$bbnDisplayNames, { .selections$bbnNames <- input$bbnDisplayNames print(.selections$bbnNames) }) - + observeEvent(input$bbnDisplayEdges, { .selections$bbnEdges <- input$bbnDisplayEdges }) - - + + observeEvent(input$layer1Slider, { showModal( modalDialog({ @@ -410,8 +410,8 @@ server <- function(input, output, session) { sliderInput("ssgr", "Steady state growth rate", -0.1, 0.1,.resistanceScores[6], step=0.01), sliderInput("l1PressSD", "Pressure Std Dev", 0.1, 1.0, .resistanceScores[7], step=0.01) ) - }, - title='Layer 1 controls', + }, + title='Layer 1 controls', footer=tagList( modalButton("Cancel"), actionButton("modalOK", "OK") @@ -419,7 +419,7 @@ server <- function(input, output, session) { size='s') ) }) - + observeEvent(input$modalOK, { print('Modal ok pressed') @@ -430,7 +430,7 @@ server <- function(input, output, session) { .resistanceScores['ins'] <<- -input$l1VL .resistanceScores['ssgr'] <<- input$ssgr .resistanceScores['pressSD'] <<- input$l1PressSD - + print('Running calc') #.likelihoods$p_ba <<- calcLikelihood(1, .selections$pressStatus) #.likelihoods$ba_os <<- calcLikelihood(2, .selections$pressStatus) @@ -439,31 +439,31 @@ server <- function(input, output, session) { removeModal() }) - - + + output$nodeTable <- DT::renderDataTable( - + modelList[[.selections$model]]$nodes, selection = 'single',options = list(searching = TRUE, pageLength = 10, editable=TRUE),server = TRUE, escape = FALSE,rownames= TRUE ) output$edgeTable <- DT::renderDataTable( - + modelList[[.selections$model]]$edges, selection = 'single',options = list(searching = TRUE, pageLength = 10, editable=TRUE),server = TRUE, escape = FALSE,rownames= TRUE ) - + getLabel <- function(value) { sign <- ifelse(value<0, "-", "+") idx <- min(which((abs(value)>=thresholds)==TRUE)) return(paste0(sign, impLabels[idx])) } - + makeBbnGraph <- function(model) { nodes <- model$nodes - + if (.selections$bbnEdges) {labels <- sapply(model$edges$values, getLabel)} else {labels <- rep("", nrow(model$edges))} - + edges <- data.frame( id = rownames(model$edges), from=match(model$edges$input, nodes$code), @@ -472,13 +472,13 @@ server <- function(input, output, session) { label=labels, arrows="to", stringsAsFactors=FALSE - ) + ) if (.selections$bbnNames) {labels <- nodes$name} else {labels <- nodes$code} - + nodeSpacing <- ifelse(.selections$bbnNames, 600, 150) - + palette <- brewer.pal(length(legends), "RdYlGn") - + nodes <- data.frame( id = rownames(nodes), label = labels, @@ -487,81 +487,81 @@ server <- function(input, output, session) { color = palette[as.integer(nodes$layer)], code = nodes$code, stringsAsFactors=FALSE - ) - + ) + edges <- edges[(abs(edges$values)>=.selections$bbnImpact),] - + nodeNet <- nodes[(nodes$code %in% .selections$pressStatus$code[.selections$pressStatus$status %in% c('On')]),] - + #save(nodes, edges, nodeNet, file = 'tmp.RData') - + if (nrow(nodeNet)>0) { - #do pressures + #do pressures edgeNet <- edges[edges$from %in% nodeNet$id, ] idx = 1 repeat { nodesToAdd <- nodes[nodes$id %in% edgeNet$to, ] nodesToAdd <- nodesToAdd[!(nodesToAdd$id %in% nodeNet$id),] - + edgesToAdd <- edges[edges$from %in% nodesToAdd$id, ] edgesToAdd <- edgesToAdd[!(edgesToAdd$id %in% edgeNet$id),] - + idx <- idx + 1 if ((idx>20) || ((nrow(nodesToAdd)==0) && (nrow(edgesToAdd)==0))) break nodeNet <- rbind(nodeNet, nodesToAdd) edgeNet <- rbind(edgeNet, edgesToAdd) - + } #until finished } else edgeNet <- edges - + legendDF <- data.frame( id = 1:length(legends), label = legends, color = palette, stringsAsFactors = FALSE ) - - visNetwork(nodeNet, edgeNet, width = "100%", main='Bayesian Belief Network', submain=input$modelSelect) %>% + + visNetwork(nodeNet, edgeNet, width = "100%", main='Bayesian Belief Network', submain=input$modelSelect) %>% visExport() %>% visLegend(useGroups=FALSE, addNodes=legendDF) %>% visHierarchicalLayout(nodeSpacing=nodeSpacing, direction='LR') %>% visOptions(highlightNearest = TRUE) #%>% #visInteraction(navigationButtons = TRUE, dragNodes = TRUE, dragView = TRUE, zoomView = TRUE) } - + output$bbnGraphPlot <- renderVisNetwork({ makeBbnGraph(modelList[[.selections$model]]) }) - + #observe({ # visNetworkProxy("bbnGraphPlot") %>% # visStabilize(iterations=10) #}) - + getModelName <- function() { paste0('data/', input$modelSelect, '.xlsx') } - + genPlot <- function(boxPlot, title) { if (nrow(boxPlot)>0) { - + palette <- brewer.pal(length(legends), "RdYlGn") #print(palette) - + colours <- palette[as.integer(boxPlot$Group)] - + #print(paste('Box plot, colours', nrow(boxPlot), length(colours))) #cat(colours) xform <- list(categoryorder = "array", categoryarray = boxPlot[,1], zerolinewidth=10) - # + # plot_ly(boxPlot, x = boxPlot[,1], y = ~Range, color = colours, colors = palette, type = "box") %>% layout(xaxis = xform, showlegend=FALSE, title=title) - - } + + } } - + prepPlot <- function(code="ba", name="Bio-Assemblage") { if (!is.null(.likelihoods$p_es)) { inScope <- startsWith(.likelihoods$p_es$code, code) @@ -571,22 +571,22 @@ server <- function(input, output, session) { genPlot(thisPlot, title) } } - + output$layer1 <- renderPlotly({ prepPlot("ba", "Bio-Assemblage") }) - + output$layer2 <- renderPlotly({ prepPlot("op", "Output Processes") }) - + output$layer3 <- renderPlotly({ prepPlot("es", "Ecosystem Services") }) - - + + export <- function(model) { - + #Get the network graph l1 <- orca(prepPlot("ba", "Bio-Assemblage"), 'tmp/layer1.png') l2 <- orca(prepPlot("op", "Output Processes"),'tmp/layer2.png') @@ -594,27 +594,27 @@ server <- function(input, output, session) { #Save pressure list, confidence levels, node and edge tables in xlsx l <- list( - pressures = .selections$pressStatus, + pressures = .selections$pressStatus, nodes = model$nodes, edges = model$edges, settings = as.data.frame(cbind(names(.resistanceScores), .resistanceScores), stringsAsFactors=FALSE) ) - + xl <- write.xlsx(l, 'tmp/dataset.xlsx') - + print('saving xlsx file export tmp/dataset.xlsx') - + zipFile <- zipr(paste0('tmp/MESO-', format(Sys.time(), "%m%d_%H%M"), '.zip'), c('tmp/layer1.png', 'tmp/layer2.png', 'tmp/layer3.png', 'tmp/dataset.xlsx')) - + print(paste('zip file complete', zipFile)) - + return(zipFile) } - - - - - + + + + + output$linkBackgroundData <- downloadHandler( filename = getModelName(), content = function(file) { @@ -622,7 +622,7 @@ server <- function(input, output, session) { }, contentType = "application/xlsx" ) - + output$download <-downloadHandler( filename = paste0('MESO-', format(Sys.time(), "%m%d_%H%M"), '.zip'), content = function(file) { @@ -632,7 +632,7 @@ server <- function(input, output, session) { contentType = "application/zip" ) - + } shinyApp(ui, server)