Danyang-gun is a county in the northeast region of Chungcheongbuk-do Province, South Korea, with a population of approximately 29,970. It is located in the range of 128°13′ to 128°39′E and 36°47′ to 37°09′N, and has an area of 780.67 km² (Figure 1). It is well-known for its scenic surroundings, including the Sobaek Mountain range and Namhan River. Sobaek Mountain is the second highest mountain in South Korea (elevation = 1,439 m) and the Namhan River flows for 23.7 km from northeast to southwest along the Sobaek Mountains. Only 11.2% of the county is cultivable, and 83.7% is mountainous. Because of this rugged terrain, both settlements and urban areas are developing in the hills and valleys. The annual mean precipitation is 1,113 mm and the annual mean temperature is 11.5°C; the highest and lowest temperatures are 17.5 and 6.6°C, respectively (KMA, 2021).

Danyang-gun is composed of various lithological units and strata, as well as complex and diverse geological structures. It consists of Precambrian base rock, Paleozoic sedimentary rock, Mesozoic sedimentary rock, and igneous rock (Figure 2). The Precambrian rock is located in the eastern study area and coincides with Sobaek Mountain. This rock has undergone granitization after regional metamorphism; it is divided into granitic and migmatitic gneisses (Won and Lee, 1967). The sedimentary rock includes unknown age quartzite, Paleozoic clastic sedimentary rock, and carbonate rock. The quartzite covers carbonate rock on the northwest side with an unconformity and is located at the western end of the study area. However, the sequence of formation is unclear because there is no direct contact with other formations (Won and Lee, 1967). The clastic sedimentary rock is composed of Cambrian quartzite and slate; it generally shows a strike of N30°E or N45°E. The carbonate rock is Cambrian–Ordovician and courses in the NE and NW directions (Aum et al., 2019). The NE carbonate rock consists of limestone, dolomitic limestone, dolomite, and banded limestone. The Mesozoic sedimentary rock covers this carbonate rock with a clinounconformity. This formation is mostly composed of clastic sedimentary rock such as shale, sandstone, and conglomerate; layers containing anthracite have also been identified.

A fault exists in the northern part of the most recent Mesozoic formation, and carbonate rock from the NW direction is distributed to the north of the study area according to this fault. Most carbonate rock from the NW direction is composed of limestone and dolomite; several types of clastic sedimentary rock of unknown age are also present. Mesozoic rock is divided into sedimentary and igneous rock. The sedimentary rock is distributed in the NE direction, as described above. The igneous rock intruded in the Cretaceous period; it includes biotite granite, quartz porphyry, and granite porphyry. The biotite granite, which is widely distributed in the south, is in contact with sedimentary rock; this forms a contact metamorphic zone. There are faults in the NE and NW directions in the study area. The faults in the NE direction cross the center of the study area, and the geology on both sides is clearly distinguished by these faults. The faults in the NW direction cut the sedimentary formations with an NE strike in an almost perpendicular direction (Won and Lee, 1967).

Various geological, geochemical, and topographical factors are associated with IRC values. Following a literature review and assessment of the available data, as well as application of the FR method, 10 effective factors were identified for IRC modeling (Table 1). These factors included lithology; fault density; mean soil calcium oxide (CaO), copper (Cu), lead (Pb), and ferric oxide (Fe₂O₃) concentrations; elevation; slope; valley depth; and TWI. Importantly, the FR values reflect probabilistic spatial relationships of dependent variables (IRC values, obtained from field measurements) with the various classes of each independent variable (“radon factors”). The FR values can be calculated as follows: F R = N r / T r N p / T p (13) where N r is the number of training samples in each subclass of IRC effective factors, T r denotes the total number of training samples, N p is the is the number of pixels of each sub-class of the effective factor, and T p indicates the total pixels of the study area. An FR value >1 indicates a high correlation between radon level and a particular factor, an FR value <1 indicates a low correlation, and an FR value of 1 indicates a moderate correlation (Al-Abadi et al., 2016).

To identify relationships among the included effective factors, multicollinearity analysis was performed based on the variance inflation factor ( V I F ) and tolerance ( T O L ) (Arabameri et al., 2021c). Importantly, some factors were found to exert a negative influence on the predictive capacity of the model. Such variables were removed from the model to increase its prediction accuracy (Miraki et al., 2019). The relative importance and predictive abilities of the various factors were determined using the information gain ratio ( I G R ). This is an entropy-based method that only considers variables associated with occurrence of an event (Bui et al., 2018). A higher IGR value indicates that a factor has greater impact on the model predictions (Panahi et al., 2021).

