一類多重積分蒙特卡羅近似求解及其局部加權回歸擬合

許昌林，舒洪銘

(1.北方民族大學 數學與信息科學學院，寧夏 銀川 750021;2.北方民族大學 寧夏智能信息與大數據處理重點實驗室，寧夏 銀川 750021)

0 引言

對于多重積分的計算問題，若利用傳統的將多重積分化為累次積分方法進行求解，則有時會受到如積分區域的形狀、維數、被積函數的復雜程度等因素的影響，使多重積分求解變得復雜甚至無法得到解。為此，許多學者從概率論的角度出發，尋求多重積分的概率求解方法，如通過積分區域的投影、依概率收斂、控制收斂定理以及勒貝格積分控制收斂定理等方法進行多重積分的求解［1-4］。還有學者從數值計算的角度出發，通過隨機模擬法尋求多重積分的近似求解，如隨機均勻數法［5］、蒙特卡羅方法［6］、利用抽樣技術改進重積分的蒙特卡羅計算［7］、平均值估計計算法［8］、數論網格法［9］等。

多重積分在社會各個領域有重要應用，如利用地震多重積分對儲層位置和厚度進行識別與表示，從而對地震阻抗反演進行建模研究［10-11］，以反映地震剖面上的地質層位和構造特征;在計及光伏發電相關性的研究方面，利用多重積分法和Gram-Charlier 級數解決計及光伏出力相關性概率潮流的輸入和輸出變量之間非線性關系帶來的計算復雜性［12］;在稠密稀疏N 體問題解析函數近似計算方面，利用多重積分表示粒子間相互作用徑向分布函數的解析表達式，從而得到稀疏區域到稠密球域的分析解，以計算蛋白質的徑向函數［13］;非線性弱奇異多重積分不等式中未知參數的估計和應用［14-15］。

本文中，首先，對文獻［3］提出的閉區域上一類復雜n 重積分在n→∞時的極限問題，通過構造適當的概率分布，利用辛欽大數定律和依概率收斂理論，對該n 重積分的極限問題給出了新的證明方法和思路，并該證明方法簡潔易懂;其次，利用隨機變量函數的數學期望與重積分之間的關系，對n 重積分進行離散化處理，通過蒙特卡羅法對n 重積分進行模擬仿真，并在蒙特卡羅法近似計算結果的基礎上利用局部加權回歸法對計算結果進行擬合，利用R 軟件給出計算過程的可視化，進行近似計算和曲線擬合;然后，結合模擬計算和文獻［3-4］的結論對文獻［3］給出的一種n 重積分極限中的參數取值范圍進行修正，將文獻［3-4］給出的固定區域［0，1］×［0，1］×…×［0，1］上的n 重積分極限的結論推廣至一般區域［0，u］×［0，u］×…×［0，u］上，并利用本文所提出的證明方法對該結論進行了證明;最后，結合蒙特卡羅算法對一般區域［0，u］×［0，u］×…×［0，u］上n 重積分進行模擬計算，并利用局部加權回歸法對其進行擬合，從而進一步驗證推廣結論的合理性。

1 基于概率分布的一類多重積分極限問題的求解

文獻［3］應用依概率收斂和控制收斂定理及其相關結論，對式(1)的n 重積分的極限進行證明，證明過程復雜且抽象，不容易理解。本文中借助概率論中隨機變量的概率分布，給出一種新的證明方法。

針對式(1)，下面給出具體證明過程。

為了證明和表示方便，令

并將W 劃分為

由于

為了求解式(2)中的解，根據概率知識，假設連續性隨機變量序列{Xi}(i=1，2，3，…，n)獨立且同分布于上的均勻分布，則隨機變量序列{Xi}的聯合密度函數為:

基于上述概率分布以及連續性隨機變量在某區域上的概率與積分的關系，首先考慮式(2)在積分區域上的解。

其次，考慮式(2)在區域S 上的積分情況。

由式(4)和式(5)知:

另外，由于隨機變量Xi與sin Xi的期望存在，且分別為

故由辛欽大數定律以及依概率收斂可知［16］，對任意ε＞0 有

因此，對于式(6)由極限的夾迫性定理［17］可得

綜上，式(1)得證。

2 一類多重積分近似求解的蒙特卡羅法及其局部加權回歸擬合

2.1 n 重積分離散化

為了給出式(1)n 重積分在n→∞時的近似求解過程，首先基于多維連續性隨機變量函數的數學期望與重積分之間的關系，將式(1)n 重積分離散化處理。由式(3)知，式(1)n 重積分表示為隨機變量{Xi}(i=1，2，3，…，n)的函數的數學期望，即

