基于K均值算法的X射線熒光光譜法檢驗(yàn)藥用鋁塑包裝片的研究

劉金坤，李春宇，呂航，李飛，姜紅，滿吉

(1.中國人民公安大學(xué) 偵查學(xué)院，北京 100038；2.北京華儀宏盛技術(shù)有限公司，北京 100123)

藥用鋁塑包裝片是由鋁箔和聚氯乙烯(PVC)塑料泡罩組成的藥品包裝物，廣泛應(yīng)用于膠囊藥品和片劑藥品的包裝[1]。在犯罪現(xiàn)場勘查過程中，藥用鋁塑包裝片是一種常見物證，通過檢驗(yàn)其PVC塑料的類別，可縮小偵查范圍，為偵查和審判提供線索和依據(jù)。

法庭科學(xué)中，X射線熒光光譜法是塑料檢驗(yàn)的常用方法[2-3]。本文收集到30種不同品牌、不同廠家的藥用鋁塑包裝片，用X射線熒光光譜儀測定每個(gè)樣品的元素含量。首先按照元素是否存在對樣品分組，對不確定組別的樣品用K均值算法聚類，運(yùn)用誤差平方和法(SSE)確定最佳聚類簇?cái)?shù)K，再用輪廓系數(shù)評估K值的有效性，得到了比較準(zhǔn)確的分類結(jié)果。

1 實(shí)驗(yàn)部分

1.1 材料與儀器

不同品牌不同廠家的藥用鋁塑包裝片聚氯乙烯(PVC)30個(gè)(部分樣品見表1)；乙醇，分析純。

表1 藥用鋁塑包裝片樣品表Table 1 Medicinal aluminum-plastic packaging samples chart

X-MET8000 X射線熒光光譜儀(XRF)；大面積SDD高分辨率半導(dǎo)體探測器。

1.2 實(shí)驗(yàn)方法

藥用鋁塑包裝片表面沾有灰塵、油漬等雜質(zhì)，需要對樣品進(jìn)行簡單清洗，以保證測量數(shù)據(jù)的準(zhǔn)確。在實(shí)驗(yàn)測定之前，使用酒精棉對樣品進(jìn)行擦拭，并將處理過的樣品晾干。用X射線熒光光譜儀分別測定樣品的元素含量，每種樣品實(shí)驗(yàn)3次,取平均值。

2 結(jié)果與討論

2.1 XRF分析

X射線熒光光譜儀檢測結(jié)果見表2。

由表2可知，X射線熒光光譜測定的元素大多來源于PVC制作過程中的填料。其中，Cl元素含量最高，這是因?yàn)樗幱娩X塑包裝片的主要成分是聚氯乙烯；Sn元素來自有機(jī)錫類穩(wěn)定劑，如PVC制備過程中，為防止其熱分解，通常會加入二甲基氧化錫穩(wěn)定劑；V元素來自釩的氧化物，有催化劑的作用；Cu元素主要來自塑料加工過程中的填料硫酸銅，具有良好的凝聚性，可以去除雜質(zhì)；Zn、Ti元素主要來自常用金屬氧化物，有著色劑、增白劑的作用；Ba元素主要來自硫酸鋇,可以提高塑料熱穩(wěn)定性，具有一定的潤滑性[4-7]。通過分析元素指標(biāo)的來源，根據(jù)樣品中元素是否存在進(jìn)行分組,結(jié)果見圖1。

表2 藥用鋁塑包裝片樣品X射線熒光分析結(jié)果(μg/g)Table 2 Analysis of medicinal aluminum-plastic packaging samples by XRF

由圖1可知，按照是否還有Ba元素進(jìn)行分組時(shí)，19#樣品單獨(dú)分為一組，其他樣品分為另外一組。同理，將剩下的樣品按照是否含有Mn、Ti、Zn、V、Cu元素依次分組，最后將分組結(jié)果附上識別標(biāo)簽，初步將30個(gè)樣品分為13組。其中，No.1，2，3，6，8僅含有一個(gè)樣品，得到了準(zhǔn)確區(qū)分；No.12，13含有兩個(gè)樣品，用Pearson相關(guān)系數(shù)判斷其相關(guān)性，可以將兩個(gè)樣品區(qū)分；No.4，5，7，9，10，11含有3個(gè)及以上樣品，可根據(jù)K均值聚類法進(jìn)行分類。

