Sustainable Nitrogen Management in Rice Farming: Spatial Patterns of Nitrogen Availability and Implications for Community-Level Practices

Abstract

Sustainable nitrogen management is crucial for long-term food security and environmental protection in rice farming systems. However, the spatial patterns of nitrogen availability at the community level remain poorly understood, hindering the development of effective sustainable management strategies. This study introduces a novel application of spatial autoregressive analysis to investigate available nitrogen distribution in paddy soils across a rice farming community in Kyoto, Japan. Soil samples from 61 plots, including organically farmed ones, were analyzed for available nitrogen and various physicochemical properties. Contrary to the hypothesis of high variability between adjacent plots, significant positive spatial autocorrelation in available nitrogen was observed, revealing previously unrecognized community-level patterns. The spatial Durbin model outperformed traditional regression approaches and revealed complex spatial interactions in soil properties. Water-soluble organic carbon and humus content showed strong but opposing effects, with a positive direct impact but negative spatial interaction, suggesting topography-driven accumulation processes. Water-soluble nitrogen exhibited reverse patterns with negative direct effects but positive spatial interaction, indicating potential nutrient transport through water movement. These findings highlight the importance of considering both direct and indirect spatial effects in understanding soil fertility patterns, challenging the conventional plot-by-plot management approach. This methodological advancement provides new perspectives for more effective, community-scale soil management strategies in rice farming systems. Moreover, it demonstrates an innovative approach to maximizing the value of outsourced soil analysis data, providing a model for more comprehensive utilization of such data in agricultural research. By enabling more targeted and efficient nitrogen management practices that consider both plot-level processes and landscape-scale interactions, this study potentially contributes to the development of more sustainable and resilient rice production systems.

1. Introduction

Rice cultivation is a cornerstone of global food security, feeding more than half of the world’s population [1]. However, the sustainability of rice production faces significant challenges, particularly in terms of nitrogen management. While the application of inorganic nitrogen fertilizers has significantly contributed to increased rice yields, their low use efficiency poses serious environmental concerns, including greenhouse gas emissions and eutrophication of surrounding waters [2,3].

Crucially, the nitrogen use efficiency of inorganic fertilizer is usually low; rice plants can take up only 30–50% of the inorganic nitrogen applied, substantially relying on soil fertility for nitrogen nutrition [4,5]. This underscores the critical importance of soil fertility and its management in sustainable rice production systems. The capacity of soils to supply nitrogen to rice plants, often evaluated as available nitrogen, is therefore crucial for sustainable rice production across various management practices [6,7].

These challenges have led to growing interest in alternative management practices, including organic rice farming and precision agriculture, as potential pathways toward more sustainable production systems [8,9]. However, the effectiveness of these practices depends on a thorough understanding of soil nitrogen dynamics at relevant spatial scales where management decisions often occur.

Rice farming communities are characterized by a patchy landscape created by bunds and terraces [10,11]. This landscape structure, combined with variations in management practices among farmers, potentially creates complex spatial patterns of soil fertility [12,13,14,15,16]. While discontinuities created by bunds and terraces could lead to significant variations in soil properties between adjacent fields, the actual spatial patterns of nitrogen availability at the community level remain unexplored.

Recent advancements in spatial statistics, particularly spatial autoregressive models, offer new opportunities to analyze complex spatial patterns in agricultural systems [17,18]. While these models have been successfully applied to analyze plant disease incidence at both field [19,20] and community scales [21], their application to soil fertility studies at the community level remains novel.

This study introduces a novel application of spatial autoregressive analysis to investigate available nitrogen distribution in paddy soils across a rice farming community in Kyoto, Japan. We hypothesize that (1) soil fertility, as indicated by available nitrogen, shows distinct spatial patterns across rice fields within a community, and (2) these patterns are associated with soil variables affected by both inherent soil properties and agronomic practices, which can be better understood through spatial autoregressive analysis.

By employing spatial autoregressive analysis, we seek to characterize the spatial distribution of available nitrogen at the community level and identify key factors influencing its variability. This approach aims to inform more sustainable nitrogen management practices and support community-level agricultural planning in rice production systems. Furthermore, it demonstrates an innovative approach to maximizing the value of outsourced soil analysis data, potentially leading to more informed and effective soil management strategies.

2. Materials and Methods

2.1. Study Site and Sampling Design

