一種病變眼底圖像的配準(zhǔn)方法

魏麗芳 潘 林 余 輪

（福州大學(xué)計(jì)算機(jī)圖象圖形研究所，福州 350002）

引言

眼底圖像配準(zhǔn)技術(shù)在輔助眼科診斷和治療過(guò)程中有著廣泛的應(yīng)用［1-2］。通過(guò)配準(zhǔn)眼底圖像，醫(yī)生能夠更好地診斷和檢測(cè)各種與眼底相關(guān)的疾病，如糖尿病、青光眼、黃斑部變性等，圖1給出兩幅病變眼底圖像示例。

近年來(lái)，各種眼底圖像配準(zhǔn)技術(shù)已得到廣泛研究，目前眼底圖像配準(zhǔn)方法主要分為基于灰度的方法和基于特征的方法［3］?；趫D像灰度的配準(zhǔn)方法直接利用圖像的灰度值來(lái)確定配準(zhǔn)的空間變換，其充分利用圖像中所包含的信息；但基于灰度的方法容易陷入一個(gè)局部極小值，而且變換模型的搜索空間是巨大的。基于特征的方法是以待配準(zhǔn)圖像對(duì)之間存在的對(duì)應(yīng)特征點(diǎn)（如血管分支和交叉點(diǎn)）為基礎(chǔ)，通過(guò)測(cè)量對(duì)應(yīng)點(diǎn)之間的相似性最值完成。鑒于兩種方法的特點(diǎn)，基于特征的眼底圖像配準(zhǔn)算法得到越來(lái)越多的應(yīng)用。Can等提出層次匹配方法［2］，基于血管分支和交叉點(diǎn)，對(duì)于血管分支和交叉點(diǎn)的提取利用了血管追蹤和特征抽取的模式。Zana等利用形態(tài)學(xué)方法提取血管樹，將血管樹分支和交叉點(diǎn)（特征點(diǎn)）用鄰近血管的方向進(jìn)行標(biāo)注，實(shí)現(xiàn)每對(duì)匹配的特征點(diǎn)之間的角度不變性［4］。但這類方法仍然依賴于眼底血管結(jié)構(gòu)的提取。文獻(xiàn)［3］將適合曲面配準(zhǔn)的迭代最近點(diǎn)iterative closest point，ICP）算法引入眼底圖像配準(zhǔn)，提出雙重引導(dǎo)的迭代最近點(diǎn)算法（dual-bootstrap iterative closest point，DB-ICP）。該算法對(duì)于眼底的曲面結(jié)構(gòu)可以獲得精確的配準(zhǔn)，在較好初值的前提下，算法具有較好的收斂性。但DB-ICP秉承了ICP算法耗時(shí)大的缺點(diǎn)，而且對(duì)于低質(zhì)量的眼底圖像，也很難配準(zhǔn)成功。文獻(xiàn)［5］則將尺度不變性特征變換（scale invariable feature transformation，SIFT）引入眼底圖像配準(zhǔn)，但使用的隨機(jī)抽樣一致性算法（random sample consensus，RANSAC）算法進(jìn)行初始估計(jì)，它在去除誤匹配同時(shí)也去除了大量正確匹配特征。圖1為糖尿病誘發(fā)的眼底病變示例。

圖1 不同病變的眼底圖像。（a）增殖性糖尿病視網(wǎng)膜五期玻璃體病變；（b）非增殖性糖尿病視網(wǎng)膜三期出血灶病變Fig.1 Different pathological retinal images.（a）proliferative diabetic retinopathy Ⅴ vitreous gel pathologic；（b）nonproliferative diabetic retinopathy Ⅲ hemorrhagic foci pathologic

針對(duì)病變眼底圖像血管結(jié)構(gòu)模糊同時(shí)也具有復(fù)雜的球體解剖學(xué)結(jié)構(gòu)和光學(xué)圖像系統(tǒng)的特點(diǎn)，提出了一種基于不變特征的眼底圖像配準(zhǔn)方法。首先提取SIFT特征作為匹配的特征點(diǎn)，提出了雙邊“或”的best-bin-first（BBF）特征匹配算法對(duì)特征點(diǎn)進(jìn)行初始匹配，通過(guò)SIFT特征點(diǎn)的方向特征和空間幾何特性去除誤匹配，精化匹配特征對(duì)，并保留足夠的正確匹配對(duì)；根據(jù)獲得的匹配特征通過(guò)迭代加權(quán)最小二乘法和M-估計(jì)進(jìn)行參數(shù)估計(jì)，實(shí)現(xiàn)病變眼底圖像的精確配準(zhǔn)。