從而根據數學期望的性質，將式(7)進一步離散化為

其中，式(8)右端隨機變量{Xi}(i=1，2，3，…，n)獨立同分布于上的離散均勻分布。

2.2 基于蒙特卡羅法的近似模擬算法

根據式(3)、式(7)和式(8)，得到計算式(1)n 重積分在n→∞時的蒙特卡羅算法，具體表述如算法1。

算法1:式(1)中n 重積分近似求解的蒙特卡羅算法

利用算法1，通過R 軟件，對不同的m 值(m=2 000，m=5 000)，當n 重積分的重數n 從100 以步長100 增加至10 000 時，式(1)的n 重積分的近似值AV 與n→∞時極限值的仿真散點圖如圖1 所示。通過圖1 可以發現，隨著重積分重數n 的不斷增大，近似值AV 越來越接近n→∞時的極限值LV，并在極限值LV 附近隨機波動。由此可以看出，式(8)的離散化結果能更好地反映n 重積分當n→∞時的極限過程。

圖1 不同m 值下式(1)n 重積分n 增加時AV 和LV 的變化趨勢

2.3 基于局部加權回歸的曲線擬合

圖1 仿真結果是根據蒙特卡洛算法隨機模擬計算得到。從圖1 中可以看出，當重積分重數n 增加時，近似值AV 在其極限值LV 附近具有一定的波動性。為了能更好地反映式(1)的n 重積分當n→∞時的近似值AV 與其極限值LV 的逼近過程，在圖1 相應條件下的散點圖(數據)基礎上，利用非參數回歸的方法對其進行擬合。其中，局部加權回歸(Lowess)是一種非參數回歸方法［18-20］，能很好地處理這種問題，它主要是把樣本劃分成一個個小區間，并對區間中的樣本進行多項式擬合，不斷重復這個過程得到在不同區間的加權回歸曲線，最后再把這些回歸曲線的中心連在一起合成完整的回歸曲線，具體過程如下:①決定擬合點的數量和位置;②以擬合點為中心，確定k 個最接近的點;③通過權重函數計算這k 個點的權重;④通過加權線性回歸進行多項式擬合，確定k 個平滑點;⑤對所有擬合點重復以上步驟;⑥將所有平滑點連接起來得出擬合曲線。

基于上述Lowess 擬合過程，分別對圖1 中m=2 000 和m=5 000 的條件下，重積分重數n 從100 以步長100 增加至10 000 時的隨機模擬計算結果進一步進行Lowess 擬合，擬合結果如圖2 所示。從圖2可以看出，對于不同m 取值，當重積分重數n 不斷增加時，Lowess 擬合曲線LWRV 逐步逼近重積分的極限值LV。由此看出，Lowess 擬合曲線能較好地反映出n 重積分收斂的可視化過程。

圖2 不同m 值下(1)式n 重積分n 增加時AV、LV 與LWRV 的變化趨勢

通過蒙特卡羅法和局部加權回歸擬合對n 重積分的收斂過程進行可視化，一方面可以驗證重積分的收斂性，另一方面可以利用蒙特卡羅法對復雜重積分進行近似求解，這為求解多重積分又提供了一種新的思路和方法。為此，進一步對文獻［3-4］中給出的n 重積分極限及其相關結論進行研究和推廣，并采用蒙特卡羅法對其進行近似計算，同時利用局部加權回歸擬合對其進行模擬。

3 一類推廣的多重積分極限的概率方法求解與模擬計算

3.1 已有結論的修正與模擬計算

根據上述式(1)的n 重積分極限的概率求解過程，類似地可應用于對文獻［3］給出的如式(9)所示的n 重積分極限的證明。

其中:文獻［3］針對式(9)中參數p，q 的要求是q＞p 且q≠-1。顯然，當(1+p)(1+q)＜0 時，式(9)中的結論不成立，所以結合文獻［3-4］如定理1 所示的結論，參數p，q 的要求修正為p，q＞0 且均為常數。

定理1［3-4］設f(x)與g(x)都是［0，1］上的實值連續函數，若存在某個常數c＞0，有不等式0≤f(x)＜cg(x)成立，則有

為了對式(9)的n 重積分極限過程進行模擬計算，類似式(7)、(8)的離散化過程，式(9)的n 重積分離散化為

