李福娟，曹麗歆，王鑫平，李光梅，周 青，孫培艷＊

（國(guó)家海洋局1.北海環(huán)境監(jiān)測(cè)中心；2.海洋溢油鑒別與損害評(píng)估技術(shù)重點(diǎn)實(shí)驗(yàn)室，山東 青島266033）

在海洋石油勘探開發(fā)過(guò)程中，極易因?yàn)樵囉?、運(yùn)輸、儲(chǔ)存以及其他原因造成突發(fā)性溢油事件。自1980年代以來(lái)，溢油事件呈上升趨勢(shì)，幾乎每年都發(fā)生由于井噴、漏油以及原油運(yùn)輸船舶的撞船、沉船等造成的溢油事件，這些事件會(huì)造成事故海域及流域的嚴(yán)重污染，而且還能間接危及陸地生物、人類健康和自然資源。對(duì)已發(fā)生的溢油污染事故，需要及時(shí)準(zhǔn)確地確定污染源，以便實(shí)施應(yīng)急措施、確定責(zé)任方、解決糾紛[1－2]。因此辨別溢油的來(lái)源便成為一個(gè)亟待解決的問(wèn)題。

目前，在所有的方法技術(shù)中，氣相色譜法、氣相色譜－質(zhì)譜聯(lián)用是應(yīng)用最廣泛的。氣相色譜分析技術(shù)尤其是毛細(xì)管柱氣相色譜－質(zhì)譜聯(lián)用儀技術(shù)，已經(jīng)得到較大提高，能夠較好的分析石油中烴類物質(zhì)。

近年來(lái)，越來(lái)越多的化學(xué)計(jì)量學(xué)方法應(yīng)用到實(shí)驗(yàn)數(shù)據(jù)處理中來(lái)，并且已經(jīng)變得非常重要。由于油品組成非常復(fù)雜，很難實(shí)現(xiàn)對(duì)所有樣品的信息進(jìn)行分析比較。大多數(shù)化學(xué)計(jì)量學(xué)方法只能從所獲得的數(shù)據(jù)中提取最能代表原油特征的信息加以利用。這樣做往往會(huì)丟失一些有用數(shù)據(jù)，存在著一定的局限性[3]。因此，可以在化學(xué)計(jì)量學(xué)分析中使用整個(gè)色譜數(shù)據(jù)來(lái)解決這個(gè)問(wèn)題，然后通過(guò)聚類分析衡量不同樣品的相似性。但是，實(shí)驗(yàn)過(guò)程中條件的稍稍改變就會(huì)引起保留時(shí)間的漂移，此方法對(duì)在時(shí)間軸上相近信息的提前或推后出現(xiàn)會(huì)產(chǎn)生較大的偏差，必須要對(duì)整個(gè)色譜圖的保留時(shí)間進(jìn)行校正。近年來(lái)，分析工作者提出了一系列檢測(cè)和消除譜峰漂移的方法，例如，偏線性擬合方法（Partial linear fit，PLF）[4－5]、偽主成分回歸 （Pseudo－principal component regression，PPCR）[6]、自動(dòng)峰疊合算法（Automatic peak alignment algorithm，APAA）[7]、動(dòng)態(tài)時(shí)間校正（Dynamic time warping，DTW）和相關(guān)系數(shù)優(yōu)化校正（Correlation optimized warping，COW）[8－14]等。其中在校正時(shí)間軌道、色譜圖和光譜圖方面，相關(guān)系數(shù)優(yōu)化校正因其耗時(shí)少，效果好而引起了人們的廣泛關(guān)注。目前，有多種聚類分析方法，根據(jù)計(jì)算方法的不同分為：歐氏距離法、馬氏距離法、類平均法、重心法、最長(zhǎng)距離法、最短距離法、密度估計(jì)法、Ward最小方差法、主成分分析和系統(tǒng)聚類分析等[15－18]。

本文采用氣相色譜和氣相色譜－質(zhì)譜聯(lián)用2種方法進(jìn)行油樣分析，針對(duì)色譜組分保留時(shí)間的漂移問(wèn)題，采用COW方法對(duì)色譜圖的保留時(shí)間校正，大大提高了采集點(diǎn)的重合性。然后通過(guò)歐氏距離對(duì)校正后的新色譜數(shù)據(jù)進(jìn)行聚類分析，達(dá)到了預(yù)期的分類效果。這種方法充分利用全指紋譜圖的信息，使分析結(jié)果更為可靠，為溢油來(lái)源的判斷提供了一種輔助鑒別方法。

