王元忠，趙艷麗，張 霽，金 航

云南省農(nóng)業(yè)科學院藥用植物研究所，云南 昆明 650200

近紅外光譜信息篩選在瑪咖產(chǎn)地鑒別中的應用

王元忠，趙艷麗，張 霽，金 航*

云南省農(nóng)業(yè)科學院藥用植物研究所，云南 昆明 650200

食藥植物瑪咖富含多種營養(yǎng)成分，極具藥用價值。采用近紅外漫反射光譜，對采自秘魯及云南共139份瑪咖樣品進行產(chǎn)地鑒別。采用多元信號校正結(jié)合二階導數(shù)和Norris平滑預處理光譜，利用光譜標準偏差初步選擇光譜波段(7 500～4 061 cm-1)，結(jié)合主成分-馬氏距離(principal component analysis-mahalanobis distance，PCA-MD)篩選出適宜的主成分數(shù)為5?；谒Y選的光譜波段及主成分數(shù)，利用“模群迭代奇異樣本診斷”方法剔除2個異常樣品后，分別采用競爭自適應重加權(quán)法(competitive adaptive reweighted sampling，CARS)、蒙特卡洛-無信息變量消除法(monte carlo-uninformative variable elimination，MC-UVE)、遺傳算法(genetic algorithm，GA)和子窗口重排(subwindow permutation analysis，SPA)四種方法篩選光譜變量信息，利用模型集群分析(model population analysis，MPA)思想對所篩選的光譜變量信息進行評價。結(jié)果顯示，RMSECV(SPA)>RMSECV(CARS)>RMSECV(MC-UVE)>RMSECV(GA)，分別為2.14，2.05，2.02，1.98，光譜變量數(shù)分別為250，240，250和70。采用偏最小二乘判別分析法(partial least squares discriminant analysis，PLS-DA)對四種方法篩選的光譜變量建立判別模型，隨機選擇97份樣品作為建模集，其余40份樣品作為驗證集。通過R2，RMSEC和RMSEP分析可知，R2： GA>MC-UVE>CARS>SPA，RMSEC和RMSEP： GA

瑪咖； 近紅外光譜； 鑒別； 光譜信息篩選； 模型集群分析

引 言

瑪咖(LepidiummeyeniiWalp.)為十字花科(Cruciferae)獨荇菜屬(Lepidium)植物，常用名Maca(音譯“瑪咖”)，又名Maka，Maca-maca，Peruvian ginseng，Maino，Ayakwillku，Ayakchichira等，原產(chǎn)于海拔3 500～4 500 m的南美洲安第斯山區(qū)，現(xiàn)主要分布于秘魯中部的Puna生態(tài)區(qū)和秘魯東南部城市Puno[1-2]?，斂K根可藥食兩用，極具營養(yǎng)價值和經(jīng)濟價值，富含蛋白質(zhì)、人體必需氨基酸、礦物質(zhì)和多種維生素，在安第斯山區(qū)已有悠久的使用歷史。研究發(fā)現(xiàn)瑪咖中除含有多種營養(yǎng)成分外[3]，還含有瑪咖酰胺、瑪咖烯、芥子油苷及其衍生物、甾醇等多種活性物質(zhì)[4-6]?，F(xiàn)代藥理學研究表明瑪咖具有抗氧化[7]、增強生育力[8]、改善性功能[9]、抗疲勞[10]、調(diào)節(jié)內(nèi)分泌[11]、增強免疫力[12]、抗壓力[13]、緩解更年期綜合癥[14]等多種功效。近年來，瑪咖因其獨特功效而倍受關(guān)注。我國云南、新疆、吉林、西藏和湖南等省(區(qū))均有引種栽培，其中以云南省種植規(guī)模最大、品種最多。據(jù)統(tǒng)計2012年，云南全省瑪咖種植面積達22 000 hm2，已成為云南省高寒山區(qū)的重要經(jīng)濟作物[15]。然而，不同產(chǎn)地種植的瑪咖其營養(yǎng)價值和藥用價值存在一定差異，且隨著市場需求逐年增加，市售瑪咖品質(zhì)差異較大，粉末樣品摻雜現(xiàn)象突出。因此，不同產(chǎn)地瑪咖樣品的鑒別，為后續(xù)指導瑪咖GAP種植、品質(zhì)及安全評價具有重要的作用。

