李爽 向瑋

(北京林業(yè)大學(xué)，北京,100083)

生長(zhǎng)輪分析作為研究樹木生長(zhǎng)與環(huán)境因子相互作用的定量分析方法，在不同空間尺度下的氣候重建中得到了廣泛應(yīng)用[1-2]。此外，生長(zhǎng)輪本身指示著林分生長(zhǎng)動(dòng)態(tài)，可用于模擬林分生物量[3-4]、林木枯損[5]、碳蓄積[6]，一直是生態(tài)學(xué)研究的關(guān)注熱點(diǎn)。隨著研究手段的進(jìn)步，生長(zhǎng)輪分析在考古學(xué)[7]、地質(zhì)學(xué)[8]、災(zāi)害重建[9-10]等領(lǐng)域也取得了豐富的研究成果。生長(zhǎng)輪寬度作為樹木年代學(xué)中最常用的指標(biāo)，其準(zhǔn)確測(cè)量是樹木年代學(xué)研究的前提，直接影響了年表構(gòu)建和氣候重建的可靠性[11]，也對(duì)測(cè)量方法提出了更高的要求。

目前，生長(zhǎng)輪寬度測(cè)量多使用WinDENDRO或LinTAB等商業(yè)化系統(tǒng)，這類軟硬件高度集成的專用工具為研究者提供了成熟、穩(wěn)定、高精度的分析平臺(tái)，但售價(jià)昂貴，使用門檻較高[12-13]。隨著計(jì)算機(jī)技術(shù)的快速發(fā)展，數(shù)字圖像處理技術(shù)為生長(zhǎng)輪寬度測(cè)量帶來(lái)了全新的思路和啟發(fā)[14-15]。生長(zhǎng)輪數(shù)字圖像的使用在保證測(cè)量精度的前提下，不僅提高了測(cè)量效率[16-17]，也顯著降低了軟硬件成本[12，18-19]。對(duì)于早晚材過渡不明顯的樹種，圖像增強(qiáng)技術(shù)可以提高人眼的辨識(shí)度[20]，而圖像拼接技術(shù)也為L(zhǎng)inTAB難以測(cè)量的大體積圓盤樣品提供了解決方案[21-22]。

R作為一款開源軟件，具有強(qiáng)大的統(tǒng)計(jì)和繪圖功能，可在Windows和Mac等主流操作系統(tǒng)運(yùn)行。相比其他軟件，R的開源特性使得研究者可以根據(jù)自身需求進(jìn)行二次開發(fā)，具有極大的靈活性[23]。近年來(lái)，國(guó)內(nèi)外研究者開發(fā)了一系列用于生長(zhǎng)輪分析的R軟件包[17，24-26]。這些R包的出現(xiàn)為樹木年代學(xué)提供了一個(gè)全新的數(shù)據(jù)分析平臺(tái)，實(shí)現(xiàn)了交叉定年、去趨勢(shì)、年表構(gòu)建等功能[27]，是對(duì)COFECHA和ARSTAN等傳統(tǒng)軟件的良好補(bǔ)充，推動(dòng)了樹輪年代學(xué)研究的發(fā)展。

盡管R軟件有許多優(yōu)點(diǎn)，但命令行的操作方式對(duì)使用者編程水平提出了較高要求，這也成為R軟件推廣普及的主要障礙。目前國(guó)內(nèi)少有介紹R軟件在樹木年代學(xué)中應(yīng)用的文獻(xiàn)，本研究目的是以實(shí)例呈現(xiàn)在R軟件中使用MtreeRing測(cè)量生長(zhǎng)輪寬度的流程，對(duì)其圖像處理算法和誤差校正原理進(jìn)行簡(jiǎn)單介紹，提出需要注意的問題和使用細(xì)節(jié)，以期為研究者提供參考。

1 MtreeRing安裝

MtreeRing的主要用途是從掃描圖像中獲取生長(zhǎng)輪寬度序列，其R函數(shù)提供了以下功能(見表1)：①讀取生長(zhǎng)輪掃描圖像；②使用邊緣檢測(cè)算法自動(dòng)識(shí)別生長(zhǎng)輪邊界；③手動(dòng)刪除或添加生長(zhǎng)輪邊界；④根據(jù)邊界坐標(biāo)計(jì)算和校正生長(zhǎng)輪寬度。該軟件包在安裝目錄提供了2張分辨率為1 200 dpi的圖像，其中，001.png用于演示生長(zhǎng)輪寬度測(cè)量，而missing_pith.png用于演示髓心丟失樣品的寬度校正，本研究也將使用這兩張圖片來(lái)演示MtreeRing的功能。

