Linking soils and human health: geospatial analysis of ground-sampled soil data in relation to community-level podoconiosis data in North West Cameroon.

BACKGROUND: Podoconiosis is a form of leg swelling, which arises when individuals are exposed over time to red clay soil formed from alkaline volcanic rock. The exact causal agent of the disease is unknown. This study investigates associations between podoconiosis disease data and ground-sampled soil data from North West Cameroon.
METHODS: The mineralogy and elemental concentrations were measured in the soil samples and the data were spatially interpolated. Mean soil values were calculated from a 3 km buffer region around the prevalence data points to perform statistical analysis. Analysis included Spearman's rho correlation, binary logistic regression and principal component analysis (PCA).
RESULTS: Six elements, barium, beryllium, potassium, rubidium, strontium and thallium, as well as two minerals, potassium feldspar and quartz, were identified as statistically related to podoconiosis. PCA did not show distinct separation between the spatial locations with or without recorded cases of podoconiosis, indicating that other factors such as shoe-wearing behaviour and genetics may significantly influence podoconiosis occurrence and prevalence in North West Cameroon.
CONCLUSION: Several soil variables were statistically significantly related to podoconiosis. To further the current study, future investigations will look at the inflammatory pathway response of cells after exposure to these variables.

Introduction

Podoconiosis is a non-infectious geochemical disease, affecting individuals exposed to red clay soil derived from alkaline volcanic rock.[1] The disease causes asymmetrical bilateral swelling of the lower legs. The exact causal agent in the soil is unknown, but has previously been linked to elements such as aluminium,[2,3] zirconium, beryllium[4] and silica,[5] as well as minerals such as quartz.[5] However, recent studies have identified that a genetic aspect of the disease exists that affects individual susceptibility.[6]

The majority of studies investigating associations between environmental variables and podoconiosis data were located in Ethiopia,[5,7-10] with a smaller number of studies taking place in Cameroon[11-13] and Rwanda.[14] A previous study of significance focused on podoconiosis prevalence in Cameroon, which successfully employed a large-scale approach, mapping the environmental influence on prevalence using country-wide digital soil maps at a spatial resolution of 250 m.[11] These maps, and specifically the clay and silt fraction content contained therein, were used to provide inputs for ensemble models to model the prevalence of podoconiosis. Results indicated that where the silt fraction exceeded 25%, the probability of podoconiosis increased; conversely, where the clay fraction exceeded 40% in the top soil, the probability of podoconiosis decreased.[11] Another notable example of a large-scale study used samples from an online database representing bedrock from five regions of Africa known to be associated with podoconiosis. The study also employed a sixth region with the Hawaiian islands as a control.[15] In the study, weight percentage data for oxides within the samples from across the six studied regions were analysed using a combination of principal component analysis (PCA), discriminant function analysis and analysis of variance, with results suggesting that a unique alkaline- and silicon-rich geochemistry underlies regions associated with podoconiosis.[15] Large-scale studies such as these have proved useful in identifying broad soil and geochemical type associations with podoconiosis over large areas; however, they cannot discriminate small-scale variations at local levels. Smaller scale podoconiosis studies have previously investigated the associations between detailed ground-sampled soil data and community-level podoconiosis data in northern Ethiopia[5] and in Kenya.[2] The Ethiopia study identified several statistically significant variables linked to prevalence, including quartz, mica and smectite content levels.[5] The Kenya study identified, through multivariate modelling, and while adjusting for the frequency of participants wearing shoes, that iron was significantly associated with podoconiosis.[2] Additionally, when the model was controlled for iron, aluminium concentrations became significant.[2]

Establishing consistent associations between environmental variables and podoconiosis has been problematic and performance throughout the literature has been variable, particularly when analysed at different spatial scales.[2,5,7-9,15] Previous large-scale podoconiosis studies have indicated that climate variables can be important predictors of podoconiosis[7,11]; however, in the current study, the focus will remain solely on high resolution soil data due to the limited climatic variability across the study area.

