王文輝，李 艷，雷 坤，張 萌，魏 明

（1.中國(guó)環(huán)境科學(xué)研究院，北京 100012；2.沈陽(yáng)理工大學(xué)，沈陽(yáng) 110159；3.復(fù)旦大學(xué)大氣科學(xué)研究院，上海 200438；4.中國(guó)科學(xué)院西北生態(tài)環(huán)境資源研究院，蘭州 730000）

0 引 言

水質(zhì)監(jiān)測(cè)通過(guò)對(duì)水污染物含量的定量分析，客觀評(píng)價(jià)一個(gè)地區(qū)水污染程度，為水污染防治和水環(huán)境管理提供決策服務(wù)［1］。水體中TN 的豐歉與水體是否富營(yíng)養(yǎng)化、水生生物能否平衡密切相關(guān)，是衡量水環(huán)境質(zhì)量的重要指標(biāo)［2］，長(zhǎng)期以來(lái)，生活污水、含氮工業(yè)廢水的大量排放及耕地化肥的施用，引起地表水水體TN 含量升高，導(dǎo)致溶解氧（Dissolved Oxygen，DO）濃度急劇下降，水生生物平衡被破壞，水環(huán)境質(zhì)量一度下降到無(wú)法滿足飲用和灌溉的需求。遼河流域作為我國(guó)最大工業(yè)基地所在地，接納了遼寧省境內(nèi)大多數(shù)冶金、石化、藥物和印染等重污染行業(yè)的廢水，水污染狀況嚴(yán)峻。自1997年以來(lái)，遼河流域的水污染治理經(jīng)歷了從單一的控源截污到流域的污染綜合治理［3］，“十一五”、“十二五”期間遼河流域以氨氮（NH3-N）作為水污染防治的重點(diǎn)關(guān)注指標(biāo)之一，參與水環(huán)境質(zhì)量的監(jiān)測(cè)與評(píng)價(jià)，污染程度有所緩解，但通過(guò)水質(zhì)監(jiān)測(cè)數(shù)據(jù)來(lái)看，水質(zhì)達(dá)標(biāo)狀況不穩(wěn)定，水污染控制程度仍不樂(lè)觀：2013-2018年五年間，劣Ⅴ類斷面占比增長(zhǎng)16.7%，2019年又降低11.8%?！笆濉庇痔岢鲞M(jìn)一步的要求，將TN也列入水環(huán)境監(jiān)測(cè)的重要指標(biāo)。

對(duì)水污染物的監(jiān)測(cè)通常以人工現(xiàn)場(chǎng)采樣并進(jìn)行實(shí)驗(yàn)分析的方式，或在線自動(dòng)分析儀器來(lái)完成［4］。近年來(lái)隨著空間探測(cè)技術(shù)的革新，水環(huán)境遙感技術(shù)相比于傳統(tǒng)的手段，具備監(jiān)測(cè)范圍廣、速度快、成本低及便于進(jìn)行動(dòng)態(tài)監(jiān)測(cè)的優(yōu)勢(shì)［5］，而且還能在一定程度上解決歷史數(shù)據(jù)難獲取的問(wèn)題，因此受到眾多研究人員的青睞。早期水環(huán)境遙感在研究區(qū)域上較多地居于海洋、湖庫(kù)等寬闊水域，隨著傳感器成像空間分辨率的提高，逐漸拓展到小面積水域的河流；研究?jī)?nèi)容從水域的定性識(shí)別［6］發(fā)展到水質(zhì)參數(shù)的定量估計(jì)［7，8］；研究方法經(jīng)歷了從分析法［9，10］到經(jīng)驗(yàn)法［11］，再到半經(jīng)驗(yàn)-半分析法［12，13］的耦合。