This study was conducted in Yosano town, northern Kyoto, Japan. The town’s landscape is predominantly forested, with approximately 880 ha of rice fields situated on an alluvial fan formed by the Noda River. These rice fields occupy narrow, flat areas and gentle slopes. Archaeological evidence suggests that rice cultivation in this region dates back more than 1200 years. In recent years, the Yosano municipal government has been promoting organic farming, leading many farmers to apply organic fertilizers to their rice fields. Volunteer farmers were recruited through a public call by the municipal government, including users of Mamekko fertilizer (an organic fertilizer produced locally from rice bran, soy pulp, and fish waste). A total of 61 paddy plots were selected across the town to represent the diversity of management practices and geographical locations (Figure 1).

2.2. Soil Sampling and Analysis

Soil samples were collected in March 2020, prior to the puddling process for the rice-growing season. From each plot, five to nine sub-samples were collected and mixed to create one composite sample representative of the plot. All samples were air-dried at room temperature for approximately two weeks, then sieved through a 2 mm mesh to remove stones and plant residues.

The soil samples were divided into three portions for different analyses:

(1) Available nitrogen: One-third of each sample was examined by the Japan Soil Association (Tokyo, Japan). Each sample was incubated at 30 °C under water-saturated conditions for four weeks. A colorimeter was used to measure the amount of ammonium and nitrate nitrogen mineralized during incubation using the indophenol and diazotization methods, respectively. The sum of ammonium and nitrate nitrogen was determined as available nitrogen.

(2) Chemical and physical properties: A third of each sample was analyzed by the Tokachi Federation of Agricultural Cooperatives (Obihiro, Japan). Available phosphate was extracted using the Truog method and measured using an autoanalyzer. The phosphorus absorption coefficient was determined by measuring the phosphate concentration in the soil solution using the autoanalyzer 0.5 h after the soil was diluted in ammonium phosphate solution. Exchangeable potassium, calcium, and magnesium were extracted in a 1 M ammonium acetate solution and measured using atomic-absorption spectrophotometry. Cation exchange capacity was determined based on all measured cations after extraction in potassium chloride solution. Bulk density was estimated by weighing a 10 mL container filled with each sample. Soil type (alluvial or volcanic), texture (clay loam, sandy loam, or loam), and humus content (rich or poor) were examined visually.

(3) Other soil variables: The remaining third of each sample was examined for pH (1:2.5 method), electrical conductivity (1:5 method), total carbon and nitrogen (dry combustion method using a CHNS elemental analyzer, Vario EL cube, Elementar, Langenselbold, Germany), and water-soluble nitrogen and organic carbon (using a total organic carbon analyzer, TOC-Vwp, Shimadzu, Kyoto, Japan).

A questionnaire was sent to volunteer farmers to obtain information on fertilizer applications (Mamekko, chemical, or organic). A summary of soil analyses is presented in

Table 1. Soil variables and their mean values. Measured variables are shown with their units, while visually examined variables or collected information are expressed as dummy variables. Calculated values are unitless.

Variable	Unit	Mean
available nitrogen	mg 100 g−1	15.84
pH	-	5.59
electric conductivity	mS cm−1	47.07
available phosphate	mg 100 g−1	11.34
phosphate absorption coefficient	-	528.7
exchangeable potassium	mg 100 g−1	10.65
exchangeable calcium	mg 100 g−1	87.52
exchangeable magnesium	mg 100 g−1	41.72
cation exchange capacity	me 100 g−1	11.33
cation saturation	%	48.06
calcium saturation	%	28.51
potassium absorption ratio	-	0.97
calcium-magnesium ratio	-	2.17
magnesium-potassium ratio	-	10.57
total carbon (C)	% w	2.23
total nitrogen (N)	% w	0.21
CN ratio	-	10.63
Water-soluble nitrogen	mg 100 g−1	8.75
Water-soluble organic carbon	mg 100 g−1	102.7
bulk density	g cm−3	0.97
dummy variable of alluvial	1 = alluvial, 0 = volcanic	0.93
dummy variable of clay loam	1 = clay loam, 0 = sandy loam or loam	0.07
dummy variable of sandy loam	1 = sandy loam, 0 = clay loam or loam	0.15
dummy variable of rich in humus	1 = rich in humus, 0 = poor in humus	0.15
dummy variable of Mamekko fertilizer * application	1 = Mamekko, 0 = chemical or organic	0.69
dummy variable of chemical fertilizer application	1 = chemical, 0 = Mamekko or organic	0.18

