2.1. 자료수집

분석에 활용한 자료는 환경부·국립환경과학원의 “수생태 계 건강성 조사 및 평가(‘08 ~ ‘13년)”와 “하천 수생태계 현 황 조사 및 건강성 평가(‘14년)”, 한강수계관리위원회의 “기 후변화가 수생태계에 미치는 영향과 대응전략(‘10 ~ ‘12)”에 서의 조사결과 중 일부였다. 이중 환경부·국립환경과학원의 자료는 동 사업의 저서성 대형무척추동물 조사지점 중 환 경부 물환경측정망의 수질측정망 지점과 일치하는 지점에 해당하는 것만 선별하였다. 수질측정자료가 수반되지 않는 지점의 생물조사자료는 생물종의 지표성 분석은 물론 생물 지수와 수질과의 상관성 분석에도 활용할 수 없으므로 분 석에서 제외하였다. 수질과 생물이 조사된 지점에서의 연중 또는 연간 반복된 조사를 포괄한 총 표본단위(sampling unit)는 7,086개였다.

2.2. 생물출현도의 확률분포

하천의 유기오염(organic pollution)을 대표하는 주요 항목 인 5일 BOD (5-day biochemical oxygen demand, BOD₅) 농도를 변수로 한 생물종의 출현 특성을 확률적으로 분석하 였다. BOD₅ 농도를 x, BOD₅ 농도의 급구간(class interval) 에 대한 일련번호를 i라 할 때, i 급구간의 평균 출현개체 수로 표현한 확률질량함수(probability mass function)는 식 (1)과 같다. 본 연구에서는 이 식에서 도출된 함수 값을 이 용하여 생물종이 서식하는 BOD₅ 농도의 중심 및 산포와 관련된 통계량(statistic)을 산출하였다.

PMF_i : BOD₅ 농도 i 급구간의 확률질량함수(probability mass function)

A _i : BOD₅ 농도 i 급구간의 평균 출현개체수(mean number of individuals)

Δx_i : BOD₅ 농도 i 급구간(class interval)

환경요인의 구배에 따라 생물종이 대수정규분포를 하고, 대수적인 환경요인의 범위가 생물종의 분포에 등가의 영향 을 미친다고 가정할 경우에는 생물종의 분포해석에 있어 환경요인의 실제 값보다는 대수 값을 사용하는 것이 더욱 적합할 수 있다. 오수생물계열(saprobic series)에 대응하는 환경요인의 범위가 이에 해당할 수 있다. 본 연구에서는 BOD₅ 농도의 2진로그(binary logarithm) 값을 사용한 급구간 으로 오수생물계열별 생물종의 분포를 검토하였다(식 (2)).

λ : BOD₅ 농도의 2진로그(binary logarithm) 값 ( λ = log 2 x )

Δλ_i : BOD₅ 농도의 2진로그 값의 i 급구간

이때 Δx_i가 x_i+1/x_i = 2인 2배수 급구간(octave)이라면, Δ λ i = λ i + 1 − λ i = log 2 x i + 1 − log 2 x = log 2 x i + 1 / x i = log 2 2 = 1 이 되기 때문 에, 생물종의 출현도에 대한 확률질량함수는 특정 급구간의 평균 출현개체수(A_i)를 각 급구간별 평균 출현개체수를 모 두 합한 값(ΣA_i)으로 나눈 비율이 된다(식 (3)).

식 (1)과 식 (2)에서 급구간을 미분한 확률밀도함수 (probability density function, PDF)를 각각 f(x)와 f(λ)라 하고 각각 대수정규분포(lognormal distribution)와 정규분 포(normal distribution)를 한다고 가정하면 이에 해당하는 확률분포식과 모수(parameter)는 Table 1과 같다. 확률분포 식의 모수 값은 i 급구간에 대한 누적이산확률(cumulative discrete probability) 값과 누적분포함수(cumulative distribution function, CDF) 값의 평균제곱근오차(root mean square error, RMSE)가 최소가 되도록 시행착오를 통해 도출하였 다(식 (4)).

P_i: i 급구간의 이산확률(PMF_iΔx_i 또는 PMF_iΔλ_i)

F_i: i 급구간의 누적분포함수(cumulative distribution function)

n : 급구간의 총수

