2.1 자료 수집

본 연구에 이용된 자료는 한강물환경연구소의 “기후변화가 수생태계에 미치는 영향과 대응전략(’10~’12)”, 환경부의 “환경생태유량 산정기준 연구 및 시범 산정(’14)”, 환경부⋅국립환경과학원의 “수생태계 건강성 조사 및 평가(’12~’13)”, “하천 수생태계 현황 조사 및 건강성 평가(’14~’21)”, “수생태계 참조하천 선정 및 활용방안 마련 연구(’18)”, “4대강 보 개방에 따른 수생태계 변화 조사(’18~’20)”, 국립환경과학원의 “생물측정망 모니터링 및 평가기법 개발연구(’17)”의 지점 또는 정점 조사 결과의 일부였으며 총 표본단위는 23,957개였다. 이 중 환경부 물환경측정망의 수질측정자료가 없는 지점의 생물조사자료는 생물과 수질과의 상관성 분석을 위해 사용할 수 없으므로 분석에서 제외하였다. 수질과 생물이 조사된 지점에서의 연중 또는 연간 반복된 조사를 포괄한 총 조사지점은 총 6,664개였다. 각 조사지점은 민통선 이남에서 제주도, 서해안 접경부터 동해안 접경, 해발고도 0 m부터 888 m까지 전국적인 범위를 포괄한 것이었다. 각 조사단위의 자료는 하루살이속 약충의 출현유무와 해당 조사지점의 수질에 대한 것이며 성충에 대한 조사결과는 포함되지 않은 것이다.

하루살이속 중 동양하루살이와 사할린하루살이는 약충 단계에서는 구분이 쉽지 않아 본 연구에 사용된 자료의 대부분은 동정 단계에서 사할린하루살이가 동양하루살이에 포함되어 기록된 것들이고 두 종은 동소서식하는 것으로 알려져 있으므로(^{Hwang et al., 2013}), 이를 동양하루살이-사할린하루살이 군으로 묶어서 분석하였다.

본 연구에서 고려한 화학적 환경요인은 하천의 유기오염을 대표하는 주요 항목인 5일 BOD (5-day biochemical oxygen demand, BOD₅) 농도와, 총인(Total phosphorus, T-P) 농도 및 총부유물질(Total suspended solids, TSS) 농도였다. 사용된 수질자료는 조사지점 중 환경부 물환경측정망의 수질측정망과 일치하는 지점에 해당하는 것만 선별하여 사용하였다.

2.2 조사자료

저서성 대형무척추동물 자료는 수심이 얕은 곳에서는 Surber net 또는 D-frame net을 사용하였고 이들을 사용 할 수 없는 수심이 깊은 곳에서는 Dredge, Ekman grab 또는 Ponar grab으로 정량채집되어 동정 및 계수된 개체수 밀도였다.

2.3 통계분석

2.3.2 확률분포모형 선정

본 연구에서 적용한 확률분포모형은 ^{Kong and Kang (2023)}의 연구에서 적용된 바와 같이 역치값을 포함한 2 인자 지수분포(Exponential distribution), 2인자 절단분포[정규분포(Normal distribution), 로지스틱분포(Logistic distribution), 굼벨분포(Gumbel distribution)], 역치값을 포함한 3 인자 분포 [대수정규분포(Lognormal distribution), 베이블분포(Weibull distribution), 감마분포(Gamma distribution)], 역치값을 포함한 4 인자 베타분포(Beta distribution)로 8개 유형이었으며, 이 중 4개의 모형(지수분포, 대수정규분포, 베이블분포, 감마분포)에 대해서는 부적편포(negatively skewed distribution)의 해석을 위한 역분포 모형(reversed distribution type)이 추가되었다(Table 1).

Table 1. Probability density function (PDF), cumulative distribution function (CDF), mode, mean, and variance according to the environmental variables ()

Exponential distribution (Type 1.1)		Type 1.1.1	Type 1.1.2
		,	,
	CDF	,
	PDF
	Mode	0
	Mean
	Variance
Reversed exponential distribution (Type 1.2)		Type 1.2.1 ,	Type 1.2.2
	CDF	,
	PDF
	Mode
	Mean
	Variance
Normal distribution (Type 2)		Type 2.1	Type 2.2
	CDF
	PDF
	Mode
	Mean	,
	Variance
Lognormal distribution (Type 3.1)		Type 3.1.1	Type 3.1.2
		,	,
	CDF
	PDF
Lognormal distribution (Type 3.1)	Mode
	Mean	,
	Variance
Reversed lognormal distribution (Type 3.2)		Type 3.2.1 ,	Type 3.2.2
	CDF
Reversed lognormal distribution (Type 3.2)	PDF
	Mode
	Mean	,
	Variance
Logistic distribution (Type 4)		Type 4.1	Type 4.2
	CDF
	PDF
	Mode
	Mean
	Variance	,
