基于分位數回歸方法的函數型數據在線控制圖

張秀玲，訾雪旻

（天津職業(yè)技術師范大學理學院，天津 300222）

利用函數型數據刻畫產品的某些特性，進而分析和監(jiān)控此類函數型數據是統計過程控制領域研究的熱點問題之一。目前，很多研究都是基于傳統最小二乘估計的方法監(jiān)控模型的回歸系數，從而建立有效控制圖。然而，相對于僅關注響應變量在解釋變量已知情況下條件期望的變化，在實際應用中監(jiān)控響應變量條件中位數或其他條件分位數的變化越來越受到學術界的重視。如2001年，Abrevaya[1]分析影響低體重新生兒因素問題；2001年，Bassett等[2]評估共同基金的投資類型；2005—2006年，Machado等[3-4]研究勞動力市場、工資結構分布等。假定樣本均值是位置參數的最優(yōu)估計的關鍵取決于假設觀測值來自一般的正態(tài)分布。如果觀測值不是來自同一個分布，如來自不同方差的正態(tài)分布，采用中位數做估計明顯比用均值做估計效果更好。自1978年Koenker等[5]提出的分位數回歸以來，出現了很多關于分位數回歸的研究，然而對于構造其應用統計過程控制中相應的控制圖這一問題研究相對較少。本文針對上述問題，基于分位數回歸方法結合多元指數加權移動平均控制圖（MEWMA），給出一種新的函數型數據在線監(jiān)控控制圖，并通過數值模擬說明該控制圖的有效性。

1 監(jiān)控線性函數型數據的MEWMA控制圖

對一般線性函數型模型和2007年Zou等[6]基于最小二乘估計方法建立的MEWMA控制圖進行描述。假設隨著時間收集第j個時刻的隨機樣本，則有觀測值（Xj，Yj），其Yj=（y1j，y2j，…，ynjj）′是nj維響應變量，Xj是nj×p（nj＞p）矩陣。當統計過程可控時，假設基本模型為：

式中：β=（β1，β2，…，βp）′是p維系數向量；誤差項εj=（ε1j，ε2j，…，εnjj）′都是獨立同分布于均值為0且協方差陣為σ2I的nj維多元正態(tài)隨機向量。

在不失一般性的前提下，假設Xj的形式是（1，Xj*），其Xj*是正交于1，1是nj維常數向量且所有分量都是1。否則也可以通過合適的變換得到這種形式。通常nj都是相等的（記作n），假設對不同時刻j的解釋變量Xj是固定的（記作X）。

將隨著時間收集的第j個時刻的隨機樣本記作{（xi，yij），i=1，2，…，n}。當統計過程可控時，假設響應變量和解釋變量的關系為：

Zou等[6]用MEWMA控制圖同時監(jiān)控模型（2）的截距β0、斜率β1和標準差σ。根據模型（1）定義：

式中：Φ-1（·）為標準正態(tài)累積分布函數的反函數；F（·；v）為自由度v的卡方分布函數即（σ））′是（p+1）維隨機向量。當統計過程可控時，Zj服從均值為0且協方差陣為∑的多元正態(tài)分布，其∑=將EWMA控制圖統計量定義為：

式中：W0為（p+1）維初始向量；λ為光滑參數且0＜λ≤1。若：

控制圖就會報警。

式中：L為控制限且L＞0。

Zou等[6]通過比較失控平均運行長度說明MEWMA控制圖比Kim等[7]3 個EWMA結合的控制圖能更有效地探測截距、斜率和標準差發(fā)生的漂移。

2 分位數回歸方法

式中：τ為下分位數且τ∈（0，1）。

式（9）是xi的線性函數，即：

式中：xi′=（xi1，xi2，…，xip）；β（τ）=（β1（τ），β2（τ），…，βp（τ））′為p維τ分位數系數向量且β（τ）取決于τ。對第j個時刻的{yij，i=1，2，…，n；j=1，2，…}，可用簡單優(yōu)化問題的解將yij的經驗條件分位數函數定義為：

式中：u為觀測值與其估計值的殘差。

在模型（2）中通過對yij與xi′β的殘差的估計損失和函數進行最小化求解可得β（τ）。假設R（τ）是最小化估計損失和函數，即：

式（13）根據式（12）可寫成：

根據Koenker等[8-9]，式（11）于τ的左導數?？梢酝ㄟ^單純形法、內點法和平滑法等對式（14）進行求解

3 基于分位數回歸方法的MEWMA控制圖

對一般線性函數型模型基于分位數回歸方法，結合MEWMA控制圖，給出一種新的MEWMA控制圖在線監(jiān)控回歸系數分位數。通過τ=0.5即中位數說明基于分位數回歸方法建立MEWMA控制圖。

根據式（2）和式（9），則yij的條件分位數函數也可以寫成：

式中：Qεij（τ）為εij的τ分位數函數。

