# LUNTIAN: An Agent-Based Model of an Industrial Tree Plantation for Promoting Sustainable Harvesting in the Philippines

Zenith Arnejo, Benoit Gaudou, Mehdi Saqalli, Nathaniel Bantayan

► To cite this version:

Zenith Arnejo, Benoit Gaudou, Mehdi Saqalli, Nathaniel Bantayan. LUNTIAN: An Agent-Based Model of an Industrial Tree Plantation for Promoting Sustainable Harvesting in the Philippines. Forests, 2025, 10.3390/f16081293. hal-05203822

HAL Id: hal-05203822

https://hal.science/hal-05203822v1

Submitted on 1 Oct 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
Article

# LUNTIAN: An Agent-Based Model of an Industrial Tree Plantation for Promoting Sustainable Harvesting in the Philippines

Zenith Arnejo<sup>1,2,*</sup> ID, Benoit Gaudou<sup>1</sup> ID, Mehdi Saqalli<sup>3</sup> ID and Nathaniel Bantayan<sup>4</sup> ID

## Abstract

Industrial tree plantations (ITPs) are increasingly recognized as a sustainable response to deforestation and the decline in native wood resources in the Philippines. This study presents LUNTIAN (Labor, University, Timber Investment, and Agent-based Nexus), an agent-based model that simulates an experimental ITP operation within a mountain forest managed by University of the Philippines Los Baños. The model integrates biophysical processes—such as tree growth, hydrology, and stand dynamics—with socio-economic components such as investment decision making based on risk preferences, employment allocation influenced by local labor availability, and informal harvesting behavior driven by job scarcity. These are complemented by institutional enforcement mechanisms such as forest patrolling, reflecting the complex interplay between financial incentives and rule compliance. To assess the model's validity, its outputs were compared to those of the 3PG forest growth model, with results demonstrating alignment in growth trends and spatial distributions, thereby supporting LUNTIAN's potential to represent key ecological dynamics. Sensitivity analysis identified investor earnings share and community member count as significant factors influencing net earnings and management costs. Parameter calibration using the Non-dominated Sorting Genetic Algorithm yielded an optimal configuration that ensured profitability for resource managers, investors, and community-hired laborers while minimizing unauthorized independent harvesting. Notably, even with continuous harvesting during a 17-year rotation, the final tree population increased by 55%. These findings illustrate the potential of LUNTIAN to support the exploration of sustainable ITP management strategies in the Philippines by offering a robust framework for analyzing complex social-ecological interactions.

**Keywords:** agent-based modeling; industrial tree plantation; sustainability; social-ecological systems; Philippines

## 1. Introduction

Industrial tree plantations (ITPs) in the Philippines represent a strategic response to decades of deforestation and the decline in the native wood industry [1]. Historically, the country was home to over 20 million hectares of dipterocarp forests, but aggressive post-war logging and agricultural expansion rapidly diminished this vast resource, leading to

<figure><img src="image_2.png" /></figure>

check for
updates

Academic Editor: Stefano Grigolato

Received: 25 June 2025

Revised: 30 July 2025

Accepted: 4 August 2025

Published: 8 August 2025

**Citation:** Arnejo, Z.; Gaudou, B.; Saqalli, M.; Bantayan, N. LUNTIAN: An Agent-Based Model of an Industrial Tree Plantation for Promoting Sustainable Harvesting in the Philippines. *Forests* 2025, 16, 1293. https://doi.org/10.3390/f16081293

**Copyright:** © 2025 by the authors.
Licensee MDPI, Basel, Switzerland.
This article is an open access article
distributed under the terms and
conditions of the Creative Commons
Attribution (CC BY) license
(https://creativecommons.org/
licenses/by/4.0/).

<footnote n="1">1 UMR 5505, Institut de Recherche en Informatique de Toulouse (IRIT), Université Toulouse Capitole, 31000 Toulouse, France; benoit.gaudou@ut-capitole.fr</footnote>

<footnote n="2">2 Institute of Computer Science, University of the Philippines, Los Baños 4031, Philippines</footnote>

<footnote n="3">3 UMR 5602, Géographie de l'Environnement (GEODE), Centre National de la Recherche Scientifique (CNRS), Université Toulouse Jean Jaures, 31000 Toulouse, France; mehdi.saqalli@univ-tlse2.fr</footnote>

<footnote n="4">4 Institute of Renewable Natural Resources, University of the Philippines, Los Baños 4031, Philippines; ncbantayan@up.edu.ph</footnote>

<footnote>* Correspondence: zoarnejo@up.edu.ph</footnote>
significant shifts in the wood industry. Once flourishing during a "golden era" in the 1980s, the sector later faced strict regulatory restrictions and bans on natural forest extraction, prompting a pivot toward ITPs as a sustainable alternative [2].

Established under Presidential Decree No. 705, ITPs were designed to promote refor-
estation through long-term leases and various incentives, including low lease fees, rental
exemptions, and tax breaks for qualified developers [3]. This policy-driven shift was fur-
ther supported by initiatives from agencies such as the Philippine Council for Agriculture,
Aquatic and Natural Resources Research and Development, which advanced ITPs as a mod-
ern approach to rehabilitating public timberlands. Additionally, government measures such
as DENR Administrative Orders Nos. 01 and 41, Executive Order No. 725, and Proclama-
tion No. 2351 aimed to streamline lease agreements, reforest degraded lands, and convert
understocked areas into productive plantations, ensuring the long-term sustainability of
the country's wood industry [2].

Despite policy innovations, the ITP sector still faces persistent socio-economic, en-
vironmental, and management issues such as inequitable land tenure, resource scarcity,
limited labor opportunities, and complex plantation ecosystem management [4]. These
issues—along with unbalanced cost-benefit sharing and unplanned land use—continue to
hindер sustainable development, while environmental degradation adds further pressure
onto reforestation efforts. In addition, management challenges stemming from inadequate
design, poor policy alignment, and misdirected incentives often result in non-technical
failures, such as conflicts, political backslashes, and labor issues. Labor intensity remains a
critical factor, and the economic viability of slower-growing native species is frequently
compromised when pitted against fast-growing alternatives and even illegally harvested
natural forest wood [4]. As the Philippines strives to balance economic development with
environmental sustainability, understanding the multifaceted challenges of ITPs is essential
to shaping future policies and practices. While ITPs are sometimes critiqued for potentially
undermining environmental values, they also offer substantial social benefits—ranging
from increased employment to improved rural incomes. These benefits are particularly vital
as the nation continues to address persistent timber shortages through a mix of domestic
innovation and reliance on imported resources.

To address these multifaceted challenges, we ask the following: *Can industrial tree plantations be designed and managed in ways that balance economic viability, ecological restoration, and social inclusion?* This question demands a system-level approach that integrates human and ecological dynamics over time. Viewing an ITP as a social-ecological system (SES) allows us to consider the feedback among stakeholder decisions, forest conditions, and institutional mechanisms. In this context, a simulation-based approach—particularly agent-based modeling—is well suited for testing management scenarios under realistic constraints by representing diverse stakeholders (e.g., investors, laborers, and managers) as autonomous agents whose behaviors are shaped by rules, preferences, and dynamic environmental conditions. To this end, we introduce LUNTIAN, an agent-based model developed to simulate an experimental ITP within a mountain forest managed by the University of the Philippines Los Baños Land Grant Management Office. Framed as an SES, LUNTIAN seeks to strike a balance between two primary objectives: preserving forest cover and ensuring the long-term profitability of ITP operations for both investors and local communities.

Unlike established forest models such as 3PG [5] and FORMIND [6], which concentrate almost exclusively on biophysical processes such as tree growth, competition, and stand dynamics, LUNTIAN places socio-economic dynamics at the core of its simulation. It actively models labor decisions, informal harvesting behavior, and investment processes—critical human factors that traditional forestry simulations often treat as external or overlook entirely. While spatial forest-landscape models such as LANDIS [7] have made strides in
incorporating management-driven change by simulating disturbance regimes and policy interventions, they typically treat socio-economic factors as external inputs or future integration points rather than intrinsic system components. LUNTIAN distinguishes itself by employing an agent-based architecture that directly simulates the decisions, behaviors, and feedback loops of key human actors—including investors, laborers, and enforcement agents—within the forest ecosystem. This design embeds socio-economic processes as internal, dynamic elements of the system, enabling LUNTIAN to evaluate sustainability in a more holistic and integrated manner across ecological, economic, and social dimensions.

The remainder of the paper is organized as follows: In Section 2, we detail the methodological framework used to develop and assess the LUNTIAN model, including the case-study context, model design, and assessment procedures. In Section 3, we present the simulation results, covering forest dynamics, model benchmarking, sensitivity analysis, and scenario-based optimization. In Section 4, we discuss the broader implications for sustainable ITP management and identify key limitations and directions for future research.

## 2. Materials and Methods

In this section, we detail the methodological framework used to develop and assess the LUNTIAN agent-based model. We begin by describing the case study area and the social-ecological context that informed the model design. Next, we justify the use of an agent-based modeling (ABM) approach, highlighting its capacity to capture and simulate complex interactions in managed forest systems. We then outline the model's structure following the ODD (Overview, Design Concepts, and Details) protocol, including agent definitions, scheduling logic, and key submodels. Finally, we explain the methods used for model validation, sensitivity analysis, and scenario exploration. The model code is openly available in a public repository, with the link provided in Supplementary Materials.

### 2.1. *Description of the Case Study*

#### 2.1.1. Study Site

The study site is located inside the University of the Philippines Laguna Land Grant (UP LLG)—at 14°23′ N, 121°29′ E. It is one of the two land grants in the southern portion of the Sierra Madre Mountain Range in the Luzon Island of the Philippines, awarded to UP through Republic Act 3990 dated 18 June 1964 and confirmed through Special Patent No. 532 and Original Certificate of Title P-919 issued 8 August 1965 [8]. Figure 1A shows the location of the industrial tree plantation (ITP) area in relation to Laguna de Bay, while Figure 1B shows its location within the entire UP LLG site.

<figure><img src="image_3.png" /><figcaption>**Figure 1.** Location of the case study site.</figcaption></figure>

