Stream Chemistry and Forest Recovery Assessment and Prediction Modeling in Coal Mine Affected Watersheds in Kentucky, USA

1. Introduction

Water is one of the most essential natural resources on the earth but due to human activities including land use and climate change impacts the water quantity and quality and hence creates the water scarcity problems and threatens aquatic biodiversity (Vörösmarty et. al., 2010). Coal mining is one of the dominant factors to degrade water quality in the Appalachia region where hundreds of headwater streams impacted (Palmer et al., 2010; Holzman, 2011; Palmer and Hondula, 2014; Evans et al., 2021). Commercial coal production started in Kentucky even before Kentucky was recognized as a state. In 1790, 18 tons (20 short tons) of coal was produced in the first commercial coal operation in what is now Lee County. The demand for more energy led to an increase in coal production so that peak levels amount exceeded 162 million tons in 1990 [1]. Kentucky was the fifth-largest coal-producing state in 2016 [2] with 42 million tons of coal production (EEC, 2018). Seventy-nine percent of Kentucky's net electricity generation was coal-fired in 2017 [2] although production has declined consistently since 1990.

Even though coal-generated electricity is less expensive than other energy sources, its impact on the environment and humans is under debate [3,4]. The Surface Mining Control and Reclamation Act (SMCRA) of 1977 created two programs: one for regulating active coal mines and a second for reclaiming abandoned mine lands [5]. Many studies, however, have shown that reclamation practices do not produce very effective results [3,4,6,7].

Major coal mining impacts generally observed in Kentucky and Central Appalachia are natural land cover loss, hydrological pattern changes, valley fill, acid drainage, and water quality degradation (Figure 1; [3]). Mountaintop coal mining, which is widely used in the Appalachian region, involves blasting rocks, clearing forests, and removing soil to reach coal reserves [8,9]. For example, the landscape was greatly modified by smelting and mining activities in Copperhill, Tennessee which is situated in the Blue Ridge Mountains at the convergence of Northern Georgia, Western North Carolina, and Southern Tennessee which was began in 1854 and landscape appeared desert-like during 1900s (Baes and McLaughlin, 1984; Newman 2020). Central Appalachia has the highest earth movement rate in the United States with each surface mine generating large quantities of spoil that are typically translocated to stream valleys close to mining areas [10]. Mountaintop removal coal mining requires forest clearing, which in some cases may be visible on the 1: 3,000,000 scale from National Agriculture Imagery Program (NAIP) aerial Images. Mountaintop mining poses a potential threat to intact characteristics of Appalachian forests [11]. Interior forest loss (change in interior forest to forest edge) in southern Appalachia between 1991 and 2002 was 1.75 to five times higher than direct forest loss (loss from edges), which may have been the result of mountaintop mining [12]. Land cover alterations from active mine to reclaimed mines in central Appalachia almost entirely occurred in forested areas [3].

The Forest Reclamation Approach (FRA) provides a set of guidelines that are based on research conducted over several decades to promote the regrowth of forests on reclaimed mine lands which requires using loose soil [7]. This could potentially lead to elevated concentration of ions and metals; moreover, the effect of FRA on water quality needs further investigation [9]. The potential for total dissolved solids-source control practices that incorporate FRA may improve mine-water quality [7]. The effect of FRA on water quality may require additional data on the spatial extent of FRA locations and reclamation age.

Reclaimed forests and areas can be detected with remote sensing techniques [13,14,15] and assessed in relation to water quality. [16] found that, after seven years of monitoring, water quality improvement was more obvious in sub-watersheds that were heavily affected by past mining activities and reclaimed by reforestation than in lands with abandoned mines. These findings demonstrate that good reclamation practices can have a positive influence on water quality over time.

Besides forest loss, surface coal mining may have a significant effect on soil hydrological properties [3]. Surface mining may have negative effects on watershed hydrologic characteristics [6,17,18]. Valley fills can bury headwater streams under tens to hundreds of meters of spoil [8,9], which minimizes the drainage capacity of a valley (Figure 2).

[19] researched hydrologic differences in a reforested area and a grass covered area that were reclaimed 14 years prior to the study with similar processes. The reforested area had a higher infiltration rate and the grassy area had more surface flow paths; however, the water quality in these two areas was not documented.

Figure 3 displays how surface mining affects the natural land cover and the appearance of land cover during recovery.

The Clean Water Act Amendment of 1972 (CWA) established standards for regulating pollutant discharges into surface water [20]. However, groundwater and surface water from mined areas may be contaminated with numerous toxic solutes such as sulfate (SO₄^2-), iron (Fe), aluminum (Al), and selenium (Se) [21,22]. Increased conductivity, pH changes, and elevated dissolved ion levels are commonly found in streams near surface mines [4]. A relationship exists between water quality degradation and coal mining at various scale of watersheds, e.g., elevated SO₄ ^2-, alkalinity, conductivity, Ca²⁺, Mg²⁺, Mn²⁺, Al³⁺, and Fe^2+,3+ [4,7,23,24,25,26].

[27] found that conductivity fluctuated seasonally and was highest in summer and lowest in spring, which may be due to the dilution effect of water in the watershed. Similarly, in consecutive seasons, conductivity was two or more times higher in mined watersheds than unmined watersheds [23, 27, Table 1). These findings suggest that precipitation or seasonal conditions may affect the results; however, differences between unmined and mined areas are stable if water samples are collected under the same conditions, e.g., at the same season and a clear day without rain.

Coal mining typically releases pyrite (Fe₂S), which forms in association with coal [28,29,30,31,32]. When pyrite comes in contact with water and oxygen, it is oxidized by autotrophic bacteria, leading to acid mine drainage (Equation 1; [33]):

FeS₂ + 3.75O₂ + 3.5H₂O → Fe(OH)₃ + 2SO₄^2- + 4H⁺

Eq. 1

Carbonate minerals, e.g., calcite (CaCO₃) and dolomite (CaMg(CO₃)₂), can neutralize the acidity (Equations 2 and 3; [34]).

2CaCO_{3 +} 2H⁺ → Ca²⁺ + 2HCO₃^-

Eq. 2

CaMg(CO₃)_{2 +} 2H⁺ → Ca²⁺ + Mg^{2 +} + 2(HCO₃)^-

Eq. 3

Using the equations above, in addition to previous studies, SO₄^2-, alkalinity, conductivity, Ca²⁺, Mg²⁺, Mn²⁺, Al³⁺, and Fe^2+,3+ values in streams near coal mines can be described as coal mine-related stream chemistry (CMRSC).

