基于Matlab的中藥主成分分析數(shù)學(xué)實驗教學(xué)案例

董鴿 閔建中 陳立范 王宏杰

摘?要：本文探討醫(yī)學(xué)院校數(shù)學(xué)實驗課教學(xué)內(nèi)容與專業(yè)相結(jié)合的教學(xué)案例。以中藥專業(yè)為例，介紹了中藥專業(yè)中常用的主成分分析的統(tǒng)計方法，對10批三黃片樣品的成分含量進行主成分分析，運用Matlab7.0軟件進行計算，給出各主成分與原始成分的線性表達式和因子載荷陣，分析各主成分與原始成分的關(guān)聯(lián)程度以及各主成分在原始成分中所產(chǎn)生的作用，最后運用綜合評價函數(shù)計算10批樣品的F值，根據(jù)F值大小進行排序，對10批樣品進行綜合評價。

關(guān)鍵詞：數(shù)學(xué)實驗;案例分析;Matlab;主成分分析

中圖分類號：G4?????文獻標(biāo)識碼：A??????doi：10.19311/j.cnki.1672-3198.2021.34.068

高等數(shù)學(xué)是許多醫(yī)學(xué)專業(yè)的必修課程，但理論性強，與專業(yè)知識銜接不夠緊密。數(shù)學(xué)實驗課程是聯(lián)系高等數(shù)學(xué)理論知識與專業(yè)實際應(yīng)用的橋梁，提高學(xué)生對高等數(shù)學(xué)知識的應(yīng)用能力和計算機技術(shù)應(yīng)用能力，培養(yǎng)學(xué)生的科研精神、創(chuàng)新意識和實際操作能力。在醫(yī)學(xué)院校中開始數(shù)學(xué)實驗課，課程內(nèi)容的選擇非常重要，不同的專業(yè)對數(shù)學(xué)知識的需求也不盡相同，因此數(shù)學(xué)實驗的內(nèi)容也要跟著進行相應(yīng)調(diào)整。以中藥專業(yè)為例，我們在基礎(chǔ)知識的基礎(chǔ)上，選擇了與中藥數(shù)理統(tǒng)計的案例，由于篇幅限制，我們僅舉例主成分分析模塊的數(shù)學(xué)實驗教學(xué)案例。

1?主成分分析的步驟

主成分分析是一種降維數(shù)據(jù)處理的方法，在人臉數(shù)據(jù)識別、基因數(shù)據(jù)列分析、食品成分分析、中藥成分分析等方面有著重要應(yīng)用。主成分分析的具體步驟如下：

（1）設(shè)原始數(shù)據(jù)的樣本個數(shù)為n，有m個觀測指標(biāo)X₁，X₂，…，X_m，其中X_j=（x_1j，x_2j，…，x_nj）^T，j=1，…，m，記矩陣X=（X₁，X₂，…，X_m），欲尋找可以概括m個觀測指標(biāo)綜合信息的綜合指標(biāo)Z₁，Z₂，…，Z_m，滿足線性組合。

Z_i=a_1iX₁+a_2iX₂+…+a_miX_m，i=1，…，m，（1）

其中：

a_1i，a_2i，…，a_mi為常數(shù)，i=1，…，m。

（2）利用公式。

（6）主成分的個數(shù)可以按照以下方法確定。

方法一：若主成分的特征值λ_i≥1，則保留該主成分Z_i，否則就去掉。

方法二：若前p個主成分的累計貢獻率達到某一特定的值時（例如，以大于90%），則保留前p個主成分。

（7）設(shè)已經(jīng)確定的主成分為Z₁，…，Z_p，選擇每個主成分的貢獻率c_i作為權(quán)重，構(gòu)造綜合評價函數(shù)。

F=c₁Z₁+c₂Z₂+…+c_pZ_p（3）

利用公式（1）得到這p個主成分的得分，再代入公式（3）求得每個樣品的F值。一般地，F(xiàn)值越大，表明該樣品的綜合評價效果越好，當(dāng)然還要根據(jù)各主成分的專業(yè)意義而定。

2?三黃片成分含量的主成分分析案例

三黃片具有清熱解毒，瀉火通便的功效。文獻利用HPLC 法分析了測定了不同廠家10 批三黃片的大黃素、大黃酚、鹽酸小檗堿和黃芩苷的成分含量，運用主成分分析和聚類分析的方法對10批樣品的成分含量進行分析評價，為三黃片質(zhì)量控制提供了依據(jù)。在運用主成分分析法時使用的是SPSS 24.0軟件，只給出了主成分特征值和貢獻率表格以及主成分得分圖。本文將采用Matlab7.0軟件對文獻中10批樣品的成分含量進行主成分分析，作為數(shù)學(xué)實驗課教學(xué)的案例，給出各主成分與原始成分的線性表達式和因子載荷陣，分析各主成分與原始成分的關(guān)聯(lián)程度以及各主成分在原始成分中所產(chǎn)生的作用，并運用綜合評價函數(shù)計算10批樣品的F值，根據(jù)F值大小進行排序，對10批樣品進行綜合評價，反映藥品質(zhì)量的高低。