對(duì)TN 遙感監(jiān)測(cè)的研究以經(jīng)驗(yàn)法和半分析-半經(jīng)驗(yàn)法為主。經(jīng)驗(yàn)法類似于“黑箱”，在一定程度上忽略光在大氣與水體之間的傳輸過(guò)程，基于TN 濃度與遙感反射率之間的相關(guān)關(guān)系，以統(tǒng)計(jì)學(xué)的方法完成對(duì)TN 濃度的反演，無(wú)需大量的長(zhǎng)時(shí)序?qū)崪y(cè)數(shù)據(jù)，具有易操作的性質(zhì)，這種方法的核心在于遙感數(shù)據(jù)源與敏感特征選取，不乏有大量學(xué)者對(duì)此進(jìn)行深入的研究，如Torbick等［14］利用TM 影像對(duì)密歇根州42 個(gè)湖泊的TN 濃度反演取得了較好的結(jié)果；雷坤等［15］利用另一種遙感數(shù)據(jù)源（CBERS-1-CCD），建立了太湖表層TN 濃度的反演模型（MRE=12.2%）；何報(bào)寅等［16］、楊柳等［17］則以ETM+影像，分別反演武漢東湖和溫榆河TN 濃度，取得了較好的結(jié)果；而王學(xué)軍等［18］通過(guò)單波段、多波段組合和主成分分析的方式選擇敏感特征，基于TM 影像反演了太湖TN 濃度。隨著高光譜技術(shù)的發(fā)展，研究人員又嘗試?yán)没趯?shí)測(cè)光譜的半分析-半經(jīng)驗(yàn)法，以期取得更準(zhǔn)確的TN反演結(jié)果，如鞏彩蘭等［19］通過(guò)分析實(shí)測(cè)高光譜歸一化值、波段比值和微分值與TN 濃度的相關(guān)性，以656.53 nm 反射率歸一化值與880 nm 反射率歸一化一階微分值反演了黃浦江TN 濃度；徐良將等［20］同樣以實(shí)測(cè)光譜，通過(guò)微分和波段比值的方式，發(fā)現(xiàn)455 nm 反射率微分值與1 015/528 nm 反射率比值為TN 濃度最敏感變量，并得到了R2=0.839 的TN 反演模型。近幾年面對(duì)TN反演準(zhǔn)確度較低的問(wèn)題，以機(jī)器學(xué)習(xí)為基礎(chǔ)的非參數(shù)算法無(wú)需確定目標(biāo)函數(shù)的形式，具有更高的靈活性與泛化性，表現(xiàn)略勝一籌，王建平等［21］以人工神經(jīng)網(wǎng)絡(luò)模型對(duì)鄱陽(yáng)湖TN 濃度的反演研究、王云霞［22］以最小二乘支持向量機(jī)法對(duì)遼寧清河水庫(kù)的TN反演，都取得了較好的結(jié)果。

大流域尺度河流水體水質(zhì)遙感監(jiān)測(cè)受遙感成像空間分辨率、成像幅寬等諸多因素的制約，有待進(jìn)一步認(rèn)識(shí)。在充分了解遙感水質(zhì)監(jiān)測(cè)優(yōu)勢(shì)和發(fā)展現(xiàn)狀的前提下，針對(duì)遼河流域干流及部分主要支流，采用最大坡度下降算法選取純水體像元，利用GF1-WFV 遙感影像反演TN 濃度，及時(shí)掌握TN 在流域尺度上的時(shí)間、空間變化狀況，為大流域尺度河流水體的TN 反演研究積累一定的經(jīng)驗(yàn)參考和數(shù)據(jù)，為遼河流域水生態(tài)環(huán)境精細(xì)化管理與保護(hù)提供科學(xué)依據(jù)。

1 研究區(qū)域概況與數(shù)據(jù)采集

1.1 研究區(qū)域概況

遼河流域（40°31＇～45°17＇N，116°54＇～125°32＇E）位于遼寧省境內(nèi)，流域面積約21.9 萬(wàn)km2，流域內(nèi)有遼河、渾河和太子河三大水系。流域覆蓋遼寧省鐵嶺、沈陽(yáng)、撫順等8 個(gè)地市全境，及黑山、北鎮(zhèn)和彰武等4縣（市）。由于地處東北老工業(yè)基地，長(zhǎng)期以來(lái)形成了以煤炭與石油開(kāi)采、焦化和鋼鐵等重工業(yè)為主的工業(yè)結(jié)構(gòu)，導(dǎo)致污染排放強(qiáng)度高，致使遼河流域污染狀況愈加嚴(yán)峻。自“九五”被列為國(guó)家重點(diǎn)治理的流域至今，遼河流域干流河段COD 劣V 類已基本被消除，部分區(qū)域生境有所恢復(fù)，水質(zhì)退化的趨勢(shì)已經(jīng)得到基本的控制，但水質(zhì)達(dá)標(biāo)狀況依舊不穩(wěn)定：2016年考核斷面水質(zhì)達(dá)標(biāo)率為80%，2017年同比增長(zhǎng)3.7%，到2018年又同比下降27.4%；水質(zhì)優(yōu)良比例（Ⅰ～Ⅲ類）從2016年的32%增長(zhǎng)到2017年的46.9%，2018年又下降了7.3%；喪失水體使用功能（劣V類）占比總體呈上升趨勢(shì)。