Topographic variables may play role in stream chemistry even though they have been used rarely by coal mine-water quality researchers. However, some land cover change (urban, agriculture) studies utilized topographic variables as contributing factors to stream chemistry [35,36]. [36] found mean watershed elevation and slope significantly correlated with conductivity. Contrarily, [35] did not find any significant relationship among watershed slope, drainage density, and conductivity. [25] did not find any relation between elevation and stream chemistry in coal mine sites.

Various studies have documented a negative correlation between coal mining and the health of local ecology, stream habitat, community structure, and ecosystem functions [37,38,39]. Prior studies have focused on water quality degradation [4,7,23,25] or land change due to coal mining [3,12] whereas few studies have focused on CMRSC on mined watersheds along with land change [26,39,40], specifically vegetation cover change in eastern Kentucky. Among those studies that have been conducted in this area, [26] found a high correlation between mining percentage and SO₄^-2 and conductivity while [40] created a general linear model between land use (surface mining, residential development, and underground mining) and specific conductivity. Future comprehensive research should focus on the interrelationship among land, vegetation, and water quality changes over the years after disturbances caused by active mining and reclamation [26]. Moreover, data are needed to determine how the spatial and temporal extent of surface coal mining affects the watershed based on the vegetation cover change, reclamation age, and topographic factors. The aim of the current research was to examine how coal mining affects water quality and temporal vegetation cover change between 1986-2017 in mined watersheds on a regional scale.

2. Materials and Methods

2.1. Study site

This research was conducted in eastern Kentucky, U.S.A. (Figure 4. Eastern U.S. coal reserves and research area). The research area covers Johnson Creek, Troublesome Creek, and Quicksand Hydrologic Unit Code 10 (HUC10) watersheds in Magoffin, Knott, Perry, and Breathitt Counties of Kentucky (Figure 5). The research area was approximately 1,768 sq. km, including 58 stream reach watersheds, which were sub-watersheds of HUC10 watersheds (Figure 5). Stream reach watersheds ranged from 2.65 km² to 16.94 km². This study defined unmined watersheds as those that were less than 5% mined [40,41]. Similarly, since fertilizer usage would have affected conductivity, watersheds that did not have combined developed and agricultural lands that exceeded 5% were selected in this study to control the agriculture and urban effects on CMRSC. The research area contained various watersheds that were 2.5% to 90% mined and had a reclamation age range of 4 years to 29 years.

2.1. Data Collection, Preparation and Analysis

2.1.1. Water Samples Data Collection and Analysis Method

National Hydrography Dataset (NHD) 100K stream shapefile (USGS 2007-2014), HUC10 watershed polygons (USGS 2007-2014), and stream reach watersheds (Kentucky Division of Mine Permits 2004) were used to determine watershed boundaries and water sample locations. Three field trips were conducted to collect water samples: 15 May 2017 (spring), 21-24 July 2017 (summer), and 21-25 October 2017 (fall). Weather was stable and there was no surface runoff for any of the sampling days. During the first two trips, 9 (spring) and 14 (summer) water samples were collected at drainage exit points. Any water sample collected from the watersheds’ exit points represented the watersheds according to their geological structures and land cover types [26]. Under dry (no precipitation) conditions, a perennial stream in a watershed is recharged by groundwater and all stream water exits through the lowest elevation point within the watershed. Three water samples were collected from each exit point at five-meter intervals towards the headstream.

This study followed Kentucky ambient/watershed water quality monitoring procedures for collecting water samples [42]. Water samples were collected in high density polyethylene (HDPE) bottles that were cooled to 4°C. The samples collected in the spring and summer were sent to the West Virginia University National Research Center for Coal and Energy Water Analysis Lab for analysis. The samples were filtered at the laboratory through 11-micron filter paper (Whatman 1001-125, Pennsylvania, USA) to separate the sediments from the water. The samples then were analyzed for SO4²-, alkalinity, conductivity, Ca²⁺, Mg²⁺, Mn²+, Al³+, and Fe^2+,Fe3+. In the fall, the third water sampling analysis was conducted for in situ conductivity measurements with a larger sample size (n = 58) using a Hydrolab Quanta Water Probe (OTT Hydromet, Kempten, Germany). The conductivity levels were measured at the drainage exit points during the sampling. Data from the spring and summer sampling were used for bivariate correlation analysis. Data from the fall visit, which had the largest sample size, were used in the regression models. Figure 6 is the schematic of water quality data preparation and analysis steps.

2.1.2. Vegetation Cover Change Data Collection and Analysis Method

We used ArcGIS Desktop 10.5 for all GIS analyses [43]. Landsat and National Aerial Imagery Program (NAIP) images, and Digital Orthophoto Quadrangle images acquired between 1986 and 2017 (NASA Landsat Program 1986-2017; Kentucky Geography Network 1990-2016). Additionally, the National Land Cover Dataset (NLCD) was downloaded from Multi-Resolution Land Characteristics Consortium (MRLC 1992-2011). These images were compiled for the mined and unmined boundaries to derive the vegetation cover data and accuracy assessment. Spectral reflectance values of near infrared (NIR) wavebands are commonly used in combination with other spectral bands to detect different land cover classes [44,45,46,47]. Satellite images are particularly beneficial for detecting land cover conversion from natural vegetation cover (e.g., deforestation; [48,50–54].

[55,56] compared four different vegetation indices (Normalized Difference Vegetation Index [NDVI], Normalized Difference Moisture Index [NDMI], Normalized Burn Ratio [NBR], and Tasseled Cap Green-Blue [TC G-B]) for known land covers (barren, vegetation), to delineate coal mine disturbed areas; NDVI provided the most accurate results (Equation 4).

NDVI = (NIR - RED)/ (NIR + RED)

Eq. 4

We followed this methodology to derive the NDVI maps between 1986 and 2017 for vegetation cover change analysis. NDVI measures vegetation greenness by computing the proportion of visible red and NIR spectral reflectance [57,58].

We studied unmined (mining percentage ≤5%) and mined sites covered with undisturbed forest, reclaimed land (reclaimed forest, reclaimed woods, and reclaimed grass/pasture land), and barren ground. We extracted NDVI maps from Landsat 5 and 8 clear images (<10% cloud coverage), which were taken between 1986 and 2014.

Using the Image Analysis tool of ArcMap (ArcGIS Desktop) 10.5 and Equation 4, time series NDVI maps were derived from each Landsat images acquired in 1986, 1990, 1994, 1999, 2002, 2006, 2010, 2014, and 2017 since RGB-NDVI produces the most accurate results [59]. The NDVI data were further grouped into five classes: 1) barren (active mine); 2) reclaimed grassland; 3) reclaimed woodland (<50% forest cover); 4) reclaimed forest (>50% forest cover); and 5) undisturbed forest using the [3,59] framework. Assigning years to the barren pixels with con tool and running cell statistics tool created a composite map (Figure 7). The resulting map displayed the location of barren (active mine) areas with the latest mining year. Afterwards, we compared barren areas with 2017 vegetation classifications. For instance, we compared an area that was mined in 1994 and then later converted to reclaimed forest and was reclaimed woods in 2017. Consequently, temporal NDVI class change within a specific area was considered to be vegetation cover change for that specific area and time. For instance, forest that was converted to barren ground during the coal mining activity in a particular year was expected to be converted to reclaimed forest after several years. Vegetation cover change was estimated based on vegetation improvement from no vegetation (barren) to final vegetation cover type (reclaimed grass, reclaimed woods or reclaimed forest). We evaluated the classified pixels from National Land Cover Dataset (NLCD), Digital Orthophoto Quadrangle (DOQ), and NAIP imageries. If more than half of a pixel contained trees, they were classified as forest, if less than half they were classified as woods regardless of the remained class or classes in a pixel.