* Mamekko is an organic fertilizer produced from a mixture of rice bran, soy pulp, and fish waste.

.

2.3. Statistical Analysis

Categorical variables such as soil type, texture, humus content, and fertilizer application were transformed into dummy variables. For soil type, i = 1 if a soil sample was classified as alluvial; otherwise, i = 0 (alluvial dummy). For texture, i = 1 if a soil sample was classified as clay loam; otherwise, i = 0 (clay loam dummy). Similarly, i = 1 if a soil sample was classified as sandy loam; otherwise, i = 0 (sandy loam dummy). For humus content, i = 1 if a soil sample was evaluated as rich; otherwise, i = 0 (rich-in-humus dummy). For fertilizer application, i = 1 if i was applied with Mamekko fertilizer; otherwise, i = 0 (Mamekko fertilizer dummy), and i = 1 if i was applied with chemical fertilizer and i = 0 otherwise (chemical fertilizer dummy).

All variables, including dummy variables, were scaled using the “preProcess” command of the “caret” package in R v.4.0.4. Available nitrogen was regressed on other variables, and parameters were estimated using the ordinary least squares method using Equation (1). Parameter estimation was repeated for different combinations of explanatory variables, and the best combination was selected using Akaike’s information criterion (AIC) and the “step” command in R. Some explanatory variables were further excluded to avoid multicollinearity when the variance inflation factor (VIF) was greater than 10 using the “vif” command of the “car” package in R.

A spatial weights matrix was constructed to quantify the spatial relationships between sampling plots to account for potential spatial dependencies in the data. The current study evaluated three different approaches for constructing spatial weight matrices: radius distance, k-nearest neighbor, and inverse distance methods. The radius distance method, which defines neighboring plots based on a fixed distance threshold, resulted in isolated plots without neighbors. The k-nearest neighbor method, while ensuring each plot had a specified number of neighbors, created asymmetric spatial relationships that were theoretically difficult to justify.

The inverse distance method was ultimately adopted, as it assigns weights that decay continuously with distance. Different power values (1–4) for distance decay and distance thresholds (0.8–2.0 km) were explored within this framework. Although statistical criteria (AIC) suggested optimal values at power = 1 and threshold = 2.0 km (AIC = 105.20), power = 3 and threshold = 0.8 km (AIC = 111.59) were selected based on theoretical and practical considerations. This specification assigns weights as wij = 1/dij³ for plots i and j separated by distance dij when dij ≤ 0.8 km and wij = 0 when dij > 0.8 km. The weights were row-standardized to ensure the sum of weights for each plot equals one.

Using this spatial weight matrix, spatial autocorrelation in available nitrogen was tested using Moran’s I and Geary’s C statistics. Lagrange Multiplier (LM) tests and their robust versions were conducted for both spatial error and spatial lag dependencies, following the decision framework [22]. Based on these diagnostics, four regression models were constructed and compared. The first model is the ordinary least squares (OLS) regression obtained above, expressed as:

$$y=X\beta + \epsilon$$

(1)

where y is the vector of available nitrogen, X is the matrix of explanatory variables, β is the vector of regression parameters, and ε is the vector of independent and identically distributed error terms. The second model is the spatial error model (SEM), expressed as:

$$y=X\beta + u (u=\lambda Wu+\epsilon ) (\epsilon ~\left(0, {\sigma }^2{I}_n\right) )$$

(2)

where u is the spatially autocorrelated error term, λ is the spatial error parameter indicating the strength of spatial autocorrelation in the error terms, W is the spatial weights matrix representing the spatial relationships between observations, and ε is the vector of normally distributed error terms with variance σ²I_n where I_n is a unit matrix. The third model is the spatial lag model (SLM), expressed as:

$$y=\rho Wy+X\beta + \epsilon$$

(3)

where ρ is the spatial lag parameter representing the strength of spatial dependence in the dependent variable, and Wy is the spatially lagged dependent variable. The fourth model is the spatial Durbin model (SDM), expressed as:

$$y=\rho Wy+X\beta + WX\gamma +\epsilon$$

(4)

where γ is the vector of parameters capturing the spatial spillover effects of the explanatory variables, and WX represents the spatially lagged explanatory variables. All parameters were estimated using the maximum likelihood method and the “errorsarlm”, “lagsarlm”, and “spatialreg” commands of the “spdep” package in R. Model selection was based on the AIC values and the results of spatial dependence tests. For the selected model, direct effects (the impact of a variable in a plot on that plot’s available nitrogen) and indirect effects (the impact on neighboring plots’ available nitrogen) were calculated using the ‘impacts’ function in the ‘spdep’ package.