1 特征提取與匹配

1.1 特征點(diǎn)的提取

圖像不變特征（如，角點(diǎn)、SIFT、SURF）和描述符的不斷發(fā)展為圖像精確配準(zhǔn)提供了有利的前提［6-9］。在之前的研究中已經(jīng)證明SIFT算法在眼底圖像配準(zhǔn)中的有效性［5］。SIFT特征檢測(cè)算法詳見文獻(xiàn)［8-9］，是一種基于尺度空間的、對(duì)圖像縮放、旋轉(zhuǎn)甚至仿射變換保持不變性的圖像局部特征描述算子，該算法首先在尺度空間進(jìn)行特征檢測(cè)，并確定關(guān)鍵點(diǎn)的位置和關(guān)鍵點(diǎn)所處的尺度，然后使用關(guān)鍵點(diǎn)鄰域梯度的主方向作為該點(diǎn)的方向特征，以實(shí)現(xiàn)算子對(duì)尺度和方向的無(wú)關(guān)性。其對(duì)旋轉(zhuǎn)、尺度縮放、亮度變化保持不變性，對(duì)視角變化、仿射變換、噪聲也保持一定程度的穩(wěn)定性，并對(duì)每一個(gè)關(guān)鍵點(diǎn)產(chǎn)生一個(gè)128維的特征描述符以作為匹配依據(jù)。

1.2 特征匹配

為了提高特征正確匹配率，同時(shí)確保具有足夠的匹配對(duì)應(yīng)用于圖像的兩兩配準(zhǔn)，提出一種雙邊“或”的Best-Bin-First（BBF）算法建立兩個(gè)圖像特征點(diǎn)之間的相應(yīng)性，利用特征點(diǎn)的方向信息和空間幾何信息去除誤匹配，盡可能地保留正確的匹配。

根據(jù)特征點(diǎn)的鄰域信息建立高維特征描述符，通過(guò)最近鄰比值方法，由逼近高維空間的最近鄰搜索的BBF算法找到待匹配圖像中對(duì)應(yīng)的特征匹配點(diǎn)對(duì)［9-10］；BBF算法通過(guò)尋找最近鄰和次近鄰，同時(shí)對(duì)最近鄰與次近鄰的比值設(shè)定閾值T（文中T取0.5）。

設(shè)圖像I1和I2的SIFT描述符分別是DES1和DES2，計(jì)算相應(yīng)的距離［10］：

式中，·代表向量的點(diǎn)積。這個(gè)集合包含了des與I2中所有描述符之間的距離。dess對(duì)應(yīng)于Dis的最大元素，表示des的最近鄰。如果最近鄰距離與次近鄰距離的比值小于T，則說(shuō)明構(gòu)成最近鄰的兩點(diǎn)是匹配的，否則就不匹配。該算法可有效達(dá)到穩(wěn)定的匹配：在匹配對(duì)中，正確的匹配對(duì)應(yīng)的最短距離應(yīng)該比錯(cuò)誤匹配對(duì)應(yīng)的的最短距離短得多，而在高維特征空間中相似的距離內(nèi)可能還存在著一些其他的錯(cuò)誤匹配，利用次近鄰匹配作為這部分空間錯(cuò)誤匹配密度的估計(jì)。

以上匹配過(guò)程是單邊匹配的，然而在單邊匹配時(shí)，存在一些遺漏的正確匹配。對(duì)于具有曲面結(jié)構(gòu)的眼底圖像，正確的匹配對(duì)越多且分布越廣泛（重疊區(qū)域內(nèi)），獲得的圖像對(duì)之間的變換關(guān)系越精確。基于保留更多正確匹配對(duì)的考慮，建立雙邊的匹配，即找出I1到I2之間匹配和I2到I1匹配，建立I1和I2之間的完全匹配，進(jìn)一步解決了漏匹配問題。雙邊“或”的BBF算法在獲得更多正確匹配的同時(shí)，也獲得了更多的誤匹配。