Vegetation cover change was evaluated from the latest barren land to the current (2017) vegetation cover. The Zonal histogram tool of ArcMap 10.5 produced counts of the classified pixels (barren, grass, woods, forest) in each watershed.

Figure 7. Active mining areas by years between 1986 and 2017.

Figure 7. Active mining areas by years between 1986 and 2017.

We conducted an accuracy assessment to evaluate and validate the land cover or NDVI classes derived from the satellite images using ancillary or secondary or in situ (ground truth) data [60]. NAIP, digital orthophoto quadrangle images, and NLCD were used to assess classification accuracy. Randomly selected (a total of 200 points, 50 for each class) stratified points were used for data from 1990, 2002, and 2014 [3]. We used NLCD (1992) and Digital Orthophoto Quadrangle (DOQ) (1985) images for accuracy assessment of data from 1990, NLCD (2004), DOQ (2000-2001), and NAIP (2004) images for data from 2002; NLCD (2011), and NAIP (2014) for data from 2014. The criteria used for separating reclaimed woods from reclaimed forest was visual inspection of the abundance of trees in a pixel. If a pixel contained more than 50% trees, it was classified as reclaimed forest. The accuracy assessment results provided overall accuracy of 79.5%, 80.5%, and 89.5% for 1990, 2002, and 2014, respectively. We obtained the following Kappa () for the same time periods: 67.6%, 73.1%, and 83.1%, respectively. Overall accuracy exceeded the target threshold accuracy (80%) as assumed high accuracy by [3]. The overall accuracy and Kappa statistics suggest that the vegetation cover change data were valid and usable for further statistical analysis. Accuracy assessment results are reported in Tables 3a, 3b, and 3c.

Table 3a. Accuracy assessment for land cover classes for 1990.

		Producer
	Classification	Barren	Grass	Woods	Forest	Row Total
User	Barren	21	1	0	2	24
	Grass	0	39	5	14	58
	Woods	0	2	11	4	17
	Forest	1	8	4	88	101
	Column Total	22	50	20	108	200
	Omission Error (%)	4.5	22	45.00	18.51
	Producer Accuracy (%)	95.5	78.0	55.0	81.5
	Commission Error (%)	12.5	48.7	35.3	12.9
	User Accuracy (%)	87.5	51.3	64.7	87.1
	Overall Accuracy (%) =	79.5	Kappa (%) =		67.6

Table 3a. Accuracy assessment for land cover classes for 1990.

		Producer
	Classification	Barren	Grass	Woods	Forest	Row Total
User	Barren	21	1	0	2	24
	Grass	0	39	5	14	58
	Woods	0	2	11	4	17
	Forest	1	8	4	88	101
	Column Total	22	50	20	108	200
	Omission Error (%)	4.5	22	45.00	18.51
	Producer Accuracy (%)	95.5	78.0	55.0	81.5
	Commission Error (%)	12.5	48.7	35.3	12.9
	User Accuracy (%)	87.5	51.3	64.7	87.1
	Overall Accuracy (%) =	79.5	Kappa (%) =		67.6

Table 3.b. Accuracy assessment for land cover classes for 2002.

		Producer
	Classification	Barren	Grass	Woods	Forest	Row Total
User	Barren	42	0	0	2	44
	Grass	1	41	5	19	66
	Woods	0	1	19	7	27
	Forest	0	1	3	59	63
	Column Total	43	43	27	87	200
	Omission Error (%)	2.3	4.6	29.60	32.1
	Producer Accuracy (%)	97.7	95.4	70.4	67.9
	Commission Error (%)	4.5	37.9	29.6	6.3
	User Accuracy (%)	95.5	62.1	70.4	93.7
	Overall Accuracy (%)	=80.5	Kappa (%)	=73.1

Table 3.b. Accuracy assessment for land cover classes for 2002.

		Producer
	Classification	Barren	Grass	Woods	Forest	Row Total
User	Barren	42	0	0	2	44
	Grass	1	41	5	19	66
	Woods	0	1	19	7	27
	Forest	0	1	3	59	63
	Column Total	43	43	27	87	200
	Omission Error (%)	2.3	4.6	29.60	32.1
	Producer Accuracy (%)	97.7	95.4	70.4	67.9
	Commission Error (%)	4.5	37.9	29.6	6.3
	User Accuracy (%)	95.5	62.1	70.4	93.7
	Overall Accuracy (%)	=80.5	Kappa (%)	=73.1

Table 3.c. Accuracy assessment for land cover classes for 2014.

		Producer
	Classification	Barren	Grass	Woods	Forest	Row Total
User	Barren	5	0	0	0	5
	Grass	2	39	2	3	46
	Woods	4	0	35	5	44
	Forest	0	1	4	100	105
	Column Total	11	40	41	108	200
	Omission Error (%)	54.5	2.5	14.60	7.4
	Producer Accuracy (%)	45.5	97.5	85.4	92.6
	Commission Error (%)	0	15.2	20.4	4.8
	User Accuracy (%)	100	84.8	79.6	95.2
	Overall Accuracy (%)	=89.5	Kappa (%)	=83.1

Table 3.c. Accuracy assessment for land cover classes for 2014.

		Producer
	Classification	Barren	Grass	Woods	Forest	Row Total