2.3. 저서성 대형무척추동물의 지표성

2.4. 저서동물지수(Benthic Macroinvertebrates Index, BMI) 의 구성

현재까지 우리나라에서 사용되어 온 오수생물지수 (saprobic index)의 형태는 ^{Zelinka and Marvan (1961)}에 따 라 생물종의 출현도와 지표가중치(indicator weight value)로 가중하여 오탁지수를 평균한 것이다. 즉 ^{Yoon et al. (1992a)}이 ^{Zelinka and Marvan (1961)}의 오수생물지수를 도입하여 Yoon-Kong의 SI를 제안한 이후 ^{Won et al. (2006)}의 KSI나 ^{Kong et al. (2012)}의 KBI가 모두 같은 유 형의 오수생물지수에 근간을 둔 것이었다. 본 연구에서 제 안하는 저서동물지수(BMI) 역시 ^{Zelinka and Marvan (1961)}의 오수생물지수 값을 단지 가독성을 높이기 위해 0 ~ 100점으로 계량화한 것일 뿐 그 개념에는 차이가 없다 (식 (14)).

BMI : 저서동물지수 (Benthic Macroinvertebrates Index)

i : 지정된 지표생물종의 일련번호 (the number assigned to the species)

n : 출현한 지표생물종의 총수 (the number of species)

s_i : i 종의 오탁지수 (the saprobic value of the species i)

h_i : i 종의 상대개체수출현도 (the relative abundance of the species i)

g_i : i 종의 지표가중치 (the indicator weight value of the species i)

^{Yoon et al. (1992a)}, ^{Won et al. (2006)}, ^{Zelinka and Marvan (1961)}은 오수생물지수 산정 시 개체수 출현도(h_i) 에 출현 개체수(number of individuals)를 적용한 바 있다. 그러나 이 경우 특정 생물의 생활사에서 부화시기에는 많 은 개체수가 출현하기 때문에 같은 환경상태라도 평가치가 이들 종의 특성으로 크게 왜곡되는 사례가 발생하며, 또한 특정 종이 매우 많이 출현할 경우 다른 많은 희소종(rare species)들의 출현특성은 무시된다는 단점을 가지고 있다 (^{Kong et al., 2012}). 이와 관련하여 ^{Kong et al. (2012)}은 출 현개체수를 적용하는 것보다 상대출현도(relative abundance) 를 적용하는 것이 생물학적 하천환경평가에 대한 변별력이 더욱 크다는 것을 확인한 바 있다. 따라서 본 연구에서도 ^{Kong et al. (2012)}에 따라 출현개체수 대신 상대출현도를 적용하였다.

독일의 ^{DIN 38410 (1990)}의 오수생물지수(saprobic index) 산정 시 출현도를 7단계로 구분하여 적용하고 있으나, 본 연구에서는 출현개체수의 순위를 ^{Zelinka and Marvan (1961)}과 같이 5단계로 구분하였다. 각 종의 출현개체수 순 위를 백분위수로 변환하고(식 (15)), 이를 5 등급으로 구분 하여 상대출현도 값을 1부터 5까지 부여하였다(Table 5).

R_i : i 지표생물종의 출현개체수 순위 백분율 (the rank percentage of species i based on individual abundance)

N_i : i 지표생물종의 출현개체수 순위 (the rank of species i based on individual abundance)

N : 총 출현종수 (the total number of species)

2.5. 통계분석

BMI의 평가등급별 BOD₅ 농도의 분포는 대수정규분포 (lognormal distribution)로, BMI의 평가등급별 Shannon-Weaver 다양도지수 H'의 분포는 와이블모형으로 해석하였 다. 와이블모형은 금속 및 복합재료의 강도나 전자 및 기 계부품의 수명분포를 나타내는 데 활용되고 있으나, 생물종 의 서식처합도 분석(^{Kong and Kim, 2017}) 등 생태학적 연 구에서도 활용되고 있다.

3.1. 생물출현도의 확률분포

