典范對(duì)應(yīng)分析數(shù)據(jù)模擬物種豐度計(jì)算方法

劉志麗， 吳亞玲， 周鐵軍

(湖南農(nóng)業(yè)大學(xué)信息與智能科學(xué)技術(shù)學(xué)院，長(zhǎng)沙 410128)

1 引 言

(1)

典范對(duì)應(yīng)分析需要檢驗(yàn)每個(gè)物種j在不同樣地上的豐度標(biāo)準(zhǔn)值qij與這些樣地上環(huán)境因子是否存在顯著的線性關(guān)系.對(duì)檢驗(yàn)方法的有效性分析需要由某種分布隨機(jī)生成的環(huán)境因子的數(shù)據(jù)依據(jù)線性關(guān)系或無(wú)關(guān)系生成標(biāo)準(zhǔn)矩陣Q.現(xiàn)在的問(wèn)題是，如何由矩陣Q確定絕對(duì)頻率矩陣Y=[yij]n×p？正是因?yàn)檫@個(gè)問(wèn)題的存在，使得參考文獻(xiàn)[1]附錄D在進(jìn)行檢驗(yàn)方法的有效性分析時(shí)沒(méi)有提供檢驗(yàn)發(fā)生第二類錯(cuò)誤的仿真計(jì)算結(jié)果.因此解決由標(biāo)準(zhǔn)矩陣Q計(jì)算物種豐度矩陣Y問(wèn)題具有理論意義.當(dāng)然從式(1)可以知道，對(duì)于兩組成比例的豐度數(shù)據(jù)，其標(biāo)準(zhǔn)矩陣Q是相同的，即由Q并不能唯一確定物種豐度矩陣Y，但我們希望能夠確定其中一組豐度值.

2 等效的齊次線性方程組存在半正解問(wèn)題

由式(1)得

(2)

對(duì)j求和得

于是有

(3)

類似地，由(2)對(duì)j求和還可得

(4)

所以

(5)

從等式(3)-(5)，問(wèn)題化為要求滿足下列條件的非負(fù)不全為0的si(i=1，2，…，n)和tj(j=1，2，…，p):

(6)

(7)

(8)

用矩陣表示為

QTs=0， Qt=0， sTs=tTt.

顯然s=0， t=0能夠滿足上述方程(6)-(8)，但沒(méi)有實(shí)際意義.因此希望求方程組(6)-(8)的非負(fù)不全為0的解.確定了si與tj后，就可以根據(jù)下式計(jì)算yij，

(9)

注意到條件(8)可以通過(guò)將齊次方程(6)與(7)的非零解單位化后得到.所以問(wèn)題的關(guān)鍵是如何求齊次方程(6)與(7)的非負(fù)非零解使yij≥0.齊次線性方程組的非負(fù)非零解一般稱為半正解，如果解的全部分量均為正，則稱為正解.

參考文獻(xiàn)[2]獲得了齊次線性方程組Ax=0有正解的充分必要條件，即這個(gè)方程組存在若干個(gè)減列方程組， 它們都是有正解的極小方程組，并且這些極小方程組全部列向量的并集的秩等于A的秩.但確定極小方程組及其正解不是一件簡(jiǎn)單的事.參考文獻(xiàn)[3]對(duì)一類正線性方程組討論了它的半正解的存在性，而參考文獻(xiàn)[4]也從系數(shù)矩陣A出發(fā)給出了齊次線性方程組Ax=0有非負(fù)解的充分條件及必要條件.

通過(guò)以下方法求齊次方程(6)與(7)的半正解.顯然容易求得齊次線性方程組(6)或(7)的基礎(chǔ)解系，所以問(wèn)題轉(zhuǎn)化為如何用齊次線性方程組的一組基礎(chǔ)解系表示出它的一個(gè)半正解.參考文獻(xiàn)[5]對(duì)基礎(chǔ)解系中自由未知量取標(biāo)準(zhǔn)單位向量且只有兩個(gè)向量的條件下討論了半正解的存在性.本文研究用齊次線性方程組的任一組基礎(chǔ)解系(含兩個(gè)或三個(gè)向量)表示出它的一個(gè)半正解的條件，并給出數(shù)值計(jì)算方法.

3 由基礎(chǔ)解系表示半正解的條件