This novel study develops previous research, particularly for North West Cameroon, by aiming to increase both the spatial resolution of the investigation and detail of soil analysis. The study relates community-level podoconiosis data, collected by community health implementers (CHIs), with mineralogical and elemental data obtained via x-ray diffraction (XRD) and inductively coupled plasma mass spectrometry (ICP-MS), respectively.

Materials and Methods

Study area

This study was conducted in the North West Region of Cameroon (Figure 1) between latitudes of 5°53′ to 6°19′ and longitudes of 9°42′ to 10°18′. The North West Region of Cameroon is a podoconiosis-endemic region.[11] The main occupation in North West Cameroon is farming.

Figure 1.

Map of Africa, with two inset maps of Cameroon. The inset map shows the location of soil and community-level podoconiosis sample locations.

Data collection

Podoconiosis data collection

The disease data were collected as part of a previous study in which the collection strategy is reported.[16] The study consisted of individuals aged >18 y and those who had lived in the area for at least 10 y. The census was carried out by CHIs, who visited all the households in the selected communities, registering every individual in each household and screening those individuals for podoconiosis. Research team members re-examined the cases considered positive by CHIs. The data included global positioning system (GPS) coordinates of the community, podoconiosis prevalence data (the proportion of the community with the disease) and a binary variable representing presence or absence of the disease, which will be referred to as the occurrence data. A correction factor of 48.5 was applied to the prevalence data to account for overdiagnosis of podoconiosis during data collection by CHIs[16]; 168 sampling points were selected and utilised in this study due to their spatial contiguity with the soil sampling data points (Figure 2).

Figure 2.

Map highlighting the spatial contiguity and location of the soil and podoconiosis community-level sample locations in the North West Region of Cameroon.

Soil data collection

One hundred and fifty-two sampling sites were spaced 4.5 km apart in a gridded formation. Samples were taken at the centre of each grid square. At each soil sampling site, GPS coordinates (WGS 84, decimal degrees), elevation (m), vegetation type and fertiliser/insecticide usage were recorded. Any vegetation or rocks were removed from the surface and a rock hammer was then used to mix soil to a depth of 10 cm. Any large roots, rocks or stones >0.5 cm in diameter were removed and a trowel was used to collect one scoop (approximately 300 g) of soil, which was stored in a sample bag.

This process was then repeated twice at each sampling site, at 10 paces due west and 10 paces due east from the original sampling point. At each sampling location the completed samples were all mixed into the same sample bag.

In addition to this standard sampling approach, 10 random grid squares featured 5 extra sampling sites. This was implemented to capture the soil variability at a greater resolution; the format of this collection is shown in Figure 3.

Figure 3.

Sampling technique used for grids with five extra soil samples (a-e); distances are shown in m (not to scale).

The central blue circle represents the original soil sample and a–e represent the five additional sampling points.

However, some samples were not collected, due to inaccessibility of sample locality, or were deemed missing, resulting in 194 soil samples being collected out of the planned 202 samples (Figure 2).

Analysis of the chemical element constituents of the 194 soil samples was performed using ICP-MS. The total element content was obtained for 0.25 g subsamples of soil digested using a mixture of hydrofluoric, perchloric and nitric acids. Semiquantitative XRD was employed to derive the mineralogical content of the soil samples. A subsample of 100 of the 194 soil samples were tested using XRD. XRD and ICP-MS analysis was completed by the British Geological Survey. The soil variables measured can be found in the supplementary data Tables 1–3.

Data preparation

Spatial interpolation of soil variables

High resolution soil sampling over large areas is not always achievable, due to logistical, temporal and financial constraints. Spatial interpolation is therefore a vital tool that is frequently used for the mapping of soil properties over large spatial extents.[17] Although several authors have assessed and compared different interpolation methods for soil data modelling, no consensus exists on the best method.[17,18] Four spatial interpolation techniques commonly employed for soil mapping were examined in this study to determine the effectiveness of each to represent the spatial distribution of the soil variables. These methods included inverse distance weighting (IDW),[18] ordinary kriging (OK),[19] empirical Bayesian kriging (EBK)[20] and universal kriging (UnK).[21] These interpolation methods were carried out using ArcGIS version 10.5.1 (ESRI, Redlands, CA, USA). The ‘leave-one-out’ cross validation method was utilised to assess model performance. Prediction error statistics were used for model selection, which can be found in supplementary data table 4, and the workflow model followed can be found in supplementary data figure 1.