近紅外光譜(NIRS)對復雜體系中C—H，O—H，N—H等含氫基團振動的倍頻和合頻吸收極為敏感，全光譜具有同時反映物料的物化信息等特點。但近紅外全光譜變量中，很多是無效信息、甚至是干擾變量，引入此類變量建立判別模型增加了模型的復雜性，同時降低了模型的判別精度，進而對樣本的鑒別和變量的計算產(chǎn)生重要影響[16]。為建立較為準確的判別模型，需對光譜信息進行選擇及優(yōu)化，使所篩選的光譜信息具有全面性和客觀性。收集了139份不同產(chǎn)區(qū)種植的瑪咖樣品，擬采用競爭自適應重加權(quán)法(CARS)[17]、蒙特卡洛-無信息變量消除法(MC-UVE)[18-19]、遺傳算法(GA)[20]和子窗口重排(SPA)[21]篩選瑪咖樣品的近紅外光譜變量后，利用模型集群分析(MPA)[22]思想對所篩選的光譜變量信息進行評價，并確定篩選的光譜波數(shù)，進一步采用偏最小二乘判別分析法(PLS-DA)建立各方法的判別模型，用于預測驗證樣品。該方法旨在為探析構(gòu)建近紅外光譜判別模型的新方法提供參考，通過集群分析結(jié)合PLS-DA建立不同產(chǎn)地瑪咖樣品的判別模型，以期為近紅外光譜鑒別和評價中藥材產(chǎn)地、品種等提供依據(jù)。

1 實驗部分

1.1 材料

2014年采集15個不同產(chǎn)地種植的藥用植物瑪咖139份，經(jīng)鑒定，樣品來源見表1。

Table 1 Source and number of materials

1.2 儀器

AntarisⅡ近紅外光譜儀(美國賽默飛世爾科技，配備InGaAs檢測器的漫反射模塊，Result 2.1軟件采集光譜圖、TQ 8.6軟件定性分析)； DFT-100型中藥粉碎機(浙江溫嶺市林大機械有限公司)，80目不銹鋼篩盤(北京，中西泰安)； SIMCA-P+11.0軟件(瑞典，UMETRICS)； MATLAB R2010a分析軟件，代碼來源于http://code.google.com/p/carspls/和http://www.mathworks.cn/。

1.3 樣品制備

樣品采集后用自來水洗凈根莖部位，再用蒸餾水沖洗3次，將其切成薄片，置于干凈的白紙上在常溫下陰干。用中藥粉碎機粉碎，過80目篩，樣品不少于20.0 g，存放于自封袋中，備用。

1.4 近紅外光譜采集

將樣品置于平衡箱中平衡2 h，使樣品水分含量在10%～12%。準確稱取平衡后的樣品20.0 g，充分混勻，置于近紅外采樣杯中，壓緊。在儀器預熱2 h后，使用Result 2.1軟件漫反射模塊采集近紅外光譜圖。測試條件： 掃描次數(shù)： 64次，分辨率： 4 cm-1，掃描范圍： 10 000～4 000 cm-1。樣品平行測定3次，取平均光譜。隨機選擇99份樣品作為訓練集，其余40份樣品作為驗證集。采用TQ 8.6軟件對光譜進行優(yōu)化處理，消除基線漂移和噪音，結(jié)果見圖1和圖2。

Fig.1 Original NIRS spectra of Maca

Fig.2 Second derivative spectra of NIRS of Maca

2 數(shù)據(jù)處理

2.1 近紅外光譜信息初步篩選[23]

利用TQ 8.6軟件，使用MSC+SD+ND(13∶3)方法初步優(yōu)化所采集的光譜，對不同產(chǎn)地樣品類型分別賦值： 香格里拉1#“1”，香格里拉2#“2”，香格里拉3#“3”至轎子雪山13#“13”，大山包14#“14”，秘魯15#“15”。通過光譜標準偏差選擇光譜波段，利用主成分-馬氏距離(PCA-MD)初步建立分類模型，并對主成分數(shù)進行選擇，結(jié)果見圖3和表2。