1 理論與算法

COW是由Nielsen在1998年提出的一種分段優(yōu)化數(shù)據(jù)的算法[19]。在運(yùn)算過(guò)程中，向量的端點(diǎn)固定不動(dòng)，根據(jù)松弛參數(shù)（t）將向量分成相同的段數(shù)。從最后一段開始同參照向量進(jìn)行比照、校正，在松弛參數(shù)正負(fù)范圍內(nèi)進(jìn)行優(yōu)化，得到最大相關(guān)系數(shù)的一組數(shù)據(jù)向量，然后再在此基礎(chǔ)上對(duì)第二段數(shù)據(jù)進(jìn)行優(yōu)化，依此類推，得到一組最優(yōu)的重組數(shù)據(jù)向量。當(dāng)抽樣向量的時(shí)間點(diǎn)數(shù)和參照向量的點(diǎn)數(shù)不相同時(shí)，就在抽樣向量?jī)?nèi)線性插入合適的點(diǎn)數(shù)得到相同段長(zhǎng)的預(yù)處理向量。相關(guān)系數(shù)的計(jì)算公式如下：

COW只需2個(gè)輸出參數(shù)進(jìn)行分段線性相關(guān)系數(shù)優(yōu)化校正[13]，這成為它的一大優(yōu)點(diǎn)。運(yùn)算過(guò)程中可以選擇一個(gè)適中的段長(zhǎng)和較小的松弛參數(shù)來(lái)補(bǔ)償色譜圖中產(chǎn)生的時(shí)間漂移。

歐氏距離是兩項(xiàng)間的差，即：每個(gè)變量值差值的平方和再平方根，目的是計(jì)算其間的整體距離，即不相似性。其公式如2所示：

其中：xik表示為第i個(gè)序列的第k個(gè)指標(biāo)的測(cè)定值；yjk為第j個(gè)序列的第k個(gè)指標(biāo)的測(cè)定值。Dij為第i個(gè)序列與第j個(gè)序列之間的歐氏距離。其具體應(yīng)用的一般算法過(guò)程如下：（1）收集特征數(shù)據(jù)并且建立模型特征表；（2）規(guī)格化特征表；（3）計(jì)算各序列間距離并產(chǎn)生一個(gè)距離向量；（4）實(shí)施聚類分析；（5）根據(jù)分類距離等級(jí)要求決定把目標(biāo)對(duì)象總體細(xì)分為幾組，否則回到第3步繼續(xù)；（6）產(chǎn)生分組結(jié)果。本文把信號(hào)點(diǎn)與點(diǎn)之間的歐氏距離的大小作為判別分類的依據(jù)。

2 實(shí)驗(yàn)過(guò)程

2.1 樣品信息及前處理方法

在某次溢油事故中采集2個(gè)溢油樣品，命名為樣品1與樣品2。在此次溢油事故中采集一個(gè)可疑油源命名為樣品3。另外選擇2個(gè)不相關(guān)油田的原油樣品進(jìn)行比較，其中樣品4～7為來(lái)自油田A的原油樣品，其中樣品4為非降解油，樣品5、6、7受到一定程度的降解；樣品8～11為油田B的原油樣品，均為重度降解原油（見表1）。

稱取油樣約0.8g，溶于正己烷，定容至10mL，離心10min后取上層油樣200μL至進(jìn)樣瓶，同時(shí)加入800μL正己烷，混合均勻后上機(jī)分析。

2.2 儀器及工作條件

儀器 島津GC2010氣相色譜儀；島津GCMS－QP2010氣相色譜－質(zhì)譜儀。毛細(xì)管色譜柱：DB－5（30m×0.32 mm×0.25μm）、DB－5MS（30m×0.25mm×0.25μm）（長(zhǎng)度×內(nèi)經(jīng)×膜厚）。

分析條件 正構(gòu)烷烴采用氣相色譜/氫火焰離子化檢測(cè)器（GC/FID）分析。毛細(xì)管色譜柱涂層為5%苯基、95%二甲基聚硅氧烷，涂層厚度為0.25μm，內(nèi)徑為0.32mm，長(zhǎng)度為30m。色譜分析條件如下：載氣：高純氦氣，1.0mL/min；進(jìn)樣方式：不分流；進(jìn)樣口溫度：290℃；檢測(cè)器溫度：300℃；升溫程序：在50℃保持2min，以6℃/min的速度升到300℃，保持16min。