하루살이류 중 출현도가 높은 Uracanthella를 대표 분류 군으로 선별하여 BOD₅ 농도에 따른 출현 특성을 분석하 였다. 식 (1)의 BOD₅ 농도의 급구간별 Uracanthella의 평 균 출현개체수를 가지고 확률질량함수를 도출한 결과 이 에 부합하는 확률밀도함수는 대수정규분포로 나타났다(Fig. 1(a), 1(b)). 즉 Uracanthella의 출현도는 BOD₅ 농도에 따라 3.6의 왜도(skewness)를 보이며 정적편포(positive skewed distribution)함에 따라 중심치 중 평균이 2.1 mg L^-1로 최빈 수인 0.8 mg L^-1나 중위수인 1.6 mg L^-1보다 크게 나타나고 있다. 즉 BOD₅ 농도의 절대농도를 기준으로 한 생물종 (Uracanthella를 대표 종으로 볼 때)의 분포는 편포하는 형 태를 보이기 때문에 Table 2와 같이 배수체계의 오수생물 계열을 적용하는 경우에는 정규분포로의 변형에 따른 분포 해석이 필요함을 시사하는 것이다.

Fig. 1. Relationship between the mean individual abundance ofUracanthellaand the BOD5concentration (a) PDF of abundance vs. BOD5(b) CDF of abundance vs. BOD5(c) PDF of abundance vs. octaves of BOD5(d) PDF of abundance vs. tail-tied octaves of BOD5; Noteλ=log2x

식 (2)의 2진로그 BOD₅ 값의 급구간을 2배수 체계(octave) 로 하고 급구간별 Uracanthella의 평균 출현개체수를 가지 고 확률질량함수를 도출한 결과 이에 부합하는 확률밀도함 수는 정규분포로 나타났다(Fig 1(c)). 이때 BOD₅ 농도가 1 mg L^-1 이하인 범위를 극빈부수성으로, 8 mg L^-1를 초과하 는 범위를 강부수성으로 통합하여 평균개체수를 적용하면 식 (2)는 식 (10)과 같은 형태가 된다. 즉 BOD₅ 농도의 2 진로그 값이 0이하인 곳과 4 초과인 것을 통합(tail-tied)하 여, 식 (5)로부터 각 BOD₅의 급구간별 오탁계량치를 구하 고, 이를 식 (11)에 대입하여 구해진 Uracanthella의 오탁지 수는 0.8이었다.

여기에서 Fig. 1은 BOD₅ 농도에 따른 생물종의 출현 특 성을 바탕으로 오탁지수를 도출하는 과정을 개념화한 것일 뿐이다. 실제 오탁지수의 산정은 모든 종이 Uracanthella와 같이 오수생물계열에 따라 정규분포하는 것으로 가정하고 식 (10)과 식 (11)에 따라 산출되었고 그 결과는 Appendix 1과 같다. 대부분의 분류군들은 오탁계량치가 특정 오수생 물계열에 단일최빈수(unimodal)의 경향을 보이고 있으나 Phagocata와 같이 일부 분류군에서는 오탁계량치의 정규성 이 나타나지 않았다. 그 이유는 아직 분류학적 완성도가 낮아 한 분류군에 생태적 특성이 다른 미지의 여러 종이 포함되어 있는 경우, 오동정 또는 다른 오류로 자료가 혼 돈된 경우, 표본크기(sample size)가 작아 아직 생태학적 정 보가 미흡한 경우 등으로 구분할 수 있다. 이는 향후 지속 적으로 정보를 축적하고 갱신하면서 보완해야 할 사항이다.

3.2. 오탁계량치(saprobic valency)

오수생물계열별 오탁계량치의 분포는 출현속성별로 다르 게 나타났다. 전체적으로 출현도가 높았던 파리류의 깔다구 과(Chironomidae)를 예로 들 때 실평균개체수를 기준으로 산정된 오탁계량치의 계열 간 차이가 가장 작았다(Fig. 2(b)). 평균출현개체수를 기준으로 산정된 오탁계량치의 분 포(Fig. 2(a))는 상대출현빈도를 기준으로 하였을 때(Fig. 2(c))에 비해 계열간 차이가 비슷한 수준이었다. 복합상대 출현도를 기준으로 산정된 오탁계량치의 계열 간 차이는 전체 출현속성 중에 가장 크게 나타났다(Fig. 2(d)).

Fig. 2. Saprobic valency of 190 indicator taxa according to each saprobic series value. Error bars indicate standard deviation.