2.2 近紅外光譜集群分析

2.2.1 奇異樣品診斷

由于光譜數(shù)據(jù)(X)和量測指標(Y)的采集和測定中，儀器的穩(wěn)定性及異常波動會使部分X或Y偏離整體分布的數(shù)據(jù)點，導致所建立的判別模型存在奇異點。為提高判別模型的預測精度，建立穩(wěn)健的判別模型，采用SIMCA-P+11.0軟件將2.1中所篩選的光譜波段(7 500～4 061 cm-1)轉(zhuǎn)置為數(shù)據(jù)格式，并利用“模群迭代奇異樣本診斷”方法通過多次循環(huán)校正，將所選擇的光譜數(shù)據(jù)(X=7 500～4 061 cm-1)與量測指標(Y=產(chǎn)區(qū)分類賦值1～15)建立判別模型，并計算模型的偏度和峰度，結(jié)果見圖4。

Table 2 Contribution of principal component

Fig.3 Spectral range selected by standard deviation

Fig. 4 Singular sample diagnosis

2.2.2 集群分析優(yōu)化光譜變量

剔除2.2.1中異常樣本1—2和12—6后，選擇主成分數(shù)為5，對光譜數(shù)據(jù)與產(chǎn)地分類指標進行判別分析。將建模集樣品和驗證集樣品按0.8的比例(即Ratio=0.8)，分別利用CARS(運行次數(shù)為N=100次，主成分數(shù)A=5，交叉驗證數(shù)K=10，統(tǒng)計每個光譜波數(shù)選中頻率)、MC-UVE、GA(MC-UVE運行次數(shù)為N=100次，主成分數(shù)A=5； GA運行次數(shù)為N=100次，Autoscaling=2，Deletion groups=5，Chromosomes=30，Variables as a maximum=30，Probability of mutation=0.01, Probability of cross-over=0.5，統(tǒng)計每個光譜波數(shù)選中頻率)和SPA(主成分數(shù)p.component=5，顯著性水平p.sig=0.01，載荷p.Q=0.3)優(yōu)化光譜波數(shù)，計算標表征變量重要性的指標，CARS，MC-UVE，GA和SPA分別為頻率(Frequency，F(xiàn))、變量穩(wěn)定指數(shù)(reliability index，RI)[17]、頻率和條件協(xié)同得分(conditional synergetic score，COSS)[22]，按降序排列，以蒙特卡洛采樣技術(shù)，采用向前推移選擇變量重要性指標引入變量(最大變量數(shù)為250)重復建模，使用蒙特卡洛交互驗證評價預測性能[22]，以RMSECV值最小或變化趨是較小時確定最佳變量數(shù)，結(jié)果見圖5。

Fig.5 RMSECV change trend chart

采用上述四種方法對所篩選的光譜波數(shù)進行分析，結(jié)果見圖6—圖10。由圖6—圖10可知，CARS法選中的頻率大于1%的光譜波數(shù)分布于少數(shù)幾個波段，其中主要分布于4 072～4 539，4 917～4 933，5 191～5 993，7 046～7 116 cm-1等波段，選中的光譜變量數(shù)為240個； GA法選中的頻率大于4%的光譜波數(shù)分布均勻，主要集中于4 315～7 254 cm-1，選中的光譜變量數(shù)為70個； MC-UVE法選中的變量穩(wěn)定指數(shù)(RI)大于3的光譜波數(shù)主要分布于4 396～7 170 cm-1，選中的光譜變量數(shù)為250個； SPA法選中的條件協(xié)同得分(COSS)大于1的光譜波數(shù)主要分布于4 072～6 981 cm-1，選中的光譜變量數(shù)為240個。

Fig.6 Spectrum wave number selected by CARS

Fig.7 Spectrum wave number selected by GA

Fig.8 Spectrum wave number selected by MC-UVE

Fig.9 Spectrum wave number selected by SPA

2.3 PLS-DA判別模型的建立