Buffer analysis of prevalence data

The podoconiosis prevalence data points, collected by CHIs, were geolocated at the centroid of each community, thereby representing a localised point source. To enable effective statistical analysis and to capture the lifestyles of participants more realistically, a larger region than a centroid position was considered. This provided more accurate models for residents who travel large distances for daily activities, for example, tending to their farms, which are often beyond village boundaries (S. Wanji, personal communication). Spatial buffer regions were created around disease data points to provide an estimation of the area residents interacted with during these activities. As there were no studies in Cameroon that investigated daily activity distances, a previous study looking at villages in Uganda was used as a precedent.[22] The Uganda study estimated that 3 km was representative of the average distance that inhabitants walked from their villages on a daily basis.[22] Although this study does not present the lifestyles of citizens in Cameroon and Uganda as being similar, they do share certain common traits, for example, in the rural areas of both countries villagers walk large distances to collect water and there is also a strong reliance on agriculture and livestock husbandry. Therefore, in the current study, 3 km buffer zones were created as a default around each prevalence data point, from which mean soil variable values could then be calculated from the interpolated soil property surfaces within the buffer areas.

Statistical analysis

Bivariate analysis

The associations between the soil and disease data were analysed using bivariate analysis.

Spearman's rho correlation analysis

Spearman's rho correlation analysis was carried out between mean soil variable values from the interpolation surface calculated using the 3 km buffer zones and the corresponding podoconiosis prevalence value. In total, 168 data points were used for the correlation analysis. Variables with p>0.05 and variables with a significant negative correlation were removed from the subsequent multivariate analysis.

Binary logistic regression analysis

Logistic regression analysis was carried out between mean soil variable values from the interpolation surface calculated using the 3 km buffer zones and the corresponding binary variable representing disease occurrence. Logistic regression moves beyond correlation analysis because, instead of using prevalence (which measures the amount of podoconiosis in relation to population), it utilises occurrence data. By analysing occurrence data, the soil variables can be examined to predict how a unit increase of the individual soil variable affects the likelihood of a community having at least one case of podoconiosis. In total, 168 data points were used for the analysis. Variables with p>0.05 and variables with OR<1 were removed from the subsequent multivariate analysis.

Multivariate analysis

PCA

PCA was utilised to reduce the dimensionality of the soil element and mineralogical data, while retaining the variation present in the original data set. This was completed by a linear transformation of the original data to create a set of new orthogonal variables called principal components.

The PCA included those soil variables identified from the bivariate analysis as statistically significant (p<0.05) with either a positive correlation or OR>1. To visualise the underlying structure of the data, a biplot was created and ellipses of 95% confidence were produced representing areas with podoconiosis (presence) and areas with no identified cases of podoconiosis (absence).

Multiple testing

Within this study, multiple testing corrections were applied due to the large number of individual tests performed; however, all (premultiple testing) significant results were included in the subsequent PCA. This approach was implemented as the aim of this study was not to state that the statistically significant variables are the ‘causal agents of podoconiosis’, but rather to identify potential variables, which could be investigated in subsequent analysis.

Results

Spatial interpolation outcome

For several variables it was not possible to perform spatial interpolation, due to their excess zero values and the inability to produce a reliable interpolation surface. Zero values occurred due to the absence of the variable or due to the sensitivity of the testing procedures employed. Variables affected by these issues included the minerals anatase, cristobalite, tridymite, amphibole, olivine, chromite, ilmenite, smectite, vermiculite and hydrobiotite.

Using prediction error statistics, interpolations associated with each technique were examined then chosen for further processing (see Figure 4 for an example interpolation surface). These included the following surfaces: 0 from IDW, 2 from UnK, 20 from OK and 42 from EBK. The prediction error statistics for these surfaces can be found in the supplementary data table 5.

Figure 4.

Empirical Bayesian kriging (EBK) interpolation surface for the element beryllium (Be). Prediction error statistics: root mean square error (RMSE)=0.927, mean rrror (ME)=−0.004, average standard error (ASE)=0.927, root mean square standardized error (RMSSE)=1.002.

