基于穩定同位素和礦物元素的牛奶產地溯源技術

王默涵，高 銘，李雅麗，柳 嘉，富宏坤，孫曉辰，段盛林,

（1.中國食品發酵工業研究院有限公司，北京 100015；2.功能主食創制與慢病營養干預北京市重點實驗室，北京 100015；3.中國合格評定國家認可中心，北京 100062）

乳業是畜牧產業發展的重要組成，我國牛乳產業雖然起步較晚，但日益興盛。牛奶是優質鈣和蛋白質的良好來源，隨著生活水平的提高和牛奶銷售量的大幅增長，保證牛奶質量是迫切所需。自2008年三聚氰胺事件后，乳品質量安全備受關注，具有特定產地、草飼散養等聲稱的高端乳制品層出不窮，但也使其更容易被混淆甚至偽造，乳源控制因而開始成為生產商和消費者的關注點，牛奶產地屬性也對確定產品價值具有重要意義[1]。

關于產地屬性問題，世界各國通過出臺法律法規或建立追溯制度來推進食品質量安全管理[2]。歐盟、美國和日本是較早開展食品追溯標準化工作的國家和地區[3]，歐盟早在1992年就已提出食品產地的概念，并在2000年要求食品的生產、加工和流通等環節建立溯源制度來明確產品來源[4]。我國在《農產品質量安全》以及《食品安全法》中也明確指出了食品農產品的追蹤與溯源的重要性，對依法建立食品安全追溯體系、保證食品的可追溯性進行了明確說明[4- 5]。

穩定同位素和礦物元素是生物天然攜帶的一種專屬標記，穩定同位素比值和礦物元素指紋信息與其生存環境和生長方式密切相關[5]，它已成功應用于葡萄酒[6]、谷物[7]、茶葉[8]、果蔬[9-10]、畜禽肉品[11-12]及水產品[13]中，是產地屬性的有效保護工具，能為食品溯源提供一種科學獨立、穩定流動的身份信息[14]，且由于受人為因素和檢測環境的影響較小而具有準確度高、靈敏度高、無污染等特點[15]。

現有的牛奶產地溯源研究主要為僅測定牛奶中的穩定同位素比值或結合同位素和十余種甚至幾十種礦物元素進行建模[5]，但若僅測定同位素，可能影響訓練模型的鑒別準確度，若納入的元素種類過多，成本和技術層面在推廣時具有一定難度。綜合考慮模型鑒別精度和成本，本研究采集黑龍江省、云南省、河北省、甘肅省、內蒙古自治區五地的牛奶樣本，以牛奶中碳、氮、氫、氧四種穩定同位素比值（δ13C、δ15N、δ2H和δ18O）和鈣、鐵、鋅三種重要礦物元素的質量濃度為指標對檢測結果進行多元統計分析，基于樣品產地建立數據模型，探究其對牛奶產地的判別能力，開發一種便于應用推廣的牛奶產地溯源技術，助力奶源控制，提振國民對國產乳品的消費信心。

1 材料與方法

1.1 材料與儀器

9 mm × 5 mm錫舟、5 mm × 3.5 mm銀囊 德國賽默飛世爾科技有限公司；鈣、鐵、鋅單元素標準液（1 mg/mL） 國家標準物質認證中心；濃硝酸（優級純） 江蘇晶銳化學股份有限公司；正己烷 天津市大茂化學試劑廠；NaN3青島雪潔助劑有限公司；實驗用水均為GB/T 6682規定的一級水。

MARS 6微波消解系統 美國CEM公司；iCAP6000電感耦合等離子體發射光譜儀、Flash EA 2000 元素分析儀、高溫裂解元素分析儀、 Delta V Advantage穩定同位素比值質譜儀 美國 Thermo Fisher公司；SL2002N電子天平 瑞士 Mettler-Toledo公司； DHG-9075A電熱干燥箱 上海一恒科技有限公司；LGJ-18真空冷凍干燥機 北京松源華興科技發展有限公司；GN-RO-40型超純水儀 北京雙峰眾邦科技發展有限公司。

