An application of Bayesian Belief Networks to assess management scenarios for aquaculture in a complex tropical lake system in Indonesia

A Bayesian Belief Network, validated using past observational data, is applied to conceptualize the ecological response of Lake Maninjau, a tropical lake ecosystem in Indonesia, to tilapia cage farms operating on the lake and to quantify its impacts to assist decision making. The model captures ecosystem services trade-offs between cage farming and native fish loss. It is used to appraise options for lake management related to the minimization of the impacts of the cage farms. The constructed model overcomes difficulties with limited data availability to illustrate the complex physical and biogeochemical interactions contributing to triggering mass fish kills due to upwelling and the loss in the production of native fish related to the operation of cage farming. The model highlights existing information gaps in the research related to the management of the farms in the study area, which is applicable to other tropical lakes in general. Model results suggest that internal phosphorous loading (IPL) should be recognized as one of the primary targets of the deep eutrophic tropical lake restoration efforts. Theoretical and practical contributions of the model and model expansions are discussed. Short- and longer-term actions to contribute to a more sustainable management are recommended and include epilimnion aeration and sediment capping.


Introduction
Inland water cage culture farms (IWCCF) operate at a range of spatial scales, from traditional households to corporate modern farms. IWCCF is rapidly growing around the globe as an answer to the increasing demand for aquatic products [1]. In Asia, cage culture is one of the most widely practised aquaculture systems in terms of both production and value [2,3]. The rapid development of cage culture in the region is often attributed to the necessity to provide fish as an affordable nutritional source and to generate employment opportunity, both in support of poverty alleviation efforts [4][5][6]. Despite their benefits, there are increasing concerns about the negative impacts of IWCCF on the aquatic environment. One of the widely identified effects is eutrophication [2,7,8]. Eutrophication leads to blooming algae and the degradation process of the algae consumes dissolved oxygen (DO). This process then leads to DO depletion causing mass fish kills (MFK) [9,10]. Other impacts are changes in macro and microbenthic population, and effects on fish biodiversity [11,12].
MFK events have been long reported from various lakes, rivers, and reservoirs in Indonesia and many other countries (see [13][14][15][16][17]) and they have caused huge financial loss [18]. This raises the need to provide a model that can predict MFK to mitigate the loss. However, data limitations can sometimes hinder the development of mechanistic models. Thus, engaging with approaches such as Bayesian Belief Networks (BBN) becomes an interesting alternative. Two papers have used BBN to predict algal blooms related to several nutrient input scenarios (see [19,20]). Both studies were conducted in temperate lakes in high income countries where monitoring data of nutrient inputs are reasonably adequate. They use a relatively low number of nodes (maximum of nine) to account for the principal drivers of algal blooms. However, we are not aware of any studies that use BBN to predict the MFK in deep tropical lakes with complex and specific ecological characteristics [21,22].
In this study, we apply BBN to generate information on the ecological response to IWCCF that can underpin appraisals of IWCCF management scenarios of Lake Maninjau, Indonesia. We specifically aim to answer the research question: to what extent can BBN be used to conceptualize the lake's complex ecological system related to IWCCF management and quantify its impacts to assist decision making? We surmise that BBN is particularly suited to investigate the impacts of IWCCF on the occurrence of MFK incidents due to upwelling processes in locations where data is scarce. The developed model can then help improve our understanding of both the economy-environment trade-offs of IWCCF expansion and the effects of IWCCF and upwelling processes on a native fish species, Gobiopterus sp. This species is the only finfish species reported by the fishers to diminish during and after major upwelling events; thus, it is chosen as the indicator of IWCCF impacts on the native biota.
This study differs from previous BBN applications [19,20] because it includes more complex relationships between physical, biological, and chemical components of the lake ecology, in addition to human activities-IWCCF management. Thus, it better illustrates the complexity of the ecosystem. Moreover, this study emphasizes the advantage of BBN in overcoming limited data, which is a common problem in the Global South. Through the development and testing of the BBN in our case study, we are able to identify important gaps in data in the studied lake and beyond that should be addressed to improve model accuracy and facilitate development of mechanistic models. We apply the BBN to investigate a range of management intervention scenarios. This provides important information to inform sustainable lake management.
Beyond the applied perspective, the paper makes several contributions to the academic literature. First, we develop an integrated model to predict the probability of MFK incidents in a complex tropical lake system. To the best of our knowledge, this is the first such modelling effort. In line with this, [23] identified that BBN application in tropical region is limited. Second, this paper contributes to the growing literature on the use of BBN to assess ecosystem service (ES) trade-offs; this paper specifically investigates trade-offs between two provisioning services: fish provided by aquaculture and a native species capture fishery. There are few BBN case studies that aim to simulate ES trade-offs especially in complex-spatial systems [23][24][25]. More importantly, this study rigorously applies a variable correlation table to build the direct acyclic diagram (DAG), which is rarely conducted by previous BBN studies [26]. The applied and methodological contributions highlighted above carry insights that may be used in different lakes across the tropics that are subject to planned or ongoing IWCCF.
The paper is organised as follows: first, a brief overview of the history of IWCCF in Lake Maninjau and of the lake's current ecological condition. An introduction of the BBN approach follows this. The third section will present the methodology applied to build and validate the model. In the fourth section, the developed model is shown and the results of several IWCCF management scenarios are presented. A concluding section then follows a discussion of model performance and scenario results.