3. Results

3.1. Distribution of Available Nitrogen

The available nitrogen across the 61 sampled plots in Yosano town showed considerable variation (Figure 2). The values ranged from 7.3 to 27.0 mg 100 g⁻¹, with a mean of 15.8 mg 100 g⁻¹. The majority of the plots (approximately 80%) fell within the range of 8 to 20 mg 100 g⁻¹, which is considered suitable for maintaining high productivity in rice production according to the guidelines set by the Ministry of Agriculture, Forestry, and Fisheries in Japan [23].

3.2. Spatial Patterns of Available Nitrogen

Visual inspection of the spatial distribution of available nitrogen (Figure 2) revealed some clustering of similar values. Notably, plots in the small valley in the middle section of the study area (Figure 1c) tended to have higher available nitrogen values (15.5 to 27.04 mg 100 g⁻¹), while plots along the riverside in the southern section (Figure 1d) generally showed lower values (7.3 to 16.1 mg 100 g⁻¹).

This visual observation was confirmed by formal spatial autocorrelation tests (

Table 2. Spatial autocorrelation and dependence tests for available nitrogen.

Test Type	Test Statistic	Statistic
Spatial autocorrelation	Moran’s I	0.415	***
	Geary’s C	0.567	***
Spatial dependence	LM test (error)	4.950	**
	LM test (lag)	5.016	**
	Robust LM test (error)	0.913
	Robust LM test (lag)	0.980
	SARMA	5.929	*

LM: Lagrange Multiplier; SARMA: Spatial AutoRegressive Moving Average. *, **, and *** indicate statistical significance at 10, 5, and 1% levels, respectively.

). Both Moran’s I (0.415, p < 0.001) and Geary’s C (0.567, p < 0.002) indicated significant positive spatial autocorrelation, suggesting that similar available nitrogen values tend to cluster spatially. The spatial dependence diagnostics showed significant LM tests for both error (4.950, p = 0.026) and lag (5.016, p = 0.025) dependencies. However, their robust versions were not significant (error: 0.913, p = 0.339; lag: 0.980, p = 0.322), indicating that simple spatial error or lag specifications might not adequately capture the spatial structure. The marginally significant SARMA test (5.929, p = 0.052) suggested the presence of complex spatial interactions.

3.3. Model Selection and Comparison

Based on the results of spatial dependence diagnostics, four regression models were constructed and compared (

Table 3. Estimates from ordinary least squares and spatial autoregressive models for available nitrogen. Explanatory variables.

	OLS		SEM		SLM		SDM
							Wx
intercept	−0.148		−0.164		−0.144	*	−0.002
electric conductivity	−0.091		−0.036		−0.084		−0.062		−0.127
exchangeable potassium	0.108		−0.024		0.073		0.015		0.325	***
exchangeable calcium	0.167		0.115		0.139		0.085		0.188
calcium saturation	−0.436	***	−0.197		−0.355	***	−0.226		−0.416	**
calcium-magnesium ratio	0.189		0.058		0.143		0.188		0.204
water-soluble organic carbon	0.868	***	1.139	***	0.864	***	1.096	***	−0.876	***
water-soluble nitrogen	−0.381	**	−0.517	***	−0.395	***	−0.437	***	0.355	**
rich in humus	0.581	**	0.646	**	0.578	***	0.864	***	−0.762	*
sandy loam	0.419	*	0.524	*	0.473	**	−0.373		0.337
bulk density	−0.250	**	−0.206	***	−0.235	***	−0.301	***	−0.139
Spatial parameters
λ			0.4048	***
ρ					0.1811		0.237	**
Model fit
AIC	119.27		111.59		115.68		111.54

λ and ρ: Spatial parameters; AIC: Akaike information criterion; *, **, and *** indicate statistical significance at 10, 5, and 1% levels, respectively.

). The AIC values demonstrated that both SDM (111.54) and SEM (111.59) provided substantially better fits than SLM (115.68) and OLS (119.27). The SDM was selected as the most appropriate model due to its ability to capture both direct effects and spatial interactions, supported by its significant spatial parameter (ρ = 0.237, p < 0.01) and its slightly lower AIC value.