The biophysical environment of the study site consists of hilly to mountainous terrain with a mix of forest cover types. A forest inventory conducted in 2019 provides descriptive data for the gray-shaded portion of the UP LLG shown in Figure 1B. This
inventory was based on ten established transects within that area and does not represent the entire land grant; no comprehensive forest inventory is currently available for the full extent of the UP LLG. Within the surveyed portion, 48 species were recorded, 45 native and 3 naturalized, with native trees constituting approximately 99% of the recorded population, underscoring the ecological integrity and semi-natural character of the inventoried area. Among the native species, members of the Dipterocarpaceae family—such as *Anisoptera thurifera* (locally known as Palosapis)—are dominant, accounting for over 34% of the native tree population. The remaining 1% comprises naturalized species, most notably *Swietenia macrophylla* (Mahogany), which is represented by three individuals in the inventory. This represents a marginal fraction of the tree population in the surveyed area, highlighting Mahogany’s limited ecological footprint. These naturalized individuals may be present due to earlier planting activities or unintentional spread.

### 2.1.2. Social Issue

The University of the Philippines Laguna Land Grant (UP LLG) has remained under the ownership of the University since its designation in the mid-20th century, though it has experienced changes in administrative management and oversight over time. Originally used for timber production under government license, it was formally transferred to the University of the Philippines in 1964 for academic and research purposes [9]. Today, the site is managed by the UP Los Baños Land Grant Management Office (UPLB LGMO) [10]. Among current initiatives, the ITP has emerged as a promising strategy for generating income through managed timber harvests while maintaining ecological integrity [8], and approximately 1411 hectares of the 3435.4-hectare site has been allocated for this project. Given its proximity to local communities, the ITP also presents an opportunity to promote employment and rural development, highlighting the importance of integrating social considerations into ITP design and management.

These social considerations are especially relevant given the evolving land-use dynamics surrounding the UP LLG. In 2020, Barangay San Antonio—the closest barangay (the smallest administrative division in the Philippines, similar to a village or district)—had a population of 11,263 individuals distributed across multiple sitios (smaller hamlets or sub-village units within a barangay) [11]. A field visit in January 2023 revealed that these sitios developed gradually over time, with the most recent, Sitio Balatkahoy, established about two decades ago near the proposed edge of the ITP area. Sitio Balatkahoy has since become a focal point of land-use contention. In 2020, the Department of Agrarian Reform (DAR) flagged the land grant as potentially eligible for redistribution under the Comprehensive Agrarian Reform Program (CARP), raising questions about the future of land designated for academic purposes [12]. A subsequent UPLB Legal Office site inspection in 2023 confirmed the presence of informal settlers, underscoring the complex challenges of undocumented occupancy and overlapping claims [13]. Once sparsely populated, the area has experienced increased in-migration driven by economic activity and shifting land use—from timber extraction to academic and plantation development.

### 2.1.3. Forest Management Choice

In this study, the forest management system is modeled as a Social–Ecological System (SES), following Ostrom’s definition [14], where the university, investors, and laborers are considered users (or actors) and institutional components are represented by the rules governing harvesting and investment. The roles of various entities—including the university as resource manager, tree growth dynamics, planting and harvesting activities, investors, and community stakeholders—are explicitly integrated to reflect the complex interactions that characterize the UP LLG.
Human resource dynamics are central to the sustainability of the ITP and are closely linked to investment decisions. The model follows the university's implementation plan: private logging is excluded, and community members are hired as formal laborers for nursery maintenance, planting, and harvesting. Employment is offered based on investment-driven demand; thus, the number of active investors directly affects labor opportunities. When labor demand is unmet or when working conditions are perceived as unfavorable, some individuals may turn to unauthorized independent harvesting. This dynamic reflects real-world tensions between job scarcity and informal forest use. Figure 2 illustrates the institutional relationships and resource utilization links among the key actors in the ITP system, with a focus on how the university mediates interactions among investors, community labor, and the forest landscape.

<figure><img src="image_4.png" /></figure>

**Figure 2.** Schematic representation of key actor interactions within the ITP system. Solid arrows indicate communication channels, while dashed arrows represent utilization or resource exchange links. The diagram highlights how the university serves as an intermediary between investors and local community members, allocating labor and coordinating plantation activities. Labor availability and forest use are shaped by these institutional relationships and influence both formal employment and potential unauthorized harvesting behavior. Icons by various contributors via Canva.com (Free Content License). Figure created by the author.

A key component of the plantation strategy involves nursery management, which supports regeneration by maintaining a steady supply of seedlings—particularly native species—for enrichment planting. Managed by the university, nurseries serve as propagation hubs where community laborers are hired to collect wild seedlings, nurture them in controlled plots, and prepare them for transplanting into investment areas. The nursery system not only sustains replanting cycles but also plays a vital role in reinforcing species selection strategies aimed at balancing ecological integrity and economic viability.

Species selection, therefore, is tightly linked to nursery operations and broader management practices. The UP LLG hosts a mix of native and exotic species, as documented in a 2019 forest inventory. In operational contexts, fast-growing exotics such as *S. macrophylla* are often favored for their quicker economic returns. However, native species—particularly from the Dipterocarpaceae family—are more ecologically beneficial, contributing to biodiversity and long-term resilience. This study incorporates both groups to reflect these trade-offs. The selection was shaped by three constraints: (1) limited availability of species-specific growth data, (2) the need for manageable and interpretable scenario comparisons, and (3) a focus on testing strategic trade-offs over replicating full silvicultural complexity. This approach reflects the university's implementation strategy and supports evidence-based exploration of sustainable plantation configurations.

## 2.2. *Justification for an Agent-Based Modeling Approach*

Social-Ecological Systems (SESs) exhibit complex and dynamic behavior, driven by the co-evolving interactions between human societies and natural environments [15,16]. In these systems, human decisions and ecological processes influence one another across
multiple spatial and temporal scales. Industrial tree plantations (ITPs) represent a specialized form of SES by combining biological growth, competitive tree dynamics, and intensive management to pursue both economic and ecological goals [17,18].

Researchers have long used simulation models to understand SESs. For instance [17], developed an individual-based growth-yield model capable of simulating decades of forest development in minutes, while [18] introduced spatially explicit models that captured forest succession and competitive dynamics. Building on these foundational efforts, agent-based models (ABMs) have emerged as powerful approaches to representing individual-level behavior in SESs, revealing system-level patterns—such as resilience and regime shifts—that arise from local interactions [19,20]. Despite the strengths of these models, researchers have seldom applied ABMs to forest-based systems—particularly within industrial plantations—even though these models are well suited for simulating the interplay of stakeholder behavior, environmental feedback, and institutional policies in ITP contexts [17–19]. For example, ref. [20] demonstrated how ABMs, when combined with land-use optimization, can assess the ecological and economic feasibility of sustainable resource use. Meanwhile, most existing forest models continue to emphasize biophysical processes while underrepresenting socio-economic factors such as investment decisions, labor dynamics, and enforcement mechanisms [17,18]. Recent applications, such as [21] on palm oil landscapes, show that ABMs can effectively capture complex decision making and policy interactions in managed systems. Furthermore, a review by [16] confirms that ABMs and system dynamic models remain underutilized in intensively managed settings such as ITPs. These gaps underscore the need for models that integrate detailed ecological processes with agent-level decision making to enable hypothesis testing, scenario exploration, and the design of more adaptive and sustainable management strategies.

Unlike aggregate models that rely on probabilistic assumptions, ABMs allow for the fine-grained representation of actor heterogeneity and dynamic role switching [22], such as laborers becoming informal harvesters in response to limited employment opportunities. These transitions are not isolated decisions—they are part of feedback loops, where institutional rules, economic incentives, and enforcement presence shape behavior over time. For instance, persistent unemployment can lead to increased unauthorized harvesting, which in turn triggers enforcement responses and alters labor dynamics in future cycles. Such interdependent processes are difficult to simulate realistically using non-agent methods without significant abstraction or loss of spatial and behavioral resolution. Capturing this localized, emergent feedback is key to understanding the full complexity of social-ecological systems.

In this study, we adopt an agent-based modeling approach to simulate ITP dynamics. ABMs offer a scientifically rigorous way to evaluate long-term management strategies under uncertain environmental and socio-economic conditions. Unlike system dynamic or analytical models, ABMs simulate autonomous actors—such as investors, laborers, and managers—within a spatially explicit environment. This bottom-up structure allows us to (1) represent behavioral diversity (e.g., passive vs. risk-seeking investors), (2) simulate feedback loops (e.g., labor competition driving informal harvesting), (3) test interventions (e.g., incentives or enforcement), and (4) observe emergent outcomes that traditional models may miss. By capturing these interdependencies, our ABM supports the systematic analysis of ecological, economic, and social outcomes. We use it to evaluate how harvesting intensity, species composition, and community engagement affect forest cover, biodiversity, and financial returns. The model also enables evidence-based decision making by providing interpretable metrics—such as projected net earnings, forest regeneration, and employment outcomes—across a range of management scenarios. Rather than relying on static assumptions or averages, we simulate and compare specific strategies to identify those that best align with sustainability goals.
## 2.3. *Description of the ABM*

In this section, we detail the ABM by using the Overview, Design Concepts, and Details (ODD) protocol [23,24], a standardized framework for describing agent-based models. The Overview Section outlines the model's purpose, entities, spatial and temporal scales, and simulation schedule. The Design Concepts Section presents the key behavioral assumptions and theoretical principles embedded in the model. Finally, the Details Section describes the implementation of the model's subcomponents and simulation scenarios.

### 2.3.1. Overview

LUNTIAN is a spatially explicit agent-based model developed to simulate an industrial tree plantation within the UP Laguna Land Grant (UP LLG). Built on the GAMA platform (version 1.9.3) [25], it integrates empirical data and stakeholder input to explore management strategies that balance ecological regeneration, financial viability, and community engagement. The model captures interactions among key actors (university, investors, laborers, and community members) and biophysical processes such as tree growth, hydrology, and land-use change.

#### Purpose

The LUNTIAN model—short for Labor, University, Timber Investment, and Agent-based Nexus—embodies an interdisciplinary approach to promoting sustainable forest recovery in the Philippine context. The name also draws from the Filipino word luntian, meaning “green,” to emphasize its commitment to ecological restoration. Grounded in the institutional goals of the UP LLG, the model provides a decision-support tool for exploring plantation management strategies, especially where field experimentation is limited. Rather than intervening in unresolved land tenure or political issues, LUNTIAN offers a neutral platform for simulating ecological, economic, and social trade-offs within a defined implementation framework.

#### Entities and State Variables

Two main agent types—biophysical and social—interact within the simulated ITP environment to represent the system's ecological and institutional dynamics. Each agent is characterized by internal state variables—such as employment status, investment history, and patrol assignment—that influence its behavior and interactions over time. These variables are informed by empirical data, the literature, and stakeholder input, and are operationalized through submodules described in later sections. Figure 3 presents the model's core entities and their interrelationships in a UML class diagram.