1.2 數(shù)據(jù)采集與預(yù)處理

1.2.1 水質(zhì)數(shù)據(jù)

選取流域面積在1 000 km2以上的支流，以谷歌地圖（Level=14，空間分辨率約為17 m）為底圖，采用人機(jī)交互的方式在支流入干前，選取水陸區(qū)分度較高、水面寬度大于1 個(gè)像元的區(qū)域，布設(shè)水樣采集點(diǎn)位；干流水樣采集點(diǎn)位的布設(shè)方法與支流相同，盡可能上中下游均布設(shè)點(diǎn)位；布設(shè)的點(diǎn)位盡量與日常監(jiān)測(cè)點(diǎn)位保持一致，以便對(duì)數(shù)據(jù)進(jìn)行補(bǔ)充。全流域布設(shè)40個(gè)水質(zhì)采樣點(diǎn)位（圖1）。于2018年4-9月，每月4-10 號(hào)采集水樣并分析水質(zhì)。

圖1 采樣點(diǎn)位布設(shè)Fig.1 Layout of sampling points

1.2.2 遙感數(shù)據(jù)

對(duì)比非商業(yè)多光譜遙感數(shù)據(jù)之間的參數(shù)，綜合考慮遙感衛(wèi)星成像的時(shí)間、空間分辨率，及成像幅寬和數(shù)據(jù)的批量易獲取性：高分一號(hào)WFV（GF1-WFV）在星下點(diǎn)成像空間分辨率（16 m）上優(yōu)于Landsat 系列（30 m，Landsat8 融合后為15 m）、MODIS（最高250 m）和HJ1A、1B（30 m），在成像幅寬（800 km）、成像時(shí)間分辨率（2 d）和數(shù)據(jù)批量易獲取性上要優(yōu)于Sentinel-2，多光譜相機(jī)成像譜段涵蓋藍(lán)光（0.45～0.52 μm）、綠光（0.52～0.59 μm）、紅光（0.63～0.69 μm）和近紅外光（0.77～0.89 μm），滿足研究的需求。選用GF1-WFV 數(shù)據(jù)作為研究遼河流域干流及主要支流TN濃度反演的遙感數(shù)據(jù)源。

從資源衛(wèi)星應(yīng)用中心網(wǎng)站（http：//www.cresda.com/CN/）篩選并獲取云量小于10%的研究區(qū)域GF1-WFV遙感影像，7月全月云量較多，不能滿足研究需求，因此保留了4-6月、8-9月的影像，共計(jì)25 景，成像時(shí)間與水樣采集時(shí)間準(zhǔn)同步（±5 d）。依次按照輻射定標(biāo)、大氣校正、正射校正的順序，基于IDL語(yǔ)言，集成了ENVI5.3 下的輻射定標(biāo)、FLAASH 大氣校正和正射校正模塊，對(duì)GF1-WFV 影像形成流程化的批量處理程序，其中正射校正所使用的DEM數(shù)據(jù)為ASTER GDEM。

2 研究方法

2.1 數(shù)據(jù)分析

2.1.1 TN濃度潛在特征構(gòu)造

GF1-WFV 影像帶寬高，波段少，但通過(guò)波段之間的數(shù)學(xué)運(yùn)算可以增強(qiáng)或屏蔽某些信號(hào)［23］，使有用的信號(hào)更為突出。在GF1-WFV 四個(gè)波段的基礎(chǔ)上，參考相關(guān)研究［24］，以雙波段、三波段、四波段組合的方式構(gòu)建了59個(gè)TN濃度潛在特征，累計(jì)63個(gè)潛在特征，計(jì)算方式見(jiàn)表1。

表1 TN濃度潛在特征Tab.1 Latent features of TN concentration

2.1.2 純水體像元與水體提取

水體的反射率通常小于其他地面物體，純水體像元（遙感影像中只包含水體的像元，以GF1-WFV 為例，1 個(gè)像元為邊長(zhǎng)16 m 的正方形）的歸一化植被指數(shù)（NDVI）值一般不大于0，為了確保最大限度地降低偶然誤差，采取最大坡降算法［25］，以采樣點(diǎn)的像元為中心，固定與其相鄰的8個(gè)像元，計(jì)算中心像元與臨近像元的坡降度［式（1）］，選取坡降最大的臨近像元作為采樣點(diǎn)的純水體像元。