其中，連續性隨機變量序列{Xi}(i=1，2，3，…，n)獨立且同分布于［0，1］上的均勻分布。基于式(10)，得到如算法2 所示的蒙特卡羅模擬算法。

算法2:式(9)n 重積分近似求解的蒙特卡羅算法

利用算法2 和Lowess 擬合過程，對不同m 值(m=2 000，m=5 000)，當n 重積分的重數n 從100 以步長100 增加至10 000 時，式(9)的n 重積分的近似值AV、n→∞時極限值與Lowess 值LWRV的仿真散點圖如圖3 所示，這里p、q 分別取2 組不同值。通過圖3 可以發現，隨著重積分重數n 的不斷增大，近似值AV 越來越接近n→∞時的極限值LV，并在極限值LV 附近隨機波動，Lowess 擬合曲線也逐步逼近重積分的極限值LV，并且能看出參數p、q 修正為p，q＞0 是合理的。

圖3 不同m，p，q 值下(9)式n 重積分n 增加時AV、LV 與LWRV 的變化趨勢

3.2 一類推廣的多重積分的極限

定理1 是在［0，1］×［0，1］×…×［0，1］區域上n 重積的極限，該結論可進一步推廣到區域［0，u］×［0，u］×…×［0，u］上，得到下列n 重積分的極限。

定理2設a(x)與b(x)為［0，u］上的實值連續函數，若存在某常數c＞0，使得0＜a(x)＜cb(x)成立，則有

證明類似地，引入概率相關知識進行證明。假設連續性隨機變量序列{Xi}(i=1，2，3，…，n)獨立且同分布于［0，u］上的均勻分布，則隨機變量序列{Xi}的聯合密度函數為

記隨機變量函數a(Xi)與b(Xi)的期望分別為

由于

再考慮式(12)在積分區域T 上的情況。由于E(a)=CE(b)，所以

綜合式(13)(14)，可得

由辛欽大數定律以及依概率收斂可知［16］，對任意ε＞0 有

所以，對式(15)由極限的夾迫性定理［17］，可得

綜上，式(11)得證。

3.3 推廣的一類多重積分極限的模擬計算

根據定理2，結合式(9)的n 重積分極限，對式(16)的n 重積分及其極限利用蒙特卡羅法進行近似求解，并利用Lowess 曲線擬合進行模擬。

首先將(16)式n 重積分極限離散化為

其中，連續性隨機變量序列{Xi}(i=1，2，3，…，n)獨立且同分布于［0，u］上的均勻分布。基于式(17)，得到如算法3 所示的蒙特卡羅模擬算法。

算法3:式(16)n 重積分近似求解的蒙特卡羅算法

為了對式(16)的n 重積分的極限過程利用算法3 和Lowess 擬合過程進行模擬計算，將積分區域中的參數u 按照0＜u＜1 和u＞1 兩種情況取值，被積函數中的參數p，q 按照p＞q＞0 和0＜p＜q 兩種情況取值，參數u，p，q 取值情況如表1 所示。當n 重積分的重數n 從100 以步長100 增加至10 000 時，式(16)n 重積分的近似值AV、n→∞時極限值與Lowess 值LWRV 的仿真散點圖分別如圖4、圖5 所示。從圖4 和圖5 顯示的結果中可以看出，當0＜u＜1 以及u＞1 時，在不同參數p，q 取值下，式(16)相應重積分的近似模擬值(AV 曲線)都能隨著n 的增大逐漸逼近其極限值(LV 曲線)，數值計算結果表明，定理2 的結論是合理的。

表1 參數u，p，q 取值與n→∞時(16)式n 重積分的極限值

圖4 當m=5 000，0＜u＜1，不同p，q 值式(16)n 重積分n 增加時AV、LV 與LWRV 的變化趨勢

圖5 當m=5 000，u＞1，不同p，q 值式(16)n 重積分n 增加時AV、LV 與LWRV 的變化趨勢

4 結論

從概率論角度出發，通過構造適當概率論分布，在利用辛欽大數定律和依概率收斂定理等相關概率知識的基礎上，對一類復雜n 重積分當n→∞時的極限問題進行了詳細證明。基于多維連續性隨機變量的數學期望與重積分之間的關系，對n 重積分進行離散化處理，并利用蒙特卡羅算法對離散化結果在n→∞時的極限過程進行近似計算和可視化。同時，利用一種非參數回歸方法——Lowess 法對n→∞時的近似計算結果進行曲線擬合，發現Lowess 擬合曲線能較好地反映n 重積分收斂的可視化過程。