Radon levels and their controlling factors vary spatially, and the selection of appropriate predictive variables is important for radon mapping accuracy. As shown in Table 1, 10 geogenic effective factors were used to model the indoor radon level. Local lithology and fault density are crucial factors affecting radon production and distribution, even in adjacent areas (Buttafuoco et al., 2010). Radon is released naturally via uranium-bearing mineral decay, such that fractures and faults provide an important route for radon migration from bedrock to the surface (Ciotoli et al., 2017). For the current study, geological and geochemical maps from the Korea Institute of Geoscience and Mineral Resources (https://www.kigam.re.kr/) were used (Figure 3). In addition to site geological characteristics, the concentrations of some chemical elements (i.e., CaO, Cu, Pb, and Fe₂O₃) remaining in minerals and soil after erosion can affect IRCs. Soil geochemistry can serve as a predictor of radon level (Ball et al., 1991; Schumann and Gundersen, 1997; Drolet et al., 2014). The effects of bedrock geochemistry on IRC are reportedly greater than those of topsoil properties, because a large portion of the topsoil tends to be removed during construction; thus, only a few centimeters remain (Appleton, 2013).

In addition to geological variables, topographical factors were considered for our indoor radon potential mapping. The data were derived from a digital elevation model with a resolution of 10 m, provided by the National Geographic Information Institute (http://www.ngii.go.kr). The data were processed by SAGA software (http://www.saga-gis.org/en/index.html) to produce slope, valley depth, and TWI maps (Figure 3). In the present study, the TWI was used as a proxy of the spatial distribution of soil moisture, and was calculated as follows: T W I = ln β tan ⁡ α (14) where β and α are the cumulative catchment area in m² and slope angle in radians, respectively (Arabameri et al., 2021b). The TWI can reflect the water transmissivity and infiltration rate at a given location. Areas with low slope angles have high TWI values, while steeper areas have low TWI values (Mattivi et al., 2019). Notably, pores saturated with water trap radon in the soil and slow its transport through soil into the atmosphere (Kellenbenz and Shakya, 2021; Shahrokhi and Kovacs, 2021). However, soil moisture content can influence radon escape from mineral matter if fewer than 30% of the soil pores are filled with water; higher soil moisture content leads to a considerable reduction in radon emanation because of decreased gas permeability (Je et al., 1999). Furthermore, large valley depth values indicate areas with low elevation and gentle slopes (Figure 3). In such areas, the infiltration rate is high; the high soil wetness and fine texture lead to low permeability, in turn causing convective radon flow and slow soil gas exhalation (Wiegand, 2001).

The generation of a radon inventory map is important for developing a machine learning-based model. In the current study, with the aim of obtaining representative samples of indoor radon levels, 1,452 dwellings were selected at random throughout the study area. Since 2008, passive IRC measurements have been conducted by National Institute of Environmental Research (NIER) using alpha-track detectors (Raduet; Radosys Ltd., Budapest, Hungary). The detectors were typically positioned in the living room, where residents spent most of their time. Each measurement period (all in winter) was 90 days in duration; the collected data were returned to NIER for analysis, and showed that the IRC value exceeded the recommended level of 148 Bqm⁻³ in 726 samples. To develop the model, samples were classified in a binary manner in terms of their IRC values. Samples with IRC >148 Bqm⁻³ were coded as 1, indicating locations with high radon levels. All remaining samples with low IRC values were coded as 0, indicating locations with low radon levels. Two classes of data (high and low radon levels) with equal numbers of samples (726) were randomly divided into training and testing subsets at the ratio of 70:30 (Kadirhodjaev et al., 2020; Panahi et al., 2021; Roy et al., 2021). The distribution of the training and testing samples is illustrated in Figure 1. To build the model, the training dataset was constructed by combining 508 samples belonging to the high and low radon level classes. Similarly, to validate the predictive accuracy of the model, the testing dataset was constructed from 218 samples that belonged to the high and low radon level classes. The training and testing datasets were then superimposed with all of the radon factors to extract their attribute characteristics. Finally, the data were transferred into MATLAB software (https://www.mathworks.com) to construct the LSTM, ELM, and RVFL models.

Model validation is critical to confirm the reliability of machine learning algorithms. Various statistical analysis methods are used to evaluate modeling accuracy. The area under the receiver operating characteristic curve (AUROC) is a useful quantitative parameter, where accurately detected events are plotted on the y-axis (i.e., sensitivity) against false detections on the x-axis (i.e., 1–specificity). AUROC can be constructed from both training and testing datasets, to yield success and prediction rates, respectively. The success rate curve represents model accuracy according to the locations of the samples; the prediction rate curve indicates the predictive power or generalizability of the model (Golkarian and Rahmati, 2018). The AUROC takes a value between 0 and 1, where values closer to 1 reflect better predictive ability (Park et al., 2017).The root mean square error (RMSE) and standard deviation (StD) are another statistical approaches used to assess the prediction accuracy of a model with n total variables, as follows: R M S E = 1 n ∑ i = 1 n ( X p r e d i c t e d − X a c t u a l ) 2 (15) S t D = 1 n ∑ i = 1 n ( X p r e d i c t e d − X ¯ p r e d i c t e d ) 2 (16) where X ¯ p r e d i c t e d is the mean value of the predicted dataset, X p r e d i c t e d and X a c t u a l indicate the predicted and actual values of the variable, respectively. An overview of the methods used for indoor radon potential mapping is provided in Figure 4.

Collinearity among effective radon factors was determined by calculating the V I F and T O L , where V I F > 10 and T O L < 0.1 indicate collinearity among predictors (Arabameri et al., 2021a). As shown in Table 2, the V I F and T O L values of the selected factors were lower than the critical values; thus, there was no collinearity among inputs. Notably, elevation had the lowest T O L (0.316) and highest V I F (3.160).

The I G R method was applied to rank the predictive capabilities of the variables; the results indicated that elevation had the strongest effect on radon-prone area mapping ( I G R = 0.61), followed by lithology (0.32), valley depth (0.30), and mean soil Cu concentration (0.29). In fact, all factors with I G R > 0 had predictive power (Table 2).

Various geogenic factors can affect radon levels; this can be quantified through statistical modeling, such as the FR method. Stronger correlations are indicated by higher FR values, while F R < 1 indicates a weak relationship between a given predictor and the IRC value. As shown in Table 3, 10 variables were used to predict areas with potentially dangerous radon levels. The results implied that an elevation of 120–242 m (FR = 3.29) had the greatest influence on the IRC. Importantly, approximately 80% of Danyang-gun is mountainous, with a shallow soil profile that mostly contains coarse fragments. This promotes soil permeability and movement of radon gas within the soil (Hauri et al., 2012). In contrast, soil gas accumulation in lowlands is high; thus, it can easily infiltrate indoor environments via the soil through openings and cracks in basement foundations. Analysis of the FR values for the slope factor showed that the highest value (2.41) was associated with the class of 0–13.5°. This finding implies that the IRC decreases in sloped areas because the released radon is rapidly diluted in outdoor air (Appleton, 2013). In terms of valley depth and TWI, the highest FR values were found in the sixth (161–370 m) and fifth (5.59–7.32) classes, respectively. These factors reflect the effects of hydrological variables (e.g., rainfall-runoff and infiltration rates) and soil moisture content on soil gas exhalation capacity; this capacity is generally diminished when soil wetness is increased (Sasaki et al., 2004; Raduła et al., 2018). Analysis of the relationship between the radon emanation rate and presence of specific uranium- and radium-containing minerals (i.e., in the host rock and remaining soils after weathering) showed that the highest values of FR were 1.47, 1.59, 1.76, and 1.48 for mean soil CaO, Cu, Pb, and Fe₂O₃ concentrations, respectively. For all factors, the FR values were >1, indicating strong correlations with radon levels at the monitoring sites. Furthermore, radon levels were high in areas where the fault density varied between 0.58 and 1.2 (FR = 2.32). Notably, fault systems located in fracture zones provide a route for radon to migrate upward from deeper sources (Han et al., 2006). Finally, in terms of lithology, the FR analysis yielded higher values, of 200.65, 32.74, and 13.71, for Cretaceous acidic dike (Kad), Cretaceous quartz porphyry (Kqp), and Cambrian quartzite and slate (CEdy) units, respectively. Generally, sedimentary, igneous, and metamorphic rocks contain variable amounts of uranium and radium, depending on the rock formation processes (Pasculli et al., 2014).

The maps generated using the LSTM, ELM, and RVFL algorithms are shown in Figure 5. The maps included five classes of radon-prone areas (very low, low, moderate, high, and very high), based on the quantile method (Khosravi et al., 2018). The percentage area of each class on each map is shown in Figure 6. The ELM model was the most accurate; it categorized 19.62, 20.64, 19.84, 20.01, and 19.88% of the study area into the very low, low, moderate, high, and very high classes. As depicted in Figure 5, high radon levels were observed in central and southwestern parts of the study area due to the distribution of sedimentary rock and unconsolidated deposits such as carbonate, shale, sandstone, conglomerate, limestone, and dolomite, all of which are rich in uranium and organic materials. These findings were consistent with the results of previous studies (Cho et al., 2015; Hwang et al., 2017; Kim and Ha, 2018; Park et al., 2019).

The reliability of the results was checked using the FR method, which revealed that most of the samples with high radon levels were from the very high and high radon potential areas. Thus, the models exhibited satisfactory performance in terms of study area classification. The AUROC values were calculated to quantitatively evaluate the predictive accuracy of each model. The AUROC values for the success rate curve analysis of the LSTM, ELM, and RVFL models were 0.81, 0.83, and 0.82, respectively. The AUROC value for the prediction rate curve analysis was 0.82 for the ELM model; the LSTM and RVFL models had lower values of 0.80 and 0.78, respectively (Figure 7). The RMSE values exhibited a similar pattern. As shown in Figure 8, analysis based on training data showed that the RMSE was lowest for the ELM model (0.152); the LSTM and RVFL models exhibited higher RMSEs of 0.163 and 0.182, respectively. Further analysis based on the testing data showed that the RMSEs of the ELM, LSTM, and RVFL models were 0.209, 0.232, and 0.0286, respectively. The standard deviation (StD) values for the ELM model (0.152 and 0.207) were lower than those for the LSTM and RVFL models, during both the training and validation phases. In summary, by comparison of the AUROC, RMSE, and StD values calculated using the training and testing datasets, all of the evaluated models had acceptable performance in terms of classifying radon-prone areas; however, the ELM model was slightly superior to the two other models.