顯然，如果基礎(chǔ)解系只含一個(gè)非零向量，除非它的分量全部為正或全部為負(fù)，否則它不能表示出一個(gè)半正解.

3.1 基礎(chǔ)解系至少含分量不全為正的兩個(gè)向量

先考慮基礎(chǔ)解系至少含分量不全為正的兩個(gè)向量情形，即它們的分量中即有正分量也有負(fù)分量.這些正負(fù)分量在兩個(gè)向量中的位置有三種情形，一是負(fù)分量個(gè)數(shù)少的向量其負(fù)分量所在位置和負(fù)分量個(gè)數(shù)多的相應(yīng)位置上的分量都是負(fù)分量；二是負(fù)分量個(gè)數(shù)少的向量其負(fù)分量所在位置和負(fù)分量個(gè)數(shù)多的相應(yīng)位置上的分量都是正分量；三是負(fù)分量個(gè)數(shù)少的向量其負(fù)分量所在位置和負(fù)分量個(gè)數(shù)多的相應(yīng)位置上的分量有一部分是負(fù)分量，有一部分是正分量.設(shè)基礎(chǔ)解系中向量α1中有k個(gè)負(fù)分量，向量α2中有l(wèi)個(gè)負(fù)分量，不妨假設(shè)k≥l.

情形1 設(shè)基礎(chǔ)解系中含有向量

α1=(-a11，…，-a1l，-a1，l+1，…，-a1k，a1，k+1，…，a1n)，

α2=(-a21，…，-a2l，a2，l+1，…，a2k，a2，k+1，…，a2n)，

其中aij≥0.記A={1，2，…，l}，B={l+1，l+2，…，k}，C={k+1，k+2，…，n}.如果k=l，則B=?.

證向量α1與α2的線性組合可以表示為

c1α1+c2α2=(…，-c1a1i-c2a2i，…，-c1a1i1+c2a2i1，…，c1a1j+c2a2j，…).

如果a1j=0，必有c1<0，c2>0，使得

-c1a1i-c2a2i>0， -c1a1i1+c2a2i1>0，c1a1j+c2a2j>0.

于是得到

-c1a1i-c2a2i>0， -c1a1i1+c2a2i1>0，c1a1j+c2a2j>0.

故得到

c1α1+c2α2>0.

情形2 設(shè)基礎(chǔ)解系中含有向量

α1=(-a11，…，-a1k，a1，k+1，…，a1，k+l，a1，k+l+1，…，a1n)，

α2=(a21，…，a2k，-a2，k+1，…，-a2，k+l，a2，k+l+1，…，a2n)，

其中aij≥0，n≥k+l.記A={1，2，…，k}，B={k+1，k+2，…，k+l}，C={k+l+1，k+l+2，…，n}.如果n=k+l，則C=?.由于

-α2=(-a21，…，-a2k，a2，k+1，…，a2，k+l，-a2，k+l+1，…，-a2n)，

于是-α2中負(fù)分量個(gè)數(shù)n-l大于等于α1中負(fù)分量個(gè)數(shù)k，即n-l≥k，且α1中負(fù)分量所在位置為A，-α2中負(fù)分量位置為A∪C?A，故就是情形1.因此得到定理2.

情形3 設(shè)α2中l(wèi)個(gè)負(fù)分量只有s個(gè)位置與α1中負(fù)分量位置相同，其余負(fù)分量位于α1的正分量位置.不妨設(shè)基礎(chǔ)解系中向量

α1=(-a11，…，-a1s，-a1，s+1，…，-a1k，a1，k+1，…，a1，k+l-s，a1，k+l-s+1，…，a1n)，

α2=(-a21，…，-a2s，a2，s+1，…，a2k，-a2，k+1，…，-a2，k+l-s，a2，k+l-s+1，…，a2n)，

其中aij≥0，s

α1=(-a11，…，-a1s，-a1，s+1，…，-a1k，a1，k+1，…，a1n)，

-α2=(a21，…，a2s，-a2，s+1，…，-a2k，a2，k+1，…，a2n).

下面考慮n>k+l-s，即D≠?.假設(shè)存在常數(shù)c1與c2，使得

c1α1+c2α2=(…，-c1a1i-c2a2i，…，-c1a1i1+c2a2i1，…，c1a1j-c2a2j，…，c1a1j1+c2a2j1，…)>0.

其中i∈A，i1∈B，j∈C，j1∈D.因此有