式中：Drop為中心像元與臨近像元的坡降度；NDVIcenter為中心像元的NDVI值；NDVIε為 與中心像元臨近像元的NDVI值，ε=1，2，3，…，8；L為中心像元與臨近像元的距離，對(duì)角線臨近像元，其余臨近像元L=16。

GF1-WFV 不具備中紅外波段，無(wú)法計(jì)算MNDWI，故以NDWI計(jì)算結(jié)果結(jié)合人機(jī)交互的方式進(jìn)行水體的提取，GF1-WFV影像NDWI計(jì)算見(jiàn)式（2）。

式中：Rrs(B2)為綠波段遙感反射率，GF1-WFV 為第2 波段；Rrs(B4)為近紅外波段遙感反射率，GF1-WFV為第4波段。

2.2 參數(shù)回歸模型

參數(shù)回歸模型對(duì)自變量和因變量之間有明確的函數(shù)式，通過(guò)分析潛在特征與TN 濃度的相關(guān)性，確定TN 濃度的敏感特征，以TN 濃度作為因變量，敏感特征作為自變量，建立一個(gè)相關(guān)性較高、且具有明確函數(shù)關(guān)系的方程，實(shí)現(xiàn)對(duì)TN 濃度的反演。但高維度的可能存在冗余特征，將這些特征全部引入模型，雖然擬合效果較理想，但冗余特征的共線性會(huì)過(guò)度稀釋有用特征的信息量，使目標(biāo)值的特征難以被解釋。逐步回歸的方法將特征逐個(gè)引入模型，根據(jù)特征對(duì)目標(biāo)值的解釋度剔除掉冗余特征，從而保留解釋度較高的特征。按照8∶2 的比例分割數(shù)據(jù)集，采用逐步線性回歸的方法，80%的數(shù)據(jù)來(lái)建立逐步回歸模型，剩余20%以作驗(yàn)證。

2.3 非參數(shù)回歸模型

非參數(shù)回歸是隨著近年機(jī)器學(xué)習(xí)的興起而衍生的一類算法，這類算法根據(jù)特征集預(yù)測(cè)目標(biāo)值時(shí)，無(wú)需假設(shè)目標(biāo)函數(shù)的形式，具有更高的靈活性與泛化性能。采用隨機(jī)森林［26］（Ran‐dom Forest，RF）和極端梯度提升［27］（Extreme Gradient Boosting，XGBoost）兩種非參算法探索遼河流域干流及主要支流TN 濃度反演最優(yōu)方式。同樣按照8∶2 的比例分割數(shù)據(jù)集，80%的數(shù)據(jù)作非參數(shù)算法的訓(xùn)練集，其余的20%驗(yàn)證算法的預(yù)測(cè)準(zhǔn)確度。

2.3.1 隨機(jī)森林算法

隨機(jī)森林算法是一種以決策樹(shù)為基學(xué)習(xí)器的有監(jiān)督集成學(xué)習(xí)算法，面對(duì)回歸分析的問(wèn)題時(shí)，取多個(gè)決策樹(shù)預(yù)測(cè)結(jié)果的均值作為目標(biāo)值，具有泛化能力突出、易于實(shí)現(xiàn)、抗噪能力強(qiáng)、計(jì)算速度快等優(yōu)勢(shì)，有利于遼河流域干流及主要支流TN 濃度的預(yù)測(cè)。RF 在訓(xùn)練的過(guò)程中使用bootstrapping 抽樣方法，約有36.8% 的數(shù)據(jù)不會(huì)參與RF 的訓(xùn)練：對(duì)于數(shù)據(jù)集D={Y|X1，…，Xn}，RF先使用這部分?jǐn)?shù)據(jù)計(jì)算預(yù)測(cè)錯(cuò)誤率，隨機(jī)打亂某一特征Xa的排序后，再次計(jì)算預(yù)測(cè)錯(cuò)誤率，取兩次錯(cuò)誤率之差量化特征Xa的權(quán)重，若RF 構(gòu)建了β個(gè)基學(xué)習(xí)器，則特征Xa的權(quán)重為β個(gè)基學(xué)習(xí)器前后兩次預(yù)測(cè)錯(cuò)誤率之差的均值，因此RF也具備特征選擇的能力。

2.3.2 極端梯度提升算法