3.4. Effects of Soil Properties on Available Nitrogen

3.4.1. Organic Matter and Nitrogen Components

The SDM revealed complex relationships between soil properties and available nitrogen, with several variables showing significant but contrasting direct and spatial interaction effects (Table 3).

Water-soluble organic carbon demonstrated the strongest but opposing effects among all variables. It showed a strong positive direct effect (coefficient = 1.096, p < 0.001), indicating that higher organic carbon content is associated with increased available nitrogen in the same location. However, its negative spatial interaction (coefficient = −0.876, p < 0.001) suggests that higher organic carbon levels in the same location are associated with lower available nitrogen across spatial boundaries in neighboring locations.

Water-soluble nitrogen exhibited a similar pattern of opposing effects but in the reverse direction. A negative direct effect (coefficient = −0.437, p < 0.001) indicates that higher water-soluble nitrogen is associated with lower available nitrogen in the same location, while a positive spatial interaction (coefficient = 0.355, p < 0.01) indicates that higher water-soluble nitrogen levels in the same location are associated with higher available nitrogen across spatial boundaries in the neighboring locations.

Rich-in-humus content showed a strong positive direct effect (coefficient = 0.864, p < 0.001) but a negative spatial interaction (coefficient = −0.762, p < 0.05), suggesting that while humus content directly promotes nitrogen availability locally, it might influence nitrogen dynamics differently at broader spatial scales. This pattern of opposing direct and spatial effects parallels the pattern observed for water-soluble organic content.

3.4.2. Soil Chemical Properties

Several chemical properties showed significant effects primarily through spatial interactions.

Exchangeable potassium demonstrated a significant positive spatial effect (coefficient = 0.325, p < 0.001) without a significant direct effect, suggesting that potassium availability in neighboring areas might influence nitrogen availability through indirect mechanisms.

Calcium saturation showed a significant negative spatial interaction (coefficient = −0.416, p < 0.05), indicating that higher calcium saturation in surrounding areas is associated with reduced nitrogen availability in the focal plots.

3.4.3. Physical Properties

Among physical properties, bulk density was the only variable showing a significant effect, with a negative direct impact (coefficient = −0.301, p < 0.001) but no significant spatial interaction. This suggests that soil compaction primarily affects nitrogen availability through local processes rather than spatial interactions.

3.4.4. Total Effects and Model Performance

The direct, indirect, and total effects estimated from the spatial Durbin model allow for a comprehensive assessment of each variable’s influence on nitrogen availability (

Table 4. Direct, indirect, and total effects from spatial Durbin model for available nitrogen.

Explanatory Variables	Direct	Indirect	Total
electric conductivity	−0.083	−0.163	−0.247
exchangeable potassium	0.066	0.380	0.446
exchangeable calcium	0.117	0.241	0.357
calcium saturation	−0.298	−0.543	−0.841
calcium-magnesium ratio	0.227	0.288	0.515
water-soluble organic carbon	1.001	−0.713	0.288
water-soluble nitrogen	−0.398	0.291	−0.107
rich in humus	0.778	−0.645	0.134
sandy loam	−0.335	0.288	−0.047
bulk density	−0.333	−0.243	−0.577

). The total effects, combining direct and indirect (spatial spillover) effects, revealed that water-soluble organic carbon maintained a positive overall influence (0.288) despite its opposing direct and spatial effects. Rich-in-humus content showed a modest positive total effect (0.134), while bulk density exhibited a consistently negative influence (−0.577). Water-soluble nitrogen demonstrated a small negative total effect (−0.107). Other variables included in the model (electric conductivity, exchangeable calcium, calcium-magnesium ratio, and sandy loam) did not show statistically significant associations with available nitrogen (Table 3). The dummy variables for fertilizer application types (Mamekko and chemical fertilizers) were not retained in the final model, suggesting that the type of fertilizer applied did not have a significant direct effect on available nitrogen when accounting for other soil properties and spatial autocorrelation.

Direct effects represent the impact of a change in an explanatory variable on the dependent variable in the same location. Indirect effects represent the impact on neighboring locations.

4. Discussion

4.1. Spatial Autocorrelation of Available Nitrogen

This study revealed significant spatial autocorrelation in available nitrogen across the rice farming community in Yosano town (Table 2), contrary to our initial hypothesis of high variability between adjacent plots. This finding suggests that soil fertility, as indicated by available nitrogen, tends to be similar among contiguous plots despite the physical barriers created by bunds and terrace boundaries. The positive spatial autocorrelation could be attributed to shared landscape features, similar management practices among neighboring farmers, or underlying geological factors influencing soil properties across adjacent plots.