圖1 藥用鋁塑包裝片樣品分組圖Fig.1 Medicinal aluminum-plastic packaging samples sub-group chart

2.2 Pearson相關(guān)系數(shù)及假設(shè)檢驗(yàn)

Pearson相關(guān)性分析是指對兩個(gè)或多個(gè)具備相關(guān)性特征元素進(jìn)行分析，通過相關(guān)系數(shù)反映出特征元素間的相關(guān)關(guān)系[8]，其表達(dá)式如下：

(1)

式中l(wèi)XX——X的離均差平方和；

lYY——Y的離均差平方和；

lXY——X、Y間的離均差積和。

根據(jù)式(1)r值可判斷元素間的相關(guān)性，當(dāng)|r|為0.00～0.19時(shí)，元素相關(guān)性極低；r=0.20～0.39時(shí)元素低度相關(guān)；0.40～0.69時(shí)中度相關(guān)；0.70～0.89時(shí)高度相關(guān)；0.90～1.00時(shí)相關(guān)性極高。

確定樣品的Pearson相關(guān)系數(shù)后，通過假設(shè)檢驗(yàn)來判斷相關(guān)系數(shù)的有效性。假設(shè)檢驗(yàn)首先提出假設(shè)，無關(guān)假設(shè)為H0，相關(guān)假設(shè)為H1。在統(tǒng)計(jì)學(xué)中，顯著性水平α通常設(shè)為0.05，當(dāng)取得r值的概率0

由表3可知，No.12，13組的Pearson相關(guān)系數(shù)都接近于1，且P值遠(yuǎn)小于0.05，相關(guān)性很強(qiáng)。通過比對鋁塑包裝片樣品表1可知，No.12組的10#和15#樣品、No.13組的9#和28#樣品分別屬于相同品牌、不同廠家的鋁塑包裝片，表明同一品牌的不用廠家生產(chǎn)鋁塑包裝片的材質(zhì)差異較小。

表3 Pearson相關(guān)系數(shù)與假設(shè)檢驗(yàn)Table 3 Pearson correlation coefficients and hypothesis tests

2.3 K均值算法

K均值算法的思想是首先選定一個(gè)K值和K個(gè)初始類簇中心點(diǎn)，將樣品分別歸到離自己最近的簇中，然后重新計(jì)算每個(gè)簇的中心點(diǎn)，通過不斷迭代，當(dāng)達(dá)到規(guī)定的迭代次數(shù)或者類簇中心點(diǎn)最小時(shí)，聚類完成[9]。通常情況下，K值的選擇有一定的不確定性，本文嘗試用誤差平方和法SSE來尋找K值[10]，SSE的表達(dá)式如下：

(2)

式中K——聚類數(shù)量；

p——聚類樣品；

mk——k個(gè)聚類的中心點(diǎn)。

由式(2)可知，隨著K值增大，每個(gè)聚類簇的聚合程度隨著增加，SSE的值慢慢減小；當(dāng)K值接近真實(shí)聚類數(shù)時(shí)，再增加K所得到的聚合程度效果會迅速變小，SSE的下降幅度會驟減；當(dāng)K值繼續(xù)增大時(shí)，SSE的下降幅度趨于平緩，那么最先趨于平緩的點(diǎn)就是合適的K值。

K均值聚類時(shí)選定No.11類藥用鋁塑包裝片樣品集(見表4)，從Python語言的sklearn工具包中調(diào)用K-Means模塊，用Pycharm Community Edition實(shí)現(xiàn)代碼運(yùn)行及數(shù)據(jù)分析[11]。

表4 No.11組藥用鋁塑包裝片樣品數(shù)據(jù)(μg/g)Table 4 No.11 medicinal aluminum-plastic packaging samples data