首創(chuàng)于2016年的XGBoost 也是一種集成學(xué)習(xí)算法，其基學(xué)習(xí)器可以是樹(shù)模型，也可以是線性模型。與RF 不同的是，XGBoost采用了boosting 集成框架，使學(xué)習(xí)器根據(jù)當(dāng)前的損失殘差，以梯度提升的方式加入另一個(gè)學(xué)習(xí)器，降低損失殘差，通過(guò)不斷地迭代達(dá)到更高的預(yù)測(cè)準(zhǔn)確率，屬于RF模型的拓展。

2.4 結(jié)果驗(yàn)證

參數(shù)、非參數(shù)兩類回歸模型反演遼河流域干流及部分支流TN 濃度的準(zhǔn)確度采用決定系數(shù)（R2）、均方誤差（Root Mean Square Error，RMSE）、平均相對(duì)誤差（Average Relative Error，MAE）評(píng)價(jià)，計(jì)算公式見(jiàn)式（3）～（5）。

式中：xi為驗(yàn)證數(shù)據(jù)集中第i個(gè)樣本的實(shí)測(cè)值（或?qū)崪y(cè)值對(duì)數(shù)），i=1，2，3，…，n；xi，p為驗(yàn)證數(shù)據(jù)集中第i個(gè)樣本的對(duì)應(yīng)的預(yù)測(cè)值（或預(yù)測(cè)值對(duì)數(shù)），i=1，2，3，…，n。

3 結(jié)果與討論

3.1 水質(zhì)變化趨勢(shì)與大氣校正結(jié)果

從水樣檢測(cè)結(jié)果來(lái)看，遼河流域干流及部分支流TN 濃度較高。240 組TN 濃度檢測(cè)結(jié)果中超地表水Ⅴ類標(biāo)準(zhǔn)（＞2 mg/L）的占比76.67%，其中4月、9月劣Ⅴ類占比均在80%以上，8月最低為65%。在時(shí)間上也存在較大的差異，4-9月變異系數(shù)（CV）在0.67～1.14 之間，其中6-8月3 個(gè)月TN 濃度均值均小于4 mg/L，低于4-5月、9月［圖2（a）］，由于遼河流域年內(nèi)降水量分布差異較大［28］，降水集中在6-8月，故而徑流量較大，污染物濃度低。

抽取30%的點(diǎn)位，觀察大氣校正后的光譜曲線，如圖2（b）所示，藍(lán)光（B）區(qū)普遍呈吸收狀態(tài)，隨著波長(zhǎng)的增加，約81.1%的點(diǎn)位在綠光（G）區(qū)出現(xiàn)明顯的反射，紅光（R）與近紅外（NIR）區(qū)反射率逐漸降低，光譜曲線形態(tài)整體呈“倒U”形，反射率均小于10%；其余點(diǎn)位峰值由綠光區(qū)向紅光區(qū)偏移，反射率值高于10%。根據(jù)已有的研究［29］：天然水體的反射率普遍居低（小于10%），在綠光區(qū)會(huì)存在明顯的反射峰，隨著波長(zhǎng)的增加，反射率與波長(zhǎng)成反比；當(dāng)水體渾濁或含沙量較高時(shí)，反射峰會(huì)從綠光區(qū)移向紅光區(qū)，反射率甚至?xí)哂?0%，這種現(xiàn)象通常出現(xiàn)在河流下游區(qū)域。大氣校正結(jié)果接近準(zhǔn)確。

圖2 TN濃度（4-9月）與大氣校正后的光譜曲線Fig.2 TN concentration（April to September）and Water spectrum curve after atmospheric correction

3.2 基于參數(shù)回歸模型的TN濃度反演

為探究遼河流域干流及部分支流TN 的反演是否可以按季節(jié)、水期分布，將時(shí)間按照春季（4-5月）、夏季（6-8月）、平水期（5-6月）、豐水期（8-9月）進(jìn)行劃分。對(duì)各月、各季度、各水期TN 濃度的正態(tài)分布性進(jìn)行了檢驗(yàn)，均不服從正態(tài)分布（K-S 檢驗(yàn)Sig＜P=0.05），故采用Spearman 秩相關(guān)系數(shù)來(lái)描述TN 濃度與63個(gè)潛在特征的相關(guān)性，以相關(guān)顯著的特征參與參數(shù)回歸模型的建立。相關(guān)性計(jì)算結(jié)果如圖3所示，4-6月、8-9月、春季、夏季、平水期、豐水期分別與TN 濃度相關(guān)性顯著的特征累計(jì)為54、26、30、15、5、37、31、41和16個(gè)。