User	Barren	5	0	0	0	5
	Grass	2	39	2	3	46
	Woods	4	0	35	5	44
	Forest	0	1	4	100	105
	Column Total	11	40	41	108	200
	Omission Error (%)	54.5	2.5	14.60	7.4
	Producer Accuracy (%)	45.5	97.5	85.4	92.6
	Commission Error (%)	0	15.2	20.4	4.8
	User Accuracy (%)	100	84.8	79.6	95.2
	Overall Accuracy (%)	=89.5	Kappa (%)	=83.1

NDVI-based land cover classes provided an opportunity to compute average reclamation age (Equation 5 and 6); mining percentage (Equation 7); reclaimed forest percentage (Equation 8); and reclaimed woods percentage (Equation 9) for each watershed area studied. Many of the watersheds were mined multiple times in different years; therefore, average reclamation year for the watersheds was calculated by taking the weighted average of the initial reclamation year (latest mining year) and pixel counts (Equation 5).

$$Avg.Rec.Year=\frac{1986*\left(BarrenPixelCount\right)+\dots \dots .+2017*\left(BarrenPixelCount\right)}{TotalBarrenPixelCount}$$

Eq. 5

Subtracting the average reclamation year from 2017 (Equation 6) provided the reclamation age (years).

$$AverageReclamationAge=2017-avgrecyear$$

Eq. 6

Mining percentage was defined as ratio of sum of active years’ pixel values to total pixels counted for a watershed, multiplied by 100:

$$MiningPercentage\left(\mathrm{\%}\right)=\frac{\left(MinedPixelsof1986\right)+\dots .+\left(MinedPixelsof2017\right)}{TotalCountofPixelsforawatershed}*100$$

Eq. 7

Equation 8. was used to calculate reclaimed forest percentage:

$$ReclaimedForest(\mathrm{\%})=\frac{ReclaimedForestPixelCount(1986\dots 2017)}{TotalCountofPixelsforawatershed}*100$$

Eq. 8

Reclaimed woods percentage was calculated as follows:

$$ReclaimedWoods(\mathrm{\%})=\frac{ReclaimedWoodsPixelCount(1986..2017)}{TotalCountofPixelsforawatershed}*100$$

Eq. 9

Four variables were derived for further analysis: (1) Mined (mining percentage); (2) Reclamation Age (average reclamation age); (3) Reclaimed Forest; and (4) Reclaimed Woods.

2.1.3. Topographic Data

KYAPED DEM (5 ft.) aggregated to 30 m. horizontal spatial resolution was used to derive topographic variables. Then, mean elevation and mean slope for each watershed were calculated with the zonal statistics tool. Drainage density (density of streams in length in a watershed) was computed using the following equation [61].

$$D_d=(\sum L)/A$$

Eq. 10

• D_d: Drainage density
• ∑L: Total length of streams within the watershed
• A: Area of the watershed

The NHD 100K stream shapefile was extracted with the identity and intersect tools and all streams were assigned to their watershed using the spatial join tool. Next, total stream length was calculated with the summary statistics tool for each watershed. Then, drainage density was obtained for each watershed using Equation 10.

Hydrologic Soil Group maps were created with the Soil Data Viewer [62] tool. Hydrologic Soil Group maps comprise soil infiltration rates with seven ordinal categories.

Hydrologic Soil Group Classification

Hydrologic soil groups are based on estimates of runoff potential. Soils are assigned to one of four groups according to the rate of water infiltration when the soils are not protected by vegetation, are thoroughly wet, and receive precipitation from long-duration storms.

The soils in the United States are assigned to four groups (A, B, C, and D) and three dual classes (A/D, B/D, and C/D). The groups are defined as follows:

Group A. Soils having a high infiltration rate (low runoff potential)

Group B. Soils having a moderate infiltration rate

Group C. Soils having a slow infiltration rate

Group D. Soils having a very slow infiltration rate (high runoff potential)

If a soil is assigned to a dual hydrologic group (A/D, B/D, or C/D), the first letter is for drained areas and the second is for undrained areas. Only the soils that in their natural condition are in group D are assigned to dual classes. Source: USDA NRCS, 2015

These categories were translated to numbers. Based on categories and their representative scores, mean infiltration was calculated for each watershed with zonal statistics.

2.2. Empirical Models

We used mined, reclaimed woods, reclaimed forest, and reclamation age derived from Landsat images as vegetation cover and mined (mining percentage) variables (Table 4) to develop three regression models as specified in Equations 11, 12, and 13. All topographic variables added to the models as control variables (elevation, slope, drainage density, infiltration). Infiltration was an ordinal variable whereas rest of the variables comprise continuous data only.

Conductivity =β₀+ β_mined X_mined + β_{reclamation age}X_{reclamation age} + β_elevationX_elevation + β_slope X_slope+ β_{drainage density} X_{drainage density +} β_infiltrationX_infiltration + e

Eq. 1

First, we aimed to predict conductivity without vegetation cover effect with all other independent and control variables.

Conductivity =β₀+ β_minedX_mined + β_{reclamation age}X_{reclamation age} + β_elevationX_elevation + β_slope X_slope+ β_{drainage density} X_{drainage. density +} β_infiltrationX_infiltration + β_{reclaimed forest}X_{reclaimed forest} + β_{reclaimed woods}X_{reclaimed woods} + e

Eq. 12

Second, we added the vegetation cover variables to the previous model to predict conductivity and determine how vegetation recovery affects conductivity.

Reclaimed Forest =β₀+ β_mined X_mined + β_{reclamation age}X_{reclamation age} + β_elevationX_elevation + β_slope X_slope+ β_{drainage density} X_{drainage. density +} β_infiltrationX_infiltration + e

Eq. 13

Finally, we created an empirical model to predict the vegetation recovery with independent variables.

Where,

β₀ = Constant; β = Unstandardized beta coefficients for independent variables, and Conductivity = dependent variable (Eq. 12), and Reclaimed Forest = dependent variable (Eq. 13) and -e = error term.

2.3. Statistical Analyses

Standard deviation, skewness, kurtosis, standard errors of skewness, and kurtosis values were calculated for each variable as part of the descriptive statistics. Kurtosis values between -2 and 2 were accepted as an indication of approximately normal distribution [63,64]. Skewness between -1 and 1 was accepted as a normal distribution [65]. In addition, histogram plots with normal distribution curves corroborated descriptive statistics visually. Multicollinearity was tested with tolerance and variance inflation factor (VIF); all independent variables had tolerance values larger than 0.1 and VIF values smaller than 5.0 indicating absence of multicollinearity [66]. Durbin-Watson value (1.5 < d < 2.5) [67] assured absence of auto-correlation. In addition, scatterplots of residual distribution, P-P plot, and standardized residual-standardized predicted values plot were computed to validate homoscedasticity. The data were analyzed using SPSS v 22. This study accepted results with p≤0.05 as significant in the models.