어떤 종이 특정 계열을 선호하는 경우 상대출현빈도는 물론 평균출현개체수도 높아지는 경우가 많으므로 이를 곱 한 속성인 복합상대출현도에 의한 오탁계량치는 해당 계열 에서 상대적으로 더욱 커지는 것으로 이해된다(^{Kong and Kim, 2016}).

3.3. 저서동물지수(BMI)

3.3.1. 평균오탁지수(Average saprobic value)의 적용

네 가지 출현속성별로 평균오탁지수를 기반으로 산정된 BMI 값과 BOD₅ 농도와의 관계는 Fig. 3과 같다.

Fig. 3. Relationship between the concentration of BOD5and Benthic Macroinvertebrates Index (BMI) using average saprobic values.

양 극단의 오수생물계열을 선호하는 종일지라도 해당 종 은 일반적으로 여러 계열에 걸쳐 출현할 수 있기 때문에 오탁계량치로 가중평균한 오탁지수는 오수생물계열의 양 극단치(0 또는 4)보다는 중심쪽으로 치우치는 경향을 보이 게 된다. 특히 실평균개체수를 기준으로 산정된 평균오탁지 수를 적용한 BMI 값은 BOD₅ 농도에 따른 변화도가 상대 적으로 적은데, 이는 오수생물계열에 따른 오탁계량치가 큰 차이를 보이지 않은 점(Fig. 2(b))과 관련이 있다.

이러한 이유로 평균오탁지수를 기반으로 산정된 BMI 값 은 BOD₅ 농도에 따른 변화율이 완만하여 BOD₅ 농도가 매 우 낮거나 높은 구간에서 나타나는 값은 식 (14)의 극단 값(100 또는 0)보다는 중심치로 치우치는 경향을 보인다. 복합상대출현도는 이러한 중심화 경향이 가장 작게 나타나 도록 하는 변수로서 이를 기반으로 한 BMI 값은 다른 변 수의 BOD₅ 농도와의 상관성이 약간 낮지만 그 변별력은 가장 크다고 할 수 있다.

3.3.2. 대표오탁지수(Representative integer of saprobic value) 의 적용

출현속성별로 대표오탁지수를 기반으로 산정된 BMI 값 은 평균오탁지수를 적용하였을 때에 비해 BOD₅ 농도에 대한 상관도가 약간 낮았다(Fig. 4). 대표오탁지수는 출현 도가 가장 큰 오수생물계열치인 정수(integer)가 되기 때문 에 오염에 대한 분포범위가 넓어도 그 지표성이 특정치에 제한된다는 점에서 수질판정의 변별력이 떨어지는 것으로 보인다.

Fig. 4. Relationship between the concentration of BOD5and Benthic Macroinvertebrates Index (BMI) using representative integer saprobic values.

한 예로서 Appendix 1에서 4번 Cipangopaludina (논우렁 이속)는 오수생물계열의 모든 계열에서 출현하며 특히 오 염된 곳에 치우치는 분포를 보이지만 대표치만 고려하면 실제는 내성도가 높은 종임에도 불구하고 빈수부성 분류군 으로 처리되는데 이를 적용할 경수 불합리한 평가결과를 초래할 수 있다.

복합상대출현도를 기준으로 산정된 평균오탁지수를 적용 한 BMI 값은 BOD₅ 농도에 따른 변별력과 상관도 양 측면 에서 적정하였으므로(Fig. 3(d)), 이를 기준으로 저서성 대 형무척추동물 총 190분류군의 오탁계량치(V)를 산정하고, 평균오탁지수(s)와 지표가중치(g)를 결정하였다(Appendix 1). BMI 값의 계산은 종 수준을 기반으로 한 것이므로, 오탁 지수와 지표가중치가 속(genus) 수준 이상의 상위범주로 설 정된 경우에는 각 종의 오탁지수와 지표가중치는 해당 상 위 범주의 값을 따라야 한다.

3.4. 기존 KSI와 BMI 결과 비교

복합상대출현도를 기준으로 산정된 BMI 값(Fig. 5(a))은 기존의 ^{Won et al. (2006)}의 KSI (Fig, 5(b))에 비해 BOD₅ 농도에 대한 상관도가 약간 더 높았다. 또한 BOD₅가 고농 도인 구간에서 KSI는 BMI에 비해 산포 정도가 크게 나타 나 고농도의 BOD₅ 농도로 대표되는 유기오염에 대한 생물 학적 평가에는 적용성이 상대적으로 낮은 경향을 보였다.