圖3 潛在特征與TN濃度的Spearman秩相關(guān)性Fig.3 Spearman rank correlation between latent featuresand TN concentration

“3σ”法對(duì)異常數(shù)據(jù)分析后，步進(jìn)條件設(shè)置為0.15＜Per‐cent_F＜0.25，參數(shù)回歸模型的建立與驗(yàn)證結(jié)果見(jiàn)表2，4-6月、9月宜逐月反演，8月宜劃分為豐水期來(lái)反演，所對(duì)應(yīng)特征組合分別為（X41）、（X41，X1，X11）、（X41，X36）、（X53）和（X53，X49，X17）。隨著敏感特征的逐步增加，R2不斷遞增，F(xiàn)在遞減，但并非最高R2所對(duì)應(yīng)的驗(yàn)證R2也最高，可能是冗余特征稀釋了敏感特征的信息，產(chǎn)生了過(guò)擬合的情況。綜合考慮驗(yàn)證R2、驗(yàn)證RMSE和驗(yàn)證MAE，最終以X41 作為4月和平水期TN 濃度的敏感特征；X53 為8、9月和夏季TN 濃度敏感特征；X41、X1 和X11共同作為5月TN 濃度敏感特征；X41和X36共同作為6月TN 濃度敏感特征；X41、X10、X4 共同作為春季TN 濃度敏感特征；X53、X49 和X17 共同作為豐水期TN 濃度敏感特征。同時(shí)從驗(yàn)證R2來(lái)看：夏季參數(shù)回歸模型驗(yàn)證結(jié)果最好，5、4、8月、平水期和9月次之，6月、豐水期和春季居后；但從驗(yàn)證RMSE和MAE來(lái)看，夏季和平水期均要高于6月，因此6月TN 濃度不宜劃分為夏季來(lái)反演；8月驗(yàn)證RMSE雖然小于夏季，但與豐水期相比，驗(yàn)證RMSE與MAE較高，且在模型建立時(shí)R2也不及豐水期，因此8月TN濃度宜劃分為豐水期來(lái)反演，4-6月、9月逐月反演。

表2 TN濃度反演的參數(shù)回歸模型Tab.2 Parameter regression model for TN concentration inversion

3.3 基于非參數(shù)回歸模型的TN濃度反演

非參數(shù)回歸模型能夠很好地處理高維特征數(shù)據(jù)，先將63個(gè)特征全部引入RF 的訓(xùn)練，超參數(shù)的尋優(yōu)采用3 折交叉驗(yàn)證（20%的驗(yàn)證數(shù)據(jù)）與網(wǎng)格搜索法。從預(yù)測(cè)結(jié)果的準(zhǔn)確度上來(lái)看（圖4），非參算法反演TN濃度宜逐月進(jìn)行。圖4中曲線為R2，柱狀為特征權(quán)重，隨著將特征按照權(quán)重高低依次加入RF，R2不斷波動(dòng)，總體呈上升趨勢(shì)，到曲線尾部逐漸趨于平直?？梢园l(fā)現(xiàn)：并不是所有的特征都參與到RF 的訓(xùn)練才使得R2最高，因此非參算法反演TN濃度時(shí)，先要對(duì)特征進(jìn)行篩選。

圖4 RF反演TN濃度學(xué)習(xí)曲線Fig.4 The learning curve for RF inversion of TN concentration

為了使特征的可解釋性更高，進(jìn)一步壓縮特征維度，取RF驗(yàn)證集最高R2對(duì)應(yīng)的特征組合與Spearman 秩相關(guān)性顯著的特征組合的交集XS=Xrf∩XSpearman，在XGBoost 中進(jìn)行逐月反演，超參的尋優(yōu)同樣采用3折交叉驗(yàn)證（20%的驗(yàn)證數(shù)據(jù)）與網(wǎng)格搜索，從圖5可以看出：整體的反演效果相比RF略有提升，MAE分布在7.27%～26.99%之間，5月最佳，8月R2較低，TN 濃度大于10 mg/L時(shí)，反演效果較差。

圖5 XGBoost反演TN濃度的結(jié)果Fig.5 XGBoost inversion results of TN concentration