用戶首先需要安裝R軟件，其安裝包可以從官方網(wǎng)站(https://cran.r-project.org/)免費(fèi)獲取，本研究使用的軟件版本為3.6.1。由于MtreeRing的圖形界面是基于Shiny框架開發(fā)的，因此推薦用戶使用RStudio運(yùn)行該R包，以獲得更好的體驗(yàn)。關(guān)于R軟件及RStudio的安裝細(xì)節(jié)，可參閱《R語(yǔ)言實(shí)戰(zhàn)》一書[23]。

啟動(dòng)R軟件后，在命令提示符“>”后輸入install.packages(“MtreeRing”)命令，在彈出窗口選擇CRAN鏡像站點(diǎn)下載安裝包，安裝完成后輸入library(MtreeRing)命令將MtreeRing加載到當(dāng)前環(huán)境。本研究使用1.4.2版本(https://CRAN.R-project.org/package=MtreeRing)。

表1 MtreeRing軟件包中的R函數(shù)簡(jiǎn)介

2 生長(zhǎng)輪寬度測(cè)量

外業(yè)工作采集的樣品首先應(yīng)進(jìn)行膠裝、打磨、掃描等預(yù)處理，具體步驟可參考羅春旺[28]與孟盛旺[19]的研究。MtreeRing中一次完整的寬度測(cè)量由圖像讀取、邊界檢測(cè)、寬度計(jì)算與修正、數(shù)據(jù)導(dǎo)出4部分組成。

2.1 圖像讀取

使用library命令加載MtreeRing后，用戶可以使用ring_read函數(shù)將圖像文件讀入R軟件，目前支持的文件格式包括png、jpg、tiff、bmp，具體命令如下：

img.path<-system.file(“001.png”,package=“MtreeRing”)；

t1<-ring_read(img=img.path,dpi=1200,plot=TRUE)。

第一行命令system.file將從MtreeRing的安裝目錄獲取圖像001.png的完整路徑，并以字符串格式賦值到對(duì)象img.path。第二行命令ring_read完成圖像讀取和繪圖，并將圖像數(shù)據(jù)保存到了名為t1的對(duì)象。其中，參數(shù)img指定了圖像文件路徑；dpi為表示圖像分辨率的數(shù)值型變量，應(yīng)根據(jù)研究需要選擇，通常不低于600 dpi；plot為邏輯值，表示打開新的R圖形窗口并繪出生長(zhǎng)輪圖像。為加快讀取和繪圖速度，圖像不應(yīng)含有無(wú)關(guān)內(nèi)容，如粘貼樣芯的木質(zhì)凹槽；讀圖時(shí)，樹皮和髓心應(yīng)分別放置在R圖形窗口的左右兩側(cè)，用戶可以使用參數(shù)rotate旋轉(zhuǎn)原始圖像。

2.2 路徑創(chuàng)建與邊界檢測(cè)

讀圖后，用戶可以使用ring_detect函數(shù)創(chuàng)建測(cè)量路徑并沿路徑識(shí)別生長(zhǎng)輪邊界。該函數(shù)是MtreeRing的核心函數(shù)，用戶可以用鼠標(biāo)交互方式在R圖形窗口中創(chuàng)建路徑，并以3種可選方法自動(dòng)識(shí)別邊界，ring_detect的重要參數(shù)見表2。本研究將使用分水嶺算法對(duì)圖像中的生長(zhǎng)輪邊界進(jìn)行識(shí)別，具體代碼如下：

t2<-ring_detect(ring.data=t1,seg=3,method=‘watershed’,border.color=‘green’)。

該命令將開啟一個(gè)或若干個(gè)R圖形窗口，并顯示檢測(cè)到的生長(zhǎng)輪邊界(圖1)。

由于R圖形窗口不支持滾動(dòng)縮放，因此當(dāng)樣品的寬度序列較長(zhǎng)時(shí)，使用參數(shù)seg將原圖像分割為等距離的若干段(圖1中為3段)，以獲得更好的顯示效果。ring_detect自動(dòng)創(chuàng)建了一條從圖像中心穿過的虛線路徑(單路徑模式)，分水嶺算法檢測(cè)到的生長(zhǎng)輪邊界以綠色突出顯示，并標(biāo)有生長(zhǎng)輪編號(hào)。該函數(shù)的圖像處理參數(shù)(未在表2列出)適用于常見針葉樹種，熟悉圖像處理技術(shù)的用戶可以修改這些參數(shù)以適應(yīng)特定解剖學(xué)特征的樣品。MtreeRing所使用的形態(tài)學(xué)運(yùn)算及分水嶺分割、Canny邊緣檢測(cè)等算法的細(xì)節(jié)可以參考Soille et al.[29]和史景寧等[17]文獻(xiàn)。

實(shí)際采樣中，樣芯常出現(xiàn)髓心丟失的現(xiàn)象，導(dǎo)致靠近髓心的生長(zhǎng)輪以一定角度傾斜(圖2)，此時(shí)若直接使用邊界之間的水平距離會(huì)導(dǎo)致測(cè)量結(jié)果偏大。這種情況下，用戶可以將ring_detect的參數(shù)incline設(shè)置為TRUE，函數(shù)會(huì)將兩條相互平行的路徑同時(shí)添加到圖像(雙路徑模式，見圖2)，隨后的測(cè)量將在兩條路徑上同時(shí)進(jìn)行。

表2 ring_detect函數(shù)的重要參數(shù)

2.3 生長(zhǎng)輪邊界修正

當(dāng)生長(zhǎng)輪圖像存在噪聲(如樹脂管道、物理?yè)p傷)或邊界不清晰時(shí)，自動(dòng)識(shí)別會(huì)出現(xiàn)過度檢測(cè)或錯(cuò)過邊界的現(xiàn)象。此時(shí)用戶可以使用ring_modify函數(shù)刪除或添加新的邊界點(diǎn)。以圖3為例，編號(hào)為21和23的生長(zhǎng)輪邊界位置錯(cuò)誤，用戶需要?jiǎng)h除這兩條邊界并使用鼠標(biāo)單擊添加2個(gè)新邊界，具體代碼如下：