Spearman's rho correlation outcome

From correlation analysis of the 3 km buffer zones around the disease data coordinates, it was identified that the mean soil values for the chemical element variables barium (Ba), beryllium (Be), potassium (K), rubidium (Rb), strontium (Sr) and thallium (Tl), as well as those of the mineral variable potassium feldspar (K-feldspar), were positively correlated with the prevalence of podoconiosis (Table 1).

Table 1.

Correlation between the mean soil variable values extracted from the 3 km buffer zones and the corresponding podoconiosis prevalence data. Spearman's rho correlation coefficient, 95% CIs and p-value. Sample size 168. Only variables with p-value<0.05 are included in the table. Variables are shown by their abbreviated names.

Variable	Correlation coefficient	95% CIs	p
Elements
As	−0.200	−0.343 to −0.049	0.009
Ba*	0.230	0.079 to 0.370	0.003
Be*	0.194	0.043 to 0.337	0.012
Cd	−0.202	−0.345 to −0.051	0.008
Ce	−0.191	−0.334 to −0.04	0.013
Cr	−0.205	−0.347 to −0.053	0.008
Fe	−0.154	−0.300 to −0.002	0.046
Gd	−0.176	−0.320 to −0.024	0.022
Hf	−0.157	−0.302 to −0.005	0.042
K*	0.182	0.03 to 0.326	0.018
LOI[a]	−0.212	−0.354 to −0.061	0.006
Lu	−0.154	−0.299 to −0.001	0.047
Mo	−0.162	−0.307 to −0.01	0.036
Nb	−0.176	−0.320 to −0.024	0.022
Nd	−0.173	−0.317 to −0.021	0.025
P	−0.266	−0.403 to −0.117	0.0001
Pr	−0.170	−0.314 to −0.018	0.028
Rb*	0.161	0.009 to 0.306	0.037
Sm	−0.173	−0.307 to −0.021	0.025
Sr*	0.155	0.003 to 0.301	0.044
Ta	−0.193	−0.336 to −0.042	0.012
Ti	−0.202	−0.344 to −0.05	0.009
Tl*	0.192	0.041 to 0.335	0.013
Tm	−0.170	−0.314 to −0.018	0.028
V	−0.252	−0.391 to −0.103	0.001
Zr	−0.160	−0.305 to −0.008	0.038
Minerals
Goethite	−0.152	−0.297 to 0.000	0.049
K-feldspar*	0.154	0.002 to 0.300	0.046

LOI=loss on ignition represents organic matter content.

*Variables with a statistically significant positive correlation coefficient.

Binary logistic regression outcome

From binary logistic regression, the mean soil values identified from the 3 km buffer zones around the disease data coordinates for both Ba and quartz had a significant OR>1 (Table 2). This suggests that with a unit increase of these variables, the likelihood that an area has at least one case of podoconiosis increases.

Table 2.

Binary logistic regression between the mean soil values extracted from the 3 km buffer zones and the corresponding binary occurrence data of podoconiosis. The coefficient, standard error, z-value, p-value, OR and 95% CIs. Sample size 168. Only variables with p-value <0.05 are included in the table. Variables are shown by their abbreviated names.

Variables	Coefficient	Std. error	Z-value	p	OR	95% CI
Elements
Ba*	0.002	0.001	2.28	0.017	1.002	1.000 to 1.003
LOI[a]	−0.184	0.062	−2.98	0.002	0.832	0.737 to 0.939
P	−0.002	0.001	−2.28	0.02	0.999	0.997 to 1.000
Th	−0.027	0.013	−2.05	0.039	0.973	0.948 to 0.999
V	−0.016	0.006	−2.5	0.011	0.985	0.973 to 0.997
Minerals
Amorphous	−0.054	0.026	−2.06	0.037	0.947	0.900 to 0.998
Quartz*	0.091	0.04	2.29	0.018	1.095	1.013 to 1.184

LOI=loss on ignition represents organic matter content.

Variables with a statistically significant OR > 1.

PCA outcome

PCA was utilised to identify any underlying patterns in the data by creating a biplot of the scores and variable loadings for the principal components PC1 and PC2 (Figure 5).