1.2 實驗方法

1.2.1 樣品采集 本研究采集的牛奶樣本均為荷斯坦乳牛常乳，于2020年7月～8月分別取自我國黑龍江省、內蒙古自治區、河北省、甘肅省、云南省牧場，各省份樣本均為同一牧場不同牛區或不同產奶時間的奶罐樣，即不同胎次、泌乳期的多頭泌乳牛的混合乳，以減少其他因素干擾，每個牧場采集5個獨立樣本。各樣品均添加0.03% NaN3用于防腐，編號、登記后立即運回，直至檢測前始終儲藏于-20 ℃。采樣信息統計如表1所示。

表1 采樣信息統計Table 1 Sampling information statistics

1.2.2 穩定同位素分析

1.2.2.1 樣品前處理 參考WIJENAYAKE等[16]的方法提取牛奶中的酪蛋白進行同位素分析。將牛奶樣品冷凍干燥24 h后搗碎，準確稱取12 g樣品于50 mL離心管中，加入10 mL超純水稀釋，再加入5 mL正己烷，渦旋振蕩5 min，靜置30 min后棄去上層溶液，用0.1 mol/L稀鹽酸將離心管下層溶液調節至pH4.6，以8000 r/min的條件離心10 min，離心后取下層沉淀，加入5 mL超純水，重復用5 mL超純水洗滌2～3次，水洗后的沉淀冷凍干燥24 h后待測。

1.2.2.2δ13C和δ15N值的測定 使用元素分析儀聯用穩定同位素比值質譜儀測定樣品的δ13C和δ15N值。元素分析儀左反應管用氧化鉻填充，測樣前升溫至650 ℃；右反應管用銅絲填充，測樣前升溫至980 ℃，色譜柱溫度為50 ℃；氦氣流速為100 mL/min。分別測定參考物質IAEA-CH-6（δ13C=10.4‰±0.2‰）和USGS-34（δ15N=-1.8‰±0.2‰）的δ13C值和δ15N值，確定儀器處于穩定狀態后，稱取0.15～0.25 mg樣品粉末于錫舟中檢測。

1.2.2.3δ2H值和δ18O值的測定 使用高溫裂解元素分析儀聯用穩定同位素比值質譜儀測定樣品的δ2H值和δ18O值。高溫裂解元素分析儀反應管用玻璃碳填充，測定δ18O值前反應管階段升溫至1250 ℃，色譜柱溫度為60 ℃，氦氣流速為120 mL/min；測定δ2H值前反應管階段升溫至1400 ℃，色譜柱溫度為85 ℃，氦氣流速為120 mL/min。分別測定參考物質NBS-22（δ2H=-120.0‰±1.0‰）和USGS-34（δ18O=-27.9‰±0.6‰）的δ2H值和δ18O值，確定儀器處于穩定狀態后，稱取0.25～0.35 mg樣品粉末于銀紙杯中檢測。

1.2.2.4 同位素比值的計算 采用國際通用的相對測量法表示物質的同位素比值，碳、氫、氧、氮元素的同位素特征均以國際標準物質為參照，結果以δ（‰）表示：

式中：X為重同位素的豐度，即13C或18O 或2H或15N；R為重同位素與輕同位素的豐度比，即13C/12C 或18O/16O或2H/1H或15N/14N；R樣品為被測樣品的同位素比值；R標準為國際標準物質的同位素比值。