我們采用文獻的表1中來自不同廠家的10批三黃片的大黃素、大黃酚、大黃素和大黃酚總量、鹽酸小檗堿、黃芩苷五個成分含量測定數(shù)據(jù)作為原始數(shù)據(jù)，構(gòu)成10行5列的矩陣，記作X。利用公式（2）對該矩陣進行標(biāo)準(zhǔn)化處理，得標(biāo)準(zhǔn)化之后的矩陣，記作SX，并求出其相關(guān)系數(shù)矩陣C.Matlab程序如下：

X=[0.61 1.28 1.89 4.02 14.33;

0.57 1.35 1.92 4.60 13.89;

1.03 1.56 2.59 5.85 17.53;

0.98 1.51 2.49 5.50 18.02;

0.50 1.96 2.46 5.26 15.97;

0.51 1.89 2.40 5.43 16.24;

0.64 1.77 2.41 4.85 20.21;

0.63 1.69 2.32 5.13 21.23;

0.46 1.36 1.82 9.07 16.05;

0.42 1.29 1.71 8.89 15.77];

m=size（X，1）;

n=size（X，2）;

SX=zeros（m，n）;

for i=1：b

SX（：，i）=（X（：，i）-mean（X（：，i）））/std（X（：，i））;

end

接著計算矩陣C的特征值與特征向量，Matlab程序為：

C=corrcoef（SX）;

[V，D]=eig（C）;

xlswrite（'D'，D）

xlswrite（'V'，V）

程序中的D為矩陣C的特征值，V為對應(yīng)的特征值向量。矩陣的特征值具體為λ₁=2.6938，λ₂=1.0168，λ₃=0.8805，λ₄=0.4089，λ₅=1.11e-16。

計算主成分的累計貢獻率，Matlab程序如下：

for j=1：n

DD（j，1）=D（n+1-j，n+1-j）;

end

for i=1：n

DD（i，2）=DD（i，1）/sum（DD（：，1））;

DD1（i）=DD（i，1）;

DD（i，3）=sum（DD1（1：i））/sum（DD（：，1））;

end

T=0.9;

for k=1：n

if DD（k，3） >=T

num=k;

break;

end

end

m1=size（V，1）;

PV=zeros（m1，num）;

for j=1：num

PV（：，j）=V（：，n+1-j）;

end

num

PV

DD

將上述程序運行后得到num=3，這表明前3個主成分的累計貢獻率大于90%。DD的第1、2、3列分別為特征值、貢獻率、累計貢獻率。由DD的第3列可以看到，前3個主成分的累計貢獻率已經(jīng)達到91.82%，因此保留了前3個主成分信息。PV運行的結(jié)果表示矩陣C特征值λ₁～λ₃對應(yīng)的三個單位正交特征向量，將其代入公式（1）得到前3個主成分的線性表達式為：

其中變量X₁～X₅分別為大黃素、大黃酚、大黃素和大黃酚總量、鹽酸小檗堿以及黃岑苷的含量。可以看出主成分1中大黃素、大黃酚、大黃素和大黃酚總量及黃岑苷的系數(shù)都大于0，與主成分1正相關(guān)，鹽酸小檗堿的系數(shù)小于0，與主成分1負(fù)相關(guān);主成分2只有與大黃素與主成分2正相關(guān)，其余都負(fù)相關(guān);主成分3除了與大黃酚負(fù)相關(guān)，其余都正相關(guān)。

運行下列運算程序可得到因子載荷矩陣Q：

Q=zeros（m1，num）;

for i=1：n

for j=1：num

Q（i，j）=PV（i，j）*sqrt（DD（j，1））;

end

end

Q

因子載荷矩陣Q的運行結(jié)果見表1，前3個主成分與原始成分的因子載荷陣如表1所示。主成分1中鹽酸小檗堿的因子載荷小于0，對主成分1產(chǎn)生負(fù)向影響，其余都大于0，大黃素、大黃酚、大黃素和大黃酚總量、黃芩苷這四個原始成分對主成分1產(chǎn)生正向影響。大黃素和大黃酚總量的因子載荷為0.9686，為高度相關(guān)，說明與主成分1的聯(lián)系密切程度非常高。大黃素、大黃酚、黃芩苷與主成分1顯著相關(guān)，說明與主成分1的聯(lián)系密切度高。鹽酸小檗堿的因子載荷為-0.5906，為顯著負(fù)相關(guān)。可見主成分1越大時大黃素、大黃酚和黃岑苷的含量也越高，而鹽酸小檗堿的含量越低。