采用SIMCA-P+11.0軟件對2.2.2中四種方法優(yōu)化的光譜波數(shù)及其分類賦值建立判別模型，剔除奇異樣品后，訓練集為97份，驗證集為40份，分類賦值見2.2.1。以GA法為例，建立判別模型，用于預測驗證集樣品。采用決定系數(shù)(R2)、校準均方根誤差(RMSEC)和預測均方根誤差(RMSEP)三個參數(shù)評估預測模型的準確性[23]，結(jié)果見圖10和表3。

Fig.10 3D plot of PLS-DA by GA analysis

Table 3 Results of PLS-DA models for discrimination of different geographical Maca samples by GA, MC-UVE, CARS and SPA analysis

續(xù)表3

11-4111110.5720.30311UI11.5730.405111110.8470.108111110.8270.12211-9111111.0000.000111111.4070.288111110.4730.373111110.6530.24612-2121211.8600.099121212.4210.297121211.9260.052121212.1020.07212-5121212.4780.338121212.4320.30612UI12.5690.40212DE12.7960.56313-2131313.1500.106131313.2440.172131312.5470.320131313.0680.04813-6131312.6860.222131313.2230.15813UI12.3070.490131312.6690.23413-9131312.6330.260131313.2280.161131312.6680.235131312.8670.09414-2141413.5980.284141414.2670.18914DE13.2380.539141413.4540.38614-4141414.0000.000141413.7910.148141413.6320.26014UI13.3640.45014-8141414.1290.091141413.5330.330141413.5250.336141413.7220.19715-1151514.8850.082151514.9680.023151514.6440.252151514.8910.07715-5151514.5720.302151515.2650.187151515.0630.045151515.1280.09115-8151514.8020.140151515.1600.113151514.5430.323151514.5860.293正確率/%95.0092.5090.0085.00R20.99620.99380.98900.9864RMSEC0.180.230.310.34RMSEP0.280.350.470.52

Note： AC: actual class; CC: calculated class; Ypre: predicted value; Ydev: deviation; UI: unidentified; DE: discriminant error

3 結(jié)果與討論

3.1 近紅外光譜初步篩選分析

由圖3和表2可知，15個不同產(chǎn)區(qū)139份樣品建立的PCA-MD分類模型，光譜波段選擇7 500～4 061 cm-1，主成分數(shù)為5時，提取的光譜信息較為充分，提取的總光譜貢獻率為98.192 7%，選擇的7 500～4 061 cm-1光譜波段提取貢獻率為98.9995%，且隨主成分數(shù)的增加，光譜信息數(shù)據(jù)量無明顯變化。因此，選擇主成分數(shù)為5較為適宜。

3.2 近紅外光譜集群分析

由圖4分析可知，共剔除奇異點(outlier)數(shù)2個，即樣品1—2和12—6，計算出模型的偏度和峰度分別為-0.31和-0.52。結(jié)果表明，樣品1—2和12—6可能在所篩選的光譜波段與其余樣品存在較大差異，致使樣品出現(xiàn)奇異。通過模型偏度和峰度分析可知，在所選擇的光譜波段范圍內(nèi)，分析樣品較為集中，說明該光譜波段對模型的建立具有重要影響。

由圖5可知，采用四種方法對光譜波數(shù)進行篩選，RMSECV(SPA)>RMSECV(CARS)>RMSECV(MC-UVE)>RMSECV(GA)，分別為2.14，2.05，2.02，1.98； 且變量個數(shù)分別為250，240，250和70。分析結(jié)果顯示： 四種集群分析方法中，通過綜合評價選擇的光譜波數(shù)和RMSECV值，其中GA法最優(yōu)，CARS法和MC-UVE法次之，SPA法稍差。

由圖6—圖10分析可知，不同產(chǎn)區(qū)采集的瑪咖樣品在近紅外光譜分析中，其光譜差異波段主要集中在4 300～7 200 cm-1之間，且光譜波數(shù)具有不連續(xù)性，四種方法篩選瑪咖樣品近紅外光譜波數(shù)均具有一定的代表性，其中GA法篩選的光譜波數(shù)較少，且分布范圍較廣，具有較好的代表性。

3.3 PLS-DA判別模型分析