為進(jìn)一步去除初步建立的匹配特征對(duì)中的誤匹配，首先利用SIFT特征旋轉(zhuǎn)不變性，根據(jù)匹配特征的方向一致性初步去除誤匹配。然后通過(guò)建立對(duì)應(yīng)特征對(duì)之間的空間幾何關(guān)系進(jìn)行一致性檢測(cè)，有效去除誤匹配。匹配點(diǎn)一致性檢測(cè)運(yùn)用類似于文獻(xiàn)［11］中的原理，代替計(jì)算對(duì)應(yīng)線段比作為它的檢測(cè)數(shù)據(jù)。通過(guò)對(duì)應(yīng)點(diǎn)空間幾何關(guān)系-空間斜率和距離進(jìn)行一致性檢測(cè)，去除誤匹配。它基于這樣的事實(shí)：對(duì)于所有對(duì)應(yīng)特征點(diǎn)對(duì)，它們的空間斜率和空間距離應(yīng)該具有一致性。

當(dāng)所有點(diǎn)對(duì)都是物理對(duì)應(yīng)點(diǎn)時(shí)，對(duì)應(yīng)的空間斜率和空間距離應(yīng)是一樣的。一致性檢測(cè)首先對(duì)所有線段比值數(shù)據(jù)進(jìn)行分級(jí)聚類，找出一致性數(shù)值。然后通過(guò)空間斜率和空間距離數(shù)據(jù)處理，迭代刪除偽匹配對(duì)?？臻g幾何關(guān)系的一致性檢測(cè)算法詳細(xì)描述如下：

步驟1：計(jì)算所有對(duì)應(yīng)特征點(diǎn)對(duì)的空間斜率p_theta和空間距離p_dis

式中，（xi，yi）和（xj，yj）分別為I1和I2中相對(duì)應(yīng)的匹配特征點(diǎn)對(duì)。

步驟2：對(duì)得到的兩組數(shù)據(jù)分別進(jìn)行分級(jí)聚類，使聚類結(jié)果的類別數(shù)少且類內(nèi)距離足夠小。聚類結(jié)果若有某一類樣本個(gè)數(shù)大于3，且大于其它類的樣本個(gè)數(shù)，則該類樣本值稱為一致性數(shù)值。

步驟3：統(tǒng)計(jì)其余樣本的非一致性數(shù)據(jù)個(gè)數(shù)和一致性數(shù)據(jù)的方差。

步驟4：將非一致性數(shù)據(jù)個(gè)數(shù)值最大的樣本集作為候選樣本，然后將方差最大的樣本集刪除。

步驟5：再次統(tǒng)計(jì)剩余樣本的非一致性數(shù)據(jù)個(gè)數(shù)和一致性數(shù)據(jù)的方差。

步驟6：若非一致性數(shù)據(jù)的樣本集數(shù)仍大于3，且至少有一個(gè)樣本集統(tǒng)計(jì)的值大于0，則轉(zhuǎn)到步驟4；否則，輸出保留的一致性數(shù)據(jù)樣本集，兩種幾何關(guān)系一致性檢測(cè)后保留的樣本集中的數(shù)據(jù)對(duì)應(yīng)的點(diǎn)對(duì)即為對(duì)應(yīng)的匹配對(duì)。

圖2為特征匹配結(jié)果。

圖2 特征匹配。（a）初始匹配；（b）最終匹配Fig.2 Feature matching.（a）initial match；（b）final match

2 變換模型及參數(shù)估計(jì)

眼底結(jié)構(gòu)通常近似于剛體，而眼底表面是一個(gè)空間的曲面，擁有所有曲面的通性，即能以一定數(shù)學(xué)表達(dá)式來(lái)表達(dá)，并反映其空間形態(tài)的本質(zhì)特征［12］。文獻(xiàn)［13-14］中認(rèn)為眼底是一個(gè)二次曲面，而二次曲面包括了橢球面、球面、拋物面及雙曲面。眼底的二次曲面模型，即多項(xiàng)式模型是在相機(jī)的弱透視模型和對(duì)眼底曲面結(jié)構(gòu)的近似，因此以二次曲面的數(shù)學(xué)模型近似眼底表面是合理的。點(diǎn)p=（x，y）T和q=（x’，y’）T分別是圖像I1和I2中的相應(yīng)匹配點(diǎn)，兩幅圖像的空間二次曲面模型的變化關(guān)系可表示為［2］