Weibull distribution (Type 5.1)		Type 5.1.1	Type 5.1.2
		,	,
	CDF	,
	PDF
	Mode
	Mean	: upper incomplete gamma function
	Variance
		Type 5.2.1 ,	Type 5.2.2
	CDF	,
Reversed Weibull distribution (Type 5.2)	PDF
	Mode
	Mean	: upper incomplete gamma function
	Variance
Gamma distribution (Type 6.1)		Type 6.1.1	Type 6.1.2
		,	,
	CDF
	PDF
	Mode
	Mean	: upper incomplete gamma function
	Variance
Reversed gamma distribution (Type 6.2)		Type 6.2.1 ,	Type 6.2.2
	CDF
	PDF
	Mode
	Mean	: upper incomplete gamma function
	Variance
Beta distribution (Type 7)		Type 7.1	Type 7.2
		,
	CDF
	PDF
Beta distribution (Type 7)	Mode	,  at  at  at	,  at  at  at
	Mean	: hypergeometric function
	Variance
Gumbel distribution (Type 8)		Type 8.1	Type 8.2
	CDF
	PDF
	Mode
	Mean	: Euler’s constant Ei: Exponential integral	: Euler’s constant
	Variance

수질에 대한 생물의 분포는 요인의 범위가 또는 인 조건에서 더미분포를 가진 유형(Fig. 1a)과 역치값을 가지고 또는 의 분포를 갖는 유형(Fig. 1b)으로 구분할 수 있다.

Fig. 1. Distribution type of aquatic organisms according to environmental factors; (a) with dummy, (b) shifted.

Fig. 1a의 유형을 제외하고는 확률변수의 정의역에서 일반적으로 알려진 PDF와 누적분포함수(cumulative distribution function, CDF) 및 통계량 함수를 적용할 수 있으나, 더미(dummy) 구간을 가진 Fig. 1a의 유형에 대해서는 확률변수의 정의역이 0 이상인 조건에 부합되도록 PDF와 CDF를 변형하고 통계량 함수를 유도하였다(Table 1).

3.1 BOD₅

BOD₅에 대한 각 확률분포모형의 CDF와 각 분류군의 출현도에서 산정된 CMF 간의 NRMSE는 Table 2와 같다. 동양하루살이-사할린하루살이 군, 무늬하루살이, 가는무늬하루살이 모두 최빈값을 중심으로 정적편포하였다. 동양하루살이-사할린하루살이는 대수정규분포의 적합도가 가장 높았으며 무늬하루살이는 굼벨분포의 적합도가 가장 높았으며 가는무늬하루살이는 역베이블분포의 적합도가 가장 높았다. 또한 이 세 분류군 모두 지수분포모형 및 역지수분포모형을 제외한 다른 모형들의 적합도 역시 전반적으로 높았다.

Fig. 2는 각 분류군의 BOD₅에 대한 PMF와 이에 대하여 최적의 적합도를 보이는 확률분포모형의 PDF 및 CDF를 도시한 것이다. BOD₅에 대한 각 분류군의 출현도는 분류군별로 차이가 컸으며 BOD₅의 범위에 따라 표본단위의 수도 차이가 컸기 때문에 각 급구간의 범위는 분류군의 출현특성에 따라 달리 설정하였다.

Fig. 2. Probability mass function based on the composite frequency of each taxon according to BOD5gradient, probability density function, and probability distribution function of the best-fit model; (a)Ephemera orientalis-sachalinensis, (b)Ephemera strigata, (c)Ephemera separigata, (d) all taxa.

동양하루살이-사할린하루살이 군은 최빈값인 0.9 mg/L 이상의 높은 농도의 BOD₅의 구간에서도 높은 출현도를 보였으며 무늬하루살이는 약 7.4 mg/L까지 출현하였지만 최빈값인 0.5 mg/L 이상에서 지속적으로 출현도가 감소하였고 가는무늬하루살이는 최빈값인 0.0 mg/L 이상에서 지속적으로 출현도가 뚜렷하게 감소하였다(Fig. 2, Table 3).

HSI에 따른 각 분류군의 BOD₅ 범위는 Table 4와 같다. HSI 값이 1.0에 해당하는 범위를 최적 적응범위라 본다면 Fig. 3의 50% 범위의 하한값과 상한값은 생물의 BOD₅에 대한 최적 적응범위의 경계값이라 볼 수 있다. 이에 따라 BOD₅에 대한 동양하루살이-사할린하루살이 군의 최적 적응범위는 0.4~1.9 mg/L, 무늬하루살이는 0.4~0.7 mg/L, 가는무늬하루살이는 0.0~0.4 mg/L로 분석되었다(Table 4). 무늬하루살이와 동양하루살이-사할린하루살이는 BOD₅에 대한 최적 적응범위에서 지위가 일부 중복되며 가는무늬하루살이는 다른 하루살이과 분류군과 지위가 뚜렷하게 구별됨을 확인할 수 있었다.