The clustering of higher available nitrogen values in the small valley and lower values along the riverside suggests that topography plays a crucial role in shaping soil fertility patterns (Figure 2). This aligns with previous studies that have found associations between available nitrogen and landscape position [24]. These findings highlight the importance of considering spatial dependencies when studying soil fertility at the community level, as traditional plot-by-plot analyses may overlook these broader patterns.

4.2. Spatial Interaction in Available Nitrogen

The spatial Durbin model revealed complex relationships between soil properties and available nitrogen, characterized by distinct direct and indirect spatial effects (Table 3 and Table 4). This complexity is particularly evident in the contrasting patterns observed for several key variables.

Water-soluble organic carbon showed strong but opposing effects, with a positive direct impact but a negative spatial interaction. This pattern suggests that while organic carbon enhances nitrogen availability locally, its influence across spatial boundaries may involve different mechanisms. The positive association between water-soluble organic carbon and available nitrogen at the local scale underscores the important influence of easily decomposable organic matter in nitrogen mineralization. This fraction of organic matter, originating from plant residues and organic fertilizers, contains nitrogen-rich compounds that are readily accessible to soil microbes [25]. The negative association across spatial boundaries may indicate topography-driven accumulation of easily decomposable organic matter in lower-lying plots from surrounding areas.

The humus content exhibited similar contrasting patterns, with its positive local association reflecting the contribution of humus, particularly its fulvic acid fraction, to the pool of easily decomposable organic substances that support nitrogen mineralization [26]. The spatial pattern suggests that a similar topographical process to that found in water-soluble organic carbon, the accumulation of organic substances in lower-lying plots from surrounding areas, may influence humus distribution across the landscape. Given that the movement and accumulation of organic matter through topographical processes requires considerable time [27], these spatial patterns of water-soluble organic carbon and humus content may reflect long-term landscape-scale processes rather than short-term dynamics. Higher water content in lower-lying plots may have contributed to the accumulation of organic matter through the slow decomposition rate over the long term [28]. This temporal aspect of organic matter redistribution may have important implications for understanding nitrogen availability patterns across agricultural landscapes.

Water-soluble nitrogen also showed contrasting patterns but in the opposite direction, with a negative direct effect and a positive spatial interaction. The negative local association might reflect nitrogen losses from plots: higher water-soluble nitrogen content could indicate greater potential for nitrogen leaching, leading to reduced available nitrogen [29]. The positive spatial interaction could potentially be explained by the movement of leached nitrogen through subsurface water flow in connected paddy field systems, though this interpretation requires further investigation. While hydrological connectivity between plots exists in paddy systems, the specific mechanisms linking water-soluble nitrogen movement to available nitrogen patterns across plots remain to be elucidated.

Among chemical properties, exchangeable potassium exhibited significant spatial effects without direct effects, while calcium saturation showed negative spatial interaction. These spatial patterns likely reflect the movement of nutrients through the connected paddy field system, where topography and irrigation networks may influence nutrient distribution. Although the precise mechanisms behind these spatial associations between chemical properties and available nitrogen remain to be elucidated, the patterns suggest the importance of considering nutrient interactions at both local and landscape scales [30,31,32].

The soil’s physical properties, particularly bulk density, showed only direct negative effects, suggesting that soil physical conditions primarily influence nitrogen availability through local processes such as the loosening of soils by plowing rather than spatial interactions. This localized effect of bulk density aligns with the understanding that soil structural properties are primarily managed and modified at the individual plot level through farming practices.

4.3. Implication for Sustainable Agriculture

Our findings have several implications for sustainable soil fertility management in rice farming communities. The detected spatial patterns suggest the importance of a multiscale approach to soil fertility management, considering both plot-level processes and landscape-scale interactions [28,30,33]. At the plot level, the positive effect of water-soluble organic carbon and humus content highlights the potential benefits of organic amendments and practices that enhance soil organic matter, which could improve nitrogen use efficiency and reduce reliance on inorganic fertilizers, contributing to more sustainable rice production. However, the negative spatial interactions of these organic components indicate that topographical features should be considered when planning such management practices, particularly in preventing excessive accumulation of organic matter in lower-lying areas.