此外，通過(guò)對(duì)比，非參數(shù)回歸模型與參數(shù)回歸模型類似，特征X41，即（B1+B3）/B2均為4-6月逐月、春季、平水期TN濃度的敏感特征；在特征X41 的基礎(chǔ)上加入B4 后，生成（B1+B3+B4）/B2即特征X51，均為8-9月逐月、夏季、豐水期TN濃度的敏感特征。從每月水樣采集期間水文監(jiān)測(cè)站泥沙含量來(lái)看［圖6（a）］，8、9月日均泥沙濃度分別為0.044 91、0.054 44 kg/m3，明顯高于4-6月［0.022 34 kg/m3，0.025 39 kg/m3］；從遙感反射率來(lái)看［圖6（b）］，B3 波段（R）到B4 波段（NIR），8、9月反射率開(kāi)始增高，而4-6月均呈下降趨勢(shì)。由于8、9月泥沙濃度的升高，水體較為渾濁，引起水體光譜反射率的變化，反射峰值移動(dòng)到NIR［23，29］，而水體中的氮元素會(huì)以泥沙作為吸附載體［30］，故NIR 波段反射率對(duì)8-9月TN濃度的反演起到重要的作用。

圖6 采樣期間撫順（二）站日均泥沙濃度和遙感反射率Fig.6 During the sampling period，the daily average sediment concentration and remote sensing reflectivity of Fushun（2）Station

3.4 模型適宜性評(píng)價(jià)及誤差分析

對(duì)比兩類算法反演結(jié)果表明：使用XGBoost 模型的反演效果較好，適合遼河流域干流及部分主要支流TN 濃度反演。由于水體的光學(xué)特性較為復(fù)雜，影響因素眾多［31-33］，可以說(shuō)明遼河流域干流及部分支流TN 濃度與遙感反射率（或組合）之間不具有明確的線性關(guān)系，非參數(shù)回歸模型類似于“黑箱”，無(wú)需確定的函數(shù)關(guān)系，且對(duì)復(fù)雜的高維特征也能夠較好地應(yīng)對(duì)。

XGBoost 模型反演遼河流域干流及部分支流TN 濃度的實(shí)測(cè)值與預(yù)測(cè)值之間存在7.27%～26.99%的相對(duì)平均誤差?？赡苡梢韵略蛞穑?/p>

（1）流域尺度遙感數(shù)據(jù)的局限性。亞米級(jí)分辨率的遙感影像對(duì)于細(xì)小河流的觀測(cè)效果固然更好，以GF2-PMS 影像為例，23 km 的成像幅寬對(duì)方圓21.9 萬(wàn)km2的遼河流域而言，經(jīng)濟(jì)適用性較低。國(guó)產(chǎn)GF1-WFV 成像空間分辨率比廣泛使用的Landsat 系列等要高，盡管在采樣時(shí)盡可能選取在河道中央，但對(duì)于河寬不足16 m的河段，單個(gè)像元內(nèi)混合了河流及兩岸其他物體的信息，其記錄的反射率值仍可能不屬于真實(shí)水體。

（2）大氣校正的系統(tǒng)誤差。FLAASH 大氣校正模型能最大程度上消除大氣對(duì)地表真實(shí)情況的影響，盡可能達(dá)到逼真的效果，但校正過(guò)程中所使用的經(jīng)驗(yàn)參數(shù)無(wú)法消除衛(wèi)星成像整個(gè)過(guò)程中的系統(tǒng)誤差，影像所記錄的地表信息與實(shí)際不完全相符。

（3）水質(zhì)監(jiān)測(cè)數(shù)據(jù)與遙感數(shù)據(jù)的時(shí)間吻合度。水質(zhì)與遙感反射率都是水體的瞬時(shí)信息，水樣采集時(shí)間與衛(wèi)星成像時(shí)間雖控制在±5 d 以內(nèi)，但這段時(shí)間內(nèi)采樣點(diǎn)的TN 濃度是存在波動(dòng)性的，總會(huì)存在差異，客觀上與遙感反射率不完全吻合。

3.5 TN濃度的時(shí)空分布

對(duì)三條水系的反演結(jié)果進(jìn)行局部展示（圖7），從時(shí)間上來(lái)看，遼河流域干流及部分支流TN 濃度在4-6月、8月期間呈逐月遞減趨勢(shì)，9月又開(kāi)始增加，豐水期低于平水期與枯水期；從空間上來(lái)看，太子河本溪段下游TN 濃度在0.43～8.02 mg/L，明顯優(yōu)于清河入干前（0.96～8.17 mg/L）和渾河沈陽(yáng)城區(qū)段下游（1.13～9.88 mg/L）。