Fig. 3. Habitat suitability index ofEphemeraspecies for biochemical oxygen demand (BOD5) gradient; (a)Ephemera orientalis-sachalinensis, (b)Ephemera strigata, (c)Ephemera separigata, (d) all taxa.

서식처적합도(Table 4)와 중앙값과 표준편차(Table 3)를  고려할 때 동양하루살이-사할린하루살이 군은 낮은 농도부터 높은 농도의 BOD₅까지 넓게 적응된 광부수성(eurysaprobic) 이었으며, 무늬하루살이는 낮은 농도의 BOD₅에 적응한 빈부수성(oligosaprobic), 가는무늬하루살이는 낮은 농도의 BOD₅에 좁게 적응된 협빈부수성(stenooligosaprobic)으로 평가된다.

동양하루살이의 출현에 대하여 ^{Park and Park (1999)}은 동양하루살이가 BOD₅의 농도가 4.3 mg/L에서도 출현한다고 보고하였다. ^{Kong et al. (2018)}은 BOD₅의 농도를 변수로 한 생물종의 출현 특성을 확률적으로 분석하였다. 이때 산출된 오탁치는 동양하루살이는 0.8, 사할린하루살이는 0.7, 무늬하루살이는 0.1, 가는무늬하루살이는 0.1로 제시되었다. 이는 다른 부수성 생물의 오탁치에 비해 낮은 정도로 비교적 좋은 수질의 환경에 적응되어 있음을 확인 할 수 있다. 또한 환경부의 환경정책기본법 시행령 제 2조 별표1 환경기준(개정 2022.12.06.)의 하천수 수질환경기준에서 BOD₅의 농도가 1.0 mg/L 이하인 하천을 수질이 매우 좋다고 평가하였으며 5.0 mg/L 이하까지는 보통의 하천이라고 평가하였다. 이를 통해 동양하루살이-사할린하루살이 군이 다른 하루살이속의 집단보다 수질이 비교적 좋지 않은 환경에 대하여까지 적응되어 있음을 확인할 수 있다.

3.2 T-P

T-P에 대한 CDF와 CMF 간의 NRMSE는 Table 5와 같다. 최빈값을 중심으로 정적편포를 보이는 동양하루살이-사할린하루살이 군의 분포에는 굼벨분포모형의 적합도가 가장 높았으며 무늬하루살이 및 가는무늬하루살이의 분포는 최빈값을 중심으로 약한 정적편포를 보이며 모두 대수정규분포모형의 적합도가 가장 높았고 역대수정지수분포모형(reversed exponential distribution model)의 적합도를 제외한 다른 모형들의 적합도 또한 전반적으로 높았다.

Fig. 4는 각 분류군의 유속에 대한 PMF와 이에 대하여 최적의 적합도를 보이는 확률분포모형의 PDF 및 CDF를 도시한 것이다. T-P에 대한 각 분류군의 출현도는 분류군별로 차이가 컸으며 T-P의 범위에 따라 표본단위의 수도 차이가 컸기 때문에 각 급구간의 범위는 분류군의 출현특성에 따라 달리 설정하였다.

Fig. 4. Probability mass function based on the composite frequency of each taxon according to total phosphorus (T-P) gradient, probability density function, and probability distribution function of the best-fit model; (a)Ephemera orientalis-sachalinensis, (b)Ephemera strigata, (c)Ephemera separigata, (d) all taxa.

동양하루살이-사할린하루살이 군은 낮은 농도의 T-P부터 높은 농도의 T-P 까지 전체적으로 출현하였으며 출현도는 최빈값인 0.052 mg/L 이상으로 갈수록 뚜렷하게 감소하였다. 무늬하루살이의 출현도는 최빈값인 0.012 mg/L 이상에서 뚜렷하게 감소하는 경향을 보였다. 가는무늬하루살이의 출현도는 최빈값인 0.000 mg/L 이상에서 감소하는 경향을 보였다(Fig. 4, Table 6).

T-P에 대한 하루살이속의 HSI 값을 산출한 결과, 동양하루살이-사할린하루살이 군의 최적 적응범위는 0.024~0.086 mg/L, 무늬하루살이는 0.007~0.018 mg/L, 가는무늬하루살이는 0.000~0.015 mg/L로 분석되었다(Fig. 5, Table 7). 동양하루살이-사할린하루살이 군의 T-P에 대한 최적 적응범위와 무늬하루살이 및 가는무늬하루살이의 최적 적응범위는 뚜렷하게 구분되었다. 무늬하루살이와 가는무늬하루살이의 최적 적응범위에서는 지위가 중복되지만 가는무늬하루살이가 적응범위가 낮은 농도의 T-P에서 더 크게 나타났다.