主成分2中大黃素的因子載荷大于0，對主成分2產(chǎn)生正向影響，其余都小于0，大黃酚、大黃素和大黃酚總量、鹽酸小檗堿、黃芩苷這四個原始成分對主成分2產(chǎn)生負(fù)向影響。大黃素的因子載荷為0.7144，與主成分2顯著相關(guān)，說明與主成分1的聯(lián)系密切度高。大黃酚的因子載荷為-0.5971，與主成分2顯著負(fù)相關(guān)?？梢娭鞒煞?越大時大黃素的含量越高，大黃酚、鹽酸小檗堿和黃岑苷的含量越低，而且大黃素和大黃酚總量越小。

主成分3中大黃素、大黃素和大黃酚總量、鹽酸小檗堿和黃岑苷的因子載荷都大于0，這三個原始成分對主成分3產(chǎn)生正向影響，大黃酚的因子載荷小于于0，對主成分3產(chǎn)生負(fù)向影響。鹽酸小檗堿、黃芩苷的因子載荷分別是0.6726、0.5458，與主成分3顯著相關(guān)，說明與主成分3的聯(lián)系密切度高?？梢娭鞒煞?越大時大黃素、鹽酸小檗堿和黃岑苷的含量越高，大黃酚的含量越低，但大黃素和大黃酚總量仍然越大。

根據(jù)公式（1-3）計算綜合評價函數(shù)F， 并進行排序， Matlab程序如下：

score=SX*PV;

for i=1：m

F（i，1）=score（i，：）*DD（1：3，2）;

F（i，2）=i;

end

Fscore=[score， F];

disp（'各主成分得分、F值及排序（按第4列的F值進行降序排列，前3列為個各成分得分，第5列為三黃片的品種）。

F_sort=sortrows（Fscore，-4）

運行程序后可得到不同生產(chǎn)廠家的各成分的F值，并進行排序，結(jié)果見表2。三黃片10批樣品從高到低的排序的編號為：S3、S4、S7、S8、S5、S6、S1、S2、S9、S10?？梢钥闯鯯3、S4的F值比較接近，S7和S8的F值接近，S5與S6的F值接近，S1、S2的F值接近，S9、S10的F值接近，這與文獻聚類分析的結(jié)果一致。S3-S8的F值都大于0，說明樣品的品質(zhì)較高。最后可以根據(jù)樣品的編號再確定對應(yīng)實際廠家所生產(chǎn)的三黃片的質(zhì)量。

3?結(jié)束語

數(shù)學(xué)實驗是聯(lián)系數(shù)學(xué)知識和醫(yī)學(xué)問題的橋梁，通過數(shù)學(xué)軟件把數(shù)學(xué)基礎(chǔ)理論知識、計算機技術(shù)和醫(yī)學(xué)實際問題有效的結(jié)合，能夠激發(fā)學(xué)生的學(xué)習(xí)興趣，培養(yǎng)學(xué)生的應(yīng)用能力，如何更好地開展醫(yī)學(xué)院校的數(shù)學(xué)實驗課程教學(xué)是我們將不斷探索的問題。

參考文獻

[1]董鴿，徐方勤.淺談與計算機專業(yè)相結(jié)合的高等數(shù)學(xué)教學(xué)改革[J].科教文匯，2018，（4）：67-69.

[2]Hui Zou，Trevor Hastie，Robert Tibshirani.Sparse Principal Component Analysis[J].Journal of Computational and Graphical Statistics，2006，15（2）：265-286.

[3]Hancock P，Burton A，Bruce V.Face processing：human perception and principal components analysis[J].Memory and Cognition，1996，（24）：26-40.

[4]Carly L.Clayman，Satish M.Srinivasan，Raghvinder S.Sangwan，K-means Clustering and Principal Components Analysis of Microarray Data of L1000 Landmark Genes[J].Procedia Computer Science，2020（168）：97-104.

[5]張婷，劉慧琴，郭勤衛(wèi)，等.十六份辣椒材料游離氨基酸組成的主成分分析與聚類分析[J].浙江農(nóng)業(yè)學(xué)報，2021，33（4）：640-650.

[6]崔新剛，王穎瑩，王新勝，等.基于主成分分析和聚類分析的不同廠家三黃片質(zhì)量分析[J].中國民族民間醫(yī)藥，2020，29（16）：36-40，43.

[7]史周華，張雪飛.中醫(yī)藥統(tǒng)計學(xué)[M].北京：科學(xué)出版社，2009.

基金項目：2020年上海健康醫(yī)學(xué)院師資百人庫計劃項目“統(tǒng)計方法在中草藥分析中的應(yīng)用”（B3020020311084）。

作者簡介：董鴿（1980-），女，漢族，江蘇泰州人，理學(xué)博士，上海健康醫(yī)學(xué)院副教授，主要從事數(shù)學(xué)教學(xué)研究。