其中θm，n是2×6參數(shù)矩陣

M估計(jì)是一種魯棒的參數(shù)估計(jì)方法，也可以消除誤匹配帶來(lái)的外點(diǎn)，從而進(jìn)一步提高估算精度。對(duì)θ的M估計(jì)算法具體如式（7）和式（8）所示。

式中，k表示待配準(zhǔn)的兩圖像具有k對(duì)匹配的特征點(diǎn)對(duì)，ρ是具有魯棒性損失函數(shù)。進(jìn)行M估計(jì)的一個(gè)重要技術(shù)就是根據(jù)實(shí)際情況合理選擇ρ。常用的三種ρ函數(shù)的性能如圖3所示［16］，圖中橫坐標(biāo)為歸一化的誤差u，縱坐標(biāo)為ρ的權(quán)重值。根據(jù)魯棒估計(jì)理論［15］，σ可通過(guò)式（9）計(jì)算：

圖3 不同損失函數(shù)權(quán)重值示意Fig.3 Plots of the weight functions for the robust loss functions

由圖3可見，Beaton-Tukey函數(shù)反應(yīng)最靈敏，收斂最快，可以有效去除最多的外點(diǎn)，故采用這個(gè)函數(shù)。

對(duì)θ的估計(jì)通過(guò)迭代加權(quán)最小二乘法（iteratively-reweighted least-squares，IRLS）實(shí)現(xiàn)。通常情況下，IRLS可以很快實(shí)現(xiàn)收斂。

在M估計(jì)中，權(quán)值函數(shù)求解方法為：

由式（10）可知，Beaton-Tukey函數(shù)對(duì)應(yīng)的權(quán)值函數(shù)為

IRLS迭代的初始矩陣取自仿射變換矩陣。

3 算法評(píng)估

源圖像由福建省附屬第一醫(yī)院眼科提供，圖像大小為3 072像素×2 048像素。目標(biāo)圖像共28組，需要說(shuō)明的是，所討論方法中一致性檢測(cè)的實(shí)現(xiàn)是基于同組圖像尺度變化較小的前提條件。為定量評(píng)價(jià)本方法的性能，運(yùn)用正確匹配特征的數(shù)量和圖像配準(zhǔn)精度進(jìn)行評(píng)估。

匹配特征對(duì)的數(shù)量和分布，是決定參數(shù)估計(jì)準(zhǔn)確性的重要因素。把所采用的雙邊“或”匹配方法與單邊方法及雙邊“與”方法進(jìn)行評(píng)測(cè)。匹配特征的數(shù)量評(píng)測(cè)準(zhǔn)則定義如下：令θ表示圖像I1和I2之間的映射關(guān)系，點(diǎn)p=（x，y）T和q=（x’，y’）T分別是圖像I1和I2中的相應(yīng)匹配點(diǎn)，理想情況下q=θX（p）。實(shí)際上，當(dāng)p和q在空間位置上足夠接近時(shí)認(rèn)為匹配關(guān)系存在，即‖q-θX（p）‖≤t，則認(rèn)為p和q是一對(duì)正確匹配。實(shí)驗(yàn)中選擇距離門限t=1.5。

圖像配準(zhǔn)精度是衡量圖像配準(zhǔn)算法性能的有效方法。這里采用均方根誤差（root of mean square error，RMSE）進(jìn)行配準(zhǔn)精度的測(cè)量：

式中，N表示最終得到的正確匹配對(duì)的數(shù)量。

在檢驗(yàn)所提出方法在特征點(diǎn)匹配及配準(zhǔn)精度等方面的優(yōu)勢(shì)時(shí)，對(duì)病變眼底圖像進(jìn)行實(shí)驗(yàn)。首先通過(guò)病變眼底圖像的配準(zhǔn)結(jié)果，說(shuō)明本方法的有效性；然后利用匹配特征的數(shù)量評(píng)估指標(biāo)，對(duì)采用雙邊“或”匹配方法與單邊方法的文獻(xiàn)［5］和文獻(xiàn)［9］及雙邊“與”方法的文獻(xiàn)［10］進(jìn)行評(píng)測(cè)；同時(shí)計(jì)算配準(zhǔn)成功圖像對(duì)的均方根誤差評(píng)估本方法，并與文獻(xiàn)［5］方法所獲得的配準(zhǔn)精度進(jìn)行對(duì)比，定量描述本方法的優(yōu)勢(shì)。

