This commit is contained in:
2022-04-07 10:21:06 +01:00
4 changed files with 489 additions and 554 deletions

View File

@@ -57,10 +57,13 @@ parseScenario <- function(press, prefix = "p") {
ncol = 3,
dimnames = list(NULL, c("growth", "confidence", "layer"))
)
for (col in 2:ncol(press)) {
coefs[col - 1, ] <- as.numeric(split(press[1, col]))[match(c("growth", "confidence", "layer"), states)]
}
press[is.na(press)] <- 0
if (sum(duplicated(pressNames)) > 0) {
cat("Duplicated pressure node names found")
print(pressNodes[duplicated(pressNames)])
@@ -85,7 +88,6 @@ getInitial <- function(string, letter) {
}
split <- function(cell) {
params <- unlist(strsplit(cell, ","))
values <- rep(0, length(states))
@@ -118,7 +120,6 @@ getOutNodes <- function(codes, codeList) {
}
buildGraph <- function(model, desc) {
# model contains the following
# node table, edge table
@@ -126,7 +127,6 @@ buildGraph <- function(model, desc) {
# 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 = "")
@@ -148,19 +148,12 @@ buildGraph <- function(model, desc) {
c(model$nodes$growth[nodeRef], model$edges$values[rows]),
c("(Intercept)", model$edges$input[rows])
)
#str(coefVal)
outDist[[idx]] <- list(coef = coefVal, sd = model$nodes$confidence[nodeRef])
}
print("Saving model prior to network modelling")
modelDefn <- paste0(inputText, edges)
save(modelDefn, file="buildGraph.RData")
#print("about to build network")
#print(paste0(inputText, edges))
net <- model2network(paste0(inputText, edges), debug = FALSE)
@@ -176,15 +169,12 @@ buildGraph <- function(model, desc) {
allDists <- as.list(setNames(c(inDist, outDist), c(inputNodes, outNodes)))
#print(allDists)
cfit <- custom.fit(net, allDists)
cat("about to calculate sample distributions")
#print(outNodes)
print("about to calculate sample distributions")
sampleDists <- cpdist(cfit, nodes = outNodes, evidence = TRUE, n = 10000, method = "lw")
summDists <- summary(sampleDists)
#stdDev <- sd(sampleDists)
print("sample distribution build successful")
@@ -205,7 +195,6 @@ buildGraph <- function(model, desc) {
getValidNodes <- function(mapping, prevOutputs, prefix) {
# Find row id for input nodes, internal and published
inputNodes <- mapping[2:nrow(mapping), 1]
@@ -272,10 +261,6 @@ getCode <- function(name, nodeDF) {
}
getValidEdges <- function(mapping, nodeDF, prevEdge = NULL, prefix) {
#utils::str(nodeDF)
#save(mapping, nodeDF, prevEdge, prefix, file="validEdges.RData")
edgeCols <- c("inputNode", "outputNode", "impact")
edgeM <- matrix(data = NA, nrow = 0, ncol = length(edgeCols), dimnames = list(NULL, edgeCols))
@@ -284,10 +269,11 @@ getValidEdges <- function(mapping, nodeDF, prevEdge = NULL, prefix) {
count <- 0
for (row in 2:nrow(mapping)) {
if (!is.na(mapping[row, col])) {
edgeM <- rbind(edgeM,
c(getCode(mapping[row, 1], nodeDF),
edgeM <- rbind(
edgeM,
c(
getCode(mapping[row, 1], nodeDF),
getCode(colnames(mapping)[col], nodeDF),
split(mapping[row, col])[match("impact", states)]
)
@@ -297,6 +283,7 @@ getValidEdges <- function(mapping, nodeDF, prevEdge = NULL, prefix) {
# if (count == 0) print(paste("No edges found for output", colnames(mapping)[col]))
}
}
if (is.null(prevEdge)) {
return(data.frame(
input = edgeM[, "inputNode"],
@@ -321,8 +308,6 @@ parseMapping <- function(mapping, prevOutputs, prefix) {
nodeDF <- getValidNodes(mapping, prevOutputs$nodes, prefix)
edgeDF <- getValidEdges(mapping, nodeDF, prevEdge = prevOutputs$edges, prefix)
#save(nodeDF, edgeDF, file="mapping.RData")
return(list(
# New structure
nodes = nodeDF,
@@ -331,19 +316,16 @@ 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)
print(sheets)
if (sum(sheets == refs) == length(refs)) {
# read all mapping tables
@@ -360,7 +342,6 @@ parseSheet <- function(fName) {
legend = legend
)
)
} else {
print(paste("Sheets found include", mappings[sheets]))
cat("Missing sheets are:")

187
app.R
View File

@@ -8,7 +8,6 @@ modules::import(shinyBS)
modules::import(bnlearn)
modules::import(visNetwork)
modules::import(RColorBrewer)
modules::import(plotly)
modules::import(openxlsx)
modules::import(zip)
modules::import(DT)
@@ -16,7 +15,6 @@ modules::import(plyr)
modules::import(magrittr)
parser <- modules::use("Parses.R")
rw <- modules::use("reWeight.R")
@@ -31,8 +29,8 @@ impLabels <- c("Very High", "High", "Medium", "Low", "Very Low")
ui <- dashboardPage(
dashboardHeader(title = "JNCC MESO online",
dashboardHeader(
title = "JNCC MESO online",
tags$li(
id = "dropdownHelp",
class = "dropdown",
@@ -80,7 +78,8 @@ ui <- dashboardPage(
)
),
dashboardSidebar(
sidebarMenu(id = "tabs",
sidebarMenu(
id = "tabs",
menuItem("Introduction", tabName = "1", icon = icon("arrow-down")),
menuItem("Pressure Test", tabName = "2", icon = icon("arrow-down")),
menuItem("Bayesian Network", tabName = "3", icon = icon("atom")),
@@ -127,8 +126,10 @@ ui <- dashboardPage(
tags$p(
style = "font-size: 12pt",
"Impact of pressures are as defined in ",
tags$a(href = "https://www.marlin.ac.uk/sensitivity/sensitivity_rationale",
"the Marine Evidence based Sensitivity Assessment (MarESA).", target = "_BLANK")
tags$a(
href = "https://www.marlin.ac.uk/sensitivity/sensitivity_rationale",
"the Marine Evidence based Sensitivity Assessment (MarESA).", target = "_BLANK"
)
),
tags$p(
style = "margin-top: 150px; font-size: 12pt",
@@ -145,7 +146,8 @@ ui <- dashboardPage(
"Copyright Notice: All images, logos and sources are property and copyright of their respected owners"
)
),
tabItem(tabName = "2", h2("Impact Distribution"),
tabItem(
tabName = "2", h2("Impact Distribution"),
fluidRow(
column(
width = 6,
@@ -169,13 +171,14 @@ ui <- dashboardPage(
p("Download results as Excel workbook")
)
),
plotlyOutput("layer1", height = "270px") %>% withSpinner(),
plotly::plotlyOutput("layer1", height = "270px") %>% withSpinner(),
h4("Effect on Ecosystem Processes"),
plotlyOutput("layer2", height = "270px") %>% withSpinner(),
plotly::plotlyOutput("layer2", height = "270px") %>% withSpinner(),
h4("Effect on Ecosystem Services"),
plotlyOutput("layer3", height = "270px") %>% withSpinner()
plotly::plotlyOutput("layer3", height = "270px") %>% withSpinner()
),
tabItem(tabName = "3",h2("Bayesian Network"),
tabItem(
tabName = "3", h2("Bayesian Network"),
fluidPage(
p("Graphical output of the Bayesian Network. Note: The graph will only draw if pressures are applied!"),
fluidRow(
@@ -224,7 +227,6 @@ server <- function(input, output, session) {
palette <- c("firebrick", "coral", "rosybrown", "tan", "salmon", "olivedrab", "seagreen", "aquamarine", "darkcyan", "dodgerblue", "steelblue", "royalblue")
models <- NULL
pressures <- NULL
@@ -262,11 +264,21 @@ server <- function(input, output, session) {
)
getImpact <- function(v) {
if ((v == "INS") || (v == "IV")) return(.resistanceScores[1])
if ((v == "HR") || (v == "III")) return(.resistanceScores[2])
if ((v == "MR") || (v == "II")) return(.resistanceScores[3])
if ((v == "LR") || (v == "I")) return(.resistanceScores[4])
if (v == "NR") return(.resistanceScores[5])
if ((v == "INS") || (v == "IV")) {
return(.resistanceScores[1])
}
if ((v == "HR") || (v == "III")) {
return(.resistanceScores[2])
}
if ((v == "MR") || (v == "II")) {
return(.resistanceScores[3])
}
if ((v == "LR") || (v == "I")) {
return(.resistanceScores[4])
}
if (v == "NR") {
return(.resistanceScores[5])
}
as.numeric(v)
}
@@ -277,7 +289,8 @@ server <- function(input, output, session) {
# save(newNameMap, file="nameMap.RData")
stripStr <- function(nodeStr) {
nodeStr %>% stringr::str_replace_all("\\.", "") %>%
nodeStr %>%
stringr::str_replace_all("\\.", "") %>%
stringr::str_replace_all(" ", "") %>%
stringr::str_replace_all("\\(", "") %>%
stringr::str_replace_all("\\)", "") %>%
@@ -286,24 +299,18 @@ server <- function(input, output, session) {
}
setNewNames <- function(wb, habName) {
#habName <- substr(fileList[idx], 1, (nchar(fileList[idx])-5))
print(habName)
possNames <- newNameMap %>%
dplyr::filter(hab == habName) %>%
dplyr::mutate(node = stripStr(node))
newNodes <- wb$p_es$nodes %>% dplyr::mutate(node = stripStr(name))
print(possNames$node)
print(newNodes$node)
newNames <- apply(newNodes, 1, function(row) {
id <- match(row["node"], possNames$node)
print(paste(id, row["node"]))
possNames$newname[id]
})
print(newNames)
wb$p_es$nodes$name <- newNames
return(wb)
}
@@ -318,15 +325,12 @@ server <- function(input, output, session) {
print(paste("attempting to load", paste0(dataStorage, fileList[idx])))
wb <- parser$parseSheet(paste0(dataStorage, fileList[idx]))
#print(tmp)
wb$p_es$edges$values <- sapply(wb$p_es$edges$impact, getImpact)
if (!is.null(wb)) {
habName <- substr(fileList[idx], 1, (nchar(fileList[idx]) - 5)) %>%
stringr::str_replace_all("_", " ")
print(habName)
wb2 <- setNewNames(wb, habName)
@@ -334,11 +338,10 @@ server <- function(input, output, session) {
models <<- c(models, habName)
print(paste("Model file successfully loaded", fileList[idx]))
#save(tmp, file = "tmp.RData")
cnt <- cnt + 1
}
}
#save(modelList, file="models.RData")
updateSelectInput(session, "modelSelect", choices = models)
return(modelList)
}
@@ -346,15 +349,8 @@ server <- function(input, output, session) {
# parse on load sheets in the input sheet folder - replace with R Data
modelList <- getAvailableModels()
save(modelList, file="model.RData")
#print(load("modelList.RData"))
calcLikelihood <- function(layer, pressStatus, forPlotly) {
isolate({
modelList[[.selections$model]]$p_es$edges$values <<- sapply(modelList[[.selections$model]]$p_es$edges$impact, getImpact)
modelList[[.selections$model]]$p_es$nodes$growth <<- .resistanceScores["ssgr"]
modelList[[.selections$model]]$p_es$nodes$confidence <<- .resistanceScores["pressSD"]
@@ -380,22 +376,12 @@ server <- function(input, output, session) {
print(names(thisModel))
# Now do it in stages with one assessment per stage
thisModel$p_es$nodes$confidence <- 0.1 * thisModel$p_es$nodes$confidence
#save(pressStatus, thisModel, file="beforeWeight.RData")
if (sum(pressStatus$status == "On") > 0) {
thisModel$p_es <- rw$reWeightModel(thisModel$p_es, pressStatus)
} # else nothing to do
#save(pressStatus, thisModel, file="afterWeight.RData")
thisNet <- parser$buildGraph(thisModel$p_es, desc = list(inputCode = "p", outputCodes = c("ba", "op", "es")))
sampleDists <- cpdist(
@@ -438,7 +424,6 @@ server <- function(input, output, session) {
stringsAsFactors = FALSE
))
} else {
return(data.frame(
name = thisModel$p_es$nodes$name,
code = thisModel$p_es$nodes$code,
@@ -449,14 +434,12 @@ server <- function(input, output, session) {
maxes = apply(sampleDists, 2, max),
stringsAsFactors = FALSE
))
}
}
observeEvent(input$modelSelect, {
.selections$model <<- match(input$modelSelect, models)
#.selections$runOnce <<- TRUE
})
observeEvent(reactiveValuesToList(input), {
@@ -478,7 +461,6 @@ server <- function(input, output, session) {
.selections$pressStatus <<- newStatus
}
}
})
@@ -500,7 +482,6 @@ server <- function(input, output, session) {
)
# This assumes all pressures are the same...
setPressures(pressures)
btnList <- apply(pressures, 1, makeRadioButtons)
}
@@ -517,13 +498,13 @@ server <- function(input, output, session) {
observeEvent(input$bbnDisplayEdges, {
.selections$bbnEdges <- input$bbnDisplayEdges
})
observeEvent(input$layer1Slider, {
showModal(
modalDialog({
modalDialog(
{
tagList(
sliderInput("l1VL", "Insensitive", 0.01, 0.2, abs(.resistanceScores[1]), step = 0.01),
sliderInput("l1L", "Low Sensitivity/High resistance", 0.15, 0.5, abs(.resistanceScores[2]), step = 0.01),
@@ -539,13 +520,12 @@ server <- function(input, output, session) {
modalButton("Cancel"),
actionButton("modalOK", "OK")
),
size = "s")
size = "s"
)
)
})
observeEvent(input$modalOK, {
.resistanceScores["nr"] <<- -input$l1VH
.resistanceScores["lr"] <<- -input$l1H
.resistanceScores["mr"] <<- -input$l1M
@@ -557,7 +537,6 @@ server <- function(input, output, session) {
.likelihoods$p_es <<- calcLikelihood(0, .selections$pressStatus, TRUE)
removeModal()
})
@@ -626,8 +605,6 @@ server <- function(input, output, session) {
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
edgeNet <- edges[edges$from %in% nodeNet$id, ]
@@ -643,14 +620,11 @@ server <- function(input, output, session) {
if ((idx > 20) || ((nrow(nodesToAdd) == 0) && (nrow(edgesToAdd) == 0))) break
nodeNet <- rbind(nodeNet, nodesToAdd)
edgeNet <- rbind(edgeNet, edgesToAdd)
} # until finished
} else {
edgeNet <- edges
}
print(paste(nrow(model$legend), length(palette)))
legendDF <- data.frame(
id = 1:nrow(model$legend),
label = model$legend,
@@ -670,35 +644,22 @@ server <- function(input, output, session) {
makeBbnGraph(modelList[[.selections$model]])
})
#observe({
# visNetworkProxy("bbnGraphPlot") %>%
# visStabilize(iterations = 10)
#})
getModelName <- function() {
paste0("data/", input$modelSelect, ".xlsx")
}
genPlot <- function(boxPlot, title, paletteLength) {
if (nrow(boxPlot) > 0) {
#print(paste('Palette length', paletteLength))
#palette <- brewer.pal(paletteLength, "Set3")
#palette <- c("red", "sienna3", "plum2", "rosybrown4", "sandybrown", "yellow", "seashell3", "palegreen", "springgreen4", "steelblue", "azure")
names(palette) <- 1:length(palette)
#print(paste("Box plot, colours", nrow(boxPlot), length(colours)))
#cat(colours)
xform <- list(categoryorder = "array",
xform <- list(
categoryorder = "array",
categoryarray = boxPlot[, 1],
zerolinewidth = 10)
#
plot_ly(boxPlot, x = boxPlot[,1], y = ~Range, color = as.character(boxPlot$Group), colors = palette, type = "box") %>%
layout(xaxis = xform, yaxis=list(dtick=0.25, range=c(-1.25, 1.25)), showlegend = FALSE, title = title)
zerolinewidth = 10
)
plotly::plot_ly(boxPlot, x = boxPlot[, 1], y = ~Range, color = as.character(boxPlot$Group), colors = palette, type = "box") %>%
plotly::layout(xaxis = xform, yaxis = list(dtick = 0.25, range = c(-1.25, 1.25)), showlegend = FALSE, title = title)
}
}
@@ -714,29 +675,33 @@ server <- function(input, output, session) {
}
}
output$layer1 <- renderPlotly({
output$layer1 <- plotly::renderPlotly({
prepPlot("ba", "Functional Groups")
})
output$layer2 <- renderPlotly({
output$layer2 <- plotly::renderPlotly({
prepPlot("op", "Ecosystem Processes")
})
output$layer3 <- renderPlotly({
output$layer3 <- plotly::renderPlotly({
prepPlot("es", "Ecosystem Services")
})
isAbsolutePath = function( path ){
if( path == "~" )
return(TRUE);
if( grepl("^~/", path) )
return(TRUE);
if( grepl("^.:(/|\\\\)", path) )
return(TRUE);
if( grepl("^(/|\\\\)", path) )
return(TRUE);
return(FALSE);
isAbsolutePath <- function(path) {
if (path == "~") {
return(TRUE)
}
if (grepl("^~/", path)) {
return(TRUE)
}
if (grepl("^.:(/|\\\\)", path)) {
return(TRUE)
}
if (grepl("^(/|\\\\)", path)) {
return(TRUE)
}
return(FALSE)
}
output$linkBackgroundData <- downloadHandler(
@@ -748,16 +713,19 @@ server <- function(input, output, session) {
)
makeLikelihoods <- function() {
likeliTab <- as.data.frame(
cbind(
.likelihoods$p_es, codeVal = sapply(
.likelihoods$p_es,
codeVal = sapply(
.likelihoods$p_es$code, function(str) {
if (startsWith(str, 'p')) as.numeric(substring(str, 2, nchar(str)))
else as.numeric(substring(str, 3, nchar(str)))
if (startsWith(str, "p")) {
as.numeric(substring(str, 2, nchar(str)))
} else {
as.numeric(substring(str, 3, nchar(str)))
}
)),
}
)
),
stringsAsFactors = FALSE
)
@@ -780,7 +748,6 @@ server <- function(input, output, session) {
max = likeliTab$range[elementRow + 6]
)
outputTab <- rbind(outputTab, tabRow)
}
likelihoods <- data.frame(
@@ -798,10 +765,10 @@ server <- function(input, output, session) {
}
output$download <- downloadHandler(
filename = function() { paste0("MESO-", format(Sys.time(), "%m%d_%H%M"), ".xlsx") },
filename = function() {
paste0("MESO-", format(Sys.time(), "%m%d_%H%M"), ".xlsx")
},
content = function(file) {
showModal(
modalDialog(
fluidRow(
@@ -817,8 +784,6 @@ server <- function(input, output, session) {
dir.create(tmp)
setwd(tmp)
l <- list(
pressures = .selections$pressStatus,
nodes = modelList[[.selections$model]]$p_es$nodes,
@@ -828,12 +793,8 @@ server <- function(input, output, session) {
)
xl <- write.xlsx(l, "dataset.xlsx")
#zipFile <- zipr(file, c("dataset.xlsx"))
file.copy("dataset.xlsx", file)
#print(paste("zip file complete", zipFile))
setwd(oldDir)
unlink(tmp)
@@ -841,8 +802,6 @@ server <- function(input, output, session) {
},
contentType = "application/xlsx"
)
}
shinyApp(ui, server)

View File

@@ -1,16 +1,21 @@
# R script to upload the existing spreadsheets and homologise them
library(magrittr)
modules::import(magrittr)
fList <- list.files("data", pattern = "*.xlsx")
# Objective to create data tables with
linkCheck <- function(nodeType, nodeString, nodeStringCheck) {
nodeString <- stringr::str_replace_all(nodeString, "\\.", " ")
res <- sapply(nodeString, match, nodeStringCheck$Nodes) %>% is.na() %>% which()
res <- sapply(nodeString, match, nodeStringCheck$Nodes) %>%
is.na() %>%
which()
if (length(res) > 0) print(paste("Clean up error found in", nodeType, "mapping at", names(res)))
}
getNodeVals <- function(nodeStr) {
params <- stringr::str_split(nodeStr, ",") %>% unlist() %>% trimws()
params <- stringr::str_split(nodeStr, ",") %>%
unlist() %>%
trimws()
paramVals <- stringr::str_split(params, "=")
vals <- c()
lapply(paramVals, function(l) {
@@ -32,7 +37,9 @@ getNodeVals <- function(nodeStr) {
sheetNames <- c("TestScenario", "Map_P_BA", "Map_BA_OP", "Map_OP_ES", "Legend")
cleanNames <- function(namVec) {
stringr::str_replace_all(namVec, "\\.", " ") %>% trimws() %>% tolower()
stringr::str_replace_all(namVec, "\\.", " ") %>%
trimws() %>%
tolower()
}
nodeTable <- tibble::tibble()
@@ -40,7 +47,7 @@ nodeTable <- tibble::tibble()
for (wbIdx in 1:length(fList)) {
wb <- openxlsx::loadWorkbook(paste0("data/", fList[wbIdx]))
hab <- stringr::str_split(fList[wbIdx], "\\.")[[1]][1]
#get pressure names
# Drop the time column no use at all....
sheet <- openxlsx::readWorkbook(wb, sheet = sheetNames[1])[, -1]
@@ -67,6 +74,7 @@ for (wbIdx in 1:length(fList)) {
# linkCheck("bioassemblages", ba, ba_check)
sheet <- openxlsx::readWorkbook(wb, sheet = sheetNames[4])[, -1]
op_check <- na.omit(sheet[, 1:2])
sheet2 <- na.omit(sheet[, -c(1, 2)])
@@ -76,6 +84,7 @@ for (wbIdx in 1:length(fList)) {
# linkCheck("outputprocesses", op, op_check)
legend <- openxlsx::readWorkbook(wb, sheet = sheetNames[5])
nodeType <- c(
@@ -85,8 +94,6 @@ for (wbIdx in 1:length(fList)) {
rep("ecosystemservice", length(es))
)
res <- t(sapply(es_nodes[1, ], getNodeVals)) %>% as.data.frame()
names(res) <- cleanNames(names(res))
res <- res %>% mutate(nodeName = names(res))
@@ -98,14 +105,9 @@ for (wbIdx in 1:length(fList)) {
res
)
)
}
mapNewNames <- function() {
newNameMap <- openxlsx::read.xlsx("MBA_MESO_Nodes.xlsx") %>%
dplyr::select(hab, nodeType, Suggestion, node, newname)
save(newNameMap, file="nameMap.RData")
}

View File

@@ -1,9 +1,7 @@
modules::import(magrittr)
reWeightLayer <- function(nestedLayerTib, fudge = 1) {
for (idx in 1:nrow(nestedLayerTib)) {
#print(nestedLayerTib$data[idx])
thisData <- nestedLayerTib$data[idx][[1]]
# Calculate the overall depletion rate
@@ -12,45 +10,43 @@ reWeightLayer <- function(nestedLayerTib, fudge=1) {
survived <- 1
grown <- 1
for (depIdx in 1:nrow(thisData)) {
if (thisData$values[depIdx]<0) survived <- survived * (1 + thisData$values[depIdx]) else
if (thisData$values[depIdx] < 0) {
survived <- survived * (1 + thisData$values[depIdx])
} else {
grown <- (1 - thisData$values[depIdx]) * grown
}
}
# Update the edge weightings to reflect the combined depletion on the BA from each of the edges
effDepRate <- survived - 1
effGrowthRate <- 1 - grown
#print(effDepRate)
if (sum(thisData$values)==0) newValues <- rep(0, length(thisData$values)) else
if (sum(thisData$values) == 0) {
newValues <- rep(0, length(thisData$values))
} else {
newValues <- round(thisData$values / sum(thisData$values) * (effDepRate + effGrowthRate), digits = 3)
#print(paste(idx, paste(newValues, collapse=",")))
}
nestedLayerTib$data[idx][[1]]$values <- newValues / fudge
}
return(nestedLayerTib %>% tidyr::unnest(cols = c(data)))
}
assignWeights <- function(
edgesTib,
incode,
outcode,
value) {
assignWeights <- function(edgesTib, incode, outcode, value) {
for (idx in 1:length(incode)) {
ref <- intersect(which(edgesTib$input == incode[idx]),
which(edgesTib$output == outcode[idx]))
utils::str(ref)
ref <- intersect(
which(edgesTib$input == incode[idx]),
which(edgesTib$output == outcode[idx])
)
if (length(ref) > 1) stop("Error has occurred with multiple edges between two nodes")
print(paste(ref, edgesTib$values[ref], value[idx]))
edgesTib$values[ref] <- value[idx]
#Set the appropriate values
}
return(edgesTib)
}
reWeightModel <- function(thisNet, pressStatus) {
print("About to recalc p - ba")
# what is the depletion factor for each of the pressures applied to the BA?
@@ -60,9 +56,6 @@ reWeightModel <- function(thisNet, pressStatus) {
dplyr::left_join(thisNet$edges, by = c("code" = "input")) %>%
dplyr::mutate(values = values * 0.9)
print("before")
print(sum(p_on$values))
p_on <- p_on %>%
dplyr::rename(presscode = code) %>%
dplyr::rename(ba_code = output) %>%
@@ -72,7 +65,6 @@ reWeightModel <- function(thisNet, pressStatus) {
newP <- reWeightLayer(p_on, fudge = 1)
print("About to recalc ba - op")
# Repeat for the linkage between ba and op
@@ -88,6 +80,7 @@ reWeightModel <- function(thisNet, pressStatus) {
newBA <- reWeightLayer(ba_impacted, fudge = 4)
print("About to recalc op - es")
# Repeat for the linkage between op and es
@@ -103,10 +96,10 @@ reWeightModel <- function(thisNet, pressStatus) {
newOP <- reWeightLayer(op_impacted, fudge = 2)
# Check for any more links through the system
print("About to recalc es - es")
ess <- unique(newOP$es_code)
es_impacted <- thisNet$nodes %>%
dplyr::filter(code %in% ess) %>%
@@ -125,8 +118,8 @@ reWeightModel <- function(thisNet, pressStatus) {
thisNet$edges <- assignWeights(thisNet$edges, incode, outcode, value)
print("exitting reweighting process")
return(thisNet)
}