模型（2）中εij均是獨立同分布于均值為0且方差為σ2的正態(tài)分布，則Qεij（0.5）=0；再根據式（10）可得β（τ）=β。假設εij的分布函數記作F（x），其概率密度函數記作f（x）。因為Qεij（0.5）=0，則F（x）的中位數為0，εij在中位數的概率密度為f（0）。根據Koenker等[5]定理，若σ2和f（0）相比滿足條件[2f（0）]-1＜σ，則最小絕對誤差估計（LAE）的漸近方差比最小二乘估計（LSE）的方差??；則LAE估計即回歸中位數（（0.5））比LSE回歸系數更有效。

當統計過程可控時，根據Bassett等[10]的定理可知的極限分布是正態(tài)分布；其均值為0且協方差陣為是來自F（x）的隨機樣本的樣本中位數的漸近方差即ω=[2f（0）]-1。

對模型（2）在線監(jiān)控p個回歸系數中位數，給出新的MEWMA控制圖。

式中：Zj為p維隨機向量。當統計過程可控時，Zj服從均值為0且協方差陣為（Q*）-1的多元正態(tài)分布，其Q*=nQ。將MEWMA控制圖統計量定義為：

式中：W0為p維初始向量。

若

控制圖就會報警。

4 模擬研究控制圖表現

根據Qi等[11-15]理論，通過平均運行長度（ARL）研究本文推薦的MEWMA控制圖監(jiān)控表現。在可控模型中τ=0.5、0.9，n=4，p=2，λ=0.2，β0（0.5）=3，β1（0.5）=2，β0（0.9）=3，β1（0.9）=2，σ2=1和xi=2，4，6，8。在參數已知的前提下，考慮其MEWMA控制圖的可控平均運行長度（IC ARL）大約為200。通過數值模擬得到失控平均運行長度（OC ARL），分析其MEWMA控制圖監(jiān)控表現。模型（2）中截距和斜率中位數及0.9分位數分別發(fā)生漂移時MEWMA控制圖的ARL如表1所示。

表1 模型（2）中截距和斜率中位數及0.9分位數分別發(fā)生漂移時MEWMA控制圖的ARL

表1中，其MEWMA 控制圖的IC ARL 為199.8，則將第I 類型錯誤概率控制在0.5%左右。在此基礎上截距和斜率中位數分別發(fā)生漂移，即β（00.5）漂移到β0（0.5）+δ1σ 和 β（10.5）漂移到 β（10.5）+δ2σ；其0.9 分位數發(fā)生同樣的漂移。截距或斜率中位數及0.9 分位數發(fā)生的漂移越大，其OC ARL 越小，則其MEWMA控制圖探測漂移的速度越快。它們分別發(fā)生0.15、0.02小漂移的 OC ARL 分別是 84.5、87.1、40.3 和 65.8；分別發(fā)生 0.3、0.05 中等漂移的 OC ARL 分別是 24、19.7、9.9 和14.7；則其MEWMA 控制圖可以有效地探測截距或斜率中位數及0.9 分位數發(fā)生的中小漂移。分別發(fā)生0.8、0.15 大漂移的OC ARL 都在2 步以內，則其MEWMA 控制圖可快速探測到截距或斜率中位數及0.9 分位數發(fā)生的大漂移。分別發(fā)生0.1、0.01 小漂移的 OC ARL 分別是 124.8、138.9、67.5 和 114.4；使用漸近控制限而不是精確控制限，會導致其MEWMA 控制圖探測漂移的速度推遲。模型（2）中截距和斜率中位數及0.9 分位數同時發(fā)生漂移時MEWMA 控制圖的OC ARL 如表2所示。

表2 模型（2）中截距和斜率中位數及0.9 分位數同時發(fā)生漂移時MEWMA 控制圖的OC ARL

表2中，截距和斜率中位數及0.9 分位數同時發(fā)生漂移時，二者漂移同時變大或者任一個發(fā)生的漂移固定，另一個發(fā)生的漂移越大，其OC ARL 越小，則其MEWMA 控制圖探測漂移的速度越快。它們同時發(fā)生0.1 和 0.01 小漂移的 OC ARL 是 71.5 和 38.6，同時發(fā)生 0.25 和 0.025 中等漂移的 OC ARL 是 10.7 和 5.8；則其MEWMA 控制圖可以有效地探測截距和斜率中位數及0.9 分位數同時發(fā)生的中小漂移。同時，發(fā)生0.4和0.045 大漂移的OC ARL 是2.3 和1.7，則其MEWMA控制圖可快速探測到截距和斜率中位數及0.9 分位數同時發(fā)生的大漂移甚至OC ARL 在2 步以內。

5 結 語

本文基于MEWMA 控制圖對統計過程的可控和失控情況進行研究，研究結果表明：其可以有效地同時在線監(jiān)控函數型模型的截距和斜率中位數及0.9 分位數并能夠探測到它們發(fā)生的中小漂移。今后將研究基于分位數回歸方法建立累積和（CUSUM）等其他控制圖。