3. Results and Discussion

3.1. Descriptive statistics

This research covers mined areas between 1986 and 2017; however, considering Kentucky’s long history of coal production, mined (mining percentage) areas were expected to be more than what Table 5 displays. The mined without reclaimed forests (RF) variable shows similar values to the mined area variable values because forest recovery was low as a percentage of a watershed. Average RF (reclaimed forest) percentage was 8.21%, which means only 21% of the mined area of a watershed converted to reclaimed forest between 1986 and 2017 in consideration of average values. Mean reclaimed woods was 17.42%, which was more than double the amount of reclaimed forest. Numerous studies found a link between macroinvertebrate community structures and physicochemical parameters in streams [22,25,41,49]. For example, [41,49] reported conductivity of greater than 300 µS/cm negatively affects the aquatic environment. Moreover, conductivity >500 µS/cm indicates ecological impairment and decreased biological diversity [22,25]. The results of our study indicate that in 42 of 58 (72%) watersheds in the study area, conductivity was higher than 300 µS/cm and in 32 (55%) of them, conductivity was over 500 µS/cm. Therefore, impairment and negative impact is expected in aquatic ecology due to change in natural biogeochemical processes in more than half of the watersheds although we did not perform any biological assessment.

3.2. Bivariate Correlations

Among all variables, reclaimed forest and drainage density variables were transformed using square root transformation to gain normal distribution. Linearity between independent and dependent variables were checked by using scatter plots prior to running bivariate correlations and Ordinary Least Square (OLS) regression models. Pearson bivariate correlation matrixes for all sampling dates and variables were employed to find significant correlations among variables (Table 6a, 6b, 6c).

There was no linear or curvilinear relationship between reclaimed forest and conductivity variables; however, the remaining scatter plots displayed linear relations between independent and dependent variables. Subsequently, another variable was created to test reclaimed forest influence on conductivity, i.e., mined without reclaimed forest (mined without RF or mined w/o RF). This variable represented the numerical difference between mined and reclaimed forest. We expected a potential influence on reclaimed forest from mined and reclamation age variables. While checking the data assumptions, a curvilinear relationship was discovered between mined and reclaimed forest variables. However, log transformation of the mined variable yielded linearity between the variables. We did not find any correlation among vegetation cover, mining percentage, conductivity and Al³⁺, Fe^2+,Fe³⁺, Mn²⁺. We concluded that either metals were not dissolved and filtered out during the filtering process or Al³⁺, Fe^2+,Fe³⁺, Mn² levels were not sufficient to detect the correlations.

Comparisons between independent and dependent variables are displayed visually with maps and symbols on Figure 8, Figure 9, Figure 10, Figure 11, Figure 12, and Figure 13. LogMined (Independent V.) versus normalized reclamation forest percentage (Dependent V.) values (n=58, R=0.45, p< 0.01). Bivariate correlation matrix results were consistence with linearity validation results. Linear relationships that were observed in scatterplots were significant at the p≤0.05 or p≤0.01 level (Table 6a, 6b, 6c). Bivariate correlation analysis demonstrated that conductivity was strongly related to CMRSC variables [39,68], except for Fe^2+,Fe3+, Al³⁺, and Mn²⁺. An in-depth literature review also suggested conductivity was highly associated with surface coal mining in eastern Kentucky and West Virginia [26,39,41]. Based on this information, we decided to use conductivity as a representative variable for CMRSC variables in multivariate regression models.

We did not find significant correlation between topographic variables and conductivity (fall) except infiltration. There was a strong correlation between infiltration and conductivity.

Table 6.a. Pearson correlation matrix for spring field study results (n=9).

	Mined (%)	Mined w/o RF (%)	Reclaimed Forest (%)	Reclaimed Woods (%)	Reclamation Age (Years)	Alkalinity (mg/L)	Conductivity (µS/cm)	SO42- (mg/L)	Ca2+ (mg/L)	Mg2+ (mg/L)
Mined (%)	1	0.980**	0.886**	0.920**	0.255	0.751*	0.761*	0.749*	0.776*	0.791*
Mined w/o RF (%)	0.980**	1	0.775*	0.940**	0.155	0.785*	0.794*	0.775*	0.797*	0.837**
Reclaimed Forest (%)	0.886**	0.775*	1	0.726*	0.443	0.553	0.563	0.568	0.601	0.558
Reclaimed Woods (%)	0.920**	0.940**	0.726*	1	0.338	0.907**	0.884**	0.843**	0.903**	0.914**
Reclamation Age (Years)	0.255	0.155	0.443	0.338	1	0.241	0.158	0.182	0.218	0.159
Alkalinity (mg/L)	0.751*	0.785*	0.553	0.907**	0.241	1	0.985**	0.897**	0.993**	0.962**
Conductivity (µS/cm)	0.761*	0.794*	0.563	0.884**	0.158	0.985**	1	0.933**	0.987**	0.969**
SO4 (mg/L)	0.749*	0.775*	0.568	0.843**	0.182	0.897**	0.933**	1	0.921**	0.964**
Ca (mg/L)	0.776*	0.797*	0.601	0.903**	0.218	0.993**	0.987**	0.921**	1	0.974**
Mg (mg/L)	0.791*	0.837**	0.558	0.914**	0.159	0.962**	0.969**	0.964**	0.974**	1

*. Correlation is significant at the 0.05 level (2-tailed). ** Correlation is significant at the 0.01 level (2 tailed).

Table 6.a. Pearson correlation matrix for spring field study results (n=9).

	Mined (%)	Mined w/o RF (%)	Reclaimed Forest (%)	Reclaimed Woods (%)	Reclamation Age (Years)	Alkalinity (mg/L)	Conductivity (µS/cm)	SO42- (mg/L)	Ca2+ (mg/L)	Mg2+ (mg/L)