4 結(jié)果

對(duì)28組實(shí)驗(yàn)數(shù)據(jù)進(jìn)行實(shí)驗(yàn)，其中的27組實(shí)驗(yàn)數(shù)據(jù)配準(zhǔn)成功。圖4給出了兩組不同病變程度眼底圖像的配準(zhǔn)結(jié)果示例，圖4（a）的圖像受玻璃體病變影響使圖像模糊、渾濁，尤其病灶區(qū)域血管結(jié)構(gòu)已無(wú)法提取，圖4（b）的圖像受出血灶影響也很難提取完整的血管結(jié)構(gòu)。從兩組圖像的配準(zhǔn)結(jié)果中，可以看出配準(zhǔn)后的圖像保留了有效信息，病變區(qū)域及模糊的血管結(jié)構(gòu)均實(shí)現(xiàn)了良好的細(xì)節(jié)對(duì)齊，無(wú)明顯的錯(cuò)誤對(duì)齊現(xiàn)象，能夠滿足眼底的醫(yī)學(xué)診斷要求。

在分析正確匹配特征數(shù)量的實(shí)驗(yàn)中，分析本方法對(duì)保留正確匹配對(duì)的影響。表1給出了圖4中兩組示例圖像組產(chǎn)生的匹配特征情況，同時(shí)比較文獻(xiàn)［5］、［9］采用單邊方法和文獻(xiàn)［10］采用雙邊“與”方法取得的正確匹配特征的數(shù)量。從表1可以看出，本文方法保留了足夠的正確匹配對(duì)，與雙邊“與”方法相比能夠獲得更多的正確匹配對(duì)，與單邊方法相比也有一定的優(yōu)勢(shì)。

圖4 病變眼底圖像配準(zhǔn)結(jié)果。（a）示例1；（b）示例2Fig.4 Pathological retinal image registration results.（a）example 1；（b）example 2

表1 特征匹配方法比較Tab.1 Comparison of feature matching method

利用均方根誤差RMSE分析本研究配準(zhǔn)成功的27組圖像的配準(zhǔn)精度。同時(shí)與文獻(xiàn)［5］采用一般的特征點(diǎn)單邊匹配方法結(jié)合常用的RANSAC算法（RANdom SAmple Consensus，即隨機(jī)抽樣一致性算法），去除誤匹配點(diǎn)的基于不變特征配準(zhǔn)方法行進(jìn)對(duì)比實(shí)驗(yàn)，其中，經(jīng)實(shí)驗(yàn)統(tǒng)計(jì)在28組圖像序列中，本方法配準(zhǔn)成功27組病變眼底圖像，文獻(xiàn)［5］配準(zhǔn)成功25組。從圖5中可看出本方法對(duì)每組病變眼底圖像的RMSE值均小于1，并浮動(dòng)于0.5左右，遠(yuǎn)達(dá)到了亞像素級(jí)的配準(zhǔn)精度。從兩種方法的RMSE值的對(duì)比可以看出，本方法中每組圖像數(shù)據(jù)的RMSE值均小于文獻(xiàn)［5］方法的RMSE值，說(shuō)明本方法優(yōu)于文獻(xiàn)［5］方法的配準(zhǔn)精度。

圖5 配準(zhǔn)精度對(duì)比Fig.5 Comparison of registration accuracy

5 討論和結(jié)論

從實(shí)驗(yàn)結(jié)果可見，采用雙邊“或”的特征匹配和利用空間幾何特性的一致性檢測(cè)去除誤匹配的配準(zhǔn)方法的優(yōu)越性。同時(shí)說(shuō)明在基于特征的配準(zhǔn)方法中，準(zhǔn)確與豐富的匹配特征在改善和提高配準(zhǔn)性能、實(shí)現(xiàn)精確配準(zhǔn)方面的重要性。因此，采用合適的去除誤匹配和保留足夠的匹配特征對(duì)的方法，對(duì)于基于特征的配準(zhǔn)方法提高配準(zhǔn)精度是有效可行的。本研究的創(chuàng)新點(diǎn)在于合理地運(yùn)用了特征的雙邊匹配策略，在基于方向特征去除誤匹配的前提下加入了空間幾何特性的一致性檢測(cè)，有助于保留足夠的正確匹配特征對(duì)，為提高參數(shù)估計(jì)的精度奠定了基礎(chǔ)。