表1 11個(gè)樣品的信息表Table 1 Information of 11samples

甾、萜烷類生物標(biāo)志化合物均采用氣相色譜/質(zhì)譜（GC/MS）分析。毛細(xì)管色譜柱涂層為5%苯基、95%二甲基聚硅氧烷，涂層厚度為0.25μm，內(nèi)徑為0.25mm，長(zhǎng)度為30m。色譜分析條件如下：載氣：高純氦氣，1.0 mL/min；進(jìn)樣方式：不分流；進(jìn)樣口溫度：290℃；接口溫度：280℃；離子源溫度：230℃；升溫程序：在50℃保持2min，以6℃/min的速度升到300℃，保持16min。

3 分析結(jié)果

3.1 數(shù)據(jù)優(yōu)化

采用COW算法對(duì)氣相色譜數(shù)據(jù)進(jìn)行校正，結(jié)果如圖1所示，圖1（a）為原始?xì)庀嗌V圖的一部分，2個(gè)色譜圖在時(shí)間軸上的相近數(shù)值之間產(chǎn)生較大的漂移。圖1（b）為相應(yīng)數(shù)據(jù)校正后的氣相色譜圖。比較這2張色譜圖，發(fā)現(xiàn)色譜峰的重合性大大提高，并且計(jì)算其相關(guān)系數(shù)，校正前為0.383 1，校正后為0.850 0，達(dá)到了理想的校正效果。

3.2 聚類分析

3.2.1 氣相色譜數(shù)據(jù)聚類結(jié)果 根據(jù)歐氏距離公式對(duì)COW校正后的氣相色譜數(shù)據(jù)進(jìn)行聚類分析。COW校正前、后色譜數(shù)據(jù)的聚類分析樹狀圖如圖2～3所示。

圖1 COW樣品校正前a、后b的部分氣相色譜圖Fig.1 Data before and after the application of the COW algorithm

圖2 基于未優(yōu)化的色譜數(shù)據(jù)的聚類結(jié)果Fig.2 Clustering result based on unoptimized chromatographic data

圖3 基于COW優(yōu)化的色譜數(shù)據(jù)的聚類結(jié)果Fig.3 Clustering result based on optimized chromatographic data

基于原始色譜數(shù)據(jù)的聚類結(jié)果為：11個(gè)樣品分為4組，溢油樣品1、2與可疑油源3為一組；來(lái)自油田A的樣品4、5、6、7分為2組，其中樣品4單獨(dú)為一組，樣品5、6、7為一組；來(lái)自油田B的樣品8、9、10、11為一組?；趦?yōu)化后的色譜數(shù)據(jù)的聚類結(jié)果為：11個(gè)樣品分為3組，樣品4與溢油樣品1、2和可疑油源3歸為一類；來(lái)自油田A的樣品5、6、7分為一組；來(lái)自油田B的樣品8、9、10、11為一組。對(duì)照樣品的色譜圖（見圖4）分析，溢油樣品1、2與可疑油源3為非降解油（以樣品1為代表），樣品4為非降解油，樣品5、6、7為降解油（以樣品5為代表），而樣品8、9、10、11為嚴(yán)重降解油（以樣品8為代表）。結(jié)合聚類分析結(jié)果與樣品的色譜圖特征來(lái)看，兩者的結(jié)果是一致的。因此，該聚類方法對(duì)于降解和非降解油來(lái)說(shuō)是行之有效的，是完全可以進(jìn)行分類的。但對(duì)受到風(fēng)化影響的油品和為了準(zhǔn)確查找溢油來(lái)源，還需要借助于質(zhì)譜數(shù)據(jù)甾萜烷生物標(biāo)志物的分析。

3.2.2 氣相色譜－質(zhì)譜數(shù)據(jù)聚類結(jié)果 甾萜烷生物標(biāo)志物代表樣品本身的性質(zhì)，且在環(huán)境中幾乎不受風(fēng)化的影響，因此，在質(zhì)譜數(shù)據(jù)中選擇甾萜數(shù)據(jù)進(jìn)行聚類分析。根據(jù)歐氏距離公式分別對(duì)COW校正前后的氣相色譜數(shù)據(jù)m/z191（萜烷）和 m/z217（甾烷）進(jìn)行聚類分析。COW校正前后的聚類分析樹狀圖如圖5～8所示。