The biophysical entities, i.e., river, road, climate, and soil agents, define the environmental landscape and influence ecological processes such as tree growth, hydrology, and site productivity. The tree entities, central components of the model, are categorized into native and exotic types, where native trees generally follow growth dynamics characteristic of the Dipterocarpaceae family, while exotic trees are modeled after *S. macrophylla*. Trees are distributed across plots, each representing one hectare of land; plot attributes are shaped by interpolated soil and climate variables—including pH, rainfall, temperature, and evapotranspiration—as well as proximity to roads and rivers.

The social entities are the human and institutional actors involved in plantation management. At the core is the university agent, representing the UPLB LGMO, which coordinates operations such as nursery supervision, labor deployment, planting, harvesting, salary disbursement, and financial monitoring. Nursery plots are managed by the university to produce seedlings for enrichment planting. The market agent supports economic dynamics by supplying external pricing information for timber and labor. While the
university handles hiring directly, wages are benchmarked against market rates. Due to the absence of reliable time-series data in the Philippine context, both timber and labor prices are assumed to be constant throughout the simulation.

<figure><img src="image_5.png" /><figcaption>**Figure 3.** UML class diagram of the LUNTIAN model.</figcaption></figure>

Other key social agents perform operational, economic, and enforcement roles within the ITP system. Investor agents finance plantation activities and exhibit different levels of risk tolerance (i.e., risk averse, risk neutral, or risk seeking), and their behavior adapts to investment outcomes, shifting between passive, active, or disengaged states. Labor agents, drawn from the local community or university, execute ground-level tasks such as planting, nursery maintenance, and harvesting, but under conditions of limited employment or low wages, some of them may transition into informal roles as independent harvesters; to mitigate such unauthorized activities, special police agents patrol designated zones and enforce institutional regulations. The broader community member agents represent the surrounding population and can dynamically switch roles—becoming formal workers, informal harvesters, or passive observers—depending on contextual factors such as employment availability and perceived benefits.

## Process Overview and Scheduling

The LUNTIAN model proceeds in monthly time steps. One complete simulation experiment executes a full harvest rotation—approximately 180 cycles, equivalent to 15 years for a 15-year harvest interval. Figure 4 presents an overview of the processes executed during each simulation cycle.

At initialization, the model constructs the biophysical environment, representing the landscape as a grid of 1-hectare parcels, which serve as the fundamental spatial units in the model. Each parcel receives structural and environmental attributes—including soil properties and baseline climate conditions—based on reference data. Using these inputs, the model calculates the initial amount of available water per parcel. As the simulation progresses, it updates temperature, precipitation, and available water monthly to reflect seasonal patterns.

Following environmental updates, the model sequentially executes five core submodels that simulate the ecological, economic, and institutional dynamics of the ITP system. The stand growth submodel updates tree development through diameter growth, fruiting, recruitment, and mortality based on species-specific parameters and monthly climate data. Forest operations—including harvesting, enrichment planting, and nursery management—are handled by the forest operation submodel, with activities restricted to investment-active plots, while nursery work occurs on selected parcels near mother
trees. The investment process submodel identifies viable plots and simulates investor decisions according to risk preferences and projected returns. The employment dynamics submodel allocates labor from the university and nearby communities for plantation tasks and models informal harvesting behavior that may arise from unemployment or limited job access. Finally, the ITP policing submodel simulates enforcement through patrolling agents who monitor designated areas, detect unauthorized harvesting, and apply penalties using proximity-based rules. Together, these submodels enable the model to capture the system's multi-dimensional interactions at a monthly scale.

<figure><img src="image_6.png" /><figcaption>**Figure 4.** Scheduling of the monthly processes simulated by LUNTIAN.</figcaption></figure>

### 2.3.2. Design Concepts

Table 1 summarizes how the LUNTIAN model addresses the design concepts from [26], indicating whether each concept is included and describing how it is implemented or represented within the model.

<table><caption>**Table 1.** Summary of model design concepts as addressed in the LUNTIAN model.</caption><thead><tr><td>Design Concept</td><td>Included in LUNTIAN?</td><td>How It Is Addressed in LUNTIAN</td></tr></thead><tbody><tr><td>Basic Principles</td><td>Yes</td><td>The model is grounded in the social-ecological system (SES) framework [14], combining agent-based modeling (ABM) with sustainability science. At the system level, it integrates concepts from resilience theory [27] (e.g., adaptive feedback between labor markets and forest conditions), sustainability triad analysis [28] (balancing ecological, economic, and social pillars), and economic decision theory [29] (investor behavior via risk-based heuristics). At the submodel level, tree dynamics follow individual-based growth theory [30], while labor decisions are shaped by rational choice theory [31] and informal economy insights. Enforcement behavior is modeled on deterrence theory [31], with spatial patrolling reflecting real-world law enforcement strategies.</td></tr><tr><td>Emergence</td><td>Yes</td><td>Key outcomes—such as forest cover change, unauthorized harvesting, investor disengagement, and labor market shifts—emerge from interactions among agents and their environment. While individual behaviors (e.g., role switching or investment withdrawal) are programmed, the broader patterns are not predetermined. They result from how agents adapt to shifting employment prospects, profitability signals, and institutional constraints—capturing a form of weak emergence [32], where complex system-level dynamics arise from simple decision rules applied across diverse contexts.</td></tr><tr><td>Adaptation</td><td>Yes</td><td>Agents modify behavior based on experience (e.g., investors go passive and communities shift roles).</td></tr></tbody></table>
<table><caption>**Table 1.** *Cont.*</caption><thead><tr><td>Design Concept</td><td>Included in LUNTIAN?</td><td>How It Is Addressed in LUNTIAN</td></tr></thead><tbody><tr><td>Objectives</td><td>Yes</td><td>Each agent has explicit goals: profit, livelihood, enforcement, or sustainability.</td></tr><tr><td>Learning</td><td>Partial/No</td><td>There is no formal learning; behavioral shifts mimic adaptive responses without strategy evolution.</td></tr><tr><td>Prediction</td><td>Yes</td><td>Investors estimate future plot profitability by using risk-based expectations.</td></tr><tr><td>Sensing</td><td>Yes</td><td>Agents perceive relevant, localized information (e.g., profitability, job status, and illegal activity).</td></tr><tr><td>Interaction</td><td>Yes</td><td>Indirect interactions through employment, harvesting pressure, and enforcement feedback.</td></tr><tr><td>Stochasticity</td><td>Yes</td><td>Randomness in profit outcomes, detection chances, and agent decisions reflects uncertainty.</td></tr><tr><td>Collectives</td><td>No</td><td>Agents act individually; no group structures or emergent collectives are modeled.</td></tr><tr><td>Observation</td><td>Yes</td><td>Key indicators (e.g., forest status, employment, and harvest levels) are tracked for analysis.</td></tr></tbody></table>

### 2.3.3. Details

#### Initialization

To create the biophysical aspect of the model, the climate, soil, road, river, plot, and tree entities are initialized at the start of the simulation. The values for climate, soil, and road remain static throughout the simulation because they represent fixed baseline conditions: soil pH is not expected to change significantly over the 17-year rotation period; roads are treated as permanent infrastructure with fixed spatial locations; and monthly climate inputs (e.g., precipitation, temperature, and evapotranspiration) are based on historical climatological averages rather than dynamic weather events, reflecting long-term environmental patterns rather than short-term variability. Dynamic climate projections were not included in this version to isolate social-ecological interactions without introducing additional uncertainty; future versions may incorporate climate change scenarios as data availability improves. The climate entities were initialized using raster files representing temperature, precipitation, and evapotranspiration values characteristic of the study region. Temperature and precipitation data were taken from 2010–2018 records provided by WorldClim, based on CRU-TS 4.06 [33] and downscaled using WorldClim 2.1 [34], while for evapotranspiration, we used the global map of monthly reference evapotranspiration at approximately 1 km resolution developed by [35]. The soil entity attributes, including pH values, were derived from maps provided by the Bureau of Soils and Water Management [36], while hydrologic soil group classifications were based on data from [37]. Finally, the road and plot entities were constructed from vector files extracted from topographic and elevation raster data of the case study site provided by UPLB LGMO; these datasets are not publicly available.

Building on this initialization, the tree entities were generated based on field inventory data and additional modeling assumptions. The individual attributes of the tree agents were derived from a 2019 tree inventory conducted in the northern part of the UP LLG site. In the absence of a dedicated inventory for the ITP section and due to the lack of publicly available, nationally aggregated data on tree densities in logged secondary forests in the Philippines, we based the initial population estimate on the modeling parameter used by [38], whose work focused on simulating tree growth and carbon accumulation in the Philippine context using an assumed density of 400 trees per hectare. The original inventory contained 2565 trees, and to approximate the density across the entire ITP area (1411 hectares), we
replicated this dataset 50 times, yielding a total of 128,250 trees—roughly one-quarter of the estimated tree density per hectare—and randomly positioned them within the ITP area, as illustrated in Figure 5. This relatively low figure reflects the logged condition of the area and the presence of mixed land uses, such as agroforestry, roads, and zones with intermittent human habitation.

<figure><img src="image_7.png" /><figcaption>**Figure 5.** Location of the source data and the case study area within the LLG site. The section labeled "Actual Tree Inventory" indicates where the tree inventory was conducted, while the "ITP Area" denotes the region where the inventory data was applied.</figcaption></figure>

Overall, the initial model environment consists of a landscape comprising 1369 parcels—excluding those intersected by rivers—each representing a 1-hectare plot of land (see Figure 6). Across this landscape, 127,652 tree stands are randomly distributed, with some parcels also containing roads designated for timber transport.

<figure><img src="image_8.png" /><figcaption>**Figure 6.** State of the environment at the start of the simulation. Green dots indicate the locations of native trees, while yellow dots indicate the locations of exotic trees.</figcaption></figure>

The model also includes social agents: the market, labor agents, the university agent, investor agents, and special police agents. Each of them was assigned a random location outside the environment during initialization, and the number of labor, investor, and special police agents was determined by parameter values associated with each agent type.
Table 2 lists the parameter variables along with their descriptions and their respective recommended value ranges.