Mined (%)	1	0.980**	0.886**	0.920**	0.255	0.751*	0.761*	0.749*	0.776*	0.791*
Mined w/o RF (%)	0.980**	1	0.775*	0.940**	0.155	0.785*	0.794*	0.775*	0.797*	0.837**
Reclaimed Forest (%)	0.886**	0.775*	1	0.726*	0.443	0.553	0.563	0.568	0.601	0.558
Reclaimed Woods (%)	0.920**	0.940**	0.726*	1	0.338	0.907**	0.884**	0.843**	0.903**	0.914**
Reclamation Age (Years)	0.255	0.155	0.443	0.338	1	0.241	0.158	0.182	0.218	0.159
Alkalinity (mg/L)	0.751*	0.785*	0.553	0.907**	0.241	1	0.985**	0.897**	0.993**	0.962**
Conductivity (µS/cm)	0.761*	0.794*	0.563	0.884**	0.158	0.985**	1	0.933**	0.987**	0.969**
SO4 (mg/L)	0.749*	0.775*	0.568	0.843**	0.182	0.897**	0.933**	1	0.921**	0.964**
Ca (mg/L)	0.776*	0.797*	0.601	0.903**	0.218	0.993**	0.987**	0.921**	1	0.974**
Mg (mg/L)	0.791*	0.837**	0.558	0.914**	0.159	0.962**	0.969**	0.964**	0.974**	1

*. Correlation is significant at the 0.05 level (2-tailed). ** Correlation is significant at the 0.01 level (2 tailed).

Table 6.b. Pearson correlation matrix for summer field study results (n=14).

	Mined (%)	Mined w/o RF (%)	Reclaimed Forest (%)	Reclaimed Woods (%)	Reclamation Age (Years)	Alkalinity (mg/L)	Conductivity (µS/cm)	SO42- (mg/L)	Ca2+ (mg/L)	Mg2+ (mg/L)
Mined (%)	1	0.940**	0.125	0.616**	-0.340	0.467*	0.734**	0.715**	0.769**	0.727**
Mined w/o RF (%)	0.940**	1	-0.221	0.495*	-0.627**	0.614**	0.838**	0.787**	0.867**	0.833**
Reclaimed Forest (%)	0.125	-0.221	1	0.320	0.847**	-0.450	-0.336	-0.243	-0.323	-0.342
Reclaimed Woods (%)	0.616**	0.495*	0.320	1	0.017	0.294	0.508*	0.485*	0.483*	0.504*
Reclamation Age (Years)	-0.340	-0.627**	0.847**	0.017	1	-0.534*	-0.603*	-0.514*	-0.610*	-0.599*
Alkalinity (mg/L)	0.467*	0.614**	-0.450	0.294	-0.534*	1	0.705**	0.604*	0.666**	0.760**
Conductivity (µS/cm)	0.715**	0.787**	-0.243	0.485*	-0.514*	0.604*	0.980**	1	0.972**	0.960**
SO4 (mg/L)	0.734**	0.838**	-0.336	0.508*	-0.603*	0.705**	1	0.980**	0.988**	0.992**
Ca (mg/L)	0.769**	0.867**	-0.323	0.483*	-0.610*	0.666**	0.988**	0.972**	1	0.967**
Mg (mg/L)	0.727**	0.833**	-0.342	0.504*	-0.599*	0.760**	0.992**	0.960**	0.967**	1

*. Correlation is significant at the 0.05 level (2-tailed). ** Correlation is significant at the 0.01 level (2-tailed).

Table 6.b. Pearson correlation matrix for summer field study results (n=14).

	Mined (%)	Mined w/o RF (%)	Reclaimed Forest (%)	Reclaimed Woods (%)	Reclamation Age (Years)	Alkalinity (mg/L)	Conductivity (µS/cm)	SO42- (mg/L)	Ca2+ (mg/L)	Mg2+ (mg/L)
Mined (%)	1	0.940**	0.125	0.616**	-0.340	0.467*	0.734**	0.715**	0.769**	0.727**
Mined w/o RF (%)	0.940**	1	-0.221	0.495*	-0.627**	0.614**	0.838**	0.787**	0.867**	0.833**
Reclaimed Forest (%)	0.125	-0.221	1	0.320	0.847**	-0.450	-0.336	-0.243	-0.323	-0.342
Reclaimed Woods (%)	0.616**	0.495*	0.320	1	0.017	0.294	0.508*	0.485*	0.483*	0.504*
Reclamation Age (Years)	-0.340	-0.627**	0.847**	0.017	1	-0.534*	-0.603*	-0.514*	-0.610*	-0.599*
Alkalinity (mg/L)	0.467*	0.614**	-0.450	0.294	-0.534*	1	0.705**	0.604*	0.666**	0.760**
Conductivity (µS/cm)	0.715**	0.787**	-0.243	0.485*	-0.514*	0.604*	0.980**	1	0.972**	0.960**
SO4 (mg/L)	0.734**	0.838**	-0.336	0.508*	-0.603*	0.705**	1	0.980**	0.988**	0.992**
Ca (mg/L)	0.769**	0.867**	-0.323	0.483*	-0.610*	0.666**	0.988**	0.972**	1	0.967**
Mg (mg/L)	0.727**	0.833**	-0.342	0.504*	-0.599*	0.760**	0.992**	0.960**	0.967**	1

*. Correlation is significant at the 0.05 level (2-tailed). ** Correlation is significant at the 0.01 level (2-tailed).

Table 6.c. Pearson correlation matrix for fall field study results (n=58).

	Mined (%)	Reclaimed Forest (%)	Reclaimed Woods (%)	Reclamation Age (year)	Conductivity (µS/cm)	Log Mined	Mined w/o RF (%)	Infiltration
Mined (%)	1	0.223	0.831**	-0.496**	0.863**	0.908**	0.965**	-0.678**
Reclaimed Forest (%)	0.223	1	0.246	0.500**	-0.072	0.451**	-0.032	-0.048
Reclaimed Woods (%)	0.831**	0.246	1	-0.329*	0.720**	0.796**	0.799**	-0.591**
Reclamation Age (year)	-0.496**	0.500**	-0.329*	1	-0.672**	0.294*	0.636**	-0.396**
Conductivity (µS/cm)	0.863**	-0.072	0.720**	-0.672**	1	0.726*	0.904**	-0.676**
Log Mined	0.908**	0.451**	0.796**	0.294*	0.726**	1	0.823**	0.580**
Mined w/o RF (%)	0.965**	-0.032	0.799**	0.636**	0.904**	0.823**	1	-0.681**
Infiltration	-0.678**	-0.048	-0.591**	0.396**	-0.676**	-0.580**	-0.681**	1

*. Correlation is significant at the 0.05 level (2-tailed). ** Correlation is significant at level 0.01.

Table 6.c. Pearson correlation matrix for fall field study results (n=58).

	Mined (%)	Reclaimed Forest (%)	Reclaimed Woods (%)	Reclamation Age (year)	Conductivity (µS/cm)	Log Mined	Mined w/o RF (%)	Infiltration
