我虛構的例子是這樣的：
 |
| 圖片由Gemini 生成，示意用。 |
有位北漂族偶而會開車回屏東老家坐坐，有時1季回去10幾次或20幾次或30幾次不等。假設從2025/1開始記錄，創建一個excel檔，第一欄為停車案例序，第二欄為停車位置。
在1/3回屏東，車停家門口 (表示事件發生，註記為1)，第1筆觀察值註記為1(第一欄)和1(第二欄); 在1/8回屏東，車停較遠處的農田小路 (表示事件未發生，註記為0)，走路約3-5分鐘吧，第2筆觀察值註記為2(第一欄)和0(第二欄)，長期監測下來的資料長這樣：
No event
1 1
2 0
3 0
4 0
5 1
6 0
......
依此類推。只要車停家門口，對面的老太婆就會出來看一下。若看到車停在路邊，她壓力值就會上升(她厭惡車停此會造成她家乖兒回來車無法直接駛入庭院停車，就算倒車退一次喬個方向再駛入也不行！)；若後續有回來，但車未停門口，壓力值就會依時間而消退（注意：假設此例1/3回來和1/8回來隔了5天壓力不隨時間消退，餘類推）。
此例想知道若每次回家車都停家門口和回來偶而停家門口，看她老人家壓力值曲線波動情形，並估算出閾值以衡量這支指標「車停家門口」表現是否顯著過高 (定義: 老太婆出面叫人把車開走表示車真的太常停家門口了，此時壓力值超過閾值)。
統計圖如下：
統計圖要表達的資訊已完美呈現就不再額外陳述。因我的指標為1/0有或無，故也可稱為Bernoulli EWMA方法(EWMA, Exponentially Weighted Moving Average)，關於Bernoulli分布，可見先前寫過的文章。在此處，y軸的正式名稱為「最近事件發生機率的指數加權平滑估計值」(如果你想貼在論文期刊，則可用EWMA statistic或EWMA estimate of event probability)，它是最近事件發生機率，最近兩字很重要，因當下機率值仍舊有考量近期事件發生有無(上圖註腳2已說明)，是有累積性的，而這些事件發生有或無有經統計學上的「平滑方法」處理，其機制是將新資料與過去平均值進行加權混合，使指標不會因單一 0/1 事件而產生劇烈波動。EWMA的核心精神就是最近資料對當下的影響最重要，愈久的資料影響愈小所以權重愈小。當 lambda 較小(如 0.05)時，曲線會非常平滑；lambda 較大(如 0.5)時，平滑程度較低且對新資料反應較敏感，實務上常使用 0.2。至於要精確到哪個值，有公式可用。假設P0=0.13, UCL=0.46,如果想連續3事件發生就呈現顯著惡化訊息，經公式推導可得知 lambda ≥ 0.145(為理論值, 但實務上會再保守點，往上拉些)。
而閾值估算過程簡易說明(我沒有精算，總之概算就是這樣)如下：
 |
