Introduction
Trematode parasites present a major global issue for the livestock industry. Specifically, fascioliasis, caused in Europe by the liver fluke Fasciola hepatica (Linnaeus, 1758), was estimated to cost the European agricultural industry €634.49 million per year in 2020 (Charlier et al., Reference Charlier, Rinaldi, Musella, Ploeger, Chartier, Vineer, Hinney, Von Samson-himmelstjerna, Băcescu, Mickiewicz, Mateus, Martinez-valladares, Quealy, Azaizeh, Sekovska, Akkari, Petkevicius, Hektoen, Höglund, Morgan, Bartley and Claerebout2020). Meanwhile, the rumen fluke Calicophoron daubneyi (Dinnik, 1962) has become more prominent in Europe over the past decades and is now widespread (Mage et al., Reference Mage, Bourgne, Toullieu, Rondelaud and Dreyfuss2002; Jones et al., Reference Jones, Brophy, Mitchell and Williams2017a). Furthermore, fascioliasis is also a human zoonotic health issue, though estimates of the number of infected individuals vary broadly between 2.4 and 50 million (Nyindo and Lukambagire, Reference Nyindo and Lukambagire2015).
Diagnosis and treatment of trematode diseases face serious challenges. There is no gold standard diagnostic test, whilst limited anthelmintic options and widespread anthelmintic hamper effective control (Fairweather et al., Reference Fairweather, Brennan, Hanna, Robinson and Skuce2020). Considering this, preventing infection via non-chemical strategies will be essential in controlling the spread and impact of these parasites and of the development of widespread anthelmintic resistance. Such strategies will use land and grazing management practices to target the parasite’s intermediate snail host by either disrupting or eliminating suitable habitat or reducing contact between the snail and grazing livestock at key periods. To successfully implement such strategies, it is imperative that the ecology of the intermediate snail host is well understood.
Galba truncatula (O. F. Müller, 1774) is the main intermediate snail host for F. hepatica and C. daubneyi in Europe (Caron et al., Reference Caron, Martens, Lempereur, Saegerman and Losson2014; Smith et al., Reference Smith, Morgan and Jones2024). The presence of G. truncatula is strongly associated with water-saturated environments. As a result, fasciolosis risk models have largely centred on climatic suitability, particularly moisture and temperature, to infer transmission risk (Rapsch et al., Reference Rapsch, Dahinden, Heinzmann, Torgerson, Braun, Deplazes, Hurni, Bär and Knubben-schweizer2008; Olsen et al., Reference Olsen, Frankena, Bødker, Toft, Thamsborg, Enemark and Halasa2011; Caminade et al., Reference Caminade, Van Dijk, Baylis and Williams2015; Baggenstos et al., Reference Baggenstos, Dahinden, Torgerson, Bär, Rapsch and Knubben-schweizer2016; Innocent et al., Reference Innocent, Gilbert, Jones, McLeod, Gunn, McKendrick and Albon2017; Roessler et al., Reference Roessler, Oehm, Knubben-schweizer and Groll2022). However, a recent systematic literature review has highlighted that soil properties determine habitat suitability for G. truncatula (Smith et al., Reference Smith, Morgan and Jones2024). Several studies noted that high levels of clay within the soil seem to be preferable for G. truncatula (Walton, Reference Walton1918; Charlier et al., Reference Charlier, Soenen, De Roec, Hantson, Ducheyne, Van Coillie, De Wulf, Hendrickx and Vercruysse2014b; Dmitrović et al., Reference Dmitrović, Savić and Pešić2016; Vignoles et al., Reference Vignoles, Dreyfuss and Rondelaud2022). Additionally, the general consensus is that G. truncatula is mostly absent from peatland habitats (Walton, Reference Walton1918; Peters, Reference Peters1938; Heppleston, Reference Heppleston1972; Moens, Reference Moens1981; Dreyfuss et al., Reference Dreyfuss, Vignoles and Rondelaud2018a), albeit these peatland habitats are rarely defined with accuracy in studies. However, the recent study of Nijhuis et al. (Reference Nijhuis, Van Doorn, Holstege, Aalberts, Wagenaar, Antonis and Van Schaik2026) found higher F. hepatica antibody positivity in bulk milk tank samples from Dutch dairy farms on peat soils, demonstrating an urgent need to further investigate and clarify the relationship between G. truncatula and soil organic matter (OM).
Water and soil pH have been extensively studied in relation to G. truncatula habitats, although most research has focused on water pH, despite soil pH being a better indicator of long-term environmental conditions (Hill and Neal, Reference Hill and Neal1997; Ji-dong et al., Reference Ji-dong, Bing-jie, Yong-chun, Ai-jun, Yun-wang, Xian-ju, Hui and Hong-bo2012). There is currently no consensus regarding the species’ pH preferences: some suggest a preference for neutral habitats (Stubbings et al., Reference Stubbings, Skuce and Williams2022) or a positive association with higher pH levels (Charlier et al., Reference Charlier, Soenen, De Roec, Hantson, Ducheyne, Van Coillie, De Wulf, Hendrickx and Vercruysse2014b), while others indicate that certain populations may be adapted to more acidic conditions (Bargues et al., Reference Bargues, Angles, Coello, Artigas, Renata Funatsu, Fernando Cuervo, Buchon and Mas-coma2021) or show limited preference (Moens, Reference Moens1981). Beyond physicochemical properties, habitat structure may also influence suitability by affecting water depth, flow rate and shade. Although there is general agreement that G. truncatula favours water-saturated soil environments with slow-moving water or damp mud surfaces (Dreyfuss et al., Reference Dreyfuss, Vignoles, Rondelaud and Cabaret2015), it remains unclear which habitat types are most suitable. Waterlogged field areas may provide favourable conditions due to reduced water flow, yet the literature review by Smith et al. (Reference Smith, Morgan and Jones2024) found that watercourses, including drainage ditches, streams and rivers, were the most frequently reported habitats.
An additional under-explored area of research is the role of other snail species, which could influence the presence of G. truncatula through predation or competition, or act as alternative intermediate hosts for the parasite F. hepatica (Relf et al., Reference Relf, Good, McCarthy and De Waal2009; Correa et al., Reference Correa, De Meeûs, Dreyfuss, Rondelaud and Hurtrez-boussès2017). Such species could potentially decrease the size of risk areas by limiting the presence of G. truncatula in otherwise unsuitable habitats, or they could expand risk areas by acting as a source of infection in habitats where G. truncatula is absent.
This study aims to investigate the influence of soil properties and habitat type on the likelihood of G. truncatula presence and analyse how interactions between these habitat properties also influence the suitability of habitats for this species. It will also investigate whether other small snail species have a potential impact on the presence or absence of G. truncatula in farmland habitats.
Methodology
Farm visits
Sixteen Welsh sheep farms, selected to capture a broad variety of farming systems and geographies across Wales, were recruited to a wider project assessing F. hepatica control on Welsh farms as described by Rees et al. (Reference Rees, Evans, Davis, Brophy, Morphew, Williams, Williams and Jones2026). To survey G. truncatula presence, these farms were visited either between the 22nd of April and the 31st of May 2024 or the 7th of May and 19th of June 2025, with farms being categorized based on the year in which sampling took place (termed year 1 and year 2 for 2024 and 2025 respectively). The farm visits were undertaken during these periods because late spring is optimal for G. truncatula detection in the project’s locality based on the extensive experiences of the project team (Jones et al., Reference Jones, Williams, Dalesman and Brophy2015, Reference Jones, Williams, Dalesman, Ayodeji, Thomas and Brophy2017b, Reference Jones, Brophy, Davis, Davies, Emberson, Rees Stevens and Williams2018, Reference Jones, Davis, Nalepa-grajcar, Woodruff, Williams, Brophy and Jones2024a). At this time of year, temperatures are consistently high enough for snails to emerge from aestivation and remain highly active, while the risk of prolonged drought conditions that could suppress activity is still low. Each farm was visited once except for farm 4, which had 2 sites that were visited on separate days (labelled as Farm 4A and Farm 4B). As such, there were a total of 17 farm sites.
During visits, water-saturated habitats on pastureland that were accessible to livestock were sampled. These habitats were either identified as perceived fluke risk areas by the farmers during a previous mapping exercise (Rees et al., Reference Rees, Evans, Davis, Brophy, Morphew, Williams, Williams and Jones2026) or were further identified by farmers or researchers during the survey. The locations of habitats were recorded using Global Positioning System coordinates, and habitats were classified as pastureland, watercourse or pond based on which was the predominant characteristic of the habitat. It was also noted whether the habitat also contained any of the other 2 habitat features. For instance, a pastureland habitat with small watercourses running through it would be classed primarily as a pastureland but noted as containing both pastureland and watercourse habitats. The definitions for the 3 habitat types are given in Figure 1. Other habitat features such as poaching, water flow direction and vegetation were also noted.
Table showing the definitions of each of the habitat types used to classify habitats in this study, along with an image example of each type.
Figure 1 Long description
The figure is a table with two columns. The first column is titled “Definition”. The cells in this column give the definition of each of the habitat types. The second column is titled “Example Photographs”. The cells in this column each have a photograph which is an example of the habitat type described in the corresponding cell in the “Definition” column. Below the image is a line of text which states the farm and habitat number that the image is a photograph of. In the next row down there is text in the Definition column which reads as follows: “Watercourse: A habitat feature consisting of a channel, stream or ditch area with flowing or standing water or damp mud. Or A low-lying area without a well-defined ditch that still contains flowing water. A pastureland may be classified as containing a watercourse if it includes any area where water or damp mud forms channels, small connected patches of water, depressions in pastureland that direct water flow when present.” The corresponding image in the “Example Photographs” column shows a shallow muddy channel near the edge of a field with a small area of standing water. The top of the image has a line of hedges and fence posts. Below the image is a line of text which reads “Farm 1, Habitat 8”. In the next row down there is text in the Definition column which reads as follows: “Pond: An area of water at least 2.5m across which contained either standing water, or mostly standing water with running water around inlets/outlets.Or A depression at least 2.5m across surrounded by a relatively steep, defined bank which has damp mud at the bottom.” The corresponding image in the “Example Photographs” column shows a roughly circular body of water in the centre of the image. In the background (the top of the image, behind the body of water), there is a line of trees. In the foreground (the bottom of the image, in front of the body of water) is a wide area of muddy ground with very little vegetation. Below the image, a line of text reads “Farm 3, Habitat 2”. In the next row down there is text in the Definition column which reads as follows: “Pastureland An area of damp mud which does not contain a bank around the edge. It may include areas of standing water smaller than 2.5m wide. Very often poached and may also contain smaller areas that fit into the above categories.” The corresponding image in the “Example Photographs” column shows an area of a field covered in rushes with damp looking soil between the plants. Below the image, a line of text reads “Farm 4A, Habitat 11”.
The boundaries of the habitat were then established, which would define the search areas. For pastureland habitats, this was judged based on whether there was water visible after pushing down on the soil by foot. This method also defined the outer edge of pond habitats, with the inner edge being defined as the transition point between mud and water. This principle also applied to watercourses, such as rivers, to ensure only the shore was sampled. However, in the case of watercourses, an upper and lower point of sampling was also defined along a stretch of habitat. These points were usually defined where there was no visible mud (e.g. the bed became too rocky), too dry (no presence of water underfoot) or if the watercourse was not accessible beyond a certain point (due to fences, dense vegetation, etc.). Habitats which had been identified as perceived fluke risk areas by farmers but that were too dry at the time of surveying (containing no water or damp soil as judged by the above test) were surveyed but were ultimately excluded from the results, as this could be a contributing factor to the absence of snails. As such, all habitats included in the survey were deemed to contain sufficient moisture to support G. truncatula populations at the time of the visit.
After this, around 40 mL of water was collected from each habitat to be used for water pH testing. The habitats were then searched by eye for G. truncatula and other snail species for 10–20 min depending on the size of the habitat to standardize search effort. Habitats in which snails were found were marked as ‘Positive’ for snail presence, while those where snails were not found were marked as ‘Negative’. Specimens of snails believed to be G. truncatula and those of other snails were placed in 2 separate 10 mL bijou tubes (Sarstedt, Nümbrecht, Germany) and taken back to Aberystwyth University. Snail species were then identified using the guidelines of Macan (1977), with the DNA of identified G. truncatula specimens extracted using the Qiagen blood and tissue kit (Qiagen, Hilden, Germany). G. truncatula DNA extracts were pooled for each habitat, and PCR was undertaken following Jones et al. (Reference Jones, Brophy, Davis, Davies, Emberson, Rees Stevens and Williams2018) to specifically amplify a 288 bp amplicon of the G. truncatula ITS2 gene to confirm species identity within each habitat.
Three 0–50 cm soil cores were collected from each habitat and were placed in a single Ziploc bag. Coring locations within each habitat were chosen to cover as much of the habitat as possible while being equidistant from each other with distances being judged by eye. The 0–50 cm cores were used to provide a stable, profile-integrated measure of longer-term soil conditions, while reducing the influence of short-term surface variation caused by recent management or climatic conditions. In pond habitats, the soil cores were taken from around the edge of the body of water. In watercourses, if no deep or running water was present, the cores were taken from the centre of the habitat. Otherwise, they were taken from the bank within 30 cm of the water’s edge where possible. In pastureland habitats, locations were picked based on the longest axis of the habitat so that each core represented one third of the total habitat. Note that due to the irregular shape of habitats, this ‘longest axis’ was not always a straight line. For instance, if the pastureland was a curved or crescent shape, the samples would be taken in the centre of the habitat in a pattern that reflected this curve.
Soil and water processing
To measure water pH, the flacon tube containing the water was shaken for 30 sec and left to rest for 5 min. To calculate soil pH, 10 g of homogenized soil was mixed with 25 mL of deionized water for 30 sec and left to rest for 5 min. In both cases, a FiveEasy pH/mV meter FP20 (Mettler Toledo, Toronto, ON, Canada) was then placed within the slurry, and the pH reading was taken after 30 sec. Two pH readings were taken per sample, and if the readings were more than 0.1 pH apart, a third reading was taken. The mean averages of these readings were then calculated.
The remaining soil was air-dried, weighed, ground using a pestle and mortar and then passed through a 2 mm sieve before being re-weighed to calculate the gravel content. Portions of this sieved soil were then used to measure the soil moisture content, organic content and particle size.
Soil moisture and organic content
To calculate soil moisture content, 10 g of sieved soil was added to a pre-weighed crucible and then placed in an oven at 100 °C for 24 h and then re-weighed. OM content was then calculated by loss on ignition. The oven-dried soil samples were then placed in a furnace at 405 °C for 16 h and then re-weighed once again.
Soil particle size
Particle size analysis was undertaken via the pipette sampling technique, following Rowell (Reference Rowell1994). OM was removed from samples with 6% (v/v) hydrogen peroxide before neutralization with 5% (w/v) sodium hexametaphosphate. The contents were then transferred to 1 L cylinders, which were then filled to the 1 L mark with deionized water before being shaken. Soil particles were allowed to settle, with samples taken at appropriate time points following Stokes’ law. Once samples were extracted, they were oven-dried at 105 °C for 16 h and weighed. Percentages were then derived from total sample weight.
Geographic Information System (GIS) data
Further soil data were attained for each habitat from a GIS dataset provided by the UK Soil Observatory (UKSO) data (UK Soil Observatory, 2014). The following UKSO data were extracted for each habitat in QGIS (QGIS Team, 2026) and included for analysis: Soil Carbon Content, Soil pH, Soil Group (listed as simple soil texture on the website), Soil Texture and Soil Depth. The sources for the UKSO data are listed on their website (UK Soil Observatory, 2014).
Statistical analysis
All statistical analyses were performed using RStudio (R Core Team, 2024), specifically the following packages:
ggplot2 (Wickham et al., Reference Wickham, Chang, Henry, Pedersen, Takahashi, Wilke, Woo, Yutani, Dunnington and Van den Brand2026a), mgcv (Wood, Reference Wood2025), ggtern (Hamilton, Reference Hamilton2025), car (Fox et al., Reference Fox, Weisberg, Price, Adler, Bates, Baud-bovy, Bolker, Ellison, Firth, Friendly, Gorjanc, Graves, Heiberger, Krivitsky, Laboissiere, Maechler, Monette, Murdoch, Nilsson, Ogle, Proctor, Ripley, Short, Venables, Walker, Winsemius and Zeileis2026), pROC (Robin et al., Reference Robin, Turck, Hainard, Tiberti, Lisacek, Sanchez, Müller, Siegert, Doering and Billings2025), glmmTMB (Brooks et al., Reference Brooks, Bolker, Kristensen, Maechler, Magnusson, Skaug, Nielsen, Berg and Van Bentham2026), visreg (Breheny, Reference Breheny2025), plotly (Sievert et al., Reference Sievert, Parmer, Hocking, Chamberlain, Ram, Corvellec, Despouy and Brüggemann2026), DHARMa (Hartig, Reference Hartig2024), emmeans (Lenth et al., Reference Lenth, Piaskowsk, Banfai, Bolker, Buerkner, Giné-vázquez, Hervé, Jung, Love, Miguez, Riebl and Singmann2025), performance (Lüdecke et al., Reference Lüdecke, Makowski, Ben-shachar, Patil, Waggoner, Wiernik, Thériault, Arel-bundock, Jullum, Bacher and Luchman2026), patchwork (Pedersen, Reference Pedersen2025), magrittr (Bache and Wickham, Reference Bache and Wickham2025), cutpointr (Thiele, Reference Thiele2018), multcompView (Graves et al., Reference Graves, Piepho, Selzer and Dorai-raj2026) and dplyr (Wickham et al., Reference Wickham, François, Henry, Müller, Vaughan and Pbc2026b).
Data transformation and clarification
As OM content of habitats was not normally distributed (Shapiro–Wilk values: W = 0.69; p < 0.001), this variable was log-transformed prior to the analysis. To examine associations between soil particle size and G. truncatula presence, isometric log-ratio transformations (ILRs) were used rather than raw percentage values. This approach addresses the constant-sum constraint, where soil percentages are not independent (as an increase in one necessitates a decrease in another). Given that prior research identifies clay as a key indicator of G. truncatula habitat suitability, clay was designated as the primary numerator within the ILRs: ILR 1: clay vs (silt plus sand) (representing the influence of the finest particles) and ILR 2: silt vs sand (representing the balance of the remaining coarser fractions).
Regarding the GIS data, Soil Type, Soil Group and Soil Depth were already presented as named categories, but the Soil pH and Carbon levels were numerical values. However, these numbers were actually categories, which represented the approximate pH and Carbon Content for that area according to the UKSO. As such, these 2 variables were included in the analysis using the ‘factor’ function.
Generalized additive model
Binomial generalized additive models (GAMs) were used to identify factors associated with water-saturated habitats positive for G. truncatula. GAM was selected as some factors, such as soil pH, seemed to have a non-linear relationship with G. truncatula’s presence. In all models, Farm was included as a random effect, while Year was included as a fixed effect (2 levels: 1 and 2) to account for unmeasured temporal variation between survey periods. Initially, univariable models were created to measure the impact of each variable in isolation. A full multivariable model was initially constructed containing all variables analysed using the univariable GAMs. Specifically, these included log OM (log percentage OM content), soil pH, water pH, the log ratio of clay to smaller particles (ILR1), the log ratio of sand to clay particles assuming clay is kept constant (ILR2), soil moisture content, the presence of a pond, pastureland or watercourse element with the habitat and the overall habitat type.
This was followed by a backward elimination process to identify the best-fit model, defined by the lowest Akaike information criterion (AIC) value (Burnham and Anderson, Reference Burnham and Anderson2002). Biologically plausible interaction effects were included alongside the individually tested variables. These included paired interactions between each soil property variable, as well as interactions between habitat type and each soil property. Variables were removed sequentially based on p-values until the GAM contained only significant predictors (p < 0.05). This model was further refined using a stepwise approach, where previously excluded variables and interaction effects were re-introduced to aim to minimize the AIC while ensuring all included terms remained significant (p < 0.05).
The validity of the final GAM model was checked using tests for concurvity, basis dimension adequacy and residual diagnostics, overdispersion, multicollinearity among parametric terms (VIF), correlation among the raw predictors, correlation between smooth contributions and parametric terms and the model predictive performance (ROC/AUC). Furthermore, soil pH was modelled as a smooth term in the GAM and exhibited a clear non-linear threshold-type response. Youden’s J statistic was then used to derive a soil pH threshold that best explained the presence or absence of G. truncatula in habitats that maximizes the sum of sensitivity and specificity.
Water pH, UKSO soil pH and habitat type were significant in univariable GAMs but were excluded from the final multivariable GAM because they were not retained in the lowest-AIC model. These variables were therefore investigated further to assess potential confounding with more strongly correlated environmental predictors. Water pH largely reflected the effect of soil pH and was not examined further. Box-and-whisker plots were used to compare UKSO soil pH categories with measured soil pH values and UKSO carbon content with measured soil OM. To explore habitat-type effects, GLMMs were fitted using the glmmTMB package (Brooks et al., Reference Brooks, Bolker, Kristensen, Maechler, Magnusson, Skaug, Nielsen, Berg and Van Bentham2026), with ILR1, pH and log OM as response variables, habitat type and year as fixed effects, and field nested within farm as a random intercept. Estimated marginal means were calculated using the emmeans package (Lenth et al., Reference Lenth, Piaskowsk, Banfai, Bolker, Buerkner, Giné-vázquez, Hervé, Jung, Love, Miguez, Riebl and Singmann2025) with Tukey-adjusted pairwise comparisons.
Results
Habitat characteristics
A total of 87 habitats from 9 sites were sampled in year 1, and 59 habitats from 8 sites were sampled in year 2. Galba truncatula snails were found in 64/146 (43.83%) of habitats, with species identity confirmed for each habitat via PCR. Snails of the genus Ampullaceana (Linnaeus, 1758) were found in 7/146 (4.79%) habitats, 4 of which were ponds, Potamopyrgus antipodarum (Gray, 1843) were found in 8/146 (5.48%) and Succinea (Draparnaud, 1801) were found in 5/146 (3.42%). A summary of how many habitats there were per farm, and how many were positive for G. truncatula and other snails, is given in Table 1.
Summary of the location (by county), total number of surveyed water-saturated habitats and the number of habitats in which Galba truncatula or other species (specifically A. balthica, Succinea sp. and Potamopyrgus antipodarum) were found within each farm
Table 1 Long description
Table showing the presence of G. truncatula and other freshwater snails in water-saturated habitats on 16 farms across Wales. Columns show farm, county, number of habitats surveyed, number and percentage positive for G. truncatula, and number and percentage positive for other snails (Ampullaceana balthica, Succinea spp., and Potamopyrgus antipodarum). A total of 146 habitats were surveyed. G. truncatula was detected in 64 habitats (43.8%), while other snails were detected in 19 habitats (13.0%). G. truncatula occurrence ranged from 0% (Farm 15) to 80% (Farm 12), with the highest numbers of positive habitats on Farm 6 (11/15) and Farm 7 (8/11). Other snails were absent from surveyed habitats on six farms and most common on Farm 16 (4/9; 44.4%) and Farm 13 (2/3; 66.7%).
Percentage = the percentage of surveyed habitats where these snails were found.
Most habitats only contained a single habitat feature, with only 14/146 (8.97%) having multiple habitat features. As shown in Table 2, Watercourses were the most common feature found in habitats (82/146 of habitats, 56.12%) and the most common overall habitat categorization (71/146, 48.63%). Most habitats in which G. truncatula was found also contained a watercourse (41/64, 64.06%). Galba truncatula was only found in a minority of habitats containing pastureland (23/40, 57.5%) and was also only found in a minority of habitats containing ponds (3/15, 20%).
A full summary of the minimum, maximum, mean and standard deviation for each soil factor across all habitats and those which were positive or negative for snail presence is shown in Table 2, while Figure 2 shows the ratio of clay, silt and sand across positive and negative habitats.
Soil texture triangle showing the proportion of clay, sand and silt in habitats where Galba truncatula was present or absent.
Figure 2 Long description
The figure is a triangle diagram which has the labels Sand (%), Silt (%) and Clay (%) at each of the corners. The percentage of silt runs along the bottom, the percentage of Clay runs from the bottom right corner to the top corner and the percentage sand runs from the top corner to the bottom left. Within the triangle are small blue triangles (representing habitats where G. truncatula was found) and red diamonds (representing habitats where G. truncatula was not found. The position of these shapes represents the percentage clay, silt and sand found within the soil of that specific habitat. Most shapes for both sets of habitats are clustered towards the centre. Habitats where G. truncatula was found are clustered slightly towards having a higher proportion of clay, while habitats where G. truncatula was not found are clustered slightly towards having a higher proportion of sand.
Summary statistics of soil properties of surveyed habitats and their types
Table 2 Long description
Table summarizing environmental characteristics of 146 habitats overall and grouped by the presence or absence of Galba truncatula. The upper section reports minimum, maximum, mean and standard deviation for soil pH, water pH, soil moisture, organic matter, log organic content, and gravel, clay, silt and sand content. Habitats positive for G. truncatula generally had slightly higher soil pH (mean 5.95 vs 5.60), higher clay content (33.07% vs 27.94%), and lower soil moisture (2.00% vs 3.45%), organic matter (10.64% vs 20.78%), and sand content (34.88% vs 39.70%) than negative habitats. Water pH, silt and gravel content were similar between groups. The lower section summarizes habitat categories. Of 146 habitats, 64 (43.84%) were G. truncatula positive and 82 (56.16%) negative. Watercourses had the highest proportion of positive habitats (50%), followed by pastureland (36.7%) and ponds (20%).
Containing = habitat contained that feature, Overall type = the habitat overall was classed as this type. Percentage = percentage of total habitats of that type or which contain that feature. Positive and negative habitats mean that G. truncatula was or was not found within a given habitat.
General additive models
Results of the univariable GAM analysis showed that the habitat type classification of watercourse, the log proportion of clay vs (silt plus sand) (ILR1), and both soil and water pH had a significant positive association with G. truncatula presence (p = 0.01, 0.04, <0.01 and 0.03 respectively), while the log of OM had a significant negative association (p < 0.01). In all non-categorical cases, the effects were significant for these factors regardless of whether data were smoothed or unsmoothed. Galba truncatula was also significantly more likely to be present in habitats within the third and first highest UKSO pH categories (pH 6.31, p = 0.01 and pH 6.04, p = 0.02, respectively). Full results are shown in Table 3. The UKSO data on Carbon Content, Soil Group and Soil Texture produced unstable parameter estimates in GAMs, characterized by extremely large coefficients and standard errors, indicative of quasi-complete separation arising from sparse category frequencies. As these variables could not be reliably estimated or interpreted, they were excluded from subsequent multivariable analyses and are not listed in the univariate analysis results.
Univariate GAM results showing the impact of measured factors on the predicted probability of snail presence
Table 3 Long description
Table summarizing univariable Generalized Additive Model (GAM) results for factors associated with Galba truncatula presence. Columns report fixed-effect estimates, standard errors and p-values, and where applicable, smoothed term effective degrees of freedom (edf), reference degrees of freedom (Ref.df) and p-values. Variables include habitat type, habitat features, soil texture, soil and water chemistry, soil moisture, the presence of other snail species, and GIS-derived soil depth and soil pH. Watercourse habitats had significantly higher predicted probability of G. truncatula presence than pastureland (estimate 1.04, p = 0.01). Significant positive associations were found for ILR1, soil pH, water pH, and GIS soil pH classes 6.04 and 6.31. Log organic matter showed a significant negative association. Soil moisture showed a significant non-linear effect, while gravel, ILR2, habitat features, other snail species and GIS soil depth were not significant.
Note that ‘Type’ refers to what the habitat was classified as overall, while ‘X Present’ means that habitat feature was present within the habitat (i.e. Type: Watercourse means that the habitat overall was classed as a watercourse, while Watercourse Present means whether a watercourse was present in the habitat). ILR1 refers to the log ratio of clay vs (silt plus sand), while ILR2 refers to the proportion of silt vs sand while keeping clay constant. GIS Soil Depth and GIS Soil pH refer to data from the UK Soil Observatory. Std. error = standard error, edf = effective degrees of freedom, Ref.df = reference degrees of freedom.
The multivariable GAM with the lowest AIC value showed that soil pH had a significant non-linear relationship (p = 0.02) with G. truncatula presence (see Figure 3A), log OM and the presence of ponds had a significant negative association with G. truncatula presence (p = 0.03 and 0.02 respectively) with G. truncatula presence (Figure 3B and D respectively), while the proportion of clay vs (silt plus sand) in the soil had a significant (p < 0.01) positive association (Figure 3C). The interaction between soil pH and log OM also had a significant non-linear association (p = 0.03), with OM having a larger impact at higher pH levels as shown by Figure 4. Habitats with high log OM commonly had low pH values, and Pearson correlation showed that soil OM was significantly negatively correlated with soil pH (r = −0.499, p < 0.001) and explained around 25% of soil pH variation (r2 = 0.249). No other interaction effects were found to be significant in multivariable GAM construction and thus were not included in the final model. Some variables were significant during univariable GAMs but were subsequently removed during the backward elimination process discussed above and thus were not present in the final multivariable GAM. The removed variables were habitat type, water pH and UKSO soil pH. The final multivariable GAM model is shown in Table 4.
Visualization of the associations between habitat factors and Galba truncatula presence in the best-fit multivariable GAM model created. Error bars and blue shaded areas around regression lines represent 95% confidence intervals. The small dashes above and below the graphs mark habitats in which snails were and were not found, respectively. Visualizations of predicted presence of G. truncatula vs model fixed/smoothed factors were generated by the best-fit model when one holds the non-focal variables constant. (A) Soil pH, (B) log soil organic matter, (C) log ratio of clay vs (silt plus sand) and (D) presence of a pond in the habitat.
Figure 3 Long description
Visualization of the associations between habitat factors and the predicted probability of Galba truncatula presence (as defined by being found during the habitat surveys) in the best fit multivariable GAM model. There are blue shaded areas around regression lines for Figure 3A, 3B and 3C which represent 95% confidence intervals. For 3D, this is represented by error bars. The small dashes above and below the graphs mark habitats in which snails were and were not found respectively. Visualizations of predicted presence of G. truncatula versus model fixed/smoothed factors were generated by the best-fit model when one holds the non-focal variables constant. Figure 3A is a graph showing the impact of soil pH on the predicted probability of G. truncatula presence within the surveyed habitats according to the best fit Generalized Additive Model. The predicted probability is plotted on the y-axis (going from 0 to 1) while the soil pH is plotted on the x-axis (going from 4.34 to 7.19). The blue regression line starts at a near zero predicted probability, but this rises rapidly at a pH level of around 5, to reach a predicted probability of around 0.8 at a pH level of around 5.5. The regression line then increases more steadily, reaching a predicted probability of nearly 1 at pH levels of around 6.3. The small marks representing habitats where G. truncatula was found tend to be clustered towards higher pH levels than those where G. truncatula was not found. Figure 3B is a graph showing the impact of Organic Matter on the predicted probability of G. truncatula presence within the surveyed habitats according to the best fit Generalized Additive Model. The predicted probability is plotted on the y-axis (going from 0 to 1). The percentage organic matters plotted on the x axis going from 0% to 87.95%. However, the scale is not linear as the model used the log organic matter percentage, but the values have been re-converted back into the actual organic matter (%) values. The regression line begins at over 90% predicted probability and then steadily decreases as Organic Matter percentage increases, dropping to around 0.5 at 87.95% organic matter. The small marks representing habitats where G. truncatula was found tend to be clustered towards lower organic matter percentage levels than those where G. truncatula was not found. Figure 3C is a regression graph showing the impact of the log ratio of clay vs silt+sand on the predicted probability of G. truncatula presence within the surveyed habitats according to the best fit Generalized Additive Model. The predicted probability is plotted on the y-axis (going from 0 to 1) while the log ratio of clay vs silt + sand is plotted on the x-axis (going from 0 to 4, though this is a log scale rather than linear). The blue regression line starts at around a 0.35% predicted probability and steadily rises as the log ratio of clay increases, reaching over 90% when the log ratio is near 4. The small marks representing habitats where G. truncatula was found tend to be clustered towards higher log clay ratios than those where G. truncatula was not found. The box and whisker plot has the predicted probability of G. truncatula presence plotted on the Y-axis from 0 to 1. the x-axis is labelled “Pond Presence”, with two pieces of text above this labelled Absent and Present. Above each of these labels are blue rectangles with a vertical line sticking out the top and bottom (the error bars) and another horizontal line cutting through each box. The horizontal line represents the mean percentage probability of G. truncatula presence in each of these groups. The box for the Absent group is much higher than the box for the Present group.
Interaction between log soil organic matter and soil pH on the predicted probability of Galba truncatula presence. Red areas represent circumstances where the predicted probability of G. truncatula presence is lower, while lighter yellow areas are where the predicted probability is higher. White areas represent a lack of data (i.e. no habitats with log % organic matter and soil pH in that range). The blue contour lines show the gradient of predicted probability of G. truncatula presence; for example, habitats with values between the lines marked 0.8 and 0.9 have an 80–90% predicted probability of G. truncatula presence while holding all other model variables constant.
Figure 4 Long description
This graph shows the relationship between soil organic matter (log-transformed) and soil pH on the predicted probability of G. truncatula presence. Soil organic matter is shown on the y-axis as percentages (0–87.95%), converted from the log scale, while soil pH is plotted on the x-axis. Red areas indicate a low predicted probability of G. truncatula presence, while lighter yellow areas indicate a higher probability. White areas represent combinations of soil pH and organic matter for which no data were available. Blue contour lines show predicted probability gradients; e.g., values between the 0.8 and 0.9 contours indicate an 80–90% predicted probability, with all other model variables held constant. Low pH and high organic matter are associated with lower predicted probabilities, whereas higher pH and lower organic matter are associated with higher probabilities. Above approximately pH 5, the contours become more horizontal, indicating that pH has a stronger effect at lower pH values.
Results of the final general additive model, identifying associations between factors and Galba truncatula presence in a water-saturated habitat
Table 4 Long description
Table summarizing the final Generalized Additive Model for predicting Galba truncatula presence in water-saturated habitats. Columns report parameter estimates, standard errors, effective degrees of freedom (edf), reference degrees of freedom (Ref.df), and p-values. Significant predictors included a negative association with log soil organic matter (estimate −0.84, p = 0.03), a negative association with pond presence (estimate −2.07, p = 0.02), and a positive association with ILR1, representing the proportion of clay relative to silt and sand (estimate 2.14, p < 0.01). Soil pH had a significant smoothed effect (edf 4.09, p = 0.02), and the smoothed interaction between soil pH and log organic matter was also significant (p = 0.03). Farm was included as a significant random effect (edf 9.21, p = 0.01). Year was not a significant predictor (p = 0.31).
Log organic matter = the natural log of the % organic content, Pond = the presence of a pond within the habitat, ILR1 = the proportion of clay vs (silt plus sand) in the soil, Soil pH = the average soil pH measured for that habitat. Std. error = standard error, edf = effective degrees of freedom, Ref.df = reference degrees of freedom.
Post hoc checks for adequacy, residual diagnostics, overdispersion, multicollinearity, correlations among raw predictors, correlations between smooth contributions and parametric terms, and predictive performance were all within acceptable limits. Concurvity results were relatively high (<0.7) for both the interaction effect and soil pH. Further analysis showed that most of the concurvity between these terms was due to their inherent interaction, which is expected. When one of the terms was removed, the concurvity score for the other dropped below an acceptable threshold of 0.5. Above this threshold, concurvity values are often seen as a cause for concern (Lai et al., Reference Lai, Tang, Li, Zhang and Mao2024). Application of Youden’s J statistic indicated that the likelihood of G. truncatula presence decreased in habitats with a soil pH below 5.51 (sensitivity = 0.78, specificity = 0.46).
The full post hoc analysis of the final multivariable GAM model used for this paper can be found in Supplementary document 1. Specifically, concurvity results are found in Table S1, convergence diagnostics checks based on the GAM check function are found in Table S2, correlation between original predictors is found in Table S3, correlation between smooth contributions and parametric terms is found in Table S4 and random effect by Farm is found in Table S5. Other post hoc test results such as overdispersion and model predictive performance are typed in-text.
UKSO and habitat type analysis
The box and whisker plot created to compare the UKSO and sampled pH data (Figure 5) demonstrated that the habitats with the highest measured pH tended to be in the highest UKSO soil category. However, the category names did not closely match the soil pH measured from habitats, and most of the categories showed extremely high variation. As for carbon content, there seemed to be almost no observable correlation beyond a large plurality of results being grouped into one UKSO category (62.8). This is likely because the UKSO data were derived from 1 km2 areas, which could contain multiple habitats (especially if they were in the same field). Given the discrepancies between categories and measured values and this large variation, it was determined that any statistical analysis to determine how much variation in measured soil pH was explained by the UKSO categories was not needed.
Box and whisker plots comparing how the values of soil properties measured from samples taken from habitats compare with how those habitats are categorized in the UK Soil Observatory maps. (A) Soil pH and (B) soil organic matter/carbon content. G. truncatula presence refers to whether G. truncatula was found within a given habitat. The numbers along the bottom are not on a continuous scale, but instead are categories in which all habitats within a 1 km2 area are grouped.
Figure 5 Long description
Both Figure 6A and Figure 6B have points marked for each habitat within the box and whisker plots. The categories according to the UK Soil Observatory maps are referred to as “GIS category” Points marked with a circle represent habitats where G. truncatula was found while points marked with an X represent habitats where G. truncatula was not found. Figure 6A shows the measured soil pH along the y-axis ranging from 4.34 to 7.19, with the GIS soil pH categories along the x-axis. In Figure 6B, the measured soil organic matter was on a scale of 0 to 100%, with the GIS soil carbon categories along the x-axis. In both cases, the graphs show no true pattern between the GIS categories and the measured soil pH and organic matter. The only real exception to this was the highest GIS soil pH category, which did have the highest measured pH values. However, the other GIS soil pH categories did not follow this trend.
Estimated marginal means derived from the generalized linear mixed models indicated that several soil properties differed among habitat types after accounting for variation at the farm and field levels (Figure 6). Watercourses had significantly higher soil pH compared to pastureland habitats (estimate = −0.23, p < 0.01), whereas no significant differences were detected between pasturelands and ponds or between ponds and watercourses (Figure 6A). Log-transformed soil OM was significantly higher in pastureland habitats than in watercourses (estimate = 0.77, p < 0.01) and ponds (estimate = 0.53, p = 0.04), while ponds and watercourses did not differ significantly (Figure 6B). The log ratio of clay relative to larger soil particles (silt + sand) was significantly greater in pasturelands than in watercourses (estimate = 0.21, p < 0.01), with no significant differences between ponds and watercourses or between ponds and pasturelands (Figure 6C). Full pairwise comparison results are provided in Table S6 in Supplementary document 1.
Differences in predicted (A) soil pH, (B) log organic matter content and (C) log proportion of clay vs (silt plus sand) by habitat type, calculated by generalized linear mixed models accounting for farm- and field-level clustering. In all cases, significant differences were found between pastureland and watercourses (p < 0.01 for all cases), and pastureland vs pond also had significantly different levels of predicted log organic matter (p = 0.04). The black diamonds are the estimated marginal means within each habitat type, while the blue bars are the 95% confidence intervals for those estimated means. The letters within each figure show whether pairs of model-estimated marginal means differ significantly, after multiple-comparison correction. Compact letter displays indicate groups that do not differ significantly based on Sidak-adjusted pairwise comparisons of estimated marginal means.
Figure 6 Long description
The three bar graphs have the predicted variables plotted along the x axis, while the y axis has three labels: Watercourse, Pond and Pastureland. A black diamond shows the estimated marginal mean of the predicted value of each of the factors, while a blue bar on either side shows the 95% confidence interval of that mean. Going by estimated marginal means, watercourses were found to have the highest predicted soil pH, lowest predicted log organic matter and lowest log ratio of clay vs (silt + sand). Pastureland habitats showed the opposite trends and ponds were an intermediate result with the largest confidence intervals.
Discussion
This study aimed to investigate the impact of soil properties and habitat type on the presence of G. truncatula in water-saturated habitats on farms. Overall, our results show that soil properties are strongly associated with G. truncatula presence; specifically, habitats with soils that have pH above a threshold of 5.51, are low in OM or have a high proportion of clay are significantly more likely to be inhabited by G. truncatula. Furthermore, analysis showed a significant interaction between soil pH and OM, with OM content more strongly associated with G. truncatula presence in habitats with a soil pH above around 5.5. Results also indicated that habitat type may play a significant role in suitability for G. truncatula, as the presence of ponds in habitats had a significant negative association with the presence of G. truncatula, while watercourse habitat soils had significantly lower OM and higher pH levels, factors found to be preferred by G. truncatula. These results provide a crucial step in determining which habitats are likely to be risk areas for fluke transmission in livestock and may provide a framework for future modelling of these areas and the likelihood of fluke transmission at farm-level.
Currently, many fascioliasis risk models have focused primarily on climatic factors rather than soil properties (Olsen et al., Reference Olsen, Frankena, Bødker, Toft, Thamsborg, Enemark and Halasa2011; Caminade et al., Reference Caminade, Van Dijk, Baylis and Williams2015; Baggenstos et al., Reference Baggenstos, Dahinden, Torgerson, Bär, Rapsch and Knubben-schweizer2016; Innocent et al., Reference Innocent, Gilbert, Jones, McLeod, Gunn, McKendrick and Albon2017). Although soil type has sometimes been included to represent drainage or waterlogging potential (Baggenstos et al., Reference Baggenstos, Dahinden, Torgerson, Bär, Rapsch and Knubben-schweizer2016), the direct influence of soil composition on snail presence has received little attention. Where soil variables are considered, they are typically modelled using linear or log-linear approaches (Charlier et al., Reference Charlier, Bennema, Caron, Counotte, Ducheyne, Hendrickx and Vercruysse2011; Novobilský et al., Reference Novobilský, Novák, Björkman and Höglund2012), which cannot easily capture non-linear relationships or interactions among edaphic factors, or have used simple characterizations of soil types (Rapsch et al., Reference Rapsch, Dahinden, Heinzmann, Torgerson, Braun, Deplazes, Hurni, Bär and Knubben-schweizer2008). General additive models have only recently been used to investigate these complex patterns, including in spatio-temporal fluke risk models incorporating soil particle composition and pH (Priddle et al., Reference Priddle, Swindells, Elliott, Ralph, Duwalage, Roberts and Mengersen2026). However, the use of such models to direct farm-scale analysis of G. truncatula presence in the UK is novel and has allowed us to move beyond simple linear assumptions to reveal complex habitat drivers.
The positive association between clay and G. truncatula is expected given prior research (Walton, Reference Walton1918; Charlier et al., Reference Charlier, Vercruysse, Morgan, Van Dijk and Williams2014a; Knubben-Schweizer and Torgerson, Reference Knubben-schweizer and Torgerson2015; Dmitrović et al., Reference Dmitrović, Savić and Pešić2016; Vignoles et al., Reference Vignoles, Dreyfuss and Rondelaud2022). The likely primary cause is that the small size of clay soil particles compared to sand or silt promotes poor drainage and greater water retention (Tsubo et al., Reference Tsubo, Fukai, Basnayake, To, Bouman and Harnpichitvitaya2007; Szejba, Reference Szejba2020). This causes long-lasting areas of pooled shallow standing water or damp mud, which could also be beneficial for the snails during periods of drought. Clay-rich soils may also facilitate burrowing during aestivation (Vignoles et al., Reference Vignoles, Dreyfuss and Rondelaud2022) and act as a buffer against extreme pH conditions, helping maintain pH above the minimum threshold for G. truncatula (Weaver et al., Reference Weaver, Kissel, Chen, West, Adkins, Rickman and Luvall2004). However, while both clay content and soil pH were significant in the final multivariable model, their interaction was not, suggesting largely independent effects on snail presence. Clay particles are also more readily transported downstream in faster-flowing environments (Church, Reference Church2006), which may explain the lower clay content observed in watercourse habitats. Despite this, G. truncatula was more frequently present in watercourses, suggesting that the benefits of these habitats may outweigh their lower clay content, while clay-rich pastureland may compensate for otherwise suboptimal conditions through enhanced water retention.
Regarding the impact of pH on G. truncatula, prior studies had tended to show a minimum pH tolerance limit of around 5.0 (Smith et al., Reference Smith, Morgan and Jones2024), which roughly aligns our data and suggests that the probability of G. truncatula presence is only markedly affected by soil pH levels below approximately 5.51. This was further evidenced by no snails being found in habitats with a soil pH below 5.16 in our study (see Table 2). However, our results primarily reflect the lower half of the pH tolerance range for G. truncatula, and it remains unclear whether similar patterns hold in soils with higher pH. Gastropods are particularly sensitive to changes in pH because their shells are made from calcium carbonate and low pH environments can result in reduced shell growth, increased shell erosion and deformities (Marshall et al., Reference Marshall, Abdelhady, Wah, Mustapha, Gӧdeke, De Silva and Hall-spencer2019).
Potential methods to decrease soil pH include the targeted application of elemental sulphur (Zhao et al., Reference Zhao, Gupta, Degryse and McLaughlin2017) or ammonium-based fertilizers (Johnson et al., Reference Johnson, Fisher, Huang, Vetanovetz and Argo2010). Such practices in theory could reduce fluke risk from an area if pH could be lowered below the G. truncatula tolerance threshold. However, soil pH is an important indicator of soil and grassland health as it influences nutrient availability (Raniolo et al., Reference Raniolo, Maretto, Benedetti Del Rio, Cournut, Cremilleux, Nowak, Michaud, Lind, Concheri, Stevanato, Squartini, Ramanzin and Sturaro2023; Najdenko et al., Reference Najdenko, Lorenz, Dittert and Olfs2024), and as such practices such as liming are commonly used to reduce soil acidity and improve soil conditions (Wenyika et al., Reference Wenyika, Enesi, Gorim and Dyck2025). Consequently, liming highly acidic land may inadvertently increase fluke transmission risk by creating more favourable conditions for G. truncatula. However, this has not been experimentally demonstrated. Anecdotal historical accounts from Wales reported lime application being blamed for fascioliasis outbreaks in sheep (Walton, Reference Walton1918), which supports this hypothesis.
An alternative to chemical acidification of land is the creation and preservation of peatland, which tends to have acidic conditions hostile to G. truncatula (Paleckiene et al., Reference Paleckiene, Navikaite and Slinksiene2021; Pastukhov et al., Reference Pastukhov, Kovaleva and Kaverin2022; Jandl et al., Reference Jandl, Negassa, Eckhardt and Leinweber2023). This aversion to peatland has been noted by several prior studies (Walton, Reference Walton1918; Peters, Reference Peters1938; Heppleston, Reference Heppleston1972; Moens, Reference Moens1981; Dreyfuss et al., Reference Dreyfuss, Vignoles and Rondelaud2018a) and is supported by the negative associations of G. truncatula presence with high OM and low soil pH habitats in this study. However, although peatland restoration is a policy across many European countries, introducing large areas of peatland onto farmland may be impractical. Farmers frequently express concerns that rewetting may increase fluke risk (NatureScot, 2026), and there are also practical limitations for farmers, such as the reduction in usable land area, which can cause them a financial loss (Niemi et al., Reference Niemi, Mattila and Seppälä2024).
The aversion of G. truncatula to peatland likely explains part of the reason for the negative association between OM content and G. truncatula presence. Peatlands form when underlying hydrological conditions promote prolonged waterlogging. This creates anaerobic conditions, incomplete decomposition and the accumulation of OM and compounds such as humic and fulvic acids, which contribute to lower soil pH (Gondar et al., Reference Gondar, Lopez, Fiol, Antelo and Arce2005), which in many cases may be lower than G. truncatula’s minimum soil pH threshold. Such environments also have low levels of dissolved oxygen and key nutrients, both of which may also contribute to low survivability of G. truncatula (Heppleston, Reference Heppleston1972; Jigyasu and Singh, Reference Jigyasu and Singh2010). Indeed, G. truncatula has been shown to be able to quickly colonize peatland habitats following fertilizer application, possibly due to the rapid algal growth which occurred as a result (Comerton and Houghton, Reference Comerton and Houghton1975). As such, the negative association of G. truncatula and OM at very high levels (e.g. 50% + OM) is likely due to these habitats being highly acidic, low in oxygen and nutrient-poor.
That said, our results demonstrated a negative linear relationship between soil OM and G. truncatula presence that extends beyond peatland soils across managed, improved pastures. Furthermore, the significant interaction effect of OM and soil pH shows that the association between OM and G. truncatula presence strengthens as soil pH increases. This finding suggests that while soil pH and OM are linked, neither factor alone fully explains the observed patterns, highlighting the need to consider additional ecological or environmental mechanisms.
Waterlogged improved pastures with high OM may be more prone than lower OM waterlogged pastures to hypoxic or anaerobic conditions, as microbial decomposition of high levels of OM increases oxygen demand (Siedt et al., Reference Siedt, Teggers, Linnemann, Schäffer and Van Dongen2023). High-OM areas may also be associated with denser vegetation or accumulations of decaying plant material, such as trampled grass or leaves (Billings, Reference Billings2006; Laganière et al., Reference Laganière, Augusto, Hatten and Spielvogel2022). These features could reduce habitat suitability by shading or smothering the soil surface, thereby limiting algal photosynthesis and food availability for G. truncatula (Dreyfuss et al., Reference Dreyfuss, Vignoles, Rondelaud and Cabaret2015; Roessler et al., Reference Roessler, Oehm, Knubben-schweizer and Groll2022; Jones et al., Reference Jones, Davis, Nalepa-grajcar, Woodruff, Williams, Brophy and Jones2024a),
The direct impact of increased OM on algal communities is complex. Increased OM has sometimes been negatively correlated with algal abundance (Wei et al., Reference Wei, Zhao, Chen, Sun, Zhang and Chen2023), and in other cases there was a positive correlation (Chang et al., Reference Chang, Gao, Zamyadi, Wang and Li2025). Furthermore, the decomposition of some plant matter can result in the release of algicidal compounds, which can inhibit cyanobacteria (Świerk et al., Reference Świerk, Celewicz, Krzyżaniak, Antoszewski, Stuper-szablewska, Szablewski, Kurasiak-popowska, Kosiada, Stoyneva-gärtner and Krawiec2025). As such, it may be the case that increased soil OM actually decreases the availability of G. truncatula’s food sources. However, our understanding of the dietary preferences of G. truncatula is limited beyond its general preference to unicellular algae (Kendall, Reference Kendall1953; Dreyfuss et al., Reference Dreyfuss, Vignoles, Rondelaud and Cabaret2015), and no thorough investigation has yet compared the diet of G. truncatula in the field with the algal communities present in its potential habitats. As such, it is unclear whether peatland algal communities differ from those typically consumed by G. truncatula.
Another complicating factor is the role of calcium-organic compounds. G. truncatula, like most snails, has a high physiological demand for calcium, with calcium concentrations in water positively associated with the snail’s presence (Olkeba et al., Reference Olkeba, Boets, Mereta, Yeshigeta, Akessa, Ambelu and Goethals2020; Şahin and Zeybek, Reference Şahin and Zeybek2020). Calcium can promote the accumulation and stabilization of soil OM (Rowley et al., Reference Rowley, Grand and Verrecchia2018, Reference Rowley, Grand, Spangenberg and Verrecchia2021), while high-OM soils may bind calcium within organo-mineral complexes (Rowley et al., Reference Rowley, Nico, Bone, Marcus, Pegoraro, Castanha, Kang, Bhattacharyya, Torn and Peña2023, Reference Rowley, Pena, Marcus, Porras, Pegoraro, Zosso, Ofiti, Wiesenberg, Schmidt, Torn and Nico2025). Given the positive association between calcium availability and G. truncatula presence, the observed negative relationship between soil OM and snail occurrence is therefore unexpected. One possible explanation is that, in OM-rich soils, a substantial proportion of calcium may be bound within stable calcium–organic or calcium (bi)carbonate–organic complexes that reduce its bioavailability despite high total calcium concentrations (Qafoku et al., Reference Qafoku, Battu, Varga, Marcus, O’callahan, Zhao, Mergelsberg, Kew, Loring, Qafoku and Leichty2023). The extent of this effect likely depends on OM type, as more humified materials bind calcium more strongly than less humified materials such as fresh manure (Masmoudi et al., Reference Masmoudi, Magdich, Rigane, Medhioub, Rebai and Ammar2020). As a result, although higher OM levels in water-saturated environments may correlate with increased total calcium concentrations, a substantial proportion of this calcium may occur in forms that are less accessible to G. truncatula. Nevertheless, this explanation remains speculative and highlights the need for further investigation into how calcium is stored and partitioned in potential G. truncatula habitats with varying quantities and types of soil OM.
Although habitat type was not found to be a significant factor in the main multivariable GAM, our results show that habitat type can significantly influence soil conditions relevant to G. truncatula suitability, even after accounting for variation between farms and fields. Given that watercourses were more commonly positive for G. truncatula, this indicates that some drainage systems may pose a risk to livestock. However, this may depend on the type of drainage system used. Open systems such as drainage ditches are also commonly reported as habitats for G. truncatula across multiple studies (Smith et al., Reference Smith, Morgan and Jones2024), but subsurface drainage systems have previously been linked to declining populations of G. truncatula (Dreyfuss et al., Reference Dreyfuss, Vignoles and Rondelaud2018b). Even if watercourses are more likely to be inhabited by G. truncatula, whether they are a risk area also depends on whether they are accessible to livestock. From the observation of notable habitat features and information given by farmers, 96% of watercourses had signs of frequent use by livestock (e.g. primary water source, hoofprints, wool), indicating that contact with these habitats is likely on many farms in Wales. Watercourses may also act as a refuge for the snail in areas that are otherwise unsuitable, as our results show these habitats tended to have lower OM and higher soil pH compared to surrounding pastureland. As such, in fields which have high OM and low pH soils, watercourses may act as a conduit of infection in an area that was otherwise low risk.
This may also explain the apparent contradiction between our findings and those of Nijhuis et al. (Reference Nijhuis, Van Doorn, Holstege, Aalberts, Wagenaar, Antonis and Van Schaik2026), where liver fluke infections in Dutch dairy herds were positively associated with peat soils. In the extensive peat meadows of the Netherlands, G. truncatula snails are predominantly found in or nearby watercourses, with these habitats demonstrated to be the key driver for F. hepatica transmission on some farms (Gaasenbeek et al., Reference Gaasenbeek, Over, Noorman and De Leeuw1992). Based on our results, these watercourse habitats may have lower OM content than the surrounding pasture, making them more suitable for G. truncatula despite the wider peatland context. However, managing these habitats presents clear challenges, as fencing watercourses is often impractical or costly, particularly when they serve as essential water sources, and field-level land management aimed at reducing snail habitat may be undermined if G. truncatula persists within adjacent watercourses.
Apart from changing the hydrology and number of suitable habitats, other strategies to reduce fluke transmission may rely on changing the properties of the soil. Drastic decreases in soil pH may be an impractical solution, as would directly reducing the proportion of clay in soil. OM enrichment of soil may be more feasible. However, this is reliant on OM having a direct negative impact on the suitability of habitats for G. truncatula, and currently the mechanisms behind the negative association are not well understood. Provided that OM does directly impact habitat suitability, soil OM content has numerous other benefits such as the detoxification of various pollutants, increasing soil fertility and increased carbon storage (Gerke, Reference Gerke2022). Current UK farm support schemes, including the Sustainable Farming Scheme in Wales (Welsh Government, 2025) and the Environmental Land Management scheme in England (Department for Environment, 2023), provide direct financial incentives for building soil carbon, potentially supplemented by emerging private markets for carbon credits in future (Jones et al., Reference Jones, Paige, Smith, Worth, Betts and Stafford2024b).
Consequently, increasing OM through regenerative grazing, the addition of organic fertilizer or soil management may have further benefits in suppressing snail populations without requiring full peatland restoration or the associated loss of agricultural productivity. However, the specific mechanisms underlying the apparent detrimental impact of increased soil OM on G. truncatula remain to be fully resolved. Given this uncertainty, it is also important to consider whether other soil amendments with known effects on soil chemistry and structure could influence snail habitat suitability. As an example, biochar is known to increase aeration and improve drainage in fine-textured clay soils (Castellini et al., Reference Castellini, Giglio, Niedda, Palumbo and Ventrella2015; Obia et al., Reference Obia, Mulder, Hale, Nurida and Cornelissen2018) and is also reported to suppress soil-borne pathogens and pests (Iacomino et al., Reference Iacomino, Idbella, Laudonia, Vinale and Bonanomi2022). However, care will be needed when implementing soil amendments, or any land management strategies generally, to ensure no unintended detrimental effects are seen either on fluke transmission or generally on pasture productivity and farmland biodiversity.
The success of such strategies targeting G. truncatula may also be influenced by the presence of other snail species acting as alternative intermediate hosts. Other snail species were found during this study, albeit only in a small number of habitats. These included A. balthica, a snail known to be capable of acting as an intermediate snail host for F. hepatica. The importance of A. balthica as an intermediate snail host is disputed, as its very low susceptibility and minimal cercarial output suggest that it plays little to no epidemiological role in fluke transmission (Caron et al., Reference Caron, Lasri and Losson2007; Rondelaud et al., Reference Rondelaud, Titi, Vignoles, Mekroud and Dreyfuss2014; Ciuca et al., Reference Ciuca, Maurelli, Bosco, Hammami, Vitiello, Sengupta, Stensgaard and Rinaldi2026). Nevertheless, it has been suggested that A. balthica snails can inhabit environments that are more acidic than G. truncatula can tolerate, explaining why liver fluke is seen on an acidic hill farm in Ireland (Relf et al., Reference Relf, Good, McCarthy and De Waal2009). In this study, the average soil pH of habitat in which A. balthica was found was lower than for G. truncatula (see Supplementary Table S7). However, the overall range of soil pH in which G. truncatula was found extended to lower pH levels than for A. balthica, and based on our findings, it is feasible that the drainage ditches where A. balthica was recorded by Relf et al. (Reference Relf, Good, McCarthy and De Waal2009) were less acidic than the surrounding farm environment.
However, one clear niche for A. balthica on surveyed farms was ponds, habitats where G. truncatula was less likely to be found, probably because A. balthica is a more aquatic snail compared to the amphibious G. truncatula and, as such, deeper water would be favourable. If A. balthica is a more important intermediate host than currently assumed, such habitats could play a greater role in F. hepatica transmission. However, ponds are declining on UK farms (Smith et al., Reference Smith, Clarke, Weldon and Robson2022; Robson and Von Benzaon, Reference Robson and Von Benzaon2025). To counter this, there is renewed interest in pond creation, including scrapes, a very shallow temporary pond which may be more suited to G. truncatula (Cuthill, Reference Cuthill2020) through current farm support schemes, which may shift their importance as fluke transmission sites in future (Jones et al., Reference Jones, Rees, Brophy, Williams, Williams, Morphew, Davis and Evans2024c). Larger ponds, if capable of supporting more complex communities, might increase competition with and predation of G. truncatula, both through resident biodiversity and by attracting mobile animals such as snail-feeding birds. In particular, establishing wader scrapes on high OM or highly acidic soils may enhance local biodiversity and attract snail-feeding birds to a farm generally, without increasing the risk of fluke transmission. Evidence for this is lacking, but the information on snail habits from this and other studies could help optimize the design of such actions.
These potential ecological and management implications should, however, be interpreted within the geographical and methodological scope of the study. Sampling was restricted to water-saturated habitats in Wales, and the observed associations may not apply unchanged to other regions, farming systems, soil types or G. truncatula populations. Snail presence was assessed using standardized short-duration searches, while soil properties were measured from 0 to 50 cm cores and therefore represent the wider soil profile rather than the immediate surface layer occupied by G. truncatula. Further studies across a broader geographical range, incorporating repeated snail surveys and depth-resolved soil sampling, are therefore required to assess the general applicability of these relationships. These limitations are particularly relevant when considering the wider use of farm-level habitat models to remotely map and predict fluke infection risk. Despite our study’s success in modelling G. truncatula habitats at farm level, challenges remain in being able to utilize such models in future to remotely map and model fluke infection risk areas. Although some soil pH classes derived from the UKSO were significantly associated with G. truncatula presence in the univariable GAM analyses, indicating some predictive value, these pH data were generally inaccurate when compared to field-measured pH, reflecting its coarse categorization and low spatial resolution. Notably, soil pH was the only UKSO variable that showed a significant effect; other GIS-derived soil factors were either non-significant or could not be reliably utilized in models. Consequently, the present results do not support substituting UKSO soil classifications for direct field measurements, and more fine-scale assessment of soil characteristics is likely required for accurate habitat suitability evaluation. Crucially, soil and habitat heterogeneity on the farms surveyed is at a much finer scale than in currently available soil maps, and as such within-farm ground assessments are a useful method to predict likely foci for G. truncatula populations.
In future, such fine-scale assessments could potentially be achieved using drone- or satellite-based remote sensing, as multispectral and hyperspectral imagery have been shown to estimate soil pH (Webb et al., Reference Webb, Barnes, Powell and Jones2021), OM content (Sevilla et al., Reference Sevilla, Moreno, Ortega and Ballesteros2025) and aspects of soil texture (Song et al., Reference Song, Gao, Song, Li, Chen, Li, Zhang and Cai2023) with far higher spatial precision than coarse national datasets. Moreover, these sensing technologies may also facilitate high-resolution detection and monitoring of water-saturated habitats on farmland (De Roeck et al., Reference De Roeck, Van Coillie, De Wulf, Soenen, Charlier, Vercruysse, Hantson, Ducheyne and Hendrickx2014; Fernandez-Camacho et al., Reference Fernandez-camacho, Barja, Revilla, Ore, Alccacontor-muñoz, Morales, Tanabe, Carrasco-escobar and Cabada2026), providing a more comprehensive remote-assessment framework for modelling G. truncatula presence and supporting the development of next-generation, field-level prediction models of fluke risk.
Commonly, farmers are advised to ‘manage or avoid grazing wet areas’ to prevent liver fluke infection of their livestock, as this is thought to be the primary factor in determining the likelihood of G. truncatula presence and F. hepatica transmission. This study is further evidence that such assessments are much too simplistic, as the physicochemical properties of soil also play a significant role in determining the suitability of potential habitats. Currently, soil factors are often excluded from models of fluke infection. However, given their importance in determining the probability of the presence of the main intermediate snail host, soil factors should be included to make these models more precise. At this stage, the underlying mechanisms behind the impact of these soil factors are not completely understood, especially in the case of soil OM. Going forward, further research is needed to examine the interactions between soil OM, habitat vegetation, pH, soil particle size, algal communities and mineral availability (especially calcium) in the context of habitat suitability for G. truncatula. The complexities and current uncertainties of how these relationships impact G. truncatula also highlight the importance of multivariable models that account for non-linear associations. Once these underlying interaction effects are elucidated, the inclusion of these factors in models will, in turn, allow for the creation of more robust strategies to help prevent fluke infection.
Supplementary material
The supplementary material for this article can be found at https://doi.org/10.1017/S0031182026102820.
Data availability statement
The data used in this paper, along with the R files used to perform statistical analysis, can be found at this web address: https://doi.org/10.20391/4146ad2f-64df-4a32-8249-a3ba21f83ccd.
Acknowledgements
We would like to thank the many farmers who allowed us to collect data on their land and provided their time and effort to complete surveys and help guide us across their property. We would also like to thank other staff at the Department of Life Sciences of Aberystwyth University who provided occasional assistance.
Author contributions
This study was conceptualized (ideas; formulation or evolution of overarching research goals and aims) by C.S. and R.J., and the methodology was developed by C.S., R.J., W.S. and S.G. Investigation (conducting a research and investigation process, specifically performing the experiments, or data/evidence collection) was performed by C.S., R.J., S.E. and C.D., and formal analysis (application of statistical, mathematical, computational or other formal techniques to analyse or synthesize study data) was performed by C.S. and R.J. Writing was primarily performed by C.S., with editing and review performed by all authors. Supervision (oversight and leadership) was conducted by R.S. and E.R., while funding acquisition was obtained by R.J., E.M., R.M., G.R., P.B., M.W. and H.W.
Financial support
This research was funded through the UKRI FoodBioSystems DTP (Grant Ref: BB/T008776/1) and the Biotechnology and Biological Sciences Research Council Precision solutions for fasciolosis control in sheep (Grant Ref: BB/X01746X/1) projects.
Competing interests
N/A.
Ethical standards
Not applicable.
Open access