1.2.3 礦物元素分析 分析方法及儀器操作條件均參考GB 5009.268-2016《食品中多元素的測定》[17]。準確稱取各牛奶樣品0.50 g于微波消解內罐中，加入5 mL濃硝酸，旋緊罐蓋放置于通風櫥中酸解過夜，次日放入微波消解儀中消解，冷卻后取出，緩慢打開罐蓋排氣并用少量水沖洗內蓋，將消解罐放在控溫電熱板上于100 ℃加熱30 min趕酸排氣，取出恢復至室溫后用超純水定容至50 mL，混勻備用，同時做空白試驗。采用電感耦合等離子體發射光譜儀，通過外標法定量測定樣品中的鈣、鐵、鋅三種礦物元素的質量濃度。

1.3 數據處理

本研究中的所有指標均進行了至少3次的平行測定，結果表示為平均值±標準差，數據處理和作圖采用IBM SPSS Statistics 25。各組間通過單因素方差分析和鄧肯多重比較確定顯著性差異（P＜0.05）；所有數據標準化后，通過聚類分析（系統聚類，Ward法，歐氏平方距離）、主成分分析和判別分析（Fisher 線性判別分析，“留一法”交叉驗證）區分每個牛奶樣品的地理來源。

2 結果與分析

2.1 穩定同位素分析

穩定同位素是理想的產地溯源參數，因為其在動物組織中的分布水平與動物的飲食攝入和地理環境密切相關[18]。其中，δ13C主要與奶牛的食物成分有關，δ15N主要反映飼料種植及奶牛生存環境的土壤情況，δ2H和δ18O則與當地的水源、空氣等環境因素及氣候條件密切相關[5]。結果顯示，來自黑龍江省、內蒙古自治區、河北省、甘肅省、云南省等5個產地的牛奶樣品，其酪蛋白的δ13C、δ15N、δ2H和δ18O值地域差異性顯著，結果如表2所示。

表2 五地牛奶樣品穩定同位素比值測定結果（‰）Table 2 The stable isotope ratio in milk samples from five region(‰)

由于C3植物和C4植物光合作用的途徑不同，可能導致碳同位素比值存在差異[5]。牧草、小麥和大豆等C3植物的δ13C值范圍通常在-24‰至-34‰之間，而玉米和高粱等C4植物的δ13C范圍則通常介于-10‰至-23‰之間[19]。經奶牛攝入和消化后，這種特征可被濃縮至奶牛組織中，因此乳汁通常帶有飼料類型或比例特征[20]。MOLKENTIN等[21]的研究表明，無玉米喂養的奶牛乳汁δ13C的極限值為-23.5‰，而無任何C4植物喂養的奶牛乳汁中的δ13C值更負，草飼奶牛乳汁的δ13C值顯著低于精飼奶牛（如以玉米、高粱、青貯等為原料經精加工而成的飼料）。本研究中的δ13C值均介于-18‰至-20‰之間，表明5個采樣牧場均表現為缺乏新鮮牧草、高谷物水平的精飼料喂養模式，顯著高于同類研究中愛爾蘭、新西蘭、澳大利亞等具有豐富牧草資源國家的酪蛋白δ13C值，這些國家草飼奶牛乳汁的δ13C值通常低于-24‰[22]。具體δ13C值表現為黑龍江省和河北省＞甘肅省＞云南省和內蒙古自治區（P＜0.05），呈現東部、低海拔地區（黑龍江和河北）顯著高于西部、高海拔地區（內蒙古、甘肅、云南）的趨勢（P＜0.05）。由于采樣牧場均為以玉米、棉籽、青貯等為主要原料的精飼喂養，因此在飼料類型大致相同的情況下的差異可能與不同地域飼料作物的光合作用有關，但尚無δ13C值與海拔之間相關趨勢的直接研究。

氮同位素比值可用于區分飼料作物的施肥方式和固氮作物種類在奶牛飲食中的比例[19]。合成肥料的δ15N值通常介于-4‰至4‰之間，而有機氮肥（包括糞便和堆肥）的δ15N值范圍可寬至0～30‰[23]，部分固氮作物（如豆科植物）的δ15N值則可能接近0‰[24]。本研究中牛奶酪蛋白具有相似的δ15N 值，均在2‰至4‰之間，且內蒙古自治區和河北省＞黑龍江省和云南省＞甘肅省（P＜0.05），均顯著低于愛爾蘭、新西蘭、澳大利亞等國家[22]，這可能是由于奶牛飲食中合成氮肥的使用比例較大而有機氮肥比例較小，且飼料中豆粕等固氮作物的比例較高，差異主要反映在耕作及喂養模式上，與經緯度、海拔的相關性較小。