Site description
Lake Maninjau (Fig 1) is one of the big volcanic lakes in Sumatra Island, Indonesia. Its total surface area is 97.37 km 2 , and its catchment area is approximately 132.6 km 2 with a full water volume of 10 million m 3 [27]. Its geographical positions are 100 o 08'5384" E to 100 o 14'02.39" E and 0 o 14' 52.50" S to 0 o 24' 12.17" S.
It is selected as the case study as it represents many Indonesian lakes formed by a volcanic and or tectonic processes in complex tropical systems [28]. Because of its geological formation, the lake has a specific character. It is naturally rich in sulphur, which may be deposited in sediment [29] and drives its ecological process [10,29]. Therefore, it has a natural sulphur purification activity locally named as tubo belerang [30]. Thus, these type of lakes are generally considered to be inappropriate for IWCCF practices from an ecological perspective [31].
The sulphur cycle interrelated with organic material is the primary key to balance this lake's ecosystem and should be preserved to prevent severe consequences such as MFK. The balance is produced by an oxidation process depending on oxygen level and phototropic sulphur bacteria activities. Sulphide, which is a reductive form of useful sulphate or elemental sulphur, actually can be easily oxidized to elemental sulphur if there is oxygen or the mentioned

PLOS ONE
bacterial activity [29]. In sufficient DO condition, the sulphide will react with reactive iron ions (Fe 2+ ) creating FeS compounds. FeS will the react to form pyrite sulphide (FeS 2 ), a stable structure which can sink to the sediment. However, in the condition when Fe 2+ is limited, more hydrogen sulphide (H 2 S) is unbounded and accumulated in pore waters [32]. Further, reactive Fe (Fe 2+ ) is oxidised into Fe 3+ , which will bind phosphate ions (PO 4 3-) to form FePO 4 compounds. However, [10] and [29] elucidated that high H 2 S can be produced in the hypolimnion if there is a high rate of organic material degradation. In this circumstance, limitation of free reactive Fe will boost reduction of Fe 3+ to Fe 2+ to form FeS, but it will unbind the phosphate ions with the help of sulphate reduction microbes. As a result, the phosphate ions are released back into the hypolimnion, which is then diffused to the metalimnion and epilimnion, where the phosphate ions contribute to eutrophication.
As for many lakes in Indonesia (and in the Global South in general), there is limited (time series) data available for Lake Maninjau concerning important biogeochemical indicators and general lake condition. The Research Centre for Limnology, Indonesian Institute of Sciences (RCL-IIOS), has collected data on selected limnological aspects over the years. However, the data was sporadically collected due to financial constraints and is therefore patchy [33]. Because these data gaps exist, a mechanistic model illustrating the connection between IWCCF management and its results on the water quality improvement has not yet been developed. We consider a BBN to be a suitable alternative to fill this gap.
Biologically, the lake is the habitat for 18 types of finfish, which become the target of capture fishery [ However, IWCCF operation in the lake has also been reported to be the cause of a severe decrease in water quality, which is indicated by eutrophication and more frequent MFK incidents [15,[39][40][41]]. An MFK taking place in 2018 was estimated to result in a financial loss in excess of 8.7 million USD [42].
Excessive nutrient input from accumulated feed waste correlated with the use of DO for decomposition is further related to the DO depletion. This process is associated with widening anoxic layer which becomes the catalyst for more frequent MFK events. The thermocline was observed in the depth of 4-20 m [43]. They recognized that the anoxic layer was found at this depth. As a consequence, a strong wind can trigger the mixing and bring the chemical substances in the anoxic layer to surface layer. The combination of DO deprivation and H 2 S released during the upwelling process are lethal to the fish [10]. Furthermore, data collected by the authors in 2019 also shows that Gobiopterus sp. could not be found for months after the MFK incidents [44]. The implication of this situation is income loss for local fishers and forgoing one of the affordable protein sources for the locals.

Bayesian Belief Networks (BBN)
A BBN is a graphical interface presenting a joint probability of various variables which can link actions to outcomes [45][46][47][48]. Mostly, it is used in Decision Support System (DSS) implementing Bayes' rule on probability (Eq 1) [49]. Increasing applications in environmental risk assessment is also documented [50-52]. Other applications as a tool for understanding complex (socio-ecological) systems, habitat suitability, management evaluation, and modelling ES are well reported [23, 25, 53-57]. Further, BBN are also increasingly applied in environmental studies to connect management actions and environmental impacts in complex environments using interdisciplinary knowledge (i.e., [57-59]).

P HjE
Note: P(H) = the probability of event H will be in a particular state P(EjH) = the conditional probability of E given H P(HjE) = the updated probability of event H will be in a particular state given that evidence E is happening The use of BBN as a modelling tool in interdisciplinary research is related to its flexibility to effectively combine various multidisciplinary perspectives [48]. Further, BBN as a modelling tool allows connecting a data learning process and expert knowledge. This flexibility makes it very useful to address research questions considering existing gaps in empirical data [46,60].
BBN may encourage communication between scientists and decision-makers, who may have a limited background in the specific science knowledge [61]. Another merit is the possibility to relate BBNs to a participatory decision-making process [62, 63] to combine both empirical data and expert knowledge to address challenges associated with data and evidence gaps [45,64]. However, [65] advised to carefully validate the constructed model if it is built based on expert domains solely.
Despite the advantages of BBN, one of the limitations is that they cannot incorporate and model feedback loops [66]. Another limitation is the difficulty to capture dynamic processes [23,60]. Further details about the advantages and limitations of the model are discussed by [24].

Ethical clearance
The research has been granted an ethical clearance from the ethic committee of the School of Geosciences, University of Edinburgh.

Conceptual model development
The key variables influencing upwelling events and MFK were determined by referring to peer-reviewed literature (S1 Table in S1 File). This approach was also used by [67]. Thirty-five selected variables comprising physical, chemical, biological, and anthropogenic aspects were selected. An initial version of a variable correlation table was created to avoid over-complexity. The collaborators from RCL-IIOS then reviewed the

BBN construction and model parameterisation
The directed acyclic graph (DAG) was used to construct the model in Netica software (Norsys Software Corp, 2019). The states were set based on literature described in Table 1 and limited to not more than five states to keep the model simple. In cases where there were more than five states, such as for the wind speed node and for wind direction, these were reduced to a maximum of five states (see [65]). A loop problem occurred when defining the network. The loop comprised light intensity, water transparency, DO epilimnion, and Chlo-a epilimnion. Actually, light intensity is directly correlated to photosynthesis capacity and reproductive ability of the algae before it is further related to the production of oxygen. As a solution to the problem, light intensity was connected directly to DO epilimnion instead of Chlo-a epilimnion to represent photosynthesis process. A mix of data and expert domains were engaged to fill Conditional Probability Tables (CPT) referring to the work of [64]. The data to set prior probabilities were obtained from Meanwhile, most daily climatic data (rainfall intensity, wind speed, wind direction, cloudy) were collected from Meteorology, Climatology, and Geophysical Agency (BMKG). The data presented in high-frequency measurement data (measured every 2 to 10 minutes) were processed into the most frequent occurrence data (i.e., wind direction, wind speed) and maximum value data (rainfall intensity).
Other CPTs (i.e., H 2 S, DO-metalimnion, phosphate released from sediment) were set by using expert knowledge referring to peer-reviewed publications and reports from Government Institutions. They were calibrated after conducting repeated sensitivity to parameters tests (to mixing, MFK, and Gobiopterus disappearance). Schmidt Stability Index (SSI), which shows the lake's resistance to mechanical mixing [68] was calculated by using R-lake analyser from highfrequency measurement data from RCL-IIOS (2014-2017). The trend of SSI related to light intensity was used as the basis to fill the CPT. Detail of the method used to fill CPT of each node is presented in Table 1 and S3 Table in S1 File.

Model evaluation
Sensitivity to findings. Sensitivity to findings helps to identify errors in either the network or the CPTs [46, 64]. First, the analysis was conducted to MFK and Gobiopterus disappearance node. The analysis results show that mixing and epilimnion zone after mixing nodes are the most influencing factors for both output modes. These results are consistent with the model as these two nodes are the direct nodes connected to the output nodes. Next, the sensitivity to findings was applied to find the most influenced factors to mixing and epilimnion zone after mixing. The sensitivity to findings is illustrated in Mutual Information Value. Mutual Information Value is the subtraction of entropy/ H(X) and the probability of distribution of the entropy given prior condition (H (XjY)). This value is used to rank the most influencing factors to the output nodes, please refer to [64] for details.
Model accuracy evaluation. A goodness of fit procedure using past field observation data of MFK events was used to measure model accuracy [46, 60, 64]. Three indicators used to show the model accuracy is the error rate, true skill statistics, and Area Under Curve (AUC). Unfortunately, there is no record data of the abundance of Gobiopterus; thus, model accuracy evaluation cannot be performed to 'Gobiopterus disappearance' node.
The past observation data include wind speed, wind direction data, and the dates when MFK events happened. The data were obtained from mixed data from RCL-IIOS online monitoring device and observation data from Agam Fisheries Agency. The data from the RCL-IIOS were recorded in a high-frequency observation from mid of July to early August 2019, December 2019, and January-February 2020 (n = 8,331 data sets). These data were then processed into frequent wind speed and wind direction data resulting in 100 observation points. Meanwhile, among 20 MFK event dates only nine dates can be used, because there were no climatological data for the other 11 dates. The data is then inputted to the model to obtain the prediction of MFK probabilities, coded as predicted MFK. Meanwhile, the actual MFK events were coded as observed MFK. The predicted MFK probability was coded as Yes = 1 or No = 0 using various cut-off points (1% to 96%). Then the predicted probability distribution is compared with actual MFK event occurrence to calculate the sensitivity and specificity of the model. The obtained sensitivity and 1-specificity were plotted to form a receiving operating curve (ROC) (see [67,79]). A nonparametric ROC was generated by using Stata 15 software following the methodology proposed by [80]. A Non-parametric ROC curve was selected as it does not require assumptions of a normal distribution of specificity and sensitivity [81]. The area under the ROC curve (AUC) is used as the indicator of model accuracy.
The next step is the calculation of an error rate following [65]. First, an optimum cut-off point should be obtained to code the projected probabilities into MFK predicted occurrence (Yes = 1 and No = 0). The optimum cut-off point was determined by engaging the minimum P-value approach [82,83]. The optimum cut-off value is the cut point resulting in maximum χ 2 value. Then, using this cut-off point, a 2x2 cross table of predicted MFK and observed MFK was produced. Finally, the error rate is calculated by dividing the incorrect prediction with total observation (n = 100). The true skill statistic (TSS) was also calculated at the obtained cut-off point. This value of TSS was obtained using the equation suggested in [84,85]: TSS = Sensitivity +specifity-1.

Management scenarios
Several management scenarios are presented. These scenarios are grouped as short-term (immediate response) and longer-term (up to 25 years) management actions. The short-term management options include with and without complete removal of IWCCF to test the implications of complete removal of IWCCF such as suggested by the environmentalists, reduction IWCCF of up to 6,000 cages (proposed by the Act of Bupati Agam No. 22 /2009), and reduction of internal phosphorous loading (IPL) (i.e., by sediment capping) and aeration recommended by aquaculture practitioners. All the short-term scenarios are generated in the model setting representing the current lake condition (i.e., high H 2 S and high IPL).
For longer-term management options, a continuation of current IWCFF situation and present lake condition as the result of no management action, complete removal of IWCCF, and reduction IWCCF up to 6,000 will be the proposed actions. A timeframe of maximum 25 years was chosen as this corresponds to the water retention time of Lake Maninjau [86]. The best management scenarios are chosen based on their ability to reduce the probabilities of MFK and Gobiopterus disappearance, while also considering social implications of the scenarios.

Constructed BBN
The built labelled box style of the model comprises physical, chemical, biological, and human dimension (input) nodes related to IWCCF and lake ecosystem is shown in Fig 2A. An alternative form of model presentation, the belief bar style (the model completed with the CPT tables), is shown in Fig 2B. The first three nodes are regarded as intermediate nodes, and the human dimension is further referred to as input nodes. Details of the variables can be accessed in Table 1. Interconnections of these parameters result in two targeted outputs (IWCCF Fish MFK and Gobiopterus disappearance). These two outputs represent the criteria of interest for informing decision-makers about trade-offs between IWCCF and the abundance of Gobiopterus. Fig 3 presents the sensitivity to findings (presented as Mutual Information Value) of four most influential parameters to changes in the nodes "mixing" and "epilimnion zone after mixing". The Figure also shows the results of sensitivity to findings to the output nodes. As it has been explained in Sensitivity to Findings Section, "mixing" and "epilimnion zone after mixing" nodes are the most influential factors of both output nodes. Thus, these two intermediate nodes reflect the influence of other parameters towards the output nodes. In the Figure, the horizontal axis shows the mutual information value, which is a measurement of influence of one node on other nodes [64]. The sensitivity to findings of changes in other factors are not presented here due to space limitations.

Results of sensitivity to findings
It can be seen in the Figure that wind speed, water current velocity, SSI, and rainfall intensity are the most affecting parameters. Thus, it can be seen that physical factors especially wind speed, are the principal drivers of MFK and Gobiopterus disappearance events. Further, it is also suggested that the output nodes are insensitive to the current lake condition, controlling external phosphate loading by regulating the number of cages, stocking density, and feed management (input nodes). Please refer to discussion section to read the detail explanation.

Results of model accuracy evaluation
An optimum cut-off value is presented to balance between the sensitivity (true positive) and specificity (true negative) results. An optimum cut-off value of 26% was calculated ( Table 2). A high cut-off value (e.g. 50%) will increase the possibility of true negative results but will reduce the possibility of true positive results. On the contrary, if the cut-off point is too low (e.g. 5%), it will increase true positive results but reduce the true negative results. Imbalance between sensitivity and specificity will eventually raise the model error rates.
The non-parametric ROC is presented in Fig 4, showing the AUC is 0.9. The AUC value shows that the constructed BBN model is useful [87] with moderate to excellent discrimination [79,88]. The computed error of the model at the optimum cut-off point is 4% consisting of four wrong predictions among 100 tested data points ( Table 3). The information in the table shows that there were two predicted MFKs that were not observed in actual observations. At Table 2. Sensitivity, specificity, calculated χ 2 , and P value of several cut-off points.  the same time, there were also two MFK events that were not predicted by the model. This result is relatively low compared to other BBN models (i.e., 9.62% for wildfire occurrence [67]). However, the categorization of high and low error rate is somewhat different for each subject of study. The calculated TSS is 0.796. The TSS value has a range from -1 to 1, with the value below zero indicating that model performance is based on chance alone (random) [84,85]. This means that the model can produce strong prediction of MFK with the available data and the correct prediction are not based on random incidence; thus, it can be considered to be statistically robust.  Table 4, and the details can be seen in S4 Table in S1 File. The scenarios are categorized as the short-term impacts of the available management options (Scenario 1a-1f) and the longer-term effects (Scenario 2a-2d). The value in the table shows predicted probability which is translated to yes (the events happen) and no (the events do not happen) by using the calculated optimum cut-off value (see Model Accuracy Evaluation Section). As an example: for Scenario 1a, the predicted probability of MFK in good weather condition is 1.8-17%; this means that no MFK event is predicted as the value is lower than 26%. At the same time, in bad weather condition, the predicted probability is 32-85.5%, which means it is highly likely that MFK will occur.  The predicted probability is shown in two contrasting weather conditions. The first one (good weather condition characterized by dry season, calm wind, no heavy rain) represents the optimum condition when MFK and Gobiopterus disappearance is the least likely to occur. On the contrary, the second condition (poor weather condition associated with wet season, gale-storm wind, heavy rain) represents the highest possibility of MFK and Gobiopterus disappearance to happen.  The short-term scenarios reflect the expected probabilities of occurrence of both MFK events and Gobiopterus disappearance within the near future after the chosen management actions are imposed. The management options can be seen in the management Scenarios Section. The longer-term effect scenarios project the impacts of the management to the expected probabilities of both output nodes up to 25 years. This time frame is chosen because it is the water retention time of the lake water [86].

Results of predicted probability of MFK and Gobiopterus disappearance for several management scenarios
We set the 'PO 4 released from sediment' and 'Anoxic layer' to the lower state to represent the expected condition in longer term after the reduction and eradication of the cages making the fulfilment of the assumption become vital.
There are two assumptions used to run scenario 2b (longer-term effects of complete removal of IWCCF) and scenario 2c (longer-term effects of keeping 6,000 IWCCF). The first assumption is that there are no other factors inhibiting lake self-purification (i.e., continuing IPL and absence of stable macrophytes communities) (see [89,90]). The second assumption is that the rate of external organic sediment run-off (i.e., domestic waste) is unchanged.
The results show that in the short-term, removing IWCCF will not have a significant impact on Gobiopterus disappearance, although the losses due to MFK can be mitigated (Scenario 1a vs 1b). Meanwhile, keeping 6,000 IWCCF offers similar results as the business-as-usual scenario (Scenario 1a vs 1c). More importantly, combating the IPL is expected to improve the lake condition and may reduce the probability of MFK and Gobiopterus disappearance considerably (Scenario 1a vs 1d). Overall, the suggested best short-term management option is by conducting aeration alone or combine it with a reduction of IPL (Scenario 1e and 1f).
For the longer-term scenarios, removing all IWCCF may become the best option to avoid MFK events and Gobiopterus disappearance as long as the assumptions above are fulfilled (Scenario 2b). Further, limiting the number of IWCCF to a moderate amount (e.g. 6,000) may cut probability of MFK and Gobiopterus disappearance by more than half, if it is complemented with a reduction of IPL and better feeding management (Scenario 2a vs 2d). Without performing these efforts, a less significant improvement can result from reducing the cages to the 6,000 units (Scenario 2a vs 2c).

Discussion
The constructed BBN model can provide a relatively accurate prediction of MFK when tested with only two explanatory variables (wind speed and wind direction). This finding indicates that the developed model can perform well with limited data availability. The BBN's ability to incorporate both data learning and expert knowledge is the principal benefit in supporting this feature. The model, therefore, can be useful to assist decision-making in the study area.
Further, the model results show that weather conditions, especially wind speed, can be used as a good predictor of the outputs. Therefore, in the current lake condition, forecasting wind speed is crucial to minimize financial loss during MFK by harvesting the fish from cages before stormy events. However, some loss is inevitable in such a strategy, because not all fish will have reached the most profitable size. Interestingly, some local IWCCF farmers have adopted this adaptation strategy in recent years [44]. They harvest the farmed fish one or two days before strong wind is predicted by the weather forecast, and they reduce stocking rates and production during the rainy season.
Our previous study [44] has correlated the farmers' strategy with the prediction of more frequently extreme weather such as La Niña which likely triggers more upwelling in the study area. Thus, it gives us an insight that developing the model to take into account the climate change impacts in the future is necessary aligning with the work of [91]. In the developed model, the output nodes are not sensitive to the input nodes/human dimension nodes. In BBN, the sensitivity of output nodes towards input nodes can be reduced by a large number of intermediate nodes [23,65]. Thus, [65] advised to limit the number of nodes to not exceed 40. The model uses 35 nodes in line with their advice. This would suggest that model over-complexity is unlikely to be the cause of the lack of output nodes sensitivity to the input nodes. One possible explanation is that high rates of IPL and H 2 S mimic the current condition of the lake (see Materials and Methods Section). This shows that even if the IWCCF was eradicated, it will not immediately result in water quality improvement. In line with this finding, [89,90,92] elucidated that chemical and biological resistance factors (i.e., IPL) can delay water quality improvement in shallow lakes. IPL can endure for 10-15 years after the external load is reduced [93,94]. However, [95] mentioned that deep stratified lakes typically have shorter longevity of internal sediment P released than the shallow lakes. Thus, limiting IPL is one of the main actions that should be conducted to achieve desired results in terms of reduced MFK frequency and presence of Gobiopterus. Further, it shows that longer-term result of complete removal of IWCCF will be achieved if other factors inhibiting lake self-purification such as IPL is controlled.
Thus, the necessity to recognize IPL as one of the primary causes of water quality deterioration and the needs of restoration actions to combat it becomes the main lesson learnt from the application of the model. International literature recognizing the importance of acknowledging IPL and the significance to combat it in the deep tropical lakes is scarce. Until recently, IPL control is not included in the priority list in most deep lake restoration programmes in tropical countries (See [44, [96][97][98] for literature on national restoration plans of several lakes in India, Thailand, and Indonesia). At the same time, most related research are case studies in shallow temperate lakes (e.g., [95,[99][100][101]) and in shallow subtropical and tropical lakes (e.g., [102,103]).
In the short-term, conducting aeration may significantly reduce the probability of MFK and Gobiopterus disappearance (Scenario 1e). This action commonly practiced in pond aquaculture by infusing oxygen in the water column is aimed at reducing the probability of fish death (see [104,105]). In 2019, research related to the application of nano-bubble technology epilimnion aeration was conducted in Lake Maninjau. The research was done by a collaborative research team of RCL-IIOS and The Ministry of Environment and Forestry (see [44]). Further, combining epilimnion aeration and other efforts to control IPL is considered as the best option to reduce the probability of MFK and Gobiopterus disappearance events (Scenario 1f). The aeration can be placed in the IWCCF area. In the future, development of spatial BBN is necessary to advise which areas should be prioritized to be aerated.
The significance of spatial BBN has been demonstrated by [106][107][108]. They use the spatial model to map the prioritize management intervention area, to incorporate stakeholders' perspectives, and to map landscape vulnerability and concluded that spatial BBN is a practical tool to achieve their objectives.
Although total cage removal can make the possibility of MFK and Gobiopterus disappearance more unlikely in the longer-term, other considerations such as poverty alleviation and significant potential to create social conflict (see [15,42]) should also be thoroughly assessed. It was reported that conflicts of interest between water users developed into riots in 2001 [103]. Thus, a balance between ecology and economy should be achieved considering these circumstances. This is particularly important if the impacts of actions such as cage removal are not certain to alleviate MFK. One alternative that can be proposed from this model is the integration of efforts to reduce IPL while maintaining a moderate amount of IWCCF (Scenario 2d). Various alternatives of techniques to reduce IPL are available. Two promising methods that can be applied are sediment interface aeration, and sediment capping [109][110][111].
To provide a more comprehensive assessment of both ecology and economy impacts of the scenarios, an endpoint such as a utility node defined in terms of effects on income/livelihood from IWCCF can be a useful amendment to the current model. Towards this end, deriving high quality time series data about the economic aspects of IWCCF and local residents' livelihood will be of high priority for future research. This would allow a more detailed quantification of trade-offs between IWCCF and native fisheries in order to inform management decisions [23]. There is a clear necessity of quantifying ES trade-offs using specific modelling tools such as BBN [27]. In the study area, information of quantification of ES trade-off is very useful to help the decision makers to allocate their limited financial budget on relatively more sustainable management actions. However, monetary valuation of trade-offs between IWCCF and native fisheries is beyond the scope of this paper.
The main limitation of the developed BBN model is the difficulty to integrate dynamic observations to provide the exact point in time when the IPL starts taking place. In this case, mechanistic models are required to complement the BBN model. Unfortunately, for the time being, such models are not available yet in the study area, and their development may be challenging and time consuming due to the complexity of this tropical lake ecosystem. Some studies have demonstrated the use of BBN to complement other models [58, [112][113][114]. These studies highlight that BBN generate results that were more user-friendly and more easily understood by end users and decision makers; nevertheless, mechanistic models can provide more detail and thus also help refining inputs for the BBN.
In this case study, the BBN points to a need for an improved understanding information gaps which is applicable to the lake management in Indonesia and the Global South in general. In short, the information gaps found by application of the constructed model are: 1. Absence of high-quality time series data on the limnological aspects of the lake to create a mechanistic model to complement the built BBN model as the built BBN model has several aforementioned limitations.
2. Lack of observational time series data on the impacts of phosphate reduction on water quality improvement which leads to poor understanding of the underlying mechanism that may delay the effects of phosphorous reduction in tropical lakes.
3. Lack of high-quality time series data about the economic aspects of IWCCF and local residents' livelihood.
4. Need to conduct an economic valuation of the management scenarios to assists the decision-making process.
Other studies aligning with these findings are for example [96] voicing the importance to conduct economic valuation research to support lake management in India. Further, [115] identified that inadequacy of basic limnological monitoring data hindered effective sciencebased policy lake management in the Philippines.

Conclusion
The constructed BBN model can be used to capture the complexity of the ecological system in a deep tropical lake. It was constructed by using both data and expert and stakeholder knowledge approaches. However, careful model parameterisation and validation is required when engaging in this approach. It overcomes limited data availability to overcome limitations regarding data availability that often characterise such systems and is shown to be accurate when tested with past field data observations (an error rate of 4%, TSS of 0.8, and AUC of 0.9).
The model can assist decision making by highlighting the importance of wind speed forecasts to reduce the loss of MFK. Other important factors that should be observed are related to lake stability (i.e., SSI and rainfall intensity). Further, it reflects that reducing external phosphorus loading (in this case by reducing the number of cages used for fish farming on the lake) may not give expected results within a short time and also in the longer-term, if IPL is not controlled. Thus, the findings underline the necessity of actions to combat IPL in deep tropical lake systems.
In the short-term, combining control of internal phosphate released with aeration is considered the best option to reduce the probabilities of MFK and Gobiopterus disappearance. In the longer-term, maintaining a moderate amount of the cages along with efforts to control IPL may become an alternative solution to balance the lake's ecological condition and to conserve economic benefits. Together, this may also help alleviating social conflicts in the study area and elsewhere.
The use of BBN does not undermine the importance of availability of adequate data to create a mechanistic model. Rather, it should be complemented by such model to provide more comprehensive decision making. Further, the identified data and knowledge gaps are relevant to improving lake management in the tropical system in general.