c1a1i<-c2a2i，c1a1i1c2a2j，c1a1j1>-c2a2j1.

如果c2>0，則由c1a1i<-c2a2i知c1<0，而由c1a1j>c2a2j知c1>0.如果c2<0，則由c1a1i1-c2a2j1知c1>0.因此上述不等式不能同時(shí)成立，即在此種情形不可能將α1與α2線性組合成一個(gè)非負(fù)非零向量.此時(shí)需要基礎(chǔ)解系中至少含有三個(gè)向量.

3.2 基礎(chǔ)解系至少含分量不全為正的三個(gè)向量

只需考慮三個(gè)向量中任意兩個(gè)的負(fù)分量位置不全部相同的情形.不妨設(shè)基礎(chǔ)解系中有向量

α1=(-a11，…，-a1s，-a1，s+1，…，-a1k，a1，k+1，…，a1，k+l-s，a1，k+l-s+1，…，a1n)，

α2=(-a21，…，-a2s，a2，s+1，…，a2k，-a2，k+1，…，-a2，k+l-s，a2，k+l-s+1，…，a2n)，

其中aij≥0，sk+l-s.記A={1，2，…，s}，B={s+1，s+2，…，k}，C={k+1，k+2，…，k+l-s}，D={k+l-s+1，k+l-s+2，…，n}，且A=A-∪A+，B=B-∪B+，C=C-∪C+，D=D-∪D+.向量α3=(b3j)中分量的正負(fù)性如下：

其中a3j≥0.

如果A-，A+，B-，B+，C-，C+，D-，D+都不是空集，則由c1α1+c2α2+c3α3>0可得

-c1a1i1-c2a2i1-c3a3i1>0，i1∈A-，

(10)

-c1a1i2-c2a2i2+c3a3i2>0，i2∈A+，

(11)

-c1a1i3+c2a2i3-c3a3i3>0，i3∈B-，

(12)

-c1a1i4+c2a2i4+c3a3i4>0，i4∈B+，

(13)

c1a1j1-c2a2j1-c3a3j1>0，j1∈C-，

(14)

c1a1j2-c2a2j2+c3a3j2>0，j2∈C+，

(15)

c1a1j3+c2a2j3-c3a3j3>0，j3∈D-，

(16)

c1a1j4+c2a2j4+c3a3j4>0，j4∈D+.

(17)

如果c2>0，c3>0，則由(10)推出c1<0，而由(14)推出c1>0，相互矛盾.

如果c2<0，c3>0，則由(12)推出c1<0，而由(16)推出c1>0，相互矛盾.

如果c2>0，c3<0，則由(11)推出c1<0，而由(15)推出c1>0，相互矛盾.

如果c2<0，c3<0，則由(13)推出c1<0，而由(17)推出c1>0，相互矛盾.

因此，如果A-，A+，B-，B+，C-，C+，D-，D+都不是空集，即α3的正、負(fù)分量分別位于A，B，C，D中，則不能由α1，α2，α3線性組合生成一個(gè)非負(fù)非零的解向量.所以A-，A+，B-，B+，C-，C+，D-和D+至少有一個(gè)空集，可以只考慮A-，B-，C-，D-四個(gè)集合中至少有一個(gè)是空集的情形，因?yàn)槿绻鸄+，B+，C+，D+中至少有一個(gè)空集，則-α3的負(fù)分量位置A-，B-，C-，D-四個(gè)集合中至少有一個(gè)是空集.A-，B-，C-，D-四個(gè)集合中至少有一個(gè)是空集共有10種情形.

a1i1x+a2i1y>a3i1，i1∈A-，

(18)

a1i2x+a2i2y>-a3i2，i2∈A+，

(19)

a1i3x-a2i3y>a3i3，i3∈B-，

(20)

a1i4x-a2i4y>-a3i4，i4∈B+，

(21)

a1j1x-a2j1y<-a3j1，j1∈C-，

(22)

a1j2x-a2j2y

(23)

a1j4x+a2j4y

(24)

其中(19)是顯然成立的.如果其它6個(gè)不等式有解，則可以由α1，α2，α3線性組合生成一個(gè)非負(fù)非零的解向量.記

于是

所以不等式(18)，(23)成立.又

所以(20)，(21)，(22)和(24)成立.因此正實(shí)數(shù)x，y就是方程組(18)，(20)-(24)的一組解.從而由基礎(chǔ)解系必可以線性表示出一個(gè)半正解.