δ2H值和δ18O值的結果趨勢表現出高度的一致性：甘肅省＞河北省和內蒙古自治區＞黑龍江省和云南省，該趨勢也與δ15N值表現出一定的一致性，這可能與氣候條件、距海遠近、海拔、空氣、水源等養殖環境有關[18]，但與愛爾蘭、新西蘭、澳大利亞等國家及我國的其他同類數據有較大差距[1,22]，這可能與不同的采樣時間和采樣地點有關。

可見，穩定同位素比值是鑒別養殖模式（草飼、精飼等）的有效參數，也為地理溯源提供可靠依據。

2.2 礦物元素分析

礦物元素常具有豐度效應，不同地域間土壤、水源及當地奶牛飼料中的礦物元素含量可能存在差異，故地域可能影響乳中礦物元素含量[25]。ZHOU等[26]評估了我國家畜乳中元素含量在地域間的差異，結果表明，不同地域家畜乳中元素含量不同，飲水與飼料是可能的影響元素。本研究通過單因素方差分析和鄧肯多重比較（見表3）顯示出牛奶中鈣和鋅的質量濃度在5個產地的分布差異顯著（P＜0.05），且分布規律一致，均在甘肅省最高、內蒙古自治區最低；而牛奶樣品中鐵元素的含量較低，在5個產地的分布差異不顯著（P＞0.05）。鑒于采樣時已排除個體差異的影響，此結果反映了不同養殖環境下差異化的營養狀態，表明其對區分五地區的樣品具有一定的可行性。

表3 五地牛奶樣品礦物元素質量濃度測定結果Table 3 The concentration of mineral elements in milk samples from five region

2.3 聚類分析

聚類分析是歸類劃分數據中相似性數據的統計方法[27]。使用以上四項同位素比值及三項礦物元素質量濃度的標準化得分進行聚類分析，產生的聚類樹狀圖如圖1所示。將聚類距離為5時作為分割線，五產地的牛奶樣品可被分為4類。其中甘肅省的樣品與其他省距離較遠，可被完全分離；河北省的樣品雖同屬一類，聚類效果較好，但其他產地的個別樣品也分類至這類中；云南省和黑龍江省的大部分樣品被歸為一類，表明兩產地的聚類距離相近，無法完全分離；內蒙古自治區的樣品聚類效果不好，分散在各類別之中。可見，聚類分析為產地鑒別提供了可能，但仍存在交叉現象。

圖1 不同產地牛奶的系統聚類樹狀圖Fig.1 Dendrogram of cluster analysis from different regions

2.4 主成分分析

利用穩定同位素比值和礦物元素組合參數的主成分分析（PCA）可對不同產地的牛奶進行分類[5]。將所有數據標準化后，KMO檢驗統計量大于0.6，Bartlett球形檢驗顯著性P＜0.05，表明本實驗模型適合做主成分提取[28]。PCA共提取了7個主成分，總方差解釋如表4所示，其中初始特征值大于1的主成分有3個，累積方差百分比為81.95%，可代表大部分數據特征。以前三個主成分為代表的模型為：

表4 總方差解釋Table 4 Total variance explained

利用此三個主成分的標準化得分值做散點圖，結果如圖2所示，甘肅省的樣品與其他省區的差異明顯，可進行區分，其他地區樣品仍有交叉性。

圖2 不同產地牛奶的主成分得分散點圖Fig.2 Scatter plot of principal component score from different regions

2.5 判別分析