Figure 5.

Biplot of principal component scores and factor loading on principal component 1 and 2 from the PCA analysis of soil variables, extracted from the 3 km buffer zones, which have been suggested to be positively related with podoconiosis through Spearman's rho correlation and/or binary logistic regression. PC1 eigenvalue=5.17. PC2 eigenvalue=1.37. Ellipses are created with a confidence level of 0.95.

The first two PCs account for 64.6% and 17.1% of the total variation of the data set, respectively. PC3 accounted for 8.6% (eigenvalue=0.69) of the variation of the data. Therefore, PC3 and subsequent PCs, with decreasing eigenvalues, were not investigated in this study as they explained low variation of the data set and would be unlikely to identify any patterns.

In the PCA ordination, the communities with at least one case of podoconiosis are represented by blue triangles and those communities with no recorded cases are represented by red dots.

The PC1 axis highlights a clustering of the positive eigenvectors, suggesting they are highly correlated with each other. The variables with the greatest eigenvectors on the horizontal axis and therefore with the greatest weighting on PC1 include the variables K-feldspar, Sr, K, Rb, Ba and Tl. The variables with the greatest eigenvectors on the vertical axis, and therefore with the greatest weighting on PC2, include the variables quartz and Be.

From a visual inspection of the PCA, it appears that there is no distinct separation between the communities with at least one case of podoconiosis and communities with no recorded cases. However, with the 95% confidence ellipses added, the data do show a separation, on both component axes, between the disease presence ellipses as the selected soil variables increase in value.

Discussion

Spatial interpolation

Continuous ground-sampled surface soil data were not available, and therefore, to co-locate and extract meaningful soil values for each prevalence data point, interpolation modelling was necessary using the existing soil sample data set.

In this study, and following independent validation, EBK was identified as the most suitable method for 42 of the 64 soil variables. However, it should be noted that regarding the prediction error statistics, OK performed almost as well in most cases. This finding is consistent with other studies.[23] The suggested advantage of EBK over the two other kriging methods (UnK and OK) is in its ability to account for the error introduced by estimating the semivariogram model through repeated simulations.

Buffered spatial selection

Previous studies comparing environmental and podoconiosis data have not considered the movement of people directly.[2,5] Instead, these studies either compared the soil variables beneath the disease data coordinates[5] or compared soil variables collected from the households in the study.[2] These studies are ultimately reliant upon geolocation accuracy, village boundaries and sedentary lifestyles.[2,5] However, the buffer approach implemented in this study considered the geolocation accuracy of the prevalence data points and the movement of people in this area of Cameroon within a fixed radius.

Statistical analysis outcome

The PCA and biplot revealed no discrete groupings of soil compositions with respect to communities in which podoconiosis was present or absent. This suggests that the relationship between soil and this disease is likely complex and nuanced, moreover, areas can have similar soil compositions but different occurrence values. However, the results do indicate separation between the disease occurrence ellipses as the selected soil variables increase in value and that these variables represent disease covariates. This further supports the literature and may indicate that other factors such as shoe-wearing behaviour, foot washing[24] and genetics[6] may play a significant role in podoconiosis occurrence. The genetic and behavioural influence of the disease may also indicate why correlation values for the significant variables were moderate to weak. Although statistical analysis was robust, the weak to moderate outcomes could also be attributed to inaccuracies in the CHIs’ collection of prevalence data[16] and the common occurrence of initial misdiagnosis of podoconiosis. Misdiagnosis may have occurred due to the occurrence of lymphatic filariasis in North West Cameroon[25] as oedema of the foot and lower leg is present in both diseases. During the training of CHIs, measures to combat the misdiagnosis of podoconiosis were implemented by explaining the physical differences between the two diseases.[16] Misdiagnosis is acknowledged as a potential limitation within this and related studies, and may require further interpretation of the required correction factor to account for errors.