<table><caption>**Table 2.** Model parameters to be defined per simulation execution.</caption><thead><tr><th>Parameter</th><th>Description</th><th>Recommended Values</th><th>Units</th></tr></thead><tbody><tr><td>Nursery count</td><td>The number of nurseries that will be established.</td><td>{1, 3}</td><td>Number of nurseries</td></tr><tr><td>Community members</td><td>The total number of community members available for hire.</td><td>{15, 30, 50}</td><td>Number of individuals</td></tr><tr><td>Special police</td><td>The number of special police available in the ITP.</td><td>{1, 3, 10}</td><td>Number of individuals</td></tr><tr><td>Investor count</td><td>The total number of potential investors.</td><td>{5, 15, 50}</td><td>Number of individuals</td></tr><tr><td>Investor share on earning</td><td>The percentage of overall ITP earnings allocated to the investor once an investment deal has been made.</td><td>{0.3, 0.4}</td><td>Proportion (unitless %)</td></tr><tr><td>Investment rotation years</td><td>The total number of years until the plot invested upon is harvested.</td><td>{10, 15, 20}</td><td>Years</td></tr><tr><td>Timber price<sup>1</sup> multiplier</td><td>The multiplier applied to the base harvesting timber price to determine the minimum price required for the ITP to be profitable.</td><td>{1.0, 2.0}</td><td>Unitless multiplier</td></tr></tbody></table>

## Input Data

The model does not use input data to represent time-varying processes.

## Submodels

To operationalize interactions within the ITP system, LUNTIAN integrates five inter-dependent submodels. Their design was informed by a combination of field observations, stakeholder input, and the relevant literature to ensure scientific validity and contextual relevance. Together, these simulate the ecological, operational, economic, labor, and governance processes introduced in the Process Overview and Scheduling Section. The following provides a detailed description of each submodel.

**Stand Growth Submodel:** It manages tree development by simulating the physical growth of tree agents, their fruit production, and mortality at each simulation step. This submodel is based on the framework for single-tree models proposed by [41]. The growth process involves three core components: (1) diameter increment using species-specific growth equations; (2) mortality prediction based on probabilistic assessments; and (3) recruitment, which captures processes of maturation, fruiting, and seedling survival. Tree life stages are defined by diameter at breast height (DBH), progressing from seedling to sapling, then pole, and finally adult. The equations used in this submodel were based on the established literature (see Appendix A) and thus had already been calibrated.

**Forest Operation Submodel:** It introduces human-driven activities such as harvesting, planting, and nursery management. The harvesting and planting routines are adapted from the GAMA Forest Model described in [42], incorporating techniques such as selective logging [43] and enrichment planting [44], which are executed only in plots under active investment agreements. Harvesting takes place twice—at the beginning and end of the investment cycle—while enrichment planting is triggered when either basal area or tree density falls below a defined threshold. During harvests, the university agent consults the market agent to determine roundwood prices (in board feet), after which profits are

<footnote n="1">¹ The default price-per-board-foot for harvested roundwood was based on the average national retail price of lumber in 2016 [39]. With a 7-year inflation rate of 1.0342 [40], the default price-per-board-foot for native trees is set to PHP 49.35, while that for exotic trees is PHP 45.06.</footnote>
computed and allocated between the university and investors based on a tunable profit-sharing ratio.

Nursery management is tightly integrated with planting and harvesting to support the sustainability of enrichment efforts. This process, led by nursery laborers, begins with collecting seedlings from the ITP area and cultivating them in designated nursery plots; once mature, saplings—prioritizing native species—are transferred to investment plots, and if the supply is insufficient, the university supplements it by purchasing additional seedlings to meet planting targets. This system ensures a steady pipeline of native tree species, reinforcing the plantation's long-term regenerative capacity.

Although specific silvicultural treatments such as pruning, thinning, and fertilization are not modeled as discrete actions, they are implicitly accounted for through labor and maintenance costs. These aggregated costs reflect typical plantation upkeep activities that enhance tree growth and productivity. However, the model considers only the revenues from final harvests when calculating economic returns, while intermediate income streams—such as those from pruned biomass—are excluded to maintain the focus on maximizing net returns from full harvest cycles.

**Investment Process Submodel:** It facilitates the evaluation and allocation of investment opportunities and is activated once a set number of nurseries are established to ensure that the ITP can adequately support incoming investment requests. Investable plots are then offered to potential investors categorized by risk preference: risk averse (investing with a 0.25 probability if projected gains are <50%), risk neutral (0.50 probability), and risk seeking (0.75 probability if gains are ≥ 50%), following [45]. These probabilities determine the likelihood of transitioning from an “interested-passive” state—where investors are aware but uncommitted—to an active investment state. After completing an investment, each investor evaluates the outcome; if the result is deemed unsatisfactory, the investor may transition into a “not-interested passive” state, where the probability of the investor remaining disengaged ($P(NI)$) or returning to the “interested-passive” state ($P(\text{Interested})$) is governed by their risk profile. Figure 7 illustrates these state transitions.

<figure><img src="image_9.png" /><figcaption>**Figure 7.** State-transition diagram for the investor agent.</figcaption></figure>

The design of the investment process submodel was guided by stakeholder input from a meeting with the Land Grant Management Office (LGMO), which played a key role in shaping the intended role of investors in the planned ITP implementation. Based on this meeting, the ITP is conceptualized as a program to engage potential investors—particularly university alumni—by inviting them to finance one-hectare parcels, guiding the representation of investor agents as key initiators of plantation activity.

**Employment Dynamics Submodel:** It governs the allocation of labor resources and captures transitions between formal employment and informal harvesting. This submodel operates in three phases, initial harvesting and planting, mid-cycle maintenance, and final harvesting and replanting, and community members and university workers form the labor pool.

To address labor shortfalls, the university issues open hiring calls whenever demand exceeds supply. Workers evaluate job opportunities by sensing the presence of other laborers and comparing expected wages with potential income from independent harvesting.
If a community member remains unemployed for several months, they may eventually shift to informal harvesting—particularly in fringe parcels at the plantation's border. This mechanism reflects field reports of labor competition observed at the study site and is contextualized by broader historical precedent. Ref. [1] documents how as early as the 1970s, an estimated 200,000 upland families across the country were already engaged in shifting cultivation; forestry experts at the time expressed concern that providing more efficient harvesting tools might unintentionally expand unauthorized timber extraction. These insights support the submodel's assumption that labor scarcity can trigger transitions from formal employment to informal, unsanctioned forest use.

The employment dynamics submodel, therefore, captures the evolving tension between cooperative labor systems and opportunistic resource use. The full progression of laborer engagement is illustrated in Figure 8.

<figure><img src="image_10.png" /><figcaption>**Figure 8.** State-transition diagram for the laborer agent.</figcaption></figure>

**ITP Policing Submodel:** It introduces a spatial enforcement mechanism to maintain order and discourage unauthorized activities—particularly illegal harvesting—within the plantation. Policing agents represent special police officers who patrol the ITP area, with only laborers and these agents being permitted to enter the plantation in the model; while real-world access could include external actors, such complexity is beyond the current model's scope. Patrols are repositioned every three months to simulate unpredictability, and enforcement occurs when a community member enters a police agent's observable range, triggering a "catch" event. This results in the laborer forfeiting any financial gain, becoming inactive, and accumulating a violation count. These penalties are designed to reinforce rule compliance and deter repeated infractions.

The design of this submodel was informed by a site visit in January 2023, during which the researchers were escorted by special police assigned to monitor the site. Informal interviews conducted during the visit revealed that their duties included regular patrolling, discouraging unauthorized activities, and intercepting informal harvesters; the officers also shared that patrol teams operate in multiple groups and rotate locations periodically—patterns directly reflected in the modeled agent behaviors, including spatial movement and group deployment.

## 2.4. *Model Assessment and Exploration*

Due to limited assessment data, model validation was conducted solely on the bio-physical layer. Simulation outputs from LUNTIAN were compared with those generated by 3PG2PY [46], an adaptation of the widely used 3PG forest growth model [5], standard-izing key input parameters across the models to enable direct comparison. We selected 3PG2PY because its structural and functional configuration closely aligns with that of LUNTIAN: it relies on similar input variables (stand age, tree biomass, stand volume, soil water, and canopy cover), operates on a monthly time step, and requires primarily weather data; additionally, it has been successfully applied to long-term simulations spanning up to 120 years, producing ecologically realistic growth trajectories [5]. The assessment focused on two core output variables: mean diameter at breast height (mean DBH, in cm)
and stand basal area (SBA, in m²). To account for differences in forest composition, the evaluation was performed separately for two scenarios—one involving only native species and another involving only exotic species—reflecting the original 3PG model’s design for monospecific plantations.

The assessment procedure comprised several key steps: First, the LUNTIAN experiment was run, during which the initial and final values of the observed variables—mean tree DBH and basal area per plot—were recorded. Next, 3PG2PY was executed under the same initial conditions, and the final values of the observed variables were also recorded. Finally, mean DBH and basal area per plot from both models were compared.

Two evaluation methods were employed to determine whether LUNTIAN accurately captures tree growth dynamics and reliably simulates forest stand development. The first method, Mean Absolute Percentage Error (MAPE) [47], quantifies prediction accuracy by calculating the average percentage deviation of forecast values from observed data, using the formula

$$ \frac{1}{N} \sum_{t=1}^{N} \left| \frac{A_t - F_t}{A_t} \right| $$

where $A_t$ represents the actual values generated by the 3PG2PY model (serving as the benchmark or reference model) and $F_t$ represents the forecast values produced by the LUNTIAN model at data point t, with N being the total number of data points. Forecasting accuracy was interpreted using established thresholds: MAPE < 10% indicates highly accurate predictions, 11–20% reflects good accuracy, 21–50% is considered reasonable, and values exceeding 51% indicate poor accuracy [48]. The second method, Fuzzy Kappa [49], compares observed and simulated maps using fuzzy logic to accommodate minor spatial discrepancies. The Fuzzy Kappa coefficient was calculated using the implementation available in the GAMA platform, which is based on the approach described in [50]. According to [51], the agreement levels on the interpretation scale for Kappa statistics are poor (<0.00), slight (0.00–0.20), fair (0.21–0.40), moderate (0.41–0.60), substantial (0.61–0.80), and almost perfect (0.81–1.00). To account for differences in growth drivers between LUNTIAN and 3PG—such as the latter’s reliance on radiant energy—model outputs were grouped into bins prior to comparison, with the optimal number of bins selected based on the highest Fuzzy Kappa score.