通過(guò)聚類分析，利用m/z191的數(shù)據(jù)可以將11個(gè)樣品分為3類。溢油樣品1、2聚為一類（以樣品1為代表）；來(lái)自同一油田的4、5、6、7分為一類（以樣品4為代表）；可疑油源3與來(lái)自另一油田的樣品8、9、10、11分為一類（以樣品8為代表）。溢油樣品1、2能很好的同其他2個(gè)油田的樣品區(qū)分開來(lái)，但是也未與可疑油源樣品3聚為一類，沒有達(dá)到理想的聚類效果。對(duì)m/z191的優(yōu)化后數(shù)據(jù)進(jìn)行聚類分析后發(fā)現(xiàn)，可疑油源3與溢油樣品1、2聚為一類，另外2個(gè)油田的樣品分別歸為一類，達(dá)到了理想的聚類效果。分析樣品的萜烷譜圖見圖9，從譜圖上分析也可以看出，樣品1與樣品3譜圖一致，而與樣品4、8譜圖不一致。聚類結(jié)果與譜圖分析和樣品信息完全相符，樣品數(shù)據(jù)的優(yōu)化達(dá)到了優(yōu)化聚類結(jié)果的目的。

對(duì)COW優(yōu)化前和優(yōu)化后的m/z 217的數(shù)據(jù)進(jìn)行聚類后發(fā)現(xiàn)，2組數(shù)據(jù)的聚類結(jié)果是一致的，溢油樣品1、2與可疑油源3聚為一類，另外2個(gè)油田的樣品分別聚為一類。但從聚類效果來(lái)看，優(yōu)化后數(shù)據(jù)溢油樣品1、2與可疑油源3的差距明顯縮小。聚類效果要優(yōu)于未優(yōu)化數(shù)據(jù)的聚類效果。分析樣品的甾烷譜圖見圖10，從譜圖上分析也可以看出，樣品1與樣品3譜圖一致，而與樣品4、8譜圖不一致。聚類結(jié)果與譜圖分析和樣品信息完全相符，樣品數(shù)據(jù)的優(yōu)化達(dá)到了優(yōu)化聚類結(jié)果的目的。

圖4 4種油樣的氣相色譜圖Fig.4 Gas chromatograms of 4kinds of crude oil

圖5 基于未優(yōu)化m/z191的質(zhì)譜數(shù)據(jù)的聚類結(jié)果Fig.5 Clustering result based on unoptimized mass spectrometry data data（m/z 191）data

圖6 基于優(yōu)化m/z191的質(zhì)譜數(shù)據(jù)的聚類結(jié)果Fig.6 Clustering result based on optimized mass spectrometry data data（m/z 191）data

圖7 基于未優(yōu)化m/z217的質(zhì)譜數(shù)據(jù)的聚類結(jié)果Fig.7 Clustering result based on unoptimized mass spectrometry data data（m/z 217）data

圖8 基于優(yōu)化m/z217的質(zhì)譜數(shù)據(jù)的聚類結(jié)果Fig.8 Clustering result based on optimized mass spectrometry data data（m/z 217）data

4 結(jié)論與展望

通過(guò)氣相色譜法與氣相色譜－質(zhì)譜法分析油樣的正構(gòu)烷烴與甾萜生物標(biāo)志物，然后結(jié)合相關(guān)系數(shù)優(yōu)化（COW）方法，解決了色譜圖保留時(shí)間漂移的問(wèn)題。計(jì)算校正前后2個(gè)色譜數(shù)據(jù)的相關(guān)系數(shù)，發(fā)現(xiàn)校正后數(shù)據(jù)的相關(guān)系數(shù)大大提高，具有很好的重合性。然后以歐氏距離為判據(jù)對(duì)校正后的色譜數(shù)據(jù)進(jìn)行聚類分析，通過(guò)與實(shí)際樣品信息進(jìn)行比較后發(fā)現(xiàn)，對(duì)氣相色譜數(shù)據(jù)，該方法能有效分類降解和非降解油品；對(duì)氣相色譜－質(zhì)譜的甾、萜烷數(shù)據(jù)的聚類結(jié)果則與實(shí)際溢油來(lái)源排查結(jié)果是一致的，達(dá)到了理想的聚類效果，該聚類結(jié)果不受風(fēng)化和降解的影響，更能代表油品本身的特點(diǎn)。并且歐氏距離為判據(jù)的聚類方法僅需幾秒鐘就可完成。該方法能夠充分利用色、質(zhì)譜采集信息，使分析結(jié)果更為可靠，為辨別溢油的來(lái)源建立了一種快速分類輔助鑒別方法。