Mined (%)	1	0.223	0.831**	-0.496**	0.863**	0.908**	0.965**	-0.678**
Reclaimed Forest (%)	0.223	1	0.246	0.500**	-0.072	0.451**	-0.032	-0.048
Reclaimed Woods (%)	0.831**	0.246	1	-0.329*	0.720**	0.796**	0.799**	-0.591**
Reclamation Age (year)	-0.496**	0.500**	-0.329*	1	-0.672**	0.294*	0.636**	-0.396**
Conductivity (µS/cm)	0.863**	-0.072	0.720**	-0.672**	1	0.726*	0.904**	-0.676**
Log Mined	0.908**	0.451**	0.796**	0.294*	0.726**	1	0.823**	0.580**
Mined w/o RF (%)	0.965**	-0.032	0.799**	0.636**	0.904**	0.823**	1	-0.681**
Infiltration	-0.678**	-0.048	-0.591**	0.396**	-0.676**	-0.580**	-0.681**	1

*. Correlation is significant at the 0.05 level (2-tailed). ** Correlation is significant at level 0.01.

3.3. Multivariate Regression Models Results

3.3.1. Effects of mined and reclamation age parameters on conductivity

Multivariate regression test results predicted conductivity (R²=0.825, adjusted R²=.818, F(2,55) = 129.213, p < 0.01). The regression model accounted for 81.8% of variance in conductivity prediction (Table 7a). Mined and reclamation age were significant (p < 0.01; Table 7b). Standardized beta (Β) for mined and reclamation age coefficients were 0.70 and -0.32, respectively (Table 7b), which indicate that the mined parameter was positively and strongly correlated with conductivity whereas reclamation age was significantly, but negatively, related to conductivity. Partial correlations for coefficients also confirmed that the mined parameter was more influential on conductivity than reclamation age. Infiltration was not a significant factor so it is excluded from the model. The following regression equation was predicted for conductivity:

Conductivity = 16.96*Mined-30.86*Reclamation Age+610.713

Eq. 14

We found, for this model, that the constant value was higher than [26,40], suggesting that it might have been due to inadequate control sites. However, the coefficients were comparable to those found by [26,40]. This supports the hypothesis that reclamation age has less influence on conductivity and that conductivity was best predicted with mining percentage [26,40,41]; however, none of these studies reported a correlation between reclamation age and conductivity.

Table 7.a. Regression Model A summary.

Regression Model A Summary –Conductivity Prediction
Step	R	R2	R2 adj	Δ R2	Fchg	p	df1	df2
Mining Percentage	0.863a	0.745	0.741	0.745	163.98	<0.01	1	56
Reclamation Age	0.908b	0.825	0.818	0.08	129.21	<0.01	1	55

Table 7.a. Regression Model A summary.

Regression Model A Summary –Conductivity Prediction
Step	R	R2	R2 adj	Δ R2	Fchg	p	df1	df2
Mining Percentage	0.863a	0.745	0.741	0.745	163.98	<0.01	1	56
Reclamation Age	0.908b	0.825	0.818	0.08	129.21	<0.01	1	55

Table 7.b. Regression Model A coefficients.

Coefficients for Final Model-Conductivity Prediction
Model	β	Β	t	p	Bivariate r	Partial r
(Constant)	610.732		4.303	<0.01
Mining Percentage	16.956	0.703	10.804	<0.01	0.863	0.610
Reclamation Age (year)	-30.860	-0.324	-4.979	<0.01	-0.672	-0.281

Table 7.b. Regression Model A coefficients.

Coefficients for Final Model-Conductivity Prediction
Model	β	Β	t	p	Bivariate r	Partial r
(Constant)	610.732		4.303	<0.01
Mining Percentage	16.956	0.703	10.804	<0.01	0.863	0.610
Reclamation Age (year)	-30.860	-0.324	-4.979	<0.01	-0.672	-0.281

3.3.2. Effects of vegetation cover change on conductivity

Bivariate Pearson correlation analysis between mined without RF and conductivity was very strong and positive (R=0.90; Table 6c). Multivariate regression test results indicated significant overall goodness of fit for two predictors (mined without RF and reclamation age) that significantly predict conductivity (R²=0.832, adjusted R²= 0.826, F(2,55) =136.659, p < 0.01). The regression model accounted for 82.6% of the variance in the conductivity prediction (Table 8a). The model showed mined without RF (p < 0.01) and reclamation age (p < 0.026) were significant coefficients (Table 8b). Unstandardized beta coefficients (β) were 0.799 (mined w/o reclaimed forest) and -0.164 (reclamation age). Pearson bivariate correlation results suggested a strong positive correlation between reclaimed woods, infiltration, and conductivity. Conversely, the model showed that reclaimed woods and/or infiltration was not a significant contributor for predicting conductivity. P values were higher than 0.05 and t values were smaller than 1 for reclaimed woods; therefore, reclaimed woods was excluded in both steps. Multivariate analysis yielded the following equation for predicting conductivity:

Conductivity= 19.54*(Mined-Reclaimed Forest)-15.64*Reclamation Age+429.20

Eq. 15

The model had a slightly higher adjusted R2 value (0.826) than Regression Model A (0.818). There was no difference between undisturbed forest and reclaimed forest from a land cover standpoint, suggesting that reclaimed forest areas can be assumed as not disturbed. Comparing the coefficients of mined and mined without RF parameters, the mined without RF was a better predictor of conductivity.

We were able to improve the model by removing the reclaimed forested areas. The result indicated a positive correlation between reclaimed forest and conductivity mitigation [16]. Our findings agree with [7], who showed that there is a potential for total dissolved solids-source-control practices that incorporate FRA to improve mine-water discharge quality. The vegetation data covered a 31-year period (1986-2017) since data prior to 1986 were not reliable. Measurements of the geographic region of the reforested areas for the entire 31-year period have proven challenging due to data inconsistency. It should be noted that reclaimed forest percentage range was relatively low within the given time interval (1986-2017). This might have caused the nonexistence of a direct relationship between reclaimed forest and conductivity. However, recent improvements in Landsat imagery will likely provide more reliable data. In addition, forest recovery is expected to increase in future. This experiment should be repeated to obtain a clear idea about the influence of reclaimed vegetation on conductivity and water quality.

Table 8.a. Regression Model B summary.

Regression Model B Summary – Conductivity Prediction
Step	R	R2	R2 adj	Δ R2	Fchg	p	df1	df2
Mined w/o RF (%)	0.904	0.816	0.813	0.816	249.09	<0.01	1	56
Reclamation Age (years)	0.912	0.832	0.826	0.016	5.26	<0.026	1	55

Table 8.a. Regression Model B summary.