Several of the variables identified as being significantly related to podoconiosis, namely, Ba,[9] Be,[4] K,[3,15] Sr,[4] quartz[5] and K-feldspar,[9] have previously been identified in the literature in relation to podoconiosis. Multiple studies in Ethiopia have identified Ba as less common in the soils of endemic regions compared with non-endemic regions,[4,9] contradicting the findings in the current study. Similarly, Sr was reported at higher levels in soils from towns where podoconiosis did not occur.[4] However, it has been suggested as probable that soil variables identified as significantly related to podoconiosis in other countries (e.g. Ethiopia) may not equally influence the distribution of podoconiosis everywhere else.[11] Conversely, soil-beryllium in Ethiopia was measured at relatively high concentrations in areas with podoconiosis compared with neighbouring areas where podoconiosis prevalence was low.[4] Be has been suggested to induce granuloma formation in the lymphoid tissue of humans.[26] The relationship with Be is intriguing. Be is relatively rare in rock-forming minerals (e.g. beryl and bertrandite) but may also be associated with feldspars mostly concentrated in granitic pegmatites.[27] Soils developed from these parent rocks may therefore contain Be, frequently bound to clay minerals or organic matter, and may be augmented by anthropogenic coal ash, ore-processing products or sewage sludge.[28] Be is one of the most toxic elements in the periodic table and, if inhaled or ingested, can lead to the frequently fatal lung disease, chronic beryllium disease (chronic berylliosis), and is listed as a class A Environmental Protection Agency (EPA) carcinogen.[28]

K has previously been identified as occurring in mineral particles of elephantiasis nodes, but Price and Henderson were unable to report if there was a significant difference between non-elephantiasis nodes.[3] Other studies reported an increase of K in podoconiosis-endemic soils.[9,15] However, it is inferred that K is an unlikely pathological variable at small-scale absorption as K is an essential constituent of the body.

In the Ethiopia study, quartz was significantly associated with the prevalence of podoconiosis in ground-sampled soil variables.[5] Quartz is a very common mineral, occurring in almost all rock types and, due to its resistance to weathering, in most soils. Although inhalation of very fine-grained (<10 µm) forms of quartz may cause silicosis,[29] the relationship between quartz and podoconiosis is perhaps surprising considering quartz's generally inert chemical behaviour. However, its relative hardness (7 on the Mohs scale) may result in increased skin abrasion, facilitating the entry of more toxic elements. As previously suggested, repeated exposure to quartz could worsen stratum corneum degradation and enable particles to penetrate the skin.[9]

In Ethiopia, feldspar was reported to have a statistically significant lower proportion in endemic compared with non-endemic soils.[9] However, it should be recognised that the feldspar identified in the Ethiopian study[9] could represent multiple types of feldspar other than the K-feldspar measured in the current study. K-feldspar is a common mineralogical component of igneous and metamorphic rocks in particular, but it is also found in sedimentary rocks and therefore soils derived from all three bedrock types. The association between Ba, K, Rb, Sr, Tl and K-feldspar is well established, as K forms a major component (11–14%) of this group of minerals, and Ba, Rb, Sr and Tl all readily substitute for K, as they have similar ionic radii.[30]

Finally, Rb and Tl have not been previously identified as potentially associated with podoconiosis. To the knowledge of the authors, no previous studies have examined these variables in relation to podoconiosis.

From this study it is further highlighted that elements and minerals identified as potential disease covariates or triggers of podoconiosis can vary from one place to another. It should be acknowledged that factors describing the soil physical properties, such as the texture and structure of the soil, or complex interactions between elements and minerals, could have a greater impact on the risk of podoconiosis than a single element or mineral. This potentially complex interaction could explain why there is no consensus on single elements or minerals from these studies across multiple countries. It is recommended that future studies should also investigate the physical properties of the soil in relation to podoconiosis.

Conclusions

From the analysis presented in this study, several soil variables were statistically identified as significantly related to podoconiosis. It is, however, not clear if the variables identified can be suggested as disease covariates or as causal agents of the disease. The weak to moderate strength of these relations could be impacted by other influences such as genetic variation and shoe-wearing behaviour, as well as the errors inferred from spatial interpolation. To further the current study, future investigations should consider the inflammatory pathway response of cells to exposure of these significant variables. This would further our understanding as to whether these soil variables have a potential role in the pathogenesis of podoconiosis.

Click here for additional data file.