圖7 遼河、渾河和太子河局部反演結(jié)果Fig.7 Local inversion results of Liaohe，Hunhe and Taizihe

清河在入干前與渾河沈陽(yáng)城區(qū)段下游TN 濃度在時(shí)間上的變化趨勢(shì)一致，4-5月明顯高于6、8、9 三個(gè)月。據(jù)《遼河治理攻堅(jiān)戰(zhàn)支流河目標(biāo)清單》，清河支流馬仲河、渾河上游撫順段支流東洲河、章黨河、將軍河為重污染河流，考核斷面水質(zhì)現(xiàn)狀均為劣Ⅴ類，且2018年上半年清河支流寇河進(jìn)行了河道清淤工程，這期間水質(zhì)較差，可能是支流匯流對(duì)干流的水質(zhì)產(chǎn)生了較大的影響；此外，2018年上半年，受降水影響，渾河流量同比減少47.9%，且撫順、沈陽(yáng)城區(qū)筑壩較多，使水流不暢，水質(zhì)也可能會(huì)隨著變差。

4 結(jié) 論

目前對(duì)湖庫(kù)等寬闊水域水環(huán)境參數(shù)定量遙感反演的研究相對(duì)成熟，大流域尺度河流水體因受遙感成像空間、時(shí)間分辨率與幅寬的制約，有待進(jìn)一步認(rèn)識(shí)。以遼河流域?yàn)檠芯繀^(qū)域，TN 為反演對(duì)象，利用2018年4-9月（不含7月）實(shí)測(cè)水質(zhì)數(shù)據(jù)和準(zhǔn)同步GF1-WFV 遙感影像，基于參數(shù)、非參數(shù)兩類回歸模型，對(duì)遼河流域干流及部分支流TN 濃度進(jìn)行了反演研究，得出結(jié)論：

（1）GF1-WFV 影像空間分辨率16 m，時(shí)間分辨率2 d，成像幅寬800 km，多光譜相機(jī)成像譜段涵蓋藍(lán)、綠、紅和近紅外4 個(gè)波段，大氣校正后的反射率普遍小于10%，綠光區(qū)為反射峰，隨著波長(zhǎng)的增加，反射率下降，滿足遼河流域干流及主要支流TN濃度反演的需求。但8-9月（豐水期）由于降雨沖刷導(dǎo)致泥沙量增大，反射峰朝近紅外偏移，水體中的氮元素會(huì)吸附在泥沙上，因此近紅外波段反射率在8-9月遼河流域干流及部分支流TN濃度反演中具有重要意義。

（2）遼河流域干流及部分支流TN 濃度與遙感反射率（或組合）之間線性關(guān)系較弱，使用參數(shù)回歸模型反演結(jié)果不及非參數(shù)回歸模型。非參數(shù)回歸模型中，極端梯度提升的樹(shù)模型（XG‐Boost）比隨機(jī)森林（RF）反演效果更好，4-6、8-9月逐月的R2均在0.575 以上，RMSE在0.54～1.899 之間，MAE在7.27%～26.99%之間。當(dāng)TN濃度在10 mg/L以上時(shí)，反演效果不理想。

（3）內(nèi)陸水體的光學(xué)特性較為復(fù)雜，影響因素眾多，遼河流域干流及部分支流TN 濃度反演需多個(gè)特征參與。與TN 濃度相關(guān)的潛在特征中存在非相關(guān)、低相關(guān)的冗余特征，不僅不能提高反演的準(zhǔn)確度，還難以被解釋，因此使用XGBoost模型前要對(duì)特征進(jìn)行選擇。取RF 驗(yàn)證集最高R2對(duì)應(yīng)的特征組合與Spearman 秩相關(guān)性顯著的特征組合交集（XS=Xrf∩XSpearman）的方式能夠較好地完成對(duì)特征的選擇。

（4）2018年遼河流域干流及部分支流TN 濃度在時(shí)間上存在較大的差異，4-6月明顯高于8-9月，8月達(dá)到最低（1.675～7.675 mg/L）；受清河支流馬仲河、寇河，渾河上游撫順段支流東洲河、章黨河、將軍河劣Ⅴ類水體匯流，以及渾河流量減少的影響，4-5月清河入干前（1.47～8.17 mg/L）和渾河沈陽(yáng)城區(qū)段下游（1.13～9.88 mg/L）TN 污染情況比太子河本溪段下游（0.48～8.02 mg/L）嚴(yán)重。