Regression Model B Summary – Conductivity Prediction
Step	R	R2	R2 adj	Δ R2	Fchg	p	df1	df2
Mined w/o RF (%)	0.904	0.816	0.813	0.816	249.09	<0.01	1	56
Reclamation Age (years)	0.912	0.832	0.826	0.016	5.26	<0.026	1	55

Table 8.b. Regression Model B coefficients.

Coefficients for Final Model - Conductivity Prediction
Model	β	B	t	p	Bivariate r	Partial r
(Constant)	429.20		2.85	<0.01
Mined w/o RF (%)	19.54	0.799	11.17	<0.01	0.904	0.617
Reclamation Age (years)	-15.64	-0.164	-2.29	<0.026	-0.672	-0.127

Table 8.b. Regression Model B coefficients.

Coefficients for Final Model - Conductivity Prediction
Model	β	B	t	p	Bivariate r	Partial r
(Constant)	429.20		2.85	<0.01
Mined w/o RF (%)	19.54	0.799	11.17	<0.01	0.904	0.617
Reclamation Age (years)	-15.64	-0.164	-2.29	<0.026	-0.672	-0.127

3.3.3. Effects of reclamation age and mined operation (mining percentage) on reclaimed forest

We used multivariate regression to predict the relationship between the dependent variable reclaimed forest and the independent variables of reclamation age and mining percentage. The results indicated that the model was significant (R²=0.641, adjusted R²=0.628, F(2,55)=49.02, p < 0.01). The regression model accounted for 62.8% of variance in the reclaimed forest prediction (Table 9a), and mined and reclamation age were significant coefficients (p < 0.01 for both; Table 9b). Adding the mined parameter improved the reclaimed forest prediction greatly (R² changed from 0.250 to 0.641). Standardized beta coefficients (B) showed that log mined (0.692) and reclamation age (0.654) were positive and had similar effects in the reclaimed forest prediction (Table 9b). None of the topographic variables were added to the model because no significant correlations were found between reclaimed forest percentage and topographic variables. Multivariate regression analysis yielded the following equation for the prediction of reclaimed forest:

Reclaimed Forest= [(1.833*log₁₀Mined)+(0.121*Reclamation Age)-1.953]²

Eq. 16

Table 9.a. Regression Model C summary.

Regression Model C Summary – Reclaimed Forest Prediction
Step	R	R2	R2adj	Δ R2	Fchg	p	df1	df2
Reclamation Age (years)	0.500	0.250	0.236	0.236	18.61	<0.01	1	56
Log Mined	0.800	0.641	0.628	0.392	35.65	<0.01	1	55

Table 9.a. Regression Model C summary.

Regression Model C Summary – Reclaimed Forest Prediction
Step	R	R2	R2adj	Δ R2	Fchg	p	df1	df2
Reclamation Age (years)	0.500	0.250	0.236	0.236	18.61	<0.01	1	56
Log Mined	0.800	0.641	0.628	0.392	35.65	<0.01	1	55

Table 9.b. Regression Model C coefficients.

Coefficients for Final Model Reclaimed Forest Prediction
Model	β	B	t	p	Bivariate r	Partial r
(Constant)	-1.953		-4.037	<0.01
Reclamation Age (years)	.121	0.692	8.185	<0.01	0.500	0.741
Log Mined	1.833	0.654	7.735	<0.01	0.451	0.722

Table 9.b. Regression Model C coefficients.

Coefficients for Final Model Reclaimed Forest Prediction
Model	β	B	t	p	Bivariate r	Partial r
(Constant)	-1.953		-4.037	<0.01
Reclamation Age (years)	.121	0.692	8.185	<0.01	0.500	0.741
Log Mined	1.833	0.654	7.735	<0.01	0.451	0.722

[69] reported that coal surface mine reclamation efforts in the southeastern USA are usually assessed after five years of mining activity. Conversely, the model estimates 25 years is necessary for 10% of a mined watershed to be converted to a reclaimed forest (Equation 16). These results suggest that five years may not be enough time to assess reclamation success. Many scientists have studied effective techniques to reclaim forests in areas impacted by coal mining by species composition, comparing mined and unmined sites, hydrologic properties, or best revegetation practices [7,19,70,71]. Conversely, not many studies measured the quantity of reclaimed forest in relation to reclamation age at the watershed level.

4. Conclusions

This study explored the relationships between water quality parameters, reclaimed forest percentage (dependent variables), mining percentage, reclamation age, reclaimed woods, slope, elevation, drainage density, and infiltration in stream-reach watersheds affected by coal mining. We found that reclamation age is a significant factor for predicting conductivity in reclaimed mines. Even though reclamation age was not a primary factor, Reclamation age appeared to increase the accuracy of the conductivity prediction although reclamation age was not a primary factor accounting for it increased the accuracy of the conductivity prediction. This research used more convenient and accurate methods for collecting and analyzing the data about mining percentage and reclamation age than the previous studies in the Appalachia region (states Furthermore, reclamation age and topographic factors were not evaluated on a regional scale in previous studies.

We also investigated the effects of reclaimed forest and reclaimed woods on conductivity to improve the Regression Model A, which included only mining percentage (mined) and reclamation age for conductivity prediction. We did not find a direct relation between conductivity and reclaimed forest percentage or reclaimed woods percentage; however, by removing reclaimed forest and performing correlation of mining percentage with conductivity, we were able to moderately improve the results. The Regression Model A results were consistent with earlier studies that found coal mining had a major impact on water quality degradation and forest structure.

In addition, we used reclamation age and mining percentage to estimate the quantity of reclaimed forested areas in the watershed. Reclaimed forest quality is usually measured by comparing the resemblance to original forest. We concluded that the quality of reclaimed forest is important for evaluating reclamation success; however, quantity is equally important. Quantity of a reclaimed forest coupled with quality, provides a better assessment of reclaimed areas when advanced remote sensing techniques are used. For example, high resolution aerial images or Light Detection and Ranging (LIDAR) data may be used to define the composition of reclaimed forested area with accurate delineation of mined area.

We included topographic variables (slope, elevation, drainage density, and infiltration) to take into account possible influences of topographic variations in the regression models. The results indicated no such relation except infiltration was positively correlated with conductivity (fall), but this was not significant. Studies by [25,35] found similar results.

Overall, this study suggests that conductivity is a predictable water quality indicator that is highly associated with CMRSC where agriculture and urban areas are limited. Furthermore, our assessment of vegetation cover change may provide insight into reclamation success in terms of restoring deforested areas. Our findings may help the scientific community and regulating agencies improve their understanding of water quality characteristics in watersheds affected by coal mining. This can be accomplished by refining land reclamation practices and monitoring protocols.