由3D圖10可知，秘魯15#、香格里拉2#、昭通5#、大山包14#、轎子雪山13#等能清晰分開，其余產(chǎn)地樣品雖集聚在一起，但仍能清晰區(qū)分。由表3可知，通過上述四種方法篩選光譜波數(shù)，采用PLS-DA對篩選的光譜波數(shù)與其分類賦值進行分析，并對驗證樣品進行預測。結(jié)果顯示： 四種方法篩選的光譜波數(shù)建立的判別模型預測驗證樣品正確率均大于85%，其中GA方法驗證結(jié)果較好，達95%。通過R2、RMSEC和RMSEP分析可知，GA(R2)>MC-UVE(R2)>CARS(R2)>SPA(R2)； GA(RMSEC，RMSEP)

4 結(jié) 論

通過采集不同產(chǎn)地瑪咖樣品的近紅外光譜，采用TQ 8.6軟件，應用光譜標準偏差初步選擇7 500～4 061 cm-1光譜波段，應用PCA-MD初步建立產(chǎn)地分類模型，提取光譜信息的主成分貢獻率選擇所提取的主成分數(shù)為5，進一步采用“模群迭代奇異樣本診斷”方法通過多次循環(huán)校正診斷，剔除異常樣品1—2和12—6，使用模型偏度和峰度初步判定所選光譜波段的實用性。

為使鑒定效果更好，采用MatlabR2010a分析軟件，分別利用GA，CARS，MC-UVE和SPA四種方法對初步篩選的光譜波段及其分類賦值樣品進行分析，優(yōu)化不同產(chǎn)地咖樣品的光譜波數(shù)，并通過RMSECV和篩選的光譜波數(shù)評估方法的可靠性。采用SIMCA-P+11.0軟件結(jié)合PLS-DA對四種方法篩選的光譜波數(shù)建立判別模型，利用R2，RMSEC和RMSEP評價模型的預測性能。結(jié)果顯示： 四種方法篩選的光譜波數(shù)建立的判別模型預測驗證樣品正確率均大于85%，其中GA方法驗證結(jié)果較好，達95%，GA(R2)>MC-UVE(R2)>CARS(R2)>SPA(R2)； GA(RMSEC，RMSEP)

[1] YU Long-jiang, JIN Wen-wen, WU Yuan-xi, et al(余龍江, 金文聞, 吳元喜, 等). Natural Product Research and Development(天然產(chǎn)物研究與開發(fā)), 2002, 14(5): 71.

[2] YU Long-jiang, SUN You-ping, CHENG Hua, et al(余龍江, 孫友平, 程 華, 等). Acta Botanica Boreali-Occidentalia Sinica (西北植物學報), 2004, 24(10): 1901.

[3] WANG Yi-qiang, CHEN Zhang-jing, WANG Qi-ye, et al(王義強, 陳章靖, 王啟業(yè), 等). Nonwood Forest Research(經(jīng)濟林研究), 2014, 32(2): 167.

[4] Zheng B L, He K, Rogers L, et al. Urology, 2000, 55(4): 598.

[5] Li G, Ammermann U, Quiros C F. Economic Botany, 2001, 55(2): 255.

[6] Piacente S, Carbone V, Plaza A, et al. Journal of Agricultural and Food Chemistry, 2002, 50(20): 5621.

[7] Sandoval M, Okuhama N N, Angeles F M, et al. Food Chemistry, 2002, 79(2): 207.

[8] Onyechi O, Lawrence U S E, Akuoma H O. Nutrition Research, 1999, 19(3): 443.

[9] McKay D. Alternative Medicine Review, 2004, 9(1): 4.

[10] Schroeck F R, Hollingsworth J M, Hollenbeck B K, et al. Urology, 2013, 81(6): 1177.

[11] Meissner H O, Reich-Bilinska H, Mrozikiewicz R, et al. Menopause, 2005, 12(6): 813.

[12] ZHANG Yong-zhong, YU Long-jiang, WAN Jun-mei, et al(張永忠, 余龍江, 萬軍梅, 等). Natural Product Research and Development(天然產(chǎn)物研究與開發(fā)), 2007, 19(2): 274.