將 No.11組樣品數(shù)據(jù)傳入K-Means函數(shù)模塊，設(shè)定初始K值范圍1～6，調(diào)用SSE函數(shù)inertia,并調(diào)用matplotlib畫圖模塊展示SSE折線圖，見圖2。

圖2 SSE折線分布圖Fig.2 SSE line distribution

由圖2可知，當(dāng)K值為1，2，3時(shí)，SSE的下降幅度驟減；K值為3，4時(shí)，折線走勢趨于平緩，故確定最佳K值為3。隨后，用K均值聚類算法進(jìn)一步得到樣品分類結(jié)果，見表5。其中，1表示樣品屬于一類，0表示樣品不屬于一類，分類結(jié)果附上識別標(biāo)簽。5個(gè)樣品中，22#、23#、24#被分成一類，11#和17# 單獨(dú)分為一類。

表5 No.11組樣品K均值算法分類表
Table 5 No.11 samplesK-meansalgorithm classification table

標(biāo)簽樣品編號1117222324100111210000301000

2.4 輪廓系數(shù)

聚類效果的評估方法通常有輪廓系數(shù)、蘭德系數(shù)、互信息、Homogeneity、Fowlkes-Mallows scores、Calinski-Harabaz Index等[12]，其中輪廓系數(shù)較為常用，當(dāng)我們不能確定實(shí)際聚類類別時(shí)，可以通過輪廓系數(shù)來進(jìn)一步評估[13]。單個(gè)樣品點(diǎn)Xi的輪廓系數(shù)表達(dá)式如下：

(3)

式中，a為Xi與它同類別中其他樣品的平均距離；b為Xi與最近簇中所有樣品的平均距離。

通常情況下，用輪廓系數(shù)的平均值作為整個(gè)樣品集的輪廓系數(shù)值，取值范圍為[-1,1]。當(dāng)同類樣品距離相近且不同類別樣品距離越遠(yuǎn)，輪廓系數(shù)值就會增大，分類越合理。

在Python中，將K值傳入K-Means模塊，從sklearn工具包中調(diào)用元素指標(biāo)驗(yàn)證模塊metrics，再從metrics中引用silhouette-score函數(shù)。通過運(yùn)算發(fā)現(xiàn)，當(dāng)K為2時(shí)，S值為0.400,K為3時(shí)，S值為0.443，K為4時(shí)，S值為0.130。因此，當(dāng)K值為3時(shí)分類效果最好，評價(jià)結(jié)果與聚類結(jié)果相一致，說明將SSE方法用于確定K值很有效。

根據(jù)此種方法將剩下的No.4，5，7，9，10組樣品繼續(xù)分類，得到了有效的聚類結(jié)果，見表6。其中，No.4，9組樣品分為兩類時(shí)輪廓系數(shù)S值較大，分類合理；No.5，7，10組樣品分2類時(shí)S值偏小，分3類時(shí)S值更小，故將組內(nèi)樣品分為2類較為合理。

表6 樣品K均值聚類及輪廓系數(shù)評估Table 6 Sample K mean clustering and silhouette coefficient evaluation

3 結(jié)論

采用X射線熒光光譜法結(jié)合統(tǒng)計(jì)學(xué)方法，實(shí)現(xiàn)了對藥用鋁塑包裝片的準(zhǔn)確分類。實(shí)驗(yàn)用X射線熒光光譜儀測定樣品元素含量，對樣品初步分組；對于只有一個(gè)樣品的組，不再繼續(xù)分類；有兩個(gè)樣品的組可根據(jù)Pearson相關(guān)系數(shù)分類；有3個(gè)以上樣品的組根據(jù)K均值算法分類，通過SSE方法確定K值并用輪廓系數(shù)評估分類結(jié)果，最終將藥用鋁塑包裝片樣品成功分類，達(dá)到預(yù)期的實(shí)驗(yàn)結(jié)果。但本方法仍有需要改進(jìn)的地方，比如樣品量應(yīng)足夠多、模型更加簡化等。基于此，下一步將探索其他機(jī)器學(xué)習(xí)的分類算法，建立更加簡便的模型對樣品進(jìn)行分類檢驗(yàn)。