For model exploration, we employed Sobol sensitivity analysis integrated into the GAMA platform [25]. Additionally, we used OpenMole [52,53], a comprehensive tool that offers a range of experimental methods to enhance our understanding of model behavior and capabilities. In particular, we utilized the Non-dominated Sorting Genetic Algorithm (NSGA2) for parameter optimization to answer the primary research question.

## 3. Results and Discussion

This section presents the key findings from the LUNTIAN model, organized into four main components: baseline forest dynamics, model benchmarking, sensitivity analysis, and simulation-based scenario exploration. We begin by examining how the simulated forest evolves under purely ecological processes to establish a reference trajectory. Next, we conduct a benchmarking exercise that evaluates the model’s biophysical plausibility by comparing tree growth outputs to those of the widely used 3PG2PY model. We then apply a global sensitivity analysis to identify which model parameters exert the greatest influence on economic outcomes, particularly net earnings. Finally, we perform a set of simulation experiments, including critical parameter testing and multi-objective optimization, to identify configurations that promote ecological regeneration, economic viability, and social inclusion.
Together, these results support the model's internal plausibility and its utility as a tool for simulation-based experimentation in sustainable industrial tree plantation strategies. While LUNTIAN has not yet been validated against real-world socio-economic outcomes due to limited field data, the model was designed primarily as an exploratory platform for understanding system behavior and trade-offs. Moreover, the model's focus on institutional design and agent-based social dynamics is the reason climate change projections were not yet incorporated, as the primary objective was to explore the operational feasibility and stakeholder interactions against a fixed-climate baseline.

### 3.1. *Baseline Forest Dynamics*

To illustrate the model's behavior, in this section, we examine the simulated evolution of the forest landscape under baseline conditions. In this environment-only configuration, the LUNTIAN model simulates forest development driven solely by internal ecological processes—namely, growth, fruiting, and mortality—without any external interventions such as harvesting or replanting. The simulation began with 127,652 trees randomly distributed throughout the ITP area. As shown in Figure 9, the total tree population, dominated by saplings and seedlings, peaked at approximately 6 million in year 70.

<figure><img src="image_11.png" /><figcaption>**Figure 9.** Tree growth trend over 2000 simulation cycles (approximately 165 years) in the absence of external influences. (A) Growth trend of trees categorized by species type—native versus exotic. (B) Growth trend of trees categorized by life stage—seedling, sapling, pole, and adult.</figcaption></figure>

By year 145, the environment had stabilized with around 4 million trees, with a gradual shift toward older cohorts: the number of adult and pole-stage trees continued to rise, while younger populations held steady. These dynamics indicate that the forest had approached an ecological equilibrium by year 145, providing a structural reference point for evaluating model plausibility and informing subsequent comparisons with models such as 3PG.

### 3.2. *Model Benchmarking*

The comparison with the 3PG2PY results revealed that LUNTIAN performed comparably in terms of mean DBH and SBA for both the native-only and exotic-only setups, for final stand ages of up to 80 years. As shown in Table 3, the error for both output variables is below 20%, meaning good forecasting accuracy, establishing that LUNTIAN results are consistent with those of the 3PG2PY model. The Fuzzy Kappa analysis revealed varying levels of spatial agreement across forest types and variables. For native-only species, mean DBH yielded a value of 0.6047, falling just below the 0.61 threshold for substantial agreement and thus classified as moderate, while SBA was 0.7543, indicating substantial agreement. In exotic-only stands, both mean DBH (0.6612) and SBA (0.7680) showed substantial agreement, suggesting that the model performs slightly better in predicting the spatial structure of exotic species. The near-substantial mean DBH score for native species still reflects reasonably strong performance, though the broader use of growth equations based on the entire Dipterocarpaceae family—compared with the more uniform, mahogany-based
approach for exotics—may have introduced additional variability and reduced predictive precision for native stands.

This performance gap reflects the broader variability in growth patterns among native species, which were modeled using generalized equations based on the Dipterocarpaceae family. These equations capture average trends but may not fully account for the species-specific traits, ecological niches, or competitive strategies of individual native taxa. In contrast, the exotic species—modeled solely on *Swietenia macrophylla* (mahogany)—benefited from more detailed, species-specific parameters derived from targeted empirical studies. This difference in parameter granularity led to more consistent and predictable outputs for exotics, while native species exhibited higher variance and reduced predictive precision due to their ecological diversity and the reliance on aggregated growth assumptions.

These findings indicate that even when accounting for spatial distribution, the performance of LUNTIAN is comparable to that of the 3PG2PY model.

<table><caption>**Table 3.** Assessment results.</caption><thead><tr><th></th><th></th><th>MAPE</th><th>Fuzzy Kappa</th></tr></thead><tbody><tr><td rowspan="2">Native Only</td><td>Mean DBH</td><td>8.4914%</td><td>0.6047</td></tr><tr><td>SBA</td><td>10.5313%</td><td>0.7543</td></tr><tr><td rowspan="2">Exotic Only</td><td>Mean DBH</td><td>8.4148%</td><td>0.6612</td></tr><tr><td>SBA</td><td>14.4600%</td><td>0.7680</td></tr></tbody></table>

To assess spatial agreement, Figure 10 presents a set of spatial maps derived from both the LUNTIAN simulation and the 3PG2PY outputs. Specifically, the figure includes three types of maps for each variable: the simulated output from LUNTIAN, the observed reference from 3PG2PY, and a comparison map. Since Fuzzy Kappa is specifically designed for evaluating similarity between categorical or class-based raster data [54], we discretized the continuous model outputs (specifically, mean DBH and stand basal area) into bins to ensure compatibility with the metric; this is a well-established preprocessing step in environmental modeling and is frequently employed to improve model interpretability and computational efficiency in approaches such as Bayesian Networks [55]. To avoid arbitrary thresholds, we tested multiple binning schemes and selected the configurations that maximized Fuzzy Kappa scores for each variable; the resulting bin counts—11 for native-only mean DBH, 15 for native-only SBA, 19 for exotic-only mean DBH, and 16 for exotic-only SBA—represent these optimized settings. This approach introduces a tolerance range that groups similar values, mitigating sensitivity to minor numerical variation and enabling more meaningful spatial comparisons.

The fuzzy map comparison highlights areas of spatial agreement and disagreement between the simulated and observed outputs. Specifically, white areas correspond to matched values between the LUNTIAN and 3PG2PY maps, gray areas indicate mismatches beyond the defined bin tolerance, and black areas fall outside the plantation boundary and are excluded from the comparison. While these maps suggest broad discrepancies and localized clustering—particularly near river networks—such patterns remain speculative. Further quantitative analysis and statistical validation are required to confirm the significance of these spatial structures and understand their ecological or model-driven origins.
<figure><img src="image_12.png" /><figcaption>**Figure 10.** Fuzzy Kappa maps for LUNTIAN and 3PG2PY outputs. Each row shows the simulated output (top), the observed 3PG2PY output (middle), and the fuzzy comparison map (bottom). Only the fuzzy map comparison uses a consistent color scheme, where white indicates agreement between the LUNTIAN and 3PG2PY model outputs, gray represents disagreement beyond bin tolerance, and black denotes areas outside the modeled environment. The simulated and observed maps use automatically generated color gradients based on the distribution of each run and therefore do not include a fixed legend to avoid misinterpretation.</figcaption></figure>

### 3.3. *Sensitivity Analysis*

For the sensitivity analysis, we focused on overall final net earnings as the primary performance metric. While the broader evaluation of LUNTIAN considers ecological, social, and economic dimensions, we centered this section on evaluating the model's robustness in representing long-term economic viability, which constitutes this study's principal research objective. Final net earnings provide a comprehensive economic indicator, reflecting the cumulative impacts of management decisions over time, and offer a practical, interpretable basis for evaluating sensitivity to key input parameters. While LUNTIAN tracks multiple outcomes, including ecological regeneration and social labor transitions, we focused on final net earnings in the sensitivity analysis, since this metric captures the long-term economic viability of the ITP. Social and environmental metrics were integrated in the optimization phase (Section 3.4), where trade-offs among the three sustainability pillars were explored in greater depth.

We applied Saltelli's variance-based sensitivity analysis to quantify the influence of model parameters on output variability [56]. Specifically, we computed both the first-order and total-order sensitivity indices for the following parameters: community members, nurseries, investors, special police officers, the share of investors in ITP earnings, and the length of one complete harvest rotation (see Figure 11).

In Sobol sensitivity analysis, the first-order sensitivity index captures the direct effect of a single input parameter on the output variance, assuming that all other parameters remain fixed, which reflects how much a parameter alone influences the model outcome. In contrast, the total-order sensitivity index accounts for the parameter's overall contribution, including both its direct impact and all interactive effects with other parameters [57]. This distinction is especially relevant in a complex agent-based system such as LUNTIAN, where socio-economic and ecological processes are tightly coupled. Parameters often do not act independently; their effects emerge through feedback loops and interdependencies. Thus, a high total-order index relative to the first-order index suggests that a parameter plays a crucial role in shaping system behavior through interactions with others—an important consideration in understanding the emergent dynamics of social-ecological systems.
<figure><img src="image_13.png" /><figcaption>**Figure 11.** Sensitivity analysis of variables affecting final net earnings.</figcaption></figure>

For each parameter set, we executed a full simulation of the ITP—from initialization to the completion of the harvest rotation for every active investor—and Sobol analysis was performed using 420 Saltelli samples.

Based on the Sobol sensitivity analysis results, the investor count clearly exerted the strongest influence on the final net earnings, as shown by its high total-order index and the notable gap between its first-order and total-order values. This indicates that the investor count not only had a substantial direct effect but also interacted meaningfully with other parameters. Higher investor numbers generally lead to more financed plots and increased labor demand, which boosts earnings; however, this relationship is not linear—saturation effects, coordination challenges, and management overhead may reduce marginal gains as investor numbers grow.

Harvest rotation also exhibited a significant influence, with its total-order index exceeding the first-order value, suggesting strong interactive effects. Shorter rotations may increase cash flow and accelerate returns, but they typically reduce final biomass and potential yield; conversely, longer rotations can enhance ecological regeneration and total harvest volume, but they accumulate higher costs and delay revenue realization. This underscores the importance of identifying an optimal rotation period that balances short-term returns with long-term sustainability.

The community count and investor share also showed moderate-to-high influence on earnings. The former's slightly higher indices reflect its impact on labor availability and associated costs, while investor share directly shapes the distribution of profits. In contrast, the police count and nursery count had relatively lower but still measurable effects, indicating that while important for governance and regeneration, their roles are less dominant in driving overall profitability.