Fig. 5. Relationship between (a) Benthic Macroinvertebrates Index (BMI) (b) Korean Saprobic Index (KSI) and the concentration of BOD5.

3.5. 저서동물지수(BMI)의 수생태계 건강성 평가 등급

복합상대출현도 기준의 평균오탁지수에 따른 BMI 값과 BOD₅, 총부유물질(TSS), 총인(TP)의 장기 평균농도와의 다 중회귀식은 식 (16)과 같이 나타났다. 이때 BOD₅는 유기오 염을, 총부유물질은 무기오염을, 총인은 영양상태(trophic state)를 각각 대변하는 수질항목으로 고려한 것이다.

BOD ₅ : the concentration of 5-day biochemical oxygen demand (mg L^-1)

TSS : the concentration of total suspended solods (mg L^-1)

TP : the concentration of total phosphorus (mg L^-1)

BMI의 평가등급은 Table 3의 5계의 오수생물계열 체계 와 2배수의 BOD₅ 농도 기준으로 구분하였다(Table 6). 우 리나라의 하천 생활환경기준에서 총인 농도의 등급치는 일 반적인 부영양화 기준에 비해 상당히 완화되어 있는 수준 이므로 본 연구에서는 ^{OECD (1982)}의 영양상태 판정기준 (<0.01 mg L^-1: oligotrophic, 0.01 ~ 0.035 mg L^-1: mesotrophic, 0.035 ~ 0.1 mg L^-1: eutrophic)과 하천의 생활환경기준의 ‘약간 나쁨’ 등급의 경계치인 0.3 mg L^-1를 적용하여 구분 하였다.

하천 생활환경기준에서 부유물질 농도의 기준은 III 등급 까지도 변별력이 없는 25 mg L^-1의 단일 기준으로 되어 있 다. 부유물질 농도가 높고 이에 따라 하상이 세립화되면 굴파는 무리를 제외한 대부분의 저서생물이 먹이활동, 산 란, 서식에 방해를 받기 때문에 부유물질은 하천의 다른 어떤 생물보다도 저서생물에 큰 영향을 준다. 따라서 평가 등급별 부유물질 농도의 기준은 하천과 호소의 생활환경기 준에서 ‘보통’ 상태에 대한 기준치의 중간 값인 20 mg L^-1 에서 배수체계로 감하여 정하였다.

각 등급의 BMI 값은 각 수질항목의 등급치를 식 (15)의 다중회귀식에 대입하여 계산치를 도출하고 이 값을 가독성 을 높이기 위해 5점 단위로 조정하여 설정하였다.

Table 7은 이러한 일련의 과정을 통해 도출된 BMI의 평 가등급과 각 등급에서 나타나는 저서성 대형무척추동물 군 집의 다양성과 민감 생물군의 출현도에 대한 서술적 상태 를 제시한 것이다.

Fig. 6은 Table 7의 기준에 따라 7,086개 표본단위에서 도출된 BMI의 각 등급별 BOD₅ 농도의 분포와 다양도지수 의 분포를 도시한 것이다. BMI의 등급별 BOD₅ 농도의 분 포는 대수정규분포의 형태로 나타났고(Fig. 5(a)), 다양도지 수의 분포는 등급별로 편포(skewed)하는 방향과 정도가 달 랐다(Fig. 5(b)).

Fig. 6. Distribution of (a) BOD5concentration and (b) Shannon-Weaver's H' according to each class of BMI in the 7,086 sampling units during 2008 ~ 2014.

BOD₅ 농도의 분포는 A 등급으로부터 하위등급으로 갈수 록 최빈수가 고농도 쪽으로 이동하고 첨도가 낮아지며 분 산도가 커졌다. 등급별 90 백분위수는 A가 1.2 mg L^-1, B 가 2.4 mg L^-1, C가 4.1 mg L^-1, D가 5.4 mg L^-1, E가 7.2 mg L^-1였다. 예를 들어 BMI의 평가결과가 A 등급이면 수 질의 약 90%가 I등급(매우 좋음) 수준으로 볼 수 있다.