判別分析是一種有監督的數據分類方式[29]。其可利用各類樣品之間的相關性建立判別模型而對未知樣品進行分類，還可利用交叉驗證把除某一個樣品之外的其余樣品作為訓練集，建立判別函數對該樣品進行判別，從而極大地檢驗了判別模型的穩定性[30]。基于牛奶產地進行元素差異分析、聚類分析和主成分分析后可以看出，本研究的組合參數對牛奶產地溯源具有一定可行性。為進一步了解各參數對牛奶產地的判別效果，利用Fisher函數建立判別模型[31]。

本次判別分析建立了4個典型判別式函數，其中函數1的方差貢獻率為96.1%，函數2的方差貢獻率為2.8%，利用這2個判別函數分別作為橫縱坐標建立坐標系，得到5個組的散點圖，如圖3所示，可以看出各產地樣品相對集中，根據產地分為五組后各組質心也可較好地分離。

圖3 不同產地牛奶的判別分析散點圖Fig.3 Scatter plot of discriminant analysis from different regions

采用交叉驗證檢驗模型的可靠性，結果如表5所示。原始分類結果顯示：7種元素分析的原始整體判別率為100%，五產地的樣品均可被100%正確識別；“留一法”交叉驗證的分類結果中，有96%的牛奶樣本被正確分類，其中黑龍江省的樣本存在20%的誤判，其余省區的樣本均可被100%正確判別。說明本次判別分析建立的判別函數對來自5個地區的牛奶樣品具有很好的分類效果。所建立的判別函數模型為：

表5 不同產地牛奶的判別分析分類結果（%）Table 5 Classification result of discriminant analysis from different regions (%)

同類的牛奶產地溯源研究中，ZHAO等[1]針對δ13C、δ15N、δ2H、δ18O四種穩定性同位素比值建立數據模型，得到河北省、寧夏回族自治區、陜西省、內蒙古四個區域的交叉驗證準確率71.55%，河北省、天津市、江蘇省、內蒙古自治區四個區域的交叉驗證準確率92.11%；侯建波等[5]除穩定性同位素外，納入了40種礦物元素，選出As、Se、Rb、Sr、Mo、Cs、Ba、Tl、δ13C、δ15N和δ2H組合指標，對德國、西班牙、意大利、澳大利亞、新西蘭、中國和奧地利七個區域的牛奶整體判別率為 92.4%。本研究的分析方法與僅納入穩定同位素的數據模型相比判別準確率提高，同時又控制了成本投入，無需檢測過多礦物元素即可實現國內多省市區牛奶樣本的判別，是一種更適合應用推廣的方法。

3 結論

本研究分別測定黑龍江省、云南省、河北省、甘肅省、內蒙古自治區五個產地的牛奶樣品的穩定同位素及礦物元素指紋信息，建立的數據模型可對奶源地喂養模式和牛奶產地實現有效判別。結果表明：δ13C、δ15N、δ2H、δ18O在產地間的差異顯著（P＜ 0.05），這可能與地理屬性和養殖模式有關，可區分如草飼、精飼等的喂養模式，也可為地理溯源提供可靠依據；鈣、鋅在產地間的差異顯著（P＜0.05），而鐵在產地間的差異不顯著（P＞0.05），表現出不同養殖環境差異化的營養狀態，是地理溯源的有力指標。基于產地的系統聚類和主成分分析對牛奶產地溯源具有一定可行性，但仍存在部分交叉；Fisher線性判別分析則可有效判別牛奶產地，原始整體判別率為100%，交叉驗證的正確判別率為96%，其中對內蒙古自治區的正確判別率為80%，其余均為100%。本研究結果證明利用穩定同位素碳、氮、氫、氧與鈣、鐵、鋅相結合的溯源技術可對我國黑龍江省、云南省、河北省、甘肅省、內蒙古自治區五個產地的牛奶進行有效的產地判別，是產地溯源的有利技術。