t2<-ring_modify(ring.data=t2,del=c(21,23),add=TRUE)。

該命令中，參數(shù)del表示要?jiǎng)h除的生長(zhǎng)輪邊界編號(hào)，當(dāng)有多個(gè)邊界要?jiǎng)h除時(shí)，使用c函數(shù)連接各編號(hào)；參數(shù)add為邏輯值，表示是否添加新邊界。當(dāng)add=TRUE時(shí)，運(yùn)行該命令后ring_detect創(chuàng)建的R圖形窗口會(huì)被依次激活(窗口左上角名稱顯示為ACTIVE)，添加邊界的具體步驟如下：①將鼠標(biāo)移動(dòng)至被激活的R圖形窗口內(nèi)，此時(shí)光標(biāo)變?yōu)槭中?；②在待添加的邊界與路徑交點(diǎn)處單擊一次，圖像上會(huì)顯示該點(diǎn)；③如果一個(gè)圖形窗口有多個(gè)邊界待添加，重復(fù)步驟2；④添加完邊界后，單擊右鍵，在菜單欄選擇“停止”；⑤如果有多個(gè)圖形窗口(參數(shù)seg大于1)，重復(fù)前述步驟即可。當(dāng)最后一個(gè)圖形窗口被關(guān)閉后，ring_modify會(huì)重新打開這些窗口，并為所有邊界點(diǎn)重新編號(hào)，供用戶檢查。

2.4 寬度計(jì)算與結(jié)果導(dǎo)出

當(dāng)所有生長(zhǎng)輪邊界都被正確識(shí)別后，可以使用ring_calculate函數(shù)計(jì)算生長(zhǎng)輪寬度，代碼如下：

rw<-ring_calculate(ring.data=t2,seriesID=“940220”)。

寬度序列rw為數(shù)據(jù)框格式，參數(shù)seriesID為數(shù)據(jù)框的列名。單路徑模式下，根據(jù)生長(zhǎng)輪邊界坐標(biāo)計(jì)算相鄰兩點(diǎn)之間的距離，并使用圖像分辨率(R)將原始單位轉(zhuǎn)換為mm，公式如下：

(1)

(2)

式中：xi和yi分別為第i個(gè)生長(zhǎng)輪邊界的橫、縱坐標(biāo)；Di為第i和i+1個(gè)生長(zhǎng)輪之間的原始距離；di為實(shí)際生長(zhǎng)輪寬度(mm)。

在雙路徑模式下，ring_calculate函數(shù)使用生長(zhǎng)輪邊界傾斜角度αi對(duì)相鄰邊界之間的水平距離Di進(jìn)行校正。傾斜角度定義為生長(zhǎng)輪邊界的切線與路徑法線之間的夾角，當(dāng)傾斜超過5°時(shí)，應(yīng)考慮對(duì)寬度序列進(jìn)行誤差校正[28]，所用公式如下：

(3)

αi=arctan(mi)；

(4)

(5)

式中：(Xi,Yi)和(xi,yi)分別為第i個(gè)生長(zhǎng)輪的邊界與上、下路徑的交點(diǎn)；mi是第i個(gè)生長(zhǎng)輪邊界切線的斜率；Di是根據(jù)三角函數(shù)關(guān)系校正后的生長(zhǎng)輪寬度。

需要指出的是，髓心丟失的樣品最內(nèi)側(cè)2～3個(gè)生長(zhǎng)輪常為圓弧狀(圖2)，且只能與一條路徑相交，因而無(wú)法使用雙路徑模式校正誤差。這種情況下，由于最內(nèi)側(cè)弧線通常為同心圓[30]，用戶可以使用pith_measure函數(shù)定位髓心，并完成寬度測(cè)量與校正。該函數(shù)帶有插圖的示例可以在GitHub下載(https://ropensci.github.io/MtreeRing/articles/pith-MtreeRing.html)。

ring_calculate所得的寬度序列可以直接在R軟件中進(jìn)行交叉定年和去趨勢(shì)，也可以使用dplR包導(dǎo)出，具體代碼如下：

library(dplR)；

write.rwl(rwl.df=rw,fname=“C:Users