Overall, these results highlight how LUNTIAN's performance is shaped by the interplay of economic, ecological, and social levers and demonstrate the need for careful calibration to achieve balanced, sustainable outcomes.

We also observed that the first-order indices for the variables nursery count and police count were slightly negative. In theory, Sobol indices are expected to lie between 0 and 1 because they represent the proportion of output variance attributable to each input parameter. However, negative values can occur, particularly when the true sensitivity of a parameter is near zero. This phenomenon is primarily attributed to numerical inaccuracies in the
estimation process. For example, the Saltelli method may yield slightly negative values for factors with negligible influence due to errors in covariance estimation—a problem that is exacerbated in stochastic simulations such as those generated with ABMs [56]. Moreover, when these methods are applied to ABMs, independent sampling and the high dimensionality of the parameter space further degrade the accuracy of the estimates, resulting in negative sensitivity values. Although some researchers set these negative values to zero for interpretability, doing so may introduce bias into the analysis [58,59]. Given that the first-order indices for these parameters remain within -0.05, we can conclude that neither variable significantly contributes to the variance in the evaluation metric (that is, the final net earnings).

## 3.4. Simulation-Based Scenario Exploration

The primary research question driving this modeling effort is the following: “Can the industrial tree plantation (ITP) be implemented sustainably, and if so, what are the optimal parameters to maximize returns?” To address this, we performed an optimization analysis aimed at identifying configurations that satisfy three key sustainability criteria—environmental, economic, and social—based on the sustainability triangle framework introduced by [60]. Environmental sustainability was evaluated by measuring the final tree population at the end of each simulation, economic sustainability by calculating investors’ and net profit, and social sustainability by examining the number of community members transitioning to independent harvesting. Therefore, in our analysis, the optimal configuration is defined as one that maximizes the tree population (ecologically sustainable), generates a positive net profit for both the university and investors (economically sustainable), and minimizes the number of community members engaging in independent harvesting (socially sustainable).

We investigated the following simulation scenarios to determine the model's behavior under varying conditions:

*   **Critical parameters:** Intended to observe the impact of key variables identified during the sensitivity analysis.

*   **Main scenario:** Focused on exploring the interplay among environmental, economic, and social factors to achieve sustainable outcomes.

### 3.4.1. Critical Parameters

#### Experimental Design

To assess which parameters most significantly influence the profitability of the ITP, we conducted a targeted analysis of four variables identified in the sensitivity analysis: (1) number of community members, (2) investor count, (3) investor profit share, and (4) harvest rotation years. Each parameter was varied individually while holding the others constant at their maximum recommended values. This method allowed us to isolate and evaluate the effect of each parameter on the model's primary economic output—net earnings.

We performed five simulation experiments:

*   **All Max:** All critical parameters set to their maximum values.

*   **Rotation Min:** Harvest rotation years set to the minimum.

*   **Investor Share Min:** Investor profit share set to the minimum.

*   **Investor Count Min:** Investor count set to the minimum.

*   **Community Members Min:** Number of community laborers set to the minimum.

#### Results and Analysis

Figure 12 summarizes the net earnings across the five simulation setups. In four scenarios—All Max, Rotation Min, Investor Count Min, and Community Members Min—the model produced negative net earnings, underscoring how small variations in key
parameters can push the ITP into unprofitability. Only the Investor Share Min scenario resulted in positive earnings, indicating that reducing the investor's profit share can effectively absorb high management costs.

These outcomes are consistent with the sensitivity analysis (Figure 11). The investor share parameter showed a strong direct influence on net earnings, as indicated by similar first- and total-order Sobol indices; in contrast, other parameters had greater interaction effects, contributing to output variability mainly through their interplay. This non-additive behavior highlights the importance of managing parameters in combination rather than isolation.

<figure><img src="image_14.png" /><figcaption>**Figure 12.** Comparison of management cost, ITP earnings, and net earnings across five parameter configurations: All Max, Rotation Min, Investor Share Min, Investor Count Min, and Community Members Min.</figcaption></figure>

### 3.4.2. Main Scenario

#### Experimental Design

We conducted an optimization experiment by using the NSGA2 algorithm implemented in OpenMole [53]. The goal was to identify a parameter configuration that satisfies environmental, economic, and social sustainability criteria. We used all the parameters specified in Table 2 at their recommended values. The outputs optimized included the following:

* Environmental sustainability: Maximize final tree population.

* Economic sustainability: Maximize net earnings and investor profit.

* Social sustainability: Minimize number of community members engaging in independent (unauthorized) harvesting.

#### Results and Analysis

During the NSGA2 optimization (run over 100 evolutions), all simulation outputs were logged across the explored parameter space. These logged results were then used to generate the pairwise correlation matrix shown in Figure 13. The matrix captures relationships between key output variables based on the full set of simulations produced during the optimization process. For instance, a strong positive correlation was observed
between tree population and management costs. This indicates that simulations with higher tree populations also tended to report higher management costs—primarily due to the increased labor and nursery demands associated with maintaining and enriching denser tree populations.

Similarly, investor profits were positively correlated with both management costs and tree population, suggesting that scenarios promoting forest growth also support economic returns, albeit with higher operational expenses. A strong negative correlation was observed between tree population and net earnings, indicating that scenarios promoting forest growth tended to incur higher labor and nursery maintenance costs, which in turn reduced net profitability—especially when not offset by proportional revenue increases.

<figure><img src="image_15.png" /><figcaption>**Figure 13.** Pairwise correlation matrix of observed variables. Asterisks indicate the strength of correlation: \* for values ≥ 0.2, \*\* for values ≥ 0.3, and \*\*\* for values ≥ 0.4.</figcaption></figure>

We also noticed a strong correlation between partner profits (community-hired laborers) and independent profits (unauthorized harvesters). This likely indicates that in scenarios with more labor opportunities and higher wages, independent harvesting may still occur—likely because not all community members are hired, as the number of available workers exceeds labor demand.

## Final Configuration and Outputs

The best-performing parameter combination obtained with the optimization process is nursery count = 3, price multiplier = 2, investor earning share (%) = 0.3, investment rotation years = 17, investor count = 50, police count = 3, and member count = 28. The outputs obtained under this configuration given a single optimized simulation are shown below (Table 4).

The optimal configuration identified by the LUNTIAN model demonstrates that sustainable ITP management is achievable through the strategic calibration of key parameters. A 55% increase in tree population over a 17-year harvest cycle—despite ongoing extraction—shows that forest regeneration can occur alongside commercial harvesting when supported by appropriate management practices and investment timing. This finding challenges the assumption that profitability and ecological resilience are inherently at odds, offering a data-driven case for integrating long-term ecological goals into plantation planning.
<table><caption>**Table 4.** Output metrics under the optimal configuration.</caption><thead><tr><td>Output Metric</td><td>Value</td></tr></thead><tbody><tr><td>trees</td><td>285,376</td></tr><tr><td>management_cost</td><td>PHP 27,923,416.49</td></tr><tr><td>investors_profit</td><td>PHP 10,699,486.66</td></tr><tr><td>net_earning</td><td>PHP 8,884,374.93</td></tr><tr><td>independent_profit</td><td>0</td></tr><tr><td>partners_profit</td><td>PHP 7,344,327.91</td></tr></tbody></table>

Economically, the configuration produced strong net earnings for both the university and investors, even with a relatively modest investor share of 30%. The 30% investor share identified in the optimal scenario demonstrates that a modest yet meaningful return is sufficient to sustain investor interest while enabling the university to cover management costs and provide stable employment for community members. In the final optimized configuration, this arrangement—combined with a 17-year harvest rotation, 3 nurseries, 50 investors, 28 community laborers, and moderate enforcement (3 special police officers)—resulted in a 55% increase in tree population, PHP 8.8 million in net earnings, and 0 instances of unauthorized harvesting. These outcomes indicate that when key parameters are carefully calibrated, it is possible to concurrently achieve ecological regeneration, financial sustainability, and social harmony.

From a practical standpoint, this configuration reflects the real-world feasibility of designing industrial tree plantation (ITP) systems that attract investment without over-reliance on high investor returns. A 30% share proves adequate not only to yield PHP 10.7 million in investor profit but also to preserve forest cover and fund essential community labor, where community members would earn over PHP 7.3 million in wages. The absence of informal harvesting further signals a healthy labor management dynamic, where equitable job access and fair compensation reduce incentives for unauthorized forest use.

In the context of the social-ecological system (SES) framework, these results affirm that the model supports positive outcomes across all three core metrics: ecological (increased forest cover), economic (viable net income and investor return), and social (harmonious labor relations and inclusive benefit sharing). LUNTIAN thus offers a practical decision-support tool for stakeholders aiming to balance long-term ecological goals with inclusive economic development.

## 4. Conclusions

In this study, we developed an agent-based model to simulate an industrial tree plantation (ITP) operation within a mountain forest managed by the University of the Philippines Los Baños Land Grant Management Office. The model integrates key biophysical processes with socio-economic dynamics, including investment decision making, labor engagement, and forest policing. The simulation results indicate that the ITP can be managed to ensure environmental, economic, and social sustainability. Notably, despite ongoing harvesting activities over a 17-year rotation, the LUNTIAN model projected a 55% increase in final tree population relative to initial levels., which suggests that carefully calibrated management practices can support forest regeneration even under commercial harvesting regimes.

The model's biophysical outputs were benchmarked against the established 3PG forest growth model, with comparisons based on mean diameter at breast height and stand basal area showing acceptable error margins. Fuzzy Kappa coefficients also demonstrated substantial spatial agreement, and a sensitivity analysis identified key parameters—particularly investor earnings share and community member count—as strong drivers of management costs and net
earnings. These results highlight how financial sustainability depends not only on ecological performance but also on the careful management of investment and labor dynamics.

Nonetheless, several limitations must be acknowledged. The behavioral rules gov-
erning social agents, including investors and laborers, are based on simplified heuristics
rather than empirical calibration or adaptive learning. While informed by the literature
and stakeholder input, these assumptions may not fully capture the nuances of real-world
decision making, especially under uncertainty. Additionally, the model treats market
variables—such as timber and labor prices—as static, thereby excluding the potential ef-
fects of economic volatility on system dynamics. Although the ecological processes were
validated, the socio-economic parameters lack similar empirical grounding. Furthermore,
the model is geographically specific to the UP LLG and may require significant adaptation
to be transferable to other ITP contexts with differing environmental, institutional, or
socio-political conditions.