BOD₅ 농도의 분포는 등급 간에 중첩되는 정도가 컸으나 다양도지수는 ESB 평가등급을 보다 명확하게 구분하는 것 으로 나타났다. 다양도지수의 분포는 상위등급으로 갈수록 최빈수가 고다양도 쪽으로 이동하고 첨도가 커지며 분산도 가 낮아졌다. 등급별 최빈수(90% 신뢰구간)는 A가 3.3(2.1 ~ 4.1), B가 2.9(1.5 ~ 3.9), C가 2.0(0.9 ~ 3.4), D가 1.6(0.6 ~ 2.9), E가 1.2(0.4 ~ 2.4)이었다. 다양도지수의 최빈수를 기준 으로 한다면 A 등급은 ^{Staub et al. (1970)}의 기준에 따라 빈부수성, B와 C는 베타중부수성(각각 빈부수성과 알파중 부수성에 가까운), D는 알파중부수성, E는 강부수성에 가 까운 알파중부수성이라 할 수 있다.

3.6. 향후 연구방향

본 연구에서 제시된 지표생물표의 190분류군은 종 수준 뿐만 아니라 속 이상의 상위분류군에서 설정된 분류군이 많기 때문에 일반적으로 하천에서 출현하는 대다수의 종이 이에 포함될 수 있다. 그러나 본 연구의 주요 대상 수역은 산간계류부터 하구에 이르는 유수역(lotic)을 대상으로 하고 있으므로 정수역(lentic)에 주로 서식하는 많은 종류들은 지 표생물군의 목록에 포함되어 있지 않다. 향후 정수역 생물 의 환경지표성에 대한 연구가 수행되어 정수역의 생물학적 환경평가를 위한 방법론을 정립할 필요가 있다.

최근 수생태계에 대한 미국과 독일의 평가체계는 생물군 집의 온전성(biological integrity)에 관한 평가요소(elements) 들을 계량화한 다중계량지수(multimetric index)에 초점을 두고 있다(^{ME, 2014}). 영국은 자연적인 상태에서 나타날 수 있는 모델기반의 BMWP score와 실제의 값을 비교하는 체계를 적용하고 있다(^{ME, 2014}). 이들 국가의 수생태계 평가체계에 대한 상세한 사항은 ^{ME (2014)}에서 파악할 수 있다. 이러한 제외국의 흐름은 과거의 단순한 생물학적 수 질평가에 그치지 않고 군집내 종 조성, 길드(guild) 조성, 풍부도 등을 총체적으로 고려한다는 점이다.

우리나라 국립환경과학원의 ‘생물측정망 조사 및 평가지 침’에서 어류의 평가지수는 이미 다중계량지수에 기초한 것이지만 부착돌말과 저서동물지수는 각각 영양상태와 부 수성을 평가하는 지수이다. 향후 저서성 대형무척추동물에 의한 수생태계 평가체계는 제외국의 흐름에 맞추어 온전성 평가체계로 전환해야 할 것으로 본다. 그러나 이 경우에도 민감종(intolerant species)의 구성에 바탕을 둔 지수나 지표 는 여전히 평가요소의 하나로 다루어져야 하기 때문에 현 재의 저서동물지수의 신뢰성을 높이기 위한 조사자료의 축 적 및 관련 연구는 지속되어야 한다.

저서동물지수의 신뢰성을 높이기 위하여 우선적으로 필 요한 사항은 환대강(Clitellata)이나 깔다구과(Chironomidae) 와 같이 동정(identification)의 문제로 상위 범주로 묶어 지 표치를 부여하는 분류군에 대하여 그간의 분류학적 성과를 평가체계에 흡수하여 활용하는 것이다. 또한 위에서 언급된 바와 같이 오수생물계열에 따른 오탁계량치가 단일최빈수 형(unimodal)이 아닌 분류군에 대해서는 향후 보다 많은 자 료를 축적하여 지표치를 보완할 필요가 있다. 끝으로 본 연구에서 제시되고 있는 190 분류군의 지표치는 주로 유수 (lotic) 환경에서 나타난 결과로부터 도출된 것이다. 정수 (lentic) 환경에서 주로 출현하는 종들에 대한 지표치에 대 한 설정은 여전히 과제로 남아있는 상태이다.