病變眼底圖像配準(zhǔn)是目前醫(yī)學(xué)眼底圖像處理和分析中的一項(xiàng)挑戰(zhàn)性課題，其研究工作的開展和研究成果的獲得將直接影響到該學(xué)科領(lǐng)域的發(fā)展和臨床眼科診斷水平的提高。在前期研究工作的基礎(chǔ)上，提出了一種基于特征的病變眼底圖像配準(zhǔn)方法。通過(guò)28組病變眼底圖像進(jìn)行實(shí)驗(yàn)仿真和分析表明，對(duì)于不同類型的病變眼底圖像該方法都具有較好的配準(zhǔn)精度，具有一定的臨床參考價(jià)值和應(yīng)用價(jià)值。在對(duì)匹配特征進(jìn)行一致性檢測(cè)時(shí)，是基于同組眼底圖像尺度變化不大的前提下進(jìn)行的，如何實(shí)現(xiàn)病變眼底圖像大尺度變化下的一致性檢測(cè)及多幅眼底圖像配準(zhǔn)是下一步研究的方向。

［1］王翠平，李新光，季亞成，等.高血壓和糖尿病患者視網(wǎng)膜病變與腦梗死的關(guān)系［J］.實(shí)用心腦肺血管病雜志，2007，15（2），152-156.

［2］Can A，Stewart C，Roysam B，et al.A feature-based，robust，hierarchical algorithm for registering pairs of images of the curved human retina［J］.IEEE Tran Pattern Ana.Mach Intel，2002，24：347-364.

［3］Stewart C，Tsai CL，Roysam B.The dual-bootstrap iterative closestpoint algorithm with application to retinal image registration［J］.IEEE Trans Med Imag，2003，22（11）：1379-1394.

［4］Zana F，Klein J.A multimodal registration algorithm of eye fundus images using vessels detection and hough transform［J］.IEEE Trans Biomed Eng，1999，18（5）：419-428.

［5］Wei LiFang，Huang LinLin，Pan，Lin，et al.The retinal image mosaic based on invariant feature and hierarchial transformation models［C］//Proceedings Image and Signal Processing，CISP.Tianjin：IEEE，2009：3463-3467.

［6］Harris C，.Stephens MJ.A combined corner and edge detector［C］//Proceedings of the 4th Alvey Vision Conference.Manchester：IEEE，1988：147-152.

［7］Bay H，Tuytelaars T，Gool LV.Surf：Speeded up robust features［C］//Proceedings of the 9th European Conference on Computer Vision.Graz：Springer，2006：404-417.

［8］Lowe D.Object Recognition from local scale-invariant features［C］//Proceedings of the International Conference on Computer Vision.Washington，DC：IEEE Computer Society，1999：1150-1157.

［9］Lowe D.Distinctive Image Features from Scale-Invariant Keypoints［J］.International Journal of Computer Vision，2004，60（2）：91-110.

［10］Chen Jian，Smith R，Tian Jie，et al.A novel registration method for retinal images based on local features［C］//Proceedings of IEEE Eng Med Biol Soc.Vancouver：IEEE，2008：2242-2245.

［11］LI H，Manjunath BS，Mitra SK.A contour-based approach to multisensor image registration［J］.IEEE Trans on Image Processing，1995，4（3）：320-334.

［12］邵婷婷，鄭穗聯(lián)，王波，等.正常個(gè)體人眼角膜空間形態(tài)數(shù)學(xué)建模路線研究［J］.眼視光學(xué)雜志，2005，7（4）：253-256.

［13］Burek H，Douthwaite WA.Mathematical models of the general corneal Surface［J］.Ophthal Physiol Opt，1993，13（1）：68-72.

［14］Richard L，George S，et al.Descriptors of corneal shape［J］.Op2tom Vis Sci，1998，75（2）：156-158.

［15］Rousseeuw PJ，Leroy AM.Robustregression and outlier detection［M］.New York：John Wiley ＆ Sons，1987：216-247.

［16］Stewart C.Robust parameter estimation in computer vision［J］.Society for Industrial and Applied Mathematics.1999，41（3）：513-537.