Environmental dynamics are also simplified, with climate inputs being based on historical averages rather than incorporating projections of climate change or extreme events. Institutional complexity is partially modeled through special police enforcement, but broader policy and governance factors—such as land tenure reform, incentive schemes, or regulatory changes—are not yet integrated. These limitations point to several promising directions for future development, particularly in enhancing the model’s realism and decision-support capabilities. Key priorities include integrating climate-change projections, simulating dynamic timber and labor markets, and incorporating adaptive learning mechanisms to better reflect how agents respond to uncertainty over time. Future work will specifically focus on using time-varying climate data to evaluate long-term system resilience, implementing learning algorithms that allow key social agents—such as laborers and investors—to adapt their behavior based on experience, and refining decision-making rules through improved empirical calibration. Together, these enhancements aim to capture more realistic socio-environmental dynamics and reinforce LUNTIAN’s utility as a tool for exploring sustainable management strategies.

As with any simulation, simplification introduces the risk of omitting critical dynam-
ics. However, models are inherently simplifications built with a specific purpose—not to
capture every detail but to illuminate key interactions and trade-offs. LUNTIAN was inten-
tionally designed as a scenario-exploration platform, aiming to capture essential system
behaviors while maintaining interpretability, and its structure prioritizes the representation
of core social-ecological processes relevant to industrial tree plantations. Future work may
expand this foundation to include currently excluded variables, such as market volatility
and broader institutional complexity.

While this study does not propose direct prescriptions for ITP implementation, the simulation results offer meaningful insights into the conditions under which plantations could achieve ecological, economic, and social sustainability. The optimal configuration identified by the model reflects a balance of parameters—such as rotation length, investor share, nursery count, and labor availability—that produced positive financial returns, increased forest cover, and eliminated unauthorized harvesting. Rather than presenting a fixed blueprint, the model provides a flexible, evidence-based platform for exploring plausible management strategies and trade-offs. As such, LUNTIAN is best viewed as a tool for scenario exploration, communication, and in vitro experimentation—supporting stakeholders in identifying potential pathways, testing assumptions, and engaging in more informed discussions around sustainable plantation development.

**Supplementary Materials:** The model code can be found at https://github.com/zarnejo-exers/LUNTIAN (accessed on 24 June 2025).
**Author Contributions:** Conceptualization, N.B. and Z.A.; methodology, Z.A., B.G. and M.S.; software, Z.A.; validation, Z.A., B.G. and M.S.; formal analysis, Z.A., B.G. and M.S.; investigation, Z.A.; resources, N.B., B.G., M.S. and Z.A.; data curation, Z.A.; writing—original draft preparation, Z.A.; writing—review and editing, Z.A., B.G., M.S. and N.B.; visualization, Z.A.; supervision, B.G., M.S. and N.B. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research study received no direct external funding

**Data Availability Statement:** The data used to generate the simulated forest are based on an unpublished dataset. Only transformed and aggregated data can be shared.

**Acknowledgments:** We would like to thank the Department of Science and Technology—Science Education Institute (DOST-SEI), Philippines, and the French Embassy to the Philippines and Micronesia for their financial support through the PhilFrance-DOST Fellowship of Z. Arnejo.

**Conflicts of Interest:** The authors declare no conflicts of interest. The funders had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript; or in the decision to publish the results.

## Abbreviations

The following abbreviations are used in this manuscript:

<table><tr><td>LUNTIAN</td><td>Labor, University, Timber Investment, and Agent-based Nexus</td></tr><tr><td>ABM</td><td>agent-based model</td></tr><tr><td>ITP</td><td>industrial tree plantation</td></tr><tr><td>SES</td><td>Social-Ecological System</td></tr><tr><td>DBH</td><td>diameter at breast height</td></tr></table>

## Appendix A. Parameter Sources and Literature References

<table><caption>Table A1. Parameter sources and literature references for modeling *Swietenia macrophylla*.</caption><thead><tr><th>Modeled Tree Behavior</th><th>Tree Parameters and/or Traits</th><th>Reference(s)</th></tr></thead><tbody><tr><td rowspan="3">Growth</td><td>DBH increment</td><td>[61]</td></tr><tr><td>DBH size classes</td><td>[62]</td></tr><tr><td>Effects of neighborhood interaction</td><td>[63]</td></tr><tr><td>Mortality</td><td>Shade tolerance and DBH</td><td>[64]</td></tr><tr><td>Recruitment</td><td>DBH</td><td>[65]</td></tr></tbody></table>

<table><caption>Table A2. Parameter sources and literature references for modeling Dipterocarp species.</caption><thead><tr><th>Modeled Tree Behavior</th><th>Tree Parameters and/or Traits</th><th>Reference(s)</th></tr></thead><tbody><tr><td rowspan="3">Growth</td><td>DBH increment</td><td>[66]</td></tr><tr><td>DBH size classes</td><td>[67]</td></tr><tr><td>Effects of neighborhood interaction</td><td>[63]</td></tr><tr><td>Mortality</td><td>Stand basal area</td><td>[66]</td></tr><tr><td>Recruitment</td><td>Stand basal area and count of new Dipterocarp trees</td><td>[66]</td></tr></tbody></table>

## References

1. Jurvélius, M. *Labor-Intensive Harvesting of Tree Plantations in the Southern Philippines: A Case Study*; Food and Agriculture Organization of the United Nations: Rome, Italy, 1997.

2. Philippine Council for Agriculture, Aquatic and Natural Resources Research and Development. Industrial Tree Plantation (ITP). 2022. Available online: https://ispweb.pcaarrd.dost.gov.ph/industrial-tree-plantation-itp/ (accessed on 17 February 2025).
3. Official Gazette of the Republic of the Philippines. Presidential Decree No. 705: Revising the Forestry Reform Code of the Philippines. 1975. Available online: https://www.officialgazette.gov.ph/1975/05/19/presidential-decree-no-705-s-1975/ (accessed on 17 February 2025).

4. Niskanen, A.; Saastamoinen, O. *Tree Plantations in the Philippines and Thailand: Economic, Social and Environmental Evaluation*; UNU World Institute for Development Economics Research: Helsinki, Finland, 1996.

5. Landsberg, J.; Waring, R. A generalised model of forest productivity using simplified concepts of radiation-use efficiency, carbon balance and partitioning. *For. Ecol. Manag.* **1997**, *95*, 209–228. [CrossRef]

6. Köhler, P.; Huth, A. The effects of tree species grouping in tropical rainforest modelling: Simulations with the individual-based model Formind. *Ecol. Model.* **1998**, *109*, 301–321. [CrossRef]

7. Mladenoff, D.J. LANDIS and forest landscape models. *Ecol. Model.* **2004**, *180*, 7–19. [CrossRef]

8. University of the Philippines. Gazette. Volume XLIII No. 7 ISSN No. 0115-7450. 2012. Available online: https://osu.up.edu.ph/wp-content/uploads/gazette/2012-DEC.pdf (accessed on 27 February 2025).

9. Supreme Court of the Philippines. Philippine Supreme Court Decision. G.R. No. 52518. August 1991. Available online: https://lawphil.net/judjuris/juri1991/aug1991/gr_1_52518_1991.html (accessed on 18 February 2025).

10. Office of the Secretary of the University (OSU). Sustainable Management and Development Plan for the UP Sierra Madre Land Grants, 1260th Meeting of the Board of Regents, Board Resolution. 24 September 2010. Available online: https://legal.uplb.edu.ph/up-policies/sustainable-management-and-development-plan-for-the-up-mierra-madre-land-grants/ (accessed on 24 June 2025).

11. Philippine Statistics Authority. Region IV-A (Calabarzon): Total Population by Province, City, Municipality and Barangay. 2020. Available online: https://psa.gov.ph/system/files/phcd/2022-12/Region%25204A.xlsx (accessed on 8 July 2021).

12. Manila Bulletin. DAR Identifies UP Property for Possible CARP Coverage. 2020. Available online: https://mb.com.ph/2020/7/8/dar-identifies-up-property-for-possible-carp-coverage (accessed on 5 June 2025).

13. UPLB Legal Office. Historic Trek: UPLB Legal Office Attorney Ventures into Laguna Land Grant. 2023. Available online: https://legal.uplb.edu.ph/updates/historic-trek-uplb-legal-office-attorney-ventures-into-laguna-land-grant/ (accessed on 5 June 2025).

14. Ostrom, E. A General Framework for Analyzing Sustainability of Social-Ecological Systems. *Science* **2009**, *325*, 419–422. [CrossRef]

15. Gonzalez-Redin, J.; Gordon, I.; Polhill, J.; Dawson, T.; Hill, R. Navigating Sustainability: Revealing Hidden Forces in Social-Ecological Systems. *Sustainability* **2024**, *16*, 1132. [CrossRef]

16. Nugroho, S.; Uehara, T. Systematic Review of Agent-Based and System Dynamics Models for Social-Ecological System Case Studies. *Systems* **2023**, *11*, 530. [CrossRef]

17. Newnham, R.M. The Development of a Stand Model for Douglas Fir. Ph.D. Thesis, University of British Columbia, Vancouver, BC, Canada, 1964. [CrossRef]

18. Liu, J.; Ashton, P.S. Individual-based simulation models for forest succession and management. For. Ecol. Manag. **1995**, *73*, 157–175. [CrossRef]

19. Bitterman, P.; Bennett, D. Constructing stability landscapes to identify alternative states in coupled social-ecological agent-based models. *Ecol. Soc.* **2016**, *21*, 21. [CrossRef]

20. Bartkowski, B.; Beckmann, M.; Drechsler, M.; Kaim, A.; Liebelt, V.; Müller, B.; Witing, F.; Strauch, M. Aligning Agent-Based Modeling With Multi-Objective Land-Use Allocation: Identification of Policy Gaps and Feasible Pathways to Biophysically Optimal Landscapes. *Front. Environ. Sci.* **2020**, *8*, 103. [CrossRef]

21. Gonzalez-Redin, J.; Polhill, J.; Dawson, T.; Hill, R.; Gordon, I. Exploring sustainable scenarios in debt-based social-ecological systems: The case for palm oil production in Indonesia. *AMBIO* **2020**, *49*, 1530–1548. [CrossRef]

22. An, L. Modeling human decisions in coupled human and natural systems: Review of agent-based models. *Ecol. Model.* **2012**, *229*, 25–36. [CrossRef]