| 修正用字：為求整數與盡早偵測到惡化訊號的結果 |
上圖中的0.00135即對應到下下圖最右邊的0.13% (做指標監測的人，請永遠記住這個數字)
z值=(UCL-p0)/sigma z=(0.543-0.18)/0.121=3, 當z值大於3，曲線下面積(機率)為多少？
P(Z>3)=0.00135
此例子中，可發現序29-序33均出現事件，有明顯群聚事件(佐證於下段)。接著可排除這些群聚事件後，用相對較穩的資料來計算基準事件率 (即上述公式中的p0)，就可以得到比上述閾值更低的值 (即惡化標準降低了，更容易發生「惡化事件」)，我沒再算，不過新值應會比0.5再低些。所以前述說我估0.5算是取折衷。
數據雖然都公式理論算出，但不否認的說，最後也都是折衷出來的。像上述p0=0.18, 平滑參數lambda＝0.18, UCL=0.543 LCU=0，以此條件下去模擬多次，看EWMA平均要多少筆資料(ARL, Average Run Length）)才會誤觸UCL, 若50筆就是假警報很多，若200筆假警報就少，若能達370筆則會到3倍sigma管制圖程度 (也可以反推每一筆資料假警報機率)。要取一個折衷到假警報不會太多，而偵測時也不會太遲鈍。
明顯群聚事件=驗證方法即我授課中曾教過的Run chart，引用的規則為下圖的規則3，Runs數要在合理區間 (細節可看此文章)。
補充說明：
記得，X軸資料均得依發生先後序排列，因為是用時間序列方法來分析。
1.其它例子：X軸為每日，紅點位置就是班上小明(特大號問題學生)該日有否惹事，藍線就是國小導師的壓力值，當藍線超過黑線，就該判定「老師敵對凶惡分數」指標表現顯著過高，該日老師的心情非常不好，隨時會爆炸。
2..其它例子：X軸為每日，紅點位置就是上班途中是否被車撞，當藍線超過黑線就該判定「上班日出車禍事件數」指標表現顯著過高，小明應擬定改善措施(e.g., 盡早去廟宇點光明燈XD)。
3.這種模型也能輕易應用於醫療場域以衡量臨床或管理指標表現是否異常，專門對付那些發生機率極低或極高的指標。第1件正常或不正常，再來第2件正常或不正常，案件依時間序一件一件來並標定為正常(0)或不正常(1)，最後就會得到一個時間序列(e.g., 0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,1,0,0,0,0,0,0,1,1,0,0,1,0,0,1,0,0,0,0....)，這類指標千萬別用傳統管制圖方法(這類圖的標準差公式)衡量表現是否正常。我用管制圖假設理論與機率論來說明，可以更清楚理解原因。
4.相較於TCPI指標庫，疾病認證指標庫其正向指標佔比大很多(e.g., 追蹤率、完治率、復健率、衛教宣導率、給予率、評估率、...)，若指標所收的案例大都有達成(達成佔比若不高，就不用認證了，沒人會願意跟錢過意不去，這不是一筆小數目)，p管制圖上每季(若監測頻率為季，以下同)表現值將落於上界線附近(甚至大多數資料點都100%)致統計假設失效(即圖上的訊號無法信任)。故可改監測單一序列，1表示負面事件發生而0表示未發生(e.g., 1:未完整 0:完整; 1:未完治 0:完治)，然後套此篇文章所介紹方法來衡量判定指標表現異常的時機和指標當下的表現機率估計值。每一指標每季的分子分母數據仍須有一張表呈現，然當講到指標表現異常時機判定和當下指標表現情形是否異常，則改用上述 No event停車案例來排列每一樣本並以此方法計算相關參數以及繪圖(我目前已將此計算方法做成SAS版本和Excel版本)。
傳統管制圖(e.g., p chart, u chart, x bar chart, c chart,...)都有一項最核心假設，即資料必須來自具有變異的隨機過程。若指標值畫出來，該率值線幾乎貼在天花板上(100%)或躺在地上(0%)，表示變異極小，資料分布呈極端偏態，此時統計假設就失效(很多碩士論文使用統計模型(e.g., 線性迴歸)，模型使用前提基本假設未逐一驗證然後都跑得好開心^___^，下結論下的很有自信，統計學老師看了會哭哭)。若從機率分布角度來看，若資料像1,1,1,1,1,1,1,1,1,0.97,1,1,1 這樣，機率分布會變成P(X=1) =0.95, 資料卡在上界，這時透過 p chart公式，可以算出LCL和UCL都落在100％，你的值只要稍稍低於100%就異常，但它實際是偽陽性，因為它只是在正常波動。
監測這種稀少事件，有種管制圖叫g chart (分布是不對稱的)，在品質工程界常被使用(我說常被使用就是常被使用，不要質疑我的Ref，我本身就是Ref)，它適合偵測指標是否改善(上界限可以無限，即由你的資料數據分布決定)，但不適合偵測是否惡化(下界限受限於0，只要小於0，一律視為0，不視為0沒有意義)，統計界其實早就知道這問題 (有沒有覺得我很強？但不知為何，醫策會就是不願找我去教課)。
 |
| Lloyd P. Provost 是管制圖大師 |
 |
| 呼應前段說明 |
基於此，故剩CUSUM chart和EWMA chart可用 (t chart在這種案例不適用，基於某些原因)。不可否認，CUSUM chart對「群聚事件」偵測非常敏感，事件接二連三(它只是形容詞)出現真的會「立即警報」給你看！而EWMA chart偵測的是平均率上升，雖然不如CUSUM chart靈敏，但可透過閾值訂定來改善。只要算出的值超過閾值，就可以很有信心的下結論：指標表現已顯著性惡化，請單位擬定改善措施因應。因為CUSUM chart我有正式應用於醫療場域，所以本篇文章使用Bernoulli EWMA概念(因我的數據為二元資料型式, 故EWMA chart不適合)評估稀少事件指標表現是否異常，所以CUSUM chart可以下去領便當了。
在QI Macros裡，有這些管制圖可用。
延伸閱讀：
沒有留言:
張貼留言