成立，其中i1∈A-，i3∈B-，i4∈B+，j1∈A-，j2∈B+，j4∈B+，則由基礎(chǔ)解系必可以線性表示出一個(gè)半正解.

證由不等式(c)可以知道a1j2a2j1-a1j1a2j2>0.根據(jù)條件(d)，必存在正實(shí)數(shù)x使得

所以得

且

從而有

與

所以(22)與(23)成立.又

所以

從而

所以(20)，(21)成立.因此正實(shí)數(shù)x，y就是方程組(18)，(20)-(24)的一組解.從而由基礎(chǔ)解系必可以線性表示出一個(gè)半正解.

4 基礎(chǔ)解系表示半正解的數(shù)值計(jì)算方法

可以將用齊次線性方程組的一組基礎(chǔ)解系表示出它的一個(gè)半正解歸結(jié)為一個(gè)二次規(guī)劃問(wèn)題.具體計(jì)算步驟如下：

第一步，求齊次線性方程組(6)的基礎(chǔ)解系ξ1，ξ2，…，ξn-k；

第二步，由基礎(chǔ)解系ξ1，ξ2，…，ξn-k確定齊次方程(6)的一個(gè)半正解

s=c1ξ1+c2ξ2+… +cn-kξn-k，

其中系數(shù)ci是如下二次規(guī)劃問(wèn)題的最優(yōu)解：

其中b是一個(gè)任意設(shè)定的正向量，lc是一個(gè)負(fù)向量，規(guī)定c=(c1，c2，…，cn-k)T的取值下限，uc是正向量，規(guī)定c的取值上限.fi是任意正數(shù).這是一個(gè)二次規(guī)劃問(wèn)題，可以通過(guò)Matlab 函數(shù)quadprog計(jì)算得到c.

s=null(Q′);

r=size(s，2);

H=eye(r);

f=rand(r，1);

lc=-5*ones(r，1);

uc=-lc;

b=0.01*rand(n，1);

c=quadprog(H，f，-s，-b，[]，[]，lc，uc);

s=s*c;

第三步，將半正解s單位化得si.

類似地可以通過(guò)上述步驟計(jì)算齊次線性方程組(7)的半正解t，單位化后得ti.

例2取矩陣

利用上述算法可以計(jì)算得方程(6)的一個(gè)半正解

與(7)的半正解

它們都是單位向量.

進(jìn)一步根據(jù)計(jì)算得到的兩個(gè)半正解s與t，由(9)式計(jì)算得物種豐度矩陣如下：

(25)

利用(1)式，可以由上述矩陣Y計(jì)算得到矩陣Q，這驗(yàn)證了提出的由Q逆算Y的方法的準(zhǔn)確性.例2中的Q是由參考文獻(xiàn)[1]中的第1～4及第6種蜘蛛在前10個(gè)樣地上數(shù)量經(jīng)(1)轉(zhuǎn)化成的標(biāo)準(zhǔn)矩陣.如果我們將矩陣Y乘以10000/38，則可以得到

不難發(fā)現(xiàn)，除了3個(gè)帶小數(shù)的數(shù)外，其它數(shù)據(jù)和原始觀察數(shù)據(jù)一致，而且這三個(gè)帶小數(shù)部分的數(shù)據(jù)其整數(shù)部分也是與相應(yīng)觀察值一致的.

5 結(jié) 論

研究典范對(duì)應(yīng)分析置換檢驗(yàn)的有效性時(shí)需要隨機(jī)生成大量數(shù)據(jù)，利用齊次線性方程組半正解存在性有效地解決了由生成數(shù)據(jù)轉(zhuǎn)換為物種豐度的計(jì)算問(wèn)題.在齊次線性方程組基礎(chǔ)解系至少具有兩個(gè)分量不全為負(fù)的解向量情形下，給出了由這兩個(gè)解向量線性表示半正解的充分條件，在不能由含兩個(gè)分量不全為負(fù)的解向量線性表示半正解時(shí)，又給出了一組含三個(gè)分量不全為負(fù)的解向量的基礎(chǔ)解系線性表示半正解的充分條件.

致謝作者非常感謝相關(guān)文獻(xiàn)對(duì)本文的啟發(fā)以及審稿專家提出的寶貴意見.