23. Grimm, V.; Berger, U.; Bastiansen, F.; Eliassen, S.; Ginot, V.; Giske, J.; Goss-Custard, J.; Grand, T.; Heinz, S.K.; Huse, G.; et al. A standard protocol for describing individual-based and agent-based models. *Ecol. Model.* **2006**, *198*, 115–126. [CrossRef]

24. Grimm, V.; Railsback, S.F. Designing, formulating, and communicating agent-based models. In *Agent-Based Models of Geographical Systems*; Springer: Dordrecht, The Netherlands, 2012; pp. 361–377. [CrossRef]

25. Taillandier, P.; Gaudou, B.; Grignard, A.; Huynh, Q.N.; Marilleau, N.; Caillou, P.; Philippon, D.; Drogoul, A. Building, composing and experimenting complex spatial models with the GAMA platform. *GeoInformatica* **2019**, *23*, 299–322. [CrossRef]

26. Grimm, V.; Berger, U.; DeAngelis, D.L.; Polhill, J.G.; Giske, J.; Railsback, S.F. The ODD protocol: A review and first update. *Ecol. Model.* **2010**, *221*, 2760–2768. [CrossRef]

27. Holling, C.S. Resilience and Stability of Ecological Systems. *Annu. Rev. Ecol. Evol. Syst.* **1973**, *4*, 1–23. [CrossRef]

28. Purvis, B.; Mao, Y.; Robinson, D. Three pillars of sustainability: In search of conceptual origins. *Sustain. Sci.* **2019**, *14*, 681–695. [CrossRef]

29. Gigerenzer, G.; Gaissmaier, W. Heuristic Decision Making. *Annu. Rev. Psychol.* **2011**, *62*, 451–482. [CrossRef]
30. Grimm, V.; Railsback, S.F. *Individual-Based Modeling and Ecology*; Princeton University Press: Princeton, NJ, USA, 2005.

31. Becker, G.S. Crime and Punishment: An Economic Approach. *J. Political Econ.* **1968**, *76*, 169–217. [CrossRef]

32. Bedau, M.A. Weak Emergence. *Philos. Perspect.* **1997**, *11*, 375–399. [CrossRef]

33. Harris, I.; Osborn, T.J.; Jones, P.; Lister, D. Version 4 of the CRU TS monthly high-resolution gridded multivariate climate dataset. *Sci. Data* **2020**, *7*, 109. [CrossRef]

34. Fick, S.E.; Hijmans, R.J. WorldClim 2: New 1-km spatial resolution climate surfaces for global land areas. *Int. J. Climatol.* **2017**, *37*, 4302–4315. [CrossRef]

35. Zomer, R.J.; Trabucco, A. *Global Aridity Index and Potential Evapotranspiration (ET0) Climate Database v2*; Figshare: London, UK, 2022. [CrossRef]

36. Bureau of Soils and Water Management. Soil Data Laguna Quezon. 2023. Available online: https://www.foi.gov.ph/ (accessed on 21 March 2021).

37. Ross, C.; Prihodko, L.; Anchang, J.; Kumar, S.; Ji, W.; Hanan, N. *Global Hydrologic Soil Groups (HYSOGs250m) for Curve Number-Based Runoff Modeling*; ORNL DAAC: Oak Ridge, TN, USA, 2018. [CrossRef]

38. Harja, D.; Salvan, G.A.R.; Mulia, R. *Growth and Aboveground Carbon of Trees Estimated Using the FORCASTREE (SExI-FS) Model: Possible Inputs for Land Restoration and Carbon Project in the Philippines*; Technical Report, Technical report prepared for Kennemer Eco-Solutions (KenEco); World Agroforestry (ICRAF) and Kennemer Foods International, Inc.: Los Baños, Philippines, 2023.

39. Forest Management Bureau. Average Domestic Prices 2016. Available online: https://forestry.denr.gov.ph/fmb_web/wp-content/uploads/2023/05/2016-Prices.pdf (accessed on 9 May 2024).

40. The World Bank. Inflation, Consumer Prices (Annual %)—Philippines. Available online: https://data.worldbank.org/indicator (accessed on 10 June 2025).

41. Vanclay, J.K. *Modelling Forest Growth and Yield: Applications to Mixed Tropical Forests*; CAB International: Wallingford, UK, 1994.

42. Arnejo, Z.; Barua, L.; Ramirez, P.J.; Tiburan Jr., C.; Bantayan, N. An Agent-Based Model of a Sustainable Forest Operation in a Theoretical Lowland Dipterocarp Forest Modeled after Mount Makiling Forest Reserve, Philippines. *Forests* **2023**, *14*, 428. [CrossRef]

43. Tomboc, C.C.; Mendoza, G.A. A Simulation Model for the Philippine Selective Logging System: Case Study of a Timber Company. *J. Environ. Manag.* **1993**, *39*, 89–100. [CrossRef]

44. Pinto, R.C.; Pinheiro, C.; Vidal, E.; Schwartz, G. Technical and financial evaluation of enrichment planting in logging gaps with the high-value species *Swietenia macrophylla* and *Handroanthus serratifolius* in the Eastern Amazon. *For. Ecol. Manag.* **2021**, *495*, 119380. [CrossRef]

45. Abidin, S.N.Z.; Jaaman, S.H.; Bakar, A.S.A. The Analysis Level of Optimism that Influence Investor's Risk Tolerance in Asset Allocation. *Sains Malays.* **2023**, *52*, 2999–3011. [CrossRef]

46. Song, X.; Song, Y. Introducing 3-PG2Py, an open-source forest growth model in Python. *Environ. Model. Softw.* **2022**, *150*, 105358. [CrossRef]

47. Kim, S.; Kim, H. A new metric of absolute percentage error for intermittent demand forecasts. *Int. J. Forecast.* **2016**, *32*, 669–679. [CrossRef]

48. Lawrence, K.D.; Klimberg, R.K.; Lawrence, S.M. *Fundamentals of Forecasting Using Excel*; Industrial Press Inc.: South Norwalk, CT, USA, 2009.

49. McHugh, M.L. Interrater reliability: The kappa statistic. *Biochem. Medica* **2012**, *22*, 276–282. [CrossRef]

50. van Vliet, J.; Hagen-Zanker, A.; Hurkens, J.; van Delden, H. A fuzzy set approach to assess the predictive accuracy of land use simulations. *Ecol. Model.* **2013**, *261*, 32–42. [CrossRef]

51. Landis, J.R.; Koch, G.G. The measurement of observer agreement for categorical data. *Biometrics* **1977**, *33*, 159–174. [CrossRef]

52. Reuillon, R.; Chuffart, F.; Leclaire, M.; Faure, T.; Dumoulin, N.; Hill, D.R. Declarative Task Delegation in OpenMOLE. In Proceedings of the International Conference on High Performance Computing and Simulation (HPCS), Caen, France, 28 June–2 July 2010; pp. 55–62.

53. Reuillon, R.; Leclaire, M.; Rey-Coyrehourcq, S. OpenMOLE, a workflow engine specifically tailored for the distributed exploration of simulation models. *Future Gener. Comput. Syst.* **2013**, *29*, 1981–1990. [CrossRef]

54. Hagen-Zanker, A. An improved Fuzzy Kappa statistic that accounts for spatial autocorrelation. *Int. J. Geogr. Inf. Sci.* **2009**, *23*, 61–73. [CrossRef]

55. Beuzen, T.; Marshall, L.; Splinter, K.D. A comparison of methods for discretizing continuous variables in Bayesian Networks. *Environ. Model. Softw.* **2018**, *108*, 61–66. [CrossRef]

56. Saltelli, A. (Ed.) *Global Sensitivity Analysis: The Primer*; John Wiley: Chichester, UK; Hoboken, NJ, USA, 2008. [CrossRef]

57. Saltelli, A.; Annoni, P.; Azzini, I.; Campolongo, F.; Ratto, M.; Tarantola, S. Variance based sensitivity analysis of model output. Design and estimator for the total sensitivity index. *Comput. Phys. Commun.* **2010**, *181*, 259–270. [CrossRef]

58. Thiele, J.C.; Kurth, W.; Grimm, V. Facilitating Parameter Estimation and Sensitivity Analysis of Agent-Based Models: A Cookbook Using NetLogo and 'R'. *J. Artif. Soc. Soc. Simul.* **2014**, *17*, 11. [CrossRef]
59. ten Broeke, G.; van Voorn, G.; Ligtenberg, A. Which Sensitivity Analysis Method Should I Use for My Agent-Based Model? *J. Artif. Soc. Soc. Simul.* **2016**, *19*, 5. [CrossRef]

60. Campbell, S.; Fainstein, S. (Eds.) *Readings in Planning Theory*; Blackwell: Oxford, UK, 2003.

61. Shono, K.; Snook, L. Growth of big-leaf mahogany (*Swietenia macrophylla*) in natural forests in Belize. *J. Trop. For. Sci.* **2006**, *18*, 61–73.

62. Coracero, E.E. Distribution and management of the invasive *Swietenia macrophylla* King (Meliaceae) at the foot of a protected area in Luzon Island, Philippines. *J. Zool. Bot. Gard.* **2023**, *4*, 637–647. [CrossRef]

63. Yue, Q.; Geng, Y.; von Gadow, K.; Fan, C.; Zhang, C.; Zhao, X. Effects of neighborhood interaction on tree growth in a temperate forest following selection harvesting. *Ecol. Indic.* **2022**, *136*, 108663. [CrossRef]

64. Nguyen, H.; Vanclay, J.; Herbohn, J.; Firn, J. Drivers of tree growth, mortality and harvest preferences in species-rich plantations for smallholders and communities in the tropics. *PLoS ONE* **2016**, *11*, e0164957. [CrossRef]

65. Grogan, J.; Landis, R.M.; Free, C.M.; Schulze, M.D.; Lentini, M.; Ashton, M.S. Big-leaf mahogany (*Swietenia macrophylla*) population dynamics and implications for sustainable management. *J. Appl. Ecol.* **2014**, *51*, 664–674. [CrossRef]

66. Nguyen, T.T. Modelling Growth and Yield of Dipterocarp Forests in Central Highlands of Vietnam. Ph.D. Thesis, Technische Universität München, München, Germany, 2009.

67. Galias, D.C.; Cuevas, V.C. The regenerating forest of Magbukún Aeta in Morong, Bataan, Philippines: A biological hotspot for protection and conservation. *Philipp. J. Syst. Biol.* **2018**, *12*, 77–102. [CrossRef]

**Disclaimer/Publisher's Note:** The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.