[13] Rubio J, Riqueros M I, Manuel G, et al. Food and Chemical Toxicology, 2006, 44: 1114.

[14] Brooks N A, Wilcox G, Walker K Z, et al. Menopause, 2008, 15(6): 1.

[15] Yang S H, Li G Z, Xue R G, et al. Agricultural Science and Technology, 2013, 14(12): 1877.

[16] RUAN Zhi-gang, LI Bin(阮治綱，李 彬). Chinese Journal of Pharmaceutical Analysis(藥物分析雜志)，2011，31(2)： 408.

[17] Li H D, Liang Y Z, Xu Q S, et al. Analytica Chimica Acta, 2009, 648(1): 77.

[18] Cai W, Li Y, Shao X. Chemometrics and Intelligent Laboratory Systems, 2008, 90(2): 188.

[19] Han Q J, Wu H L, Cai C B, et al. Analytica Chimica Acta, 2008, 612(2): 121.

[20] Leardi R. Journal of Chemometrics, 2000, 14(5-6): 643.

[21] LIANG Yi-zeng, XU Qing-song(梁逸曾, 許青松). Instrumental Analysis of Complex Systems-White, Gray and Black Analytical Systems and Their Multivariate Methods(復雜體系儀器分析-白、灰、黑分析體系及其多變量解析方法). Beijing: Chemical Industry Press(北京: 化學工業(yè)出版社), 2012. 487.

[22] Li H D, Liang Y Z, Xu Q S, et al. Journal of Chemometrics, 2009, 24(7-8): 418.

[23] ZHAO Yan-li, ZHANG Ji, YUAN Tian-jun, et al(趙艷麗, 張 霽, 袁天軍, 等). Spectroscopy and Spectral Analysis(光譜學與光譜分析), 2014, 34(7): 1831.

Study on Application of NIR Spectral Information Screening in Identification of Maca Origin

WANG Yuan-zhong, ZHAO Yan-li, ZHANG Ji, JIN Hang*

Institute of Medicinal Plants, Yunnan Academy of Agricultural Sciences, Kunming 650200, China

Medicinal and edible plant Maca is rich in various nutrients and owns great medicinal value. Based on near infrared diffuse reflectance spectra, 139 Maca samples collected from Peru and Yunnan were used to identify their geographical origins. Multiplication signal correction (MSC) coupled with second derivative (SD) and Norris derivative filter (ND) was employed in spectral pretreatment. Spectrum range (7 500～4 061 cm-1) was chosen by spectrum standard deviation. Combined with principal component analysis-mahalanobis distance (PCA-MD), the appropriate number of principal components was selected as 5. Based on the spectrum range and the number of principal components selected, two abnormal samples were eliminated by modular group iterative singular sample diagnosis method. Then, four methods were used to filter spectral variable information, competitive adaptive reweighted sampling (CARS), monte carlo-uninformative variable elimination (MC-UVE), genetic algorithm (GA) and subwindow permutation analysis (SPA). The spectral variable information filtered was evaluated by model population analysis (MPA). The results showed that RMSECV(SPA)>RMSECV(CARS)>RMSECV(MC-UVE)>RMSECV(GA), were 2.14, 2.05, 2.02, and 1.98, and the spectral variables were 250, 240, 250 and 70, respectively. According to the spectral variable filtered, partial least squares discriminant analysis (PLS-DA) was used to build the model, with random selection of 97 samples as training set, and the other 40 samples as validation set. The results showed that,R2: GA>MC-UVE>CARS>SPA, RMSEC and RMSEP: GA

Lepidium meyenii Walp.; NIR spectroscopy; Identification; Spectral information screening; Model population analysis

Jul. 11, 2014; accepted Nov. 12, 2014)

2014-07-11，

2014-11-12

國家自然科學基金項目(31460538, 81260608)和云南省自然科學基金項目(2013FD066, 2013FZ150)資助

王元忠，1981年生，云南省農(nóng)業(yè)科學院藥用植物研究所助理研究員 e-mail： yzwang1981@126.com *通訊聯(lián)系人 e-mail： jinhang2009@126.com

O657.3

A

10.3964/j.issn.1000-0593(2016)02-0394-07

*Corresponding author