圖9 溢油樣品與可疑油源樣品的m/z191比對(duì)圖Fig.9 Compared chromatogram of spilled oil and suspicious oil sources（m/z 191）

圖10 溢油樣品與可疑油源樣品的m/z217比對(duì)圖Fig.10 Compared chromatogram of spilled oil and suspicious oil sources（m/z 217）

[1]趙玉慧，孫培艷，王鑫平，等.多環(huán)芳烴指紋用于渤海采油平臺(tái)原油的鑒別[J].色譜，2008，26（1）：43－49.

[2]楊佰娟，鄭 立，張魁英，等.原油中雙環(huán)倍半萜指紋的內(nèi)標(biāo)法分析 [J].分析測(cè)試學(xué)報(bào)，2012，31（11）：1421－1425.

[3]王鑫平，孫培艷，周青，等.原油飽和烴指紋的內(nèi)標(biāo)法分析[J].分析化學(xué)，2007，8：1121－1126.

[4]Westad F，Martens H.Shift and intensity modeling in spectroscopy－general concept and applications[J].Chemon intell Lab Syst，1999，45（1－2）：361－370.

[5]Witjes H，van den Brink M，Melssen W J，et al.Automatic correction of peak shifts in Raman spectra before PLS regression[J].Chemon intell Lab Syst，2000，52（1）：105－116.

[6]Brown T R，Atoyanova R.NMR Spectral Quantitation by Principal－Component Analysis. Ⅱ.Determination of Frequency and Phase Shifts[J].J Magn Reson，Series B，1996，112（1）：32－43.

[7]Witjes H，Melssen W J A，int Zandt H J A，et al.Automatic correction for phase shifts，frequency shifts，and lineshape distortions across a series of single resonance lines in large spectral data sets[J].J Magn Reson，2000，144（1）：35－44.

[8]Reiner E，Abbey L E，Moran T F，et al.Characterization of normal human cells by pyrolysis gas chromatography mass spectrometry[J].Biomed Mass Spectrom，1979，6：491－498.

[9]Wang C P，Isenhour T L.Time－warping algorithm applied to chromatographic peak matching gas chromatography Fourier－transform infrared mass－spectrometry[J].Anal Chem，1987，59（4）：649－654.

[10]Pravdova V，Walczak B，Massart D L.A comparison of two algorithms for warping of analytical signals[J].Anal Chim Acta，2002，456（1）：77－92.

[11]Tomasi G，van den Berg F，Andersson C.Correlation optimized warping and dynamic time warping as preprocessing methods for chromatographic data[J].J Chemometr，2004，18（5）：231－241.

[12]Itakura F. Minimum prediction residual principle applied to speech recognition[J].IEEE Trans ASSP，1975，23（1）：67－72.

[13]Sakoe H，Chiba S.Dynamic－programming algorithm optimization for spoken word recognition[J].IEEE Trans ASSP，1978，26（1）：43－49.

[14]Zhang Dabao，Huan Xiaodong，F(xiàn)red E，et al.Two－dimensional correlation optimized warping algorithm for aligning GC× GC－MS data[J].Anal Chem，2008，80（8），2664－2671.

[15]周健，成浩，曾建明，等.基于近紅外的多相偏最小二乘模型組合分析實(shí)現(xiàn)茶葉原料品種鑒定與溯源的研究[J].光譜學(xué)與光譜分析，2010，30（10）：2650－2653.

[16]張靈帥，王衛(wèi)東，谷運(yùn)紅，等.近紅外光譜的主成分分析一馬氏距離聚類判別用于卷煙的真?zhèn)舞b別[J].光譜學(xué)與光譜分析，2011，31（5）：1254－1257.

[17]劉志勇.常見氣體的聚類分析[J].Fnend of Science Amateurs，2009，11（33）：17－18.

[18]劉倩，孫培艷，高振會(huì)，等.衰減全反射傅里葉變換紅外光譜技術(shù)結(jié)合模式識(shí)別進(jìn)行油品鑒別[J].光譜學(xué)與光譜分析，2010，30（3）：663－666.

[19]Nielsen N P V，Carstensen J M，Smedsgaard J.Aligning of single and multiple wavelength chromatographic profiles for chemometric data analysis using correlation optimised warping[J].J Chromatogr A，1998，805（1－2）：17－35.