Fig. 5. Habitat suitability index ofEphemeraspecies for total phosphorus (T-P) gradient; (a)Ephemera orientalis-sachalinensis, (b)Ephemera strigata, (c)Ephemera separigata, (d) all taxa.

서식처 적합도(Table 7) 및 중앙값과 표준편차(Table 6)를 고려할 때 동양하루살이-사할린하루살이 군은 높은 농도의 T-P까지 좁게 적응되어 있었으며, 무늬하루살이는 동양하루살이-사할린하루살이 군의 적응농도보다 낮은 농도의 T-P에서 좁게 적응되어 있었고 가는무늬하루살이는 낮은 농도의 T-P에서 넓게 적응되어 있는 것으로 보였다.

환경정책기본법에서 T-P의 농도가 0.020 mg/L 이하인 하천을 수질이 매우 좋다고 평가하였으며 0.200 mg/L 이하까지는 보통의 하천이라고 평가하였다. 이에 따르면 동양하루살이-사할린하루살이 군이 다른 하루살이속의 집단보다 수질이 비교적 좋지 않은 환경까지 적응하고 있음을 확인할 수 있었다. 반면에 무늬하루살이와 가는무늬하루살이의 최적 적응범위가 0.020 mg/L 이하였으므로 수질이 매우 좋은 하천에 적응하고 있음을 확인할 수 있었다.

3.3 TSS

TSS에 대한 CDF와 CMF 간의 NRMSE는 Table 8과 같다. 최빈값을 중심으로 높은 정적편포를 보이는 동양하루살이-사할린하루살이 군, 무늬하루살이 및 가는무늬하루살이의 분포는 모두 대수정규분포모형의 적합도가 가장 높았으며 역대수정지수분포모형(reversed exponential distribution model)을 제외한 다른 모형들의 적합도 역시 전반적으로 높았다.

Fig. 6은 각 분류군별로 TSS에 대한 PMF와 이에 대하여 최적의 적합도를 보이는 확률분포모형의 PDF 및 CDF를 도시한 것이다. TSS에 대한 각 분류군의 출현도는 분류군별로 차이가 컸으며 TSS의 범위에 따라 표본단위의 수도 차이가 컸기 때문에 각 급구간의 범위는 분류군의 출현특성에 따라 달리 설정하였다.

Fig. 6. Probability mass function based on the composite frequency of each taxon according to total suspended solids (TSS) gradient, probability density function, and probability distribution function of the best-fit model; (a)Ephemera orientalis-sachalinensis, (b)Ephemera strigata, (c)Ephemera separigata, (d) all taxa.

동양하루살이-사할린하루살이 군은 최빈값인 8.0 mg/L을 전후하여 넓은 구간에서 출현하였으며 142.2 mg/L의 높은 농도까지 출현하였다. 무늬하루살이는 최빈값인 0.0 mg/L 이상의 구간으로 갈수록 출현도가 감소하였다. 가는무늬하루살이는 최빈값인 1.5 mg/L 이상의 구간에서 출현도가 감소하였다(Fig. 6, Table 9).

TSS에 대한 하루살이속의 HSI 값을 산출한 결과, 동양하루살이-사할린하루살이 군의 최적 적응범위는 2.5~22.4, 무늬하루살이는 0.0~1.7, 가는무늬하루살이는 0.5~3.1로 동양하루살이-사할린하루살이 군과 무늬하루살이 및 가는무늬하루살이 분류군은 서로 뚜렷하게 지위가 구분됨을 확인할 수 있었다. 무늬하루살이와 가는무늬하루살이는 지위가 일부 중복됨을 확인할 수 있었다(Fig. 7, Table 10).

Fig. 7. Habitat suitability indexEphemeraspecies for total suspended solids (TSS) gradient; (a)Ephemera orientalis-sachalinensis, (b)Ephemera strigata, (c)Ephemera separigata, (d) all taxa.

서식처적합도(Table 10)와 중앙값과 표준편차(Table 9)를 고려할 때 동양하루살이-사할린하루살이 군은 비교적 낮은 농도의 TSS와 높은 농도의 TSS까지 넓게 적응하였으며, 무늬하루살이는 낮은 농도의 TSS에 비교적 좁게 적응하였고, 가는무늬하루살이는 낮은 농도의 TSS에 무늬하루살이에 비해 넓게 적응된 것으로 평가된다.

환경정책기본법에서 TSS의 농도가 25.0 mg/L 이하인 하천을 수질이 보통 이상이라고 평가하였다. 이는 하루살이속이 모두 보통 이상의 수질의 하천에서 적응하고 있으며 동양하루살이-사할린하루살이 군은 다른 두 종보다 비교적 수질이 안 좋은 환경에서도 서식할 수 있다는 사실을 확인할 수 있었다.