The spatial patterns observed in both water-soluble nitrogen and chemical properties (exchangeable potassium and calcium saturation) highlight the need for coordinated nutrient management across adjacent fields. This is particularly important in connected paddy field systems where water movement can transport nutrients between plots. Management strategies should account for landscape-scale nutrient movement, potentially involving adjustments to fertilizer application rates based not only on individual plot conditions but also on the plot’s position within the nutrient flow network. The influence of landscape position indicates that management approaches should be tailored to topographical context, with different strategies potentially needed for plots in valleys versus those along hillslopes or near rivers. Farmers could potentially reduce nutrient losses by synchronizing their water management practices and considering these landscape positions when timing fertilizer applications.

The localized effect of bulk density on available nitrogen emphasizes the continued importance of plot-level soil physical management. Maintaining good soil structure through appropriate tillage and the addition of organic matter remains crucial for optimizing nitrogen availability. However, these physical management practices should be integrated within the broader landscape-level strategy, considering how changes in soil structure might affect water movement and nutrient transport between plots.

4.4. Use of Outsourced Soil Analysis Data

This study highlights the potential for maximizing the value of outsourced soil analysis data. In agricultural research, it is common practice to outsource soil analyses to specialized laboratories to ensure precision and consistency. However, the resulting data are often underutilized and frequently presented as simple tables in research papers without in-depth analysis. Our application of spatial autoregressive modeling to these outsourced data demonstrates a novel approach to extracting more meaningful insights from such datasets. This method not only enhances the value of the collected data but also provides a model for other researchers to more fully leverage outsourced soil analysis results in their studies, potentially leading to more informed and effective soil management strategies across different farming systems.

4.5. Limitations

Several limitations of this study warrant consideration in interpreting its results. The data collection was limited to a single growing season, necessitating long-term studies to validate the temporal stability of observed patterns. The visual assessment of some soil properties (texture, humus content) may introduce subjective elements that could affect the reliability of results. Additionally, the specific mechanisms behind some of the observed relationships remain unclear and require further investigation.

Future research should focus on conducting multi-season studies to assess the temporal dynamics of spatial patterns in soil fertility. Integration of high-resolution spatial data (e.g., digital elevation models, remote sensing data) with soil analysis could further elucidate the factors driving spatial patterns in soil fertility. Investigating the relationship between available nitrogen patterns and actual rice yields could provide valuable insights for precision agriculture applications. Moreover, applying similar spatial analysis techniques to other types of outsourced soil data could reveal additional patterns and relationships in various agricultural contexts, potentially leading to more informed and effective soil management strategies across different farming systems.

5. Conclusions

This study applied spatial autoregressive analysis to investigate the distribution of available nitrogen in paddy soils at the rice farming community level in Yosano town, Kyoto, Japan. Our findings reveal significant spatial autocorrelation in available nitrogen across the study area, suggesting that soil fertility patterns are influenced by factors operating at scales larger than individual fields. This insight challenges the conventional approach of managing soil fertility on a plot-by-plot basis and highlights the potential benefits of community-level soil management strategies.

The spatial Darbin model identified several soil properties significantly associated with available nitrogen, characterized by distinct direct and indirect spatial effects. Notably, water-soluble organic carbon and humus content showed positive direct effects but negative spatial interactions, suggesting topography-driven accumulation processes, while water-soluble nitrogen exhibited negative direct effects but positive spatial interactions, indicating potential nutrient transport through water movement. Our methodological approach, combining outsourced soil analysis data with spatial autoregressive modeling, offers a novel perspective on soil fertility patterns in rice farming landscapes. By accounting for spatial dependencies, we were able to reveal relationships that were not apparent in traditional regression analysis.

The findings of this study, while preliminary in nature, have important implications for sustainable soil fertility management in rice farming communities. They suggest that community-level approaches to soil management, considering both plot conditions and landscape position, may be more effective than strategies focused solely on individual plots. The observed spatial patterns could inform more targeted, site-specific nutrient management strategies, potentially improving nitrogen use efficiency and reducing environmental impacts associated with excessive fertilizer use.

Although long-term studies across diverse agricultural settings are needed to validate these findings, this study contributes to a more nuanced understanding of soil fertility patterns in rice farming communities and provides a model for leveraging outsourced soil data more effectively. The insights gained from this research can inform the development of more sustainable nitrogen management practices, contribute to precision agriculture approaches, and support community-level agricultural planning in rice production systems.